ERKALE
ERKALE - DFT from Hel
 All Classes Functions Variables Friends Pages
casida.h
1 /*
2  * This file is written by Arto Sakko and Susi Lehtola, 2011
3  * Copyright (c) 2011, Arto Sakko and Susi Lehtola
4  *
5  */
6 
7 /*
8  *
9  *
10  * This file is part of
11  *
12  * E R K A L E
13  * -
14  * DFT from Hel
15  *
16  * ERKALE is written by Susi Lehtola, 2010-2011
17  * Copyright (c) 2010-2011, Susi Lehtola
18  *
19  *
20  * This program is free software; you can redistribute it and/or
21  * modify it under the terms of the GNU General Public License
22  * as published by the Free Software Foundation; either version 2
23  * of the License, or (at your option) any later version.
24  */
25 
41 #ifndef ERKALE_CASIDA
42 #define ERKALE_CASIDA
43 
44 class BasisSet;
45 class Settings;
46 
47 #include <armadillo>
48 #include <cstdio>
49 #include <cstdlib>
50 #include <cfloat>
51 
52 #include <xc.h>
53 
55 enum coupling_mode {
57  IPA,
59  RPA,
61  TDLDA
62 };
63 
68 typedef struct {
70  int i;
72  int f;
74 
76 class Casida {
78  std::vector< std::vector<states_pair_t> > pairs;
80  enum coupling_mode coupling;
81 
83  std::vector< arma::vec > f;
85  std::vector<size_t> nocc;
87  std::vector<size_t> nvirt;
88 
90  std::vector<arma::vec> E;
92  std::vector<arma::mat> C;
94  std::vector<arma::mat> P;
96  std::vector< std::vector<arma::mat> > dipmat;
97 
98  // K matrix
99  arma::mat K;
100 
102  arma::vec w_i;
104  arma::mat F_i;
105 
116  void Kcoul(const BasisSet & basis, const Settings & set);
117 
127  void Kxc(const BasisSet & bas, double tol, int x_func, int c_func);
128 
129  /* Helper routines */
130 
132  double esq(states_pair_t ip, bool ispin) const;
133 
135  double fe(states_pair_t ip, bool ispin) const;
136 
138  void form_pairs(const Settings & set, const std::vector< std::vector<double> > occs);
139 
141  arma::mat matrix_transform(bool ispin, const arma::mat & m) const;
143  arma::cx_mat matrix_transform(bool ispin, const arma::cx_mat & m) const;
144 
145 
147  void parse_coupling(const Settings & set);
148 
150  void calc_K(const Settings & set, const BasisSet & bas);
151 
153  void coulomb_fit(const BasisSet & basis, std::vector<arma::mat> & munu, arma::mat & ab_inv, const Settings & set) const;
154 
156  void solve();
157 
159  arma::mat transition(const std::vector<arma::mat> & M) const;
161  arma::mat transition(const std::vector<arma::cx_mat> & M) const;
162 
163  public:
165  Casida();
167  Casida(const Settings & set, const BasisSet & basis, const arma::vec & E, const arma::mat & C, const arma::mat & P, const std::vector<double> & occs);
169  Casida(const Settings & set, const BasisSet & basis, const arma::vec & Ea, const arma::vec & Eb, const arma::mat & Ca, const arma::mat & Cb, const arma::mat & Pa, const arma::mat & Pb, const std::vector<double> & occa, const std::vector<double> & occb);
171  ~Casida();
172 
174  arma::mat transition(const BasisSet & bas, const arma::vec & q) const;
175 
177  arma::mat transition(const BasisSet & bas, double q) const;
178 
180  arma::mat dipole_transition(const BasisSet & bas) const;
181 };
182 
183 
184 
185 #endif
Class for performing Casida calculations.
Definition: casida.h:76
void form_pairs(const Settings &set, const std::vector< std::vector< double > > occs)
Form pairs and occupations.
Definition: casida.cpp:161
arma::mat transition(const std::vector< arma::mat > &M) const
Compute the transition speeds using the matrix elements given by M.
Definition: casida.cpp:331
std::vector< arma::mat > P
Density matrices.
Definition: casida.h:94
arma::vec w_i
Eigenvalues of Casida&#39;s equation (and later on, the excitation energies)
Definition: casida.h:102
double esq(states_pair_t ip, bool ispin) const
Calculate .
Definition: casida.cpp:279
int i
Initial state.
Definition: casida.h:70
void solve()
Solve the Casida equation.
Definition: casida.cpp:288
arma::mat dipole_transition(const BasisSet &bas) const
Compute the dipole transition spectrum.
Definition: casida.cpp:399
std::vector< arma::vec > f
Occupancies of orbitals.
Definition: casida.h:83
arma::mat matrix_transform(bool ispin, const arma::mat &m) const
Transform given AO matrix to MO.
Definition: casida.cpp:153
Definition: casida.h:68
std::vector< arma::vec > E
MO energies for spin up and spin down.
Definition: casida.h:90
Casida()
Dummy constructor.
Definition: casida.cpp:45
void Kcoul(const BasisSet &basis, const Settings &set)
Definition: casida.cpp:705
int f
Final state.
Definition: casida.h:72
Settings used for a calculations.
Definition: settings.h:74
~Casida()
Destructor.
Definition: casida.cpp:150
void calc_K(const Settings &set, const BasisSet &bas)
Construct the K matrices.
Definition: casida.cpp:122
double fe(states_pair_t ip, bool ispin) const
Calculate .
Definition: casida.cpp:284
std::vector< size_t > nvirt
Number of virtual orbitals.
Definition: casida.h:87
std::vector< arma::mat > C
MO coefficient matrices for spin up and spin down.
Definition: casida.h:92
void parse_coupling(const Settings &set)
Common routines for constructors.
Definition: casida.cpp:99
Basis set.
Definition: basis.h:187
std::vector< std::vector< states_pair_t > > pairs
List of electron-hole pairs included in calculation.
Definition: casida.h:78
void coulomb_fit(const BasisSet &basis, std::vector< arma::mat > &munu, arma::mat &ab_inv, const Settings &set) const
Compute the Coulomb fitting integrals and the matrix .
Definition: casida.cpp:500
std::vector< size_t > nocc
Number of occupied orbitals.
Definition: casida.h:85
std::vector< std::vector< arma::mat > > dipmat
Dipole matrix elements: [nspin][3][norb,norb].
Definition: casida.h:96
enum coupling_mode coupling
How to couple electron-hole pairs (IPA, RPA, TDLDA)
Definition: casida.h:80
arma::mat F_i
Eigenvectors of Casida&#39;s equation.
Definition: casida.h:104
void Kxc(const BasisSet &bas, double tol, int x_func, int c_func)
Definition: casida.cpp:764