API

SimulationInfo

Contains identification information about simulation, including the location data is written to, the simulation ID, and MPI process ID, and whether this simulation started a new simulation or resumed a previous simulation.

Fields

• filepath::String: File path to where data folder lives.
• datafolder_prefix: Prefix for the data folder name.
• datafolder_name::String: The data folder name, given by $(datafolder_prefix)_$(sID).
• datafolder::String: The data folder, including filepath, given by joinpath(filepath, datafolder_name).
• pID::Int: MPI process ID, defaults to 0 if MPI not being used.
• sID::Int: Simulation ID.
• write_bins_concurrent::Bool: Whether binned data will be written to HDF5 during the simulation or held in memory until the end of the simulation.
• bin_files::Vector{Vector{UInt8}}: Represents the HDF5 files containing the binned data.
• resuming::Bool: Whether current simulation is resuming a previous simulation (true) or starting a new one (false).
• smoqy_version::VersionNumber: Version of SmoQyDQMC.jl used in simulation.

Notes

If write_bins_concurrent = true, then the elements of bin_files correspond to the HDF5 bin filenames, assuming the vector elements are converted to strings. If write_bins_concurrent = false, then the elements of the bin_files correspond to a byte vector representation of a HDF5 file containing the binned data. For small simulations that run very fast setting write_bins_concurrent = false can make sense, as it significantly reduces the frequency of file IO during the simulation. This can cause issues on some clusters with respect to overtaxing the cluster network if data is being written to file too frequently during the simulation. However, for most larger simulations it is advisable to set write_bins_concurrent = true as this significantly reduces the memory footprint of the simulation, particularly when making time-displaced correlation measurements. Also, setting write_bins_concurrent = false dramatically increases the size of the checkpoint files if checkpointing is occurring during the simulation, as the checkpoint files now need to contain all the binned data collected during the simulation.

SimulationInfo(;
    # KEYWORD ARGUMENTS
    datafolder_prefix::String,
    filepath::String = ".",
    write_bins_concurrent::Bool = true,
    sID::Int=0,
    pID::Int=0
)

Initialize and return in instance of the type SimulationInfo.

initialize_datafolder(comm::MPI.Comm, sim_info::SimulationInfo)

initialize_datafolder(sim_info::SimulationInfo)

Initialize sim_info.datafolder directory if it does not already exist. If comm::MPI.Comm is passed as the first argument, this this function will synchronize all the MPI processes, ensuring that none proceed beyond this function call until the data folder that results will be written to is successfully initialized.

model_summary(;
    simulation_info::SimulationInfo,
    β::T, Δτ::T, model_geometry::ModelGeometry,
    tight_binding_model::Union{TightBindingModel,Nothing} = nothing,
    tight_binding_model_up::Union{TightBindingModel,Nothing} = nothing,
    tight_binding_model_dn::Union{TightBindingModel,Nothing} = nothing,
    interactions::Union{Tuple,Nothing} = nothing
) where {T<:AbstractFloat}

Write model to summary to file. Note that either tight_binding_model or tight_binding_model_up and tight_binding_model_dn need to be specified.

ModelGeometry{D, T<:AbstractFloat, N}

Contains all the information defining the lattice geometry for the model in D spatial dimensions.

Comment

The bond ID associated with a bond::Bond{D} corresponds to the index associated with it into the bonds vector field.

Fields

• unit_cell::UnitCell{D,T,N}: Defines unit cell.
• lattice::Lattice{D}: Defines finite lattice extent.
• bonds::Vector{Bond{D}}: All available bond definitions in simulation, with vector indices giving the bond ID.

ModelGeometry(unit_cell::UnitCell, lattice::Lattice)

Initialize and return a ModelGeometry instance. Defines a "trivial" bond definition for each orbital in the unit cell that connects an orbital to itself.

addbond!(modelgeometry::ModelGeometry{D,T}, bond::Bond{D}) where {D, T}

Add bond definition to model_geometry, returning the bond ID i.e. the index to bond in the vector model_geometry.bonds. This method first checks that bond is not already defined. If it is this method simply returns the corresponding bond ID. If bond is not already defined, then it is appended to the vector model_geometry.bonds.

get_bond_id(model_geometry::ModelGeometry{D,T}, bond::Bond{D}) where {D, T}

Return the bond ID associated with the bond defintion bond, returning bond_id=0 if the it is not a recorded bond.

FermionPathIntegral{H<:Number, T<:Number, U<:Number, R<:AbstractFloat}

A type (mutable struct) to represent a fermion path integral. In particular, this type contains the information required to reconstruct the diagonal on-site energy matrices $V_l$ and hopping matrices $K_l$ for each imaginary time slice $l \in [1, L_\tau],$ where $\tau = \Delta\tau \cdot l$ and $\beta = \Delta\tau \cdot L_\tau.$

Types

• H<:Number: $H_l = (K_l + V_l)$ Hamiltonian matrix element type.
• T<:Number: $K_l$ kinetic energy matrix element type.
• U<:Number: $V_l$ potential energy matrix element type.
• R<:AbstractFloat: Real number type.

Fields

• β::R: Inverse temperature.
• Δτ::R: Discretization in imaginary time.
• Lτ::Int: Length of the imaginary time axis.
• N::Int: Number of orbitals in the lattice.
• neighbor_table::Matrix{Int}: Neighbor table for each pair of orbitals in the lattice connected by a hopping.
• t::Matrix{T}: Hopping amplitudes for imaginary-time slice $l$ are stored in t[:,l].
• V::Matrix{U}: The diagonal on-site energy matrices $V_l$ for imaginary-time slice $l$ are stored in V[:,l].
• K::Matrix{T}: Used to construct hopping matrix to calculate exponentiated hopping matrix if checkerboard approximation is not being used.
• Sb::H: Keeps track of total bosonic action associated with fermionic path integral.
• eigen_ws::HermitianEigenWs{T,Matrix{T},R}: For calculating eigenvalues and eigenvectors of K while avoiding dynamic memory allocations.
• u::Matrix{H}: Temporary matrix to avoid dynamic allocation when performing local updates.
• v::Matrix{H}: Temporary matrix to avoid dynamic allocation when performing local updates.

FermionPathIntegral(;
    # KEYWORD ARGUMENTS
    tight_binding_parameters::TightBindingParameters{T,R},
    β::R, Δτ::R,
    forced_complex_kinetic::Bool = false,
    forced_complex_potential::Bool = false
) where {T<:Number, R<:AbstractFloat}

Initialize an instance of FermionPathIntegral an instance of TightBindingParameters.

If forced_complex_kinetic = true, then the off-diagonal kinetic energy matrices $K_l$ are assumed to be complex, otherwise the matrix element type is inferred.

If forced_complex_potential = true, then the diagonal potential energy matrices $V_l$ are assumed to be complex, otherwise the matrix element type is inferred.

initialize_propagators(
    # ARGUMENTS
    fermion_path_integral::FermionPathIntegral;
    # KEYWORD ARGUMENTS
    symmetric::Bool,
    checkerboard::Bool
)

Initialize a propagator for each imaginary time slice, returning a vector of type Vector{<:AbstractPropagators{T,E}}.

calculate_propagators!(
    # ARGUMENTS
    B::Vector{P},
    fpi::FermionPathIntegral;
    # KEYWORD ARGUMENTS
    calculate_exp_V::Bool,
    calculate_exp_K::Bool
) where {P<:AbstractPropagator}

Calculate the propagator matrices $B_l$, given by B[l], for all imaginary time slice $\tau = \Delta\tau \cdot l.$ If calculate_exp_V = true, then calculate the diagonal exponentiated on-site energy matrices. If calculate_exp_K = true, then calculate the exponentiated hopping matrices.

calculate_propagator!(
    # ARGUMENTS
    B::P,
    fpi::FermionPathIntegral{H,T,U},
    l::Int;
    # KEYWORD ARGUMENTS
    calculate_exp_V::Bool,
    calculate_exp_K::Bool
) where {H<:Number, T<:Number, U<:Number, P<:AbstractPropagator{T,U}}

Calculate the propagator matrix $B_l$ for imaginary time slice $\tau = \Delta\tau \cdot l.$ If calculate_exp_V = true, then calculate the diagonal exponentiated on-site energy matrix. If calculate_exp_K = true, then calculate the exponentiated hopping matrix.

update_stabilization_frequency!(
    Gup::Matrix{H}, logdetGup::R, sgndetGup::H,
    Gdn::Matrix{H}, logdetGdn::R, sgndetGdn::H;
    fermion_greens_calculator_up::FermionGreensCalculator{H,R},
    fermion_greens_calculator_dn::FermionGreensCalculator{H,R},
    Bup::Vector{P}, Bdn::Vector{P}, δG::R, δθ::R, δG_max::R
) where {H<:Number, R<:Real, P<:AbstractPropagator}

If the corrected error in the Green's function matrix is too large, δG > δG_max, then increase the frequency of numerical stabilization by decrementing n_stab such that it is updated to n_stab = max(n_stab - 1, 1), and update the equal-time Green's function matrices and all related variables and types. If the frequency of stabilization is udpated, then δG and δθ are reset to zero. This method returns a tuple of the following variables:

(logdetGup, sgndetGdn, logdetGup, sgndetGdn, δG, δθ),

where updated = true if n_stab was decremented.

update_stabilization_frequency!(
    G::Matrix{H}, logdetG::R, sgndetG::H;
    fermion_greens_calculator::FermionGreensCalculator{H,R},
    B::Vector{P}, δG::R, δθ::R, δG_max::R
) where {H<:Number, R<:Real, P<:AbstractPropagator}

If the corrected error in the Green's function matrix is too large, δG > δG_max, then increase the frequency of numerical stabilization by decrementing n_stab such that it is updated to n_stab = max(n_stab - 1, 1), and update the equal-time Green's function matrices and all related variables and types. If the frequency of stabilization is udpated, then δG and δθ are reset to zero. This method returns a tuple of the following variables:

(updated, logdetG, sgndetG, δG, δθ),

where updated = true if n_stab was decremented.

TightBindingModel{T<:Number, E<:AbstractFloat, D}

Defines a tight binding model in D dimensions. Note that spin = 1 (spin = 2) corresponds to spin-up (spin-down), and spin = 0 corresponds to both spin-up and spin-down.

Types

• T<:Number: The type of the hopping energy, which can be a real or complex number.
• E<:AbstractFloat: The type of the chemical potential and on-site energy, which must be real.
• D: The number of spatial dimensions the model lives in.

Fields

• μ::E: Chemical potential.
• ϵ_mean::Vector{E}: Mean on-site energy for each orbital in the unit cell.
• ϵ_std::Vector{E}: Standard deviation of on-site energy for each orbital in the unit cell.
• t_bond_ids::Vector{Int}: The bond ID for each bond/hopping definition.
• t_bonds::Vector{Bond{D}}: Bond definition for each type of hopping in the tight binding model.
• t_mean::Vector{T}: Mean hopping energy for each type of hopping.
• t_std::Vector{E}: Standard deviation of hopping energy for each type of hopping.
• η::Vector{E}: Relative twist angle for each lattice vector direction $d \in [1, D]$ such that $\eta_d \in [0, 1).$
• expniϕ::Vector{T}: Twist angle phase $\exp(i n_d \phi_d)$ for each hopping.

TightBindingModel(;
    # KEYWORD ARGUMENTS
    model_geometry::ModelGeometry{D,E,N},
    μ::E,
    ϵ_mean::Vector{E},
    ϵ_std::Vector{E} = zeros(eltype(ϵ_mean), length(ϵ_mean)),
    t_bonds::Vector{Bond{D}} = Bond{ndims(model_geometry.unit_cell)}[],
    t_mean::Vector{T} = eltype(ϵ_mean)[],
    t_std::Vector{E} = zeros(eltype(ϵ_mean), length(t_mean)),
    η::Union{Vector{E},Nothing} = nothing
) where {T<:Number, E<:AbstractFloat, D, N}

Initialize and return an instance of TightBindingModel, also adding/recording the bond defintions t_bonds to the ModelGeometry instance model_geometry.

TightBindingParameters{T<:Number, E<:AbstractFloat}

A mutable struct containing all the parameters needed to characterize a finite tight-binding Hamiltonian for a single spin species $\sigma$ on a finite lattice with periodic boundary conditions of the form

\[\hat{H}_{0,\sigma}=-\sum_{\langle i,j\rangle}(t_{ij} \hat{c}_{\sigma,i}^{\dagger}\hat{c}_{\sigma,j}+\textrm{h.c.})+\sum_{i}(\epsilon_{i}-\mu)\hat{n}_{\sigma,i},\]

where $\hat{c}_{\sigma,i}^\dagger$ is the fermion creation operator for an electron with spin $\sigma$ on orbital $i,$ $t_{i,j}$ are the hopping energies, $\epsilon_i$ are the on-site energies for each orbital in the lattice, and $\mu$ is the chemical potential.

Fields

• μ::E: The chemical potential $\mu.$
• const ϵ::Vector{E}: A vector containing the on-site energy $\epsilon_i$ for each orbital in the lattice.
• const t::Vector{T}: The hopping energy $t_{i,j}$ associated with each pair of neighboring orbitals connected by a bond in the lattice.
• const neighbor_table::Matrix{Int}: Neighbor table containing all pairs of orbitals in the lattices connected by a bond, with a non-zero hopping energy between them.
• const bond_ids::Vector{Int}: The bond ID definitions that define the types of hopping in the lattice.
• const bond_slices::Vector{UnitRange{Int}}: Slices of neighbor_table corresponding to given bond ID i.e. the neighbors neighbor_table[:,bond_slices[i]] corresponds the bond_ids[i] bond defintion.
• const norbital::Int: Number of orbitals per unit cell.

TightBindingParameters(;
    tight_binding_model::TightBindingModel{T,E,D},
    model_geometry::ModelGeometry{D,E},
    rng::AbstractRNG
) where {T,E,D}

Initialize and return an instance of TightBindingParameters.

measure_onsite_energy(
    tight_binding_parameters::TightBindingParameters{T,E},
    G::Matrix{H}, orbital_id::Int
) where {H<:Number, T<:Number, E<:AbstractFloat}

Measure and return the on-site energy $\epsilon_\textrm{on-site} = (\epsilon - \mu)\langle \hat{n}_\sigma \rangle$ for the orbital_id in the unit cell.

measure_hopping_energy(
    tight_binding_parameters::TightBindingParameters{T,E},
    fermion_path_integral::FermionPathIntegral{H},
    G::Matrix{H}, hopping_id::Int
) where {H<:Number, T<:Number, E<:AbstractFloat}

Calculate the average hopping energy $\epsilon_{\rm hopping} = -\langle t_{l,\langle i,j \rangle} \hat{c}^\dagger_{\sigma,i} \hat{c}_{\sigma,j} + {\rm h.c.} \rangle$ for the hopping defined by the the hopping_id.

measure_bare_hopping_energy(
    tight_binding_parameters::TightBindingParameters{T,E},
    G::Matrix{H}, hopping_id::Int
) where {H<:Number, T<:Number, E<:AbstractFloat}

Calculate the average bare hopping energy $\epsilon_{\rm hopping} = -\langle t_{\langle i,j \rangle} \hat{c}^\dagger_{\sigma,i} \hat{c}_{\sigma,j} + {\rm h.c.} \rangle$ for the hopping defined by the hopping_id.

measure_hopping_amplitude(
    tight_binding_parameters::TightBindingParameters{T,E},
    fermion_path_integral::FermionPathIntegral{H},
    hopping_id::Int
) where {H<:Number, T<:Number, E<:AbstractFloat}

Calculate the average hopping amplitude for the hopping defined by the hopping_id.

measure_hopping_inversion(
    tight_binding_parameters::TightBindingParameters{T,E},
    fermion_path_integral::FermionPathIntegral{H},
    hopping_id::Int
) where {H<:Number, T<:Number, E<:AbstractFloat}

Measure the fraction of time the sign of the instantaneous modulated hopping amplitude $t_{l,(\mathbf{i},\nu),(\mathbf{j},\gamma)}$ is inverted relative to the bare hopping amplitude $t_{(\mathbf{i},\nu),(\mathbf{j},\gamma)}$, where $l$ is the imaginary time-slice index.

HubbardModel{T<:AbstractFloat}

A type to represent a, in general, multi-orbital Hubbard model.

If the type field ph_sym_form = false, then the particle-hole asymmetric form for the Hubbard interaction

\[\hat{H}_{U}=\sum_{\mathbf{i},\nu}U_{\nu,\mathbf{i}}\hat{n}_{\uparrow,\nu,\mathbf{i}}\hat{n}_{\downarrow,\nu,\mathbf{i}}\]

is used, where $\mathbf{i}$ specifies the unit cell, and $\nu$ denotes the orbital in the unit cell. In the case of a bipartite lattice with only nearest neighbor hopping, this convention results in an on-site energy corresponding to half-filling and particle-hole symmetry given by $\epsilon_{\nu,\mathbf{i}} = -U_{\nu,\mathbf{i}}/2.$

If ph_sym_form = true, then the particle-hole symmetric form for the Hubbard interaction

\[\hat{H}_{U}=\sum_{\mathbf{i},\nu}U_{\nu,\mathbf{i}}(\hat{n}_{\uparrow,\nu,\mathbf{i}}-\tfrac{1}{2})(\hat{n}_{\downarrow,\nu,\mathbf{i}}-\tfrac{1}{2})\]

is used instead. In this case, for a bipartite lattice with only nearest neighbor hopping, the on-site energy corresponding to half-filling and particle-hole symmetry is $\epsilon_{\nu,\mathbf{i}} = 0.$

Fields

• ph_sym_form::Bool: Determines whether the particle-hole symmetric form of the Hubbard interaction is used.
• U_orbital_ids::Vector{Int}: Orbital species/IDs in unit cell with finite Hubbard interaction.
• U_mean::Vector{T}: Average Hubbard interaction strength $U_\nu$ for a given orbital species in the lattice.
• U_std::Vector{T}: Standard deviation of Hubbard interaction strength $U_\nu$ for a given orbital species in the lattice.

HubbardModel(;
    # KEYWORD ARGUMENTS
    ph_sym_form::Bool,
    U_orbital::AbstractVector{Int},
    U_mean::AbstractVector{T},
    U_std::AbstractVector{T} = zero(U_mean)
) where {T<:AbstractFloat}

Initialize and return an instance of the type HubbardModel.

Keyword Arguments

