Many-body perturbation theory. More...

Classes
class	Feynman
	Class to construct Feynman diagrams, Green's functions and polarisation op. More...

class	Goldstone
	Class to construct Feynman diagrams, Green's functions and polarisation op. More...

class	RDMatrix

class	StructureRad
	Calculates Structure Radiation + Normalisation of states, using diagram method. More...

Typedefs
using	GMatrix = RDMatrix< double >

using	ComplexGMatrix = RDMatrix< std::complex< double > >

using	ComplexDouble = std::complex< double >

Enumerations
enum class	SigmaMethod { Goldstone , Feynman }

enum class	HoleParticle { exclude , include , include_k0 }
	Options for including hole-particle interaction. include mean all k; include_k0 means k=0 term only.

enum class	Screening { exclude , include , only }
	Options for including Screening.

enum class	GreenStates { both , core , excited }
	Which states to include in Green's function:

enum class	Denominators { RS , BW }
	Type of energy demoninators: Rayleigh-Schrodinger (RS), Brillouin-Wigner (BW). (not exact)

Functions
double	Lkmnij (int k, const DiracSpinor &m, const DiracSpinor &n, const DiracSpinor &i, const DiracSpinor &j, const Coulomb::QkTable &qk, const std::vector< DiracSpinor > &core, const std::vector< DiracSpinor > &excited, bool include_L4, const Angular::SixJTable &SJ, const Coulomb::LkTable *const Lk=nullptr, const std::vector< double > &fk={})
	Calculates ladder integral, L^k_mnab.

double	L1 (int k, const DiracSpinor &m, const DiracSpinor &n, const DiracSpinor &i, const DiracSpinor &j, const Coulomb::QkTable &qk, const std::vector< DiracSpinor > &excited, const Angular::SixJTable &SJ, const Coulomb::LkTable *const Lk=nullptr, const std::vector< double > &fk={})
	Ladder integral, L^k_mnij := L1_mnij + L2_mnij + L2_nmji.

double	L4 (int k, const DiracSpinor &m, const DiracSpinor &n, const DiracSpinor &i, const DiracSpinor &j, const Coulomb::QkTable &qk, const std::vector< DiracSpinor > &core, const Angular::SixJTable &SJ, const Coulomb::LkTable *const Lk, const std::vector< double > &fk)

double	L2 (int k, const DiracSpinor &m, const DiracSpinor &n, const DiracSpinor &i, const DiracSpinor &j, const Coulomb::QkTable &qk, const std::vector< DiracSpinor > &core, const std::vector< DiracSpinor > &excited, const Angular::SixJTable &SJ, const Coulomb::LkTable *const Lk=nullptr, const std::vector< double > &fk={})
	Ladder integral, L^k_mnab := L1_mnij + L2_mnij + L3_nmji, L3_mnij = L2_nmji.

void	fill_Lk_mnib (Coulomb::LkTable lk, const Coulomb::QkTable &qk, const std::vector< DiracSpinor > &excited, const std::vector< DiracSpinor > &core, const std::vector< DiracSpinor > &i_orbs, bool include_L4, const Angular::SixJTable &sjt, const Coulomb::LkTable const lk_prev=nullptr, bool print_progbar=true, const std::vector< double > &fk={})
	Fills Lk matrix.

double	L3 (int k, const DiracSpinor &m, const DiracSpinor &n, const DiracSpinor &i, const DiracSpinor &j, const Coulomb::QkTable &qk, const std::vector< DiracSpinor > &core, const std::vector< DiracSpinor > &excited, const Angular::SixJTable &SJ, const Coulomb::LkTable *const Lk=nullptr, const std::vector< double > &fk={})

template<typename Qintegrals , typename QorLintegrals >
double	de_valence (const DiracSpinor &v, const Qintegrals &qk, const QorLintegrals &lk, const std::vector< DiracSpinor > &core, const std::vector< DiracSpinor > &excited, const std::vector< double > &fk={}, const std::vector< double > &etak={})
	Calculate energy shift (either ladder, or sigma2) for valence.

template<typename Qintegrals , typename QorLintegrals >
double	de_core (const Qintegrals &qk, const QorLintegrals &lk, const std::vector< DiracSpinor > &core, const std::vector< DiracSpinor > &excited)
	Calculate energy shift (either ladder, or sigma2) for CORE.

