ok

TDHF (RPA) core polarisation above ionisation threshold(s): the open channels are solved with the outgoing-wave boundary condition, so the corrections and the corrected matrix elements are complex.

As TDHF, for a frequency above the ionisation threshold of one or more core orbitals,

\[ \en^a_+ = \en_a + \omega > 0. \]

Every X (+) channel of such an orbital is open, and describes the ejected electron. For it, the solution of the TDHF equation regular at the origin is not unique (the continuum orbital at en_+ may be added freely); the physical solution is the purely outgoing wave at large r, which is complex.

All other channels (closed orbitals, and every Y (-) partner) are bound, and are solved as in TDHF. Through dV every correction is complex: the real and imaginary parts are stored as separate sets (m_X, m_Y and m_Xi, m_Yi) and iterated together. Below every threshold the imaginary sets are exactly zero and the class reproduces TDHF.

In dV the (-) corrections enter as complex conjugates, so the Y sets store \( (\varphi_-)^* \) and the (real-linear) dV builder of TDHF applies to each set unchanged; the external field drives the real set only, and the two sets couple solely through the outgoing boundary condition.

Open channels: The ejected electron moves in the field of the residual ion: the self-interaction of one electron in the ionised orbital, \( V^a_0 = y^0_{aa} + X_a \) (direct and one-electron exchange), is moved from the HF Hamiltonian into the source. This is an exact rearrangement: the channel Hamiltonian has the correct -1/r tail, and the source is short ranged (the y^0_aa term cancels, pointwise, the 1/r tail of the a term of dV phi_a).

With the real regular and irregular solutions F_reg, F_irr of the local channel Hamiltonian at en_+ (F_reg + i F_irr is outgoing), the standing-wave solve returns \( \varphi_P \to K F_{\rm irr} \), and the outgoing solution is \( \varphi_+ = \varphi_P - iK F_{\rm reg} \).

With the non-local exchange iterated in the source, F_reg is replaced by the exchange-dressed regular solution \( F^{\rm HF} \to F_{\rm reg} + K_{\rm ex} F_{\rm irr} \) (built once per channel per omega), and

\[ K_+ = \frac{K}{1 + iK_{\rm ex}}, \]

and

\[ \varphi_+ = \varphi_P - iK_+F^{\rm HF} \to -iK_+\left(F_{\rm reg} + iF_{\rm irr}\right). \]

For a complex source, phi_P is two real standing-wave solves (solveContinuumMixedState).

Iteration: The TDHF map is linear in the corrections, and damped iteration diverges near resonances (autoionising resonances, omega ~ en_b - en_a). The self-consistency is therefore driven by Anderson mixing (anderson_coefficients) of the whole state, all corrections together with the channel amplitudes K_+ (linear in the state), which converges wherever the linear problem is non-singular.

Usage

As TDHF: solve_core (omega), then

A_phys: the ionisation amplitude of each open channel, \( A = K_+^* / \pi \), in channel_list order. |A| replaces \( |\redmatel{\en\kappa}{t}{a}| \) in cross-sections (energy-normalised continuum states); the conjugate is the conventional normalisation of the final state to incoming waves. The phase is relative to the regular solution of the local channel potential: the relative phase of two operators in the same channel is physical, absolute phases and phases between channels are not.
dV_complex (Fa, Fb): the complex \( \redmatel{a}{\delta V}{b} \); equals TDHF::dV below every threshold.

#include <TDHFcntm.hpp>

Inheritance diagram for ExternalField::TDHFcntm:

Classes
struct	Channel
	Open channel: core-orbital index, ejected-electron kappa, and its energy en = en_a + omega. More...

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

void	solve_core (double omega, int max_its=100, bool print=true) override
	Solves the (complex) TDHF equations self-consistently at omega.

void	clear () override
	Clears the corrections and the continuum channel data.

Method	method () const override final
	Returns RPA method.

std::vector< Channel >	channel_list () const
	Open channels at the omega of the last solve_core(), in A_phys order; empty before solve_core(), or if no channel is open.

std::vector< std::complex< double > >	A_phys () const
	Ionisation amplitude A = K_+^*‍/pi of every open channel (see class description), in channel_list order. Requires solve_core().

std::complex< double >	A_phys (const DiracSpinor &Fa, int kappa) const
	A (see A_phys) of the channel of core orbital Fa with ejected-electron kappa; zero if that channel is closed.

std::complex< double >	dV_complex (const DiracSpinor &Fa, const DiracSpinor &Fb, bool conj=false) const
	Reduced matrix element of the (complex) induced potential, <a\|\|dV\|\|b>, or the conjugate <a\|\|dV^dag\|\|b> if conj.

double	dV (const DiracSpinor &Fa, const DiracSpinor &Fb, bool conj=false) const override
	Real part of dV_complex, for bound Fa and Fb only.

