TBTK
RPASusceptibilityCalculator.h
1 /* Copyright 2017 Kristofer Bj√∂rnson
2  *
3  * Licensed under the Apache License, Version 2.0 (the "License");
4  * you may not use this file except in compliance with the License.
5  * You may obtain a copy of the License at
6  *
7  * http://www.apache.org/licenses/LICENSE-2.0
8  *
9  * Unless required by applicable law or agreed to in writing, software
10  * distributed under the License is distributed on an "AS IS" BASIS,
11  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
12  * See the License for the specific language governing permissions and
13  * limitations under the License.
14  */
15 
23 #ifndef COM_DAFER45_TBTK_RPA_SUSCEPTIBILITY_CALCULATOR
24 #define COM_DAFER45_TBTK_RPA_SUSCEPTIBILITY_CALCULATOR
25 
26 #include "TBTK/InteractionAmplitude.h"
29 #include "TBTK/RPA/SusceptibilityCalculator.h"
30 #include "TBTK/RPA/LindhardSusceptibilityCalculator.h"
31 
32 #include <complex>
33 
34 //#include <omp.h>
35 
36 namespace TBTK{
37 
39 public:
42  const MomentumSpaceContext &momentumSpaceContext,
44  = SusceptibilityCalculator::Algorithm::Lindhard
45  );
46 
49 
54 
58  void precompute(unsigned int numWorkers = 129);
59 
60  const MomentumSpaceContext& getMomentumSpaceContext() const;
61 
64  const std::vector<InteractionAmplitude> &interactionAmplitudes
65  );
66 
69  enum class EnergyType {Real, Imaginary, Complex};
70 
72 // void setSusceptibilityMode(SusceptibilityCalculator::Mode mode);
73 
75 // SusceptibilityCalculator::Mode getSusceptibilityMode() const;
76 
78  void setEnergyType(EnergyType energyType);
79 
81  EnergyType getEnergyType() const;
82 
85  void setEnergies(
86  const std::vector<std::complex<double>> &energies
87  );
88 
99  bool energiesAreInversionSymmetric
100  );
101 
105 
109  bool susceptibilityIsSafeFromPoles
110  );
111 
115 
117 // void setNumSummationEnergies(unsigned int numSummationEnergies);
118 
120  void saveSusceptibilities(const std::string &filename) const;
121 
123  void loadSusceptibilities(const std::string &filename);
124 private:
127 
129  IndexedDataTree<SerializableVector<std::complex<double>>> rpaChargeSusceptibilityTree;
130 
133 
135  EnergyType energyType;
136 
138  std::vector<std::complex<double>> energies;
139 public:
141  std::vector<std::complex<double>> calculateRPASusceptibility(
142  const DualIndex &kDual,
143  const std::vector<int> &orbitalIndices
144  );
145 
147  std::vector<std::complex<double>> calculateRPASusceptibility(
148  const std::vector<double> &k,
149  const std::vector<int> &orbitalIndices
150  );
151 
153  std::vector<std::complex<double>> calculateChargeRPASusceptibility(
154  const DualIndex &kDual,
155  const std::vector<int> &orbitalIndices
156  );
157 
159  std::vector<std::complex<double>> calculateChargeRPASusceptibility(
160  const std::vector<double> &k,
161  const std::vector<int> &orbitalIndices
162  );
163 
165  std::vector<std::complex<double>> calculateSpinRPASusceptibility(
166  const DualIndex &kDual,
167  const std::vector<int> &orbitalIndices
168  );
169 
171  std::vector<std::complex<double>> calculateSpinRPASusceptibility(
172  const std::vector<double> &k,
173  const std::vector<int> &orbitalIndices
174  );
175 
177  void setU(std::complex<double> U);
178 
180  void setUp(std::complex<double> Up);
181 
183  void setJ(std::complex<double> J);
184 
186  void setJp(std::complex<double> Jp);
187 private:
189  SusceptibilityCalculator *susceptibilityCalculator;
190 
192  std::vector<InteractionAmplitude> interactionAmplitudes;
193 
195  std::vector<InteractionAmplitude> interactionAmplitudesCharge;
196 
198  std::vector<InteractionAmplitude> interactionAmplitudesSpin;
199 
202  bool interactionAmplitudesAreGenerated;
203 
206  SusceptibilityCalculator &susceptibilityCalculator
207  );
208 
210  Index getSusceptibilityResultIndex(
211  const Index &kIndex,
212  const std::vector<int> &orbitalIndices
213  ) const;
214 
216  void invertMatrix(std::complex<double> *matrix, unsigned int dimensions);
217 
219  void multiplyMatrices(
220  std::complex<double> *matrix1,
221  std::complex<double> *matrix2,
222  std::complex<double> *result,
223  unsigned int dimensions
224  );
225 
227  std::vector<std::vector<std::vector<std::complex<double>>>> rpaSusceptibilityMainAlgorithm(
228  const DualIndex &kDual,
229  const std::vector<int> &orbitalIndices,
230  const std::vector<InteractionAmplitude> &interactionAmpltiudes
231  );
232 
234  std::complex<double> U, Up, J, Jp;
235 
239  void generateInteractionAmplitudes();
240 };
241 
242 inline const MomentumSpaceContext& RPASusceptibilityCalculator::getMomentumSpaceContext(
243 ) const{
244  return susceptibilityCalculator->getMomentumSpaceContext();
245 }
246 
248  const std::vector<InteractionAmplitude> &interactionAmplitudes
249 ){
250  this->interactionAmplitudes = interactionAmplitudes;
251 }
252 
253 inline Index RPASusceptibilityCalculator::getSusceptibilityResultIndex(
254  const Index &kIndex,
255  const std::vector<int> &orbitalIndices
256 ) const{
257  return Index(
258  kIndex,
259  {
260  orbitalIndices.at(0),
261  orbitalIndices.at(1),
262  orbitalIndices.at(2),
263  orbitalIndices.at(3)
264  }
265  );
266 }
267 
268 /*inline void RPASusceptibilityCalculator::setSusceptibilityMode(
269  SusceptibilityCalculator::Mode mode
270 ){
271  susceptibilityCalculator->setMode(mode);
272 }
273 
274 inline SusceptibilityCalculator::Mode RPASusceptibilityCalculator::getSusceptibilityMode(
275 ) const{
276  return susceptibilityCalculator->getMode();
277 }*/
278 
280  EnergyType energyType
281 ){
282  this->energyType = energyType;
283  switch(energyType){
284  case EnergyType::Real:
285  susceptibilityCalculator->setEnergyType(
286  SusceptibilityCalculator::EnergyType::Real
287  );
288  break;
289  case EnergyType::Imaginary:
290  susceptibilityCalculator->setEnergyType(
291  SusceptibilityCalculator::EnergyType::Imaginary
292  );
293  break;
294  case EnergyType::Complex:
295  susceptibilityCalculator->setEnergyType(
296  SusceptibilityCalculator::EnergyType::Complex
297  );
298  break;
299  default:
300  TBTKExit(
301  "RPASusceptibilityCalculator::setEnergyType()",
302  "Unknown energy type.",
303  "This should never happen, contact the developer."
304  );
305  }
306 }
307 
309 ) const{
310  return energyType;
311 }
312 
314  const std::vector<std::complex<double>> &energies
315 ){
316  susceptibilityCalculator->setEnergies(
317  energies
318  );
319 
320  this->energies = energies;
321 
322  rpaSusceptibilityTree.clear();
323  rpaChargeSusceptibilityTree.clear();
324  rpaSpinSusceptibilityTree.clear();
325 }
326 
328  bool energiesAreInversionSymmetric
329 ){
330  susceptibilityCalculator->setEnergiesAreInversionSymmetric(
331  energiesAreInversionSymmetric
332  );
333 }
334 
336 ) const{
337  return susceptibilityCalculator->getEnergiesAreInversionSymmetric();
338 }
339 
341  bool susceptibilityIsSafeFromPoles
342 ){
343  TBTKAssert(
344  susceptibilityCalculator->getAlgorithm()
345  == SusceptibilityCalculator::Algorithm::Lindhard,
346  "RPASusceptibilityCalculator::setSusceptibilityIsSafeFromPoles()",
347  "Only valid function call if the underlying algorithm is"
348  << " SusceptibilityCalculator::Algorithm::Lindhard.",
349  ""
350  );
351 
352  ((LindhardSusceptibilityCalculator*)susceptibilityCalculator)->setSusceptibilityIsSafeFromPoles(
353  susceptibilityIsSafeFromPoles
354  );
355 }
356 
358  TBTKAssert(
359  susceptibilityCalculator->getAlgorithm()
360  == SusceptibilityCalculator::Algorithm::Lindhard,
361  "RPASusceptibilityCalculator::getSusceptibilityIsSafeFromPoles()",
362  "Only valid function call if the underlying algorithm is"
363  << " SusceptibilityCalculator::Algorithm::Lindhard.",
364  ""
365  );
366 
367  return ((LindhardSusceptibilityCalculator*)susceptibilityCalculator)->getSusceptibilityIsSafeFromPoles();
368 }
369 
370 /*inline void RPASusceptibilityCalculator::setNumSummationEnergies(
371  unsigned int numSummationEnergies
372 ){
373  susceptibilityCalculator->setNumSummationEnergies(
374  numSummationEnergies
375  );
376 }*/
377 
379  const std::string &filename
380 ) const{
381  susceptibilityCalculator->saveSusceptibilities(filename);
382 }
383 
385  const std::string &filename
386 ){
387  susceptibilityCalculator->loadSusceptibilities(filename);
388 }
389 
390 inline std::vector<std::complex<double>> RPASusceptibilityCalculator::calculateRPASusceptibility(
391  const std::vector<double> &k,
392  const std::vector<int> &orbitalIndices
393 ){
395  DualIndex(
396  susceptibilityCalculator->getMomentumSpaceContext().getKIndex(k),
397  k
398  ),
399  orbitalIndices
400  );
401 }
402 
403 inline std::vector<std::complex<double>> RPASusceptibilityCalculator::calculateChargeRPASusceptibility(
404  const std::vector<double> &k,
405  const std::vector<int> &orbitalIndices
406 ){
408  DualIndex(
409  susceptibilityCalculator->getMomentumSpaceContext().getKIndex(k),
410  k
411  ),
412  orbitalIndices
413  );
414 }
415 
416 inline std::vector<std::complex<double>> RPASusceptibilityCalculator::calculateSpinRPASusceptibility(
417  const std::vector<double> &k,
418  const std::vector<int> &orbitalIndices
419 ){
421  DualIndex(
422  susceptibilityCalculator->getMomentumSpaceContext().getKIndex(k),
423  k
424  ),
425  orbitalIndices
426  );
427 }
428 
429 inline void RPASusceptibilityCalculator::setU(std::complex<double> U){
430  this->U = U;
431  rpaSusceptibilityTree.clear();
432  rpaChargeSusceptibilityTree.clear();
433  rpaSpinSusceptibilityTree.clear();
434  interactionAmplitudesAreGenerated = false;
435 }
436 
437 inline void RPASusceptibilityCalculator::setUp(std::complex<double> Up){
438  this->Up = Up;
439  rpaSusceptibilityTree.clear();
440  rpaChargeSusceptibilityTree.clear();
441  rpaSpinSusceptibilityTree.clear();
442  interactionAmplitudesAreGenerated = false;
443 }
444 
445 inline void RPASusceptibilityCalculator::setJ(std::complex<double> J){
446  this->J = J;
447  rpaSusceptibilityTree.clear();
448  rpaChargeSusceptibilityTree.clear();
449  rpaSpinSusceptibilityTree.clear();
450  interactionAmplitudesAreGenerated = false;
451 }
452 
453 inline void RPASusceptibilityCalculator::setJp(std::complex<double> Jp){
454  this->Jp = Jp;
455  rpaSusceptibilityTree.clear();
456  rpaChargeSusceptibilityTree.clear();
457  rpaSpinSusceptibilityTree.clear();
458  interactionAmplitudesAreGenerated = false;
459 }
460 
461 }; //End of namespace TBTK
462 
463 #endif
Definition: RPASusceptibilityCalculator.h:38
EnergyType
Definition: RPASusceptibilityCalculator.h:69
RPASusceptibilityCalculator * createSlave()
void setJp(std::complex< double > Jp)
Definition: RPASusceptibilityCalculator.h:453
void setJ(std::complex< double > J)
Definition: RPASusceptibilityCalculator.h:445
Definition: LindhardSusceptibilityCalculator.h:35
void setUp(std::complex< double > Up)
Definition: RPASusceptibilityCalculator.h:437
void setEnergies(const std::vector< std::complex< double >> &energies)
Definition: RPASusceptibilityCalculator.h:313
Property container for density.
void loadSusceptibilities(const std::string &filename)
Definition: RPASusceptibilityCalculator.h:384
Definition: DualIndex.h:33
Algorithm
Definition: SusceptibilityCalculator.h:49
bool getEnergiesAreInversionSymmetric() const
Definition: RPASusceptibilityCalculator.h:335
bool getSusceptibilityIsSafeFromPoles() const
Definition: RPASusceptibilityCalculator.h:357
std::vector< std::complex< double > > calculateChargeRPASusceptibility(const DualIndex &kDual, const std::vector< int > &orbitalIndices)
Definition: SusceptibilityCalculator.h:41
void setEnergiesAreInversionSymmetric(bool energiesAreInversionSymmetric)
Definition: RPASusceptibilityCalculator.h:327
Definition: IndexedDataTree.h:43
Definition: MomentumSpaceContext.h:32
RPASusceptibilityCalculator(const MomentumSpaceContext &momentumSpaceContext, SusceptibilityCalculator::Algorithm algorithm=SusceptibilityCalculator::Algorithm::Lindhard)
Flexible physical index.
Definition: Index.h:69
std::vector< std::complex< double > > calculateRPASusceptibility(const DualIndex &kDual, const std::vector< int > &orbitalIndices)
Definition: ModelFactory.h:35
void saveSusceptibilities(const std::string &filename) const
Definition: RPASusceptibilityCalculator.h:378
void setInteractionAmplitudes(const std::vector< InteractionAmplitude > &interactionAmplitudes)
Definition: RPASusceptibilityCalculator.h:247
void precompute(unsigned int numWorkers=129)
void setU(std::complex< double > U)
Definition: RPASusceptibilityCalculator.h:429
void setEnergyType(EnergyType energyType)
Definition: RPASusceptibilityCalculator.h:279
void setSusceptibilityIsSafeFromPoles(bool susceptibilityIsSafeFromPoles)
Definition: RPASusceptibilityCalculator.h:340
Serializable wrapper of std::vector.
EnergyType getEnergyType() const
Definition: RPASusceptibilityCalculator.h:308
int & at(unsigned int n)
Definition: Index.h:351
std::vector< std::complex< double > > calculateSpinRPASusceptibility(const DualIndex &kDual, const std::vector< int > &orbitalIndices)