Unitary Basis Transformations in Mixed Quantum-Classical Dynamics
Ken Miyazaki Alex Krotz Roel Tempelaar




    < g r a p h i c s >



A Bi-directional Quantum Search Algorithm
    Zain Hafeez^∗
Department of Physics, 

Purdue University, West Lafayette, USA


zhafeez@purdue.edu
Debanjan Konar^∗ and Vaneet Aggarwal
Purdue Quantum Science and Engineering Institute (PQSEI), and  

School of Industrial Eng., Purdue University, West Lafayette, USA 


dkonar@purdue.edu and  vaneet@purdue.edu

    
===========================================================================================================================================================================================================================================================================================================================




A common approach to minimizing the cost of quantum computations is by transforming a quantum system into a basis that can be optimally truncated. Here, we derive classical equations of motion subjected to similar unitary transformations, and propose their integration into mixed quantum–classical dynamics, enabling this class of methods to be applied within arbitrary bases for both the quantum and classical coordinates. To this end, canonical positions and momenta are combined into a set of complex-valued classical coordinates amenable to unitary transformations. We demonstrate the potential of the resulting approach by means of surface hopping calculations of an electronic carrier scattering onto a single impurity in the presence of phonons. Appropriate basis transformations, capturing both the localization of the impurity and the delocalization of higher-energy excitations, are shown to faithfully capture the dynamics within a fraction of the classical and quantum basis sets.






§ INTRODUCTION


The application of unitary basis transformations is common practice in quantum-mechanical computations, derivations, and analyses. Computations typically require bases to be truncated, and depending on the quantum system at hand, the effectiveness of such truncations varies with the transformed representation <cit.>. In addition, transformations enable one to construct perturbative expansions with optimal convergence properties, while also allowing quantum-mechanical equations to be cast in their most intuitive form.

This is exemplified by Bloch's theorem <cit.>, which finds widespread application in the modeling of materials. Bloch's theorem invokes a complex Fourier transform of physical basis states over a crystal lattice in order to yield a representation in reciprocal space <cit.>. This representation naturally captures the conservation of lattice momentum, and describes phenomena in terms of quasiparticle bands that can be included or excluded in order to modulate the computational cost <cit.>. Another example is provided by Redfield theory <cit.>, which is formulated within the eigenbasis of a quantum system, thereby allowing the nonadiabatic coupling between (adiabatic) eigenstates to be captured perturbatively. Here, the computational cost can be modulated by including or excluding eigenstates <cit.>. More generally, for a quantum system expressed in a physical basis (also referred to as local basis or “site” basis), an arbitrary unitary transformation takes the form

    |ξ⟩ = ∑_n u_ξ n|n⟩,

with u ∈ℂ and u^† u = 1, and where n and ξ label the physical and transformed basis states, respectively.

In classical mechanics, basis transformations take the form of canonical transformations. For a given set of canonical position and momentum coordinates, denoted {q_n} and {p_n}, respectively, the transformed canonical coordinates take the form q_ξ = q_ξ({q_n, p_n}) and p_ξ = p_ξ({q_n, p_n}). These transformed coordinates serve the role of “position” and “momentum” in the transformed Hamilton equations of motion, respectively, although they do not necessarily correspond to a physical position and momentum. All coordinates are real valued, as is required for classical trajectories. Canonical transformations offer the same benefits to classical mechanics as unitary transformations do to quantum mechanics, yet the interconnection between canonical and unitary transformations is somewhat opaque. Solidifying this interconnection is of particular relevance to mixed quantum–classical (MQC) dynamics, where electronic excitations are described quantum-mechanically and nuclear vibrations (phonons) classically.

MQC dynamics (sometimes referred to as nonadiabatic molecular dynamics) finds widespread application to the transient modeling of excited-state phenomena in molecular systems <cit.>. Virtually all implementations of MQC dynamics adopt a physical basis for the nuclear coordinates. This practice is perhaps motivated by the notion that molecular excited-state phenomena are commonly localized, with site-to-site energy transfer being mediated by local vibrations. Such is indeed captured most efficiently by a physical basis, as localization allows this basis to be truncated to only the spatial regions of interest.

