ERKALE
ERKALE - DFT from Hel
 All Classes Functions Variables Friends Pages
lmtrans.h
1 /*
2  * This source code is part of
3  *
4  * E R K A L E
5  * -
6  * HF/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 #ifndef ERKALE_LMTRANS
19 #define ERKALE_LMTRANS
20 
21 #include "../global.h"
22 #include <armadillo>
23 
24 #include "../gaunt.h"
25 #include "../lmgrid.h"
26 
28 typedef struct {
30  size_t i;
32  size_t f;
34  double q;
35 
37  std::vector<arma::cx_mat> itg;
38 } rad_int_t;
39 
41 typedef struct {
43  double q;
45  std::vector< std::vector<double> > jl;
46 } bessel_t;
47 
49 class lmtrans {
53  int lmax;
56 
62  arma::cx_mat radial_integral(size_t i, size_t f, int l, const bessel_t & bes) const;
63 
64  public:
66  lmtrans();
68  lmtrans(const arma::mat & C, const BasisSet & bas, const coords_t & cen, size_t Nrad=200, int lmax=5, int lquad=30);
70  ~lmtrans();
71 
76  rad_int_t compute_radial_integrals(size_t i, size_t f, const bessel_t & bes) const;
77 
81  bessel_t compute_bessel(double q) const;
82 
88  arma::cx_mat transition_amplitude(const rad_int_t & rad, double qx, double qy, double qz) const;
89 
96  std::vector<double> transition_velocity(size_t i, size_t f, const bessel_t & bes) const;
97 
99  void print_info() const;
100 
102  void write_prob(size_t o, const std::string & fname) const;
103 };
104 
105 
106 
107 #endif
std::vector< double > transition_velocity(size_t i, size_t f, const bessel_t &bes) const
Definition: lmtrans.cpp:180
arma::cx_mat transition_amplitude(const rad_int_t &rad, double qx, double qy, double qz) const
Definition: lmtrans.cpp:125
lmtrans()
Dummy constructor.
Definition: lmtrans.cpp:27
rad_int_t compute_radial_integrals(size_t i, size_t f, const bessel_t &bes) const
Definition: lmtrans.cpp:81
size_t f
Final state.
Definition: lmtrans.h:32
expansion_t exp
Orbital expansions.
Definition: lmtrans.h:51
~lmtrans()
Destructor.
Definition: lmtrans.cpp:44
Gaunt gaunt
Table of Gaunt coefficients.
Definition: lmtrans.h:55
void write_prob(size_t o, const std::string &fname) const
Write radial distribution of orbital o into file.
Definition: lmtrans.cpp:225
Bessel function stack.
Definition: lmtrans.h:41
double q
Momentum transfer.
Definition: lmtrans.h:34
Table of Gaunt coefficients.
Definition: gaunt.h:31
size_t i
Initial state.
Definition: lmtrans.h:30
std::vector< std::vector< double > > jl
Stack of Bessel functions, for different values l.
Definition: lmtrans.h:45
Basis set.
Definition: basis.h:187
Expansion of orbitals.
Definition: lmgrid.h:53
void print_info() const
Print information about orbitals.
Definition: lmtrans.cpp:169
int lmax
Maximum angular momentum in expansion.
Definition: lmtrans.h:53
double q
Momentum transfer.
Definition: lmtrans.h:43
Coordinates structure.
Definition: basis.h:50
bessel_t compute_bessel(double q) const
Definition: lmtrans.cpp:93
std::vector< arma::cx_mat > itg
Stack of radial integrals, for different values of l.
Definition: lmtrans.h:37
Radial integrals stack.
Definition: lmtrans.h:28
arma::cx_mat radial_integral(size_t i, size_t f, int l, const bessel_t &bes) const
Definition: lmtrans.cpp:47
Class for calculating XRS.
Definition: lmtrans.h:49