• ph_sym_form::Bool: Determines whether the particle-hole symmetric form of the Hubbard interaction is used.
• U_orbital::Vector{Int}: Orbital species/IDs in unit cell with finite Hubbard interaction.
• U_mean::Vector{T}: Average Hubbard interaction strength $U_\nu$ for a given orbital species in the lattice.
• U_std::Vector{T}: Standard deviation of Hubbard interaction strength $U_\nu$ for a given orbital species in the lattice.

HubbardParameters{T<:AbstractFloat}

Hubbard parameters for finite lattice.

Fields

• U::Vector{T}: On-site Hubbard interaction for each site with finite Hubbard interaction.
• sites::Vector{Int}: Site index associated with each finite Hubbard U interaction.
• orbital_ids::Vector{Int}: Orbital ID/species in unit cell with finite Hubbard interaction.
• ph_sym_form::Bool: Convention used for Hubbard interaction, refer to HubbardModel for more information.

HubbardParameters(;
    hubbard_model::HubbardModel{T},
    model_geometry::ModelGeometry{D,T},
    rng::AbstractRNG
) where {D, T<:AbstractFloat}

Initialize an instance of HubbardParameters.

initialize!(
    fermion_path_integral_up::FermionPathIntegral,
    fermion_path_integral_dn::FermionPathIntegral,
    hubbard_parameters::HubbardParameters
)

initialize!(
    fermion_path_integral::FermionPathIntegral,
    hubbard_parameters::HubbardParameters
)

Initialize the contribution from the Hubbard interaction to a FermionPathIntegral instance.

measure_hubbard_energy(
    hubbard_parameters::HubbardParameters{E},
    Gup::Matrix{T}, Gdn::Matrix{T},
    hubbard_id::Int
) where {T<:Number, E<:AbstractFloat}

Calculate the average Hubbard energy $U \langle \hat{n}_\uparrow \hat{n}_\downarrow \rangle$ if ph_sym_form = false and $U \langle (\hat{n}_\uparrow - \tfrac{1}{2})(\hat{n}_\downarrow - \tfrac{1}{2})\rangle$ if ph_sym_form = true for the orbital corresponding orbital_id in the unit cell.

ExtendedHubbardModel{T<:AbstractFloat}

A type to represent extended Hubbard interactions.

If the type field ph_sym_form = false then the particle-hole asymmetric form of the extended Hubbard interaction

\[\begin{align*} \hat{H}_{V} = \sum_{\mathbf{j},\mathbf{r},\nu,\eta}V_{(\mathbf{j}+\mathbf{r},\nu),(\mathbf{j},\eta)} & \hat{n}_{\mathbf{j}+\mathbf{r},\nu}\hat{n}_{\mathbf{j},\eta} \\ = \sum_{\mathbf{j},\mathbf{r},\nu,\eta}V_{(\mathbf{j}+\mathbf{r},\nu),(\mathbf{j},\eta)} & \bigg[\tfrac{1}{2}(\hat{n}_{\mathbf{j}+\mathbf{r},\nu}+\hat{n}_{\mathbf{j},\eta}-2)^{2}-1 \\ & -\hat{n}_{\mathbf{j}+\mathbf{r},\nu,\uparrow}\hat{n}_{\mathbf{j}+\mathbf{r},\nu\downarrow}-\hat{n}_{\mathbf{j},\eta,\uparrow}\hat{n}_{\mathbf{j},\eta\downarrow}+\tfrac{3}{2}\hat{n}_{\mathbf{j}+\mathbf{r},\nu}+\tfrac{3}{2}\hat{n}_{\mathbf{j},\eta}\bigg] \end{align*}\]

is used, where $\mathbf{j}$ specifies a unit cell in the lattice, $\mathbf{r}$ is a displacement in units, and $\nu$ and $\eta$ specify the orbital in a given unit cell. Here, $\hat{n}_{\mathbf{j},\eta} = (\hat{n}_{\uparrow,\mathbf{j},\eta} + \hat{n}_{\downarrow,\mathbf{j},\eta})$ is the electron number operator for orbital $\eta$ in unit cell $\mathbf{j}$ in the lattice. Therefore, $V_{(\mathbf{j}+\mathbf{r},\nu),(\mathbf{j},\eta)}$ controls the strength of the extended Hubbard interaction between orbital $\eta$ in unit cell $\mathbf{j}$ and orbital $\nu$ in unit cell $\mathbf{j}+\mathbf{r}$.

If the type field ph_sym_form = true, then the particle-hole symmetric for the extended Hubbard interaction

\[\begin{align*} \hat{H}_{V}=\sum_{\mathbf{j},\mathbf{r},\nu,\eta}V_{(\mathbf{j}+\mathbf{r},\nu),(\mathbf{j},\eta)}&(\hat{n}_{\mathbf{j}+\mathbf{r},\nu}-1)(\hat{n}_{\mathbf{j},\eta}-1) \\ = \sum_{\mathbf{j},\mathbf{r},\nu,\eta}V_{(\mathbf{j}+\mathbf{r},\nu),(\mathbf{j},\eta)}&\bigg[\tfrac{1}{2}(\hat{n}_{\mathbf{j}+\mathbf{r},\nu}+\hat{n}_{\mathbf{j},\eta}-2)^{2}+\tfrac{1}{2} \\ & -(\hat{n}_{\mathbf{j}+\mathbf{r},\nu,\uparrow}-\tfrac{1}{2})(\hat{n}_{\mathbf{j}+\mathbf{r},\nu\downarrow}-\tfrac{1}{2})-(\hat{n}_{\mathbf{j},\eta,\uparrow}-\tfrac{1}{2})(\hat{n}_{\mathbf{j},\eta\downarrow}-\tfrac{1}{2})\bigg] \end{align*}\]

is used instead.

Fields

• ph_sym_form::Bool: Whether the particle-hole symmetric form of the extended Hubbard interaction is used.
• V_bond_ids::Vector{Int}: Bond IDs specifying bond definition that separates a pair of orbitals with an extended Hubbard interaction between them.
• V_mean::Vector{T}: Average extended Hubbard interaction strength $V_{(\mathbf{j}+\mathbf{r},\nu),(\mathbf{j},\eta)}$ associated with bond definition.
• V_mean::Vector{T}: Standard deviation of extended Hubbard interaction strength $V_{(\mathbf{j}+\mathbf{r},\nu),(\mathbf{j},\eta)}$ associated with bond definition.

ExtendedHubbardModel(;
    # KEYWORD ARGUMENTS
    model_geometry::ModelGeometry{D,T},
    ph_sym_form::Bool,
    V_bond::Vector{Bond{D}},
    V_mean::Vector{T},
    V_std::Vector{T} = zero(V_mean)
) where {T<:AbstractFloat, D}

Initialize and return an instance of the type ExtendedHubbardModel.

ExtendedHubbardParameters{T<:AbstractFloat}

Extended Hubbard interaction parameters for finite lattice.

Fields

• V::Vector{T}: Extended Hubbard interaction strength for each pair neighbors in the lattice.
• neighbor_table::Matrix{Int}: Neighbor table for extended Hubbard interactions on lattice.
• bond_ids::Vector{Int}: Bond IDs used to define extended Hubbard interactions.
• ph_sym_form::Bool: Whether particle-hole symmetric form of extended Hubbard interaction was used.

ExtendedHubbardParameters(;
    # KEYWORD ARGUMENTS
    extended_hubbard_model::ExtendedHubbardModel{T},
    model_geometry::ModelGeometry{D,T},
    rng::AbstractRNG
) where {D, T<:AbstractFloat}

Initialize an instance of the ExtendedHubbardParameters type.

initialize!(
    fermion_path_integral_up::FermionPathIntegral,
    fermion_path_integral_dn::FermionPathIntegral,
    extended_hubbard_parameters::ExtendedHubbardParameters
)

initialize!(
    fermion_path_integral::FermionPathIntegral,
    extended_hubbard_parameters::ExtendedHubbardParameters
)

Initialize the contribution from the Hubbard interaction to a FermionPathIntegral instance.

measure_ext_hub_energy(
    ext_hub_params::ExtendedHubbardParameters{E},
    Gup::Matrix{T}, Gdn::Matrix{T},
    ext_hub_id::Int
) where {T<:Number, E<:AbstractFloat}

Measure the extended Hubbard interaction energy

\[V (\hat{n}_i-1)(\hat{n}_j-1)\]

if ph_sym_form = true and

\[V \hat{n}_i \hat{n}_j\]

if ph_sym_form = false for the specified EXT_HUB_ID.

abstract type AbstractHST{T<:Number, R<:AbstractFloat} end

Abstract type to represent a Hubbard-Stratonovich transformation. Here T is the effective Hubbard-Stratonovich field type, specifying whether the Hubbard-Stratonovich transformation is real or complex.

abstract type AbstractSymHST{T, R} <: AbstractHST{T, R} end

Abstract type to represent a Hubbard-Stratonovich transformation that couples to each spin species symmetrically. Here T is the effective Hubbard-Stratonovich field type, specifying whether the Hubbard-Stratonovich transformation is real or complex.

abstract type AbstractAsymHST{T, R} <: AbstractHST{T, R} end

Abstract type to represent a Hubbard-Stratonovich transformation that couples to each spin species asymmetrically. Here T is the effective Hubbard-Stratonovich field type, specifying whether the Hubbard-Stratonovich transformation is real or complex.

initialize!(
    fermion_path_integral_up::FermionPathIntegral{H},
    fermion_path_integral_dn::FermionPathIntegral{H},
    hst_parameters::AbstractHST{T}
) where {H<:Number, T<:Number}

initialize!(
    fermion_path_integral_up::FermionPathIntegral{H},
    fermion_path_integral_dn::FermionPathIntegral{H},
    hst_parameters::Tuple{Vararg{<:AbstractHST{T},N}}
) where {H<:Number, T<:Number, N}

initialize!(
    fermion_path_integral::FermionPathIntegral{H},
    hst_parameters::AbstractSymHST{T}
) where {H<:Number, T<:Number}

initialize!(
    fermion_path_integral::FermionPathIntegral{H},
    hst_parameters::Tuple{Vararg{<:AbstractSymHST{T},N}}
) where {H<:Number, T<:Number, N}

Initialize a FermionPathIntegral integral type to reflect the the current Hubbard-Stratonovich transformation type represented by hst_parameters.

local_updates!(
    # ARGUMENTS
    Gup::Matrix{H}, logdetGup::R, sgndetGup::H,
    Gdn::Matrix{H}, logdetGdn::R, sgndetGdn::H,
    hst_parameters::AbstractHST{T,R};
    # KEYWORD ARGUMENTS
    fermion_path_integral_up::FermionPathIntegral{H},
    fermion_path_integral_dn::FermionPathIntegral{H},
    fermion_greens_calculator_up::FermionGreensCalculator{H},
    fermion_greens_calculator_dn::FermionGreensCalculator{H},
    Bup::Vector{P}, Bdn::Vector{P},
    δG::R, δθ::R,  rng::AbstractRNG,
    δG_max::R = 1e-6,
    update_stabilization_frequency::Bool = true
) where {H<:Number, T<:Number, R<:Real, P<:AbstractPropagator}

Perform local updates to Hubbard-Stratonovich fields stored in hst_parameters. This method returns a tuple containing (acceptance_rate, logdetGup, sgndetGup, logdetGdn, sgndetGdn, δG, δθ).

Arguments

• Gup::Matrix{H}: Spin-up equal-time Green's function matrix.
• logdetGup::R: The log of the absolute value of the determinant of the spin-up equal-time Green's function matrix, $\log \vert \det G_\uparrow(\tau,\tau) \vert.$
• sgndetGup::H: The sign/phase of the determinant of the spin-up equal-time Green's function matrix, $\det G_\uparrow(\tau,\tau) / \vert \det G_\uparrow(\tau,\tau) \vert.$
• Gdn::Matrix{H}: Spin-down equal-time Green's function matrix.
• logdetGdn::R: The log of the absolute value of the determinant of the spin-down equal-time Green's function matrix, $\log \vert \det G_\downarrow(\tau,\tau) \vert.$
• sgndetGdn::H: The sign/phase of the determinant of the spin-down equal-time Green's function matrix, $\det G_\downarrow(\tau,\tau) / \vert \det G_\downarrow(\tau,\tau) \vert.$
• hst_parameters::AbstractHST{T,R}: Type representing Hubbard-Stratonovich transformation.

Keyword Arguments

• fermion_path_integral_up::FermionPathIntegral{H}: An instance of the FermionPathIntegral type for spin-up electrons.
• fermion_path_integral_dn::FermionPathIntegral{H}: An instance of the FermionPathIntegral type for spin-down electrons.
• fermion_greens_calculator_up::FermionGreensCalculator{H}: An instance of the FermionGreensCalculator type for the spin-up electrons.
• fermion_greens_calculator_dn::FermionGreensCalculator{H}: An instance of the FermionGreensCalculator type for the spin-down electrons.
• Bup::Vector{P}: Spin-up propagators for each imaginary time slice.
• Bdn::Vector{P}: Spin-dn propagators for each imaginary time slice.
• δG_max::R: Maximum allowed error corrected by numerical stabilization.
• δG::R: Previously recorded maximum error in the Green's function corrected by numerical stabilization.
• δθ::R: Previously recorded maximum error in the sign/phase of the determinant of the equal-time Green's function matrix corrected by numerical stabilization.
• rng::AbstractRNG: Random number generator used in method instead of global random number generator, important for reproducibility.
• update_stabilization_frequency::Bool = true: If true, allows the stabilization frequency n_stab to be dynamically adjusted.

local_updates!(
    # ARGUMENTS
    Gup::Matrix{H}, logdet,Gup::R, sgndetGup::H,
    Gdn::Matrix{H}, logdetGdn::R, sgndetGdn::H,
    hst_parameters::Tuple{Vararg{<:AbstractHST{T,R},N}};
    # KEYWORD ARGUMENTS
    fermion_path_integral_up::FermionPathIntegral{H},
    fermion_path_integral_dn::FermionPathIntegral{H},
    fermion_greens_calculator_up::FermionGreensCalculator{H},
    fermion_greens_calculator_dn::FermionGreensCalculator{H},
    Bup::Vector{P}, Bdn::Vector{P},
    δG::R, δθ::R,  rng::AbstractRNG,
    δG_max::R = 1e-6,
    update_stabilization_frequency::Bool = true
) where {H<:Number, T<:Number, R<:Real, P<:AbstractPropagator, N}

Perform local updates to Hubbard-Stratonovich fields for N different types of Hubbard-Stratonovich transformations. This method returns a tuple containing (acceptance_rates, logdetGup, sgndetGup, logdetGdn, sgndetGdn, δG, δθ). Note that acceptance_rates is a tuple returning the acceptance rate for local updates of each type of Hubbard-Stratonovich field that was sampled.

Arguments

• Gup::Matrix{H}: Spin-up equal-time Green's function matrix.
• logdetGup::R: The log of the absolute value of the determinant of the spin-up equal-time Green's function matrix, $\log \vert \det G_\uparrow(\tau,\tau) \vert.$
• sgndetGup::H: The sign/phase of the determinant of the spin-up equal-time Green's function matrix, $\det G_\uparrow(\tau,\tau) / \vert \det G_\uparrow(\tau,\tau) \vert.$
• Gdn::Matrix{H}: Spin-down equal-time Green's function matrix.
• logdetGdn::R: The log of the absolute value of the determinant of the spin-down equal-time Green's function matrix, $\log \vert \det G_\downarrow(\tau,\tau) \vert.$
• sgndetGdn::H: The sign/phase of the determinant of the spin-down equal-time Green's function matrix, $\det G_\downarrow(\tau,\tau) / \vert \det G_\downarrow(\tau,\tau) \vert.$
• hst_parameters::Tuple{Vararg{<:AbstractHST{T,R},N}}: Tuple of parameters for multiple different Hubbard-Stratonovich transformation fields that will be sampled.

Keyword Arguments

• fermion_path_integral_up::FermionPathIntegral{H}: An instance of the FermionPathIntegral type for spin-up electrons.
• fermion_path_integral_dn::FermionPathIntegral{H}: An instance of the FermionPathIntegral type for spin-down electrons.
• fermion_greens_calculator_up::FermionGreensCalculator{H}: An instance of the FermionGreensCalculator type for the spin-up electrons.
• fermion_greens_calculator_dn::FermionGreensCalculator{H}: An instance of the FermionGreensCalculator type for the spin-down electrons.
• Bup::Vector{P}: Spin-up propagators for each imaginary time slice.
• Bdn::Vector{P}: Spin-dn propagators for each imaginary time slice.
• δG_max::R: Maximum allowed error corrected by numerical stabilization.
• δG::R: Previously recorded maximum error in the Green's function corrected by numerical stabilization.
• δθ::R: Previously recorded maximum error in the sign/phase of the determinant of the equal-time Green's function matrix corrected by numerical stabilization.
• rng::AbstractRNG: Random number generator used in method instead of global random number generator, important for reproducibility.
• update_stabilization_frequency::Bool = true: If true, allows the stabilization frequency n_stab to be dynamically adjusted.

local_updates!(
    # ARGUMENTS
    G::Matrix{H}, logdetG::R, sgndetG::H,
    hst_parameters::AbstractSymHST{T,R};
    # KEYWORD ARGUMENTS
    fermion_path_integral::FermionPathIntegral{H},
    fermion_greens_calculator::FermionGreensCalculator{H},
    B::Vector{P},
    δG::R, δθ::R,  rng::AbstractRNG,
    δG_max::R = 1e-6,
    update_stabilization_frequency::Bool = true
) where {H<:Number, T<:Number, R<:Real, P<:AbstractPropagator}

Perform local updates to Hubbard-Stratonovich fields for a spin-symmetric (density channel) Hubbard-Stratonovich transformation. This method returns a tuple containing (acceptance_rate, logdetG, sgndetG, δG, δθ).

Arguments

• G::Matrix{H}: Equal-time Green's function matrix.
• logdetG::R: The log of the absolute value of the determinant of the equal-time Green's function matrix, $\log \vert \det G(\tau,\tau) \vert.$
• sgndetG::H: The sign/phase of the determinant of the equal-time Green's function matrix, $\det G(\tau,\tau) / \vert \det G(\tau,\tau) \vert.$
• hst_parameters::AbstractSymHST{T,R}: Type representing spin-symmetric Hubbard-Stratonovich transformation.

Keyword Arguments

