$U_{n, \mathbf{k}}$ Orbitals

An AbstractUnkOrbital represents the$u_{n, \mathbf{k}}$ orbitals.

Such a function is also a basis vector and a linear combination. This is where the single subtyping system of Julia gives me headaches.

UnkBasisOrbital(grid, elements, kpoint, index_band)

This is the concrete type that numerically stores the orbital as a function on a grid (either the homecell or the reciprocal lattice). It is also associated with a k-point and an index of a band.

More likely than not, one get these from reading the output of quantum espresso.

julia> lattice = make_grid(ReciprocalLattice3D, CARTESIAN_BASIS, size_to_domain((4, 4, 4)));

julia> ϕ = map(g->(g==lattice[1, 0, 0]) |> Float64, lattice);

julia> brillouin_zone = make_grid(BrillouinZone3D, CARTESIAN_BASIS, size_to_domain((3, 3, 3)));

julia> ϕ = UnkBasisOrbital(lattice, elements(ϕ), brillouin_zone[1, 0, -1], 2)
ket
grid:
    type: HomeCell3D
    domain: ((-2, 1), (-2, 1), (-2, 1))
    basis:
        ket: 1.000, 0.000, 0.000
        ket: 0.000, 1.000, 0.000
        ket: 0.000, 0.000, 1.000
kpoint:
    GridVector{BrillouinZone3D}:
        coefficients: [1, 0, -1]
    
band:
    2

kpoint(orbital)

Get the kpoint associated with the orbital.

julia> kpoint(ϕ)
GridVector{BrillouinZone3D}:
    coefficients: [1, 0, -1]

kpoint!(orbital, new_kpoint)

Set the kpoint associated with orbital to new_kpoint.

julia> kpoint!(ϕ, brillouin_zone[0, 0, -1])
ket
grid:
    type: HomeCell3D
    domain: ((-2, 1), (-2, 1), (-2, 1))
    basis:
        ket: 1.000, 0.000, 0.000
        ket: 0.000, 1.000, 0.000
        ket: 0.000, 0.000, 1.000
kpoint:
    GridVector{BrillouinZone3D}:
        coefficients: [0, 0, -1]
    
band:
    2

index_band(orbital)

Get the index of the band of the orbital.

julia> index_band(ϕ)
2

dagger(orbital)

Complex conjutate an orbital.

julia> dagger(ϕ) |> ket
false

dagger!(orbital)

Complex conjutate an orbital in place.

Given a few basis orbitals, we can construct a linear combination of them for a basic set of symbolic manipulation. A linear combination is named UnkOrbital due to my stupidity.

This also should be a subtype of OnGrid and Basis, but we are crippled by the single subtyping system.

UnkOrbital(orbital, orthonormal = true)

Construct a linear combination of orbitals from orbital as a basis vector. Set orthonormal to true for an orthogonal basis set.

julia> ϕ₁ = UnkOrbital(map(g->Float64(g==lattice[1, 0, 0]), lattice), true)
ket
coefficients:
    Number[1.0]
n_basis:
    1

julia> ϕ₂ = UnkOrbital(map(g->Float64(g==lattice[0, 1, 0]), lattice), true);

julia> ϕ₃ = UnkOrbital(map(g->Float64(g==lattice[0, 0, 1]), lattice), true); 

dagger(orbital)

Take the complex conjugate recursively through the basis vectors.

julia> dagger(ϕ₁) |> ket
false

dagger!(orbital)

Take the complex conjugate recursively through the basis vectors in place.

index_band(orbital)

Get the index of the band of the orbital.

orthonormal(orbital)

Whether the basis set is orthonormal.

add(o_1, o_2, [mutual_orthonormal = false])

Can also write o_1 + o_2, where mutual_orthonormal = false is assumed.

Add two linear combination of orbitals. The resulting linear combination has a orthogonal basis set if both linear combinations have orthogonal basis sets and they are mutually orthogonal as specified in mutual_orthonormal.

julia> ϕ_sum = ϕ₁ + ϕ₂
ket
coefficients:
    Number[1.0, 1.0]
n_basis:
    2
julia> orthonormal(ϕ_sum)
false
julia> ϕ_sum = add(ϕ₁, ϕ₂, true)
ket
coefficients:
    Number[1.0, 1.0]
n_basis:
    2

julia> orthonormal(ϕ_sum)
true

negate(orbital)

Negate the coefficients of the linear combination. Can also write -orbital.

mul(s, orbital)

Scale coefficients of the linear combination. Can also write s * orbital or orbital * s.

braket(o_1, o_2)

Recursively compute the inner product. o_1 must be a bra and o_2 must be a ket (this may be relaxed later). Can also write o_1 * o_2

Example:

julia> ϕ_sum' * ϕ₁
1.0
julia> ϕ_sum' * ϕ₃
0.0

zeros(UnkOrbital, dims...)

Create an array of empty linear combination. This enables semi-symbolic linear algebra.

julia> [ϕ₁, ϕ₂, ϕ₃]' * [0 0 1;
                        0 1 0;
                        1 0 0] * [ϕ₁, ϕ₂, ϕ₃]
1.0