Line data Source code
1 : /*
2 : * Copyright (c) 2016 Sebastian Weber, Henri Menke. All rights reserved.
3 : *
4 : * This file is part of the pairinteraction library.
5 : *
6 : * The pairinteraction library is free software: you can redistribute it and/or modify
7 : * it under the terms of the GNU Lesser General Public License as published by
8 : * the Free Software Foundation, either version 3 of the License, or
9 : * (at your option) any later version.
10 : *
11 : * The pairinteraction library is distributed in the hope that it will be useful,
12 : * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 : * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 : * GNU Lesser General Public License for more details.
15 : *
16 : * You should have received a copy of the GNU Lesser General Public License
17 : * along with the pairinteraction library. If not, see <http://www.gnu.org/licenses/>.
18 : */
19 :
20 : #ifndef MATRIXELEMENTCACHE_H
21 : #define MATRIXELEMENTCACHE_H
22 :
23 : #include "Basisnames.hpp"
24 : #include "State.hpp"
25 : #include "Wavefunction.hpp"
26 : #include "utils.hpp"
27 :
28 : #include <cereal/types/unordered_map.hpp>
29 : #include <cereal/types/unordered_set.hpp>
30 : #include <wignerSymbols/wignerSymbols-cpp.h>
31 :
32 : #include <memory>
33 : #include <sstream>
34 : #include <string>
35 : #include <tuple>
36 : #include <unordered_map>
37 :
38 : enum method_t {
39 : NUMEROV = 0,
40 : WHITTAKER = 1,
41 : };
42 :
43 : class StateOne;
44 : class StateTwo;
45 :
46 : bool selectionRulesMomentumNew(
47 : StateOne const &state1, StateOne const &state2,
48 : int q); // TODO rename, integrate into the MatrixElementCache namespace
49 : bool selectionRulesMomentumNew(StateOne const &state1, StateOne const &state2);
50 : bool selectionRulesMultipoleNew(StateOne const &state1, StateOne const &state2, int kappa, int q);
51 : bool selectionRulesMultipoleNew(StateOne const &state1, StateOne const &state2, int kappa);
52 :
53 : class MatrixElementCache {
54 : public:
55 : MatrixElementCache();
56 : MatrixElementCache(std::string const &cachedir);
57 :
58 : double getElectricDipole(StateOne const &state_row,
59 : StateOne const &state_col); // return value in GHz/(V/cm)
60 : double getElectricMultipole(StateOne const &state_row, StateOne const &state_col,
61 : int k); // return value in GHz/(V/cm)*um^(kappa-1)
62 : double getDiamagnetism(StateOne const &state_row, StateOne const &state_col,
63 : int k); // return value in GHz/(V/cm)*um
64 : double getMagneticDipole(StateOne const &state_row,
65 : StateOne const &state_col); // return value in GHz/G
66 : double
67 : getElectricMultipole(StateOne const &state_row, StateOne const &state_col, int kappa_radial,
68 : int kappa_angular); // return value in GHz/(V/cm)*um^(kappa_radial-1)
69 : double getRadial(StateOne const &state_row, StateOne const &state_col,
70 : int kappa); // return value in um^kappa
71 :
72 : void precalculateElectricMomentum(const std::vector<StateOne> &basis_one, int q);
73 : void precalculateMagneticMomentum(const std::vector<StateOne> &basis_one, int q);
74 : void precalculateDiamagnetism(const std::vector<StateOne> &basis_one, int k, int q);
75 : void precalculateMultipole(const std::vector<StateOne> &basis_one, int k);
76 : void precalculateRadial(const std::vector<StateOne> &basis_one, int k);
77 :
78 : void setDefectDB(std::string const &path);
79 : const std::string &getDefectDB() const;
80 : void setMethod(method_t const &m);
81 : void loadElectricDipoleDB(std::string const &path, std::string const &species);
82 :
83 : size_t size();
84 :
85 : private:
86 : int update();
87 : void precalculate(std::shared_ptr<const BasisnamesOne> basis_one, int kappa, int q, int kappar,
88 : bool calcMultipole, bool calcMomentum, bool calcRadial);
89 : double calcRadialElement(const QuantumDefect &qd1, int power, const QuantumDefect &qd2);
90 : void precalculate(const std::vector<StateOne> &basis_one, int kappa_angular, int q,
91 : int kappa_radial, bool calcElectricMultipole, bool calcMagneticMomentum,
92 : bool calcRadial);
93 :
94 : ////////////////////////////////////////////////////////////////////
95 : /// Keys ///////////////////////////////////////////////////////////
96 : ////////////////////////////////////////////////////////////////////
97 :
98 : struct CacheKey_cache_radial {
99 : CacheKey_cache_radial(method_t method, std::string species, int kappa, int n1, int n2,
100 : int l1, int l2, float j1, float j2);
101 : CacheKey_cache_radial();
102 : bool operator==(const CacheKey_cache_radial &rhs) const;
103 : std::string species;
104 : method_t method;
105 : int kappa;
106 : std::array<int, 2> n, l;
107 : std::array<float, 2> j;
108 :
109 : private:
110 : friend class cereal::access;
111 : template <class Archive>
112 570 : void serialize(Archive &ar, unsigned int /* version */) {
113 570 : ar &CEREAL_NVP(method) & CEREAL_NVP(species) & CEREAL_NVP(kappa) & CEREAL_NVP(n) &
114 1140 : CEREAL_NVP(l) & CEREAL_NVP(j);
115 570 : }
116 : };
117 :
118 : struct CacheKey_cache_angular {
119 : CacheKey_cache_angular(int kappa, float j1, float j2, float m1, float m2);
120 : CacheKey_cache_angular();
121 : bool operator==(const CacheKey_cache_angular &rhs) const;
122 : int kappa;
123 : std::array<float, 2> j, m;
124 : int sgn;
125 :
126 : private:
127 : friend class cereal::access;
128 : template <class Archive>
129 920 : void serialize(Archive &ar, unsigned int /* version */) {
130 920 : ar &CEREAL_NVP(kappa) & CEREAL_NVP(j) & CEREAL_NVP(m) & CEREAL_NVP(sgn);
131 920 : }
132 : };
133 :
134 : struct CacheKey_cache_reduced_commutes {
135 : CacheKey_cache_reduced_commutes(float s, int kappa, int l1, int l2, float j1, float j2);
136 : CacheKey_cache_reduced_commutes();
137 : bool operator==(const CacheKey_cache_reduced_commutes &rhs) const;
138 : float s;
139 : int kappa;
140 : std::array<int, 2> l;
141 : std::array<float, 2> j;
142 : int sgn;
143 :
144 : private:
145 : friend class cereal::access;
146 : template <class Archive>
147 120 : void serialize(Archive &ar, unsigned int /* version */) {
148 120 : ar &CEREAL_NVP(s) & CEREAL_NVP(kappa) & CEREAL_NVP(l) & CEREAL_NVP(j) & CEREAL_NVP(sgn);
149 120 : }
150 : };
151 :
152 : struct CacheKey_cache_reduced_multipole {
153 : CacheKey_cache_reduced_multipole(int kappa, int l1, int l2);
154 : CacheKey_cache_reduced_multipole();
155 : bool operator==(const CacheKey_cache_reduced_multipole &rhs) const;
156 : int kappa;
157 : std::array<int, 2> l;
158 : int sgn;
159 :
160 : private:
161 : friend class cereal::access;
162 : template <class Archive>
163 40 : void serialize(Archive &ar, unsigned int /* version */) {
164 40 : ar &CEREAL_NVP(kappa) & CEREAL_NVP(l) & CEREAL_NVP(sgn);
165 40 : }
166 : };
167 :
168 : struct CacheKeyHasher_cache_radial {
169 : std::size_t operator()(const CacheKey_cache_radial &c) const;
170 : };
171 : struct CacheKeyHasher_cache_angular {
172 : std::size_t operator()(const CacheKey_cache_angular &c) const;
173 : };
174 : struct CacheKeyHasher_cache_reduced_commutes {
175 : std::size_t operator()(const CacheKey_cache_reduced_commutes &c) const;
176 : };
177 : struct CacheKeyHasher_cache_reduced_multipole {
178 : std::size_t operator()(const CacheKey_cache_reduced_multipole &c) const;
179 : };
180 :
181 : std::unordered_map<CacheKey_cache_radial, double, CacheKeyHasher_cache_radial> cache_radial;
182 : std::unordered_map<CacheKey_cache_angular, double, CacheKeyHasher_cache_angular> cache_angular;
183 : std::unordered_map<CacheKey_cache_reduced_commutes, double,
184 : CacheKeyHasher_cache_reduced_commutes>
185 : cache_reduced_commutes_s;
186 : std::unordered_map<CacheKey_cache_reduced_commutes, double,
187 : CacheKeyHasher_cache_reduced_commutes>
188 : cache_reduced_commutes_l;
189 : std::unordered_map<CacheKey_cache_reduced_multipole, double,
190 : CacheKeyHasher_cache_reduced_multipole>
191 : cache_reduced_multipole;
192 :
193 : std::unordered_set<CacheKey_cache_radial, CacheKeyHasher_cache_radial> cache_radial_missing;
194 : std::unordered_set<CacheKey_cache_angular, CacheKeyHasher_cache_angular> cache_angular_missing;
195 : std::unordered_set<CacheKey_cache_reduced_commutes, CacheKeyHasher_cache_reduced_commutes>
196 : cache_reduced_commutes_s_missing;
197 : std::unordered_set<CacheKey_cache_reduced_commutes, CacheKeyHasher_cache_reduced_commutes>
198 : cache_reduced_commutes_l_missing;
199 : std::unordered_set<CacheKey_cache_reduced_multipole, CacheKeyHasher_cache_reduced_multipole>
200 : cache_reduced_multipole_missing;
201 :
202 : method_t method{NUMEROV};
203 : std::string defectdbname;
204 : std::string dbname;
205 : std::unique_ptr<sqlite::handle> db;
206 : std::unique_ptr<sqlite::statement> stmt;
207 : long pid_which_created_db;
208 :
209 : ////////////////////////////////////////////////////////////////////
210 : /// Method for serialization ///////////////////////////////////////
211 : ////////////////////////////////////////////////////////////////////
212 :
213 : friend class cereal::access;
214 : template <class Archive>
215 20 : void serialize(Archive &ar, unsigned int /* version */) {
216 20 : ar &CEREAL_NVP(method);
217 20 : ar &CEREAL_NVP(dbname);
218 20 : ar &CEREAL_NVP(cache_radial) & CEREAL_NVP(cache_angular) &
219 40 : CEREAL_NVP(cache_reduced_commutes_s) & CEREAL_NVP(cache_reduced_commutes_l) &
220 20 : CEREAL_NVP(cache_reduced_multipole);
221 20 : ar &CEREAL_NVP(cache_radial_missing) & CEREAL_NVP(cache_angular_missing) &
222 40 : CEREAL_NVP(cache_reduced_commutes_s_missing) &
223 40 : CEREAL_NVP(cache_reduced_commutes_l_missing) &
224 20 : CEREAL_NVP(cache_reduced_multipole_missing);
225 :
226 20 : if (Archive::is_loading::value && !dbname.empty()) {
227 : // Open database
228 20 : db = std::make_unique<sqlite::handle>(dbname);
229 20 : stmt = std::make_unique<sqlite::statement>(*db);
230 20 : pid_which_created_db = utils::get_pid();
231 :
232 : // Speed up database access
233 20 : stmt->exec(
234 : "PRAGMA synchronous = OFF"); // do not wait on write, hand off to OS and carry on
235 20 : stmt->exec("PRAGMA journal_mode = MEMORY"); // keep rollback journal in memory during
236 : // transaction
237 : }
238 20 : }
239 : };
240 :
241 : #endif
|