• fermion_path_integral::FermionPathIntegral{H}: An instance of the FermionPathIntegral.
• fermion_greens_calculator::FermionGreensCalculator{H}: An instance of the FermionGreensCalculator type.
• B::Vector{P}: Propagators for each imaginary time slice.
• δG_max::R: Maximum allowed error corrected by numerical stabilization.
• δG::R: Previously recorded maximum error in the Green's function corrected by numerical stabilization.
• δθ::R: Previously recorded maximum error in the sign/phase of the determinant of the equal-time Green's function matrix corrected by numerical stabilization.
• rng::AbstractRNG: Random number generator used in method instead of global random number generator, important for reproducibility.
• update_stabilization_frequency::Bool = true: If true, allows the stabilization frequency n_stab to be dynamically adjusted.

local_updates!(
    # ARGUMENTS
    G::Matrix{H}, logdetG::R, sgndetG::H,
    hst_parameters::Tuple{Vararg{<:AbstractSymHST{T,R},N}};
    # KEYWORD ARGUMENTS
    fermion_path_integral::FermionPathIntegral{H},
    fermion_greens_calculator::FermionGreensCalculator{H},
    B::Vector{P},
    δG::R, δθ::R,  rng::AbstractRNG,
    δG_max::R = 1e-6,
    update_stabilization_frequency::Bool = true
) where {H<:Number, T<:Number, R<:Real, P<:AbstractPropagator, N}

Perform local updates to multiple types Hubbard-Stratonovich fields for a spin-symmetric (density channel) Hubbard-Stratonovich transformation. This method returns a tuple containing (acceptance_rate, logdetG, sgndetG, δG, δθ).

Arguments

• G::Matrix{H}: Equal-time Green's function matrix.
• logdetG::R: The log of the absolute value of the determinant of the equal-time Green's function matrix, $\log \vert \det G(\tau,\tau) \vert.$
• sgndetG::H: The sign/phase of the determinant of the equal-time Green's function matrix, $\det G(\tau,\tau) / \vert \det G(\tau,\tau) \vert.$
• hst_parameters::Tuple{Vararg{<:AbstractSymHST{T,R},N}}: Tuple of parameters for multiple different spin-symmetric Hubbard-Stratonovich transformation fields that will be sampled.

Keyword Arguments

• fermion_path_integral::FermionPathIntegral{H}: An instance of the FermionPathIntegral.
• fermion_greens_calculator::FermionGreensCalculator{H}: An instance of the FermionGreensCalculator type.
• B::Vector{P}: Propagators for each imaginary time slice.
• δG_max::R: Maximum allowed error corrected by numerical stabilization.
• δG::R: Previously recorded maximum error in the Green's function corrected by numerical stabilization.
• δθ::R: Previously recorded maximum error in the sign/phase of the determinant of the equal-time Green's function matrix corrected by numerical stabilization.
• rng::AbstractRNG: Random number generator used in method instead of global random number generator, important for reproducibility.
• update_stabilization_frequency::Bool = true: If true, allows the stabilization frequency n_stab to be dynamically adjusted.

reflection_update!(
    # ARGUMENTS
    Gup::Matrix{H}, logdetGup::R, sgndetGup::H,
    Gdn::Matrix{H}, logdetGdn::R, sgndetGdn::H,
    hst_parameters::AbstractHST{T,R};
    # KEYWORD ARGUMENTS
    fermion_path_integral_up::FermionPathIntegral{H},
    fermion_path_integral_dn::FermionPathIntegral{H},
    fermion_greens_calculator_up::FermionGreensCalculator{H,R},
    fermion_greens_calculator_dn::FermionGreensCalculator{H,R},
    fermion_greens_calculator_up_alt::FermionGreensCalculator{H,R},
    fermion_greens_calculator_dn_alt::FermionGreensCalculator{H,R},
    Bup::Vector{P}, Bdn::Vector{P},
    rng::AbstractRNG
) where {H<:Number, T<:Number, R<:Real, P<:AbstractPropagator}

Perform a reflection update for a Hubbard-Stratonovich field on a randomly chosen location in the lattice for all imaginary-time slices. This function returns (accepted, logdetGup, sgndetGup, logdetGdn, sgndetGdn).

Arguments

• Gup::Matrix{H}: Spin-up equal-time Greens function matrix.
• logdetGup::R: Log of the determinant of the spin-up equal-time Greens function matrix.
• sgndetGup::H: Sign/phase of the determinant of the spin-up equal-time Greens function matrix.
• Gdn::Matrix{H}: Spin-down equal-time Greens function matrix.
• logdetGdn::R: Log of the determinant of the spin-down equal-time Greens function matrix.
• sgndetGdn::H: Sign/phase of the determinant of the spin-down equal-time Greens function matrix.
• hst_parameters::AbstractHST{T,R}: Hubbard-Stratonovich fields and associated parameters to update.

Keyword Arguments

• fermion_path_integral_up::FermionPathIntegral{H}: An instance of FermionPathIntegral type for spin-up electrons.
• fermion_path_integral_dn::FermionPathIntegral{H}: An instance of FermionPathIntegral type for spin-down electrons.
• fermion_greens_calculator_up::FermionGreensCalculator{H,R}: Contains matrix factorization information for current spin-up sector state.
• fermion_greens_calculator_dn::FermionGreensCalculator{H,R}: Contains matrix factorization information for current spin-down sector state.
• fermion_greens_calculator_up_alt::FermionGreensCalculator{H,R}: Used to calculate matrix factorizations for proposed spin-up sector state.
• fermion_greens_calculator_dn_alt::FermionGreensCalculator{H,R}: Used to calculate matrix factorizations for proposed spin-down sector state.
• Bup::Vector{P}: Spin-up propagators for each imaginary time slice.
• Bdn::Vector{P}: Spin-down propagators for each imaginary time slice.
• rng::AbstractRNG: Random number generator used in method instead of global random number generator, important for reproducibility.

reflection_update!(
    # ARGUMENTS
    G::Matrix{H}, logdetG::R, sgndetG::H,
    hst_parameters::AbstractSymHST{T,R};
    # KEYWORD ARGUMENTS
    fermion_path_integral::FermionPathIntegral{H},
    fermion_greens_calculator::FermionGreensCalculator{H,R},
    fermion_greens_calculator_alt::FermionGreensCalculator{H,R},
    B::Vector{P},
    rng::AbstractRNG
) where {H<:Number, T<:Number, R<:Real, P<:AbstractPropagator}

Perform a reflection update for a spin-symmetric (density channel) Hubbard-Stratonovich field on a randomly chosen location in the lattice for all imaginary-time slices. This function returns (accepted, logdetG, sgndetG).

Arguments

• G::Matrix{H}: The current Hubbard-Stratonovich field matrix.
• logdetG::R: Log of the determinant of the Hubbard-Stratonovich field matrix.
• sgndetG::H: Sign/phase of the determinant of the Hubbard-Stratonovich field matrix.
• hst_parameters::AbstractSymHST{T,R}: Hubbard-Stratonovich fields and associated parameters to update.

Keyword Arguments

• fermion_path_integral::FermionPathIntegral{H}: An instance of FermionPathIntegral.
• fermion_greens_calculator::FermionGreensCalculator{H,R}: Contains matrix factorization information for current state.
• fermion_greens_calculator_alt::FermionGreensCalculator{H,R}: Used to calculate matrix factorizations for proposed state.
• B::Vector{P}: Propagators for each imaginary time slice.
• rng::AbstractRNG: Random number generator used in method instead of global random number

swap_update!(
    # ARGUMENTS
    Gup::Matrix{H}, logdetGup::R, sgndetGup::H,
    Gdn::Matrix{H}, logdetGdn::R, sgndetGdn::H,
    hst_parameters::AbstractHST{T,R};
    # KEYWORD ARGUMENTS
    fermion_path_integral_up::FermionPathIntegral{H},
    fermion_path_integral_dn::FermionPathIntegral{H},
    fermion_greens_calculator_up::FermionGreensCalculator{H,R},
    fermion_greens_calculator_dn::FermionGreensCalculator{H,R},
    fermion_greens_calculator_up_alt::FermionGreensCalculator{H,R},
    fermion_greens_calculator_dn_alt::FermionGreensCalculator{H,R},
    Bup::Vector{P}, Bdn::Vector{P},
    rng::AbstractRNG
) where {H<:Number, T<:Number, R<:Real, P<:AbstractPropagator}

Perform a swap update for a Hubbard-Stratonovich field between a pair of randomly chosen locations in the lattice for all imaginary-time slices. This function returns (accepted, logdetGup, sgndetGup, logdetGdn, sgndetGdn).

Arguments

• Gup::Matrix{H}: Spin-up equal-time Greens function matrix.
• logdetGup::R: Log of the determinant of the spin-up equal-time Greens function matrix.
• sgndetGup::H: Sign/phase of the determinant of the spin-up equal-time Greens function matrix.
• Gdn::Matrix{H}: Spin-down equal-time Greens function matrix.
• logdetGdn::R: Log of the determinant of the spin-down equal-time Greens function matrix.
• sgndetGdn::H: Sign/phase of the determinant of the spin-down equal-time Greens function matrix.
• hst_parameters::AbstractHST{T,R}: Hubbard-Stratonovich fields and associated parameters to update.

Keyword Arguments

• fermion_path_integral_up::FermionPathIntegral{H}: An instance of FermionPathIntegral type for spin-up electrons.
• fermion_path_integral_dn::FermionPathIntegral{H}: An instance of FermionPathIntegral type for spin-down electrons.
• fermion_greens_calculator_up::FermionGreensCalculator{H,R}: Contains matrix factorization information for current spin-up sector state.
• fermion_greens_calculator_dn::FermionGreensCalculator{H,R}: Contains matrix factorization information for current spin-down sector state.
• fermion_greens_calculator_up_alt::FermionGreensCalculator{H,R}: Used to calculate matrix factorizations for proposed spin-up sector state.
• fermion_greens_calculator_dn_alt::FermionGreensCalculator{H,R}: Used to calculate matrix factorizations for proposed spin-down sector state.
• Bup::Vector{P}: Spin-up propagators for each imaginary time slice.
• Bdn::Vector{P}: Spin-down propagators for each imaginary time slice
• rng::AbstractRNG: Random number generator used in method instead of global random number generator, important for reproducibility.

swap_update!(
    # ARGUMENTS
    G::Matrix{H}, logdetG::R, sgndetG::H,
    hst_parameters::AbstractSymHST{T,R};
    # KEYWORD ARGUMENTS
    fermion_path_integral::FermionPathIntegral{H},
    fermion_greens_calculator::FermionGreensCalculator{H,R},
    fermion_greens_calculator_alt::FermionGreensCalculator{H,R},
    B::Vector{P},
    rng::AbstractRNG
) where {H<:Number, T<:Number, R<:Real, P<:AbstractPropagator}

Perform a swap update for a spin-symmetric (density channel) Hubbard-Stratonovich field between a pair of randomly chosen locations in the lattice for all imaginary-time slices. This function returns (accepted, logdetG, sgndetG).

Arguments

• G::Matrix{H}: The current Hubbard-Stratonovich field matrix.
• logdetG::R: Log of the determinant of the Hubbard-Stratonovich field
• sgndetG::H: Sign of the determinant of the Hubbard-Stratonovich field matrix.
• hst_parameters::AbstractSymHST{T,R}: Hubbard-Stratonovich fields and associated parameters to update.

Keyword Arguments

• fermion_path_integral::FermionPathIntegral{H}: An instance of FermionPathIntegral.
• fermion_greens_calculator::FermionGreensCalculator{H,R}: Contains matrix factorization information for current state.
• fermion_greens_calculator_alt::FermionGreensCalculator{H,R}: Used to calculate matrix factorizations for proposed state.
• B::Vector{P}: Propagators for each imaginary time slice.
• rng::AbstractRNG: Random number generator used in method instead of global random number generator, important for reproducibility.

HubbardSpinHirschHST{T,R} <: AbstractAsymHST{T,R}

This type represents a Hubbard-Stratonovich (HS) transformation for decoupling the local Hubbard interaction in the spin channel, where the introduced HS fields take on the two discrete values $s = \pm 1$. Specifically, the Hubbard interaction is decoupled as

\[e^{-\Delta\tau U\left(n_{\uparrow}-\tfrac{1}{2}\right)\left(n_{\downarrow}-\tfrac{1}{2}\right)} = \gamma\sum_{s=\pm1}e^{-\Delta\tau\alpha(n_{\uparrow}-n_{\downarrow})s},\]

where

\[\gamma=\frac{1}{2}e^{-\Delta\tau U/4}\]

and

\[\alpha = \frac{1}{\Delta\tau}\cosh^{-1}\left(e^{\Delta\tau U/2}\right).\]

Note that when $U \ge 0$ then $\alpha$ is real, whereas if $U<0$ then $\alpha$ is purely imaginary.

HubbardSpinHirschHST(;
    # KEYWORD ARGUMENTS
    hubbard_parameters::HubbardParameters{R},
    β::R, Δτ::R, rng::AbstractRNG
) where {R<:AbstractFloat}

Initialize an instance of the HubbardSpinHirschHST type.

HubbardSpinGaussHermiteHST{T,R} <: AbstractAsymHST{T,R}

This type represents a Hubbard-Stratonovich (HS) transformation for decoupling the local Hubbard interaction in the spin channel, where the introduced HS fields take on the four discrete values $s \in \{ -2, -1, +1, +2 \}.$ Note that the Hubbard interaction can be written in the spin channel as

\[U(\hat{n}_{\uparrow}-\tfrac{1}{2})(\hat{n}_{\downarrow}-\tfrac{1}{2})=-\tfrac{U}{2}(\hat{n}_{\uparrow}-\hat{n}_{\downarrow})^{2}+\tfrac{U}{4},\]

different only by a constant energy offset $U/4$ which does not matter. Therefore, we can perform a Gauss-Hermite Hubbard-Statonovich transformation in the spin channel as

\[e^{-\Delta\tau\left[-\frac{U}{2}\right](\hat{n}_{\uparrow}-\hat{n}_{\downarrow})^{2}} = \frac{1}{4}\sum_{s=\pm1,\pm2}e^{-S_{\text{GH}}(s)-\Delta\tau\hat{V}(s)}+\mathcal{O}\left((\Delta\tau U)^{4}\right),\]

where $\hat{V}(s)=\alpha\eta(s)(\hat{n}_{\uparrow}-\hat{n}_{\downarrow})$ and $\alpha = \sqrt{U/(2\Delta\tau)}$. In the above expression,

\[S_{\text{GH}}(s)=-\log\left(1+\sqrt{6}\left(1-\tfrac{2}{3}|s|\right)\right)\]

and

\[\eta(s)=\frac{s}{|s|}\sqrt{6(1-\sqrt{6})+4\sqrt{6}|s|}.\]

Note that $\alpha$ is strictly real when $U \ge 0$ and strictly imaginary when $U < 0$.

HubbardSpinGaussHermiteHST(;
    # KEYWORD ARGUMENTS
    hubbard_parameters::HubbardParameters{E},
    β::E, Δτ::E, rng::AbstractRNG
) where {E<:AbstractFloat}

Initialize an instance of the HubbardSpinGaussHermiteHST type.

HubbardDensityHirschHST{T,R} <: AbstractSymHST{T,R}

This type represents a Hubbard-Stratonovich (HS) transformation for decoupling the local Hubbard interaction in the density channel, where the introduced HS fields take on the two discrete values $s = \pm 1$. Specifically, the Hubbard interaction is decoupled as

\[e^{-\Delta\tau U\left(n_{\uparrow}-\tfrac{1}{2}\right)\left(n_{\downarrow}-\tfrac{1}{2}\right)} = \gamma\sum_{s=\pm1}e^{-\Delta\tau\alpha(n_{\uparrow}+n_{\downarrow}-1)s},\]

where

\[\gamma = \frac{1}{2}e^{\Delta\tau U/4}\]

and

\[\alpha = \frac{1}{\Delta\tau}\cosh\left(e^{-\Delta\tau U/2}\right).\]

Note that when $U \le 0$ then $\alpha$ is real, whereas if $U > 0$ then $\alpha$ is purely imaginary.

HubbardDensityHirschHST(;
    # KEYWORD ARGUMENTS
    hubbard_parameters::HubbardParameters{E},
    β::E, Δτ::E, rng::AbstractRNG
) where {E<:AbstractFloat}

Initialize an instance of the HubbardDensityHirschHST type.

HubbardDensityGaussHermiteHST{T,R} <: AbstractSymHST{T,R}

This type represents a Hubbard-Stratonovich (HS) transformation for decoupling the local Hubbard interaction in the density channel, where the introduced HS fields take on the four discrete values $s \in \{ -2, -1, +1, +2 \}.$ Note that the Hubbard interaction can be written in the density channel as

\[U(\hat{n}_{\uparrow}-\tfrac{1}{2})(\hat{n}_{\downarrow}-\tfrac{1}{2})=-\tfrac{U}{2}(\hat{n}_{\uparrow}+\hat{n}_{\downarrow}-1)^{2}-\tfrac{U}{4},\]

different only by a constant energy offset $-U/4$ which does not matter. Therefore, we can perform a Gauss-Hermite Hubbard-Statonovich transformation in the density channel as

\[e^{-\Delta\tau\left[\frac{U}{2}\right](\hat{n}_{\uparrow}+\hat{n}_{\downarrow}-1)^{2}} = \frac{1}{4}\sum_{s=\pm1,\pm2}e^{-S_{\text{GF}}(s)-\Delta\tau\hat{V}(s)}+\mathcal{O}\left((\Delta\tau U)^{4}\right)\]

where $\hat{V}(s)=\alpha\eta(s)(\hat{n}_{\uparrow}+\hat{n}_{\downarrow}-1)$ and $\alpha = \sqrt{-U/(2\Delta\tau)}$. In the above expression,

\[S_{\text{GH}}(s)=-\log\left(1+\sqrt{6}\left(1-\tfrac{2}{3}|s|\right)\right)\]

and

\[\eta(s)=\frac{s}{|s|}\sqrt{6(1-\sqrt{6})+4\sqrt{6}|s|}.\]

Note that $\alpha$ is strictly real when $U \le 0$ and strictly imaginary when $U > 0$.

HubbardDensityGaussHermiteHST(;
    # KEYWORD ARGUMENTS
    hubbard_parameters::HubbardParameters{E},
    β::E, Δτ::E, rng::AbstractRNG
) where {E<:AbstractFloat}

Initialize an instance of the HubbardDensityGaussHermiteHST type.

ExtHubSpinHirschHST{T,R} <: AbstractAsymHST{T,R}

This type represent a Hirsch style Hubbard-Stratonovich (HS) transformation used to decouple an extended Hubbard interaction, where the introduced HS fields can take on the values $s = \pm 1$. The decomposition is done using the spin-channel. Here the extended Hubbard interaction is expressed as

