High-precision calculations for one- and two-valence atomic systems
ExternalField::TDHFcntm

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 <TDHFcomplex.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_plus, const HF::HartreeFock *const hf, const DiracOperator::TensorOperator *const h_minus=nullptr)
 Constructs for operator h_plus (with optional h_minus); 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.
 
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) const
 Reduced matrix element of the (complex) induced potential, \( \redmatel{a}{\delta V}{b} \), or the conjugate \( \redmatel{a}{\delta V^\dagger}{b} \) if en_b > en_a; as TDHF::dV, complex.
 
double dV (const DiracSpinor &Fa, const DiracSpinor &Fb) const override
 Real part of dV_complex, for bound Fa and Fb (exact below every threshold, where the class acts as TDHF). For a continuum state use dV_complex.
 
TDHFcntm & operator= (const TDHFcntm &)=delete
 
 TDHFcntm (const TDHFcntm &)=default
 
double dV (const DiracSpinor &Fa, const DiracSpinor &Fb, bool conj) const
 Returns reduced matrix element \(\redmatel{a}{\delta V}{b}\), or the conjugate \(\redmatel{a}{\delta V^\dagger}{b}\) if conj=true.
 
virtual double dV (const DiracSpinor &Fa, const DiracSpinor &Fb) const override
 Returns reduced matrix element <n||dV_pm||m> (see namespace doc for dV_pm)
 
- Public Member Functions inherited from ExternalField::TDHF
 TDHF (const DiracOperator::TensorOperator *const h_plus, const HF::HartreeFock *const hf, const DiracOperator::TensorOperator *const h_minus=nullptr)
 Constructs TDHF for operator h.
 
virtual Method method () const override
 Returns RPA method.
 
double dV (const DiracSpinor &Fa, const DiracSpinor &Fb, bool conj) const
 Returns reduced matrix element \(\redmatel{a}{\delta V}{b}\), or the conjugate \(\redmatel{a}{\delta V^\dagger}{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 \(\varphi^v_\pm\) 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 \(\varphi^v_\pm\) 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
 
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
 
const DiracOperator::TensorOperator *const m_h_minus
 
bool m_eps_sqrt {false}
 
- 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.4}
 
double m_eps {1.0e-10}
 

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_plus,
const HF::HartreeFock *const  hf,
const DiracOperator::TensorOperator *const  h_minus = nullptr 
)

Constructs for operator h_plus (with optional h_minus); 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.

Parameters
omegaFrequency (atomic units); its magnitude is used.
max_itsMaximum number of iterations; 1 gives the first-order (outgoing-wave) correction.
printIf true, write convergence progress to screen.

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.

Reimplemented from ExternalField::TDHF.

◆ clear()

void ExternalField::TDHFcntm::clear ( )
overridevirtual

Clears the corrections and the continuum channel data.

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 
) const

Reduced matrix element of the (complex) induced potential, \( \redmatel{a}{\delta V}{b} \), or the conjugate \( \redmatel{a}{\delta V^\dagger}{b} \) if en_b > en_a; as TDHF::dV, complex.

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() [1/3]

double ExternalField::TDHFcntm::dV ( const DiracSpinor &  Fa,
const DiracSpinor &  Fb 
) const
overridevirtual

Real part of dV_complex, for bound Fa and Fb (exact below every threshold, where the class acts as TDHF). For a continuum state use dV_complex.

Reimplemented from ExternalField::TDHF.

◆ dV() [2/3]

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

Returns reduced matrix element \(\redmatel{a}{\delta V}{b}\), or the conjugate \(\redmatel{a}{\delta V^\dagger}{b}\) if conj=true.

◆ dV() [3/3]

double ExternalField::TDHF::dV ( const DiracSpinor &  Fn,
const DiracSpinor &  Fm 
) const
overridevirtual

Returns reduced matrix element <n||dV_pm||m> (see namespace doc for dV_pm)

Reimplemented from ExternalField::TDHF.


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