In a recent work <cit.>, henceforth referred to as Paper I, we derived a formulation of MQC dynamics fully within reciprocal space by subjecting both quantum and classical coordinates to a complex Fourier transform. This work was motivated by a surge in applications of MQC dynamics to crystalline materials <cit.>, for which excited-state phenomena are typically delocalized and driven by an exchange of lattice momentum between electronic carriers and Bloch-like phonons <cit.>. Arriving at the reciprocal-space MQC formalism required us to first combine classical positions and momenta within a single complex-valued coordinate, z_n, amenable to a complex Fourier transform. The transformed equations of motion were shown to yield dynamics formally-equivalent to solutions in the physical basis. At the same time, the reciprocal-space representation was shown to allow for basis truncations of band-like dynamics not possible within a physical basis representation, allowing the computational cost to be optimally reduced. In Paper I, reciprocal-space MQC was introduced for mean-field dynamics (sometimes referred to as Ehrenfest dynamics). In a follow-up work, henceforth referred to as Paper II <cit.>, we have extended reciprocal-space MQC for the popular fewest-switches surface hopping (FSSH) method <cit.>. This formulation has since found application in the modeling of the Floquet nonadiabatic dynamics of laser-dressed solid-state materials <cit.>.

While reciprocal-space MQC is particularly effective in describing band-like phenomena, its effectiveness deteriorates once the periodicity of the crystal lattice becomes disrupted, and lattice momentum is no longer a good quantum number. Such disruptions may take the form of defects, such as impurities, vacancies, and dislocations. In such cases, the physical basis may not provide an effective representation either, as excited states may retain substantial delocalization lengths. Instead, the optimal representation will be provided by some other basis. More generally, for any given system, out of all possible bases one should be able to find a choice that most effectively captures the dynamics. In that regard, localized molecular excitations and band-like excitations in materials span the two Fourier-related extremes, captured optimally by the physical and reciprocal bases, respectively.

In this Article, we introduce a formulation of MQC dynamics within arbitrary basis representations for both quantum and classical coordinates, opening the opportunity to transiently model excited-state phenomena within optimized bases for any given system. As in Paper I <cit.>, we combine classical positions and momenta within a complex-valued coordinate, z_n, which is subjected to an arbitrary basis transformation, similarly to the quantum basis states in Eq. <ref>. We then derive the transformed equations of motion, which are integrated in MQC dynamics. In order to demonstrate the utility of this framework, we consider a model invoking an electronic carrier scattering onto a single impurity in an otherwise pristine lattice, under Holstein-type coupling to phonons. For this model, fewest-switches surface hopping (FSSH) calculations are presented, showing the excited state dynamics to be invariant to significant truncations of the appropriately-transformed electronic and nuclear bases.

This Paper is organized as follows. In Sec. <ref> we introduce the transformed classical equations of motion, their integration in MQC dynamics, and in FSSH in particular. In Sec. <ref>, we use the resulting formalism to re-derive the equations of motion of reciprocal-space MQC and FSSH from Papers I <cit.> and II <cit.>. In Sec. <ref> we introduce the single impurity model, and present results from transformed FSSH under basis truncations. In Sec. <ref>, we present our conclusions and offer an outlook for future directions.



§ THEORY


Since unitary basis transformations of quantum systems are trivial and well-established, we instead begin by considering such transformations for classical systems, before presenting the integration of the resultant classical equations of motion within MQC dynamics.



 §.§ Transformed classical dynamics




The general idea behind transformed classical dynamics is schematically depicted in Fig. <ref>. Here, a set of canonical position and momentum coordinates is introduced, which are described by a generic Hamiltonian function of the form

    H({q_n,p_n}) = ∑_np_n^2/2m_n + V({q_n}).

Here, the first term represents kinetic energy, where m_n denotes the mass of coordinate n, and the second term represents potential energy. The time evolution of the canonical position and momentum coordinates is governed by the Hamilton equations of motion,

    q̇_n = ∂ H/∂ p_n = p_n/m_n,   ṗ_n = -∂ H/∂ q_n = -∂ V/∂ q_n.


As in Paper I <cit.> and as illustrated in Fig. <ref>, we combine the classical position and momentum within a single complex-valued coordinate for each n,

    z_n = √(m_nh_n/2)(q_n + ip_n/m_nh_n).

Here, we introduced a variable h_n whose magnitude is arbitrary, and which modulates the scaling of the position and momentum along the real and imaginary axes in the complex plane, respectively. Expressed in terms of the complex-valued coordinates, the generic Hamiltonian function is given by

    H({z_n}) = -∑_nh_n/4 ( z_n^2 -2 z_nz_n^* + z_n^*2) + V({z_n}).
