API
ModelGeometry Type and Methods
VariationalMC.ModelGeometry
— TypeModelGeometry( unit_cell::UnitCell,
lattice::Lattice,
bond::Vector{Vector{Any}} )
A type defining model geometry.
Internal Methods
VariationalMC.x
VariationalMC.y
VariationalMC.d
VariationalMC.reduce_index_2d
VariationalMC.reduce_index_1d
VariationalMC.max_dist
VariationalMC.x
— Functionx( i::Int, model_geometry::ModelGeometry )
Convenience function for obtaining the x-coordinate of a lattice site given a lattice spindex.
VariationalMC.y
— Functiony( i::Int, model_geometry::ModelGeometry )
Convenience function for obtaining the y-coordinate of a lattice site given a lattice spindex.
VariationalMC.d
— Functiond( p1::Int, p2::Int, model_geometry::ModelGeometry )
Given lattice indices i and j, returns the distances between those 2 points, accounting for the latticed edges with different boundary conditions.
VariationalMC.reduce_index_2d
— Functionreduce_index_2d( i::Int, j::Int, model_geometry::ModelGeometry )
Reduces the indices of 2 lattice sites (i,j) to irreducible indices (0,k), where k is an integer.
VariationalMC.reduce_index_1d
— Functionreduce_index_1d( i::Int, j::Int, model_geometry::ModelGeometry )
Reduces the indices of 2 lattice sites (i,j) to irreducible indices (0,k), where k is an integer.
VariationalMC.max_dist
— Functionmax_dist( N::Int, L::Int )
Obtains the maximum irreducible index given the total number of sites N and extent of the lattice L.
Parameters Types and Methods
VariationalMC.TightBindingModel
— TypeTightBindingModel( t₀::Float64,
t₁::Float64 )
A type defining a non-interacting tight binding model.
VariationalMC.SpinModel
— TypeSpinModel( J₁::Float64,
J₂::Float64,
J₃::Float64 )
A type defining a spin model.
VariationalMC.DeterminantalParameters
— TypeDeterminantalParameters( pars::Vector{AbstractString},
vals::Vector{AbstractFloat},
num_detpars::Int )
A type defining a set of variational parameters for the determinantal wavefunction.
VariationalMC.JastrowParameters
— TypeJastrowParameters( jastrow_type::String,
jpar_map::OrderedDict{Any, Any},
num_jpars::Int,
num_jpar_opts::Int )
A type defining quantities related to Jastrow variational parameters.
Internal Methods
VariationalMC.collect_parameters
— Functioncollect_parameters( determinantal_parameters::DeterminantalParameters,
jastrow_parameters::JastrowParameters )::Vector{AbstractFloat}
Concatenates all values of determinantal and Jastrow parameters into a single vector.
VariationalMC.update_parameters!
— Functionupdate_parameters!( new_vpars::AbstractVector,
determinantal_parameters::DeterminantalParameters )::Nothing
Updates variational (determinantal) parameters after Stochastic Reconfiguration.
update_parameters!( new_vpars::AbstractVector,
determinantal_parameters::DeterminantalParameters,
jastrow_parameters::JastrowParameters )::Nothing
Updates variational (determinantal) parameters after Stochastic Reconfiguration.
VariationalMC.readin_parameters
— Functionreadin_parameters( filename::String )
Parses TOML file containing initial variational parameters.
filename::String
: name of parameter summary file in TOML format.
Hamiltonian Methods
Internal Methods
VariationalMC.build_auxiliary_hamiltonian
VariationalMC.build_tight_binding_hamiltonian
VariationalMC.build_variational_hamiltonian
VariationalMC.add_pairing_symmetry!
VariationalMC.add_spin_order!
VariationalMC.add_charge_order!
VariationalMC.add_chemical_potential!
VariationalMC.diagonalize
VariationalMC.is_openshell
VariationalMC.get_variational_matrices
VariationalMC.get_tb_chem_pot
VariationalMC.build_auxiliary_hamiltonian
— Functionbuild_auxiliary_hamiltonian( tight_binding_model::TightBindingModel,
determinantal_parameters::DeterminantalParameters,
pht::Bool )
Constructs a complete Hamiltonian matrix by combining the non-interacting matrix with matrices of variational terms.
tight_binding_model::TightBindingModel
: parameters for a non-interacting tight-binding model.determinantal_parameters::DeterminantalParameters
: set of determinantal variational parameters.optimize::NamedTuple
: field of optimization flags.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.pht::Bool
: whether model is particle-hol transformed.
VariationalMC.build_tight_binding_hamiltonian
— Functionbuild_tight_binding_hamiltonian( tight_binding_model::TightBindingModel,
model_geometry::ModelGeometry,
pht::Bool )
Constructs a 2N by 2N Hamiltonian matrix where N is the number of lattice sites, given tight binding parameters t, and t'.
tight_binding_model::TightBindingModel
: parameters for a non-interacting tight-binding model.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.pht::Bool
: whether model is particle-hole transformed.
VariationalMC.build_variational_hamiltonian
— Functionbuild_variational_hamiltonian( determinantal_parameters::DeterminantalParameters,
optimize::NamedTuple,
pht::Bool )
Constructs 2N by 2N matrices to be added to the non-interacting tight binding Hamiltonian for each variational parameter, where N is the number of lattice sites. Returns a vector of the sum of matrices and a vector of individual matrix terms.
determinantal_parameters::DeterminantalParameters
: set of determinantal variational parameters.optimize::NamedTuple
: field of optimization flags.pht::Bool
: whether model is particle-hole transformed.
VariationalMC.add_pairing_symmetry!
— Functionadd_pairing_symmetry!( symmetry::String,
determinantal_parameters::DeterminantalParameters,
optimize::NamedTuple,
H_vpars::AbstractMatrix{<:Complex},
V,
model_geometry::ModelGeometry,
pht::Bool )::Nothing
Adds a pairing symmetry term to the auxiliary Hamiltonian.
symmetry::String
: type of pairing symmetry: "s", "d", "sPDW", or "dPDWdeterminantal_parameters::DeterminantalParameters
: set of determinantal variational parameters.optimize::NamedTuple
: field of optimization flags.H_vpars::Vector{Any}
: vector of variational Hamiltonian matrices.V::Vector{Any}
: vector of variational operators.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.pht::Bool
: whether model is particle-hole transformed.
VariationalMC.add_spin_order!
— Functionadd_spin_order!( order::String,
determinantal_parameters::DeterminantalParameters,
optimize::NamedTuple,
H_vpars::Vector{Any},
V::Vector{Any},
model_geometry::ModelGeometry,
pht::Bool )::Nothing
Adds a spin order term to the auxiliary Hamiltonian.
order::String
: type of spin order: "spin-x", "spin-z", or "site-dependent"determinantal_parameters::DeterminantalParameters
: set of determinantal variational parameters.optimize::NamedTuple
: field of optimization flags.H_vpars::Vector{Any}
: vector of variational Hamiltonian matrices.V::Vector{Any}
: vector of variational operators.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.pht::Bool
: whether model if particle-hole transformed.
VariationalMC.add_charge_order!
— Functionadd_charge_order!( order::String,
determinantal_parameters::DeterminantalParameters,
optimize::NamedTuple,
H_vpars::Vector{Any},
V::Vector{Any},
model_geometry::ModelGeometry,
pht::Bool )::Nothing
Adds a charge order term to the auxiliary Hamiltonian.
order::String
: type of spin order: "density wave" or "site-dependent"determinantal_parameters::DeterminantalParameters
: set of determinantal variational parameters.optimize::NamedTuple
: field of optimization flags.H_vpars::Vector{Any}
: vector of variational Hamiltonian matrices.V::Vector{Any}
: vector of variational operators.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.pht::Bool
: whether model if particle-hole transformed.
VariationalMC.add_chemical_potential!
— Functionadd_chemical_potential!( determinantal_parameters::DeterminantalParameters,
optimize::NamedTuple,
H_vpars::Vector{Any},
V::Vector{Any},
model_geometry::ModelGeometry,
pht::Bool )::Nothing
Adds a chemical potential term to the auxiliary Hamiltonian.
determinantal_parameters::DeterminantalParameters
: set of determinantal variational parameters.optimize::NamedTuple
: field of optimization flags.H_vpars::Vector{Any}
: vector of variational Hamiltonian matrices.V::Vector{Any}
: vector of variational operators.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.pht::Bool
: whether model if particle-hole transformed.
VariationalMC.diagonalize
— Functiondiagonalize( H::Matrix{ComplexF64} )::Tuple{Vector{Float64}, Matrix{ComplexF64}}
Returns all eigenenergies and all eigenstates of the mean-field Hamiltonian, the latter being stored in the columns of a matrix. Convention: H(diag) = U⁺HU.
H::AbstractMatrix{<:Complex}
: auxiliary Hamiltonian.
VariationalMC.is_openshell
— Functionis_openshell( ε::Vector{Float64},
Ne::Int )::Bool
Checks whether a energy configuration is open shell.
VariationalMC.get_variational_matrices
— Functionget_variational_matrices( V::Vector{Any},
U_int::Matrix{ComplexF64},
ε::Vector{Float64},
model_geometry::ModelGeometry )::Vector{Any}
Returns variational parameter matrices Aₖ from the corresponding Vₖ. Computes Qₖ = (U⁺VₖU)(ην) / (εη - ε_ν), for η > Nₚ and ν ≤ Nₚ and is 0 otherwise (η and ν run from 1 to 2N).
V::Vector{Any}
: vector of variational operators.U_int::Matrix{ComplexF64}
: matrix which diagonalizes the auxiliary Hamiltonian.ε::Vector{Float64}
: initial energies.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.
VariationalMC.get_tb_chem_pot
— Functionget_tb_chem_pot( Ne::Int64,
tight_binding_model::TightBindingModel,
model_geometry::ModelGeometry )::Float64
For a tight-binding model that has not been particle-hole transformed, returns the chemical potential.
Ne::Int64
: total number of electrons.tight_binding_model::TightBindingModel
: parameters for a non-interacting tight-binding model.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.
DeterminantalWavefunction Type and Methods
VariationalMC.DeterminantalWavefunction
— TypeDeterminantalWavefunction( W::Matrix{ComplexF64},
D::Matrix{ComplexF64},
M::Matrix{ComplexF64}
U_int::Matrix{ComplexF64},
A::Vector{Any},
ε::Vector{Float64},
pconfig::Vector{Int64} )
A type defining quantities related to a determinantal wavefunction.
VariationalMC.get_determinantal_wavefunction
— Functionget_determinantal_wavefunction( tight_binding_model::TightBindingModel,
determinantal_parameters::DeterminantalParameters,
Ne::Int64,
nup::Int64,
ndn::Int64,
model_geometry::ModelGeometry,
rng::Xoshiro)::DeterminantalWavefunction
Constructs a variational wavefunction based on parameters given by the tight binding model and determinantal parameters. Returns an instances of the DeterminantalWavefunction type. If no particle configuration is specified, a random configuration will be generated.
tight_binding_model::TightBindingModel
: parameters for a non-interacting tight-binding model.determinantal_parameters::DeterminantalParameters
: set of determinantal variational parameters.optimize::NamedTuple
: field of optimization flags.Ne::Int
: total number of electrons.nup::Int
: number of spin-up electrons.ndn::Int
: number of spin-down electrons.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.rng::Xoshiro
: random number.pconfig::Union{Nothing, Vector{Int}}=nothing
: current particle configuration.
Jastrow Types and Methods
VariationalMC.JastrowFactor
— TypeJastrowFactor( Tvec_f::Vector{Float64},
Tvec_b::Vector{Float64} )
A type defining quantities related to a Jastrow factor.
VariationalMC.get_jastrow_factor
— Functionget_jastrow_factor( jastrow_parameters::JastrowParameters,
detwf::DeterminantalWavefunction,
model_geometry::ModelGeometry,
pht::Bool )::JastrowFactor
Constructs specified Jastrow factor and returns a instance of the JastrowFactor type.
jastrow_parameters::JastrowParameters
: current set of Jastrow parameters.detwf::DeterminantalWavefunction
: current variational wavefunction.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.pht::Bool
: whether model is particle-hole transformed.
Internal Methods
VariationalMC.get_fermionic_Tvec
VariationalMC.update_fermionic_Tvec!
VariationalMC.get_fermionic_jastrow_ratio
VariationalMC.map_jastrow_parameters
VariationalMC.get_fermionic_Tvec
— Functionget_fermionic_Tvec( jastrow_parameters::JastrowParameters,
detwf::DeterminantalWavefunction,
pht::Bool,
model_geometry::ModelGeometry )::Vector{Float64}
Returns T vector with entries of the form Tᵢ = ∑ⱼ vᵢⱼnᵢ(x) where vᵢⱼ are the associated Jastrow peseudopotentials and nᵢ(x) is the total electron occupation.
jastrow_parameters::JastrowParameters
: current set of Jastrow variational parameters.detwf::DeterminantalWavefunction
: current variational wavefunction.pht::Bool
: whether model is particle-hole transformed.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.
VariationalMC.update_fermionic_Tvec!
— Functionupdate_fermionic_Tvec!( markov_move::MarkovMove,
spin::Int64,
jastrow_parameters::JastrowParameters,
jastrow_factor::JastrowFactor,
model_geometry::ModelGeometry,
n_stab_T::Int64,
δT::Float64,
pht::Bool )::Nothing
Updates elements Tᵢ of the T vector after an accepted Metropolis step.
markov_move::MarkovMove
: quantities related to a Markov move.spin::Int64
: spin of the current particle.jastrow_parameters::JastrowParameters
: current set of Jastrow variational parameters.jastrow::Jastrow
: current Jastrow factor.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.n_stab_T::Int64
: frequency of T vector stabilization setps.δT::Float64
: deviation threshold for the T vector.pht::Bool
: whether model is particle-hole transformed.
VariationalMC.get_fermionic_jastrow_ratio
— Functionget_fermionic_jastrow_ratio( markov_move::MarkovMove,
jastrow_parameters::JastrowParameters
jastrow_factor::JastrowFactor,
pht::Bool,
spin::Int64,
model_geometry::ModelGeometry )
Calculates ratio J(x₂)/J(x₁) = exp[-s(Tₗ - Tₖ) + vₗₗ - vₗₖ ] of Jastrow factors for particle configurations which differ by a single particle hopping from site 'k' (configuration 'x₁') to site 'l' (configuration 'x₂') using the corresponding T vectors Tₖ and Tₗ, rsepctively.
markov_move::MarkovMove
: quantities related to a Markov move.jastrow_parameters::JastrowParameters
: current set of Jastrow variational parameters.jastrow_factor::JastrowFactor
: current Jastrow factor.pht::Bool
: whether model is particle-hole transformed.spin::Int64
: spin of the current particle.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.
get_fermionic_jastrow_ratio( k::Int64,
l::Int64,
jastrow_parameters::JastrowParameters,
jastrow_factor::JastrowFactor,
pht::Bool,
spin::Int64,
model_geometry::ModelGeometry )
Calculates ratio J(x₂)/J(x₁) = exp[-s(Tₗ - Tₖ) + vₗₗ - vₗₖ ] of Jastrow factors for particle configurations which differ by a single particle hopping from site 'k' (configuration 'x₁') to site 'l' (configuration 'x₂') using the corresponding T vectors Tₖ and Tₗ, rsepctively.
k::Int64
: initial site of the current particle.l::Int64
: final site of the current particle.jastrow_parameters::JastrowParameters
: current set of Jastrow variational parameters.jastrow_factor
::JastrowFactor: current Jastrow factor.pht
::Bool: whether model is particle-hole transformed.spin
::Int64: spin of the current particle.model_geometry
::ModelGeometry: contains unit cell and lattice quantities.
VariationalMC.map_jastrow_parameters
— Functionmap_jastrow_parameters( model_geometry::ModelGeometry,
rng::Xoshiro )::OrderedDict{Any, Any}
Generates a dictionary of irreducible indices k which reference a tuple consisting of a vector of lattice index pairs (i,j) which generate k, and Jastrow parameters vᵢⱼ. The parameter corresponding to the largest k is automatically initialized to 0, all others are randomly initialized.
model_geometry::ModelGeometry
: contains unit cell and lattice quantities.rng::Xoshiro
: random number generator.
map_jastrow_parameters( model_geometry::ModelGeometry,
vpar_dict::Dict{Symbol, Any} )::OrderedDict{Any, Any}
Generates a dictionary of irreducible indices k which reference a tuple consisting of a vector of lattice index pairs (i,j) which generate k, and Jastrow parameters vᵢⱼ. The parameter corresponding to the largest k is automatically initialized to 0, all others are read in from file.
model_geometry::ModelGeometry
: contains unit cell and lattice qunatities.vpar_dict::Dict{Symbol, Any}
: dictionary of variational parameters from file.
Markov Methods
VariationalMC.local_fermion_update!
— Functionlocal_fermion_update!( detwf::DeterminantalWavefunction,
Ne::Int,
model_geometry::ModelGeometry,
pht::Bool,
δW::Float64,
δT::Float64,
rng::Xoshiro )
Performs a local update to the electronic sector, with a certain number of equilibration steps.
detwf::DeterminantalWavefunction
: current determinantal variational wavefunction.Ne::Int
: total number of electrons.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.n_stab_W
: frequency of Green's function stabilization steps.δW::Float64
: error threshold for the Green's function.rng::Xoshiro
: random number generator.
local_fermion_update!( detwf::DeterminantalWavefunction,
jastrow_factor::JastrowFactor,
jastrow_parameters::JastrowParameters,
Ne::Int,
model_geometry::ModelGeometry,
pht::Bool,
δW::Float64,
δT::Float64,
rng::Xoshiro )
Performs a local update to the electronic sector, with a certain number of equilibration steps.
detwf::DeterminantalWavefunction
: current determinantal variational wavefunction.jastrow_parameters::JastrowParameters
: current set of Jastrow variational parameters.jastrow::Jastrow
: current Jastrow factor.Ne::Int
: total number of electrons.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.n_stab_W
: frequency of Green's function stabilization steps.n_stab_T
: frequency of T vector stabilization steps.δW::Float64
: error threshold for the Green's function.rng::Xoshiro
: random number generator.
Internal Methods
VariationalMC.metropolis_step
— Functionmetropolis_step( detwf::DeterminantalWavefunction,
Ne::Int,
n_stab_W::Int64,
δW::Float64,
model_geometry::ModelGeometry,
rng::Xoshiro )::String
Proposes a particle to hop to a random neighboring site, and then accepts or rejects using the Metropolis algorithm.
detwf::DeterminantalWavefunction
: current determinantal variational wavefunction.Ne::Int
: total number of electrons.n_stab_W::Int64
: frequency of Green's function stabilization steps.δW::Float64
: error threshold for the Green's function.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.rng::Xoshiro
: random number generator.
metropolis_step( detwf::DeterminantalWavefunction,
jastrow_factor::JastrowFactor,
jastrow_parameters::JastrowParameters,
Ne::Int,
n_stab_W::Int64,
n_stab_T::Int64,
δW::Float64,
δT::Float64,
model_geometry::ModelGeometry,
pht::Bool,
rng::Xoshiro )::String
Proposes a particle to hop to a random neighboring site, and then accepts or rejects using the Metropolis algorithm.
detwf::DeterminantalWavefunction
: current variational wavefunction.jastrow_factor::JastrowFactor
: current Jastrow factor.jastrow_parameters::JastrowParameters
: current set of Jastrow parameters.Ne::Int
: total number of electrons.n_stab_W
: frequency of Green's function stabilization steps.n_stab_T
: frequency of T vector stabilization steps.δW::Float64
: error threshold for the Green's function.δT::Float64
: error threshold for the T vector.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.pht::Bool
: whether model is particle-hole transformed.rng::Xoshiro
: random number generator.
ParticleConfiguration Types and Methods
VariationalMC.get_particle_numbers
— Functionget_particle_numbers( density::Float64 )::NTuple{4, Int64}
Given a particle density, returns the total number of particles Np, number of spin-up particles Npu, number of spin-down particles Npd, and total number of electrons Ne.
density::Float64
: desired electronic density.
VariationalMC.get_particle_density
— Functionget_particle_density( nup::Int64,
ndn::Int64 )::Tuple{Float64, Int64, Int64}
Given the number of spin-up electrons nup, and number of spin-down electrons ndn, returns the particle density, total number of particles Np, and the total number of electrons Ne.
nup::Int
: number of spin-up electrons.ndn::Int
: number of spin-down electrons.
Internal Types and Methods
VariationalMC.MarkovMove
VariationalMC.propose_random_move
VariationalMC.hop!
VariationalMC.generate_initial_fermion_configuration
VariationalMC.get_onsite_fermion_occupation
VariationalMC.get_spindex_type
VariationalMC.get_index_from_spindex
VariationalMC.get_spindices_from_index
VariationalMC.get_linked_spindex
VariationalMC.MarkovMove
— TypeMarkovMove( particle::Int,
k::Int,
l::Int,
possible::Bool )
A type defining a Markov move.
VariationalMC.propose_random_move
— Functionpropose_random_move( Ne::Int64,
pconfig::Vector{Int64},
model_geometry::ModelGeometry,
rng::Xoshiro )::MarkovMove
Proposes randomly hopping or exchanging a particle from some intial site 'k' to a neighboring site 'l' and returns an instance of the MarkovMove type.
Ne::Int64
: total number of electrons.pconfig::Vector{Int64}
: current particle configuration.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.rng::Xoshiro
: random number generator.
VariationalMC.hop!
— Functionhop!( markov_move::MarkovMove,
pconfig::Vector{Int} )::Nothing
If proposed hopping move is accepted, updates the particle positions.
markov_move::MarkovMove
: quantities related to a Markov move.pconfig::Vector{Int}
: current particle configuration.
VariationalMC.generate_initial_fermion_configuration
— Functiongenerate_initial_fermion_configuration( nup::Int64,
ndn::Int64,
model_geometry::ModelGeometry,
rng::Xoshiro )::Vector{Int64}
Generates a random initial configuration of spin-up and spin-down fermions. The first N elements correspond to spin-up and the last N correspond to spin-down. Occupation is denoted by a positive integer corresponding to that particle's creation operator label.
nup::Int64
: number of spin-up electrons.ndn::Int64
: number of spin-down electrons.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.rng::Xoshiro
: random number generator.
VariationalMC.get_onsite_fermion_occupation
— Functionget_onsite_fermion_occupation( site::Int,
pconfig::Vector{Int} )::Tuple{Int64, Int64, Int64}
Returns the number of spin-up and spin-down electrons occupying a real lattice site 'i'.
site::Int
: lattice site.pconfig::Vector{Int}
: current particle configuration.
VariationalMC.get_spindex_type
— Functionget_spindex_type( spindex::Int,
model_geometry::ModelGeometry )::Int
Returns the spin species at a given spindex.
spindex::Int
: spin index.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.
VariationalMC.get_index_from_spindex
— Functionget_index_from_spindex( spindex::Int,
model_geometry::ModelGeometry )::Int
Returns the lattice site i for a given spindex.
spindex::Int
: spin index.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.
VariationalMC.get_spindices_from_index
— Functionget_spindices_from_index( index::Int,
model_geometry::ModelGeometry )::Tuple{Int64, Int64}
Returns spin-up and spin-down indices from a given site index.
index::Int
: lattice site index.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.
VariationalMC.get_linked_spindex
— Functionget_linked_spindex( i::Int,
N::Int )::Int
Returns an index in the spin-down sector, given an index in the spin-up sector.
- i::Int: lattice index in the spin-up sector.
- N::Int: total number of lattice sites.
Optimizer Methods
VariationalMC.stochastic_reconfiguration!
— Functionstochastic_reconfiguration!( measurement_container::NamedTuple,
determinantal_parameters::DeterminantalParameters,
η::Float64,
dt::Float64,
bin::Int,
bin_size::Int64 )::Nothing
Updates variational parameters through stochastic optimization.
measurement_container::NamedTuple
: container where measurements are stored.determinantal_parameters::DeterminantalParameters
: set of determinantal variational parameters.η::Float64
: optimization stabilization factor.dt::Float64
: optimization rate.bin::Int
: current bin number.bin_size::Int64
: length of the current bin.
stochastic_reconfiguration!( measurement_container::NamedTuple,
determinantal_parameters::DeterminantalParameters,
jastrow_parameters::JastrowParameters,
η::Float64,
dt::Float64,
bin::Int,
bin_size::Int64 )::Nothing
Updates variational parameters through stochastic optimization.
measurement_container::NamedTuple
: container where measurements are stored.determinantal_parameters::DeterminantalParameters
: set of determinantal variational parameters.jastrow_parameters::JastrowParameters
: set of Jastrow variational parameters.η::Float64
: optimization stabilization factor.dt::Float64
: optimization rate.bin::Int
: current bin number.bin_size::Int64
: length of the current bin.
Internal Methods
VariationalMC.get_Δk
— Functionget_Δk( optimize::NamedTuple,
determinantal_parameters::DeterminantalParameters,
detwf::DeterminantalWavefunction,
model_geometry::ModelGeometry,
Ne::Int )::Vector{Float64}
Calculates the local logarithmic derivative Δₖ(x) = ∂lnΨ(x)/∂αₖ, with respect to the kth variational parameter αₖ, in the determinantal part of the wavefunction. Returns a vector of derivatives.
optimize::NamedTuple
: field of optimization flags.determinantal_parameters::DeterminantalParameters
: set of determinantal variational parameters.detwf::DeterminantalWavefunction
: current variational wavefunction.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.Ne::Int
: total number of electrons.
get_Δk( optimize::NamedTuple,
determinantal_parameters::DeterminantalParameters,
jastrow_parameters::JastrowParameters,
detwf::DeterminantalWavefunction,
model_geometry::ModelGeometry,
Ne::Int,
pht::Bool )::Vector{Float64}
Calculates the local logarithmic derivative Δₖ(x) = ∂lnΨ(x)/∂αₖ, with respect to the kth variational parameter αₖ, in the determinantal part of the wavefunction. Returns a vector of derivatives.
optimize::NamedTuple
: field of optimization flags.determinantal_parameters::DeterminantalParameters
: set of determinantal variational parameters.jastrow_parameters::JastrowParameters
: current set of Jastrow variational parameters.detwf::DeterminantalWavefunction
: current variational wavefunction.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.Ne::Int
: total number of electrons.pht::Bool
: whether model is particle-hole transformed.
VariationalMC.get_covariance_matrix
— Functionget_covariance_matrix( measurement_container::NamedTuple,
opt_bin_size::Int64 )
Calculates the covariance matrix S, for Stochastic Reconfiguration, with elements Skk' = <ΔkΔk'> - <Δk><Δk'>.
measurement_container::NamedTuple
: container where measurements are stored.opt_bin_size::Int64
: length of the current bin.
VariationalMC.get_force_vector
— Functionget_force_vector( measurement_container::NamedTuple,
opt_bin_size::Int64 )
Generates the force vector f, for Stochastic Reconfiguration with elements fk = <Δk><H> - <Δ_kH>.
measurement_container::NamedTuple
: container where measurements are stored.opt_bin_size::Int64
: length of the current bin.
Measurement Methods
Intitialize Measurements
VariationalMC.initialize_measurement_container
— Functioninitialize_measurement_container( N_opts::Int,
opt_bin_size::Int,
N_bins::Int,
bin_size::Int,
determinantal_parameters::DeterminantalParameters,
model_geometry::ModelGeometry )::NamedTuple
Initializes a set of dictionaries containing generic arrays for storing measurements. Each dictionary in the container has [keys => values] of the form: [observable name => (local value(s), current bin value(s))]
N_opts::Int
: number of optimization updates.opt_bin_size::Int
: length of an optimization bin.N_bins::Int
: number of simulation bins.bin_size::Int
: length of a simulation bin.determinantal_parameters::DeterminantalParameters
: set of determinantal variational parameters.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.
initialize_measurement_container( N_opts::Int,
opt_bin_size::Int,
N_bins::Int,
bin_size::Int,
determinantal_parameters::DeterminantalParameters,
jastrow_parameters::JastrowParameters,
model_geometry::ModelGeometry )::NamedTuple
Initializes a set of dictionaries containing generic arrays for storing measurements. Each dictionary in the container has [keys => values] of the form: [observable name => (local value(s), current bin value(s))]
N_opts::Int
: number of optimization updates.opt_bin_size::Int
: length of an optimization bin.N_bins::Int
: number of simulation bins.bin_size::Int
: length of a simulation bin.determinantal_parameters::DeterminantalParameters
: set of determinantal variational parameters.jastrow_parameters::JastrowParameters
: set of Jastrow variational parameters.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.
VariationalMC.initialize_measurement_directories
— Functioninitialize_measurement_directories( simulation_info::SimulationInfo,
measurement_container::NamedTuple )::Nothing
Creates file directories and for storing measurements.
simulation_info::SimulationInfo
: contains datafolder information.measurement_container::NamedTuple
: container where measurements are contained.
Scalar Measurements
Internal Methods
VariationalMC.get_local_energy
VariationalMC.get_local_kinetic_energy
VariationalMC.get_local_hubbard_energy
VariationalMC.get_double_occ
VariationalMC.get_n
VariationalMC.get_local_energy
— Functionget_local_energy( detwf::DeterminantalWavefunction,
tight_binding_model::TightBindingModel,
model_geometry::ModelGeometry,
Ne::Int64,
pht::Bool )::Float64
Calculates the local variational energy per site for a Hubbard model (without a Jastrow factor).
detwf::DeterminantalWavefunction
: current variational wavefunction.tight_binding_model::TightBindingModel
: parameters for a non-interacting tight-binding model.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.Ne::Int64
: total number of electrons.pht::Bool
: whether model is particle-hole transformed.
get_local_energy( detwf::DeterminantalWavefunction,
tight_binding_model::TightBindingModel,
jastrow_parameters::JastrowParameters,
jastrow_factor::JastrowParameters,
model_geometry::ModelGeometry
Ne::Int64,
pht::Bool )::Nothing
Calculates the local variational energy per site for a Hubbard model.
detwf::DeterminantalWavefunction
: current variational wavefunction.tight_binding_model::TightBindingModel
: parameters for a non-interacting tight-binding model.jastrow_parameters::JastrowParameters
: current set of Jastrow variational parameters.jastrow_factor::JastrowFactor
: current Jastrow factor.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.Ne::Int64
: total number of electrons.pht::Bool
: whether model is particle-hole transformed.
VariationalMC.get_local_kinetic_energy
— Functionget_local_kinetic_energy( detwf::DeterminantalWavefunction,
tight_binding_model::TightBindingModel,
model_geometry::ModelGeometry,
Ne::Int64
pht::Bool )::Float64
Calculates the local electronic kinetic energy.
detwf::DeterminantalWavefunction
: current variational wavefunction.tight_binding_model::TightBindingModel
: parameters for a non-interacting tight-binding model.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.Ne::Int64
: total number of electrons.pht::Bool
: whether model is particle-hole transformed.
get_local_kinetic_energy( detwf::DeterminantalWavefunction,
tight_binding_model::TightBindingModel,
jastrow_parameters::JastrowParameters,
jastrow_factor::JastrowFactor,
model_geometry::ModelGeometry,
Ne::Int64,
pht::Bool )::Float64
Calculates the local electronic kinetic energy.
detwf::DeterminantalWavefunction
: current variational wavefunction.tight_binding_model::TightBindingModel
: parameters for a non-interacting tight-binding model.jastrow_parameters::JastrowParameters
: current set of Jastrow variational parameters.jastrow_factor::JastrowFactor
: current Jastrow factor.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.Ne::Int64
: total number of electrons.pht::Bool
: whether model is particle-hole transformed.
VariationalMC.get_local_hubbard_energy
— Functionget_local_hubbard_energy( U::Float64,
detwf::DeterminantalWavefunction,
model_geometry::ModelGeometry,
pht::Bool )::Float64
Calculates the energy due to onsite Hubbard repulsion.
U::Float64
: Hubbard repulsion.detwf::DeterminantalWavefunction
: current variational wavefunction.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.pht::Bool
: whether model is particle-hole transformed.
VariationalMC.get_double_occ
— Functionget_double_occ( detwf::DeterminantalParameters,
model_geometry::ModelGeometry,
pht::Bool )::Float64
Calculates the double occupancy.
detwf::DeterminantalWavefunction
: current variational wavefunction.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.pht::Bool
: whether model is particle-hole transformed.
VariationalMC.get_n
— Functionget_n( detwf::DeterminantalWavefunction,
model_geometry::ModelGeometry )::Float64
Calculate the local density.
detwf::DeterminantalWavefunction
: current variational wavefunction.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.
Optimization Measurements
Internal Methods
VariationalMC.measure_parameters!
VariationalMC.measure_Δk!
VariationalMC.measure_ΔkΔkp!
VariationalMC.measure_ΔkE!
VariationalMC.measure_parameters!
— Functionmeasure_parameters!( measurement_container::NamedTuple,
determinantal_parameters::DeterminantalParameters )::Nothing
Measures all variational (determinantal) parameters and then writes them to the measurement container.
measurement_container::NamedTuple
: container where measurements are stored.determinantal_parameters::DeterminantalParameters
: current set of determinantal variational parameters.
measure_parameters!( measurement_container::NamedTuple,
determinantal_parameters::DeterminantalParameters,
jastrow_parameters::JastrowParameters )::Nothing
Measures all initialized variational parameters and then writes them to the measurment container. The first 'p' are determinantal parameters and the rest are Jastrow parameters.
measurement_container::NamedTuple
: container where measurements are stored.determinantal_parameters::DeterminantalParameters
: current set of determinantal variational parameters.jastrow_parameters::JastrowParameters
: current set of Jastrow variational parameters.
VariationalMC.measure_Δk!
— Functionmeasure_Δk!( measurement_container::NamedTuple,
detwf::DeterminantalWavefunction,
determinantal_parameters::DeterminantalParameters,
model_geometry::ModelGeometry,
Ne::Int64 )::Nothing
Measures logarithmic derivatives for all variational parameters and then writes them to the measurement container. The first 'p' are derivatives of determinantal parameters.
measurement_container::NamedTuple
: container where measurements are stored.detwf::DeterminantalWavefunction
: current variational wavefunction.determinantal_parameters::DeterminantalParameters
: current set of determinantal variational parameters.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.Ne::Int64
: total number of electrons.
measure_Δk!( measurement_container::NamedTuple,
detwf::DeterminantalWavefunction,
determinantal_parameters::DeterminantalParameters,
jastrow_parameters::JastrowParameters,
model_geometry::ModelGeometry,
Ne::Int64,
pht::Bool )::Nothing
Measures logarithmic derivatives for all variational parameters and then writes them to the measurement container. The first 'p' are derivatives of determinantal parameters and the rest are derivatives of Jastrow parameters.
measurement_container::NamedTuple
: container where measurements are stored.detwf::DeterminantalWavefunction
: current variational wavefunction.determinantal_parameters::DeterminantalParameters
: current set of determinantal variational parameters.jastrow_parameters::JastrowParameters
: current set of Jastrow variational parameters.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.Ne::Int64
: total number of electrons.pht::Bool
: whether model is particle-hole transformed.
VariationalMC.measure_ΔkΔkp!
— Functionmeasure_ΔkΔkp!( measurement_container::NamedTuple,
detwf::DeterminantalWavefunction,
determinantal_parameters::DeterminantalParameters,
model_geometry::ModelGeometry,
Ne::Int64 )::Nothing
Measures the product of variational derivatives with other variational derivatives and then writes them to the measurement container.
measurement_container::NamedTuple
: container where measurements are stored.detwf::DeterminantalWavefunction
: current variational wavefunction.determinantal_parameters::DeterminantalParameters
: current variational determinantal parameters.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.Ne::Int64
: total number of electrons.
measure_ΔkΔkp!( measurement_container::NamedTuple,
detwf::DeterminantalWavefunction,
determinantal_parameters::DeterminantalParameters,
jastrow_parameters::JastrowParameters,
model_geometry::ModelGeometry,
Ne::Int64,
pht::Bool )::Nothing
Measures the product of variational derivatives with other variational derivatives and then writes them to the measurement container.
measurement_container::NamedTuple
: container where measurements are stored.detwf::DeterminantalWavefunction
: current variational wavefunction.determinantal_parameters::DeterminantalParameters
: current set of determinantal variational parameters.jastrow_parameters::JastrowParameters
: current set of Jastrow variational parameters.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.Ne::Int64
: total number of electrons.pht::Bool
: whether model is particle-hole transformed.
VariationalMC.measure_ΔkE!
— Functionmeasure_ΔkE!( measurement_container::NamedTuple,
detwf::DeterminantalWavefunction,
tight_binding_model::TightBindingModel,
determinantal_parameters::DeterminantalParameters,
model_geometry::ModelGeometry,
Ne::Int64,
pht::Bool )::Nothing
Measures the product of variational derivatives with the local energy and then writes them to the measurement container.
measurement_container::NamedTuple
: container where measurements are stored.detwf::DeterminantalWavefunction
: current variational wavefunction.tight_binding_model::TightBindingModel
: parameter for a non-interacting tight-binding model.determinantal_parameters::DeterminantalParameters
: current set of determinantal variational parameters.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.Ne::Int64
: total number of electrons.pht::Bool
: whether model is particle-hole transformed.
measure_ΔkE!( measurement_container::NamedTuple,
detwf::DeterminantalWavefunction,
tight_binding_model::TightBindingModel,
determinantal_parameters::DeterminantalParameters,
jastrow_parameters::JastrowParameters,
model_geometry::ModelGeometry,
Ne::Int64,
pht::Bool )::Nothing
Measures the product of variational derivatives with the local energy and then writes them to the measurement container.
measurement_container::NamedTuple
: container where measurements are stored.detwf::DeterminantalWavefunction
: current variational wavefunction.tight_binding_model::TightBindingModel
: parameters for a non-interacting tight-binding model.determinantal_parameters::DeterminantalParameters
: current set of determinantal variational parameters.jastrow_parameters::JastrowParameters
: current set of Jastrow variational parameters.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.Ne::Int64
: total number of electrons.pht::Bool
: whether model is particle-hole transformed.
Simulation Measurements
Internal Methods
VariationalMC.measure_local_energy!
— Functionmeasure_local_energy!( measurement_container::NamedTuple,
detwf::DeterminantalWavefunction,
tight_binding_model::TightBindingModel,
model_geometry::ModelGeometry,
Ne::Int64,
pht::Bool )::Nothing
Measures the total local energy for a Hubbard model and writes to the measurement container.
measurement_container::NamedTuple
: container where measurements are stored.detwf::DeterminantalWavefunction
: current variational wavefunction.tight_binding_model::TightBindingModel
: parameters for a non-interacting tight-binding model.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.Ne::Int64
: total number of electrons.pht::Bool
: whether model is particle-hole transformed.
measure_local_energy!( measurement_container::NamedTuple,
detwf::DeterminantalWavefunction,
tight_binding_model::TightBindingModel,
jastrow_parameters::JastrowParameters,
jastrow_factor::JastrowFactor,
model_geometry::ModelGeometry,
Ne::Int64,
pht::Bool )::Nothing
Measures the total local energy for a Hubbard model and writes to the measurement container.
measurement_container::NamedTuple
: container where measurements are stored.detwf::DeterminantalWavefunction
: current variational wavefunction.tight_binding_model::TightBindingModel
: parameters for a non-interacting tight-binding model.jastrow_parameters::JastrowParameters
: current set of Jastrow variational parameters.jastrow_factor::JastrowFactor
: current Jastrow factor.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.Ne::Int64
: total number of electrons.pht::Bool
: whether model is particle-hole transformed.
VariationalMC.measure_double_occ!
— Functionmeasure_double_occ!( measurement_container::NamedTuple,
detwf::DeterminantalWavefunction,
model_geometry::ModelGeometry,
pht::Bool )::Nothing
Measure the average double occupancy ⟨D⟩ = N⁻¹ ∑ᵢ ⟨nᵢ↑nᵢ↓⟩.
measurement_container::NamedTuple
: container where measurements are stored.detwf::DeterminantalWavefunction
: current variational wavefunction.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.pht::Bool
: whether model is particle-hole transformed.
VariationalMC.measure_n!
— Functionmeasure_n!( measurement_container::NamedTuple,
detwf::DeterminantalWavefunction,
model_geometry::ModelGeometry )::Nothing
Measure the local particle density ⟨n⟩.
measurement_container::NamedTuple
: container where measurements are stored.detwf::DeterminantalWavefunction
: current variational wavefunction.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.
Correlation Measurements
Internal Methods
VariationalMC.measure_density_correlation!
VariationalMC.measure_spin_correlation!
VariationalMC.get_site_dependent_n
VariationalMC.get_site_dependent_s
VariationalMC.measure_density_correlation!
— Functionmeasure_density_correlation!( measurement_container::NameTuple,
detwf::DeterminantalWavefunction,
model_geometry::ModelGeometry )::Nothing
Measures the equal-time density-density correlation function.
measurement_container::NamedTuple
: container where measurements are stroed.detwf::DeterminantalWavefunction
: current variational wavefunction.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.
VariationalMC.measure_spin_correlation!
— Functionmeasure_spin_correlation!( measurement_container::NamedTuple,
detwf::DeterminantalWavefunction,
model_geometry::ModelGeometry )::Nothing
Measures the equal-time spin-spin correlation function.
measurement_container::NamedTuple
: container where measurements are stroed.detwf::DeterminantalWavefunction
: current variational wavefunction.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.
VariationalMC.get_site_dependent_n
— Functionget_site_dependent_n( detwf::DeterminantalWavefunction,
model_geometry::ModelGeometry )::Vector{Int}
Calculates the density on each lattice site.
detwf::DeterminantalWavefunction
: current variational wavefunction.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.
VariationalMC.get_site_dependent_s
— Functionget_site_dependent_s( detwf::DeterminantalParameters
model_geometry::ModelGeometry )::Vector{Int}
Calculates the spin on each lattice site.
detwf::DeterminantalParameters
: current variational wavefunction.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.
Make Measurements
VariationalMC.make_measurements!
— Functionmake_measurements!( measurement_container::NamedTuple,
detwf::DeterminantalWavefunction,
tight_binding_model::TightBindingModel,
determinantal_parameters::DeterminantalParameters,
optimize::NamedTuple,
model_geometry::ModelGeometry,
Ne::Int64,
pht::Bool )::Nothing
Measures global, optimization, and simulation observables.
measurement_container::NamedTuple
: container where measurements are stored.detwf::DeterminantalWavefunction
: current determinantal wavefunction.tight_binding_model::TightBindingModel
: non-interacting tight-binding model.determinantal_parameters::DeterminantalParameters
: set of determinantal variational parameters.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.Ne::Int64
: number of electrons.pht::Bool
: whether or not model is particle-hole transformed.
make_measurements!( measurement_container::NamedTuple,
detwf::DeterminantalWavefunction,
tight_binding_model::TightBindingModel,
determinantal_parameters::DeterminantalParameters,
jastrow_parameters::JastrowParameters,
jastrow_factor::JastrowFactor,
model_geometry::ModelGeometry,
Ne::Int64,
pht::Bool )::Nothing
Measures global, optimization, and simulation observables.
measurement_container::NamedTuple
: container where measurements are stored.detwf::DeterminantalWavefunction
: current determinantal wavefunction.tight_binding_model::TightBindingModel
: non-interacting tight-binding model.determinantal_parameters::DeterminantalParameters
: set of determinantal variational parameters.jastrow_parameters::JastrowParameters
: current set of Jastrow variational parameters.jastrow_factor::JastrowFactor
: current Jastrow factor.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.Ne::Int64
: number of electrons.pht::Bool
: whether or not model is particle-hole transformed.
make_measurements!( measurement_container::NamedTuple,
detwf::DeterminantalWavefunction,
tight_binding_model::TightBindingModel,
model_geometry::ModelGeometry,
Ne::Int64,
pht::Bool )
Measures simulation observables.
measurement_container::NamedTuple
: container where measurements are stored.detwf::DeterminantalWavefunction
: current determinantal wavefunction.tight_binding_model::TightBindingModel
: non-interacting tight-binding model.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.Ne::Int64
: number of electrons.pht::Bool
: whether or not model is particle-hoel transformed.
make_measurements!( measurement_container::NamedTuple,
detwf::DeterminantalWavefunction,
tight_binding_model::TightBindingModel,
jastrow_parameters::JastrowParameters,
jastrow_factor::JastrowFactor,
model_geometry::ModelGeometry,
Ne::Int64,
pht::Bool )
Measures simulation observables.
measurement_container::NamedTuple
: container where measurements are stored.detwf::DeterminantalWavefunction
: current determinantal wavefunction.tight_binding_model::TightBindingModel
: non-interacting tight-binding model.jastrow_parameters::JastrowParameters
: current set of Jastrow variational parameters.jastrow_factor::JastrowFactor
: current Jastrow factor.model_geometry::ModelGeometry
: contains unit cell and lattice quantities.Ne::Int64
: number of electrons.pht::Bool
: whether or not model is particle-hole transformed.
Internal Methods
VariationalMC.reset_measurements!
— Functionreset_measurements!( measurements::Dict{String, Any} )
Resets value of a dictionary (measurement container) to zero.
Write Measurements
VariationalMC.write_measurements!
— Functionwrite_measurements!( bin::Int,
step::Int,
measurement_container::NamedTuple,
simulation_info::SimulationInfo )::Nothing
Writes optimization and simulation measurements in the current bin to file. Files created are in HDF5 format.
bin::Int
: current bin.step::Int
: current step within the bin.measurement_container::NamedTuple
: container where measurements are stored.simulation_info
: contains datafolder names.
write_measurements!( measurement_container::NamedTuple,
energy_bin::Vector{Any},
dblocc_bin::Vector{Any},
param_bin::Vector{Any} )::Nothing
DEBUG version of the write_measurements!() method. Will write individuallly binned energies, double occupancy, and parameters to specified vectors.
measurement_container::NamedTuple
: container where measurements are stored.energy_bin::Vector{Float64}
: externally specified vector for storing energy measurements.dblocc_bin::Vector{Float64}
: externally specified vector for storing double occupancy measurements.param_bin::Vector{Any}
: externally specified vector for storing parameter measurements.