Universal Steps in Quantum Dynamics with Time-Dependent Potential-Energy Surfaces: Beyond the Born-Oppenheimer Picture

It was recently shown [G. Albareda, et al., Phys. Rev. Lett. 113, 083003 (2014)] that within the conditional decomposition approach to the coupled electron-nuclear dynamics, the electron-nuclear wave function can be exactly decomposed into an ensemble of nuclear wavepackets effectively governed by nuclear conditional time-dependent potential-energy surfaces (C-TDPESs). Employing a one-dimensional model system we show that for strong nonadiabatic couplings the nuclear C-TDPESs exhibit steps that bridge piecewise adiabatic Born-Oppenheimer PESs. The nature of these steps is identified as an effect of electron-nuclear correlation. Furthermore, a direct comparison with similar discontinuities recently reported in the context of the exact factorization framework allows us to draw conclusions about the universality of these discontinuities, viz. they are inherent to all nonadiabatic nuclear dynamics approaches based on (exact) time-dependent potential energy surfaces.


I. INTRODUCTION
The description of the correlated electron-nuclear dynamics remains a formidable challenge in condensedmatter physics and theoretical chemistry [1][2][3][4][5][6]. Relying on the Born-Huang expansion of the full molecular wave function, a majority of approaches that provide a numerically accurate description of the so-called "nonadiabatic" processes require the propagation of a set of many-body nuclear wavepackets on a coupled set of Born-Oppenheimer potential-energy surfaces (BOPESs) [7][8][9][10][11][12][13][14]. Whenever electron-nuclear coherence effects are unimportant, efficient mixed quantum-classical propagation techniques can be used [15][16][17][18][19][20][21]. However, in order to account for strong correlations, the access to quantum features of the nuclear motion such as wavepacket spreading, interferences, tunneling or splitting is crucial. Hence, a reliable description of molecular dynamics becomes very expensive due to the calculation of the full set of BOPESs, a (time-independent) problem that scales exponentially with both the nuclear and electronic degrees of freedom [22][23][24].
Towards an alternative description of the coupled electron-nuclear dynamics beyond the Born-Oppenheimer (BO) picture (i.e. the Born-Huang ex- * albareda@ub.edu † aliabedik@gmail.com ‡ ita@zurich.ibm.com § angel.rubio@mpsd.mpg.de pansion of the molecular wave function), two different formally exact frameworks, viz. the exact factorization (EF) [25] and the conditional decomposition (CD) [26], have been recently proposed. These two alternative approaches are to be added to the already existing frameworks that avoid the BO picture of nonadiabatic dynamics [27][28][29]. In the EF and the CD approaches the nuclear dynamics takes place on time-dependent potential energy surfaces. This is in contrast with the BO picture where the nuclear wavepacket, with contributions on several BOPESs, undergoes transitions in the regions of strong nonadiabatic coupling (NAC) (see Fig. 1). The EF and the CD approaches are, however, grounded on very different mathematical frameworks. In the EF approach, the electron-nuclear wave function is written as a direct product of electronic and nuclear wave functions. The resulting nuclear equations of motion arise from Frenkel's stationary action principle and do involve the integration of the electronic degrees of freedom [25]. Alternatively, in the CD approach electrons and nuclei are treated on an equal mathematical footing and hence their description is kept at the full configuration level (see Fig. 1). After a partial time-dependent coordinate transformation is applied to the time-dependent Schrödinger equation (TDSE), the nuclear (electronic) coordinates are propagated along with the electronic (nuclear) wave function such that the new coordinates remain located where the full molecular wave function has a significant amplitude [26]. A major distinctive characteristic of the CD approach is the fact that it avoids the integration of the electronic arXiv:1512.08531v3 [physics.chem-ph] 29 Nov 2016 degrees of freedom. It constitutes in this respect a safe starting point to draw new connections between other formally exact frameworks. In this work we aim at investigating the generic features of the exact conditional (nuclear) time-dependent potential-energy surfaces (C-TDPESs) that arise in the CD approach in the presence of strong non-adiabatic couplings. This is a necessary study towards the use of the CD method in real scenarios and will serve as the basis for future development of practical approximations to the nuclear C-TDPESs. A major result is that, besides the radically different natures of the CD and EF frameworks, the exact C-TDPESs exhibit discontinuous steps connecting different BOPESs analogous to a paradigmatic feature of the effective timedependent potential that governs the nuclear dynamics in the EF framework [25]. By establishing a formal connection between the CD and the EF frameworks we will elaborate on the universality of these discontinuities. We will be showing that the paradigm shift associated with the transition from the many BOPESs and NACs to the (single) time-dependent potential energy surface entails a discontinuous behavior of the resulting effective potentials during the branching of the nuclear probability density in nonadiabatic electronic transitions.
This paper is organized as follows. In Section II, we introduce the CD approach to electron-nuclear correlated dynamics and describe its main mathematical ingredients. In Section III we will provide some physical intuition on the role played by the C-TDPESs in connection with the BOPESs of the BO picture. Numerical simulations are thereafter, in Section IV, used to investigate the main features of the exact C-TDPESs in the presence of strong nonadiabatic couplings. In view of the results, in Section V we will establish a formal connection between the CD and the EF approaches and elaborate on the universality of the discontinuous behavior of the resulting effective potentials in fully time-dependent approaches to electron-nuclear coupled dynamics. Conclusions will be drawn in Section VI.