template<typename T >
bool	equal (const RDMatrix< T > &lhs, const RDMatrix< T > &rhs)
	Checks if two matrix's are equal (to within parts in 10^12)

template<typename T >
double	max_element (const RDMatrix< T > &a)
	returns maximum element (by abs)

template<typename T >
double	max_delta (const RDMatrix< T > &a, const RDMatrix< T > &b)
	returns maximum difference (abs) between two matrixs

template<typename T >
double	max_epsilon (const RDMatrix< T > &a, const RDMatrix< T > &b)
	returns maximum relative diference [aij-bij/(aij+bij)] (abs) between two matrices

std::pair< std::vector< DiracSpinor >, std::vector< DiracSpinor > >	split_basis (const std::vector< DiracSpinor > &basis, double E_Fermi, int min_n_core=1, int max_n_excited=999)
	Splits the basis into the core (holes) and excited states.

double	e_bar (int kappa_v, const std::vector< DiracSpinor > &excited)
	Returns energy of first state in excited that matches given kappa.

bool	Sk_vwxy_SR (int k, const DiracSpinor &v, const DiracSpinor &w, const DiracSpinor &x, const DiracSpinor &y)
	Selection rule for Sk_vwxy (differs from Qk_vwxy due to parity)

int	number_below_Fermi (const DiracSpinor &i, const DiracSpinor &j, const DiracSpinor &k, const DiracSpinor &l, double eFermi)
	Returns number of orbitals that are below Fermi level. Used for Qk selection.

std::pair< int, int >	k_minmax_S (const DiracSpinor &v, const DiracSpinor &w, const DiracSpinor &x, const DiracSpinor &y)
	Minimum/maximum k allowed by selectrion rules for Sk_vwxy. Cannot +=2.

std::pair< int, int >	k_minmax_S (int twojv, int twojw, int twojx, int twojy)

double	Sk_vwxy (int k, const DiracSpinor &v, const DiracSpinor &w, const DiracSpinor &x, const DiracSpinor &y, const Coulomb::QkTable &qk, const std::vector< DiracSpinor > &core, const std::vector< DiracSpinor > &excited, const Angular::SixJTable &SixJ, Denominators denominators=Denominators::BW)
	Reduced two-body Sigma (2nd order correlation) operator. Sum of 6 diagrams.

template<class CoulombIntegral >
double	Sigma_vw (const DiracSpinor &v, const DiracSpinor &w, const CoulombIntegral &qk, const std::vector< DiracSpinor > &core, const std::vector< DiracSpinor > &excited, int max_l_internal=99, std::optional< double > ev=std::nullopt)
	Matrix element of 1-body Sigma (2nd-order correlation) operator; de_v = <v\|Sigma\|v>.

Variables
const auto	vroot = [](auto x) { return std::sqrt(x); }

Detailed Description

Many-body perturbation theory.

Function Documentation

◆ de_core()

template<typename Qintegrals , typename QorLintegrals >

double MBPT::de_core	(	const Qintegrals &	qk,
		const QorLintegrals &	lk,
		const std::vector< DiracSpinor > &	core,
		const std::vector< DiracSpinor > &	excited
	)

Calculate energy shift (either ladder, or sigma2) for CORE.

lk may be regular Coulomb integrals [in which case this returns MBPT(2) correction], or Ladder diagrams [in which case this returns the ladder diagram correction]

◆ de_valence()

template<typename Qintegrals , typename QorLintegrals >

double MBPT::de_valence	(	const DiracSpinor &	v,
		const Qintegrals &	qk,
		const QorLintegrals &	lk,
		const std::vector< DiracSpinor > &	core,
		const std::vector< DiracSpinor > &	excited,
		const std::vector< double > &	fk = `{}`,
		const std::vector< double > &	etak = `{}`
	)

Calculate energy shift (either ladder, or sigma2) for valence.

lk may be regular Coulomb integrals [in which case this returns MBPT(2) correction], or Ladder diagrams [in which case this returns the ladder diagram correction]

◆ L1()

double MBPT::L1	(	int	k,
		const DiracSpinor &	m,
		const DiracSpinor &	n,
		const DiracSpinor &	i,
		const DiracSpinor &	j,
		const Coulomb::QkTable &	qk,
		const std::vector< DiracSpinor > &	excited,
		const Angular::SixJTable &	SJ,
		const Coulomb::LkTable *const	Lk = `nullptr`,
		const std::vector< double > &	fk = `{}`
	)