The corresponding Hamilton equations take the form

    ż_n = -i ∂ H/∂ z_n^*.


In Paper I <cit.>, the classical coordinates were assumed to represent identical harmonic modes, meaning that m_n = m and

    V({q_n}) = 1/2 m ω^2 ∑_n q_n^2.

The harmonic frequency ω then took the place of h_n in Eq. <ref>. This choice of h_n maps harmonic motion onto a perfect circle in the complex plane, as shown in Fig. <ref>. As a result, the Hamiltonian function simplifies to H({z_n}) = ω∑_n z_n z_n^* and the associated time evolution follows as ż_n = -i ω  z_n. In that case, z_n may be interpreted as the eigenvalue associated with the coherent state of the harmonic oscillator at coordinate n <cit.>. When h_n ≠ω, harmonic motion is mapped onto an ellipse instead. In the present analysis, no assumption is made about the modes, and by leaving h_n unspecified we derive equations of motion applicable under generic potentials, including anharmonicities.

With the positions and momenta combined into a set of complex coordinates, arbitrary unitary basis transformations can be applied following

    z_ξ = ∑_n U_ξ n z_n,

as depicted in Fig. <ref>. Here, U ∈ℂ and U^† U = 1. Eq. <ref> is the classical equivalent of Eq. <ref>. As can be easily verified, the Hamilton equations expressed in terms of the transformed complex-valued coordinates follow as

    ż_ξ = -i  Hz_ξ^*.
 The transformed Hamiltonian function, on the other hand, takes the slightly more complicated form

    H({z_ξ}) =    -1/4 ∑_ξ,ξ'( h̃^*_ξξ' z_ξ z_ξ' -2 h_ξξ' z_ξ z_ξ'^*
        + h̃_ξξ' z_ξ^* z_ξ'^*) + V({z_ξ}),

where h_ξξ'≡∑_n U^*_ξ n h_n U_ξ' n and h̃_ξξ'≡∑_n U_ξ n h_n U_ξ' n.

Eq. <ref> provides a generic treatment of any classical system by appropriate construction of the potential energy in terms of the complex-valued coordinates, V({z_n}). Eq. <ref> is formally equivalent to Eq. <ref>, and the associated equations of motion given by Eq. <ref> allow the system's time evolution to be determined within an arbitrary basis. Notably, the transformed complex-valued coordinates z_ξ can be decomposed into real-valued canonical coordinates q_ξ and p_ξ, for which the Hamilton equations of motion can be derived straightforwardly. This establishes the connection between the unitary transformation applied in Eq. <ref> and canonical transformations of the form q_ξ = q_ξ({q_n, p_n}) and p_ξ = p_ξ({q_n, p_n}). In the following, however, we will minimize such analyses, and instead resort to the complex-valued coordinates, which provide a simpler and general framework for describing classical dynamics. This has the added benefit of putting classical basis transformations on the same footing as quantum-mechanical basis transformations.



 §.§ Mixed quantum–classical dynamics


We now proceed to integrate the transformed classical coordinates within MQC dynamics. This is not intended as a rigorous introduction into MQC dynamics, for which we refer the interested reader to excellent sources in the literature <cit.>.

MQC dynamics relies on the subdivision of a system of interest into a quantum subsystem and a classical subsystem, where the quantum subsystem is commonly taken to represent the electronic states while the classical subsystem is reserved for nuclear coordinates. Resorting to canonical coordinates within the physical basis, the total Hamiltonian takes the form

    Ĥ_tot({q_n, p_n}) = Ĥ_q + Ĥ_q-c({q_n}) + H_c({q_n, p_n}).

Here, Ĥ_q is the Hamiltonian operator of the quantum subsystem, H_c({q_n, p_n}) is the Hamiltonian function of the classical subsystem, and Ĥ_q-c({q_n}) is the operator governing interactions between the quantum and classical subsystems, which involves a parametric dependence on the physical position coordinates. Adopting transformed coordinates, this yields

    Ĥ_tot({z_ξ}) = Ĥ_q + Ĥ_q-c({z_ξ}) + H_c({z_ξ}).