\[V(\hat{n}_i-1)(\hat{n}_j) = V\sum_{\sigma,\sigma'} (\hat{n}_{\sigma,i}-\frac{1}{2})\hat{n}_{\sigma',j}-\frac{1}{2}).\]

Then each of the four terms on the right is decouple using a HS transformation of the form

\[e^{-\Delta\tau V\left(n_{\sigma}-\tfrac{1}{2}\right)\left(n_{\sigma'}-\tfrac{1}{2}\right)} = \gamma\sum_{s=\pm1}e^{-\Delta\tau\alpha(n_{\sigma}-n_{\sigma'})s},\]

where

\[\gamma=\frac{1}{2}e^{-\Delta\tau V/4}\]

and

\[\alpha = \frac{1}{\Delta\tau}\cosh^{-1}\left(e^{\Delta\tau V/2}\right).\]

Note that when $V \ge 0$ then $\alpha$ is real, whereas if $V<0$ then $\alpha$ is purely imaginary.

ExtHubSpinHirschHST(;
    # KEYWORD ARGUMENTS
    extended_hubbard_parameters::ExtendedHubbardParameters{R},
    β::R, Δτ::R, rng::AbstractRNG
) where {R<:AbstractFloat}

Initialize an instance of the ExtHubSpinHirschHST type.

ExtHubDensityGaussHermiteHST{T,R} <: AbstractSymHST{T,R}

This type represents a Hubbard-Stratonovich (HS) transformation for decoupling the extended Hubbard interaction, where the introduced HS fields take on the four discrete values $s \in \{ -2, -1, +1, +2 \}.$

Specifically, we perform the Gauss-Hermite Hubbard-Stratonovich transformation

\[e^{-\Delta\tau\left[\tfrac{V}{2}\right](\hat{n}_{\mathbf{i}}+\hat{n}_{\mathbf{j}}-2)^{2}} = = \frac{1}{4}\sum_{s=\pm1,\pm2}e^{-S_{\text{GH}}(s)-\Delta\tau\hat{V}(s)}+\mathcal{O}\left(\left[\tfrac{\Delta\tau V}{2}\right]^{4}\right)\]

where $\hat{V}(s) = \alpha\eta(s)(\hat{n}_{\mathbf{i}}+\hat{n}_{\mathbf{j}}-2)$ and $\alpha=\sqrt{\frac{-V}{2\Delta\tau}}$. In the above expression,

\[S_{\text{GH}}(s)=-\log\left(1+\sqrt{6}\left(1-\tfrac{2}{3}|s|\right)\right)\]

and

\[\eta(s)=\frac{s}{|s|}\sqrt{6(1-\sqrt{6})+4\sqrt{6}|s|}.\]

Note that $\alpha$ is strictly real when $V \le 0$ and strictly imaginary when $V > 0$.

ExtHubDensityGaussHermiteHST(;
    # KEYWORD ARGUMENTS
    extended_hubbard_parameters::ExtendedHubbardParameters{R},
    β::R, Δτ::R, rng::AbstractRNG
) where {R<:AbstractFloat}

Initialize an instance of the ExtHubDensityGaussHermiteHST type.

init_renormalized_hubbard_parameters(;
    # KEYWORD ARGUMENTS
    hubbard_parameters::HubbardParameters{E},
    hst_parameters::ExtHubDensityGaussHermiteHST{T,E},
    model_geometry::ModelGeometry{D,E}
) where {D, T<:Number, E<:AbstractFloat}

Returns a new instance of the HubbardParameters type with the Hubbard interactions renormalized based on the ExtHubDensityGaussHermiteHST definition. Refer to the definition of the ExtendedHubbardModel to see where this renormalization of the local Hubbard interaction comes from.

Note that either both the local and extended Hubbard interactions need to be initialized using the particle-hole symmetric or asymmetric form for the interaction (as specified by ph_sym_form keyword argument), and cannot use opposite conventions. Additionally, the HubbardModel definition used to create the hubbard_parameters instance of the HubbardParameters passed to this function must initialize a Hubbard interaction on each type of orbital species/ID appearing in an extended Hubbard interaction, even if this means initializing the local Hubbard interaction to $U = 0$.

ElectronPhononModel{T<:Number, E<:AbstractFloat, D}

Defines an electron-phonon model.

Fields

• phonon_modes::Vector{PhononModes{E,D}}: A vector of PhononMode definitions.
• phonon_dispersions::Vector{PhononDispersion{E,D}}: A vector of PhononDispersion definitions.
• holstein_couplings_up::Vector{HolsteinCoupling{E,D}}: A vector of HolsteinCoupling definitions for spin-up.
• holstein_couplings_dn::Vector{HolsteinCoupling{E,D}}: A vector of HolsteinCoupling definitions for spin-down.
• ssh_couplings_up::Vector{SSHCoupling{T,E,D}}: A vector of SSHCoupling definitions for spin-up.
• ssh_couplings_dn::Vector{SSHCoupling{T,E,D}}: A vector of SSHCoupling definitions for spin-down.

ElectronPhononModel(;
    # KEYWORD ARGUMENTS
    model_geometry::ModelGeometry{D,E},
    tight_binding_model::Union{TightBindingModel{T,E,D}, Nothing} = nothing,
    tight_binding_model_up::Union{TightBindingModel{T,E,D}, Nothing} = nothing,
    tight_binding_model_dn::Union{TightBindingModel{T,E,D}, Nothing} = nothing
) where {T<:Number, E<:AbstractFloat, D}

Initialize and return a null (empty) instance of ElectronPhononModel. Note that either tight_binding_model or tight_binding_model_up and tight_binding_model_dn needs to be specified.

PhononMode{E<:AbstractFloat}

Defines a phonon mode $\nu$ at location \mathbf{r}_\nu in the unit cell. Specifically, it defines the phonon Hamiltonian terms

\[\hat{H}_{{\rm ph}} = \sum_{\mathbf{i}} \left[ \frac{1}{2} M_{\mathbf{i},\nu}\Omega_{\mathbf{i},\nu}^{2}\hat{X}_{\mathbf{i},\nu}^{2} + \frac{1}{12}M_{\mathbf{i},\nu}\Omega_{4,\mathbf{i},\nu}^{2}\hat{X}_{\mathbf{i},\nu}^{4} + \frac{1}{2M_{\mathbf{i},\nu}}\hat{P}_{\mathbf{i},\nu}^{2} \right],\]

where the sum runs over unit cell $\mathbf{i}$, $\mathbf{r}_\nu$ denotes the location of the phonon mode in the unit cell, $M_{\mathbf{i},\nu}$ is the phonon mass M, $\Omega_{\mathbf{i},\nu}$ is the phonon frequency that is distributed according to a normal distribution with mean Ω_mean and standard deviation Ω_std. Lastly, $\Omega_{4,\mathbf{i},\nu}$ is the anharmonic coefficient, and is distributed according to a normal distribution with mean Ω4_mean and standard deviation Ω4_std.

Fields

• basis_vec::SVector{D,E}: Location $\mathbf{r}_\nu$ of phonon mode in unit cell.
• M::E:: The phonon mass $M_{\mathbf{i},\nu}.$
• Ω_mean::E: Mean of normal distribution the phonon frequency $\Omega_{\mathbf{i},\nu}$ is sampled from.
• Ω_std::E: Standard deviation of normal distribution the phonon frequency $\Omega_{\mathbf{i},\nu}$ is sampled from.
• Ω4_mean::E: Mean of normal distribution the anharmonic coefficient $\Omega_{4,\mathbf{i},\nu}$ is sampled from.
• Ω4_std::E: Standard deviation of normal distribution the anharmonic coefficient $\Omega_{4,\mathbf{i},\nu}$ is sampled from.

PhononMode(;
    # KEYWORD ARGUMENTS
    basis_vec::AbstractVector{E},
    Ω_mean::E,
    Ω_std::E = 0.,
    M::E = 1.,
    Ω4_mean::E = 0.0,
    Ω4_std::E = 0.0,
) where {E<:AbstractFloat}

Initialize and return a instance of PhononMode.

HolsteinCoupling{E<:AbstractFloat, D}

Defines a Holstein coupling between a specified phonon mode and orbital density. Specifically, if ph_sym_form = true then a the particle-hole symmetric form of the Holstein coupling given by

\[\begin{align*} H = \sum_{\mathbf{i}} \Big[ & (\alpha_{\mathbf{i},(\mathbf{r},\kappa,\nu)} \hat{X}_{\mathbf{i},\nu} + \alpha_{3,\mathbf{i},(\mathbf{r},\kappa,\nu)} \hat{X}^3_{\mathbf{i},\nu}) \ (\hat{n}_{\sigma,\mathbf{i}+\mathbf{r},\kappa}-\tfrac{1}{2})\\ & + (\alpha_{2,\mathbf{i},(\mathbf{r},\kappa,\nu)} \hat{X}^2_{\mathbf{i},\nu} + \alpha_{4,\mathbf{i},(\mathbf{r},\kappa,\nu)} \hat{X}^4_{\mathbf{i},\nu}) \ \hat{n}_{\sigma,\mathbf{i}+\mathbf{r},\kappa} \Big] \end{align*},\]

is used, whereas if ph_sym_form = false Holstein coupling is given by

\[\begin{align*} H = \sum_{\mathbf{i}} \Big[ & (\alpha_{\mathbf{i},(\mathbf{r},\kappa,\nu)} \hat{X}_{\mathbf{i},\nu} + \alpha_{3,\mathbf{i},(\mathbf{r},\kappa,\nu)} \hat{X}^3_{\mathbf{i},\nu}) \ \hat{n}_{\sigma,\mathbf{i}+\mathbf{r},\kappa}\\ & + (\alpha_{2,\mathbf{i},(\mathbf{r},\kappa,\nu)} \hat{X}^2_{\mathbf{i},\nu} + \alpha_{4,\mathbf{i},(\mathbf{r},\kappa,\nu)} \hat{X}^4_{\mathbf{i},\nu}) \ \hat{n}_{\sigma,\mathbf{i}+\mathbf{r},\kappa} \Big] \end{align*}.\]

In the above, $\sigma$ specifies the sum, and the sum over $\mathbf{i}$ runs over unit cells in the lattice. In the above $\nu$ and $\kappa$ specify the phonon mode orbital species IDs respectively, and $\mathbf{r}$ is a static displacement in unit cells.

Fields

• ph_sym_form::Bool: If particle-hole symmetric form is used for Holstein coupling.
• phonon_id::Int: The ID $\nu$ specifying phonon mode getting coupled to.
• orbital_id::Int: The ID $\kappa$ specifying orbital species the phonon mode getting coupled to.
• displacement::SVector{D,Int}: Static displacement $r$ in unit cells in the direction of each lattice vector.
• α_mean::E: Mean of the linear Holstein coupling coefficient $\alpha_{1,\mathbf{i},(\mathbf{r},\kappa,\nu)}.$
• α_std::E: Standard deviation of the linear Holstein coupling coefficient $\alpha_{1,\mathbf{i},(\mathbf{r},\kappa,\nu)}.$
• α2_mean::E: Mean of the squared Holstein coupling coefficient $\alpha_{2,\mathbf{i},(\mathbf{r},\kappa,\nu)}.$
• α2_std::E: Standard deviation of the squared Holstein coupling coefficient $\alpha_{2,\mathbf{i},(\mathbf{r},\kappa,\nu)}.$
• α3_mean::E: Mean of the cubic Holstein coupling coefficient $\alpha_{3,\mathbf{i},(\mathbf{r},\kappa,\nu)}.$
• α3_std::E: Standard deviation of the cubic Holstein coupling coefficient $\alpha_{3,\mathbf{i},(\mathbf{r},\kappa,\nu)}.$
• α4_mean::E: Mean of the quartic Holstein coupling coefficient $\alpha_{4,\mathbf{i},(\mathbf{r},\kappa,\nu)}.$
• α4_std::E: Standard deviation of the quartic Holstein coupling coefficient $\alpha_{4,\mathbf{i},(\mathbf{r},\kappa,\nu)}.$

Comment

Note that the initial orbital bond.orbital[1] must match the orbital species associated with phonon mode PhononMode getting coupled to.

HolsteinCoupling(;
    # KEYWORD ARGUMENTS
    model_geometry::ModelGeometry{D,E},
    phonon_id::Int,
    orbital_id::Int,
    displacement::AbstractVector{Int},
    α_mean::E,
    α_std::E  = 0.0,
    α2_mean::E = 0.0,
    α2_std::E = 0.0,
    α3_mean::E = 0.0,
    α3_std::E = 0.0,
    α4_mean::E = 0.0,
    α4_std::E = 0.0,
    ph_sym_form::Bool = true
) where {D, E<:AbstractFloat}

Initialize and return a instance of HolsteinCoupling.

SSHCoupling{T<:Number, E<:AbstractFloat, D}

Defines a Su-Schrieffer-Heeger (SSH) coupling between a pair of phonon modes. Specifically, it defines the SSH interaction term

\[\hat{H}_{{\rm ssh}} = -\sum_{\sigma,\mathbf{i}} \left[ t_{\mathbf{i},(\mathbf{r},\kappa,\nu)} - \left( \sum_{n=1}^{4}\alpha_{n,\mathbf{i},(\mathbf{r},\kappa',\nu')} \left( \hat{X}_{\mathbf{i}+\mathbf{r},\kappa'} - \hat{X}_{\mathbf{i},\nu'}\right)^{n}\right) \right] \left( \hat{c}_{\sigma,\mathbf{i}+\mathbf{r},\kappa}^{\dagger}\hat{c}_{\sigma,\mathbf{i},\nu}+{\rm h.c.} \right),\]

where $\sigma$ specifies the sum, and the sum over $\mathbf{i}$ runs over unit cells in the lattice. In the above $\nu$ and $\kappa$ IDs specify orbital species in the unit cell, and $\kappa'$ and $\nu'$ IDs specify the phonon modes getting coupled to. Finally, $\mathbf{r}$ is a static displacement in unit cells in the direction of each lattice vector. In that above expression $t_{\mathbf{i},(\mathbf{r},\kappa,\nu)}$ is the bare hopping amplitude, which is not specified here.

Fields

• phonon_ids::NTuple{2,Int}: Pair of phonon modes getting coupled together.
• bond::Bond{D}: Bond separating the two orbitals getting coupled to, which are separated by $\mathbf{r} + (\mathbf{r}_\kappa - \mathbf{r}_\nu)$.
• bond_id::Int: Bond ID associated with the bond field.
• α_mean::T: Mean of the linear SSH coupling constant $\alpha_{1,\mathbf{i},(\mathbf{r},\kappa,\nu)}.$
• α_std::T: Standard deviation of the linear SSH coupling constant $\alpha_{1,\mathbf{i},(\mathbf{r},\kappa,\nu)}.$
• α2_mean::T: Mean of the quadratic SSH coupling constant $\alpha_{2,\mathbf{i},(\mathbf{r},\kappa,\nu)}.$
• α2_std::T: Standard deviation of the quadratic SSH coupling constant $\alpha_{2,\mathbf{i},(\mathbf{r},\kappa,\nu)}.$
• α3_mean::T: Mean of the cubic SSH coupling constant $\alpha_{3,\mathbf{i},(\mathbf{r},\kappa,\nu)}.$
• α3_std::T: Standard deviation of the cubic SSH coupling constant $\alpha_{3,\mathbf{i},(\mathbf{r},\kappa,\nu)}.$
• α4_mean::T: Mean of the quartic SSH coupling constant $\alpha_{4,\mathbf{i},(\mathbf{r},\kappa,\nu)}.$
• α4_std::T: Standard deviation of the quartic SSH coupling constant $\alpha_{4,\mathbf{i},(\mathbf{r},\kappa,\nu)}.$
• expniϕ::T: Twisted boundary conditions phase factor.

SSHCoupling(;
    # KEYWORD ARGUMENTS
    model_geometry::ModelGeometry{D,E},
    tight_binding_model::TightBindingModel{T,E,D},
    phonon_ids::NTuple{2,Int},
    bond::Bond{D},
    α_mean::Union{T,E},
    α_std::E  = 0.0,
    α2_mean::Union{T,E} = 0.0,
    α2_std::E = 0.0,
    α3_mean::Union{T,E} = 0.0,
    α3_std::E = 0.0,
    α4_mean::Union{T,E} = 0.0,
    α4_std::E = 0.0
) where {D, T<:Number, E<:AbstractFloat}

Initialize and return a instance of SSHCoupling.

PhononDispersion{E<:AbstractFloat, D}

Defines a dispersive phonon coupling between phonon modes. Specifically, it defines the dispersive phonon term

\[\hat{H}_{{\rm disp}} = \sum_{\mathbf{i}} \left( \frac{M_{\mathbf{i}+\mathbf{r},\kappa}M_{\mathbf{i},\nu}}{M_{\mathbf{i}+\mathbf{r},\kappa}+M_{\mathbf{i},\nu}} \right) \bigg[ \Omega_{\mathbf{i},(\mathbf{r},\kappa,\nu)}^{2}\Big(\hat{X}_{\mathbf{i}+\mathbf{r},\kappa}-\hat{X}_{\mathbf{i},\nu}\Big)^{2} +\frac{1}{12}\Omega_{4,\mathbf{i},(\mathbf{r},\kappa,\nu)}^{2}\Big(\hat{X}_{\mathbf{i}+\mathbf{r},\kappa}-\hat{X}_{\mathbf{i},\nu}\Big)^{4} \bigg]\]

where the sum over $\mathbf{i}$ runs over unit cells in the lattice. In the above $\nu$ and $\kappa$ IDs specify the phonon modes in the unit cell, and $\mathbf{r}$ is a static displacement in unit cells.

Fields

• phonon_ids::NTuple{2,Int}: ID's for pair of phonon modes getting coupled together.
• displacement::SVector{D,Int}: Static displacement $\mathbf{r}$ in unit cells separating the two phonon modes getting coupled.
• Ω_mean::E: Mean dispersive phonon frequency $\Omega_{\mathbf{i},(\mathbf{r},\kappa,\nu)}.$
• Ω_std::E: Standard deviation of dispersive phonon frequency $\Omega_{\mathbf{i},(\mathbf{r},\kappa,\nu)}.$
• Ω4_mean::E: Mean quartic dispersive phonon coefficient $\Omega_{4,\mathbf{i},(\mathbf{r},\kappa,\nu)}.$
• Ω4_std::E: Standard deviation of quartic dispersive phonon coefficient $\Omega_{4,\mathbf{i},(\mathbf{r},\kappa,\nu)}.$

