3 #include "HF/HartreeFock.hpp"
4 #include "MBPT/CorrelationPotential.hpp"
5 #include "Maths/Grid.hpp"
6 #include "Physics/AtomData.hpp"
7 #include "Physics/NuclearPotentials.hpp"
8 #include "Physics/PhysConst_constants.hpp"
9 #include "Physics/RadPot.hpp"
10 #include "Wavefunction/BSplineBasis.hpp"
11 #include "Wavefunction/DiracSpinor.hpp"
44 const std::string &run_label =
"");
48 const std::string &run_label =
"");
54 std::shared_ptr<const Grid> rgrid;
57 std::string m_run_label;
61 std::vector<DiracSpinor> m_valence{};
63 std::vector<DiracSpinor> m_basis{};
65 std::vector<DiracSpinor> m_spectrum{};
67 std::vector<double> m_vnuc{};
69 std::optional<HF::HartreeFock> m_HF{std::nullopt};
71 std::optional<MBPT::CorrelationPotential> m_Sigma{};
73 std::string m_core_string =
"";
74 std::string m_aboveFermi_core_string =
"";
76 std::vector<CI::PsiJPi> m_CIwfs{};
83 std::shared_ptr<const Grid>
grid_sptr()
const {
return rgrid; };
86 double alpha()
const {
return m_alpha; }
91 return (m_alpha * m_alpha / PhysConst::alpha2) - 1.0;
97 int Znuc()
const {
return m_nucleus.z(); }
99 int Anuc()
const {
return m_nucleus.a(); }
101 double get_rrms()
const {
return m_nucleus.r_rms(); }
104 const std::vector<DiracSpinor> &
core()
const {
105 static const auto empty = std::vector<DiracSpinor>{};
106 return m_HF ? m_HF->core() : empty;
110 const std::vector<DiracSpinor> &
valence()
const {
return m_valence; }
111 std::vector<DiracSpinor> &
valence() {
return m_valence; }
115 const std::vector<DiracSpinor> &
basis()
const {
return m_basis; }
116 std::vector<DiracSpinor> &
basis() {
return m_basis; }
119 const std::vector<DiracSpinor> &
spectrum()
const {
return m_spectrum; }
120 std::vector<DiracSpinor> &
spectrum() {
return m_spectrum; }
122 const std::vector<CI::PsiJPi> &CIwfs()
const {
return m_CIwfs; }
124 const CI::PsiJPi *CIwf(
int J,
int parity)
const {
125 for (
const auto &ci_wf : m_CIwfs) {
126 if (ci_wf.twoJ() == 2 * J && ci_wf.parity() == parity)
134 const std::vector<double> &
vnuc()
const {
return m_vnuc; }
142 std::vector<double>
vlocal(
int l = 0)
const;
146 std::vector<double>
Hmag(
int l = 0)
const;
153 const MBPT::CorrelationPotential *
Sigma()
const {
154 return m_Sigma ? &*m_Sigma :
nullptr;
156 MBPT::CorrelationPotential *
Sigma() {
return m_Sigma ? &*m_Sigma :
nullptr; }
188 ", Z=" + std::to_string(m_nucleus.z()) +
189 " A=" + std::to_string(m_nucleus.a());
215 void printValence(
const std::vector<DiracSpinor> &tmp_orbitals = {})
const;
218 void printBasis(
const std::vector<DiracSpinor> &the_basis)
const;
222 bool isInValence(
int n,
int k)
const;
234 void set_HF(
const std::string &method =
"HartreeFock",
235 const double x_Breit = 0.0,
const std::string &in_core =
"",
236 double eps_HF = 1.0e-13,
bool print =
true);
242 void solve_core(
const std::string &method,
const double x_Breit = 0.0,
243 const std::string &in_core =
"",
double eps_HF = 1.0e-13,
248 const bool print =
true);
256 const std::vector<double> &fit_energies);
260 const std::vector<double> &x_spd,
261 bool do_readwrite =
true,
bool print =
true);
268 void formBasis(
const SplineBasis::Parameters ¶ms);
271 void formSpectrum(
const SplineBasis::Parameters ¶ms);
274 void formSigma(
int nmin_core = 1,
int nmin_core_F = 1,
double r0 = 1.0e-4,
275 double rmax = 30.0,
int stride = 4,
bool each_valence =
false,
276 bool include_G =
false,
bool include_Breit =
false,
277 int n_max_breit = 0,
const std::vector<double> &lambdas = {},
278 const std::vector<double> &fk = {},
279 const std::vector<double> &etak = {},
280 const std::string &in_fname =
"",
281 const std::string &out_fname =
"",
bool FeynmanQ =
false,
282 bool ScreeningQ =
false,
bool holeParticleQ =
false,
283 int lmax = 6,
double omre = -0.2,
double w0 = 0.01,
285 const std::optional<IO::InputBlock> &ek = std::nullopt);
289 void copySigma(
const MBPT::CorrelationPotential *
const Sigma) {
290 if (
Sigma !=
nullptr)
300 m_HF->vnuc() = v_new;
323 double H0ab_impl(
const DiracSpinor &Fa, std::vector<double> dga,
324 const DiracSpinor &Fb, std::vector<double> dgb)
const;
327 std::vector<DiracSpinor> determineCore(
const std::string &str_core_in);
328 bool isInAboveFermiCore(
int n,
int k)
const;
Stores the CI Solutions for given J and parity (only two-electron).
Definition: CSF.hpp:62
Stores radial Dirac spinor: F_nk = (f, g)
Definition: DiracSpinor.hpp:41
Holds grid, including type + Jacobian (dr/du)
Definition: Grid.hpp:31
Solves relativistic Hartree-Fock equations for core and valence. Optionally includes Breit and QED ef...
Definition: HartreeFock.hpp:70
Stores set of nuclear parameters (all radii in fm)
Definition: NuclearPotentials.hpp:20
Class holds Flambaum-Ginges QED Radiative Potential.
Definition: RadPot.hpp:15
Stores Wavefunction (set of valence orbitals, grid, HF etc.)
Definition: Wavefunction.hpp:36
void fitSigma_hfBrueckner(const std::string &valence_list, const std::vector< double > &fit_energies)
First, fits Sigma to energies, then forms fitted Brueckner orbitals.
Definition: Wavefunction.cpp:614
void solve_valence(const std::string &in_valence_str="", const bool print=true)
Performs hartree-Fock procedure for valence: note: poplulates valnece.
Definition: Wavefunction.cpp:185
void formBasis(const SplineBasis::Parameters ¶ms)
Calculates + populates basis [see BSplineBasis].
Definition: Wavefunction.cpp:467
void solve_core(bool print=true)
Performs hartree-Fock procedure for core.
Definition: Wavefunction.cpp:164
void printValence(const std::vector< DiracSpinor > &tmp_orbitals={}) const
Prints table of valence orbitals + energies etc.
Definition: Wavefunction.cpp:406
const std::vector< double > & vnuc() const
Nuclear potential. Only provide const version, since HF and WF version of vnuc must be kept in sync.
Definition: Wavefunction.hpp:134
std::shared_ptr< const Grid > grid_sptr() const
Copy of shared_ptr to grid [shared resource] - used when we want to construct a new object that share...
Definition: Wavefunction.hpp:83
double coreEnergyHF() const
Calculates HF core energy (doesn't include magnetic QED?)
Definition: Wavefunction.cpp:177
std::tuple< double, double > lminmax_core_range(int l, double eps=0.0) const
Returns [min,max] r values for which the core density (given l) is larger than cutoff (= eps*max_valu...
Definition: Wavefunction.cpp:348
const DiracSpinor * getState(int n, int k) const
Finds requested state; returns nullptr if not found.
Definition: Wavefunction.cpp:289
int Znuc() const
Nuclear charge, Z.
Definition: Wavefunction.hpp:97
double energy_gap() const
Energy gap between lowest valence + highest core state: e(v) - e(c) [should be positive].
Definition: Wavefunction.cpp:335
const std::vector< DiracSpinor > & spectrum() const
Sprectrum: like basis, but includes correlations.
Definition: Wavefunction.hpp:119
const std::vector< DiracSpinor > & valence() const
Valence orbitals (HF or Brueckner orbitals)
Definition: Wavefunction.hpp:110
std::string identity() const
Atomic symbol, including core ionisation degree and run_label.
Definition: Wavefunction.cpp:686
std::string coreConfiguration() const
Returns full core configuration.
Definition: Wavefunction.hpp:178
int Zion() const
Effective charge (for core) = Z-N_core.
Definition: Wavefunction.hpp:209
std::vector< double > vlocal(int l=0) const
Local part of potential, e.g., Vl = Vnuc + Vdir + Vrad_el(l) - can be l-dependent....
Definition: Wavefunction.cpp:667
std::vector< double > Hmag(int l=0) const
QED Magnetic form factor. May return empty vector. Not typically l-dependent, but may be in future....
Definition: Wavefunction.cpp:672
std::vector< double > coreDensity() const
Calculates rho(r) = sum_c psi^2(r) for core states, c={n,k,m}.
Definition: Wavefunction.cpp:455
int ion_degree(int num_val) const
0 for neutral, 1 for singly-ionised etc.
Definition: Wavefunction.hpp:201
int Anuc() const
Nuclear mass number, A.
Definition: Wavefunction.hpp:99
const Grid & grid() const
Returns a const reference to the radial grid.
Definition: Wavefunction.hpp:80
std::string atomicSymbol() const
e.g., "Cs"
Definition: Wavefunction.hpp:193
void radiativePotential(QED::RadPot::Scale s, double rcut, double scale_rN, const std::vector< double > &x_spd, bool do_readwrite=true, bool print=true)
OLD: deprecated.
Definition: Wavefunction.cpp:219
std::string ion_symbol(int num_val) const
I for neutral, II for singly-ionised etc.
Definition: Wavefunction.hpp:204
double get_rrms() const
Nuclear rms charge radii, in fm (femptometres)
Definition: Wavefunction.hpp:101
const std::vector< DiracSpinor > & basis() const
Basis, eigenstates of HF potential. Used for MBPT. Includes Breit and QED (if they are included),...
Definition: Wavefunction.hpp:115
std::string coreConfiguration_nice() const
Returns core configuration, in nice output notation.
Definition: Wavefunction.hpp:181
void formSpectrum(const SplineBasis::Parameters ¶ms)
Calculates + populates Spectrum [see BSplineBasis].
Definition: Wavefunction.cpp:486
void formSigma(int nmin_core=1, int nmin_core_F=1, double r0=1.0e-4, double rmax=30.0, int stride=4, bool each_valence=false, bool include_G=false, bool include_Breit=false, int n_max_breit=0, const std::vector< double > &lambdas={}, const std::vector< double > &fk={}, const std::vector< double > &etak={}, const std::string &in_fname="", const std::string &out_fname="", bool FeynmanQ=false, bool ScreeningQ=false, bool holeParticleQ=false, int lmax=6, double omre=-0.2, double w0=0.01, double wratio=1.5, const std::optional< IO::InputBlock > &ek=std::nullopt)
Forms + stores correlation potential Sigma.
Definition: Wavefunction.cpp:506
const QED::RadPot * vrad() const
Pointer to QED radiative potnential. May be nullptr.
Definition: Wavefunction.hpp:149
void printCore() const
Prints table of core orbitals + energies etc.
Definition: Wavefunction.cpp:374
bool isInCore(int n, int k) const
Check if a state is in the core (or valence) list.
Definition: Wavefunction.cpp:260
std::string atom() const
String of atom info (e.g., "Cs, Z=55, A=133")
Definition: Wavefunction.hpp:186
const std::vector< DiracSpinor > & core() const
Core orbitals (frozen HF core)
Definition: Wavefunction.hpp:104
double alpha() const
Local value of fine-structure constant.
Definition: Wavefunction.hpp:86
double H0ab(const DiracSpinor &Fa, const DiracSpinor &Fb) const
Returns <a|H|b> for Hamiltonian H (inludes Rad.pot, NOT sigma, Breit, or exchange!...
Definition: Wavefunction.cpp:717
const Nuclear::Nucleus & nucleus() const
Returns Nuclear::nucleus object (contains nuc. parameters)
Definition: Wavefunction.hpp:95
void update_Vnuc(const std::vector< double > &v_new)
Allows extra potential to be added to Vnuc (updates both in Wavefunction.
Definition: Wavefunction.hpp:296
void set_HF(const std::string &method="HartreeFock", const double x_Breit=0.0, const std::string &in_core="", double eps_HF=1.0e-13, bool print=true)
Initialises HF object and populates core orbitals (does not solve HF equations)
Definition: Wavefunction.cpp:122
const MBPT::CorrelationPotential * Sigma() const
Returns ptr to (const) Correlation Potential, Sigma.
Definition: Wavefunction.hpp:153
void hartreeFockBrueckner(const bool print=true)
Forms Bruckner valence orbitals: (H_hf + Sigma)|nk> = e|nk>. Replaces existing valence states.
Definition: Wavefunction.cpp:604
int Ncore() const
Number of electrons in the core.
Definition: Wavefunction.cpp:678
double FermiLevel() const
Returns energy location of the "Fermi Level", - energy half way between core/valence....
Definition: Wavefunction.cpp:316
void printBasis(const std::vector< DiracSpinor > &the_basis) const
Prints table of Basis/Spectrum orbitals, compares to HF orbitals.
Definition: Wavefunction.cpp:432
const HF::HartreeFock * vHF() const
Returns ptr to Hartree Fock (class)
Definition: Wavefunction.hpp:137
double dalpha2() const
Variation in alpha^2 : x = (alpha/alpha_0)^2 - 1.
Definition: Wavefunction.hpp:89
std::string niceCoreOutput(const std::string &full_core)
Given a full electron config., returns nicer format by recognising nobel gas.
Definition: AtomData.cpp:198
std::string atomicSymbol(int Z)
e.g., 55 -> "Cs"
Definition: AtomData.cpp:64
Constucts of spinor/orbital basis using B-splines (DKB/Reno/Derevianko-Beloy method)
Definition: BSplineBasis.cpp:20
std::string int_to_roman(int a)
Converts integer, a, to Roman Numerals. Assumed that |a|<=4000.
Definition: String.hpp:282
Parmaters used to create Grid.
Definition: Grid.hpp:9
Scale factors for Uehling, high, low, magnetic, Wickman-Kroll.
Definition: RadPot.hpp:20