Effective nonlinear Ehrenfest hybrid quantum-classical dynamics

The definition of a consistent evolution equation for statistical hybrid quantum-classical systems is still an open problem. In this paper we analyze the case of Ehrenfest dynamics on systems defined by a probability density and identify the relations of the non-linearity of the dynamics with the obstructions to define a consistent dynamics for the first quantum moment of the distribution. This first quantum moment represents the physical states as a family of classically-parametrized density matrices $\hat \rho(\xi)$, for $\xi$ a classical point; and it is the most common representation of hybrid systems in the literature. Due to this obstruction, we consider higher order quantum moments, and argue that only a finite number of them are physically measurable. Because of this, we propose an effective solution for the hybrid dynamics problem based on approximating the distribution by those moments and representing the states by them.


Introduction
It is an accepted fact that the most accurate physical theories describing Nature are quantum. Nonetheless, it is also well known that ab-initio full quantum theories are not useful from the practical point of view, because of their complexity. One of the choices is then to approximate the full quantum model by a simpler one where as many degrees of freedom as possible are modelled as classical systems. The paradigmatic example of this situation is the model of a molecule. Based on the application we have two options: either to consider all the degrees of freedom as classical (this is the case in most frequent molecular dynamics simulation methods), or to create a hybrid quantum-classical model where only the valence electrons, responsible of the chemical properties of the molecule, are modelled as quantum objects, the rest of degrees of freedom being represented by classical variables. In this second case, we define a more accurate description of the chemical properties, with a much simpler model than the full-quantum one. Nonetheless, the definition of an accurate hybrid dynamical system is not an easy task. There exist several approaches to define hybrid dynamical models for molecular systems (see, for instance, [1] for a recent review). We can classify them based on different properties, but, attending only to the definition of the dynamics, we may consider, among others: • those based on hybrid dynamics on the space of hybrid states ( [2,3,4,5,6,7,8]), • those which are algorithmic (see [9,10]) • or others obtained as suitable limit equations of the full-quantum dynamics ( [11,12,13,14,15]. If we enlarge our scope and consider other hybrid dynamical models, we can also find those considering the problem of measurements of quantum systems with classical devices [16,17] and other type of systems (see [18,19,20,21] and references therein).
One of the most common choices is the Ehrenfest dynamics, where it is simple to track the accuracy of the hybrid model with respect to the original full quantum one within a semiclassical description [22]. Ehrenfest equations represent a dynamical system defined on the cartesian product of the classical and quantum phase spaces M C × M Q . This hybrid phase space contains the hybrid pure states, i.e., those states where the classical and the quantum degrees of freedom are completely determined. Our group proved [5] that this system admits a Hamiltonian description with a suitable hybrid Poisson bracket and a hybrid Hamiltonian function, which combines the classical and the quantum energies. Nonetheless, molecular systems do not usually have completely determined dynamical states, since initial conditions are impossible to fix. Therefore, a statistical description appears to be the most reasonable one. In [5] we also used the Hamiltonian structure of Ehrenfest equations to define a statistical model having Ehrenfest equations as dynamics for the microstates. In this statistical description, the state of the system is defined as a probability distribution on the hybrid phase-space following a Liouville equation. This is a consistent mechanism to define a statistical mechanical system associated to a well defined thermodynamics (see [23]).
In the last years, several relevant applications of this scheme have been published ( [6,7]) but also some limitations have appeared. Among them, the most relevant one may be the difficulty to write an entropy function and the corresponding notion of canonical ensemble [24] or, in more general terms, of an equilibrium thermodynamics with a well defined temperature.
Our goal in this paper is a careful analysis of the origin of these difficulties, and a discussion of possible solutions and their implications. From the mathematical point of view, the problems are associated with the incompatibility of the notion of hybrid state as a probability density on the phase space and the definition of hybrid entropy, which require an alternative notion of state. Instead of the probability density described above, entropy can be formulated by considering its first quantum moment, which defines the type of hybrid state which has been used extensively since the early eighties (see [25]) and is used in many of the references presented above. But the challenge is to define a consistent master equation for this object, since the Liouville equation of the probability density does not preserve the first moment of the density. We will see this property in detail, and define a series of equations for the different moments which are equivalent to the Liouville equation for the probability distribution. Nonetheless, depending on the physical situation, not all the moments are necessary to define the physically relevant averages and their evolutions. For those situations where only a finite number of them are physically relevant, just a finite system of equations characterize the evolution of the quantum moments and an effective nonlinear master equation can be defined. This is the main result of the paper.
It is important to remark that the properties of the Hamiltonian structure of the microstate dynamics, such as having two independent symplectic or Poisson structures (as Ehrenfest dynamics), is the true origin of these properties. The particular Hamiltonian function, and hence the detailed properties of the dynamics are not that relevant. Hence, most of the conclusions of our paper shall be true in any Hamiltonian system sharing independent structures for the classical and quantum subsystems.
The scheme of the paper is as follows. First, in Section 2 we will summarize the construction of the Hamiltonian structure for Ehrenfest equations and its statistical extension. We will discuss that its Hamiltonian nature implies that only a finite number of quantum moments of the hybrid distribution are physically measurable. Then, Section 3 analyzes the properties of the different definitions of the hybrid statistical state and justifies why Liouville equation can not be restricted to any of the quantum moments. The reason for that is the nonlinearity introduced in Ehrenfest equations by the classical degrees of freedom. Nonetheless, we can encode the full Liouville equation in a series of coupled differential equations of the full set of moments. But since only a finite number of moments are physically meaningful, we can approximate the behavior of Liouville equation by a finite number of those differential equations, suitably adapted. This is the main result of the paper. Finally, Section 4 summarizes the main contributions of the paper and discusses future lines of research based on them.
2 Ehrenfest dynamics: from pure states to distributions

Ehrenfest equations on pure states
Let us consider first the problem of Ehrenfest dynamics for pure states. We assume then that the state of the hybrid system is characterized by a pair (ξ, ρ ψ ), where ξ = (⃗ q, ⃗ p) ∈ M C is a point on the classical phase space and ρ ψ specifies the pure state of a quantum system by a one-dimensional projector determining a point in the complex projective space PH of a certain Hilbert space H. Ehrenfest equations define a set of coupled differential equations describing the evolution of these degrees of freedom asq whereĤ(ξ) represents the quantum energy as a self-adjoint family of operators on H parametrized by the classical degrees of freedom written in Darboux coordinates q k , p j for j, k = 1, · · · n for the 2n-dimensional symplectic manifold M C (see [26] for the definition of these concepts). The simplest image for this operator is the electronic Hamiltonian for the valence electrons of a molecular system which depends on the positions of the charges of the protons and inner electrons of the atoms which are modelled by the classical degrees of freedom. In [22] it is proved that this system of equations defines a dynamical system whose solutions approximate the solutions of the Schrödinger equation for a full-quantum molecular Hamiltonian, the accuracy depending on the ratio of the masses of the quantum and the classical particles and the ratio of the width of the nuclear wave-packet with respect to the natural length of the problem. The classical phase space M C is assumed to be endowed with a symplectic form ω C . Hence, its functions C ∞ (M C ) can be endowed with a canonical Poisson bracket {·, ·} C . On the other hand, the quantum degrees of freedom correspond to the points of the set of pure density matrices of a quantum system (which we will denote in the following as M Q ), which being a Kähler manifold, is also endowed with a canonical Poisson structure {·, ·} Q (for details see [27,28,29]). For the sake of simplicity we will consider that the quantum system is finite dimensional, and hence that M Q is a finite dimensional Kähler manifold. From the point of view of applications, this is not a serious constraint, since numerical simulation, which is usually the final goal of the hybrid model, requires a finite dimensional system.
In conclusion, we can endow the set of functions of the hybrid phase space M C ×M Q with a Poisson structure by combining the classical and the quantum brackets (for simplicity, we will fix ℏ = 1 in the following) obtained from the symplectic forms on M C and M Q : From the results in [5], it is immediate to prove that Ehrenfest equations can be given a Hamiltonian description with respect to the hybrid Poisson bracket {·, ·} H and the Hamiltonian function where we represent as H C (ξ) the general expression of the energy of the classical degrees of freedom (in Equations (1) it would reduce to just the classical kinetic term but more general situations with arbitrary classical potentials can also be considered); andÎ represents the identity operator. Thus, the integral curves of the Hamiltonian vector field correspond to the solutions of Ehrenfest equations. The dynamics preserves the purity of the quantum states but it is non-linear in general, because of the classical degrees of freedom (the term H C (ξ) in Equation (3)). For more details see [5,6].

The general case
From the definition of the hybrid Hamiltonian, we can consider the general notion of hybrid observable, representing the physical magnitudes written in terms of the classical and quantum degrees of freedom. Generally speaking, we can think at the hybrid observables as the tensor product of the C * -algebra of classical functions and the C * -algebra of quantum operators. For technical reasons, we will ask the classical functions to have compact support for integrals to be well defined. In practical applications, most frequently computer simulations, this does not seem to be a strong requirement. If we use the geometrical formulation described in [5], we can consider that the space of hybrid observables can be represented by the functions on M C × M Q which correspond to a pair of classical and quantum observables, i.e., functions of the form fÂ (ξ) := Tr ρ ψÂ (ξ) .
In [5] we proved that it is possible to define a statistical mechanical formalism for hybrid systems by introducing a measure in this set of observables. If we take a reference measure dµ QC on the hybrid phase space (for instance, the sympletic volume), the chosen measure can be written as a probability density function F QC (ξ, ρ ψ ) with respect to it, defining the average of a certain magnitudeÂ(ξ) as Naturally, the measure must be well defined and therefore the density must satisfy that F QC ≥ 0;

The first quantum moment
Notice, though, that in physical terms, only the first quantum moment of F QC is necessary to determine the average value of a hybrid observable of the form of Equation (5). Indeed, the first quantum moment of F QC determines a family of quantum operators indexed by the classical variables ξ in the form: where dµ Q is a measure on M Q (for instance, the quantum symplectic volume). Average values of hybrid observables can then be computed as where dµ C is a measure on M C (for instance, the classical symplectic volume) and then dµ QC = dµ Q dµ C . Notice that we use the term quantum moment since, considering a basis for states ρ ψ , the coordinates ofρ(ξ) correspond, precisely, to the moments of the probability distribution F QC considered as a function of the coordinates of the quantum projectors. Furthermore, from the mathematical point of view, we can relate this first quantum moment with the marginal and conditional probabilities of the original density F QC . As a bivariate hybrid measure, we can consider its classical marginal probability and the corresponding quantum conditional one: • The classical marginal probability F C (ξ) can be obtained as • The quantum conditional probability for a classical value ξ ∈ M C is a purely quantum distribution. We can consider it to be a distribution over the quantum phase space with some density F cond ξ (ρ ψ ) on the quantum phase space M Q . As such, Gleason theorem [?] ensures that there exists a well defined quantum density matrixρ ξ which is able to determine the average value of any quantum observable (for a fixed value of the classical variable ξ).
As it was explained in [24], standard probability theory allows us to write the operator ρ(ξ) as the product of these two objects, i.e., Notice that the resulting object is not normalized as a density matrix, since This operatorρ(ξ) is the usual representation of a hybrid state in the literature, since the early eighties [25]; and it is the usual choice in most of the references we presented in the Introduction. Notice, though, that from the probabilistic point of view, it just represents an element of the dual space to the set of hybrid observables of the form of Equation (5). If we consider the C * -algebra defined by the operators of that type, operatorρ(ξ) defines a state for such an algebra. From this point of view, this first quantum moment of F QC is the only relevant part of the probability distribution if we restrict the physical magnitudes to functions of the form (5), since higher moments do not contribute to any average value. Hence we can define an equivalence relation in the space of hybrid measures by the first moment: two hybrid measures represent the same physical state if and only if their first quantum moments coincide.
Furthermore, this first momentρ(ξ) is also able to capture the mutual exclusivity of hybrid events (see [24]). Indeed, notice that the points of the hybrid phase-space M C × M Q do not represent mutually exclusive events from the probabilistic point of view, as it happens with a classical manifold. While a classical system at a point ξ 1 can not be at the same time at point ξ 2 (the probabilities of one case and the other are independent hence), it is not incompatible for a quantum system to be at point ρ ψ1 and at point ρ ψ2 unless ⟨ψ 1 , ψ 2 ⟩ = 0 (probabilities are not independent then). From that point of view, representing a hybrid state with the density F QC makes very difficult to define a hybrid entropy function (which requires of a correct representation of independent events which, as we just argued, is not achieved on M C × M Q). The first moment representationρ(ξ), on the other hand, represents well the independence of hybrid events (as von Neumann entropy does for quantum systems) and allows us to define a well-behaved entropy function. From it, we can also use the MaxEnt principle to define a consistent notion of canonical ensemble and therefore a mechanism to make finite-temperature numerical simulations of hybrid systems (see [24,30] for details).

Considering dynamics: the necessity of higher order moments
From a static point of view, the picture above is perfect. For the natural set of hybrid observables, we can consider a state corresponding to a probability density F QC on the hybrid phase space M C × M Q , but only the first quantum moment is necessary to define averages of physical magnitudes. That first moment also allows us to define a consistent notion of hybrid entropy. The problem arises when we want to consider a dynamical framework for statistical averages or even to equilibritum thermodynamics.
In this paper, we are going to consider Ehrenfest dynamics as the candidate for the dynamics of the statistical microstates. It is immediate to verify that dynamical equation (1) defines also a dynamics on the set of hybrid magnitudes, which can be written by means of the Hamiltonian vector field X f H . Indeed, we can write: Clearly, a function of the form is a solution of Ehrenfest equations, is itself a solution of Equation (13) with initial condition f A (ξ, ρ ψ ; 0) = f A (ξ, ρ ψ ). Therefore, the system has solutions, and we can consider the time-dependence of statistical averages and define out-of-equilibrium statistical mechanics of hybrid systems. But it is important to realize that hybrid observables of the form (5) do not define a Poisson subalgebra for the hybrid Poisson bracket (2). Indeed, given fÂ (ξ) and fB (ξ) defined as Equation (5), in general there exists no function fD (ξ) of the form (5) This property can be proved immediately since the classical bracket {fÂ (ξ) , fB(ξ)} C defines a function which is quadratic in the quantum degrees of freedom while the quantum bracket {fÂ (ξ) , fB (ξ) } Q is still linear. Just those hybrid magnitudes which depend only on the classical or on the quantum degrees of freedom are closed under the bracket. After a straightforward computation, and choosing Darboux coordinates ξ = (⃗ q, ⃗ p) = {(q k , p j )} k,j=1,··· ,n for M C we obtain that: and where Therefore, we can write that i.e, we obtain a combination of quantum-linear functions and quantum-quadratic ones. Clearly, successive brackets will increase further the quantum-degree. From the physical point of view, this relation implies that, even if it is Hamiltonian, the dynamics on the space of observables can not be restricted to quantum-linear functions of the form of Equation (5). Nonetheless, hybrid energy is conserved, since If we average this expression with the density F QC , we can write that Notice that in this expression the requirement of the extension of the algebra is explicit, since we need the average values of pointwise products of the elements of the algebra. Only the last term belong to it. If we write it as an equation on the density F QC , the last term affects only its first moment, i.e., But the other two can not be written as a function of the first moment only, and higher order moments are required. Indeed, we can write the average value of the product of the functions fÂ (ξ) and fB (ξ) as whereρ ⊗2 (ξ) is the second quantum moment of the density F QC , defined aŝ Analogously, the average value of the product of k functions fÂ j (ξ) is obtained as whereρ represents the k-th quantum moments of the distribution F QC .
In conclusion, in order to consider Equation (13) as a dynamical system on the hybrid algebra, the algebra itself must be enlarged. As each classical bracket increases the order of the quantum degrees of freedom by one, clearly we must consider all polynomial functions in ρ ψ and arbitrary classical dependence as the minimal Poisson algebra generated by hybrid functions of the form (5). At each level, higher and higher quantum moments of the distribution are required to compute the average values of the extended functions. Indeed, if we consider the series of derivatives of order k with respect to time of the average value ⟨Â(ξ)⟩(t), we can write it as a function of different terms depending on all the momentsρ ⊗j (ξ) up to order k + 1 evaluated on different combinations ofÂ(ξ) andĤ(ξ): Hence, for a finite range of time (t 0 , t 0 + ∆t), the behavior of an average value ⟨Â(ξ)⟩(t) can be approximated to arbitrary precision by a finite number of quantum moments of the distribution F QC . Notice that the particular time scale and the precision depends on the particular observableÂ(ξ). In any case, by considering a finite time range we see that the concept of hybrid observable changes. Time introduces correlations between the classical and quantum subsystems and the nonlinear evolution makes that correlation measurable in the time-dependence of the average values. We can see that in Equation (19). If we want to consider a linearized time dependence in ⟨Â(ξ)⟩, we need to compute average values of pairs of quantum operators 1 . The time range which can be approximated by a finite number of time derivatives may depend on the observable we consider. But for any observable, there will be a certain time scale where the behavior of the system can be approximated to any precision with a sufficiently high order of quantum moments. As the application of our model will be a numerical simulation of the system at a certain time scale, this finite number of moments should be enough to characterize the state of the system. In the following sections we will learn to write the dynamics of the quantum moments and thus the dynamics of the physical macrostate for the relevant time range.

Writing a dynamics for the physical states
Let us consider now how to define an equation to write the evolution of the average value ⟨Â(ξ)⟩ given by (14) as an evolution equation on the physical state. Following [23] we can associate a Liouville equation to the probability density and define the curve on the space of probability densities which reproduces the evolution of the average values ⟨Â(ξ)⟩(t), i.e., where F QC (ξ, ρ ψ ; t) is the solution of the master equation. In [5] we proved that because the dynamics is Hamiltonian we can obtain a dynamical equation for the density F QC as the Liouville equation This equation translates the microstate dynamics to the full probability density. Nonetheless, if we are interested in the behavior of average values of physical magnitudes for a certain time scale, we know that not all the density is required, only a few of its quantum moments. Hence, we are going to study now how can we write the effect of Ehrenfest dynamics on those moments.

The initial problem: the dynamics of the first quantum moment
Let us consider now the first moment of the distribution F QC , i.e., the family of operatorsρ(ξ). We can consider Liouville equation for the density and re-write it in terms of the first moment: As the time dependence of average values can be written as a result of dynamics written on the space of operators or on the space of states, Equation (28) where the dot represents the tangent vector of each respective curve (on states or on operators). Notice that, by definition, Therefore, we can conclude that the properties of the first quantum moment will be preserved (normalization, positivity, etc), if the solution of the Liouville equation defines a curve of well-defined probability densities. As we saw in the previous section, the set of linear functions of the form f A is not closed under the hybrid Poisson bracket. This has an important impact on the Liouville equation (27). Indeed, a master equation for the density F QC which captures the dual behavior of the dynamics of the hybrid observables (which does not have linear functions as a Poisson subalgebra) can not be restricted to the first quantum moment, which is only able to capture the linear (quantum) functions of the form (5). If we must include quantum polynomial functions of orders higher than one (i.e., objects of the form f A f H ), higher quantum moments of the distribution F QC must be considered. But then the equivalence relation given by the first moment is broken: two measures having identical first quantum moments may evolve in different ways and produce different hybrid states for each value of time, depending on the higher quantum moments. Let us consider this issue in some detail in a particular example. Example 3.1. As we saw above, there exist many densities F QC with the same first quantum moment. Let us see now that the nonlinearity of Ehrenfest equations makes impossible to write Ehrenfest dynamics in termsρ(ξ) in a consistent way.
Let us now construct a particular hybrid system with Ehrenfest dynamics, with an ad-hoc Hamiltonian and density matrix. We will then find a one-parameter family of different densities with the same first momentρ(ξ) and see that the tangent vector defined by each density is different. More in particular: 1. We consider one classical degree of freedom, M C ≡ R 2 and ξ ≡ (R, P ).
2. We choose a two-level system as the quantum subsystem. Therefore, the space of projectors is: being the Pauli matrices and σ 0 the identity. The usual description of this space in quantum mechanics is the Bloch sphere. It can be proved that M Q as defined above is bijective to S 2 = {x 2 + y 2 + z 2 = 1 | x, y, z ∈ R} ⊂ R 3 , with coordinates {x = Tr[σ 1 ρ ψ ], y = Tr[σ 2 ρ ψ ], z = Tr[σ 3 ρ ψ ]}, which will be called Bloch coordinates, denoted by a subscript B.

The Hamiltonian of the system is defined aŝ
whereπ k ∈ M Q are classical-point-dependent projectors:π 1 (R, P ) = ( sin R, 0, cos R ) B andπ 2 (R, P ) = ( − sin R, 0, − cos R ) B . Note that they are orthogonal by construction. The two energy levels E k are:
Finally, we define a family of densities F θ QC in terms of the arbitrary parameter θ ∈ R. Consider a decomposition of the matrixρ(ξ) as sum of projectorŝ It is immediate to prove that the density hasρ(ξ) as first moment. We can consider a family of the decompositions choosing the coordinates of one projector to be (sin θ, 0, cos θ) B , and solving the geometric problem in the Bloch sphere to find the other projector and the corresponding weights. Thus, we construct the θ-dependent family of densities F θ QC having the same first quantum moment. For each density, consider Equation (28) and represent the corresponding tangent vector as a function of the parameter θ. We present graphically the results in Fig. (1). We observe that the x and z coordinates ofρ(ξ) (in the Bloch basis) depend of the arbitrary parameter.
As expected, we verify that the dynamics of the first quantum moment depend on the parameter θ and therefore that the dynamics of the first moment is not well (a) x-coordinate, 1 2 Tr(ρ(ξ)σ 1 ). defined. Hence, we see that different density dynamics are possible for a givenρ(ξ), depending on elements which are not physically observable, a prioiri, at the (quantum) linear level. Notice that the origin of the problem is the nonlinearity of Ehrenfest dynamics, which defines a different evolution for the two different quantum states coupled to the same classical point ξ. Had we considered a pure quantum system with unitary evolution, or even an uncoupled hybrid system, the evolution would be well defined because of the linearity of the pure quantum evolution. Let us see now how the dynamics can be determined when considering higher order moments.

The dynamics of higher quantum moments
We have seen above how at the linear level, we have part of the information of the full Liouville equation, although not all. Nonetheless, Equation (28) is enough to read the behavior of the two marginal distributions F C (ξ) = Trρ(ξ) = M Q dµ Q F QC (ξ, ρ ψ ) andρ = M C dµ C (ξ)ρ(ξ). Indeed, marginalizing one of the variables makes one of the two brackets to vanish, and that produces: and These two equations coincide with those obtained from different approaches to hybrid quantum-classical dynamics, such as [25,12], where Ehrenfest equations were not considered. This is a remarkable result, since it ensures that, at least at first order, the evolution of pure classical or pure quantum magnitudes behave the same as in those other approaches even for coupled quantum-classical systems.
What about the full hybrid (i.e., non-marginalized) behavior? Can we read any useful properties from Liouville equation? Indeed, we can. We saw above that writing Ehrenfest dynamics on the physical states requires of higher order quantum moments to be able to capture the effects of the nonlinearity and the non-closeness of the Poisson algebra. Dynamics written on the first moment only is not well defined, as the nonlinearity of the dynamics makes it dependent on the choice of the decomposition of the quantum operator at each classical point. In this section, we will see that the dynamics of each moment depends on the next one; the Liouville equation being thus equivalent to an infinite series of differential equations, one for eachρ ⊗k (ξ).
Let us consider again the definition of the k-th moment given in Equation (24). By construction, they are operators on the Hilbert space k j=1 H j . We can define partial traces if we average over some of those copies, since Trρ ψ = 1 and hence where the zeroth order moment corresponds to the marginal classical density Despite of the symmetric nature of the expressions involving the points ρ ψ in the following we will consider traces of different operators, and will denote as Tr 1 the trace acting always on the first one In an analogous way, the expression is extended to arbitrary traces Tr k . Using these relations, it is immediate to re-write Equation (19) as From that expression, and using the compact support of the classical functions, we can read the time derivative of the first quantum moment aṡ whereρ(ξ) is defined by Equation (28). Analogously, we can define the time derivative of any momentρ ⊗k (ξ) as a function ofρ ⊗k (ξ) andρ ⊗k+1 (ξ).
In this way, we find a system of differential equations for all the quantum moments {ρ ⊗k (ξ)} k=0,··· , which is equivalent to the Liouville equation for the density F QC (Equation (27)).

An effective equation
As we saw above, only a finite number of quantum moments are necessary to recover the behavior of the average value of physical observables for finite time intervals. Notice, though, that when we make the approximation for a certain time range, the approximated solution of the dynamics of the hybrid observablesÂ(ξ; t) defines a different dynamics for the average values ⟨Â(ξ)⟩(t). As in each case the average values involve a finite number of quantum moments, we can define different dynamics on the space of states with the corresponding Equation (41). But there is a problem to do that: the dynamics of the k-th quantum moments, depends on the k + 1-th one and therefore it is not possible to solve the equation for a certain k without the solution to the next level. Nonetheless, as the time scale is only able to capture the kfirst moments, it is not possible to do any measurement or sequence of measurements which depend on those degrees of freedom. Hence, as the (k + 1)-th moment can not be known, we may consider an effective dynamical equation forρ ⊗k (ξ) where the value ofρ ⊗k+1 (ξ) represents the minimum knowledge on the system compatible with the physical constraints. In order to do that, we must introduce a suitable notion of entropy.
In [24] we introduce a notion of hybrid entropy for the first quantum moment ρ(ξ). We built it based on the bivariate distribution that the first moment represents, and the corresponding factorization in marginal and conditional probabilities. That factorization makes sense for any of the quantum momentsρ ⊗k (ξ), which, on the other hand, happens to be formally analogous toρ(ξ) with the only difference of being an operator on the product Hilbert space k j=1 H: • the marginal probability of the bivariate probability distribution represented bŷ ρ ⊗k (ξ) corresponds to the classical density • The corresponding quantum conditional probability P(ρ ψ | ξ) is the probability of a pure quantum system (the classical degrees of freedom ξ are fixed). As such, Gleason theorem ([?]) ensures that a well-defined density matrixρ ⊗k ξ (i.e., a trace class, normalized and non-negative self-adjoint operator on k j=1 H) must exist to represent it.
This representation allows to implement the exclusivity of hybrid events in a simpler way, since it incorporates automatically the orthogonality of the mutually exclusive quantum events (as the usual density matrix does). See [24] for a more detailed discussion.
By using this factorization, we can adapt the construction of the hybrid entropy in [24]. As the entropy for the stateρ ⊗k (ξ) must be written as the sum of the entropy of the marginal distribution (43) and the classical average of the entropy of the conditional one we obtain This function represents the entropy of the hybrid system modelled by the k-th quantum moment. Using this entropy function, we may consider searching for the stateρ ⊗k+1 (ξ) which maximizes the entropy subject to some constraints to use in the right hand side of Equation (41). The natural constraint is asking this first unknown moment to produce the last known one by trace, i.e., we will search for a k +1 moment such that whereρ ⊗k (ξ) is the moment whose dynamics is being defined. The unknown components of the moment are thus only the degrees of freedom being integrated out by the trace. As it is not possible to determine that information by any measurement, we select the state with the maximum possible entropy for those degrees of freedom to act as a source for the dynamical equation. In the following we shall represent aŝ ρ ⊗k MaxEnt (ξ) the state which maximizes the entropy with constraint (45). Notice that we are using the MaxEnt formalism from the point of view of information (see [31]), and not as it is usually done when trying to identify thermodynamical ensembles. Those thermodynamical problems can also be considered, though. The entropy function would be the same, but depending on the problem, other constraints would be considered. For instance, we may study the state which maximizes entropy (44) while keeping the normalization and the average value of the HamiltonianĤ(ξ) fixed. This would produce the candidate to model the canonical ensemble in our context, as we considered in [24] for the case of the first moment. In a similar way we may consider a microcanonical ensembleρ ⊗k M CE (ξ). All these situations will be considered in detail in a forthcoming paper.
We conclude thus that this MaxEnt formalism with constraint (45) provides us with a natural candidate to represent the stateρ ⊗k+1 (ξ) in Equation (41). As there is no physical measurement which can inform us about the state, we consider it to be the state which maximizes the uncertainty about it, which is preciselyρ ⊗k+1 MaxEnt (ξ). Therefore, we can write as effective dynamics for the k-th quantum moments which is the equation: This equation is well defined, since it only depends on the degrees of freedom of the k-th moment (given that the constraint (45) must hold for all times) and can be used to approximate the dynamics of the average values of physical magnitudes in a finite range of times for statistical systems whose microstates follow Ehrenfest dynamics.
Example 3.2. In order to provide a practical application of our framework, let us consider again the simple example considered above of a system with two dimensional classical phase space and a two level quantum system. Let us assume, for simplicity, that we just consider Equation (46) the first quantum moment. In that case, we have to identify the solution of Maximal Entropy for the case k + 1 = 2. Let us proceed. First of all, we are going to consider the factorization in marginal and conditional probabilities given by Equation (11), forρ(ξ) andρ ⊗2 ξ : By construction,ρ ξ is written as a combination of the basis {σ j } j=0,1,2,3 , and thereforê ρ ⊗2 ξ must be written as a combination of the 10-dimensional symmetrical basis { 1 2 (σ j ⊗ σ k +σ k ⊗σ j )} j≤k=0,1,2,3 . The corresponding expressions read: where we can relate these coordinates with the average values with respect to the conditional distribution F cond ξ of the coordinate functions of the points of M Q as: µ j (ρ ψ ) representing the j-th coordinate of the pure state ρ ψ with respect to the basis {σ j }. Furthermore, asρ ξ andρ ⊗2 ξ are density matrices, their trace must be normalized. This implies that In this context, constraint (45) implies that 2µ 0j = µ j , j = 1, 2, 3.
Notice, though, that from the definition ofρ ⊗2 (ξ) there are also a few relations between the coordinates, to be fulfilled. In particular • Asρ ξ must be a well defined density matrix, we know that its purity must be lower or equal to one, i.e.
• An analogous property holds forρ ⊗2 ξ : its purity must also be lower than one, i.e.,: • As µ jj is the expectation value of the square of the function µ j (ρ ψ ) on M Q we can easily verify that since the variance of the function µ j (ρ ψ ) must be positive definite for a welldefined probabilistic system. • As points ρ ψ ∈ M Q are pure states, j µ 2 j (ρ ψ ) = 1 2 and therefore • Analogously, we can verify that, from Cauchy-Schwartz inequality With these constraints, we search for the MaxEnt solution for the hybrid entropy function (44). The result defines a relation between the free variables ofρ ⊗2 ξ with the other variables, i.e. which depend on those ofρ ξ : Computing entropy involves a complicated expression of the logarithm ofρ ⊗2 (ξ) (or, analogously, ofρ ⊗2 ξ if we use the marginal-conditional factorization). For the sake of simplicity, instead computing the spectrum of the matrix, we can use Mercator series to approximate it by a polynomial of the traces of the powers of the matrix, i.e.
Keeping only the linear approximation of the entropy, we would obtain Maximizing this function in the region compatible with the constraints above, considering the symmetry between the different group of indices, can be achieved by fixing: where P(ρ ξ ) = j µ 2 j represents the purity of the stateρ ξ . This point can be considered just a simple approximation to the real MaxEnt solution, but it helps us to illustrate how our construction works. Higher order approximations can be implemented with numerical tools on the expansion above.
We can see in these equations how classical and quantum degrees of freedom evolve coupled, the marginal classical density µ 0 (ξ) = F C (ξ) depending on the quantum degrees of freedom µ j (ξ), whose evolution is also governed by µ 0 (ξ). Furthermore, from the expressions of Equations (61) and (62) and (63), it is immediate to see that the final dynamics ofρ(ξ) is non-linear.

Conclusions
In this paper we have considered the problem of dynamical statistical systems for hybrid quantum-classical systems having Ehrenfest dynamics as microstate dynamics. We have seen how dynamical evolution (or any other type of observable-generated transformation) enlarges the set of physically observable magnitudes since it introduces correlations between the magnitudes of the formÂ(ξ) which are the natural models of hybrid magnitudes. The nonlinearity introduced by the classical subsystem is responsible of that correlation, which produces a progressive enlargement of the hybrid subalgebra. This implies that, while (linear) hybrid observables depend only on the first quantum moment of the distribution, higher quantum moments are required to compute the evolution of the average value of physical observables. For a given time interval and a certain accuracy, only a finite number of those quantum moments are required to approximate the evolution of the average values of arbitrary physical magnitudes. These quantum moments are also useful to implement a consistent notion of hybrid entropy in a simple way, since they take the form of a family of density matrices indexed by the classical degrees of freedom, which encode in a simple way the exclusivity of quantum (and hence hybrid) events. In future works, we plan to analyze in detail the implications of these new entropy functions, and some relevant examples of ensembles arising from them, such as the canonical or micro-canonical ensembles. Furthermore, the stability of the resulting distributions with respect to the dynamics can now be studied in a simple way. Following that direction, the analysis of equilibrium thermodynamics of hybrid systems is a feasible objective.
We have also been able to re-write the Liouville Ehrenfest dynamics of the full probability density at the level of the quantum moments. We have obtained a series of coupled differential equations for the set of quantum moments to encode the physical content of Liouville equation. Furthermore, using the fact that only a finite number of quantum moments are required to model the system for finite time intervals, we have considered the problem of defining the dynamics with a finite number of those equations. In order to do that, we need to write the effect of the next moments on the highest one. As those higher moments can not be determined by physical measurements, we choose to represent those degrees of freedom by the state, compatible with the constraints of the problem, which maximizes the hybrid entropy at that level. With this choice, a well defined differential equation can be written for the system.
An important consequence of our results is the change it implies at the level of implementing numerical simulations of statistical hybrid systems. We have learned that an accurate description of a statistical hybrid system requires taking into account the correlation between quantum observables mediated by the nonlinear classical subsystem, i.e., the enlargement of the algebra of hybrid observables. This property incorporates into the description more quantum moments which are not required to obtain the initial average values of physical magnitudes, but that are necessary to estimate their evolution. From the point of view of numerical simulation this introduces some difficulties, since incorporating these higher order correlations into the usual simulation methods as the use of independent (pure-state) trajectories is not an easy task. While fixing values of single operators is simple, considering correlations of several magnitudes represents a new numerical challenge. The issue of knowing the statistical value of the correlations between observables is closely related to the issue of preparation of statistical quantum (in this case, hybrid) states. On the other hand, the preparation of states is very often implemented in simulations through a sampling over quantum phase space of initial conditions that reproduce the average initial values of the characterized observables. Our results show that for non linear dynamics it is not feasible to just approximate the quantum state as a bundle of trajectories with a statistical sampling over quantum phase space of initial conditions, knowing just the expectation value of operators. The main difference with von Neumann's evolution is that such master equation given in terms ofρ is linear on quantum projectors and such linearity preserves the original sampling throughout the evolution. In this hybrid case, such sampling must reproduce correlations and higher order observables up to the desired level of precision (associated with k + 1 in the truncation above), or work directly to the set of k density matrices, characterized through such measurements, to be able to implement the coupled evolution between subsequent powers of quantum projectors in such non-linear dynamics.