DynamicsMethods
NQCDynamics.DynamicsMethods
— ModuleThis module contains functions and types necessary for performing nonadiabatic molecular dynamics.
Dynamics is performed using DifferentialEquations.jl
. As such, this module is centered around the implementation of the functions necessary to integrate the dynamics.
For deterministic Hamiltonian methods, the central function is DynamicsMethods.motion!
, which is the inplace form of the function to be integrated by DifferentialEquations.jl
.
Further, methods that have discontinuities, such as surface hopping, use the callback interface provided by DifferentialEquations.jl
.
NQCDynamics.DynamicsMethods.Method
— TypeEach type of dynamics subtypes Method
which is passed to the AbstractSimulation
as a parameter to determine the type of dynamics desired.
NQCDynamics.DynamicsMethods.DynamicsVariables
— MethodDynamicsVariables(::AbstractSimulation, args...)
For each dynamics method this function is implemented to provide the variables for the dynamics in the appropriate format.
By default, DynamicsVariables
is set up for the classical case and takes sim
, v
, r
as arguments and returns a ComponentVector(v=v, r=r)
which is used as a container for the velocities and positions during classical dynamics.
NQCDynamics.DynamicsMethods.create_problem
— MethodProvides the DEProblem for each type of simulation.
NQCDynamics.DynamicsMethods.get_callbacks
— MethodSelect the default callbacks for this simulation type.
NQCDynamics.DynamicsMethods.motion!
— Functionmotion!(du, u, sim, t)
As per DifferentialEquations.jl
, this function is implemented for each method and defines the time derivatives of the DynamicalVariables
.
We require that each implementation ensures du
and u
are subtypes of DynamicalVariables
and sim
subtypes AbstractSimulation
.
NQCDynamics.DynamicsMethods.select_algorithm
— MethodChoose a default algorithm for solving the differential equation.
ClassicalMethods
NQCDynamics.DynamicsMethods.ClassicalMethods.Classical
— TypeClassical <: DynamicsMethods.Method
Type for performing classical molecular dynamics.
sim = Simulation{Classical}(Atoms(:H), Harmonic())
# output
Simulation{Classical}:
Atoms{Float64}([:H], [1], [1837.4715941070515])
Harmonic{Float64, Float64, Float64}
m: Float64 1.0
ω: Float64 1.0
r₀: Float64 0.0
dofs: Int64 1
NQCDynamics.DynamicsMethods.ClassicalMethods.Langevin
— TypeType for performing Langevin molecular dynamics.
using Unitful
sim = Simulation{Langevin}(Atoms(:H), Free(); γ=2.5, temperature=100u"K")
# output
Simulation{Langevin{Float64}}:
Atoms{Float64}([:H], [1], [1837.4715941070515])
Free(1)
NQCDynamics.DynamicsMethods.ClassicalMethods.MDEF
— Type\[dr = v dt\\ dv = -\Delta U/M dt - \Gamma v dt + \sigma \sqrt{2\Gamma} dW\]
$\Gamma$ is the friction tensor with units of inverse time. For thermal dynamics we set $\sigma = \sqrt{kT / M}$, where $T$ is the electronic temperature.
This is integrated using the BAOAB algorithm where the friction "O" step is performed in the tensor's eigenbasis. See src/dynamics/mdef_baoab.jl
for details.
NQCDynamics.DynamicsMethods.ClassicalMethods.ThermalLangevin
— TypeType for performing Langevin ring polymer molecular dynamics.
Currently there are separate types for classical and ring polymer versions of Langevin dynamics but they should be combined. The reason they are not at the moment is that they use different integration algorithms and require slightly different fields.
using Unitful
RingPolymerSimulation{ThermalLangevin}(Atoms(:H), Free(), 10; γ=0.1, temperature=25u"K")
# output
RingPolymerSimulation{ThermalLangevin{Float64}}:
Atoms{Float64}([:H], [1], [1837.4715941070515])
Free(1)
with 10 beads.
NQCDynamics.DynamicsMethods.ClassicalMethods.acceleration!
— Methodf1
in DifferentialEquations.jl
docs.
NQCDynamics.DynamicsMethods.ClassicalMethods.friction!
— Methodfriction!(g, r, sim, t)
Evaluates friction tensor
MappingVariableMethods
NQCDynamics.DynamicsMethods.MappingVariableMethods.NRPMD
— TypeNRPMD{T} <: DynamicsMethods.Method
Nonadiabatic ring polymer molecular dynamics Uses Meyer-Miller-Stock-Thoss mapping variables for electronic degrees of freedom and ring polymer formalism for nuclear degrees of freedom.
RingPolymerSimulation{NRPMD}(Atoms(:H), DoubleWell(), 10)
# output
RingPolymerSimulation{NRPMD{Float64}}:
Atoms{Float64}([:H], [1], [1837.4715941070515])
DoubleWell{Int64, Int64, Int64, Int64}
mass: Int64 1
ω: Int64 1
γ: Int64 1
Δ: Int64 1
with 10 beads.
NQCDynamics.DynamicsMethods.MappingVariableMethods.eCMM
— TypeSurfaceHoppingMethods
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods
— ModuleSurfaceHoppingMethods
Implementation for surface hopping methods.
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.AdiabaticIESH
— TypeIESH{T} <: SurfaceHopping
Independent electron surface hopping.
References
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.BCME
— TypeBCME{T} <: ClassicalMasterEquation
Extension to CME that incorporates broadening in the potential energy surfaces.
Note that we do not rescale the velocity as this is not mentioned only in the 2016 paper, not any of the later ones, so I presume they later decided not to do it and instead keep it the same as the original CME.
- Dou, Subotnik, J. Chem. Phys. 144, 024116 (2016)
- Dou, Subotnik, J. Phys. Chem. A, 24, 757-771 (2020)
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.CME
— TypeCME{T} <: ClassicalMasterEquation
Simple surface hopping method for Newns-Anderson (Anderson-Holstein) models.
- Dou, Nitzan, Subotnik, J. Chem. Phys. 142, 084110 (2015)
- Dou, Subotnik, J. Phys. Chem. A, 24, 757-771 (2020)
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.DecoherenceCorrectionEDC
— TypeDecoherenceCorrectionEDC{T}
Energy decoherence correction of Granucci and Persico in J. Chem. Phys. 126, 134114 (2007).
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.FSSH
— TypeFSSH{T} <: SurfaceHopping
Type for fewest-switches surface hopping
Simulation{FSSH}(Atoms(:H), Free())
# output
Simulation{FSSH{Float64}}:
Atoms{Float64}([:H], [1], [1837.4715941070515])
Free(1)
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.SurfaceHopping
— TypeAbstract type for all surface hopping methods.
Surface hopping methods follow the structure set out in this file. The nuclear and electronic variables are propagated by the motion!
function. The surface hopping procedure is handled by the HoppingCallback
which uses the functions check_hop!
and execute_hop!
as its condition
and affect!
.
To add a new surface hopping scheme, you must create a new struct and define methods for evaluate_hopping_probability!
, select_new_state
, and rescale_velocity!
.
See fssh.jl
for an example implementation.
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.apply_decoherence_correction!
— MethodEq. 17 in J. Chem. Phys. 126, 134114 (2007)
Modify the wavefunction coefficients in ψ
after a successful surface hop.
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.calculate_Akj
— MethodEquation 17 in Shenvi, Roy, Tully 2009. Uses equations 19 and 20.
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.calculate_a
— Methodcalculate_a(sim::AbstractSimulation{<:SurfaceHopping}, coupling::AbstractMatrix)
Equation 40 from [27].
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.calculate_b
— Methodcalculate_b(coupling::AbstractMatrix, velocity::AbstractMatrix)
Equation 41 from [27].
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.compute_overlap!
— MethodEquation 20 in Shenvi, Roy, Tully 2009.
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.evaluate_broadening
— Methodevaluate_broadening(h, μ, β, Γ)
Evaluate the convolution of the Fermi function with a Lorentzian.
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.evaluate_hopping_probability!
— Methodevaluate_hopping_probability!(sim::Simulation{<:FSSH}, u, dt)
Evaluates the probability of hopping from the current state to all other states
Implementation
σ
is Hermitan so the choiceσ[m,s]
orσ[s,m]
is irrelevant; we take the real part.- 'd' is skew-symmetric so here the indices are important.
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.evaluate_hopping_probability!
— MethodHopping probability according to equation 21 in Shenvi, Roy, Tully 2009.
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.evaluate_hopping_probability!
— MethodThis function should set the field sim.method.hopping_probability
.
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.extract_nonadiabatic_coupling
— Methodextract_nonadiabatic_coupling(coupling, new_state, old_state)
Extract the nonadiabatic coupling vector between states new_state
and old_state
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.frustrated_hop_invert_velocity!
— Methodfrustrated_hop_invert_velocity!(
sim::AbstractSimulation{<:SurfaceHopping}, velocity, d
)
Measures the component of velocity along the nonadiabatic coupling vector and inverts that component.
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.iesh_diabatic_population
— MethodCalculate the diabatic population in J. Chem. Theory Comput. 2022, 18, 4615−4626. Eqs. 12, 13 and 14 describe the steps to calculat it though the code is written to use matrix operations instead of summations to calculate the result.
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.perform_rescaling!
— Methodperform_rescaling!(
sim::AbstractSimulation{<:SurfaceHopping}, velocity, velocity_rescale, d
)
Equation 33 from [27].
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.rescale_velocity!
— Methodrescale_velocity!(sim::AbstractSimulation{<:SurfaceHopping}, u)::Bool
Rescale the velocity in the direction of the nonadiabatic coupling.
References
[27]
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.select_new_state
— MethodThis function should return the desired state determined by the probability. Should return the original state if no hop is desired.
NQCDynamics.DynamicsMethods.SurfaceHoppingMethods.unpack_states
— Methodunpack_states(sim)
Get the two states that we are hopping between.
NQCDynamics.DynamicsUtils.acceleration!
— MethodSet the acceleration due to the force from the currently occupied states. See Eq. 12 of Shenvi, Tully JCP 2009 paper.
NQCDynamics.DynamicsUtils.set_quantum_derivative!
— MethodPropagation of electronic wave function happens according to Eq. (14) in the Shenvi, Tully paper (JCP 2009)
In IESH each electron is independent so we can loop through electrons and set the derivative one at a time, in the standard way for FSSH.
EhrenfestMethods
NQCDynamics.DynamicsMethods.EhrenfestMethods.AbstractEhrenfest
— TypeAbstract type for Ehrenfest method.
NQCDynamics.DynamicsMethods.EhrenfestMethods.Ehrenfest
— TypeEhrenfest{T} <: AbstractEhrenfest
Ehrenfest molecular dynamics. Classical molecular dynamics where the force is derived by averaging contributions from multiple electronic states.
Simulation{Ehrenfest}(Atoms(:H), DoubleWell())
# output
Simulation{Ehrenfest{Float64}}:
Atoms{Float64}([:H], [1], [1837.4715941070515])
DoubleWell{Int64, Int64, Int64, Int64}
mass: Int64 1
ω: Int64 1
γ: Int64 1
Δ: Int64 1
IntegrationAlgorithms
NQCDynamics.DynamicsMethods.IntegrationAlgorithms.MInt
— TypeMInt <: OrdinaryDiffEq.OrdinaryDiffEqAlgorithm
Second order symplectic momentum integral algorithm.
Reference
NQCDynamics.DynamicsMethods.IntegrationAlgorithms.RingPolymerMInt
— TypeRingPolymerMInt <: OrdinaryDiffEq.OrdinaryDiffEqAlgorithm
Second order symplectic momentum integral algorithm applied to NRPMD.
Reference
NQCDynamics.DynamicsMethods.IntegrationAlgorithms.get_C_propagator
— MethodGet the C
propagator for the mapping variables.
NQCDynamics.DynamicsMethods.IntegrationAlgorithms.get_D_propagator
— MethodGet the D
propagator for the mapping variables.
NQCDynamics.DynamicsMethods.IntegrationAlgorithms.get_gamma
— MethodGet the Γ
variable used to calculate the nuclear propagators.
NQCDynamics.DynamicsMethods.IntegrationAlgorithms.get_mapping_nuclear_force
— MethodGet the force due to the mapping variables.
Equivalent to this but doesn't allocate: return 0.5 * (q'Eq + p'Ep) - q'Fp
NQCDynamics.DynamicsMethods.IntegrationAlgorithms.get_mapping_nuclear_force
— MethodGet the force due to the mapping variables.
NQCDynamics.DynamicsMethods.IntegrationAlgorithms.get_xi
— MethodGet the Ξ
variable used to calculate the nuclear propagators.
NQCDynamics.DynamicsMethods.IntegrationAlgorithms.set_C_propagator!
— MethodGet the C
propagator for the mapping variables.
NQCDynamics.DynamicsMethods.IntegrationAlgorithms.set_D_propagator!
— MethodGet the D
propagator for the mapping variables.
NQCDynamics.DynamicsMethods.IntegrationAlgorithms.set_gamma!
— MethodGet the Γ
variable used to calculate the nuclear propagators.
NQCDynamics.DynamicsMethods.IntegrationAlgorithms.set_xi!
— MethodGet the Ξ
variable used to calculate the nuclear propagators.
StochasticDiffEq.alg_cache
— MethodInsecting the inputs into a Cache structure.