The quantum–classical interaction term contributes to the evolution of the quantum subsystem, which is governed by the time-dependent Schrödinger equation,

    iħ |Ψ̇⟩ = (Ĥ_q + Ĥ_q-c({z_ξ})) |Ψ⟩,

where Ψ is the quantum wavefunction. It also affects the evolution of the classical subsystem, since the total Hamiltonian function to be used in the Hamilton equation for the classical coordinates (Eq. <ref>) receives potential energy contributions from the quantum–classical interaction as well as those intrinsic to the classical subsystem. That is,

    V({z_ξ}) = V_c({z_ξ}) + V_q-c({z_ξ}),

where the quantum–classical contribution is given by some expectation value of the quantum–classical Hamiltonian, V_q-c({z_ξ}) = ⟨Ĥ_q-c({z_ξ})|$⟩.

Various MQC dynamical methods differ in the way this expectation value is determined. In case of Ehrenfest dynamics, where the quantum–classical interaction is treated as a mean-field problem, the expectation value is taken with respect to the quantum wavefunction, yielding⟨Ĥ_q-c({z_ξ})|=⟩ ⟨Ψ| Ĥ_q-c({z_ξ}) | Ψ|$⟩. For FSSH, on the other hand, the expectation value is taken based on an instantaneous eigenstate α of the total Hamiltonian operator, satisfying

    (Ĥ_q + Ĥ_q-c) |α⟩ = ϵ_α|α⟩.

That is, a single “active surface”, denoted a, is chosen and the expectation value is determined as ⟨Ĥ_q-c({z_ξ})|=⟩⟨a | Ĥ_q-c({z_ξ}) | a|$⟩.

A key ingredient of FSSH is a stochastic switching of the active surface between instantaneous eigenstates. The switching from stateαto stateβis governed by the probability <cit.>
    P_a:α→β = 2( α|∂β/∂ tA_β/A_α)Δ t.

Here,Δtis the time increment for which the switching probability is evaluated, andA_αis the coefficient ofΨexpanded in the instantaneous eigenbasis, i.e.,

    |Ψ⟩ = ∑_α A_α|α⟩.

We note that the switching probability commonly features a product of momentum and the nonadiabatic coupling vectors within the physical basis <cit.>. By means of the chain rule <cit.>, we have replaced this product by an inner product ofαand∂β/∂tin Eq. <ref>, yielding a basis-independent form of the switching probability.

Upon a switch, the conservation of total (quantum plus classical) energy is reinforced by a rescaling of the physical momenta of the classical subsystem in the direction of the nonadiabatic coupling vector. Accordingly,

    p_n' = p_n - γα̃ | ∂/∂ q_n | β̃,

wherep_n'andp_nare the new and old momentum, respectively. Here, the tilde inβ̃refers to the projection of the potentially complex-valued eigenvector|β⟩onto a real-valued vector, in order to ensure that the physical momentum coordinate remains real-valued, and same forα̃. Complex values for eigenvectors may arise due to an arbitrary global phase, complex basis transformations followed by basis truncations <cit.>, and geometric phase effects <cit.>. In a recent work, we have proposed a means to perform this projection onto real-valued vectors while ensuring gauge-invariance <cit.>.
    
Within the transformed basis, the rescaling takes the form

    z_ξ' = z_ξ - iγα̃ | ∂/∂ z_ξ^* | β̃.


 §.§ Reciprocal space


In Papers I <cit.> and II <cit.>, we derived a formulation of MQC dynamics fully within reciprocal space, and applied the resulting method to pristine one-dimensional lattice models involving a single electronic carrier interacting with harmonic nuclear vibrations. In the following, we will show that the generalized equations presented in the current work recover those from Paper I <cit.> when the unitary transformation of the classical subsystem is taken to be a complex Fourier transforms over the lattice, which underlies Bloch's theorem.
    
Accordingly, we replaceξby the wavevectork(which quantifies the lattice momentum), and take

    U_kn = 1/√(N) e^ikn,

whereNis the total number of lattice sites. The associated classical coordinates are taken to represent identical harmonic modes, such thatm_n = mandV_c({q_n})is given by Eq. <ref>, withωthe mode frequency. Settingh_n = ω, the complex-valued classical coordinates take the form

    z_n = √(mω/2)(q_n + i p_n/mω).

Subjecting these coordinates to the transformation given in Eq. <ref>, we arrive at transformed coordinatesz_k. These coordinates are associated with phonons, i.e., nuclear vibrational quasiparticles with a well-defined lattice momentum.
    
