GitHub

Analyzing a solution

groundstate returns a KSESolution, a self-contained snapshot of the converged (or stopped) calculation: energies, orbital coefficients, occupations, and enough context (model, discretization, basis) to evaluate anything else on demand.

AtomicKohnSham.KSESolution — Type
KSESolution(solver::KSESolver, name::String) -> KSESolution

Construct a solution object from a Kohn–Sham solver after a call to solve!.

This structure stores the final state of a self-consistent field (SCF) computation, including convergence information, final energies, density and orbital coefficients, orbital energies, occupation numbers, Hartree potential coefficients, iteration logs, and the minimal context needed to post-process the solution.

Arguments

  • solver::KSESolver: Solver containing the final state of the SCF computation.
  • name::String: Name or label associated with the solution.

Fields

  • success::String: Final solver status. It is "SUCCESS" if the solver stopped before reaching maxiter, and "MAXITERS" if the maximum number of iterations was reached.
  • niter::Int: Number of SCF iterations performed.
  • stopping_criteria::T: Final value of the stopping criterion.
  • energies::Energies{T}: Final energy contributions stored in an Energies object.
  • D: Coefficients of the one-body reduced density matrix.
  • U: Coefficients of the Kohn–Sham orbitals in the FEM basis.
  • ϵ: Kohn–Sham orbital energies.
  • n: Orbital occupation numbers.
  • occupied: Summary of occupied orbitals, sorted by orbital energy. Each entry contains the shell label, orbital energy, and occupation number.
  • W: Coefficients of the Hartree potential.
  • logbook::LogBook: Logbook containing recorded quantities during the SCF iterations.
  • name::String: Name of the solution, used for identification or labeling.
  • context::KSEContext: Minimal context of the solved problem, containing the model, algorithm, discretization sizes, FEM basis, and integration method.
source

Printing a solution (println(sol), or just sol at the REPL) shows its name, convergence status, and occupied orbitals with their energies and occupation numbers — see the example at the end of Solving for the ground state.

Setup

This continues the running example (each manual page is self-contained, so we solve it again here):

using AtomicKohnSham
using Libxc

model = KSEModel(Z = 11, N = 11, ex = Functional(:lda_x; n_spin = 1), ec = NoFunctional(1))
mesh = expmesh(0, 500, 60; s = 1.2)
basis = P1IntLegendreBasis(mesh; ordermax = 10)
discretization = KSEDiscretization(basis, model; lh = 2, nh = 10,
    fem_integration_method = GaussLegendre(basis, 2000))
alg = ODA(tinit = 0.6, aufbau = OptimizedAufbau(max_degen = 2, tol = 1e-1), scftol = 1e-9)

sol = groundstate(model, discretization, alg; maxiter = 100)

Energies

sol.energies is an Energies holding the kinetic, nuclear attraction (Ecou), Hartree, exchange–correlation, and total energy:

AtomicKohnSham.Energies — Type
Energies{T<:Real}

Container for the energy components computed during the SCF iterations.

Fields (all of type T):

  • Etot : total energy
  • Ekin : kinetic energy
  • Ecou : Coulomb energy
  • Ehar : Hartree energy
  • Eexc : exchange–correlation energy
  • Ekincor : kinetic correlation energy (when applicable)
source
(Etot = sol.energies.Etot, Ekin = sol.energies.Ekin, Ecou = sol.energies.Ecou,
 Ehar = sol.energies.Ehar, Eexc = sol.energies.Eexc)
(Etot = -160.6282275876745, Ekin = 160.62822758864036, Ecou = -388.0065460113607, Ehar = 79.43155311453171, Eexc = -12.681462279485856)

Evaluating orbitals and the density

Radial quantities are evaluated pointwise on a vector of radii X (which should not include r = 0, a removable singularity for the quantities below):

AtomicKohnSham.eval_orbital — Function
eval_orbital(sol, n, l, X; h10=false)
eval_orbital(sol, n, l, σ, X; h10=false)
eval_orbital(sol, idx, X; h10=false)

