ERKALE
ERKALE - DFT from Hel
 All Classes Functions Variables Friends Pages
gto_fourier.h
1 /*
2  * This source code is part of
3  *
4  * E R K A L E
5  * -
6  * DFT from Hel
7  *
8  * Written by Susi Lehtola, 2010-2011
9  * Copyright (c) 2010-2011, Susi Lehtola
10  *
11  * This program is free software; you can redistribute it and/or
12  * modify it under the terms of the GNU General Public License
13  * as published by the Free Software Foundation; either version 2
14  * of the License, or (at your option) any later version.
15  */
16 
17 
18 
19 #ifndef ERKALE_FOURIER
20 #define ERKALE_FOURIER
21 
22 #include "../global.h"
23 #include <armadillo>
24 #include <complex>
25 #include <vector>
26 #include <cstddef>
27 
29 typedef struct {
31  std::complex<double> c;
33  int l;
34 } poly1d_t;
35 
37 bool operator<(const poly1d_t & lhs, const poly1d_t & rhs);
39 bool operator==(const poly1d_t & lhs, const poly1d_t & rhs);
40 
73  std::vector<poly1d_t> poly;
74 
79  FourierPoly_1D formpoly(int l, double zeta);
80 
81  public:
85  FourierPoly_1D(int l, double zeta);
88 
90  void addterm(const poly1d_t & term);
91 
93  FourierPoly_1D operator+(const FourierPoly_1D & rhs) const;
94 
96  size_t getN() const;
98  std::complex<double> getc(size_t i) const;
100  int getl(size_t i) const;
101 
103  void print() const;
104 
105  friend FourierPoly_1D operator*(std::complex<double> fac, const FourierPoly_1D & rhs);
106 };
107 
109 FourierPoly_1D operator*(std::complex<double> fac, const FourierPoly_1D & rhs);
110 
111 
112 /* Then, the full three-dimensional transform */
113 
114 
116 typedef struct {
118  std::complex<double> c;
120  int l;
122  int m;
124  int n;
126  double z;
127 } trans3d_t;
128 
130 bool operator<(const trans3d_t & lhs, const trans3d_t& rhs);
132 bool operator==(const trans3d_t & lhs, const trans3d_t& rhs);
133 
147 class GTO_Fourier {
149  std::vector<trans3d_t> trans;
150  public:
152  GTO_Fourier();
154  GTO_Fourier(int l, int m, int n, double zeta);
156  ~GTO_Fourier();
157 
159  void addterm(const trans3d_t & term);
160 
162  GTO_Fourier operator+(const GTO_Fourier & rhs) const;
164  GTO_Fourier & operator+=(const GTO_Fourier & rhs);
165 
167  std::vector<trans3d_t> get() const;
168 
170  void print() const;
171 
172  // Clean out the expansion
173  void clean();
174 
175  // Evaluate the expansion at p
176  std::complex<double> eval(double px, double py, double pz) const;
177 
178  friend GTO_Fourier operator*(std::complex<double> fac, const GTO_Fourier & rhs);
179  friend GTO_Fourier operator*(double fac, const GTO_Fourier & rhs);
180 };
181 
183 GTO_Fourier operator*(std::complex<double> fac, const GTO_Fourier & rhs);
185 GTO_Fourier operator*(double fac, const GTO_Fourier & rhs);
186 
187 class BasisSet;
188 
190 std::vector< std::vector<GTO_Fourier> > fourier_expand(const BasisSet & bas, std::vector< std::vector<size_t> > & idents);
191 
193 double eval_emd(const BasisSet & basis, const arma::cx_mat & P, const std::vector< std::vector<GTO_Fourier> > & fourier, const std::vector< std::vector<size_t> > & idents, double px, double py, double pz);
194 
195 #endif
friend FourierPoly_1D operator*(std::complex< double > fac, const FourierPoly_1D &rhs)
Multiply the polynomial with a complex factor.
Definition: gto_fourier.cpp:142
std::complex< double > c
Full polynomial is of the form .
Definition: gto_fourier.h:31
int l
Full polynomial is of the form .
Definition: gto_fourier.h:33
Computes the 1-dimensional Fourier polynomial needed for the Fourier transforms of Gaussian basis fun...
Definition: gto_fourier.h:71
~GTO_Fourier()
Destructor.
Definition: gto_fourier.cpp:223
std::vector< poly1d_t > poly
1-dimensional Fourier polynomial
Definition: gto_fourier.h:73
Compute Fourier transform of Gaussian Type Orbital.
Definition: gto_fourier.h:147
FourierPoly_1D()
Dummy constructor.
Definition: gto_fourier.cpp:31
Helper for 1D Fourier polynomials.
Definition: gto_fourier.h:29
friend GTO_Fourier operator*(std::complex< double > fac, const GTO_Fourier &rhs)
Scale Fourier transform of GTO by factor fac.
Definition: gto_fourier.cpp:281
std::vector< trans3d_t > trans
The terms of the Fourier transformed GTO.
Definition: gto_fourier.h:149
Fourier transform of GTO is of the form .
Definition: gto_fourier.h:116
void print() const
Print Fourier transform.
Definition: gto_fourier.cpp:299
GTO_Fourier()
Dummy constructor.
Definition: gto_fourier.cpp:182
void addterm(const trans3d_t &term)
Add a term in the contraction.
Definition: gto_fourier.cpp:226
GTO_Fourier operator+(const GTO_Fourier &rhs) const
Addition operator.
Definition: gto_fourier.cpp:247
FourierPoly_1D formpoly(int l, double zeta)
Definition: gto_fourier.cpp:43
double z
exp(-z*p^2)
Definition: gto_fourier.h:126
int getl(size_t i) const
Get the exponent of p in the i:th term.
Definition: gto_fourier.cpp:129
Basis set.
Definition: basis.h:187
std::complex< double > getc(size_t i) const
Get the i:th contraction coefficient.
Definition: gto_fourier.cpp:125
std::complex< double > c
Expansion coefficient.
Definition: gto_fourier.h:118
int l
px^l
Definition: gto_fourier.h:120
~FourierPoly_1D()
Destructor.
Definition: gto_fourier.cpp:87
GTO_Fourier & operator+=(const GTO_Fourier &rhs)
Addition operator.
Definition: gto_fourier.cpp:256
void addterm(const poly1d_t &term)
Add a term in the contraction.
Definition: gto_fourier.cpp:90
void print() const
Print polynomial.
Definition: gto_fourier.cpp:133
int n
pz^n
Definition: gto_fourier.h:124
size_t getN() const
Get number of terms in the polynomial.
Definition: gto_fourier.cpp:121
int m
py^m
Definition: gto_fourier.h:122
FourierPoly_1D operator+(const FourierPoly_1D &rhs) const
Addition operator.
Definition: gto_fourier.cpp:111