ok

Uses TDHF to include core-polarisation (RPA) corrections to matrix elements of an external field operator.

Solves the TDHF equations for each core orbital \( \phi_a \),

\[ (h_{\rm HF} - \en_a \mp \omega)\varphi^a_\pm = -(t_\pm + \delta V_\pm - \delta\en^a_\pm)\phi_a, \]

self-consistently to determine \( \delta V_\pm \). See the ExternalField namespace documentation for full physics description.

Each core orbital acquires both a forward ( \( \varphi^a_+ \), stored as X) and backward ( \( \varphi^a_- \), stored as Y) correction. Both contribute to \( \delta V_\pm \):

\[ \delta V_\pm \phi_i = \sum_a^{\rm core} \left[ \matel{\phi_a}{Q}{\varphi^a_+}\phi_i - \matel{\phi_a}{Q}{\phi_i}\varphi^a_+ + \matel{\varphi^a_-}{Q}{\phi_a}\phi_i - \matel{\varphi^a_-}{Q}{\phi_i}\phi_a \right]. \]

It is via \( \delta V_\pm \) that the \( e^{\pm i\omega t} \) terms are mixed.

Construction: Requires a pointer to the operator ( \( t_+ \), at \( +|\omega| \)) and a const HF::HartreeFock object. The user is responsible for updating the operator frequency externally before calling solve_core(). The backward operator \( t_- \) (at \( -|\omega| \)) is \( t_+ \) times TensorOperator::h_minus_sign(): -1 for an operator odd in the frequency (E1v), +1 otherwise.

Usage: solve_core (omega) solves the TDHF equations for a given frequency. dV (Fa, Fb) then returns the RPA correction to the reduced matrix element \( \redmatel{a}{\delta V}{b} \).

#include <TDHF.hpp>

Inheritance diagram for ExternalField::TDHF:

Public Member Functions
	TDHF (const DiracOperator::TensorOperator const h_plus, const HF::HartreeFock const hf)
	Constructs TDHF for operator h.

virtual void	solve_core (double omega, int max_its=100, bool print=true) override
	Solves TDHF equations self-consistently for core electrons at frequency omega.

void	solve_core_damped (double omega, int max_its=128, bool print=true)
	As solve_core(), by damped iteration only (no fallback).

void	solve_core_anderson (double omega, int max_its=128, bool print=true)
	As solve_core(), by Anderson mixing only (no fallback).

void	set_anderson (bool andersonQ)
	Method of solve_core(): Anderson mixing if true (default), damped iteration if false; the other is the fallback.

virtual Method	method () const override
	Returns RPA method.

virtual void	clear () override
	Clears the internal state back to pre solve_core()

double	dV (const DiracSpinor &Fa, const DiracSpinor &Fb, bool conj=false) const override
	Returns reduced matrix element <a\|\|dV\|\|b>, or the conjugate <a\|\|dV^dag\|\|b> if conj=true.

DiracSpinor	dV_rhs (int kappa_n, const DiracSpinor &Fm, bool conj=false) const override
	Returns [dV_pm * phi_m]_kappa: RHS of TDHF eq., projected onto kappa (see namespace doc)

const std::vector< DiracSpinor > &	get_dPsis (const DiracSpinor &Fc, dPsiType XorY) const
	Returns const ref to dPsi orbitals for given core orbital Fc.

const DiracSpinor &	get_dPsi_x (const DiracSpinor &Fc, dPsiType XorY, const int kappa_x) const
	Returns const ref to dPsi orbital of given kappa.

DiracSpinor	solve_dPsi (const DiracSpinor &Fv, const double omega, dPsiType XorY, const int kappa_beta, const MBPT::CorrelationPotential *const Sigma=nullptr, StateType st=StateType::ket, bool incl_dV=true) const
	Forms the mixed states phi^v_+- for valence state Fv (including core pol.): single kappa channel.

std::vector< DiracSpinor >	solve_dPsis (const DiracSpinor &Fv, const double omega, dPsiType XorY, const MBPT::CorrelationPotential *const Sigma=nullptr, StateType st=StateType::ket, bool incl_dV=true) const
	Forms the mixed states phi^v_+- for all kappa channels; see solve_dPsi.

TDHF &	operator= (const TDHF &)=delete

	TDHF (const TDHF &)=default

Public Member Functions inherited from ExternalField::CorePolarisation
double	last_eps () const
	Returns eps (convergance) of last solve_core run.