Expressed in terms of the transformed coordinates, the purely-classical potential energy contribution is given by

    V_c({z_k}) = ω/4∑_k( z_kz_-k + 2 z_kz_k^* + z_k^*z_-k^*).

As a result, the Hamiltonian to be used in the Hamilton equations reduces to

    H({z_k}) = ω∑_k z_kz_k^* + V_q-c({z_k}).

In accordance with Eq. <ref>, this yields

    ż_k = -i ω z_k - i V_q-c({z_k})z_k^*.


In Paper I <cit.>, canonical coordinates were reconstructed fromz_kfollowing

    q_k≡√(2/mω)z_k,   
        p_k≡√(2mω)z_k,

which obeys

    z_k = √(mω/2)( q_k + i p_k/mω).

Expressed in terms of such canonical coordinates, the Hamilton equations then follow as

    q̇_k = √(2/mω)ż_k = p_k/m + V_q-cp_k,
       ṗ_k = √(2mω)ż_k = -mω^2q_k - V_q-cq_k,

which is indeed in agreement with Paper I <cit.>.

Notably, from the appearance of gradient contributions to bothq̇_kandṗ_k, it can be clearly seen that while the involved canonical coordinates play the role of “position” and “momentum” within the Hamilton equations, they are not to be associated with physical positions and momenta <cit.>.
This is further reflected by the coordinate rescaling applied upon a switch in FSSH. This rescaling is obtained by taking the real and imaginary parts of Eq. <ref>, yielding

    p_k' = p_k - γα̃ | ∂/∂ q_k | β̃, 
       q_k' = q_k + γα̃ | ∂/∂ p_k | β̃,

in agreement with Paper II <cit.>.

In Papers I <cit.> and II <cit.>, application of Eq. <ref> in conjunction with a reciprocal-space description of the quantum subsystem was shown to yield results in agreement with those obtained fully within a physical basis. This is a direct consequence of the formal equivalence between the Hamilton equations expressed in the physical basis, Eq. <ref>, and the reciprocal-space variant, Eq. <ref>.

For pristine lattices, electronic wavepackets tend to localize in low-energy regions of reciprocal space, allowing a reciprocal-space quantum basis to be truncated to those regions <cit.>. Moreover, since phonons act so as to absorb or emit lattice momentum during electron–phonon scattering events, a reciprocal-space classical basis can be truncated to select regions in conjunction with quantum basis truncations <cit.>. Such truncations cannot be performed in the physical basis, since electronic carriers and phonons are both delocalized over the entire lattice. Importantly, such truncations offer radical savings of the computational cost of simulating the electron–phonon scattering dynamics.



§ APPLICATION TO AN IMPURITY MODEL


While the reciprocal-space treatment of pristine lattices presented in Papers I <cit.> and II <cit.> provides a compelling demonstration of the utility of transformed MQC dynamics, it is the generalized equations of motion presented in the present Article that allows optimal basis transformations to be applied even when the lattice is disrupted. In such cases, lattice momentum is no longer a good quantum number, compromising the effectiveness of reciprocal-space MQC dynamics and prompting the need for alternative representations in order to enable effective basis truncations. To demonstrate this, we proceed to present results for a lattice involving a single impurity.



 §.§ Model


As in Paper I <cit.>, we consider a single electronic carrier while representing the lattice by a tight-binding model. The associated purely-electronic quantum Hamiltonian is given by

    Ĥ_q = 
        -J ∑_n( ĉ_n+1^†ĉ_n + ĉ_n^†ĉ_n+1) - Δ ĉ_n̅^†ĉ_n̅.

Here,ĉ_n^†andĉ_nrepresent the annihilation and creation operators for the electronic carrier, respectively, associated with lattice siten, andJis the nearest-neighbor interaction term. Furthermore,n̅denotes the impurity site andΔis the associated energetic detuning relative to the other lattice sites. Periodic boundaries are imposed, so thatn=Ncorresponds ton=0. This Hamiltonian is solved for in order to obtain the purely-electronic eigenstates obeying,

    Ĥ_q|ϕ_i⟩ = E_i |ϕ_i⟩.