Evaluate a Kohn–Sham orbital on the FEM basis at radial points X.

This function evaluates the radial part of a Kohn–Sham orbital associated with the quantum numbers (n, l) (and spin σ if applicable) using the FEM basis stored in the solution context.

By default, the returned orbital corresponds to the physical radial wavefunction u(r) / r. If h10=true, the function instead returns the FEM representation u(r) in the space H¹₀, without division by r.

Arguments

  • sol::KSESolution: Converged Kohn–Sham solution.
  • n::Int: Principal quantum number.
  • l::Int: Orbital angular momentum quantum number.
  • σ::Int: Spin index (required for spin-polarized calculations).
  • idx::String: Shell label (e.g. "1s", "2p↑"), parsed internally.
  • X::AbstractVector: Radial evaluation points.

Keyword Arguments

  • h10::Bool: If true, return the FEM H¹₀ representation; otherwise return the physical radial orbital.

Returns

A vector containing the orbital values evaluated at points X.

source
AtomicKohnSham.eval_density — Function
density!(discretization, U, n, D)

Assemble the one-particle density matrix D from Kohn–Sham orbitals U and occupation numbers n.

The density matrix is constructed as D = ∑ₗₖσ nₗₖσ |Uₗₖσ⟩⟨Uₗₖσ|, where the sum runs over angular momentum, radial index, and spin. The result is symmetric and stored in-place in D.

source
eval_density(discretization, D, X)

Evaluate the radial electronic density at all positions in X from the density matrix D.

Returns a vector containing ρ(x) for each x ∈ X.

source
eval_density(sol, X)

Evaluate the (spin-summed) electron density on radial points X.

For spin-unpolarized calculations (nspin == 1), this returns ρ(X) computed from the density matrix stored in sol. For spin-polarized calculations (nspin == 2), this returns the total density ρ↑(X) + ρ↓(X).

Arguments

  • sol::KSESolution: Converged Kohn–Sham solution.
  • X::AbstractVector: Radial evaluation points.

Returns

A vector containing the density values evaluated at X.

source
eval_density(sol, X, σ)

Evaluate the spin density ρσ on radial points X.

This method is available for spin-polarized calculations and returns the density associated with spin channel σ.

Arguments

  • sol::KSESolution: Converged Kohn–Sham solution.
  • X::AbstractVector: Radial evaluation points.
  • σ::Int: Spin index (1-based).

Returns

A vector containing ρσ(X) evaluated at X.

source
X = [0.1, 0.5, 1.0, 2.0, 5.0]
eval_orbital(sol, "3s", X)
5-element Vector{Float64}:
 -0.8365159585491283
  0.36244422750242633
  0.020156528515160553
 -0.18125749662288823
 -0.08113931104303011
eval_density(sol, X)
5-element Vector{Float64}:
 95.63599694310363
  2.922939372897653
  0.41373420256253646
  0.010499565103770224
  0.0005242201546608217

Evaluating potentials

The individual contributions to the effective potential, and their sum, can also be evaluated pointwise — useful for plotting the potential well an orbital sits in (see Reports, logs & plots):

AtomicKohnSham.eval_hartree — Function
eval_hartree(sol::KSESolution, X::AbstractVector{<:Real})

Evaluate the radial Hartree potential at the radial points X.

The Hartree potential solves -1/r (r Vᴴ)'' = 4πρ with the finite-domain boundary lifting θ(r) = r/Rmax, giving (see the model derivation)

Vᴴ(r) = W(r) / r + N / Rmax,

where W is the FEM representation of the regular part w of the Hartree potential (stored in sol.W, solving 𝔸 W = 𝔹[𝔻]), N is the number of electrons, and Rmax is the last point of the radial mesh. As with eval_orbital, X should not contain r = 0 (the FEM representation W vanishes there, so W(r)/r needs to be evaluated as a limit).

Arguments

  • sol::KSESolution: Converged or stopped Kohn–Sham solution.
  • X::AbstractVector{<:Real}: Radial evaluation points.

Returns

