/
bse.h
159 lines (128 loc) · 4.88 KB
/
bse.h
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
/*
* Copyright 2009-2020 The VOTCA Development Team
* (http://www.votca.org)
*
* Licensed under the Apache License, Version 2.0 (the "License")
*
* You may not use this file except in compliance with the License.
* You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*
*/
#pragma once
#ifndef VOTCA_XTP_BSE_H
#define VOTCA_XTP_BSE_H
// Local VOTCA includes
#include "logger.h"
#include "orbitals.h"
#include "qmstate.h"
#include "threecenter.h"
namespace votca {
namespace xtp {
struct BSE_Population;
template <Index cqp, Index cx, Index cd, Index cd2>
class BSE_OPERATOR;
typedef BSE_OPERATOR<1, 2, 1, 0> SingletOperator_TDA;
typedef BSE_OPERATOR<1, 0, 1, 0> TripletOperator_TDA;
template <class T>
class QMFragment;
class BSE {
public:
// BSE(Logger& log, TCMatrix_gwbse& Mmn, const Eigen::MatrixXd& Hqp_in)
// : log_(log), Mmn_(Mmn), Hqp_in_(Hqp_in){};
BSE(Logger& log, TCMatrix_gwbse& Mmn) : log_(log), Mmn_(Mmn){};
struct options {
bool useTDA;
Index homo;
Index rpamin;
Index rpamax;
Index qpmin;
Index qpmax;
Index vmin;
Index cmax;
Index nmax; // number of eigenvectors to calculat
std::string davidson_correction;
std::string davidson_tolerance;
std::string davidson_update;
Index davidson_maxiter;
double min_print_weight; // minimium contribution for state to print it
bool use_Hqp_offdiag;
Index max_dyn_iter;
double dyn_tolerance;
};
void configure(const options& opt, const Eigen::VectorXd& RPAEnergies,
const Eigen::MatrixXd& Hqp_in);
void Solve_singlets(Orbitals& orb) const;
void Solve_triplets(Orbitals& orb) const;
Eigen::MatrixXd getHqp() const { return Hqp_; };
SingletOperator_TDA getSingletOperator_TDA() const;
TripletOperator_TDA getTripletOperator_TDA() const;
void Analyze_singlets(std::vector<QMFragment<BSE_Population> > fragments,
const Orbitals& orb) const;
void Analyze_triplets(std::vector<QMFragment<BSE_Population> > fragments,
const Orbitals& orb) const;
void Perturbative_DynamicalScreening(const QMStateType& type, Orbitals& orb);
private:
options opt_;
struct Interaction {
Eigen::VectorXd exchange_contrib;
Eigen::VectorXd direct_contrib;
Eigen::VectorXd qp_contrib;
};
struct ExpectationValues {
Eigen::VectorXd direct_term;
Eigen::VectorXd cross_term;
};
Logger& log_;
Index bse_vmax_;
Index bse_cmin_;
Index bse_size_;
Index bse_vtotal_;
Index bse_ctotal_;
Index max_dyn_iter_;
double dyn_tolerance_;
Eigen::VectorXd epsilon_0_inv_;
TCMatrix_gwbse& Mmn_;
Eigen::MatrixXd Hqp_;
tools::EigenSystem Solve_singlets_TDA() const;
tools::EigenSystem Solve_singlets_BTDA() const;
tools::EigenSystem Solve_triplets_TDA() const;
tools::EigenSystem Solve_triplets_BTDA() const;
void PrintWeights(const Eigen::VectorXd& weights) const;
template <typename BSE_OPERATOR>
void configureBSEOperator(BSE_OPERATOR& H) const;
template <typename BSE_OPERATOR>
tools::EigenSystem solve_hermitian(BSE_OPERATOR& h) const;
template <typename BSE_OPERATOR_ApB, typename BSE_OPERATOR_AmB>
tools::EigenSystem Solve_nonhermitian(BSE_OPERATOR_ApB& apb,
BSE_OPERATOR_AmB&) const;
template <typename BSE_OPERATOR_A, typename BSE_OPERATOR_B>
tools::EigenSystem Solve_nonhermitian_Davidson(BSE_OPERATOR_A& Aop,
BSE_OPERATOR_B& Bop) const;
void printFragInfo(const std::vector<QMFragment<BSE_Population> >& frags,
Index state) const;
void printWeights(Index i_bse, double weight) const;
void SetupDirectInteractionOperator(const Eigen::VectorXd& DFTenergies,
double energy);
Eigen::MatrixXd AdjustHqpSize(const Eigen::MatrixXd& Hqp_in,
const Eigen::VectorXd& RPAInputEnergies);
Interaction Analyze_eh_interaction(const QMStateType& type,
const Orbitals& orb) const;
template <typename BSE_OPERATOR>
ExpectationValues ExpectationValue_Operator(const QMStateType& type,
const Orbitals& orb,
const BSE_OPERATOR& H) const;
template <typename BSE_OPERATOR>
ExpectationValues ExpectationValue_Operator_State(
const QMState& state, const Orbitals& orb, const BSE_OPERATOR& H) const;
};
} // namespace xtp
} // namespace votca
#endif // VOTCA_XTP_BSE_H