PhononDispersion(;
    # KEYWORD ARGUMENTS
    model_geometry::ModelGeometry{D,E},
    phonon_ids::NTuple{2,Int},
    displacement::AbstractVector{Int},
    Ω_mean::E,
    Ω_std::E=0.0,
    Ω4_mean::E=0.0,
    Ω4_std::E=0.0
) where {E<:AbstractFloat, D}

Initialize and return a instance of PhononDispersion.

add_phonon_mode!(;
    # KEYWORD ARGUMENTS
    electron_phonon_model::ElectronPhononModel{T,E,D},
    phonon_mode::PhononMode{E,D}
) where {T<:Number, E<:AbstractFloat, D}

Add a PhononMode to an ElectronPhononModel.

add_holstein_coupling!(;
    # KEYWORD ARGUMENTS
    model_geometry::ModelGeometry{D,E},
    electron_phonon_model::ElectronPhononModel{T,E,D},
    holstein_coupling::Union{HolsteinCoupling{E,D}, Nothing} = nothing,
    holstein_coupling_up::Union{HolsteinCoupling{E,D}, Nothing} = nothing,
    holstein_coupling_dn::Union{HolsteinCoupling{E,D}, Nothing} = nothing
) where {T,E,D}

Add the HolsteinCoupling to an ElectronPhononModel. Note that either holstein_coupling or holstein_coupling_up and holstein_coupling_dn must be specified.

add_ssh_coupling!(;
    # KEYWORD ARGUMENTS
    electron_phonon_model::ElectronPhononModel{T,E,D},
    tight_binding_model::Union{TightBindingModel{T,E,D}, Nothing} = nothing,
    ssh_coupling::Union{SSHCoupling{T,E,D}, Nothing} = nothing,
    tight_binding_model_up::Union{TightBindingModel{T,E,D}, Nothing} = nothing,
    tight_binding_model_dn::Union{TightBindingModel{T,E,D}, Nothing} = nothing,
    ssh_coupling_up::Union{SSHCoupling{T,E,D}, Nothing} = nothing,
    ssh_coupling_dn::Union{SSHCoupling{T,E,D}, Nothing} = nothing
) where {T,E,D}

Add a SSHCoupling to an ElectronPhononModel. Note that either ssh_coupling and tight_binding_model or ssh_coupling_up, ssh_coupling_dn, tight_binding_model_up and tight_binding_model_dn need to be specified.

add_phonon_dispersion!(;
    # KEYWORD ARGUMENTS
    electron_phonon_model::ElectronPhononModel{T,E,D},
    phonon_dispersion::PhononDispersion{E,D},
    model_geometry::ModelGeometry{D,E}
) where {T,E,D}

Add a PhononDispersion to an ElectronPhononModel.

ElectronPhononParameters{T<:Number, E<:AbstractFloat}

Describes all parameters in the electron-phonon model.

Fields

• β::E: Inverse temperature.
• Δτ::E: Discretization in imaginary time.
• Lτ::Int: Length of imaginary time axis.
• x::Matrix{E}: Phonon fields, where each column represents the phonon fields for a given imaginary time slice.
• phonon_parameters::PhononParameters{E,D}: Refer to PhononParameters.
• holstein_parameters_up::HolsteinParameters{E}: Spin up HolsteinParameters.
• holstein_parameters_dn::HolsteinParameters{E}: Spin down HolsteinParameters.
• ssh_parameters_up::SSHParameters{T}: Spin up SSHParameters.
• ssh_parameters_dn::SSHParameters{T}: Spin down SSHParameters.
• dispersion_parameters::DispersionParameters{E}: Refer to DispersionParameters.

ElectronPhononParameters(;
    β::E, Δτ::E,
    model_geometry::ModelGeometry{D,E},
    tight_binding_parameters::Union{TightBindingParameters{T,E}, Nothing} = nothing,
    tight_binding_parameters_up::Union{TightBindingParameters{T,E}, Nothing} = nothing,
    tight_binding_parameters_dn::Union{TightBindingParameters{T,E}, Nothing} = nothing,
    electron_phonon_model::ElectronPhononModel{T,E,D},
    rng::AbstractRNG
) where {T,E,D}

Initialize and return an instance of ElectronPhononParameters.

PhononParameters{E<:AbstractFloat}

Defines the parameters for each phonon in the lattice, includes the phonon field configuration.

Fields

• nphonon::Int: Number of type of phonon modes.
• Nphonon::Int: Total number of phonon modes in finite lattice.
• M::Int: Mass of each phonon mode.
• Ω::Int: Frequency of each phonon mode.
• Ω4::Int: Quartic phonon coefficient for each phonon mode.
• basis_vecs::Vector{Vector{E}}: Basis vector for each of the nphonon types of phonon mode.`

HolsteinParameters{E<:AbstractFloat}

Defines the Holstein coupling parameters in lattice.

Fields

• nholstein::Int: The number of type of holstein couplings.
• Nholstein::Int: Total number of Holstein couplings in lattice.
• α::Vector{T}: Linear Holstein coupling.
• α2::Vector{T}: Quadratic Holstein coupling.
• α3::Vector{T}: Cubic Holstein coupling.
• α4::Vector{T}: Quartic Holstein coupling.
• ph_sym_form::Vector{Bool}: If particle-hole symmetric form is used for Holstein coupling.
• coupling_to_site::Vector{Int}: Maps each Holstein coupling in the lattice to the corresponding site in the lattice.
• coupling_to_phonon::Vector{Int}: Maps each Holstein coupling in the lattice to the corresponding phonon mode.
• phonon_to_coupling::Vector{Vector{Int}}: Maps each phonon model to correspond Holstein couplings.

SSHParameters{T<:Number}

Defines the SSH coupling parameters in lattice.

Fields

• nssh::Int: Number of types of SSH couplings.
• Nssh::Int: Number of SSH couplings in lattice.
• α::Vector{T}: Linear SSH coupling.
• α2::Vector{T}: Quadratic SSH coupling.
• α3::Vector{T}: Cubic SSH coupling.
• α4::Vector{T}: Quartic SSH coupling.`
• neighbor_table::Matrix{Int}: Neighbor table to SSH coupling.
• coupling_to_phonon::Matrix{Int}: Maps each SSH coupling onto that pair of coupled phonons.
• init_phonon_to_coupling::Vector{Vector{Int}}: Maps initial phonon mode to corresponding SSH coupling(s).
• final_phonon_to_coupling::Vector{Vector{Int}}: Maps final phonon mode to corresponding SSH coupling(s).
• hopping_to_couplings::Vector{Vector{Int}}: Maps hopping in the tight-binding model onto SSH couplings.
• coupling_to_hopping::Vector{Int}: Maps each SSH coupling onto the corresponding hopping in the tight-binding model.

DispersionParameters{E<:AbstractFloat}

Defines the dispersive phonon coupling parameters in the lattice.

Fields

• ndispersion::Int: Number of types of dispersive couplings.
• Ndispersion::Int: Number of dispersive couplings in the lattice.
• Ω::Vector{E}: Frequency of dispersive phonon coupling.
• Ω4::Vector{E}: Quartic coefficient for the phonon dispersion.
• dispersion_to_phonon::Matrix{Int}: Pair of phonon modes in lattice coupled by dispersive coupling.
• init_phonon_to_coupling::Vector{Vector{Int}}: Maps initial phonon mode to corresponding dispersive phonon coupling.
• final_phonon_to_coupling::Vector{Vector{Int}}: Maps final phonon mode to corresponding dispersive phonon coupling.

initialize!(
    fermion_path_integral_up::FermionPathIntegral{H,T},
    fermion_path_integral_dn::FermionPathIntegral{H,T},
    electron_phonon_parameters::ElectronPhononParameters{T,R}
) where {H<:Number, T<:Number, R<:AbstractFloat}

Initialize the contribution of an ElectronPhononParameters to a FermionPathIntegral.

initialize!(
    # ARGUMENTS
    fermion_path_integral::FermionPathIntegral{H,T},
    electron_phonon_parameters::ElectronPhononParameters{T,R};
    # KEYWORD ARGUMENTS
    spin::Int = +1
) where {H<:Number, T<:Number, R<:AbstractFloat}

Initialize the contribution of an ElectronPhononParameters to a FermionPathIntegral.

update!(
    fermion_path_integral_up::FermionPathIntegral{H,T},
    fermion_path_integral_dn::FermionPathIntegral{H,T},
    electron_phonon_parameters::ElectronPhononParameters{T,R},
    x′::Matrix{R},
    x::Matrix{R}
) where {H<:Number, T<:Number, R<:AbstractFloat}

Update a FermionPathIntegral to reflect a change in the phonon configuration from x to x′.

update!(
    # ARGUMENTS
    fermion_path_integral::FermionPathIntegral{H,T},
    electron_phonon_parameters::ElectronPhononParameters{T,R},
    x′::Matrix{R},
    x::Matrix{R};
    # KEYWORD ARGUMENTS
    spin::Int = +1
) where {H<:Number, T<:Number, R<:AbstractFloat}

Update a FermionPathIntegral to reflect a change in the phonon configuration from x to x′.

update!(
    # ARGUMENTS
    fermion_path_integral::FermionPathIntegral{H,T},
    electron_phonon_parameters::ElectronPhononParameters{T,R},
    x::Matrix{R},
    sgn::Int;
    # KEYWORD ARGUMENTS
    spin::Int = +1
) where {H<:Number, T<:Number, R<:AbstractFloat}

Update a FermionPathIntegral according to sgn * x.

measure_phonon_kinetic_energy(
    electron_phonon_parameters::ElectronPhononParameters{T,E},
    n::Int
) where {T<:Number, E<:AbstractFloat}

Evaluate the average phonon kinetic energy for phonon mode n. The measurement is made using the expression

\[\langle K \rangle = \frac{1}{2\Delta\tau} - \frac{M}{2}\bigg\langle\frac{(x_{l+1}-x_{l})^{2}}{\Delta\tau^{2}}\bigg\rangle. \]

measure_phonon_potential_energy(
    electron_phonon_parameters::ElectronPhononParameters{T,E},
    n::Int
) where {T<:Number, E<:AbstractFloat}

Calculate the average phonon potential energy, given by

\[U = \frac{1}{2} M \Omega^2 \langle \hat{X}^2 \rangle + \frac{1}{24} M \Omega_4^2 \langle \hat{X}^4 \rangle,\]

for phonon mode n in the unit cell.

measure_phonon_position_moment(
    electron_phonon_parameters::ElectronPhononParameters{T,E},
    n::Int, m::Int
) where {T<:Number, E<:AbstractFloat}

Measure $\langle X^m \rangle$ for phonon mode n in the unit cell.

measure_holstein_energy(
    holstein_parameters::HolsteinParameters{E},
    G::Matrix{T},
    x::Matrix{E},
    holstein_id::Int
) where {T<:Number, E<:AbstractFloat}

Calculate and return both the spin-resolved Holstein interaction energy

\[\epsilon_{{\rm hol},\sigma} = \left\langle [ \alpha \hat{X} + \alpha_2 \hat{X}^2 + \alpha_3 \hat{X}^3 + \alpha_4 \hat{X}^4 ] \left( \hat{n}_\sigma - \frac{1}{2} \right) \right\rangle,\]

corresponding to holstein_id.

measure_ssh_energy(
    ssh_parameters::SSHParameters{T},
    G::Matrix{T}, x::Matrix{E}, ssh_id::Int
) where {T<:Number, E<:AbstractFloat}

Calculate the return the SSH interaction energy

\[\epsilon_{\rm ssh} = \left\langle [\alpha \hat{X} + \alpha_2 \hat{X}^2 \alpha_3 \hat{X}^3 + \alpha_4 \hat{X}^4] (\hat{c}^\dagger_{\sigma,i} \hat{c}_{\sigma,j} + {\rm h.c.}) \right\rangle\]

for coupling definition specified by ssh_id.

measure_dispersion_energy(
    electron_phonon_parameters::ElectronPhononParameters{T,E},
    dispersion_id::Int
) where {T<:Number, E<:AbstractFloat}

Evaluate the average dispersion energy

\[\epsilon_{\rm disp} = \frac{1}{2} M_{\rm red} \Omega^2 \langle(\hat{X}_i - \hat{X}_j)^2\rangle + \frac{1}{24} M{\rm red} \Omega_4^2 \langle(\hat{X}_i - \hat{X}_j)^4\rangle,\]

where $M_{\rm red} = \frac{M_i M_j}{M_i + M_j}$ is the reduced mass, for the dispersive coupling definition specified by dispersion_id.

EFAHMCUpdater{T<:Number, E<:AbstractFloat, PFFT, PIFFT}

Defines an Exact Fourier Acceleration Hamiltonian/Hybrid Monte Carlo (EFA-HMC) update for the phonon degrees of freedom.

Fields

• Nt::Int: Number of time-steps in HMC trajectory.
• Δt::E: Average time-step size used in HMC update.
• δ::E: Time-step used in EFA-HMC update is jittered by an amount Δt = Δt * (1 + δ*(2*rand(rng)-1)).
• x::Matrix{E}: Records initial phonon configuration in position space.
• p::Matrix{E}: Conjugate momentum in HMC dynamics.
• dSdx::Matrix{E}: Stores the derivative of the action.
• Gup′::Matrix{T}: Intermediate spin-up Green's function matrix during HMC trajectory.
• Gdn′::Matrix{T}: Intermediate spin-down Green's function matrix during HMC trajectory.
• efa::ExactFourierAccelerator{E,PFFT,PIFFT}: Type to perform exact integration of equations of motion of quantum harmonic oscillator.

EFAHMCUpdater(;
    # KEYWORD ARGUMENTS
    electron_phonon_parameters::ElectronPhononParameters{T,E},
    G::Matrix{T},
    Nt::Int,
    Δt::E = π/(2*Nt),
    reg::E = 0.0,
    δ::E = 0.05
) where {T<:Number, E<:AbstractFloat}

Keyword Arguments

• electron_phonon_parameters::ElectronPhononParameters{T,E}: Defines electron-phonon model.
• G::Matrix{T}: Sample Green's function matrix.
• Nt::Int: Number of time-steps used in EFA-HMC update.
• Δt::E = π/(2*Nt): Average step size used for HMC update.
• reg::E = 0.0: Regularization used for mass in equations of motion.
• δ::E = 0.05: Amount of jitter added to time-step used in EFA-HMC update.

hmc_update!(
    # ARGUMENTS
    Gup::Matrix{H}, logdetGup::R, sgndetGup::H,
    Gdn::Matrix{H}, logdetGdn::R, sgndetGdn::H,
    electron_phonon_parameters::ElectronPhononParameters{T,R},
    hmc_updater::EFAHMCUpdater{T,R};
    # KEYWORD ARGUMENTS
    fermion_path_integral_up::FermionPathIntegral{H,T},
    fermion_path_integral_dn::FermionPathIntegral{H,T},
    fermion_greens_calculator_up::FermionGreensCalculator{H,R},
    fermion_greens_calculator_dn::FermionGreensCalculator{H,R},
    fermion_greens_calculator_up_alt::FermionGreensCalculator{H,R},
    fermion_greens_calculator_dn_alt::FermionGreensCalculator{H,R},
    Bup::Vector{P}, Bdn::Vector{P},
    δG::R, δθ::R, rng::AbstractRNG, 
    update_stabilization_frequency::Bool = false,
    δG_max::R = 1e-5,
    δG_reject::R = 1e-2,
    recenter!::Function = identity,
    Nt::Int = hmc_updater.Nt,
    Δt::R = hmc_updater.Δt,
    δ::R = hmc_updater.δ
) where {H<:Number, T<:Number, R<:Real, P<:AbstractPropagator{T}}

Perform EFA-HMC update to the phonon degrees of freedom. This method returns (accepted, logdetGup, sgndetGup, logdetGdn, sgndetGdn, δG, δθ), where accepted is a boolean field indicating whether the proposed HMC update was accepted or rejected.

Arguments

• Gup::Matrix{H}: Green's function matrix for spin up.
• logdetGup::R: Log determinant of Green's function matrix for spin up.
• sgndetGup::H: Sign of determinant of Green's function matrix for spin up.
• Gdn::Matrix{H}: Green's function matrix for spin down.
• logdetGdn::R: Log determinant of Green's function matrix for spin down.
• electron_phonon_parameters::ElectronPhononParameters{T,R}: Electron-phonon model parameters.
• hmc_updater::EFAHMCUpdater{T,R}: EFA-HMC updater.

Keyword Arguments

• fermion_path_integral_up::FermionPathIntegral{H}: Fermion path integral for spin up.
• fermion_path_integral_dn::FermionPathIntegral{H}: Fermion path integral for spin down.
• fermion_greens_calculator_up::FermionGreensCalculator{H,R}: Fermion greens calculator for spin up.
• fermion_greens_calculator_dn::FermionGreensCalculator{H,R}: Fermion greens calculator for spin down.
• fermion_greens_calculator_up_alt::FermionGreensCalculator{H,R}: Alternative fermion greens calculator for spin up.
• fermion_greens_calculator_dn_alt::FermionGreensCalculator{H,R}: Alternative fermion greens calculator for spin down.
• Bup::Vector{P}: Spin up propagators.
• Bdn::Vector{P}: Spin down propagators.
• δG::R: Numerical error in Green's function corrected by numerical stabilization.
• δθ::R: Numerical error in the phase of the determinant of the Green's function matrix corrected by numerical stabilization.
• rng::AbstractRNG: Random number generator.
• update_stabilization_frequency::Bool = false: Whether to update the stabilization frequency.
• δG_max::R = 1e-5: Maximum numerical error in Green's function corrected by numerical stabilization.
• δG_reject::R = 1e-2: Reject the update if the numerical error in Green's function corrected by numerical stabilization is greater than this value.
• Nt::Int = hmc_updater.Nt: Number of time-steps used in EFA-HMC update.
• Δt::R = hmc_updater.Δt: Average step size used for HMC update.
• δ::R = hmc_updater.δ: Amount of jitter added to time-step used in EFA-HMC update.