A vector containing the Hartree potential values evaluated at X.

source
AtomicKohnSham.eval_nuclear — Function
eval_nuclear(sol::KSESolution, X::AbstractVector{<:Real})

Evaluate the electron-nucleus attraction potential at the radial points X.

The nuclear potential is

Vₙᵤ꜀(r) = -Z / r,

where Z is the nuclear charge stored in sol.context.model. At r = 0, the value is -Inf.

Arguments

  • sol::KSESolution: Kohn–Sham solution containing the model parameters.
  • X::AbstractVector{<:Real}: Radial evaluation points.

Returns

A vector containing the nuclear potential values evaluated at X.

source
AtomicKohnSham.eval_kinetic_potential — Function
eval_kinetic_potential(l::Int, X::AbstractVector{<:Real})
eval_kinetic_potential(sol::KSESolution, l::Int, X::AbstractVector{<:Real})

Evaluate the centrifugal part of the radial kinetic operator at the radial points X.

For an angular momentum quantum number l, the centrifugal potential is

Vₖᵢₙ,l(r) = l(l + 1) / (2r²).

At r = 0, the value is 0 if l == 0, and +Inf if l > 0. The sol argument is accepted (and ignored) only so this function shares its call signature with the other potential evaluators (see eval_effective_potential).

Arguments

  • l::Int: Orbital angular momentum quantum number.
  • X::AbstractVector{<:Real}: Radial evaluation points.
  • sol::KSESolution: Optional solution argument, included for API consistency.

Returns

A vector containing the centrifugal kinetic potential values evaluated at X.

source
AtomicKohnSham.eval_vxc — Function
eval_vxc(sol::KSESolution, X::AbstractVector{<:Real})
eval_vxc(sol::KSESolution, X::AbstractVector{<:Real}, σ::Int)

Evaluate the local (LDA/LSDA) exchange–correlation potential vₓ꜀ = dεₓ꜀/dρ at the radial points X.

The density is first evaluated at X with eval_density, then the exchange–correlation functionals stored in sol.context.model are evaluated pointwise on that density, using the same evaluate_vrho! routine as during SCF assembly (assemble_exc!). Returns a vector of zeros if the model has no exchange–correlation (has_exchcorr(model) == false).

The first method is for spin-unpolarized calculations (nspin == 1); the second requires a spin index σ and is for spin-polarized calculations, returning vₓ꜀,σ(ρ↑(r), ρ↓(r)).

Arguments

  • sol::KSESolution: Converged Kohn–Sham solution.
  • X::AbstractVector{<:Real}: Radial evaluation points.
  • σ::Int: Spin index (required for spin-polarized calculations).

Returns

A vector containing the exchange–correlation potential evaluated at X.

source
AtomicKohnSham.eval_effective_potential — Function
eval_effective_potential(sol::KSESolution, l::Int, X::AbstractVector{<:Real}; σ::Int=1)

Evaluate the total effective radial potential seen by an orbital of angular momentum l (and spin σ if applicable), i.e. the potential part of the radial Kohn–Sham mean-field Hamiltonian:

Vₑff,l,σ(r) = Vₙᵤ꜀(r) + Vₖᵢₙ,l(r) + Vᴴ(r) + Vₓ꜀,σ(r),

summing eval_nuclear, eval_kinetic_potential, eval_hartree (only if model.hartree != 0), and eval_vxc (only if has_exchcorr(model)).

Arguments

  • sol::KSESolution: Converged Kohn–Sham solution.
  • l::Int: Orbital angular momentum quantum number.
  • X::AbstractVector{<:Real}: Radial evaluation points.

Keyword Arguments

  • σ::Int: Spin index, used only for spin-polarized calculations.

Returns

A vector containing the effective potential values evaluated at X.

source
eval_effective_potential(sol, 0, X)  # l = 0 channel
5-element Vector{Float64}:
 -85.7634862695105
  -7.93613432532158
  -1.9426059046507773
  -0.45333278789340203
  -0.09195067666122991

Continue to Reports, logs & plots to export sol to text and plots.