We propose a novel general approximation to transform and simplify the description of a complex fully quantized system describing the interacting light and matter. The method has some similarities to the time-dependent Born–Oppenheimer approach: we consider a quantum description of light rather than of nuclei and follow a similar separation procedure. Our approximation allows us to obtain a decoupled system for the light-excited matter and “dressed” light connected parametrically. With these equations at hand, we study how intense light as a quantum state is affected due to the back-action of the interacting matter. We discuss and demonstrate the possibility of the light-mode entanglement and nonclassical light generation during the interaction.

## MOTIVATION

Light is a unique tool to control and probe the world of matter/particles because of its fascinating properties: nowadays, technologies allow us to generate, measure, and control light radiation at an extremely precise level, both on spatial and temporal scales. This opens extraordinary possibilities to explore the quantum nature of our world at the micro/nanometer spatial scales and attosecond time scales and possibly will bring some new quantum phenomena to the macroscopic level. The quantum-mechanical description of both the matter (charged particles) and the electromagnetic radiation is universal: within conventional non-relativistic quantum theory, everything can be described by the total wave function, which is a solution of the corresponding time-dependent Schrödinger equation (TDSE) with a complete Hamiltonian including all the interactions. Since light is a quantum object,^{1} the complete and consistent treatment requires a fully quantum description for the whole interacting system including both the charged particles and the light modes. It is known that the calculation of the complex quantum problems is a very complicated issue, especially for strongly interacting, multidimensional systems. Here, we present a novel theoretical method that simplifies the problem significantly while still taking various purely quantum properties into account.

One of the most general and rigorous approaches that combine classical electrodynamics (non-relativistic for the particles) and quantum principles within the time-dependent Schrödinger equation (TDSE) employs the so-called *Pauli–Fierz Hamiltonian* (PFH) [see Ref. 2 and Eq. (1)]. This Pauli–Fierz Hamiltonian contains the description of the dynamics of the charged particles coupled with the electromagnetic fields in a fully quantum way. This unifies Maxwell’s equations and Newton’s equations at the fully quantum description. Within this theory, it is possible also to include the spin of the particles and several semi-relativistic effects.^{3,4} The possible solutions of this TDSE could open various opportunities for applications ranging from new sources of light or states of quantum matter to unique pure quantum collective effects. Recent experiments demonstrate some surprising quantum collective phenomena beyond classical or even semiclassical pictures.^{5,6}

This paper is organized as follows: we start with the TDSE description with the PFH. Next, we derive the related simplified system of the parametrically connected equations. Finally, we consider selected examples/solutions and discuss the consequences and the possible applications.

## GENERAL DESCRIPTION OF THE LIGHT–MATTER INTERACTION

The total Hamiltonian for the charged particles (numbered with indices *k* and *k*′) and the multi-mode electromagnetic field (modes numbered with index *j*) in the Schrödinger representation with Coulomb gauge reads as follows (known as Pauli–Fierz Hamiltonian; see overview and details in Ref. 2):

Here, *e*_{k} and *m*_{k} are the charges and masses of the particles with index *k*, and $\u2211kk\u2032Ukk\u2032$ includes particle–particle mutual interactions, as well as interaction of each particle with external longitudinal fields. Each mode *j* of the electromagnetic field is quantized via the term $\omega jN\u0302j$ in the Hamiltonian. This is the most general formulation, without including spin (which can be extended separately within Pauli–Fierz Hamiltonian) within a non-relativistic description of matter (particles).

## THE IDEA OF PARAMETRICAL-CONNECTION APPROXIMATION

Based on the PFH, we first consider the quantum interaction of two objects of different nature (one electron and the light field) having the total Hamiltonian of the form $H\u0302=H\u0302e+W\u0302int+H\u0302field$. We also consider that the initial state of the whole system is of the form Ψ_{0}(*x*, *q*) = *F*_{0}(*x*) ⊗ *G*_{0}(*q*), where *x* is related to the electronic coordinates and *q* is related to the field-modes coordinates, respectively. In the strong-field case, when $H\u0302e\u223cW\u0302int\u226aH\u0302field$, it is reasonable to assume that the quantum state of the light field changes relatively weakly during the interaction compared to the initial state *G*_{0}(*q*) since the interacting electron is a small perturbation to the strong field. This leads to the idea of the following (full first-order) approximation of the total solution: $\Psi (x,q,t)\u2248F(G0)(x,q,t)\u2297G(F)(q,t)$, with the error proportional to some degree of the ratio $W\u0302int/H\u0302field$. The most significant simplification in the following is that the solution for the electronic subsystem $F(G0)(x,q,t)$ is obtained from the TDSE with only local operators on *q* (there are no derivatives), i.e., it is an effectively one-object TDSE with parametrical *q*-dependence. The corresponding *q*-derivatives can be transformed to the simple *q*-dependent expressions using a non-perturbed *G*_{0}(*q*) within this first-order approximation. Thus, we transform the initial two-object TDSE to a simplified system of two equations, one for $F(G0)(x,q,t)$, where *q* is a parameter, and other for the field-state evolution *G*_{(F)}(*q*, *t*) including a small electron back-action via *F*. The detailed derivation of the system for the coherent (Gaussian) initial state of the intense light is given below. The weaker approximation Ψ(*x*, *q*, *t*) ≈ *F*(*x*, {*q*}, *t*) ⊗ *G*_{0} leads immediately after the field-state averaging of the total Hamiltonian in the high photon number limit to the known (fixed-amplitude) classical field limit.^{7}