Here and throughout,iis used to label the purely-electronic eigenstates, which are not to be confused with the instantaneous eigenstates of the total quantum Hamiltonian (cf. Eq. <ref>), which are denoted byαinstead. While the instantaneous eigenstates depend parametrically on the classical coordinates, thereby attaining a time-dependence, the purely-electronic eigenstates form a time-independent, diabatic basis. As a convention,iis taken to run with increasing eigenenergy.

Shown in Fig. <ref> are the purely-electronic eigenfunction amplitudes within the physical basis for a pristine lattice withJ = 1.0andΔ= 0and for an impurity model with the detuning was adjusted toΔ= 2.0. Here and throughout, we take parameter values to be unitless. We note, however, that, when taking the thermal energy at room temperature (293 K) as a reference, a unit of energy amounts to 25 meV. As can be seen in Fig. <ref>, for the pristine lattice, the eigenstates assume the periodic oscillatory profiles indicative of Bloch states. For the impurity model, however, the lowest-energy eigenstate (withi = 1) is largely localized on the impurity site,n̅. All other eigenstates with oddivalues (i = 3, 5, 7, …) reproduce the Bloch-like states (save for a slight and arbitrary phase shift). For eveni, the eigenstates tend to approach the Bloch-like states with increasing values ofi, although maintaining a deviation close ton̅. In the course of nonadiabatic dynamics, eigenstate population is expected to funnel from the delocalized Bloch-like states towards the localized state at lowest energy. It is this lowest-energy state that is challenging to resolve within a reciprocal space basis, as a localized state involves contributions from many reciprocal-space basis states due to the underlying Fourier relationship.

In order to find an alternative basis within which to effectively represent the nonadiabatic dynamics of the impurity model, we first observe that the previously-applied reciprocal-space representation <cit.> effectively invokes the purely-electronic eigenstates of the pristine lattice, i.e., the eigenstates ofĤ_qshown in Fig. <ref>. It should be noted that the eigenstates of the pristine lattice shown in this figure were solved for through a real-valued eigenvalue decomposition of the real-valued HamiltonianĤ_q. This yields cosine and sine solutions rather than the complex exponents commonly appearing in the Bloch formalism. However, one can trivially transform between both solutions by taking symmetric and antisymmetric combinations of degenerate eigenstates. Under Bloch-like solutions, truncations of the reciprocal-space basis were previously invoked by introducing a wavevector cutoff, denotedk_0, such that basis states having|k|>k_0were excluded <cit.>. For the real-valued solutions shown in Fig. <ref>, it is more appropriate to introduce an energy cutoff,E_c, such that basis states are excluded havingE_i>E_c. Importantly, this effectively yields identical basis truncations, since the Bloch-state energies increase monotonically withk(provided thatJ>0).

While the reciprocal-space basis may not offer much benefit to describing the impurity model, the notion of using the purely-electronic eigenbasis, which differs between the pristine lattice and the impurity model, is an interesting choice to consider. Indeed, the purely-electronic eigensolutions of the impurity model tend to simultaneously capture the extendedness and localization necessary to describe the nonadiabatic scattering onto the impurity, as shown in Fig. <ref>. Moreover, the underlying transformation also captures the localized nature of the nuclear vibration responsible for self-trapping at the impurity site, while simultaneously representing the (approximately) momentum-carrying phonons throughout the rest of the lattice. For that reason, in describing the impurity model, we will resort to the purely-electronic eigenstates as a basis for both the quantum and classical subsystems.

We should stress that this choice of basis is heuristic, and that superior bases are likely to exist. Moreover, there is no need to keep with the same basis when describing the quantum and classical subsystems. Identifying the theoretically-optimal bases is not a trivial task, however, and we reserve a thorough exploration of this topic for a follow-up study. Instead, by adopting the purely-electronic eigenbasis, we will present a proof-of-principles of the general applicability of transformed MQC dynamics, and the possibilities for basis truncations it affords. Within the purely-electronic eigenbasis representation for the impurity model, we invoke basis truncations similarly to that applied to the real-valued solutions of the pristine lattice, by introducing an energy cutoff,E_c.

As an initial condition of the quantum system, we consider the single-carrier excitation with zero lattice momentum, given by

    |Ψ_0⟩ = |k=0⟩ = 1/√(N)∑_n|n⟩.

This initial condition is representative of a tightly-bound electron–hole pair (Frenkel exciton) produced upon impulsive optical excitation (under the long-wavelength approximation).