hmc_update!(
    # ARGUMENTS
    G::Matrix{H}, logdetG::R, sgndetG::H,
    electron_phonon_parameters::ElectronPhononParameters{T,R},
    hmc_updater::EFAHMCUpdater{T,R};
    # KEYWORD ARGUMENTS
    fermion_path_integral::FermionPathIntegral{H,T},
    fermion_greens_calculator::FermionGreensCalculator{H,R},
    fermion_greens_calculator_alt::FermionGreensCalculator{H,R},
    B::Vector{P},
    δG::R, δθ::R, rng::AbstractRNG,
    update_stabilization_frequency::Bool = false,
    δG_max::R = 1e-5,
    δG_reject::R = 1e-2,
    recenter!::Function = identity,
    Nt::Int = hmc_updater.Nt,
    Δt::R = hmc_updater.Δt,
    δ::R = hmc_updater.δ
) where {H<:Number, T<:Number, R<:Real, P<:AbstractPropagator{T}}

Perform EFA-HMC update to the phonon degrees of freedom. This method returns (accepted, logdetG, sgndetG, δG, δθ), where accepted is a boolean field indicating whether the proposed HMC update was accepted or rejected.

Arguments

• G::Matrix{H}: Green's function matrix for spin up.
• logdetG::R: Log determinant of Green's function matrix for spin up.
• sgndetG::H: Sign of determinant of Green's function matrix for spin up.
• electron_phonon_parameters::ElectronPhononParameters{T,R}: Electron-phonon model parameters.
• hmc_updater::EFAHMCUpdater{T,R}: EFA-HMC updater.

Keyword Arguments

• fermion_path_integral::FermionPathIntegral{H}: Fermion path integral.
• fermion_greens_calculator::FermionGreensCalculator{H,R}: Fermion greens calculator.
• fermion_greens_calculator_alt::FermionGreensCalculator{H,R}: Alternative fermion greens calculator.
• B::Vector{P}: Spin up propagators.
• δG::R: Numerical error in Green's function corrected by numerical stabilization.
• δθ::R: Numerical error in the phase of the determinant of the Green's function matrix corrected by numerical stabilization.
• rng::AbstractRNG: Random number generator.
• update_stabilization_frequency::Bool = false: Whether to update the stabilization frequency.
• δG_max::R = 1e-5: Maximum numerical error in Green's function corrected by numerical stabilization.
• δG_reject::R = 1e-2: Reject the update if the numerical error in Green's function corrected by numerical stabilization is greater than this value.
• Nt::Int = hmc_updater.Nt: Number of time-steps used in EFA-HMC update.
• Δt::R = hmc_updater.Δt: Average step size used for HMC update.
• δ::R = hmc_updater.δ: Amount of jitter added to time-step used in EFA-HMC update.

reflection_update!(
    # ARGUMENTS
    Gup::Matrix{H}, logdetGup::R, sgndetGup::H,
    Gdn::Matrix{H}, logdetGdn::R, sgndetGdn::H,
    electron_phonon_parameters::ElectronPhononParameters{T,R};
    # KEYWORD ARGUMENTS
    fermion_path_integral_up::FermionPathIntegral{T,U},
    fermion_path_integral_dn::FermionPathIntegral{T,U},
    fermion_greens_calculator_up::FermionGreensCalculator{H,R},
    fermion_greens_calculator_dn::FermionGreensCalculator{H,R},
    fermion_greens_calculator_up_alt::FermionGreensCalculator{H,R},
    fermion_greens_calculator_dn_alt::FermionGreensCalculator{H,R},
    Bup::Vector{P}, Bdn::Vector{P}, rng::AbstractRNG,
    phonon_types = nothing
) where {H<:Number, T<:Number, U<:Number, R<:Real, P<:AbstractPropagator{T,U}}

Randomly sample a phonon mode in the lattice, and propose an update that reflects all the phonon fields associated with that phonon mode $x \rightarrow -x.$ This function returns (accepted, logdetGup, sgndetGup, logdetGdn, sgndetGdn).

Arguments

• Gup::Matrix{H}: Spin-up equal-time Greens function matrix.
• logdetGup::R: Log of the determinant of the spin-up equal-time Greens function matrix.
• sgndetGup::H: Sign/phase of the determinant of the spin-up equal-time Greens function matrix.
• Gdn::Matrix{H}: Spin-down equal-time Greens function matrix.
• logdetGdn::R: Log of the determinant of the spin-down equal-time Greens function matrix.
• sgndetGdn::H: Sign/phase of the determinant of the spin-down equal-time Greens function matrix.
• electron_phonon_parameters::ElectronPhononParameters{T,R}: Electron-phonon parameters, including the current phonon configuration.

Keyword Arguments

• fermion_path_integral_up::FermionPathIntegral{T,U}: An instance of FermionPathIntegral type for spin-up electrons.
• fermion_path_integral_dn::FermionPathIntegral{T,U}: An instance of FermionPathIntegral type for spin-down electrons.
• fermion_greens_calculator_up::FermionGreensCalculator{H,R}: Contains matrix factorization information for current spin-up sector state.
• fermion_greens_calculator_dn::FermionGreensCalculator{H,R}: Contains matrix factorization information for current spin-down sector state.
• fermion_greens_calculator_up_alt::FermionGreensCalculator{H,R}: Used to calculate matrix factorizations for proposed spin-up sector state.
• fermion_greens_calculator_dn_alt::FermionGreensCalculator{H,R}: Used to calculate matrix factorizations for proposed spin-up sector state.
• Bup::Vector{P}: Spin-up propagators for each imaginary time slice.
• Bdn::Vector{P}: Spin-down propagators for each imaginary time slice.
• rng::AbstractRNG: Random number generator used in method instead of global random number generator, important for reproducibility.
• phonon_types = nothing: Collection of phonon types in the unit cell to randomly sample a phonon mode from. If nothing then all phonon modes in the unit cell are considered.

reflection_update!(
    # ARGUMENTS
    G::Matrix{H}, logdetG::R, sgndetG::H,
    electron_phonon_parameters::ElectronPhononParameters{T,R};
    # KEYWORD ARGUMENTS
    fermion_path_integral::FermionPathIntegral{H,T},
    fermion_greens_calculator::FermionGreensCalculator{H,R},
    fermion_greens_calculator_alt::FermionGreensCalculator{H,R},
    B::Vector{P}, rng::AbstractRNG,
    phonon_types = nothing
) where {H<:Number, T<:Number, R<:Real, P<:AbstractPropagator{T}}

Randomly sample a phonon mode in the lattice, and propose an update that reflects all the phonon fields associated with that phonon mode $x \rightarrow -x.$ This function returns (accepted, logdetG, sgndetG).

Arguments

• G::Matrix{H}: equal-time Greens function matrix.
• logdetG::R: Log of the determinant of the equal-time Greens function matrix.
• sgndetG::H: Sign/phase of the determinant of the equal-time Greens function matrix.
• electron_phonon_parameters::ElectronPhonParameters{T,R}: Electron-phonon parameters, including the current phonon configuration.

Keyword Arguments

• fermion_path_integral::FermionPathIntegral{H,T}: An instance of FermionPathIntegral type.
• fermion_greens_calculator::FermionGreensCalculator{H,R}: Contains matrix factorization information for current state.
• fermion_greens_calculator_alt::FermionGreensCalculator{H,R}: Used to calculate matrix factorizations for proposed state.
• B::Vector{P}: Propagators for each imaginary time slice.
• rng::AbstractRNG: Random number generator used in method instead of global random number generator, important for reproducibility.
• phonon_types = nothing: Collection of phonon types (specified my PHONON_ID) in the unit cell to randomly sample a phonon mode from. If nothing then all phonon modes in the unit cell are considered.

swap_update!(
    # ARGUMENTS
    Gup::Matrix{H}, logdetGup::R, sgndetGup::H,
    Gdn::Matrix{H}, logdetGdn::R, sgndetGdn::H,
    electron_phonon_parameters::ElectronPhononParameters{T,R};
    # KEYWORD ARGUMENTS
    fermion_path_integral_up::FermionPathIntegral{H,T},
    fermion_path_integral_dn::FermionPathIntegral{H,T},
    fermion_greens_calculator_up::FermionGreensCalculator{H,R},
    fermion_greens_calculator_dn::FermionGreensCalculator{H,R},
    fermion_greens_calculator_up_alt::FermionGreensCalculator{H,R},
    fermion_greens_calculator_dn_alt::FermionGreensCalculator{H,R},
    Bup::Vector{P}, Bdn::Vector{P}, rng::AbstractRNG,
    phonon_id_pairs = nothing
) where {H<:Number, T<:Number, R<:Real, P<:AbstractPropagator{T}}

Randomly sample a pairs of phonon modes and exchange the phonon fields associated with the pair of phonon modes. This function returns (accepted, logdetGup, sgndetGup, logdetGdn, sgndetGdn).

Arguments

• Gup::Matrix{H}: Spin-up equal-time Greens function matrix.
• logdetGup::R: Log of the determinant of the spin-up equal-time Greens function matrix.
• sgndetGup::H: Sign/phase of the determinant of the spin-up equal-time Greens function matrix.
• Gdn::Matrix{H}: Spin-down equal-time Greens function matrix.
• logdetGdn::R: Log of the determinant of the spin-down equal-time Greens function matrix.
• sgndetGdn::H: Sign/phase of the determinant of the spin-down equal-time Greens function matrix.
• electron_phonon_parameters::ElectronPhononParameters{T,R}: Electron-phonon parameters, including the current phonon configuration.

Keyword Arguments

• fermion_path_integral_up::FermionPathIntegral{H,T}: An instance of FermionPathIntegral type for spin-up electrons.
• fermion_path_integral_dn::FermionPathIntegral{H,T}: An instance of FermionPathIntegral type for spin-down electrons.
• fermion_greens_calculator_up::FermionGreensCalculator{H,R}: Contains matrix factorization information for current spin-up sector state.
• fermion_greens_calculator_dn::FermionGreensCalculator{H,R}: Contains matrix factorization information for current spin-down sector state.
• fermion_greens_calculator_up_alt::FermionGreensCalculator{H,R}: Used to calculate matrix factorizations for proposed spin-up sector state.
• fermion_greens_calculator_dn_alt::FermionGreensCalculator{H,R}: Used to calculate matrix factorizations for proposed spin-up sector state.
• Bup::Vector{P}: Spin-up propagators for each imaginary time slice.
• Bdn::Vector{P}: Spin-down propagators for each imaginary time slice.
• rng::AbstractRNG: Random number generator used in method instead of global random number generator, important for reproducibility.
• phonon_id_pairs = nothing: Collection of phonon type pairs (specified by pairs of PHONON_ID values) in the unit cell to randomly sample a phonon modes from. If nothing then all phonon mode pairs in the unit cell are considered.

swap_update!(
    # ARGUMENTS
    G::Matrix{H}, logdetG::R, sgndetG::H,
    electron_phonon_parameters::ElectronPhononParameters{T,R};
    # KEYWORD ARGUMENTS
    fermion_path_integral::FermionPathIntegral{H,T},
    fermion_greens_calculator::FermionGreensCalculator{H,R},
    fermion_greens_calculator_alt::FermionGreensCalculator{H,R},
    B::Vector{P}, rng::AbstractRNG,
    phonon_id_pairs = nothing
) where {H<:Number, T<:Number, R<:Real, P<:AbstractPropagator{T}}

Randomly sample a pairs of phonon modes and exchange the phonon fields associated with the pair of phonon modes. This function returns (accepted, logdetG, sgndetG).

Arguments

• G::Matrix{H}: equal-time Greens function matrix.
• logdetG::R: Log of the determinant of the equal-time Greens function matrix.
• sgndetG::H: Sign/phase of the determinant of the equal-time Greens function matrix.
• electron_phonon_parameters::ElectronPhononParameters{T,R}: Electron-phonon parameters, including the current phonon configuration.

Keyword Arguments

• fermion_path_integral::FermionPathIntegral{H,T}: An instance of FermionPathIntegral type.
• fermion_greens_calculator::FermionGreensCalculator{H,R}: Contains matrix factorization information for current state.
• fermion_greens_calculator_alt::FermionGreensCalculator{H,R}: Used to calculate matrix factorizations for proposed state.
• B::Vector{P}: Propagators for each imaginary time slice.
• rng::AbstractRNG: Random number generator used in method instead of global random number generator, important for reproducibility.
• phonon_id_pairs = nothing: Collection of phonon type pairs in the unit cell to randomly sample a phonon modes from. If nothing then all phonon mode pairs in the unit cell are considered.

radial_update!(
    # ARGUMENTS
    Gup::Matrix{H}, logdetGup::R, sgndetGup::H,
    Gdn::Matrix{H}, logdetGdn::R, sgndetGdn::H,
    electron_phonon_parameters::ElectronPhononParameters{T,R};
    # KEYWORD ARGUMENTS
    fermion_path_integral_up::FermionPathIntegral{H,T},
    fermion_path_integral_dn::FermionPathIntegral{H,T},
    fermion_greens_calculator_up::FermionGreensCalculator{H,R},
    fermion_greens_calculator_dn::FermionGreensCalculator{H,R},
    fermion_greens_calculator_up_alt::FermionGreensCalculator{H,R},
    fermion_greens_calculator_dn_alt::FermionGreensCalculator{H,R},
    Bup::Vector{P}, Bdn::Vector{P}, rng::AbstractRNG,
    phonon_id::Union{Int, Nothing} = nothing,
    σ::R = 1.0
) where {H<:Number, T<:Number, R<:Real, P<:AbstractPropagator{T}}

Perform a radial update to the phonon fields, as described by Algorithm 1 in the paper arXiv:2411.18218. Specifically, the proposed update to the phonon fields $x$ is a rescaling such that $x \rightarrow e^{\gamma} x$ where $\gamma \sim N(0, \sigma/\sqrt{d})$ and $d$ is the number of phonon fields being updated.

Arguments

• Gup::Matrix{H}: Spin-up equal-time Greens function matrix.
• logdetGup::R: Log of the determinant of the spin-up equal-time Greens function matrix.
• sgndetGup::H: Sign/phase of the determinant of the spin-up equal-time Greens function matrix.
• Gdn::Matrix{H}: Spin-down equal-time Greens function matrix.
• logdetGdn::R: Log of the determinant of the spin-down equal-time Greens function matrix.
• sgndetGdn::H: Sign/phase of the determinant of the spin-down equal-time Greens function matrix.
• electron_phonon_parameters::ElectronPhononParameters{T,R}: Electron-phonon parameters, including the current phonon configuration.

Keyword Arguments

• fermion_path_integral_up::FermionPathIntegral{H,T}: An instance of FermionPathIntegral type for spin-up electrons.
• fermion_path_integral_dn::FermionPathIntegral{H,T}: An instance of FermionPathIntegral type for spin-down electrons.
• fermion_greens_calculator_up::FermionGreensCalculator{H,R}: Contains matrix factorization information for current spin-up sector state.
• fermion_greens_calculator_dn::FermionGreensCalculator{H,R}: Contains matrix factorization information for current spin-down sector state.
• fermion_greens_calculator_up_alt::FermionGreensCalculator{H,R}: Used to calculate matrix factorizations for proposed spin-up sector state.
• fermion_greens_calculator_dn_alt::FermionGreensCalculator{H,R}: Used to calculate matrix factorizations for proposed spin-up sector state.
• Bup::Vector{P}: Spin-up propagators for each imaginary time slice.
• Bdn::Vector{P}: Spin-down propagators for each imaginary time slice.
• rng::AbstractRNG: Random number generator used in method instead of global random number generator, important for reproducibility.
• phonon_id::Union{Int, Nothing} = nothing: Apply radial update to phonon fields corresponding tp specified PHONON_ID. If phonon_id = nothing, then radial update is applied to all phonon fields.
• σ::R = 1.0: Relative size of the radial update.

radial_update!(
    # ARGUMENTS
    G::Matrix{H}, logdetG::R, sgndetG::H,
    electron_phonon_parameters::ElectronPhononParameters{T,R};
    # KEYWORD ARGUMENTS
    fermion_path_integral::FermionPathIntegral{H,T},
    fermion_greens_calculator::FermionGreensCalculator{H,R},
    fermion_greens_calculator_alt::FermionGreensCalculator{H,R},
    B::Vector{P}, rng::AbstractRNG,
    phonon_id::Union{Int, Nothing} = nothing,
    σ::R = 1.0
) where {H<:Number, T<:Number, R<:Real, P<:AbstractPropagator{T}}

Perform a radial update to the phonon fields, as described by Algorithm 1 in the paper arXiv:2411.18218. Specifically, the proposed update to the phonon fields $x$ is a rescaling such that $x \rightarrow e^{\gamma} x$ where $\gamma \sim N(0, \sigma/\sqrt{d})$ and $d$ is the number of phonon fields being updated.

Arguments

• G::Matrix{H}: equal-time Greens function matrix.
• logdetG::R: Log of the determinant of the equal-time Greens function matrix.
• sgndetG::H: Sign/phase of the determinant of the equal-time Greens function matrix.
• electron_phonon_parameters::ElectronPhononParameters{T,R}: Electron-phonon parameters, including the current phonon configuration.

Keyword Arguments

• fermion_path_integral::FermionPathIntegral{H,T}: An instance of FermionPathIntegral type.
• fermion_greens_calculator::FermionGreensCalculator{H,R}: Contains matrix factorization information for current state.
• fermion_greens_calculator_alt::FermionGreensCalculator{H,R}: Used to calculate matrix factorizations for proposed state.
• B::Vector{P}: Propagators for each imaginary time slice.
• rng::AbstractRNG: Random number generator used in method instead of global random number generator, important for reproducibility.
• phonon_id::Union{Int, Nothing} = nothing: Apply radial update to phonon fields corresponding tp specified PHONON_ID. If phonon_id = nothing, then radial update is applied to all phonon fields.
• σ::R = 1.0: Relative size of the radial update.

update_chemical_potential!(
    # ARGUMENTS
    Gup::Matrix{H}, logdetGup::R, sgndetGup::H,
    Gdn::Matrix{H}, logdetGdn::R, sgndetGdn::H;
    # KEYWORD ARGUMENTS
    chemical_potential_tuner::MuTunerLogger{R,H},
    tight_binding_parameters::Union{TightBindingParameters, Nothing} = nothing,
    tight_binding_parameters_up::Union{TightBindingParameters, Nothing} = nothing,
    tight_binding_parameters_dn::Union{TightBindingParameters, Nothing} = nothing,
    fermion_path_integral_up::FermionPathIntegral{H},
    fermion_path_integral_dn::FermionPathIntegral{H},
    fermion_greens_calculator_up::FermionGreensCalculator{H},
    fermion_greens_calculator_dn::FermionGreensCalculator{H},
    Bup::Vector{P}, Bdn::Vector{P}
) where {H<:Number, R<:AbstractFloat, P<:AbstractPropagator}