II. THE CONDITIONAL DECOMPOSITION
Throughout this manuscript we use atomic units, and electronic and nuclear coordinates are collectively denoted by r = {r 1 , .., r Ne } and R = {R 1 , .., R Nn }, where N e and N n are the total number of electrons and nuclei. We have recently proved [26] that the full (nonrelativistic) electron-nuclear wave function Ψ(r, R, t) can be exactly decomposed into an ensemble of conditional nuclear wave functions defined as: Provided that the electronic trajectories {r α (t) ≡ r 1,α (t), .., r Ne,α (t)} explore the electronic support of |Ψ| 2 at any time [30], these conditional wave functions can be always used to reconstruct the full wave function [31] (or equivalently any observable). Notice, however, that since the initial conditions of a trajectory-based simulation can be generated with importance-sampling techniques, conditional decompositions can be exploited to circumvent the problem of storing and propagating a many-particle wave function [26,[32][33][34][35][36]. The conditional wave function concept is illustrated in Fig. 2 for a simple twodimensional (one electron plus one nucleus) system. Throughout this work we will use quantum trajectories r ξ,α (t) and R ν,α (t), that are defined through [37]: where S is the phase of the full electron-nuclear wave function Ψ = |Ψ|e iS . Note that the choice of quantum (hydrodynamic) trajectories is not mandatory. Alternatively, trajectory-based Monte Carlo or importancesampling techniques can be used, provided that they sample the (time-dependent) quantum-probability density. It is however useful for the purpose of this work to choose Bohmian trajectories because they do sample the support of the full wave function by definition [37,38].
For the sake of simplicity, hereafter we will omit the explicit dependence of the trajectories on time, i.e. {r α } ≡ {r α (t)}.
In the absence of time-dependent external vector fields, the system is described by the Hamiltonian: whereT e = − Ne ξ=1 ∇ 2 ξ /2 andT n = − Nn ν=1 ∇ 2 ν /2M ν are the electronic and nuclear kinetic energy operators, and W (r, R, t) = W ee (r) + W nn (R) + W en (r, R) denotes the Coulombic interactions. The wave functions ψ α (R, t) then obey the following equations of motion [26]: where the effective potentials, are named conditional time-dependent potential-energy surfaces (C-TDPESs). Each C-TDPES consists of three terms: the Coulombic potential, the kinetic correlation potential: and a pure advective potential, FIG. 1. Schematic representation of the four (exact) approaches to correlated electron-nuclear dynamics that will be discussed in this work. We will be considering the following representations of the molecular wave function: (upper-right panel) full electron-nuclear wave function, (upper-left panel) Born-Huang expansion of the electron-nuclear wave function, (bottom-left) exact factorization of the electron-nuclear wave function, and (bottom-right) conditional decomposition of the electron-nuclear wave function. These four pictures are respectively referred to as Full Picture, BO Picture, EF Picture, and CD Picture. The nuclear equations of motion that arise in each representation are mathematically very different. In particular, we distinguish between approaches where the nuclear dynamics occurs in a reduced nuclear subspace (viz. BO and EF Pictures) and approaches where it occurs in the full electron-nuclear Hilbert space (viz. Full and CD Pictures). Further, we distinguish between approaches where the (effective) nuclear potential-energy surfaces are time-independent (viz. Full and BO Pictures) and approaches where these are time-dependent (viz. EF and CD Pictures). For the sake of simplicity, we omitted any external field in the scheme.
Here we emphasis that unlike the parametrized nuclear static potentials that arise in the BO picture or the effective exact time-dependent potential in the EF approach to electron-nuclear dynamics, the C-TDPESs of the CD method in Eqs. (7)-(9) do not involve the tracing-out of the electronic degrees of freedom. The C-TDPESs do live in the full electron-nuclear configuration space, and therefore they cannot be defined (strictly speaking) as potential-energy surfaces. The C-TDPESs consists of three terms, viz. the bear Coulomb potential plus two additional fictitious potentials, i.e. the kinetic and the advective correlations. The advective correlation potential represents the rate of change of the nuclear conditional wave functions due to the time-dependence of the new electronic coordinates (the Lagrangian trajectories). Instead, the kinetic correlation potential defines a viscous/drag-like quantum effect of the electrons on the nuclear degrees of freedom.
As shown in Eq. (5), the evolution of each nuclear wave function ψ α (R, t) is, in general, non-unitary due to the complex nature of the kinetic and advective potentials. Therefore, the conditional wave functions ψ α (R, t) do not obey a normalization condition individually, but rather a collective normalization constraint given by: where the transformation D has been defined as [39]: when ∞ α=1 δ(a α − a) = 0, and it is zero otherwise. Equation (10) states that advective and kinetic correlations alone can not be a net source or sink of probability density when integrated over a closed or infinite domain.

FIG. 2.
Schematic representation of the Conditional Decomposition approach to molecular dynamics for a 2-dimensional system. The full nuclear probability-density |Ψ(r, R, t)| 2 is plotted at two different times t0 and t f , together with a conditional amplitude |ψα(R, t)| 2 (in red) for a particular trajectory {rα(t), Rα(t)}. Black arrows denote the velocity field {ṙα(t),Ṙα(t)}, and contour plots of the full nuclear wave function are also shown for clarity.
Hence, only after renormalization, the conditional amplitudes |ψ α | 2 can be associated to conditional probabilities, i.e.: Notice that this is different, e.g., from the effective (norm preserving) nuclear wave function considered in the EF framework [40]. The derivation of the CD mathematical framework corresponds to the transformation of the many-body TDSE to the partially co-moving frame in which the electronic coordinates move attached to the electronic flow and the nuclear coordinates are kept in the original inertial frame. Within the new coordinates, the electronic convective motion is described by a set of trajectories of infinitesimal fluid elements (Lagrangian trajectories), while the nuclear motion is determined by the evolution of the (nuclear) conditional wave functions in an Eulerian frame. The purpose of this partial time-dependent coordinate transformation is to propagate the electronic coordinates along with the electronic flow such that they remain located where the full molecular wave function has a significant amplitude. In this way, the amount of "dead wood", as it is commonly called in the context of configuration interaction expansions [41], can be significantly reduced. Furthermore, the possibilities of the core ideas behind the CD rely on the fact that the underlying mathematical structure admits approximations which substantially reduce the computational complexity without deteriorating the quality of results. This has been demonstrated for example in Refs. [26,32], where a very simple (zeroorder) propagation scheme based on the CD approach was able to accurately describe not only non-adiabatic electronic transitions, but also intricate nuclear quantum effects.

III. ON THE ROLE OF THE C-TDPES
To provide some physical intuition on the role played by the C-TDPESs, we further decompose V α in Eq. (6) into real and imaginary parts, i.e., where and are respectively the real and imaginary parts of the kinetic and advective correlation potentials. In Eq. (14) we have respectively defined the ξ-th components of the so-called electronic quantum potential [37,38], (18) and the current probability density j ξ = |Ψ| 2 ∇ ξ S. While the classical kinetic potential, (∇ ξ S) 2 /2, is in general a smooth function of the nuclear coordinates, we expect the contributions, Q e ξ,α and K α , to be rather "discontinuous" because of their dependence on the inverse of the probability density. Notice that these last two quantities are implicit functions of the electronic coordinates. Specifically, the electronic quantum potential Q e ξ,α accounts for changes of the curvature of the full probability density along the electronic coordinates (as a function of the nuclear positions). Alternatively, K α describes the dispersion (e.g. spreading or splitting) of the full wave function in the electronic direction. Both contributions become large in the vicinity of a node, and thus, they will be important whenever the full probability density splits apart along the electronic coordinates. The advective correlations in Eq. (16) are weighted by the electronic velocities,ṙ ξ,α , and hence, they will be significant during a fast reconfiguration of the electronic degrees of freedom.
As it will be discussed later, the contribution of the kinetic correlation potential to the full C-TDPESs is much larger than the advective correlation potential in general situations. In order to analyse the role played by the kinetic potential in the evolution of the nuclear dynamics we first discuss the connection between the C-TDPESs and the standard Born-Oppenhaimer (BO) picture. The BO electronic states, Φ If only one BO state, e.g. Φ , then the kinetic correlation potential reads: In the absence of advective correlations, Eq. (20) implies that the C-TDPESs (independently of the electronic trajectory r α ) are all real, time-independent, and equal to BOPES The situation is, however, not so trivial as we move to a nonadiabatic scenario. Consider for instance the case where there are only two BO states involved in the dynamics, i.e. Ψ(r, R, t) = C 1 (R, t)Φ R (r). Then the kinetic correlation potential reads: where we have defined the relation and By comparing Eqs. (14) with (23) and (15) with (24) we can write the following equalities for a two level system (in the absence of advective correlations): and To better understand the meaning of these quantities, in the following we evaluate and discuss three particularly relevant limits of Eqs. (25) and (26): BO (R), Case 2: |C 1 | = 1 (which implies |C 2 | = 0 and also ∇ ξ S = 0) BO (R), Case 3: ∆ 12 → 0 (which implies ∇ ξ S = 0) Case 1 and Case 2 represent the limits in which one of the BO states is dominant at a certain time and in a certain region in R-space. In these two cases, the imaginary part of the kinetic correlation potential is zero and the real part is identical with one of the BOPESs. Departing from these limiting cases, the kinetic potential becomes complex where its real part represents the bridge between the BO-like pieces of the potential. Case 3, shows that when the BOPESs come close to each others, i.e., nearby conical intersections or avoided crossings, the kinetic potential is real and it coincides with the BOPESs. Therefore, in the adiabatic limit as well as in the extreme nonadiabatic limit the kinetic correlation potential is real.

IV. NUMERICAL RESULTS
To study numerically the features of the exact C-TDPESs during nonadiabatic processes in the CD approach, we employ the model introduced by Shin and Metiu [42]. This model has become a standard one (just as Tully's models but in full configuration space) as it is simple and provides all the key ingredients to test new approaches to nonadiabatic dynamics. The model is very flexible and accepts an infinite number of parameter configurations that give rise to a number of interesting situations where electron-nuclear correlations play a crucial role.
The Shin-Metiu model consists of three ions and a single electron. Two ions are fixed at a distance L = 19.0a 0 , and the third ion and the electron are free to move in one dimension along the line joining the fixed ions (see Fig. 3). The Hamiltonian for this system readŝ where erf () represents the error function, the symbols r and R are replaced by r and R, and the coordinates of the electron and the movable nucleus are measured from the center of the two fixed ions. In this work, we choose the remaining parameters to be M = 1836a.u. and R f = 5.0a 0 , R l = 4.0a 0 , and R r = 3.1a 0 such that the first BOPES, BO , is strongly coupled to the second, BO around R ac = −2a 0 . The coupling to the remaining BOPESs is negligible.

FIG. 3.
Schematic representation of the Shin-Metiu model [42]. Two ions are fixed (in black) and a third one (in red) and an electron (in blue) are free to move in one dimension.
Let us remark that by taking such a confined system, one is exposed to strong electron-nuclear correlations that usually do not show up in scattering problems with a continuum of states. It is specially under this highly confined conditions that the existing (nonadiabatic) theories fail to describe the correlated motion of electrons and nuclei. Notice also that the parameters used in this work differ significantly from those used in Ref. [26]. More specifically, in Ref. [26] we described nuclear quantum effects such as interferences or tunneling, while in this work we delve into the splitting of the nuclear probability density.
The time-independent electronic problem is fully characterized through the so-called electronic Hamiltonian, i.e. (T e + W )Φ (j) R (r). We suppose the system to be initially excited to (2) BO and the initial nuclear wave function to be a Gaussian wavepacket with σ = 1/ √ 2.85, centered at R = −4.0a 0 . On Fig. 4 we show the first three BOPESs together with the evolution of the adiabatic populations (in the inset). Figure 5 shows the first two Born-Oppenheimer states Φ (1)

R (r) and Φ
(2) R (r) along with the evolution of three selected trajectories {r α , R α } labeled α = 1, 2, 3 that will be used to analyse the C-TDPESs. In Fig. 6 we present four time snapshots containing relevant information about the C-TDPESs as well as the conditional nuclear wave functions. For the sake of clarity we also define approximated real and imaginary components of C-TDPESs respectively as V α,app = Q e α + W α and V α,app = K α . At the initial time, due to the choice of Ψ(r, R, t 0 ), the C-TDPESs are real, α-independent, and by construction all equal to the first excited BOPES (2) BO : R (r) when ∂ r S = 0. As it is made clear in the right panels of Fig. 5, at t = 13.44fs, the trajectories (r 2 , R 2 ) and (r 3 , R 3 ), associated respectively with the conditional wave functions ψ 2 and ψ 3 , are running straight (i.e.ṙ 2,3 ≈ 0) from the support of Φ   R (r). Therefore, as can be seen from Fig. 6 (bottom panels), the C-TDPESs resemble diabatic potential-energy surfaces, coinciding with (2) BO for R < R ac and smoothly going through the avoided crossing region to follow (1) BO for R > R ac . The trajectory (r 1 , R 1 ), while staying most of the time in the support of Φ (2) R (r), is now tunneling from one Coulomb potential well to the other. The gradient of the phase, (∇ r S)| r1 , that appears in both the kinetic and advective correlation terms (also implicitly in the electronic velocityṙ 1 ), is now very large due to the fast reconfiguration of the electronic degrees of freedom. Hence V 1,app and V 1 show important differences. At the later time t = 19.32fs, the advective correlation and classical kinetic terms become once more negligible, and therefore V α ≈ V α,app + V α,app (see Fig. 7). While (r 1 , R 1 ) stays preferentially in the support of Φ   R (r). This mixture leads to the formation of the step observed for V 1 around R = 1.5a 0 , and two wiggles accompanying the C-TDPESs V 2 and V 3 . All these features indicate the nonadiabatic nature of the conditional wave functions and can be directly associated with the formation of nodes in the full probability density. Finally, at time t = 31.5 fs, the full molecular wave function has been split both along the electronic and nuclear directions. While the three conditional wave functions , ψ 1 , ψ 2 and ψ 3 , still embody contributions from Φ (1) R (r) and Φ (2) R (r), these contributions are now very well separated along the nuclear coordinates with a minimum at around R sp = 4a 0 (see also Fig. 5). For nuclear coordinates less than R sp , the support of the full probability density perfectly fits in the support of the first excited state Φ (2) R (r), while for R > R sp it mainly falls in the ground state Φ (1) R (r). As a direct consequence, the quantum electronic potential acquires a discontinuity, i.e. Q e α (R < R sp ) ≈T e Φ (2) R (r). The real parts of the C-TDPESs are therefore piecewise connecting adiabatic surfaces, i.e. V α ≈ V α,app ≈ (1) BO for R < R sp . In the transition between these two regions (at around R sp ), the conditional wave functions ψ 1 , ψ 2 and ψ 3 become linear combinations of the lowest two adiabatic states. As a result, the C-TDPESs exhibit abrupt peaks, for both real and imaginary parts, mainly originating from the potentials Q e α and K α .

V. DISCUSSION
In view of the above discussion, the characteristic features of an ensemble of C-TDPESs of the CD framework during a nonadiabatic transition are very similar to that of the gauge invariant part of the exact TD-PES [40,43,44] within the EF framework: while in the vicinity of the avoided crossing both potentials follow a diabatic shape, far from the region of avoided crossings (when the nonadiabatic couplings are negligible) they are piecewise identical with two different adiabatic BOPESs with nearly discontinuous steps in between. By discussing the link between the potentials of the EF and CD approaches, in the remaining part of the manuscript we elaborate on the origin of the discontinuities, wich seem to be a universal feature of the time-dependent potentials driving the nuclear dynamics within formally exact approaches to correlated electron-nuclear dynamics.
The CD and the EF approaches are essentially different exact mathematical formulations of the electronnuclear coupled problem. While the CD approach is a "trajectory-guided" method that treats electrons and nuclei on an equal mathematical footing [26], the main equations of the EF approach arise from Frenkel's stationary action principle and hence involve the tracing out of either the electronic or nuclear degrees of freedom [25]. Therefore, while the CD and the EF approaches can be FIG. 6. The nuclear dynamics before the splitting of the probability density is mainly characterized by C-TDPESs smoothly connecting (2) BO and (1) BO . PANEL a: Contour lines of the full electron-nuclear probability density (in white) together with three conditional wave functions , ψ1(R, t), ψ2(R, t), and ψ3(R, t) (black lines), defined along with trajectories {rα=1,2,3, Rα=1,2,3} (in cyan). In the background (in copper color scale), full (2D) dependence of the electron-nuclear potential energy W (r, R). PANELS 1,2,3 (bottom): Conditional amplitudes |ψα=1,2,3(R, t)| 2 (filled in gray) along with the first three BOPESs (in black) and the real part of the C-TDPESs V α = Wα +K α +A α (in red). For comparison we also show the potential V α,app = Q e α +Wα (blue circles). PANELS 1,2,3 (middle): Electronic quantum potential Q e α . PANELS 1,2,3 (top): Imaginary part of the C-TDPESs V α = K α + A α (red line) along with the approximated potential V α,app = K α (blue circles). formally linked together (as it will be briefly outlined in the following), there is no straightforward connection between the conditional wave function ψ α (R, t) of the CD and the effective wave functions of the EF approach.
The exact TDPES that governs the nuclear dynamics together with a vector potential within the EF framework consists of two parts [40,43,44], the gauge invariant term, gi (R, t), and a gauge dependent term (see Appendix A for a full derivation of the connection between the CD and the EF approaches). By virtue of the EF theorem, the conditional nuclear wave function can be written as a direct product of electronic and nuclear wave functions , i.e., ψ α (R, t) = Φ R (r α , t)χ(R, t). Hence, gi (R, t) can be expressed in terms of C-TDPESs and conditional nuclear wave functions , ψ α (R, t), as: As the second term, O(M −1 n ), is mostly negligible due to its dependence on the inverse of the nuclear mass, the first term on the r.h.s of Eq. (32) establishes a direct connection between the gauge invariant part of the TDPES and the C-TDPESs: gi (R, t) is an ensemble average of C-TDPESs, V α , which are integrated along the conditional amplitudes |ψ α | 2 and are weighted afterwards by the nuclear probability density, D r [|ψ α | 2 ]dr = |Ψ| 2 dr. As D r [|ψ α | 2 V α ]dr is zero (see in Appendix A the full derivation), Eq. (32) can be written as: Furthermore, the electronic quantum potentials are the main source of the discontinuities of V α as discussed before. Hence, the steps in gi (R, t) can be viewed as an ensemble average of the discontinuities of the electronic quantum potential.

VI. CONCLUSION
The branching of nuclear wavepackets at nonadiabatic transitions is a key nuclear quantum effect in photophysics and it is traditionally described within the Born-Oppenheimer picture where the dynamics is described by means of population transfer among different BOPESs. Within the exact factorization and conditional decomposition frameworks the same physics can also be successfully described by a single effective time-dependent potential that drives the nuclear dynamics along adiabatic segments of BOPESs connected by sudden transitions or jumps. Within the CD framework, an ensemble of C-TDPESs governs the dynamics of the conditional nuclear wave functions and provide us with an alternative approach for the study of nonadiabatic processes. In particular, we demonstrated that within the EF framework the discontinuous steps that connects two adiabatic pieces of the TDPES can be reproduced by an ensemble average of discontinuities of the C-TDPESs.
Our study reveals that the shape and the position of the steps are mainly dictated (both in the CD and the EF approaches) by the electronic quantum potential, Q e α , which is implicit in the electronic kinetic energy operator of the full time-dependent Schrödinger equation. The discontinuities associated to the quantum electronic potential are inevitable when a closed time-dependent form for the electron-nuclear wave function is used (as in the EF, CD or any other equivalent representation) and give rise to a universal mechanism for the explanation of the branching of the nuclear probability density in the nonadiabatic regime (beyond the Born-Oppenheimer picture). In other words, the paradigm shift associated with the transition from the many BOPESs to the single timedependent potential entails a discontinuous behavior of the resulting effective time-dependent potentials. More-over, these discontinuities represent a clear signature of the ultimate quantum nature of electron-nuclear correlations; the effective (nuclear) potentials have a direct dependence on the full wave function through the action of the electronic quantum potential, which is ultimately the source of the discontinuities of the time-dependent surfaces within the CD and EF frameworks.
It is worth mentioning that other attempts to reduce the complexity of the time-dependent (manyelectron/ion) Schrödinger equation, such as in timedependent density-functional theory, also lead to the emergence of striking discontinuities in the effective timedependent potentials that govern the reduced variables of interest. This is the case, e.g., in the field of electronelectron correlated dynamics, where the nature of such discontinuities is currently under study [45][46][47].
Acknowledgements In order to establish a formal connection between the C-TDPESs and the TDPES concepts we first write the molecular wave function as a direct product of electronic and nuclear probability amplitudes, i.e. ψ α (R, t) = Φ R (r α , t)χ(R, t).