We invoke the Holstein model in order to account for the nuclear modes driving the nonadiabatic dynamics. Accordingly, the classical subsystem is taken to consist of identical harmonic modes, such thatV_c({q_n})is given by Eq. <ref> andz_nis given by Eq. <ref>.
The operator governing the quantum–classical interactions is then given by

    Ĥ_q-c = g √(2ω^3)∑_nĉ_n^†ĉ_n q_n.

In what follows, we adopt the parameters from Fig. <ref>, and additionally setω= 0.2andg = 1.0, while adjusting the total number of lattice sites toN=30. The classical coordinatesq_nandp_nare initially and independently drawn from a Boltzmann distribution <cit.> at a temperatureT = 1.0.
We reiterate that, when taking the thermal energy at room temperature (293 K) as a reference, a unit of energy amounts to 25 meV. By the same token, a unit of time amounts to 164 fs.



 §.§ Results


In the following, we will present FSSH calculations for both the case of a pristine lattice, by settingΔ= 0, and that of a single impurity model, by settingΔ= 2.0. With the pristine lattice, we will be revisiting a system that was addressed in Papers I <cit.> and II <cit.>, but with adjusted parameters. The FSSH calculations presented were performed within the reciprocal-space basis as well as in the purely-electronic eigenbasis. For the former, real-valued eigenvector projections (cf. Sec. <ref>), necessary for momentum rescalings, were obtained following the procedure outlined in Paper II <cit.>. For the purely-electronic eigenbasis, all relevant eigenvectors are real-valued by construction, and were taken as is. Both treatments are fully consistent with our recently-proposed gauge-invariant generalization <cit.>.

The electronic populations were evaluated by first constructing the total electronic density matrix within the instantaneous eigenbasis of the total Hamiltonian operator,Ĥ_q + Ĥ_q-c. As is commonly done <cit.>, the diagonal elements of this density matrix were based on the active surfaces, while the offdiagonal elements were based on the electronic wavefunction coefficients, i.e.,

    ρ_αβ = δ_αβδ_α a + (1 - δ_αβ) A_α^* A_β,

whereδis the Kronecker delta function. This density matrix was then transformed <cit.> to the physical basis, reciprocal-space basis, and purely-electronic eigenbasis, respectively, after which corresponding populations were obtained throughP_n = ρ_nn,P_k = ρ_kk, andP_i = ρ_ii.



Shown in Fig. <ref> are results for the pristine lattice. Here, only reciprocal-space FSSH was applied, recognizing the formal equivalence with the purely-electronic eigenbasis in this limit (as discussed in Sec. <ref>). Fig. <ref> (a) presents time-dependent reciprocal-space electronic populations,P_k, obtained without any basis truncation imposed. As seen here, and as discussed in Papers I <cit.> and II <cit.>, a scattering of the electronic carrier out of thek=0initial state is observed.

Fig. <ref> (b) depicts the time-dependent zero-momentum electronic population,P_k=0, resulting from untruncated calculations, together withP_k=0calculated under increasing basis truncations. This population is seen to remain invariant under truncations of up to 50% of the reciprocal-space basis, which is consistent with the findings of Paper II <cit.>, where FSSH was applied to a pristine lattice under Holstein-type coupling between the carrier and the nuclear modes. This is rationalized by the carrier being both initiated and thermally-biased towardsk=0, with high-momentum basis states providing only small contributions to the nonadiabatic dynamics.



A similar analysis is presented in Fig. <ref>, but for the impurity model subjected to FSSH within reciprocal space and within the purely-electronic eigenbasis. Shown as a reference in Fig. <ref> (a) are time-dependent populationsP_n, obtained through FSSH in the physical basis. Here, thek = 0initial state is seen to be fully delocalized, as expected from Eq. <ref>, while nonadiabatic dynamics funnels the carrier into the impurity withn = n̅. A schematic illustration of this process is depicted in Fig. <ref> (b). Shown in Fig. <ref> (c) and (e) areP_kandP_i, respectively, as obtained through FSSH within the reciprocal-space basis and within the purely-electronic eigenbasis. While the dynamics ofP_kappears to be quite similar to that shown for the pristine lattice in Fig. <ref> (a),P_iclearly exhibits the gradual trapping of excitation in the lowest-energy eigenstate withi = 1.