double	last_its () const
	Returns its (# of iterations) of last solve_core run.

double	last_omega () const
	Returns omega (frequency) of last solve_core run.

int	rank () const
	Rank of the operator.

int	parity () const
	Parity of the operator.

bool	imagQ () const
	Returns true if the operator is imaginary.

double &	eps_target ()
	Convergance target.

double	eps_target () const
	Convergance target.

double	eta () const
	Damping factor; 0 means no damping. Must have 0 <= eta < 1.

void	set_eta (double eta)
	Set/update damping factor; 0 means no damping. Must have 0 <= eta < 1.

CorePolarisation &	operator= (const CorePolarisation &)=delete

	CorePolarisation (const CorePolarisation &)=default

Protected Member Functions
std::vector< std::vector< DiracSpinor > >	form_hFcore (const DiracOperator::TensorOperator *h) const

std::vector< std::vector< DiracSpinor > >	form_hFcore_minus (const std::vector< std::vector< DiracSpinor > > &hFcore) const

DiracSpinor	dV_rhs_sets (int kappa_n, const DiracSpinor &Fa, bool conj, const std::vector< std::vector< DiracSpinor > > &X, const std::vector< std::vector< DiracSpinor > > &Y) const

Protected Member Functions inherited from ExternalField::CorePolarisation
	CorePolarisation (const DiracOperator::TensorOperator *const h)

Protected Attributes
std::vector< std::vector< DiracSpinor > >	m_X {}

std::vector< std::vector< DiracSpinor > >	m_Y {}

std::vector< std::vector< DiracSpinor > >	m_hFcore {}

std::vector< std::vector< DiracSpinor > >	m_hFcore_minus {}

const HF::HartreeFock *const	p_hf

const std::vector< DiracSpinor >	m_core

const double	m_alpha

const HF::Breit *const	p_VBr

bool	m_eps_sqrt {false}

bool	m_use_anderson {true}

bool	m_use_anderson_ms {true}

Protected Attributes inherited from ExternalField::CorePolarisation
const DiracOperator::TensorOperator *	m_h

double	m_core_eps {1.0}

int	m_core_its {0}

double	m_core_omega {0.0}

int	m_rank

int	m_pi

bool	m_imag

double	m_eta {0.35}

double	m_eps {1.0e-15}

Constructor & Destructor Documentation

◆ TDHF()

ExternalField::TDHF::TDHF	(	const DiracOperator::TensorOperator *const	h_plus,
		const HF::HartreeFock *const	hf
	)

Constructs TDHF for operator h.

Parameters

h_plus	The operator \( t_+ \); must be set to positive frequency before each call to solve_core().
hf	HF::HartreeFock object defining the core.

Member Function Documentation

◆ solve_core()

void ExternalField::TDHF::solve_core	(	double	omega,
		int	max_its = `100`,
		bool	print = `true`
	)

overridevirtual

Solves TDHF equations self-consistently for core electrons at frequency omega.

Iterates the TDHF equations until \( \delta V_\pm \) converges to within eps_target(), or max_its iterations have been performed. Can be re-run at a different frequency without restarting from scratch.

Solved by solve_core_anderson() (default) or solve_core_damped(). If that does not converge, the solve is started again from scratch with the other method, and the better result is kept.

Parameters

omega	External-field frequency \( \omega \) in atomic units.
max_its	Maximum number of iterations. Set to 1 to get the first-order correction (no damping on first iteration).
print	If true, write convergence progress to screen.

Note: Solved at \( |\omega| \): both \( \varphi^a_+ \) and \( \varphi^a_- \) are found simultaneously (see class docs), so a single call covers both components, \( \delta V(\pm|\omega|) \); see dV().

Note: Does not update the frequency of the operator itself; for frequency-dependent operators, update the operator frequency externally before calling.

Implements ExternalField::CorePolarisation.

Reimplemented in ExternalField::TDHFbasis, and ExternalField::TDHFcntm.

◆ solve_core_damped()

void ExternalField::TDHF::solve_core_damped	(	double	omega,
		int	max_its = `128`,
		bool	print = `true`
	)

As solve_core(), by damped iteration only (no fallback).

Each iteration mixes the new corrections with the previous ones, with damping factor eta().

Note: nb: if the iterations stall, the solve is restarted from scratch with different damping (weaker, then stronger), up to a few times, each with a further max_its iterations; the best attempt is kept. eta() itself is not changed.

◆ solve_core_anderson()

void ExternalField::TDHF::solve_core_anderson	(	double	omega,
		int	max_its = `128`,
		bool	print = `true`
	)

As solve_core(), by Anderson mixing only (no fallback).

Self-consistency is driven by Anderson mixing (anderson_coefficients) of all the corrections together, X and Y. Each iteration applies the TDHF map undamped, then extrapolates over the recent iterates. The TDHF equations are linear in the corrections, so this is equivalent to truncated GMRES: there is no damping parameter, and it converges wherever the equations are non-singular.

◆ set_anderson()

void ExternalField::TDHF::set_anderson ( bool andersonQ )

inline

Method of solve_core(): Anderson mixing if true (default), damped iteration if false; the other is the fallback.

◆ method()

virtual Method ExternalField::TDHF::method ( ) const

inlineoverridevirtual

Returns RPA method.

Implements ExternalField::CorePolarisation.

Reimplemented in ExternalField::TDHFbasis, and ExternalField::TDHFcntm.

◆ clear()

void ExternalField::TDHF::clear ( )

overridevirtual

Clears the internal state back to pre solve_core()

Implements ExternalField::CorePolarisation.

Reimplemented in ExternalField::TDHFcntm.

◆ dV()

double ExternalField::TDHF::dV	(	const DiracSpinor &	Fa,
		const DiracSpinor &	Fb,
		bool	conj = `false`
	)		const

overridevirtual

Returns reduced matrix element <a||dV||b>, or the conjugate <a||dV^dag||b> if conj=true.

Implements ExternalField::CorePolarisation.

Reimplemented in ExternalField::TDHFcntm.

◆ dV_rhs()

DiracSpinor ExternalField::TDHF::dV_rhs	(	int	kappa,
		const DiracSpinor &	Fm,
		bool	conj = `false`
	)		const

overridevirtual

Returns [dV_pm * phi_m]_kappa: RHS of TDHF eq., projected onto kappa (see namespace doc)

Reimplemented from ExternalField::CorePolarisation.

◆ get_dPsis()

const std::vector< DiracSpinor > & ExternalField::TDHF::get_dPsis	(	const DiracSpinor &	Fc,
		dPsiType	XorY
	)		const

Returns const ref to dPsi orbitals for given core orbital Fc.

◆ get_dPsi_x()

const DiracSpinor & ExternalField::TDHF::get_dPsi_x	(	const DiracSpinor &	Fc,
		dPsiType	XorY,
		const int	kappa_x
	)		const

Returns const ref to dPsi orbital of given kappa.

◆ solve_dPsi()

DiracSpinor ExternalField::TDHF::solve_dPsi	(	const DiracSpinor &	Fv,
		const double	omega,
		dPsiType	XorY,
		const int	kappa_beta,
		const MBPT::CorrelationPotential *const	Sigma = `nullptr`,
		StateType	st = `StateType::ket`,
		bool	incl_dV = `true`
	)		const

Forms the mixed states phi^v_+- for valence state Fv (including core pol.): single kappa channel.

Solves

\[ (h_{\rm HF} + \Sigma - \en_v - \omega)\varphi^v_+ = -(t + \delta V - \delta\en^v)\phi_v \]

or

\[ (h_{\rm HF} + \Sigma - \en_v + \omega)\varphi^v_- = -(t^\dagger + \delta V^\dagger - \delta\en^v)\phi_v \]

Returns \( \chi_\beta \) for given kappa_beta, where

\[ X_{j,m} = (-1)^{j_\beta-m}tjs(j,k,j;-m,0,m)\chi_j \]

Parameters

Fv	Valence state \(\phi_v\).
omega	Perturbation frequency \(\omega\).
XorY	Selects X or Y solution; see dPsiType.
kappa_beta	Kappa quantum number of the target channel.
Sigma	Optional correlation potential; see MBPT::CorrelationPotential.
st	Bra or ket convention; see StateType.
incl_dV	Include the induced potential \(\delta V\) if true.

◆ solve_dPsis()

std::vector< DiracSpinor > ExternalField::TDHF::solve_dPsis	(	const DiracSpinor &	Fv,
		const double	omega,
		dPsiType	XorY,
		const MBPT::CorrelationPotential *const	Sigma = `nullptr`,
		StateType	st = `StateType::ket`,
		bool	incl_dV = `true`
	)		const

Forms the mixed states phi^v_+- for all kappa channels; see solve_dPsi.

The documentation for this class was generated from the following files:

ExternalField/TDHF.hpp
ExternalField/TDHF.cpp

ok

Public Member Functions

Protected Member Functions

Protected Attributes

Constructor & Destructor Documentation

◆ TDHF()

Member Function Documentation

◆ solve_core()

◆ solve_core_damped()

◆ solve_core_anderson()

◆ set_anderson()

◆ method()

◆ clear()

◆ dV()

◆ dV_rhs()

◆ get_dPsis()

◆ get_dPsi_x()

◆ solve_dPsi()

◆ solve_dPsis()