## DETAILED DERIVATION OF THE GENERAL APPROXIMATION

We first consider one particle—an electron interacting with an arbitrary given configuration of propagating electromagnetic plane waves. A similar problem statement was recently considered for a cavity-photon Hamiltonian in the dipole approximation.^{8} We note that our approach can be transformed for more complex systems with many particles, as well as for more specific ones (atomic units *ℏ* = *m*_{e} = |*e*| = 1 are used in what follows). We use the coordinate representation of the field operators,^{9–11} with the creation and annihilation operators $a\u0302j+=12qj\u2212\u2202\u2202qj$ and $a\u0302j=12qj+\u2202\u2202qj$ and the field dynamical coordinate *q*_{j}. We assume that the initial state (at the time *t* = 0) of each mode *j* is coherent, $Gj(qj,0)=G0j=C0j\u22c5exp\u221212(qj\u2212e\u2212i\theta j2Nj)2$, with *N*_{j} being an initial arbitrary average number of photons and *θ*_{j} being the initial phase (see the Appendix). In this case, if we apply the following transformation $\u220fjexp[\u2212i(\omega jt+\theta j)N\u0302j]$ (see also Ref. 9) to the general TDSE with Hamiltonian equation (1), we obtain the TDSE in the interaction picture with phase-dependent field operators and phase-independent mode states,

The initial state for each mode *j* of the light field in this representation is $Gj(qj,0)=C0j\u22c5exp\u221212(qj\u22122Nj)2$, the TDSE for the field states in vacuum is $iG\u0307jvac=0$, and the quantized field vector potential operator reads

