QuantumModels
NQCModels.QuantumModels
— ModuleQuantumModels
Models defined within this module subtype the QuantumModel
and provide potentials as Hermitian matrices and derivatives as arrays of Hermitian matrices.
NQCModels.QuantumModels.AltDebyeSpectralDensity
— TypeAltDebyeSpectralDensity{T} <: SpectralDensity
Standard Debye spectral density but uses an alternative discretization scheme that requires a cutoff parameter ωᵐ
.
References
Najeh Rekik, Chang-Yu Hsieh, Holly Freedman, Gabriel Hanna, J. Chem. Phys. 138, 144106 (2013)
NQCModels.QuantumModels.AnanthModelOne
— TypeAnanthModelOne(a=0.01, b=1.6, c=0.005, d=1.0)
Ananth's simple avoided crossing model (similar to Tully's first model) from J. Chem. Phys. 127, 084114 (2007).
NQCModels.QuantumModels.AnanthModelTwo
— TypeAnanthModelTwo(a=0.04, b=0.01, c=0.005, d=1.0, e=0.7, f=1.6)
Ananth's asymmetric model from J. Chem. Phys. 127, 084114 (2007).
NQCModels.QuantumModels.AndersonHolstein
— TypeNewns-Anderson model
The Anderson Impurity Model (AIM) describes a localized impurity state interacting with a continuous band of bath states. AIM is a foundamental model in condensed matter physics and quantum chemistry, introduced by P.W. Anderson in 1961. The Newns-Anderson model is a generalization of the AIM, which includes the possibility of multiple impurity states and a more complex interaction with the bath. A key advantage of using the AIM lies in its ability to yield analytical solutions for the energy level distribution and the hybridization (coupling) density, making it a powerful tool for theoretical analysis.
NQCModels.QuantumModels.BosonBath
— TypeBosonBath(density::SpectralDensity, N::Integer)
Bosonic bath with given spectral density.
Useful for sampling the bath uncoupled from the spin for spin-boson dynamics.
NQCModels.QuantumModels.DebyeSpectralDensity
— TypeDebyeSpectralDensity{T} <: SpectralDensity
Debye density as detailed in: Xin He, Jian Liu, J. Chem. Phys. 151, 024105 (2019)
NQCModels.QuantumModels.DoubleWell
— TypeDoubleWell(mass=1, ω=1, γ=1, Δ=1)
Two state double well, also called the one-dimensional spin-boson model. See: J. Chem. Phys. 150, 244102 (2019)
NQCModels.QuantumModels.ErpenbeckThoss
— Typestruct ErpenbeckThoss{T<:AbstractFloat} <: QuantumModel
1D two-state Quantum system capable of modelling a molecule adsorbed on a metal surface or a single-molecule junction.
In the two references, all of the parameters are identical except for the particle mass m
and the vertical shift c
applied to the ϵ₀ state. Both references modify the shift to ensure the quantum ground-state has an energy of 0 eV. Note that the mass m
is specified in atomic mass units (amu) not atomic units. If a value for the vertical offset c
is not explicitly provided whne constructing the model, it is automatically determined in the constructor from the Morse potential zero-point energy.
References
- PHYSICAL REVIEW B 97, 235452 (2018)
- J. Chem. Phys. 151, 191101 (2019)
NQCModels.QuantumModels.GatesHollowayElbow
— TypeGatesHollowayElbow()
Simple two state elbow potential from Gates and Holloway: Journal of Electron Spectroscopy and Related Phenomena, 64/65 (1993) 633-639
Has two diabatic states each comprised of the sum of a Morse and a repulsive potential. The coupling between them is an exponential function of z
(distance from the surface).
NQCModels.QuantumModels.MiaoSubotnik
— TypeMiaoSubotnik{T<:AbstractFloat} <: QuantumModel
Double well model with parameters matching those of Miao and Subotnik in the reference. This model should be paired with the AndersonHolstein
model to couple to the bath of metallic states.
References
- J. Chem. Phys. 150, 041711 (2019)
NQCModels.QuantumModels.OhmicSpectralDensity
— TypeOhmicSpectralDensity{T} <: SpectralDensity
Ohmic density as detailed in: Xin He, Jian Liu, J. Chem. Phys. 151, 024105 (2019)
NQCModels.QuantumModels.OuyangModelOne
— TypeOuyangModelOne(A=0.01, B=1.6, Γ=1e-4, N=10, ΔE=1.6e-2, D=1.0)
Model #1 from Ouyang and Subotnik. See also Ouyang's thesis.
NQCModels.QuantumModels.QuantumFrictionModel
— TypeQuantumFrictionModel <: QuantumModel
These models are defined identically to a typical QuantumModel
but allocate extra temporary arrays when used with NQCDynamics.jl
.
This allows for the calculation of electronic friction internally from the diabatic potential after diagonalisation and calculation of nonadiabatic couplings.
When a molecular dynamics with electronic friction simulation is set up in NQCDynamics, the QuantumFrictionModel
is paired with a FrictionEvaluationMethod
in order to calculate the electronic friction from the potential and derivative matrices.
NQCModels.QuantumModels.QuantumModel
— TypeQuantumModel <: Model
QuantumModel
s are used when a system has multiple electronic states and the dynamics of the system are propagated by a Hamiltonian in the diabatic representation. This is the case for the majority of model systems.
Implementation
QuantumModel
s should implement:
potential!(model, R)
derivative!(model, D, R)
nstates(model)
ndofs(model)
Example
In this example we create a simple 2 state, 1 dimensional quantum model MyModel
. As noted above, we implement the 4 relevant functions then evaluate the potential. Potential and Derivative functions take in positions as abstract matrices, since this is a 1D model the argument R
should be a Real
wrapped in a 1x1 matrix. It is recommended that you use the hcat() function to do this.
using StaticArrays: SMatrix
using LinearAlgebra: Hermitian
struct MyModel <: NQCModels.QuantumModels.QuantumModel end
NQCModels.nstates(::MyModel) = 2
NQCModels.ndofs(::MyModel) = 1
function NQCModels.potential!(::MyModel, V::Hermitian, R::AbstractMatrix)
V11 = R[1]
V22 = -R[1]
V12 = 1
V = Hermitian(SMatrix{2,2}(V11, V12, V12, V22))
end
function NQCModels.derivative!(::MyModel, D::Hermitian, R::AbstractMatrix)
D = Hermitian(SMatrix{2,2}(1, 0, 0, 1))
end
model = MyModel()
V = Hermitian(zeros(2,2))
NQCModels.potential!(model, V, hcat(10))
# output
2×2 Hermitian{Int64, SMatrix{2, 2, Int64, 4}}:
10 1
1 -10
NQCModels.QuantumModels.SpinBoson
— TypeSpinBoson(density::SpectralDensity, N::Integer, ϵ, Δ)
Spin boson model with N
bosons with given spectral density.
References
NQCModels.QuantumModels.ThreeStateMorse
— TypeThreeStateMorse()
Three state morse potential referred to as Model IA here: J. Chem. Phys. 150, 244102 (2019)
Models IB and IC retain the same functional form and need only a change of parameters.
NQCModels.QuantumModels.TullyModelOne
— TypeTullyModelOne(a=0.01, b=1.6, c=0.005, d=1.0)
Tully's simple avoided crossing model from J. Chem. Phys. 93, 1061 (1990).
NQCModels.QuantumModels.TullyModelThree
— TypeTullyModelThree(a=0.0006, b=0.1, c=0.9)
Tully's extended coupling with reflection model from J. Chem. Phys. 93, 1061 (1990).
NQCModels.QuantumModels.TullyModelTwo
— TypeTullyModelTwo(a=0.1, b=0.28, c=0.015, d=0.06, e=0.05)
Tully's dual avoided crossing model from J. Chem. Phys. 93, 1061 (1990).
NQCModels.QuantumModels.WideBandBath
— Typestruct: WideBandBath <: QuantumFrictionModel
An explicit impurity-bath model which couples a quantum model to a discritised bath. The bath density of states is explicitly given by the wide-band approximation and the bath is discretised using the trapezoidal rule.
It is recommended that this is the explicit bath model used in conjuction with MDEF and the friction evaluators defined in NQCCalculators.jl. For the modelling of explicit bathstate dynamics (i.e. using dynamics methods such as IESH) it is recommended that the user use the Anderson-Holstein model alongside the Shenvi-Gauss Legendre discretisation scheme.
NQCModels.QuantumModels.discretize
— MethodDiscretize a given spectral density for N oscillators. Returns frequencies and couplings.
NQCModels.derivative!
— Methodfunction NQCModels.derivative!(model::AndersonHolstein, D::AbstractMatrix{<:Hermitian}, R::AbstractMatrix)
output: nothing
Updates the derivative of the Anderson-Holstein Hamiltonian with respect to all spatial degrees of freedom to have the correct values for a given position R.
This function is multiple dispatched over the shape of derivative(model.impurity_model) as these impurities may be defined over different numbers of spatial degrees of freedom.
NQCModels.derivative!
— Methodfunction NQCModels.derivative!(model::WideBandBath, D::AbstractMatrix{<:Hermitian}, R::AbstractMatrix)
output: nothing
Updates the derivative of the Anderson-Holstein Hamiltonian with respect to all spatial degrees of freedom to have the correct values for a given position R.
This fucntion is multiple dispatched over the shape of derivative(model.model) as these sub-models may be defined over different numbers of spatial degrees of freedom.
NQCModels.potential!
— Methodpotential!(model::GatesHollowayElbow, V::Hermitian, R::AbstractMatrix)
finds the 2x2 potential matrix for a particle with an internal degree of freedom, x, and a distance from the surface, z, modelled as the sum of a repulsive and Morse potential. As it's currently implemented, the model is only defined for a single particle as the potential is only caclulated based on a single pair of x and z values.
NQCModels.potential
— Methodpotential(model::GatesHollowayElbow, R::AbstractMatrix)
finds the 2x2 potential matrix for a particle with an internal degree of freedom, x, and a distance from the surface, z, modelled as the sum of a repulsive and Morse potential. As it's currently implemented, the model is only defined for a single particle as the potential is only caclulated based on a single pair of x and z values.