Atomic coulomb repulsion

NewOperator(“AtomicU”,indices,RadialWavefunctions1,RadialWavefunctions2,options) calculates the full Coulomb interaction operator between all given atomic orbitals.

• indices : indices of atomic orbitals belonging to the radial wave-functions
• RadialWavefunctions1 : table of InterpolatingFunction objects
• RadialWavefunctions2 : optional - table of InterpolatingFunction objects - if given, the function calculates relativistic Slater integrals where RadialWavefunctions1 are the large parts and RadialWavefunctions2 the small ones
• options : possible Options are {“conserve”, bool} and {“1P-Operator”, bool}

  • “conserve” : if true only occupation conserving Slater integrals are calculated
  • “1P-Operator” : only if “conserve” is true and additional arguments are given, transforms interaction to 1-particle operator (see example)

The radial wave-functions can be read from FPLO with ReadFPLOBasisFunctions(), or calculated with CalculateAtomicBasisFunctions(). For the functions of CalculateAtomicBasisFunctions() the Slater integrals are calculated exactly on the B-spline basis of the atom, without a radial grid (see the second example). All functions must then come from the same call.

The example below shows how to create the full Coulomb interaction between multiple atomic orbitals with relativistic wave-functions. It is recommended to use CreateAtomicIndicesList(orbitallabels,{ {“Kappas”,true} }) for index creation, since it automatically calculates the quantum numbers kappa and stores them in the field ind.kappas which needs to be set when calling NewOperator(“AtomicU”,…) with relativistic wave-functions.

Example.Quanty

-- the Coulomb interaction of the electrons of Ne on its own Dirac-Fock orbitals
G, F, x, Atom = CalculateAtomicBasisFunctions(10)
ind, NF = CreateAtomicIndicesList(Atom.Labels, {{"Kappas", true}})
 
U = NewOperator("AU", NF, ind, G, F)
 
psi = NewWavefunction(NF, 0, { {string.rep("1", NF), 1} })
print(string.format("<psi|U|psi> = %.8f Hartree, Edirect + Eexchange = %.8f Hartree", psi*U*psi, Atom.Edirect + Atom.Eexchange))

Example.out

<psi|U|psi> = 54.08660432 Hartree, Edirect + Eexchange = 54.08660432 Hartree