Update the chemical potential $\mu$ in the simulation to approach the target density/filling. This method returns the new values for (logdetGup, sgndetGup, logdetGup, sgndetGup). Note that either the keyword tight_binding_parameters needs to be specified, or tight_binding_parameters_up and tight_binding_parameters_dn both need to be specified.

update_chemical_potential!(
    # ARGUMENTS
    G::Matrix{H}, logdetG::R, sgndetG::H;
    # KEYWORD ARGUMENTS
    chemical_potential_tuner::MuTunerLogger{R,H},
    tight_binding_parameters::TightBindingParameters,
    fermion_path_integral::FermionPathIntegral{H},
    fermion_greens_calculator::FermionGreensCalculator{H},
    B::Vector{P}
) where {H<:Number, R<:AbstractFloat, P<:AbstractPropagator}

Update the chemical potential $\mu$ in the simulation to approach the target density/filling. This method returns the new values for (logdetG, sgndetG).

save_density_tuning_profile(
    # ARGUMENTS
    simulation_info::SimulationInfo,
    chemical_potential_tuner::MuTunerLogger{R, H};
    # KEYWORD ARGUMENTS
    export_to_h5::Bool = true,
    export_to_csv::Bool = false,
    scientific_notation::Bool = false,
    decimals::Int = 9,
    delimiter::String = " ",
) where {R<:AbstractFloat, H<:Number}

Record the history of chemical potential and density tuning that occurred during the simulation, writing the information to an HDF5 and/or CSV file.

const GLOBAL_MEASUREMENTS = (
    "logdetGup",
    "logdetGdn",
    "sgndetGup",
    "sgndetGdn",
    "sgn",
    "action_total",
    "action_bosonic",
    "action_fermionic",
    "density",
    "density_up",
    "density_dn",
    "double_occ",
    "Nsqrd"
)

List of all the global measurements that are made.

const LOCAL_MEASUREMENTS = Base.ImmutableDict(
    "density"                  => "ORBITAL_ID",
    "density_up"               => "ORBITAL_ID",
    "density_dn"               => "ORBITAL_ID",
    "double_occ"               => "ORBITAL_ID",
    "onsite_energy"            => "ORBITAL_ID",
    "onsite_energy_up"         => "ORBITAL_ID",
    "onsite_energy_dn"         => "ORBITAL_ID",
    "bare_hopping_energy"      => "HOPPING_ID",
    "bare_hopping_energy_up"   => "HOPPING_ID",
    "bare_hopping_energy_dn"   => "HOPPING_ID",
    "hopping_energy"           => "HOPPING_ID",
    "hopping_energy_up"        => "HOPPING_ID",
    "hopping_energy_dn"        => "HOPPING_ID",
    "hopping_amplitude_up"     => "HOPPING_ID",
    "hopping_amplitude_dn"     => "HOPPING_ID",
    "hopping_amplitude"        => "HOPPING_ID",
    "hopping_inversion_up"     => "HOPPING_ID",
    "hopping_inversion_dn"     => "HOPPING_ID",
    "hopping_inversion"        => "HOPPING_ID",
    "hubbard_energy"           => "HUBBARD_ID",
    "ext_hub_energy"           => "EXT_HUB_ID",
    "phonon_kin_energy"        => "PHONON_ID",
    "phonon_pot_energy"        => "PHONON_ID",
    "X"                        => "PHONON_ID",
    "X2"                       => "PHONON_ID",
    "X3"                       => "PHONON_ID",
    "X4"                       => "PHONON_ID",
    "holstein_energy"          => "HOLSTEIN_ID",
    "holstein_energy_up"       => "HOLSTEIN_ID",
    "holstein_energy_dn"       => "HOLSTEIN_ID",
    "ssh_energy"               => "SSH_ID",
    "ssh_energy_up"            => "SSH_ID",
    "ssh_energy_dn"            => "SSH_ID",
    "dispersion_energy"        => "DISPERSION_ID"
)

List of all the local measurements than can be made, with a mapping to the corresponding type of ID each measurement is reported in terms of.

CORRELATION_FUNCTIONS = Base.ImmutableDict(
    "greens"           => "ORBITAL_ID",
    "greens_up"        => "ORBITAL_ID",
    "greens_dn"        => "ORBITAL_ID",
    "greens_tautau"    => "ORBITAL_ID",
    "greens_tautau_up" => "ORBITAL_ID",
    "greens_tautau_dn" => "ORBITAL_ID",
    "density"          => "ORBITAL_ID",
    "density_upup"     => "ORBITAL_ID",
    "density_dndn"     => "ORBITAL_ID",
    "density_updn"     => "ORBITAL_ID",
    "density_dnup"     => "ORBITAL_ID",
    "spin_x"           => "ORBITAL_ID",
    "spin_z"           => "ORBITAL_ID",
    "pair"             => "BOND_ID",
    "bond"             => "BOND_ID",
    "bond_upup"        => "BOND_ID",
    "bond_dndn"        => "BOND_ID",
    "bond_updn"        => "BOND_ID",
    "bond_dnup"        => "BOND_ID",
    "current"          => "HOPPING_ID",
    "current_upup"     => "HOPPING_ID",
    "current_dndn"     => "HOPPING_ID",
    "current_updn"     => "HOPPING_ID",
    "current_dnup"     => "HOPPING_ID",
    "phonon_greens"    => "PHONON_ID"
)

List of all the correlation functions that can be measured, along with the corresponding type of ID the correlation measurement is reported in terms of. Correlation functions are well defined in both position and momentum space.

initialize_measurement_container(
    model_geometry::ModelGeometry{D,T,N},
    β::T, Δτ::T
) where {T<:AbstractFloat, D, N}

Initialize and return a measurement container of type NamedTuple.

initialize_measurements!(
    measurement_container::NamedTuple,
    tight_binding_model_up::TightBindingModel{T,E},
    tight_binding_model_dn::TightBindingModel{T,E},
) where {T<:Number, E<:AbstractFloat}

initialize_measurements!(
    measurement_container::NamedTuple,
    tight_binding_model::TightBindingModel{T,E}
) where {T<:Number, E<:AbstractFloat}

Initialize tight-binding model related measurements.

Initialized Measurements

• onsite_energy: Refer to measure_onsite_energy.
• onsite_energy_up: Refer to measure_onsite_energy.
• onsite_energy_dn: Refer to measure_onsite_energy.
• hopping_energy: Refer to measure_hopping_energy.
• hopping_energy_up: Refer to measure_hopping_energy.
• hopping_energy_dn: Refer to measure_hopping_energy.

initialize_measurements!(
    measurement_container::NamedTuple,
    hubbard_model::HubbardModel{T}
) where {T<:AbstractFloat}

Initialize Hubbard model related measurements.

Initialized Measurements:

• hubbard_energy: Refer to measure_hopping_energy.

initialize_measurements!(
    measurement_container::NamedTuple,
    extended_hubbard_model::ExtendedHubbardModel{T}
) where {T<:AbstractFloat}

Initialize Extended Hubbard model related measurements.

Initialized Measurements:

• ext_hubbard_energy: Refer to measure_ext_hub_energy.

initialize_measurements!(
    measurement_container::NamedTuple,
    electron_phonon_model::ElectronPhononModel{T, E, D}
) where {T<:Number, E<:AbstractFloat, D}

Initialize electron-phonon model related measurements.

Initialized Measurements:

• phonon_kinetic_energy: Refer to measure_phonon_kinetic_energy.
• phonon_potential_energy: Refer to measure_phonon_potential_energy.
• X: Measure $\langle \hat{X} \rangle$, refer to measure_phonon_position_moment.
• X2: Measure $\langle \hat{X}^2 \rangle$, refer to measure_phonon_position_moment.
• X3: Measure $\langle \hat{X}^3 \rangle$, refer to measure_phonon_position_moment.
• X4: Measure $\langle \hat{X}^4 \rangle$, refer to measure_phonon_position_moment.
• holstein_energy: Refer to measure_holstein_energy.
• holstein_energy_up: Refer to measure_holstein_energy.
• holstein_energy_dn: Refer to measure_holstein_energy.
• ssh_energy: Refer to measure_ssh_energy.
• ssh_energy_up: Refer to measure_ssh_energy.
• ssh_energy_dn: Refer to measure_ssh_energy.
• dispersion_energy: Refer to measure_dispersion_energy.

initialize_correlation_measurements!(;
    measurement_container::NamedTuple,
    model_geometry::ModelGeometry{D,T,N},
    correlation::String,
    pairs::AbstractVector{NTuple{2,Int}},
    time_displaced::Bool,
    integrated::Bool = false
)  where {T<:AbstractFloat, D, N}

Initialize measurements of correlation for all ID pairs; refer to CORRELATION_FUNCTIONS for ID type associated with each correlation measurement. The name correlation must therefore also appear in [CORRELATION_FUNCTIONS]@ref. If time_displaced = true then time-displaced and integrated correlation measurements are made. If time_displaced = false and integrated = false, then just equal-time correlation measurements are made. If time_displaced = false and integrated = true, then both equal-time and integrated correlation measurements are made.

initialize_composite_correlation_measurement!(;
    measurement_container::NamedTuple,
    model_geometry::ModelGeometry{D,T,N},
    name::String,
    correlation::String,
    ids::Union{Nothing,Vector{Int}} = nothing,
    id_pairs::Union{Nothing,Vector{NTuple{2,Int}}} = nothing,
    coefficients,
    displacement_vecs = nothing,
    time_displaced::Bool,
    integrated::Bool = false
)  where {T<:AbstractFloat, D, N}

Initialize a composite correlation measurement called name based on a linear combination of local operators used in a standard correlation measurement.

If the keyword ids is passed and id_pairs = nothing, then the composite correlation function is given by

\[\begin{align*} C_{\mathbf{r}}(\tau) & = \frac{1}{N}\sum_{\mathbf{i}}\langle\hat{\Phi}_{\mathbf{i}+\mathbf{r}}^{\dagger}(\tau)\hat{\Phi}_{\mathbf{i}}^{\phantom{\dagger}}(0)\rangle \\ & = \frac{1}{N}\sum_{\eta,\nu}\sum_{\mathbf{i}}c_{\eta}^{*}c_{\nu}\langle\hat{O}_{\mathbf{i}+\mathbf{r},\eta}^{\dagger}(\tau)\hat{O}_{\mathbf{i},\nu}^{\phantom{\dagger}}(0)\rangle \\ & = \sum_{\eta,\nu}c_{\eta}^{*}c_{\nu}C_{\mathbf{r}}^{\eta,\nu}(\tau) \end{align*}\]

where the composite operator is

\[\hat{\Phi}_{\mathbf{\mathbf{r}}}(\tau)=\sum_{\nu}c_{\nu}\hat{O}_{\mathbf{r},\nu}(\tau).\]

The sum over $\mathbf{i}$ runs over all unit cells, $\mathbf{r}$ denotes a displacement in unit cells and $N$ is the number unit cells. The operator type $\hat{O}^{\nu}$ and corresponding correlation function type $C_{\mathbf{r}}^{\eta,\nu}(\tau)$ are specified by the correlation keyword, while $\nu$ corresponds to labels/IDs specified by the ids keyword argument. Lastly, the $c_\nu$ coefficients are specified using the coefficients keyword arguments. The corresponding fourier transform of this composite correlation function measurement is given by

\[S_{\mathbf{q}}(\tau)=\sum_{\eta,\nu}\sum_{\mathbf{r}}e^{-{\rm i}\mathbf{q}\cdot(\mathbf{r}+\mathbf{r}_{\eta}-\mathbf{r}_{\nu})}C_{\mathbf{r}}^{\eta,\nu}(\tau),\]

where the static vectors $\mathbf{r}_\nu$ are specified using the displacement_vecs keyword arguments. If displacement_vecs = nothing then $\mathbf{r}_\nu = 0$ for all label/ID values $\nu$.

On the other hand, if id_pairs is passed and ids = nothing, then the composite correlation function is given by

\[\begin{align*} C_{\mathbf{r}}(\tau) & = \sum_{n}c_{n}C_{\mathbf{r}}^{\eta_{n},\nu_{n}}(\tau) \\ & = \frac{1}{N}\sum_{n}\sum_{\mathbf{i}}c_{n}\langle\hat{O}_{\mathbf{i}+\mathbf{r},\eta_{n}}^{\dagger}(\tau)\hat{O}_{\mathbf{i},\nu_{n}}^{\phantom{\dagger}}(0)\rangle, \end{align*}\]

where the $n$ index runs over pairs of labels/IDs $(\nu_n, \eta_n)$ specified by the id_pairs keyword argument. Note that the order of the label/ID pair $(\nu_n, \eta_n)$ reflects how each tuple in the id_pairs vector will be interpreted. Once again, operator type $\hat{O}^{\nu_n}$ and corresponding correlation function type $C_{\mathbf{r}}^{\eta_n,\nu_n}(\tau)$ are specified by the correlation keyword. The corresponding fourier transform of this composite correlation function measurement is given by

\[S_{\mathbf{q}}(\tau)=\sum_{n}\sum_{\mathbf{r}}e^{-{\rm i}\mathbf{q}\cdot(\mathbf{r}+\mathbf{r}_{n})}C_{\mathbf{r}}^{\eta_{n},\nu_{n}}(\tau),\]

where the static displacement vectors $\mathbf{r}_n$ are specified by the displacement_vecs keyword argument. As before, if displacement_vecs = nothing, then $\mathbf{r}_n = 0$ for all $n$.

Note that the specified correlation type correlation needs to correspond to one of the keys in the global CORRELATION_FUNCTIONS dictionary, which lists all the predefined types of correlation functions that can be measured.

make_measurements!(
    # ARGUMENTS
    measurement_container::NamedTuple,
    logdetGup::E, sgndetGup::T, Gup::AbstractMatrix{T},
    Gup_ττ::AbstractMatrix{T}, Gup_τ0::AbstractMatrix{T}, Gup_0τ::AbstractMatrix{T},
    logdetGdn::E, sgndetGdn::T, Gdn::AbstractMatrix{T},
    Gdn_ττ::AbstractMatrix{T}, Gdn_τ0::AbstractMatrix{T}, Gdn_0τ::AbstractMatrix{T};
    # KEYWORD ARGUMENTS
    fermion_path_integral_up::FermionPathIntegral{T,E},
    fermion_path_integral_dn::FermionPathIntegral{T,E},
    fermion_greens_calculator_up::FermionGreensCalculator{T,E},
    fermion_greens_calculator_dn::FermionGreensCalculator{T,E},
    Bup::Vector{P}, Bdn::Vector{P},
    model_geometry::ModelGeometry{D,E,N},
    tight_binding_parameters::Union{Nothing, TightBindingParameters} = nothing,
    tight_binding_parameters_up::Union{Nothing, TightBindingParameters} = nothing,
    tight_binding_parameters_dn::Union{Nothing, TightBindingParameters} = nothing,
    coupling_parameters::Tuple,
    δG::E, δθ::E, δG_max::E = 1e-6
) where {T<:Number, E<:AbstractFloat, D, N, P<:AbstractPropagator}

Make measurements, including time-displaced correlation and zero Matsubara frequency measurements. This method also returns (logdetGup, sgndetGup, logdetGdn, sgndetGdn, δG, δθ). Note that either the keyword tight_binding_parameters needs to be specified, or tight_binding_parameters_up and tight_binding_parameters_dn both need to be specified.

make_measurements!(
    # ARGUMENTS
    measurement_container::NamedTuple,
    logdetG::E, sgndetG::T, G::AbstractMatrix{T},
    G_ττ::AbstractMatrix{T}, G_τ0::AbstractMatrix{T}, G_0τ::AbstractMatrix{T};
    # KEYWORD ARGUMENTS
    fermion_path_integral::FermionPathIntegral{T},
    fermion_greens_calculator::FermionGreensCalculator{T,E},
    B::Vector{P},
    model_geometry::ModelGeometry{D,E,N},
    tight_binding_parameters::TightBindingParameters,
    coupling_parameters::Tuple,
    δG::E, δθ::E, δG_max::E = 1e-6
) where {T<:Number, E<:AbstractFloat, D, N, P<:AbstractPropagator}

Make measurements, including time-displaced correlation and zero Matsubara frequency measurements. This method also returns (logdetG, sgndetG, δG, δθ).

write_measurements!(;
    measurement_container::NamedTuple,
    simulation_info::SimulationInfo,
    model_geometry::ModelGeometry{D, E, N},
    Δτ::E,
    bin_size::Int,
    measurement::Int = 0,
    bin::Int = measurement ÷ bin_size
) where {D, E<:AbstractFloat, N}

Write the measurements contained in measurement_container to file if update % bin_size == 0. Measurements are written to file in a binary format using the JLD2.jl package.

This function also does a few other things:

1. Normalizes all the measurements by the bin_size i.e. the number of measurements that were accumulated into the measurement container.
2. Take position space correlation function measurements and fourier transform them to momentum space.
3. Integrate relevant time-displaced correlation function measurements over imaginary time to get the corresponding zero Matsubara frequency correlation function.
4. Reset all the measurements in measurement_container to zero after the measurements are written to file.

merge_bins(
    # ARGUMENTS
    simulation_info::SimulationInfo
)

Merge the separate HDF5 files containing the binned measurements into a single HDF5 file. This is true even if the HDF5 "files" containing the binned data were held in memory during the simulation (simulation_info.write_bins_concurrent = false) instead of being actively written to file during the simulation (simulation_info.write_bins_concurrent = true).

rm_bins(
    comm::MPI.Comm,
    datafolder::String
)

rm_bins(
    datafolder::String
)

Delete the binned data stored in the datafolder directory. This function essentially deletes the directory datafolder/bins and its contents.

write_jld2_checkpoint(
    # Arguments
    comm::MPI.Comm,
    simulation_info::SimulationInfo;
    # Required Keyword Arguments
    model_geometry::ModelGeometry,
    measurement_container::NamedTuple,
    # Optional Keyword Arguments
    checkpoint_timestamp::T = 0.0,
    checkpoint_freq::T = 0.0,
    start_timestamp::T = 0.0,
    runtime_limit::T = Inf,
    error_code::Int = 13,
    # Arbitrary Keyword Arguments Written to Checkpoint
    kwargs...
) where {T<:AbstractFloat}