double	dV_re (const DiracSpinor &Fa, const DiracSpinor &Fb, bool conj=false) const
	Real part of dV_complex. Convenience only: the full complex element is calculated, so prefer dV_complex if both parts are needed.

double	dV_im (const DiracSpinor &Fa, const DiracSpinor &Fb, bool conj=false) const
	Imaginary part of dV_complex. Convenience only: the full complex element is calculated, so prefer dV_complex if both parts are needed.

TDHFcntm &	operator= (const TDHFcntm &)=delete

	TDHFcntm (const TDHFcntm &)=default

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

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.

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

Additional Inherited Members
Protected Member Functions inherited from ExternalField::TDHF
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 inherited from ExternalField::TDHF
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}

Class Documentation

◆ ExternalField::TDHFcntm::Channel

struct ExternalField::TDHFcntm::Channel

Class Members
size_t	i_core
int	kappa
double	en

Constructor & Destructor Documentation

◆ TDHFcntm()

ExternalField::TDHFcntm::TDHFcntm	(	const DiracOperator::TensorOperator *const	h,
		const HF::HartreeFock *const	hf
	)

Constructs for operator h; see TDHF::TDHF.

Member Function Documentation

◆ solve_core()

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

overridevirtual

Solves the (complex) TDHF equations self-consistently at omega.

Re-solving at the same omega warm-starts from the previous solution; the continuum data of the open channels is rebuilt only when omega changes.

Parameters

omega	Frequency (atomic units); its magnitude is used.
max_its	Maximum number of iterations; 1 gives the first-order (outgoing-wave) correction.
print	If true, write convergence progress to screen.

Reimplemented from ExternalField::TDHF.

◆ clear()

void ExternalField::TDHFcntm::clear ( )

overridevirtual

Clears the corrections and the continuum channel data.

Reimplemented from ExternalField::TDHF.

◆ method()

Method ExternalField::TDHFcntm::method ( ) const

inlinefinaloverridevirtual

Returns RPA method.

Reimplemented from ExternalField::TDHF.

◆ channel_list()

std::vector< TDHFcntm::Channel > ExternalField::TDHFcntm::channel_list ( ) const

Open channels at the omega of the last solve_core(), in A_phys order; empty before solve_core(), or if no channel is open.

◆ A_phys() [1/2]

std::vector< std::complex< double > > ExternalField::TDHFcntm::A_phys ( ) const

Ionisation amplitude A = K_+^*‍/pi of every open channel (see class description), in channel_list order. Requires solve_core().

◆ A_phys() [2/2]

std::complex< double > ExternalField::TDHFcntm::A_phys	(	const DiracSpinor &	Fa,
		int	kappa
	)		const

A (see A_phys) of the channel of core orbital Fa with ejected-electron kappa; zero if that channel is closed.

◆ dV_complex()

std::complex< double > ExternalField::TDHFcntm::dV_complex	(	const DiracSpinor &	Fa,
		const DiracSpinor &	Fb,
		bool	conj = `false`
	)		const

Reduced matrix element of the (complex) induced potential, <a||dV||b>, or the conjugate <a||dV^dag||b> if conj.

The complex version of TDHF::dV. Real, and equal to TDHF, below every threshold. For a continuum bra (en_a > 0), Fa must be the energy-normalised continuum state of hole Fb in the V^{N-1} potential (direct y^0_bb and one-electron exchange of Fb removed); the source then includes the hole term V^b_0 phi of the rearranged channel equation, so that <Fa||t||b> + dV_complex(Fa, Fb) is the outgoing-wave amplitude, |A| of A_phys in the phase reference of Fa.

◆ dV()

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

overridevirtual

Real part of dV_complex, for bound Fa and Fb only.

Exact below every threshold, where the class acts as TDHF.

Warning: Asserts (and aborts) if either state is a continuum state (en > 0): dV_complex must be used there.

Reimplemented from ExternalField::TDHF.

◆ dV_re()

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

Real part of dV_complex. Convenience only: the full complex element is calculated, so prefer dV_complex if both parts are needed.

◆ dV_im()

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

Imaginary part of dV_complex. Convenience only: the full complex element is calculated, so prefer dV_complex if both parts are needed.

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

ExternalField/TDHFcntm.hpp
ExternalField/TDHFcntm.cpp

ok

Classes

Public Member Functions

Additional Inherited Members

Class Documentation

◆ ExternalField::TDHFcntm::Channel

Constructor & Destructor Documentation

◆ TDHFcntm()

Member Function Documentation

◆ solve_core()

◆ clear()

◆ method()

◆ channel_list()

◆ A_phys() [1/2]

◆ A_phys() [2/2]

◆ dV_complex()

◆ dV()

◆ dV_re()

◆ dV_im()