where $\chi \u20d7j$ are polarization vectors (linear in this case) and $\beta j=c2\pi /\omega jV$ are parameters (with the same dimensionality as *A*_{j}) that depend on quantization volume *V*, according to the conventional plane-wave quantization procedure.^{9} In our case, when the initial mode state is a coherent state, each *β*_{j} parameter naturally connects the field mode amplitude *A*_{0 j} and the initial average photon number in this mode *N*_{j} as follows: $A0j=\beta j2Nj$, and consequently, $\beta \u223cW\u0302int/H\u0302field12$. A typical laser pulse, which is finite in space and time and consists of a number of modes, can be described in terms of an effective quantization volume *V*_{eff}. This volume is related to the sizes of the radiation focusing area and the effective interaction volume (see also Refs. 9 and 12). Based on the estimations for the number of relevant experiments with strong laser pulses interacting with gases or solids, where quantum-optical effects are measurable,^{5,6} the effective values are (*β*_{j}/*A*_{0j}) ∼ 10^{−7} or even much less. These parameters are also responsible for the quantum dispersion of the quantized radiation in the presence of charged particles (see Ref. 13). We refer to *β*_{j} as the fundamental scaling parameters for the corresponding light–matter interaction problem when *E*_{0j} ≪ *E*_{atomic} in the optical wavelength range. Due to its very small values, we may conclude from Eq. (2) that even in a strong-field case, the back-action of one electron to the light-mode wavefunction is small and depends somehow on *β*_{j}. Therefore, the state of electromagnetic field *G*(*q*_{j}, *t*) is affected as *G* ≈ *G*_{0} + *∑O*(*β*_{j}). Finding the equation for *G* is a key point of our study.

It was rigorously demonstrated that the solution of the two-component TDSE with an arbitrary Hamiltonian can be written as a product two special wavefunctions each satisfying the normalization condition (see Ref. 14 for full details). In our case, it leads to the following exact factorization: $\Psi =Fex(r\u20d7,qj,t)\u22c5Gex(qj,t)$, with the related system of equations for the *F*_{ex} and *G*_{ex} of the same complexity as the original TDSE (see also Ref. 14). Instead, based on our previous considerations regarding the dependence on parameters *β*_{j}, we consider the solution of Eq. (2) of the form $\Psi \u2248FG0(r\u20d7,qj,\beta j,t)\u22c5G(qj,\beta j,t)+\u2211O(\beta js)$. In our approximation, we have the same form of the solution with the same partial normalization conditions: $F|F=\u222b|F(r\u20d7,qj,t)|2d3r\u22611$ and $G|G=\u222b|G(qj,t)|2[dqj]\u22611$. The representation is not unique since one can construct other one by choosing $\Psi =F(r\u20d7,qj,t)ei\chi (qj,t)\u22c5G(qj,t)e\u2212i\chi (qj,t)$ with an arbitrary function *χ*(*q*_{j}, *t*) matching the initial condition at *t* = 0 (see also Ref. 14). In order to have a uniqueness of the representation for a given initial state, we require the condition $F|F\u0307=0$, which is often referred to as a static gauge of the factorization.^{7}

Using the considered form of the solution, we can write

and obtain the decoupling in the first-order approximation,

From this, it follows that in a first step of a reasonable approximation, for the electron subsystem, we may use *β*_{j}-independent states of light modes for the approximation of the vector potential operator. Thus, we can exclude the non-local operator $\u2202\u2202qj$ in the equation for the electron subsystem (since $\beta j\u2202\u2202qjF\u223c\beta j2$), having instead a simple local operator as follows: $\beta j\u2202\u2202qj\u2192{\beta jG0j\u22121\u2202\u2202qjG0j}$. Starting from Eq. (5a), we assume that *G*(*q*, *t*) does not have zeros within a time interval of consideration. This can be in the case when one starts with coherent state and considers conditions of relatively small deviations, leading, for example, to the various nonclassical/squeezed states. Otherwise the procedure can be modified in a way that instead of Eq. (5a), we can write $iGF\u0307=H\u0302G0F\u2212iG\u0307F$ with the similar subsequent related considerations. The described approach allows us to decouple the total system Eq. (2) and separate the equations for the electronic subsystem and “dressed” light subsystem via *β*_{j} parametrical connections. Finally, for the quantum light–electron interaction problem, we obtain

where $A\u20d7jpar=\chi \u20d7j{\beta jqj}cos(\kappa \u20d7jr\u20d7\u2212\omega jt\u2212\theta j)\u2212i\chi \u20d7j{\beta jG0j\u22121\u2202\u2202qjG0j}sin(\kappa \u20d7jr\u20d7\u2212\omega jt\u2212\theta j)$ is a field vector potential operator that depends on *q*_{j}-expressions and a set of scaling parameters {*β*_{j}}. In the strong-field case with coherent initial state of light, when *N*_{j} ≫ 1 (and consequently $\u27e8qj\u27e9=2Nj\u226b\u27e8\u2202\u2202qj\u27e9=0$), we can simply write $A\u20d7jpar=\chi \u20d7j{\beta jqj}cos(\kappa \u20d7jr\u20d7\u2212\omega jt\u2212\theta j)$, noting that this is nothing else than the parametric description of the classical vector potential ({*β*_{j}*q*_{j}} corresponding to the mode amplitudes in the classical description). This description, however, connects naturally both subsystems: the interacting electron subsystem *F* with field parameters {*β*_{j}*q*_{j}} and the “dressed” light field that depends on *F* via the averaging procedure of $H\u0302B=FH\u0302F$. The error of our first-order approximation in this case is proportional to some combination of $\beta j1$.

Now, let us apply our approach to some simple model cases of light-electron interaction. We consider the strong-field case when *N*_{j} ≫ 1. For simplicity, we consider also the so-called dipole approximation when the dependence of the vector potential on the coordinate is neglected, since $|\kappa \u20d7j(r\u20d7\u2212\u27e8r\u20d7\u27e90)|\u226a1$, which is common for the interaction in the optical wavelength range.

First, let us consider a non-bound electron wave packet ($U(r\u20d7,t)\u22610$), localized in some region of space with initial momentum 0. Then, the system of Eqs. (6a) and (6b) can be solved analytically as follows:

Equation (7a) describes the dynamics of the electronic wave packet in a given field vector potential $A\u20d7P=\u2211jA\u20d7jpar({\beta jqj},t)$, while Eq. (7b) gives the correction to the initial field state *G*({*q*_{j}}, 0) due to the interaction with the electron. The solution of Eq. (7b) includes terms exp(*d*_{ij}(*t*)*q*_{i}*q*_{j}) [where *d*_{ij}(*t*) are some time-dependent functions]; thus, it demonstrates entanglement between different modes during the interaction. In the case of a one-mode field, the solution of Eqs. (7a) and (7b) corresponds to the closed-form solution^{15} with a relative error ∼*β*^{2}.

Let us now consider the more complicated case, when the electron is bound in a quadratic potential *U* = *u*(*t*)**x**^{2}. This situation can be related to some model bound potential [*u*(*t*) = const] or to the problem of an electron in magnetic fields (time-independent or time-dependent). Since, in that case, the total Hamiltonian $H\u0302P$ for the interacting system is quadratic, the solution for Φ can be obtained analytically. The analysis of the solution is beyond the scope of this article; however, the main exponential part of the solution is evident,

where the coefficients *a*_{j}(*t*), *b*_{j}(*t*), and *d*_{ij}(*t*) depend on the initial states of the modes, initial state of the electron, and *u*(*t*). We would like to emphasize that the parameters of the exponential part [Eq. (8)] can be controlled by external parameters in various ways. This means that depending on the conditions and the binding potential, one can generate nonclassical states of light. In particular, when $(aj<\u221212)$, it is possible to observe amplitude-squeezed states. It is possible to obtain also multi-mode squeezed states of different types.

Finally, we would like to note our approximation in the context of multi-electron systems. The generalization is straightforward,

This system can be further used to unify various multi-electron theories (for example, band structure theory) with quantized light description for the strong-field interaction problem. The noticeable difference with the one-electron consideration is that the fundamental scaling parameter here is $\beta \u223cNe\u22c5We/Wfield12$, where *N*_{e} · *W*_{e} is the total oscillating energy of charged particles. Although, in many relevant cases, it still can be very small, in some other cases, higher order approximations may be necessary.

## CONCLUSIONS

We present an approach that opens the possibilities to simplify and further study some quantum phenomena in interaction of light with complex quantum systems. This approach is analogous to the BO-approximation for the multi-system interactions of different nature, including possible non-local interactions. The general scheme can be further improved by additional subsequent steps beyond the considered first-order approximation. Our approach can be applied to the general light–matter interaction problem and the various related particular problem. We obtain a system of equations connected parametrically, which allows us to simplify the problem significantly. We demonstrate multi-mode and light-electron entanglement during the interaction and the possibility to generate nonclassical states of light. Our approach can explain recent experiments with intense interacting light^{5,6} and shed light to some related phenomena toward strong-field quantum optics and beyond.

## ACKNOWLEDGMENTS

I.G. and S.G. gratefully acknowledge funding from the ERC consolidator grant QUEM-CHEM (Grant No. 772676).

## DATA AVAILABILITY

The data that support the findings of this study are available within the article.

### APPENDIX: COHERENT STATE DETAILS

We consider the following function:

where *q* is a parameter of the light wavefunction (light state internal “coordinate”), Θ is an arbitrary given real number (phase of the light wave, see later), *N*_{0} is an arbitrary given constant (average photon number, see later), and *ω* is an arbitrary given frequency of the light mode. This function is called quantum coherent state of the plane wave light mode in coordinate representation because of its following properties:

- It is an exact solution of the following Schrödinger equation:which describes the time-dependent states of a light plane wave with frequency(A2)$i\u2202t\Psi c(q,t)=\omega N\u0302\Psi c(q,t)=\omega 2q2\u2212\u22022\u2202q2\Psi c(q,t),$
*ω*propagating in vacuum. Note that this is a Schrödinger representation of Schrödinger equation. It is normalized: $\u222b\u2212\u221e+\u221e|\Psi c(q,t)|2dq\u22611$.

- Within the following definitions of the field operators,(A3a)$A\u20d7\u0302=x012\beta a\u0302ei\kappa z+a\u0302+e\u2212i\kappa z,$(A3b)$a\u0302=12q+\u2202\u2202q,$(A3c)$a\u0302+=12q\u2212\u2202\u2202q,$(A3d)$N\u0302=a\u0302+a\u0302+12=12q2\u2212\u22022\u2202q2,$we can obtain the following quantum-mechanically averaged values:(A3e)$\beta =c2\pi /\omega V,$(A4a)$\u3008q\u3009c=2N0cos\omega t\u2212\Theta ,$(A4b)$\u2202\u2202qc=\u2212i2N0sin\omega t\u2212\Theta ,$(A4c)$\u3008A\u20d7\u0302\u3009c=x0\beta 2N0cos\omega t\u2212\Theta \u2212\kappa z=x0A0cos\omega t\u2212\Theta \u2212\kappa z,$where(A4d)$\u3008N\u0302\u3009c=12+N0,$$\u27e8R\u0302\u27e9c=\u222b\u2212\u221e+\u221e\Psi c*(q,t)R\u0302\Psi c(q,t)dq.$

Finally, we note that $a\u0302\Psi c=e\u2212i\omega t+i\Theta N0\Psi c$, establishing a connection with the definition of the canonical coherent state in the Schrödinger representation.^{9}