write_jld2_checkpoint(
    # Arguments
    simulation_info::SimulationInfo;
    # Required Keyword Arguments
    model_geometry::ModelGeometry,
    measurement_container::NamedTuple,
    # Optional Keyword Arguments
    checkpoint_timestamp::T = 0.0,
    checkpoint_freq::T = 0.0,
    start_timestamp::T = 0.0,
    runtime_limit::T = Inf,
    error_code::Int = 13,
    # Arbitrary Keyword Arguments Written to Checkpoint
    kwargs...
) where {T<:AbstractFloat}

Checkpoint a simulation by writing a new checkpoint file if necessary The checkpoint file is a JLD2 binary file.

Arguments

• comm::MPI.Comm: (optional) MPI communicator object used to synchronize processes. Ensures all MPI processes remain synchronized.
• simulation_info::SimulationInfo: Contains datafolder and MPI process ID information.

Keyword Arguments

• checkpoint_timestamp::T = 0.0: (optional) Epoch timestamp of previously written checkpoint file.
• checkpoint_freq::T = 0.0: (optional) Frequency with with checkpoint files are written; new checkpoint is written only if this many seconds has elapsed since previous checkpoint.
• start_timestamp::T = 0.0: (optional) Epoch timestamp of the start time of the simulation.
• runtime_limit::T = Inf: (optional) Maximum runtime for simulation in seconds; if after writing a new checkpoint file the next checkpoint file that would be written in the future exceeds the runtime limit then exit the simulation.
• error_code::Int = 13: (optional) Error code used to exit simulation if the runtime limit is exceeded.
• kwargs...: Additional keyword arguments containing the information that will stored in the checkpoint file; keyword arguments can point to arbitrary Julia objects.

Notes

The default values for the checkpoint_timestamp, checkpoint_freq, start_timestamp, and runtime_limit keyword arguments result in there being no runtime limit for the simulation and a new checkpoint file being written every time this function is called.

read_jld2_checkpoint(
    simulation_info::SimulationInfo
)

Read in a checkpoint file written using the write_jld2_checkpoint function and return its contents as a dictionary. This function returns the tuple (checkpoint, checkpoint_timestamp) where checkpoint is a dictionary containing the contents of the checkpoint file and checkpoint_timestamp is the epoch timestamp corresponding to when the checkpoint file was read in. Behind the scenes, the JLD2.jl package is used to read (and write) the checkpoint files.

rm_jld2_checkpoints(
    # ARGUMENTS
    comm::MPI.Comm,
    simulation_info::SimulationInfo
)

rm_jld2_checkpoints(
    # ARGUMENTS
    simulation_info::SimulationInfo
)

Delete the JLD2 checkpoint files.

rename_complete_simulation(
    # Arguments
    comm::MPI.Comm,
    simulation_info::SimulationInfo;
    # Keyword Arguments
    delete_jld2_checkpoints::Bool = true
)

rename_complete_simulation(
    # Arguments
    simulation_info::SimulationInfo;
    # Keyword Arguments
    delete_jld2_checkpoints::Bool = true
)

When a simulation is complete, this function renames the data folder the results were written to such that the directory name now begins with "complete_", making it simpler to identify which simulations no longer need to be resumed if checkpointing is being used. This function also deletes the any checkpoint files written using the write_jld2_checkpoint function if delete_jld2_checkpoints = true.

save_simulation_info(
    # ARGUMENTS
    sim_info::SimulationInfo,
    metadata = nothing;
    # KEYWORD ARGUMENTS
    filename = @sprintf "simulation_info_sID-%d_pID-%d.toml" sim_info.sID sim_info.pID
)

Save the contents sim_info to a TOML file, and add an optional additional table to the output file based on the contents of a dictionary metadata.

process_measurements(
    # ARGUMENTS
    comm::MPI.Comm;
    # KEYWORD ARGUMENTS
    datafolder::String,
    ...
)

process_measurements(;
    # KEYWORD ARGUMENTS
    datafolder::String,
    ...
)

Process the HDF5 files containing the binned data generated by a DQMC simulation, and then write final statistics to a new HDF5 file(s). There is also functionality available that allows the final stats to then be exported to CSV files as well.

When the first passed arguments is comm::MPI.Comm, then the workflow is parallelized using MPI.

Keyword Arguments

This function has many possible keyword arguments. Below I group them together into related categories and define their meanings.

Required Keyword Arguments

• datafolder::String: Specify the directory generated by a DQMC simulation into which all the results were written.

Keyword for Controlling Workflow

• n_bins::Union{Int, Nothing} = nothing: Number of bins used to calculate statistics. If nothing then set equal to the number of data bins written to file during the simulation. Must be a factor of the number of bins written to file during the DQMC simulation.
• pIDs::Union{Int,Vector{Int}} = Int[]: Specifies for which process IDs to calculate average statistics. If pIDs = Int[], the calculate for all process IDs. If comm::MPI.Comm is passed as first function argument then pIDs must be of type Vector{Int} and not Int.
• filename_prefix::String = "stats": Start of filename for HDF5 containing final statistics. HDF5 files containing statistics for a single process ID will end with pID-$(pID).h5.
• rm_binned_data::Bool = false: Whether to delete the binned data after final statistics are computed.

Keywords for Exporting Statistics to CSV

• export_to_csv::Bool = true: Whether to export the final statistics to CSV file.
• scientific_notation::Bool = false: Whether to use scientific notation when exporting statistics to CSV file.
• decimals::Int = 9: How many decimal places to include when exporting statistics to CSV files.
• delimiter::String = " ": Delimiter used when writing CSV files.

Keyword Acting as Boolean Flags Indicating Which Statistics to Compute

• process_global_measurements::Bool = true: Whether to calculate the statistics for global measurements.
• process_local_measurements::Bool = true: Whether to calculate the statistics for local measurements.
• process_all_equal_time_measurements::Bool = true: Whether to calculate statistics for all equal-time correlation measurements.
• process_all_time_displaced_measurements::Bool = false: Whether to calculate statistics for all time-displaced correlation measurements.
• process_all_integrated_measurements::Bool = true: Whether to calculate statistics for all integrated correlation measurements.

Keyword Specifying Specific Correlation Statistics to Compute

If process_all_equal_time_measurements = false, then the keyword arguments below can be used to specify which specific equal-time correlation measurements to calculate statistics for.

• standard_equal_time::Vector{String} = String[]
• composite_integrated::Vector{String} = String[]

If process_all_time_displaced_measurements = true, then the keyword arguments below can be used to specify which specific time-displaced correlation measurements to calculate statistics for.

• standard_time_displaced::Vector{String} = String[]
• composite_time_displaced::Vector{String} = String[]

If process_all_integrated_measurements = false, then the keyword arguments below can be used to specify which specific integrated correlation measurements to calculate statistics for.

• standard_integrated::Vector{String} = String[]
• composite_integrated::Vector{String} = String[]

compute_correlation_ratio(
    comm::MPI.Comm;
    # KEYWORD ARGUMENTS
    datafolder::String,
    correlation::String,
    type::String,
    id_pairs::Vector{NTuple{2,Int}},
    id_pair_coefficients::Vector{T},
    q_point::NTuple{D,Int},
    q_neighbors::Vector{NTuple{D,Int}},
    num_bins::Int = 0,
    pIDs::Vector{Int} = Int[],
) where {D, T<:Number}

compute_correlation_ratio(;
    # KEYWORD ARGUMENTS
    datafolder::String,
    correlation::String,
    type::String,
    id_pairs::Vector{NTuple{2,Int}},
    id_pair_coefficients::Vector{T},
    q_point::NTuple{D,Int},
    q_neighbors::Vector{NTuple{D,Int}},
    num_bins::Int = 0,
    pIDs::Union{Int,Vector{Int}} = Int[]
) where {D, T<:Number}

Compute the correlation ratio at the $\mathbf{k}$-point using a linear combination of standard correlation function measurements. The linear combination of correlation functions used is defined by id_pairs and coefs. If type is "equal-time" or "time-displaced" then the equal-time correlation ratio is calculated. If type is "integrated" then the integrated correlation ratio is calculated.

compute_composite_correlation_ratio(
    comm::MPI.Comm;
    # KEYWORD ARGUMENTS
    datafolder::String,
    name::String,
    type::String,
    q_point::NTuple{D,Int},
    q_neighbors::Vector{NTuple{D,Int}},
    pIDs::Vector{Int} = Int[]
) where {D}

compute_composite_correlation_ratio(;
    # Keyword Arguments
    datafolder::String,
    name::String,
    type::String,
    q_point::NTuple{D,Int},
    q_neighbors::Vector{NTuple{D,Int}},
    num_bins::Int = 0,
    pIDs::Union{Int,Vector{Int}} = Int[]
) where {D}

compute_function_of_correlations(
    # ARGUMENTS
    comm::MPI.Comm;
    # KEYWORD ARGUMENTS
    f::Function,
    datafolder::String,
    correlations::AbstractVector{<:NamedTuple},
    num_bins::Int = 0,
    pIDs::Vector{Int} = Int[],
)

compute_function_of_correlations(;
    # KEYWORD ARGUMENTS
    f::Function,
    datafolder::String,
    correlations::AbstractVector{<:NamedTuple},
    num_bins::Int = 0,
    pIDs::Union{Int,Vector{Int}} = Int[]
)

Calculate the mean and error associated with computing a function of measured correlation functions.

The correlation measurements that are passed as arguments to the function to evaluate are specified by the vector of named tuples correlations. The keys of the named tuple used to specify the correlation measurements are given below:

• name::String: Name of correlation function.
• type::String: Specifies whether to use "EQUAL-TIME", "TIME-DISPLACED" or "INTEGRATED" correlation measurement.
• id_pair::NTuple{2,Int}: If a standard correlation measurement, then specifies the ID pair for the correlation measurement.
• R::NTuple{D,Int} xor K::NTuple{D,Int}: Specifies either displacement vector or momentum point in D dimensions.
• τ::Float64 xor l::Int: Imaginary-time or imaginary-time slice used if the type is "TIME-DISPLACED".

Keyword Arguments

• f::Function: The function to evaluate.
• datafolder::String: Output directory generated by SmoQyDQMC simulation.
• correlations::AbstractVector{NamedTuple}: Vector of named tuples specifying the correlation measurement function arguments.
• num_bins::Int = 0: Number of bins used to compute statistics. If zero then use all number of bins.
• pIDs = Int[]: Process ID's used when measuring function of correlation measurements. If is Int[], then all process ID's are used.

export_global_stats_to_csv(
    # ARGUMENTS
    comm::MPI.Comm;
    # KEYWORD ARGUMENTS
    datafolder::String,
    h5filename::String = "stats.h5",,
    csvfilename_prefix::String = "global",
    measurements::Vector{String} = String[],
    decimals::Int = 6,
    scientific_notation::Bool = false,
    delimiter::String = " "
)

export_global_stats_to_csv(;
    # KEYWORD ARGUMENTS
    datafolder::String,
    h5filename::String = "stats.h5",,
    csvfilename_prefix::String = "global",
    measurements::Vector{String} = String[],
    decimals::Int = 6,
    scientific_notation::Bool = false,
    delimiter::String = " "
)

This function writes the global measurement statistics stored in the h5filename HDF5 file found in the directory datafolder to CSV file, returning the name of the CSV file that was written. The measurements keyword argument specifies the measurements to be exported. If measurements = String[], then all measurements are exported.

export_global_bins_to_h5(
    # ARGUMENTS
    comm::MPI.Comm;
    # KEYWORD ARGUMENTS
    datafolder::String,
    filename_prefix::String = "global_bins",
    pIDs::Union{Vector{Int},Int} = Int[],
    global_measurements::Vector{String} = String[] 
)

export_global_bins_to_h5(;
    # KEYWORD ARGUMENTS
    datafolder::String,
    filename_prefix::String = "global_bins",
    pIDs::Union{Vector{Int},Int} = Int[],
    global_measurements::Vector{String} = String[]
)

Export the binned global measurements for specified process IDs pIDs to a single HDF5 file. If pIDs = Int[], then binned global measurements for all process IDs are exported. You can specify a subset of specific global measurements using the global_measurements keyword argument. If global_measurements = String[], then all global measurements are exported.

export_global_bins_to_csv(
    # ARGUMENTS
    comm::MPI.Comm;
    # KEYWORD ARGUMENTS
    datafolder::String,
    filename_prefix::String = "global_bins",
    pIDs::Union{Vector{Int},Int} = Int[],
    global_measurements::Vector{String} = String[],
    decimals::Int = 9,
    scientific_notation::Bool = false,
    delimiter::String = " "
)

export_global_bins_to_csv(;
    # KEYWORD ARGUMENTS
    datafolder::String,
    filename_prefix::String = "global_bins",
    pIDs::Union{Vector{Int},Int} = Int[],
    global_measurements::Vector{String} = String[],
    decimals::Int = 9,
    scientific_notation::Bool = false,
    delimiter::String = " "
)

Export the binned global measurements for specified process IDs pIDs to a single CSV file. If pIDs = Int[], then binned global measurements for all process IDs are exported. You can specify a subset of specific global measurements using the global_measurements keyword argument. If global_measurements = String[], then all global measurements are exported.

export_local_stats_to_csv(
    # ARGUMENTS
    comm::MPI.Comm;
    # KEYWORD ARGUMENTS
    datafolder::String,
    h5filename::String = "stats.h5",
    csvfilename_prefix::String = "local",
    measurements::Vector{String} = String[],
    decimals::Int = 6,
    scientific_notation::Bool = false,
    delimiter::String = " "
)

export_local_stats_to_csv(;
    # KEYWORD ARGUMENTS
    datafolder::String,
    h5filename::String = "stats.h5",
    csvfilename_prefix::String = "local",
    measurements::Vector{String} = String[],
    decimals::Int = 6,
    scientific_notation::Bool = false,
    delimiter::String = " "
)

This function writes the local measurement statistics stored in the h5filename HDF5 file found in the directory datafolder to CSV file, returning the name of the CSV file that was written. The measurements keyword argument specifies the measurements to be exported. If measurements = String[], then all measurements are exported.

export_local_bins_to_csv(
    # ARGUMENTS
    comm::MPI.Comm;
    # KEYWORD ARGUMENTS
    datafolder::String,
    filename_prefix::String = "local_bins",
    pIDs::Union{Vector{Int},Int} = Int[],
    local_measurements::Vector{String} = String[],
    decimals::Int = 9,
    scientific_notation::Bool = false,
    delimiter::String = " "
)

export_local_bins_to_csv(;
    # KEYWORD ARGUMENTS
    datafolder::String,
    filename_prefix::String = "local_bins",
    pIDs::Union{Vector{Int},Int} = Int[],
    local_measurements::Vector{String} = String[],
    decimals::Int = 9,
    scientific_notation::Bool = false,
    delimiter::String = " "
)

Export the binned local measurements for specified process IDs pIDs to a single CSV file. If pIDs = Int[], then binned local measurements for all process IDs are exported. You can specify a subset of specific local measurements using the local_measurements keyword argument. If local_measurements = String[], then all local measurements are exported.

export_local_bins_to_h5(
    # ARGUMENTS
    comm::MPI.Comm;
    # KEYWORD ARGUMENTS
    datafolder::String,
    filename_prefix::String = "local_bins",
    pIDs::Union{Vector{Int},Int} = Int[],
    local_measurements::Vector{String} = String[],
)

export_local_bins_to_h5(;
    # KEYWORD ARGUMENTS
    datafolder::String,
    filename_prefix::String = "local_bins",
    pIDs::Union{Vector{Int},Int} = Int[],
    local_measurements::Vector{String} = String[],
)

Export the binned local measurements for specified process IDs pIDs to a single HDF5 file. If pIDs = Int[], then binned local measurements for all process IDs are exported. You can specify a subset of specific local measurements using the local_measurements keyword argument. If local_measurements = String[], then all local measurements are exported.

export_correlation_stats_to_csv(
    # ARGUMENTS
    comm::MPI.Comm;
    # KEYWORD ARGUMENTS
    datafolder::String,
    correlation::String,
    type::String,
    space::String,
    h5filename::HDF5.File = "stats.h5",
    decimals::Int = 6,
    scientific_notation::Bool = false,
    delimiter::String = " "
)

export_correlation_stats_to_csv(;
    # KEYWORD ARGUMENTS
    datafolder::String,
    correlation::String,
    type::String,
    space::String,
    h5filename::HDF5.File = "stats.h5",
    decimals::Int = 6,
    scientific_notation::Bool = false,
    delimiter::String = " "
)

Export statistics for specified type of correlation function from HDF5 file to a CSV file.

export_correlation_bins_to_csv(
    # ARGUMENTS
    comm::MPI.Comm;
    # KEYWORD ARGUMENTS
    datafolder::String,
    correlation::String,
    type::String,
    space::String,
    pIDs::Union{Vector{Int},Int} = Int[],
    write_index_key::Bool = true,
    decimals::Int = 6,
    scientific_notation::Bool = false,
    delimiter::String = " "
)

export_correlation_bins_to_csv(;
    # KEYWORD ARGUMENTS
    datafolder::String,
    correlation::String,
    type::String,
    space::String,
    pIDs::Union{Vector{Int},Int} = Int[],
    write_index_key::Bool = true,
    decimals::Int = 6,
    scientific_notation::Bool = false,
    delimiter::String = " "
)

Export the binned data for a specified type of correlation to an CSV file living in the directory /datafolder/type/correlation/space. The type of correlation function is specified by type ∈ ("equal-time", "time-displaced", "integrated"). Where the correlation function is in position or momentum space is given by space ∈ ("momentum", "position"). The pIDs keyword specifies for which process IDs the binned correlation data is exported. If pIDs = Int[], then binned local measurements for all process IDs are exported. If write_index_key = true, then another CSV file is written to the /datafolder/type/correlation/space directory which provides a key on how to interpret the INDEX column appearing in the CSV file containing the binned data.

export_correlation_bins_to_h5(
    # ARGUMENTS
    comm::MPI.Comm;
    # KEYWORD ARGUMENTS
    datafolder::String,
    correlation::String,
    type::String,
    space::String,
    pIDs::Union{Vector{Int},Int} = Int[]
)

export_correlation_bins_to_h5(;
    # KEYWORD ARGUMENTS
    datafolder::String,
    correlation::String,
    type::String,
    space::String,
    pIDs::Union{Vector{Int},Int} = Int[]
)

Export the binned data for a specified type of correlation to an HDF5 file living in the directory /datafolder/type/correlation/space. The type of correlation function is specified by type ∈ ("equal-time", "time-displaced", "integrated"). Where the correlation function is in position or momentum space is given by space ∈ ("momentum", "position"). The pIDs keyword specifies for which process IDs the binned correlation data is exported. If pIDs = Int[], then binned local measurements for all process IDs are exported.