We will proceed to assess to what extent truncations of the reciprocal-space basis and the purely-electronic eigenbasis allows capturing of both the delocalization of the initial excitation as well as the transient localization at the impurity site. To this end, we evaluateP_k=0as well asP_n̅. Fig. <ref> (d) and (f) show both populations as obtained through FSSH under truncations of the reciprocal-space basis and of the purely-electronic basis, respectively. As can be seen here, both truncated bases perform well in describingP_k=0. A markedly-worse performance is found forP_n̅under reciprocal-space basis truncations, which is a direct consequence of reciprocal space providing a poor representation of local phenomena. Due to the Fourier relationship, accounting for the local population at siten̅requires the complete reciprocal-space basis to be included, and truncations of this basis yield a proportional drop in this site population.
    
The purely-electronic eigenbasis, however, is found to perform remarkably well in describingP_n̅under truncations of up to 50%. This is the result of the lowest-energy eigenstate capturing most of the local population at siten̅. As such, the truncated purely-electronic eigenbasis is shown to outperform the truncated reciprocal-space basis in governing the funneling of electronic population into the impurity, as we anticipated based on Fig. <ref>.



§ CONCLUSIONS AND OUTLOOK


In summary, we have derived a formulation of MQC dynamics within arbitrary bases for the classical and quantum coordinates. This allows any given system to be optimally treated by finding preferred bases which can be efficiently truncated while retaining good agreement with untruncated calculations. Such is demonstrated by our application of FSSH to a system involving an electronic carrier scattering onto a single impurity in an otherwise pristine lattice while interacting with phonons. The funneling from an initial, delocalized state towards a low-energy state localized at the impurity was faithfully captured when transforming the quantum and classical coordinates into the purely-electronic eigenbasis, followed by significant basis truncations. A similar level of truncation within the reciprocal-space basis instead yielded significant distortions of the impurity population.

The encouraging results obtained for the impurity model offer promising prospects for the realistic application of MQC dynamics, and in particular FSSH, to materials involving disruptions of the crystal lattice. In recent years, FSSH has found applications to materials involving an impurity <cit.>, defects <cit.>, and wrinkling <cit.>, as well as materials where the lattice symmetry was disrupted by nearby molecules <cit.>. While no basis transformations were invoked in these studies, our present analysis suggest them to offer significant gains in performance and scalability when combined with efficient basis truncations.
    
We re-emphasize that the purely-electronic eigenbasis adopted for the impurity model was chosen based entirely on heuristic arguments and, in spite of its good performance, may not be the theoretically-optimal choice. Finding the optimal choice for a given system is nontrivial, and as such we reserve such inquiry for future research. Here, it should be noted that the bases adopted for the quantum and classical subsystems need not be the same, as was adopted in the present Article. Rather both bases can be independently optimized, resorting to the full special unitary group for each <cit.>.

Importantly, the transformed equations of motion presented in this Article are applicable to the full suite of MQC methods, beyond FSSH. An example of particular interest is Ehrenfest dynamics, which is unrivaled in terms of simplicity and computational affordability. These favorable attributes notwithstanding, Ehrenfest dynamics is known to suffer from over-thermalization in the asymptotic time limit <cit.>. In addition to contributing to the method's inaccuracy, this deteriorates the effectiveness of basis truncations, as quantum excitations tend to delocalize over an excess number of basis states <cit.>. We should emphasize that the intrinsic accuracy of MQC modeling is not addressed in the present Article. In Paper II <cit.>, we found the results generated by FSSH for a lattice model to show good agreement with full-quantum reference calculations, yet much remains to be learned about the performance of FSSH in materials modeling. In that regard, a recently-proposed coherent generalization of FSSH is noteworthy <cit.>, which may offer particular advantages in capturing the coherent dynamics prevalent in materials.
    
Lastly, we note that the transformed classical dynamics derived in Sec. <ref> is completely general. Its ability to incorporate mode anharmonicities, although not taken advantage of in our current demonstration of the impurity model, enables broad application to molecular dynamics simulations, with or without accompanying quantum modeling. Subjecting a given classical system to unitary basis transformations, the same way quantum-mechanical systems are treated, may provide a straightforward route to optimally representing its dynamics by a truncated set of coordinates.


This material is based upon work supported by the National Science Foundation under Grant No. 2145433. K.M. gratefully acknowledges support from the Mark A. Ratner Postdoctoral Fellowship and the Northwestern University International Institute for Nanotechnology (IIN).