Ladder integral, L^k_mnij := L1_mnij + L2_mnij + L2_nmji.

L1^k_mnij = sum_{rs,ul} A^{kul}_mnrsij * Q^u_mnrs * (Q+L)^l_rsij / (e_ij - e_rs)

A^{kul}_mnrsij = (-1)^{m+n+r+s+i+j+1} * [k] * {m,i,k;l,u,r} * {n,j,k;l,u,s}

◆ L2()

double MBPT::L2	(	int	k,
		const DiracSpinor &	m,
		const DiracSpinor &	n,
		const DiracSpinor &	i,
		const DiracSpinor &	j,
		const Coulomb::QkTable &	qk,
		const std::vector< DiracSpinor > &	core,
		const std::vector< DiracSpinor > &	excited,
		const Angular::SixJTable &	SJ,
		const Coulomb::LkTable *const	Lk = `nullptr`,
		const std::vector< double > &	fk = `{}`
	)

Ladder integral, L^k_mnab := L1_mnij + L2_mnij + L3_nmji, L3_mnij = L2_nmji.

L2^k_mnij = sum_{rc,ul} (-1)^{k+u+l+1} A^{klu}_mjcrin Q^u_cnir * (Q+L)^l_mrcj / (e_cj - e_mr)

◆ Lkmnij()

double MBPT::Lkmnij	(	int	k,
		const DiracSpinor &	m,
		const DiracSpinor &	n,
		const DiracSpinor &	i,
		const DiracSpinor &	j,
		const Coulomb::QkTable &	qk,
		const std::vector< DiracSpinor > &	core,
		const std::vector< DiracSpinor > &	excited,
		bool	include_L4,
		const Angular::SixJTable &	SJ,
		const Coulomb::LkTable *const	Lk = `nullptr`,
		const std::vector< double > &	fk = `{}`
	)

Calculates ladder integral, L^k_mnab.

Lk pointer is pointer to previous iteration of Lk. fk is optional vector of screening factors

◆ Sigma_vw()

template<class CoulombIntegral >

double MBPT::Sigma_vw	(	const DiracSpinor &	v,
		const DiracSpinor &	w,
		const CoulombIntegral &	qk,
		const std::vector< DiracSpinor > &	core,
		const std::vector< DiracSpinor > &	excited,
		int	max_l_internal = `99`,
		std::optional< double >	ev = `std::nullopt`
	)

Matrix element of 1-body Sigma (2nd-order correlation) operator; de_v = <v|Sigma|v>.

Matrix element of 1-body Sigma (2nd-order correlation) operator; de_v = <v|Sigma|v>. qk (CoulombIntegral) may be YkTable or QkTable.

ev is optional energy at which Sigma is calculated.
If ev is not given, uses 0.5*(ev+ew)
max_l_internal is largest angular momentum l to include in summation in internal lines; only used for tests.
qk (CoulombIntegral) may be YkTable or QkTable.

◆ Sk_vwxy()

double MBPT::Sk_vwxy	(	int	k,
		const DiracSpinor &	v,
		const DiracSpinor &	w,
		const DiracSpinor &	x,
		const DiracSpinor &	y,
		const Coulomb::QkTable &	qk,
		const std::vector< DiracSpinor > &	core,
		const std::vector< DiracSpinor > &	excited,
		const Angular::SixJTable &	SixJ,
		Denominators	denominators = `Denominators::BW`
	)

Reduced two-body Sigma (2nd order correlation) operator. Sum of 6 diagrams.

Note: these have fewer symmetries to Q^k; S_vwxy = S_xyvw

◆ split_basis()

std::pair< std::vector< DiracSpinor >, std::vector< DiracSpinor > > MBPT::split_basis	(	const std::vector< DiracSpinor > &	basis,
		double	E_Fermi,
		int	min_n_core = `1`,
		int	max_n_excited = `999`
	)

Splits the basis into the core (holes) and excited states.

States with energy below E_Fermi are considered core/holes; only core states with n>=min_n_core, and excited states with n<=max_n_excited are kept.

Classes

Typedefs

Enumerations

Functions

Variables