Time-dependent Dyson orbital theory

Although time-dependent density functional theory (TDDFT) has become the tool of choice for realtime propagation of the electron density r(t) of N-electron systems, it also encounters problems in this application. The first problem is the neglect of memory effects stemming from the, in TDDFT virtually unavoidable, adiabatic approximation, the second problem is the reliable evaluation of the probabilities P(t) of multiple photoinduced ionization, while the third problem (which TDDFT shares with other approaches) is the reliable description of continuum states of the electrons ejected in the process of ionization. In this paper time-dependent Dyson orbital theory (TDDOT) is proposed. Exact TDDOT equations of motion (EOMs) for time-dependent Dyson orbitals are derived, which are linear differential equations with just static, feasible potentials of the electron–electron interaction. No adiabatic approximation is used, which formally resolves the first TDDFT problem. TDDOT offers formally exact expressions for the complete evolution in time of the wavefunction of the outgoing electron. This leads to the correlated probability of single ionization P(t) as well as the probabilities of no ionization ( P(t)) and multiple ionization of n electrons, P(t), which formally solves the second problem of TDDFT. For two-electron systems a proper description of the required continuum states appears to be rather straightforward, and both P(t) and P(t) can be calculated. Because of the exact formulation, TDDOT is expected to reproduce a notorious memory effect, the ‘‘knee structure’’ of the non-sequential double ionization of the He atom.


I. Introduction
Real-time propagation of the electron density r N (x 1 ,t) provides valuable information about dynamics of N-electron systems in strong rapidly changing electric fields.It is applied to simulate multiphoton ionization of atoms, [1][2][3][4] while the related power spectrum |m(o)| 2 of the dipole moment m i ðtÞ ¼ Ð r i r N x 1 ; t ð Þdx 1 describes high-harmonic generation in atomic systems. 4,5lectronic dynamics with fixed nuclei serves as a first stage of the simulation of similar processes in molecules. 4,6Coupling to nuclear dynamics provides not only more complete simulation of the electronic dynamics but also the description of photoinduced molecular dissociation 7 and isomerisation. 8or real-time propagation of r N (x 1 ,t) a dynamical orbital approach, time-dependent density functional theory (TDDFT) 3,[9][10][11][12][13][14][15] has become the tool of choice.In TDDFT r N (x 1 ,t) is expressed as the sum of the densities of the occupied Kohn-Sham (KS) orbitals f j (x 1 ,t), which are obtained from the TDDFT equations of motion (EOMs) The driving force of evolution in (2) is the external potential v ext (x 1 ,t), which includes the contribution from the applied electric field, while the third term in the figure brackets is the Hartree potential.7][18] At present, direct construction of the complicated functional v xc [r N ;C N (0),F N (0)] appears to be not feasible.Then, in order to make the propagation with the EOMs (2) possible, the adiabatic approximation is employed, in which the xc potential is evaluated instantaneously as the ground-state functional v 0 xc [{f j (t)}] of the density r N (t) or the orbitals {f j (t)} at the time t v xc [r N ;C N (0),F N (0)](x 1 ,t) E v 0 xc [{f j (t)}](x 1 ,t).
Since not even for the ground state the exact functional v 0 xc [{f j }] is known, it itself is, usually, replaced with an approximate functional v 0,appr xc ({f j }) of some DFT approximation.The resultant adiabatic EOMs are coupled non-linear differential equations with respect to the time-dependent KS orbitals {f j (t)}, which are solved with one of the standard or especially designed techniques of numerical integration. 19he first problem of standard TDDFT is the lack of the so-called memory effects.The memory effect is the influence on the dynamics of an electronic system of the dependence of v xc [r N ;C N (0),F N (0)](x 1 ,t 1 ) evaluated at t 1 on the densities r N (t) obtained at earlier times t o t 1 .What is called the adiabatic approximation (4), is essentially ignoring this dependence on t o t 1 and just using r(t 1 ).This can seriously affect the quality of simulated electronic dynamics. 3,18,20,21The paradigmatic memory effect is the famous ''helium knee'' indicating an enhancement of the non-sequential double ionization of the He atom by several orders of magnitude over what sequential ionization models would predict. 3,4A second problem TDDFT is facing here is the reliable evaluation of the correlated timedependent probabilities P n (t) of multiple ionizations, which in the case of He reduces to estimation of the probabilities P 1 (t) and P 2 (t) of single and double ionization, respectively.Usually, these probabilities are estimated in TDDFT in a rather crude independent-particle approximation.Note, that time-dependent Hartree-Fock (TDHF) theory, which has been also applied to real-time simulation of electronic dynamics, [22][23][24][25][26] has these two problems as well.Then, because of this, standard TDDFT and TDHF approaches fail to reproduce the knee structure of P 2 (t) as a function of the intensity of an applied field. 3,4There is also a third problem with the application of TDDFT to photoionization, which it shares with all other theoretical approaches to photoionization.This is the important problem of reliable representation of continuum states of ejected electrons.][30][31][32] In spite of all these problems, useful information about electronic dynamics can be obtained already in a rather crude single-active-electron (SAE) approximation, 1,2,6 which provides a considerable simplification of the EOMs (4).SAE employs frozen static Hartree and xc potentials which turn (4) into uncoupled linear differential eqn (5) with respect to the time-dependent KS orbitals {f j (t)}, so that in SAE individual KS orbitals are propagated independently.
In this paper time-dependent Dyson orbitals are introduced and it is shown that these orbitals provide both a theoretically appealing description of the time-evolution of the one-electron state of the electron that is ejected and a practical means of simulating this process.The time-dependent Dyson orbital theory (TDDOT) formally resolves the memory and multiple ionization problems (first and second of the above mentioned TDDFT problems) in the general N-electron case.In the special two-electron case of the He atom the formalism becomes particularly simple and allows investigation of truncation errors (in the limitation of summations over the ion states) and progressive uncoupling approximations in the spirit of the SAE approximation, e.g.uncoupling the propagation of a specific DO, e.g. the DO d 0 corresponding to the ion ground state C NÀ1 0 (which at t E 0 is practically the 1s orbital) from other DOs.7][38][39][40][41][42][43] We note that these calculations were concerned with the static Dyson orbitals, the amplitudes of the overlap between the groundstate N-electron wave function C N and the cationic ground-and excited-state (N À 1)-electron wave functions C (NÀ1) They can be obtained from one-electron equations (the ''orbital Dyson equations''), which feature the energy-dependent selfenergy, but they can also be calculated straightforwardly from good wavefunctions for the neutral starting system and the ionized system according to eqn (6).Katriel and Davidson 44 derived energy-independent but coupled equations for the static Dyson orbitals.These static Dyson orbitals are actually very close to the occupied Hartree-Fock and Kohn-Sham orbitals. 45We stress that we are dealing in this paper with time-dependent one-electron states, which only at t = 0 coincide with the familiar Dyson orbitals.In Section II the timedependent DOs d j (x 1 ,t) are introduced, with the sum of the DO densities giving the time-dependent electron density r N (t).The EOMs of TDDOT are derived with the projection and partial integration of the time-dependent Schro ¨dinger equation for an N-electron wave function C N (t), which is expanded in terms of the eigenstates of the time-independent (N À 1)-electron Hamiltonian and the time-dependent DOs d j (x 1 ,t).In spite of being exact, these EOMs (similar to the SAE EOMs (5)) contain just static potentials of the electron-electron interaction, and this solves the memory problem mentioned above.In the case of a time-independent Hamiltonian the EOMs describe stationary temporal oscillations of the DOs with frequencies equal to the corresponding vertical ionization potentials (VIPs).
The TDDOT EOMs reduce in this case to the static one-electron equations for the ground-state DOs. 44In Section III application of TDDOT to real-time simulation of electronic dynamics is discussed and the matrix form of the TDDOT EOMs is presented.
The matrix elements of their static electron-electron potentials can be precalculated and later used at all times during propagation, which offers an efficient propagation technique.A particularly promising application of TDDOT is the simulation of the dynamics of photoinduced ionization.In the general N-electron case TDDOT offers a method for the calculation of the time-dependent correlated probability of single ionization P 1 (t) as well as the probabilities P n (t) of multiple ionization of n electrons, the second problem mentioned above.In the case of a two-electron system (the He atom) the calculation of P 1 (t) and P 2 (t) appears to be particularly straightforward.The cation states of the DO expansion reduce to hydrogenic orbitals, which could be extended with the wellknown hydrogenic continuum states of the He 2+ ion to obtain a set of cation states that affords the simulation of double ionisation.Then, TDDOT is expected to successfully reproduce the ''helium knee'' in the experimental ion yield curves. 46

II. EOMs of TDDOT
The main ingredient of TDDOT proposed in this paper, the time-dependent DO d j (x 1 ,t), is introduced in analogy to the stationary case 39,44,45 as the amplitude of the overlap between the time-dependent state C N (t) and C ðNÀ1Þ j e iE ðNÀ1Þ j t , the latter describes temporal oscillations of an eigenstate of the time-independent (N À 1)-electron Hamiltonian with C (NÀ1) j being the stationary cation state with the energy E (NÀ1) j .In (7) the wave function C N (t) is the solution of the timedependent Schro ¨dinger equation for the N-electron system i @C N x 1 ; x 2 ; . . .; with as initial state at t = 0 the time-independent ground state C N 0 with the energy E N 0 .The Hamiltonian H ˆN(t) in ( 9) contains the time-dependent external potentials v ext (x l ,t) together with the kinetic operators in its first sum, while its third sum combines the electron-electron interaction operators.We derive the TDDOT EOMs for the DOs (7) from the N-electron Schro ¨dinger eqn (9) with the expansion of C N (t) in terms of the stationary cation states C ðNÀ1Þ k e ÀiE ðNÀ1Þ k t and the DOs d k (x 1 ,t) At t = 0 this is the familiar time-independent ''Dyson expansion'' of the ground state wavefunction in ion states times Dyson amplitudes. 44When t -N the propagation of the wavefunction in time will lead to a superposition of ion states times the oneelectron state (orbital) describing the emitted electron, which is by definition the time-dependent DO d k (x 1 ,t).The probability to find the system (atom or molecule) ionized at time t to the stationary where n m is the norm of the DO d m .One often considers the probability of finding the system in ion state C ðNÀ1Þ m x 1 ; x 2 ; . . .; When we are just interested in the appearance of ion state C NÀ1 m we have to integrate over all possible one-electron states {c k (x 1 ,t)} which gives the same result as ( 12) where completeness of the set {c k (x 1 ,t)} has been used.The time-dependence of the d m (x 1 ,t) and their relative amplitudes, and thus the relative probabilities, are determined by the shape and magnitude of the applied perturbing field.In the so-called This journal is © the Owner Societies 2016 ''sudden approximation'' the probability of ionization to ion state C NÀ1 m is taken proportional to the norm of the corresponding Dyson orbital at t = 0, 37 This approximation is applied with high-energy ionization, but it would be possible, according to eqn (12), to obtain probabilities without such approximation, as a function of time and of the nature of the applied field, if we could describe the timedynamics of the time-dependent Dyson orbital we have defined.Suppose we have a system that is reasonably well described by an independent particle model both in the ground state and in the ionized state (a good Hartree-Fock molecule and simple orbital ionization to a frozen orbital (Koopmans) ion state or relaxed orbital ion state).In that case the t = 0 DO are just the occupied orbitals in the ground state determinant, and in the correlated wavefunction they will still be approximately equal to those orbitals. 45Then the time evolution of the time-dependent DO would describe the photoionization process as the evolution of a one-electron state from a ground state orbital, to the fully propagated orthogonalized continuum state d m (x 1 ,t -N) for the emitted electron.We will consider in the Discussion section how to treat multiple (or no) ionization probabilities.
In order to develop solvable equations for the time evolution of the DO {d k (x 1 ,t)} we introduce the following partitioning of the Hamiltonian (10) x 2 ; x 2 ; . . .
Substituting expansion (11) into the left hand side, these equations can be written in the form of EOMs for the DOs d j (t) with static potentials for the electron-electron interaction where D j (t) and D jk (t) are the time-dependent scalar functions coming from interaction of Dv ext (t) with r (NÀ1) and In ( 18)-( 20) r (NÀ1) j is the electron density of the cation state From ( 7) and ( 11) follows, that the sum of the densities of the DOs d j (t) yields the total electron density r N (x 1 ,t) of the wave function C N (t) Note the principal difference between the derivation of the EOMs in the proposed TDDOT and the TDDFT EOMs (2). 9 The Kohn-Sham (KS) orbitals of TDDFT represent the corresponding non-interacting system, which is introduced through the action functional.Then, the stationary action principle is applied 9 and the orbital EOMs are obtained with the functional differentiation of the action with respect to the time-dependent density.In contrast to this, the time-dependent DOs are defined in (7) through the wave functions of the time-dependent N-electron and time-independent (N À 1)-electron interacting systems, so that no auxiliary non-interacting system is required.
Then, the EOMs of TDDOT are obtained directly in ( 9)-( 18) from the corresponding projected and partially integrated N-electron time-dependent Schro ¨dinger equation, in which the wave function is expanded according to (11) and the Hamiltonian is partitioned according to (15).In the special case of the time-dependent Schro ¨dinger equation with a time-independent Hamiltonian its solution C N (t) exhibits temporal oscillations with the exponential prefactor only Then, from the DO definition (7) follows that d j (t) acquires in this case the simple complex exponential for the oscillations in time The frequency I j is equal to the vertical ionization potential (VIP), is the conventional stationary DO (6).Insertion of (26) in eqn (18), in which the scalar functions D j (t) and D jk (t) are absent because the Hamiltonian is timeindependent in this case, yields the stationary equations for the static DOs d j (x 1 ) 44 ( ) This provides an energy-independent potential instead of the energy-dependent self-energy operator in the traditional Dyson equation.It does couple, as do the time-dependent eqn (18), the different DOs by the coupling terms with r (NÀ1) jk .The EOMs (18) are the key equations of TDDOT and they constitute the main result of the paper.Their remarkable feature is that they contain just static potentials of the electronelectron interaction and they could be classified as coupled homogeneous linear partial differential equations.Contrary to the TDDFT EOMs (3), no adiabatic approximation is employed.Thus, the memory problem of standard TDDFT mentioned in the Introduction is not present in TDDOT.We also note that the equations not only hold for primary ion states (resulting from orbital ionizations), as in the case of (the SAE to) the timedependent KS orbital eqn ( 5), but are just as well applicable for Dyson orbitals corresponding to satellite ion states, resulting from an orbital ionization plus an excitation.The eqn (18) are completely general, the DOs that result are only determined by the applied field (the potential v ext (x,t)).An intriguing possibility is to study the relative probabilities of ionization to specific ion states as a function of the frequency of the applied radiation, which may sweep through resonance with the subsequent ion states, or alternatively to study the relative probabilities of the first set of ionizations (valence and subvalence) as a function of the energy of the radiation, e.g.ranging from UPS energies to typical XPS energies.A prerequisite is of course an accurate solution method (see Section III).
One may also consider simplification of the full set of eqn (18).For instance, one could work towards a type of SAE approximation by first simplifying the present TDDOT EOMs by replacing in eqn (18) for a particular d j (t), all other DOs d k (t), k a j with the time-independent DOs d k (x 1 )e iIkt , which correspond to the initial time-independent potential (retaining of course the exponential time factor of the stationary orbital).With this, eqn (18) are simplified to which are uncoupled inhomogeneous linear partial differential differential equations with the last sum representing the free terms.We may also call them single active electron (SAE) equations, like ( 5), but note that they are for the time-dependent Dyson orbital and incorporate the Coulomb potential of the density of the corresponding ion state C NÀ1 j , not the Coulomb potential of the ground state density r N (x).They also retain time-dependent potentials v ext (x 1 ,t) and D j (t).Possible application of the TDDOT EOMs to real-time simulation of the electronic dynamics will be discussed in the next section.

III. Discussion
For applications, propagation of the one-electron functions {d k (x 1 ,t)} on a 3D grid would give complete flexibility to the evolving Dyson orbitals.It is also possible to operate with the matrix form of the proposed TDDOT EOMs, which is obtained with the expansion of the time-dependent DOs d j (t) in a convenient time-independent orthonormal basis {f n } where C nj (t) are the elements of the DO evolution matrix C(t).This introduces the well known problem of the representation of a continuum wavefunction, into which the DOs evolve in the course of time, with a suitable basis.Solutions for continuum electron wavefunctions have been obtained with mixed basis sets of atom based functions (Gaussians or Slater-type orbitals) and plane waves, [28][29][30][31][32] and rather successfully, for small molecules, with very flexible B-spline representation of the radial part of spherical harmonic expansions. 6,27Inserting the expansion (29) in the TDDOT EOMs (18) and integrating it with f m *, one obtains their matrix form This journal is © the Owner Societies 2016 Here, h mn (t) are the matrix elements of the time-dependent oneelectron operator W j mn are the elements of the static potential of r (NÀ1) and W jk mn are those of the static potential of r (NÀ1) With ( 23) and ( 29), the matrix form r N mn (t) of the time-dependent density reads r N mn ðtÞ ¼ The form of the linear matrix EOMs (30) provides a perspective on the application of TDDOT to real-time simulation of electronic dynamics.This simulation includes the following stages: (1) At the preliminary stage the N-electron ground state C N 0 and a chosen number M 4 N of the ground and excited stationary (N À 1)-electron states C (NÀ1) j and their energies E N 0 and E (NÀ1) j are obtained with a standard ab initio configuration interaction (CI) or coupled cluster (CC) calculation.Also, the orbitals to be used in the propagation eqn (30) are to be determined.The actual basis to be used may be the {f n } of the chosen static basis set, that has for instance been used in the initio calculations, or it may be a basis of linear combinations of the primitive basis that diagonalizes a convenient one-electron hamiltonian (Hartree-Fock or Kohn-Sham), plus some set of continuum functions.The quantities obtained are used to calculate the initial evolution matrix C(0), the cation electron density r (NÀ1) j and transition densities r (NÀ1) jk as well as the matrix elements W j mn and W jk mn of the static electron-electron interaction potential.
(2) At the real-time propagation stage the evolution matrix C(t) is propagated through a number of time steps t n with a numerical solution of the matrix TDDOT EOMs (30).In order to accomplish this, the one-electron elements h mn (t n ), the phase , and the scalar functions D j (t n ), D jk (t n ) are evaluated at the time t n , while the pre-calculated matrix elements W j mn and W jk mn are used in (30).With this, the present TDDOT offers a viable computational scheme.Indeed, the expensive ab initio approach is employed only in the preliminary stage.Also, during the bottleneck real-time simulation of electronic dynamics only relatively inexpensive orbital propagation is performed, which requires re-calculation of simple one-electron integrals with the external potential, while the same pre-calculated matrix elements of the static electronelectron interaction potentials are used at all time steps t n .Since the EOMs (30) are linear equations for the evolution matrix elements C nj (t), one can expect a rather stable propagation according to these equations.Obviously, the exact TDDOT requires the full expansion (11) with an infinite number of time-dependent DOs d j (t).Of course, the application of the TDDOT scheme would be particularly profitable, if for a reasonably accurate simulation only a fairly small number M of the DOs would have to be retained.
The TDDOT EOMs (30) can be further simplified with the SAE approximation of eqn (28).Inserting the expansion ( 29) in (28), one obtains the following inhomogeneous TDDOT-SAE matrix equations where the last sum represents their free terms.Because the generic TDDOT-SAE eqn ( 28) are uncoupled, the matrix TDDOT-SAE EOMS (35) are partitioned into uncoupled sub-sets of equations for the independent propagation of the expansion coefficients of individual DOs d j (t).
Coupled to nuclear dynamics, the matrix TDDOT EOMs could be applied to simulate photoinduced molecular dissociation and isomerisation.In particular, such simulation of the induced dissociation of the H 2 molecule and isomerisation of the diimide N 2 H 2 molecule could be compared with the corresponding TDHF and TDDFT simulations of ref. 7 and 8. Note, that in the case of a two-electron system X (for instance, X = H 2 ) the TDDOT EOMs (30) are especially simple.Indeed, in this case the static cation state C (NÀ1) j in the DO expansion (11) reduces to the ground-state wave function (orbital) c j (x 1 ) of the one-electron cation X + .Then, the cation density r (NÀ1) j in the matrix element W j mn of (32) reduces to the density of the orbital c j (x 1 ), while the transition density r (NÀ1) jk in the matrix element W jk mn of (33) turns to the product of the corresponding orbitals.
The very structure of TDDOT provides an especially interesting perspective on its application to real-time simulation of the electron dynamics of ionization induced by a strong rapidly changing electric field.Indeed, the expansion (11) can be rewritten as a weighted sum of normalized functions Y N j (t) which are the products of the normalized DO where n j (t) is the DO norm The expansion (36)  In its turn, the normalized DO d j x 1 ; t ð Þ ffiffiffiffiffiffiffiffiffi ffi n j ðtÞ p can be interpreted as the state of the electron that is being ejected during ionization to the ion state C (NÀ1) j .In the course of the simulated ionization, the DOs experience gradual delocalization, with diminishing norm for the DOs starting out at t = 0 as practically Hartree-Fock or Kohn-Sham orbitals (we mentioned the identity of the DOs at t = 0 with the orbitals if the wavefunction is approximated as a determinant).For satellites the DOs start out with very low norm (norm zero within the determinantal approximation 45 ), and will acquire increasing norm.For long-time simulation the inclusion of sufficient continuum states of various energies in the basis set expansion (29) (the third problem mentioned in the Introduction) will be mandatory.We have mentioned the alternative option of representation of the TDDOT EOMs (18) on a threedimensional grid in 3D physical space with the fully numerical solution of these equations.
Note, that simulation of multiple ionization will require also the inclusion of unbound cation states in the DO expansion (11).The present TDDOT allows in the first place evaluation of the correlated probability P 1 (t) of a single ionization of an N-electron system at the time t.The expression for P 1 (t) reads 4 where the entire physical space R is subdivided into a chosen finite volume V and the outside space O, R = V + O.The electrons in the region O are assumed to be ejected, while those inside the volume V are considered to be bound.Inserting the expansion ( 11) in (39), one obtains the TDDOT expression for P 1 (t) The cation states C (NÀ1) j in ( 40) are available from the above mentioned stage (1) of the calculations, while the DOs d(t) are obtained with their propagation at the stage (2).If the cation states are continuum states, i.e. bound states of the +2 ion times a free electron wavefunction, corresponding to two-fold ionization (and similar for multiple ionization), they would at least have one electron in a continuum state with negligible amplitude in the volume V, so the integrals over V will yield zero.So with V a limited volume, sufficient to practically fully contain the bound ion states, the summation over cation states C NÀ1 j reduces to one over bound cation states.This is also implied when we are considering the probability of single ionization, where the cation remains behind in a bound state.Orthogonality between the bound cation states is obtained by integration over the volume V, and the probabilities of single ionization (P 1 ), no ionization (P 0 ), and multiple ionization (P 41 ) become where j runs over the set of bound singly ionized (+1) ion states.Note that the outer space O in eqn (41) has to be the space complementary to a space V in which the bound states are fully contained, at least to the extent that the integrals over the bound states in eqn (40) on account of the completeness of the set of free electron states {c k }.In the traditional sudden approximation this is already evaluated at time t = 0, i.e. for the initial Dyson orbital, and provides an estimate of the relative intensities of the various peaks in the spectrum in the case of very high energy radiation.Taking into consideration the different types DOs, this affords the differentiation between primary ionizations and satellites.A primary ionization is an ''orbital ionization'' characterized by an ion wavefunction which can be reasonably approximated by a single determinant with a single orbital removed and the others possibly relaxed but recognizable.Its DO resembles closely the Hartree-Fock orbital that is vacated (the Hartree-Fock orbitals are the Dyson orbitals in a frozen-orbital approximation to the ion states).For satellites (primary ion states plus an excitation which are close in energy to some other primary ion state) the Dyson orbital at t = 0 has, for correlated wavefunctions, the same shape as the orbital corresponding to the main peak to which it is a satellite, but with (much) reduced amplitude (the DOs are not normalized), see ref. 45 for discussion.It would be zero in the frozen orbital determinantal approximation for ground state and ion (there are no satellites in the Koopmans' approximation) but will be different from zero in a relaxed orbital and in particular in a correlated treatment of the ion state.The present formalism affords a study of these features in a photoionization spectrum either without approximation or with limitations of the number of ion states included in the coupling (eqn (18)) or in the free terms of the uncoupled treatment (eqn (28)).Also the effect of more or less sophistication in the treatment of the ion states (Koopmans' approximation, SCF with orbital relaxation, correlated) can be investigated.The present TDDOT formalism also covers the simulation of the modern very short and strong laser pulse experiments which allow for multiple ionization.Note, that the total number N e of electrons ejected from the volume V because of all kinds of single and multiple ionizations is given by the integral of the electron density r(x 1 ,t) over the outer region O. 4 Using the probabilities P n (t) of ejection of n electrons, we have for the expectation value of N e N e h i ¼ which is the integral over the outer space of the standard expression of the electron density written as a sum of squares of Dyson orbitals, now time-dependent.The electron density in the outer space can be larger than 1 if the summation over j in (42) includes ion states which are continuum states, i.e. ion states that describe a doubly (or higher) ionized bound ion and one (or more) outgoing electrons.The norms of the Dyson orbitals belonging to such continuum ion states give the probability to find doubly (or higher) ionized ions and two (or more) electrons in the outer space.Then, from ( 40) and ( 42) follows that TDDOT allows also the evaluation of the number N 41 of electrons ejected in multiple ionizations The summation over j (in fact, an integration) runs over all the unbound +1 ion states, i.e. states that have a bound multiply charged ion (n-fold ionized atom or molecule, n = 2. ..(NÀ 1)) and n À 1 other free electrons (the n-th outgoing electron is the one already described by d j (x 1 ,t)).We can also focus on one individual multiple ionization yield P n , with n one of the numbers in the range 2. ..(NÀ 1).TDDOT at least formally resolves for the general N-electron case the problem of the evaluation of such an individual multiple ionization yield, one of the problems mentioned in the Introduction.For such a specific probability P n , n 4 1 we need to take those terms in the sum over j 2 B in eqn (43) which refer to states with (n À 1) free electrons and a bound atomic or molecular ion state of (N À n) electrons, with charge +n.We call this the set (+n) of unbound ion states and write With this differentiation of the summation over the unbound states in eqn ( 43) into subsets we formally solve the problem of obtaining individual P n s.
In the special case of two-electron systems eqn (43) reduces to the expression for the individual probability P 2 (and P 0 can be obtained from the sum rule P 2 n¼0 P n ðtÞ ¼ 1).As was mentioned earlier in this Section, the cation states C (NÀ1) j reduce in this case to the one-electron states (orbitals) of the cation, c j (x 1 ).In the basis set expansion for the DO we would then need, as always, functions that can represent continuum wave functions, e.g. by inclusion of unbound orbitals in the basis set expansion (29), which should be Coulomb waves in the +1/r field of the cation.Grid-based propagation of the DOs would obviate this basis set problem, of course.For the multiple ionization a second set of continuum states are needed, namely the continuum states of the cation have to be included in the expansion (11).In e.g. the He case these are the well known unbound states of the one-electron hydrogenic system, 47,48 i.e. the Coulomb waves in the +2/r field of the bare nucleus.This will be sufficient for the rigorous description of electron ejection in the process of non-sequential double ionization.Since there is no adiabatic approximation, there is not a memory problem (the first problem mentioned in the Introduction for applications of TDDFT) and TDDOT is expected to correctly reproduce the ''helium knee'' structure.The defining formula (7) for the time-dependent DOs d j (x 1 ,t) leads to the master equations of motion for the DOs, the TDDOT EOMs (18), which describe the evolution in the course of time of the one-electron state of the electron being ejected from the initial static Dyson orbital (very similar to a Hartree-Fock or Kohn-Sham orbital) into the plane wave like state of the final ionized electron.Remarkably, the equations of motion only require static electron-electron interaction potentials, but the price to be paid is coupling of the various one-electron equations, and therefore accurate solutions of the different ion states that are included in the coupled equations are required.Importantly, the present TDDOT fully covers the memory effects and does not resort to the adiabatic approximation, thus resolving the first of the mentioned problems of conventional TDDFT.

IV. Conclusions
The master eqn (18) have been cast into the convenient form (30) of the EOMs for the evolution matrix C(t) in a basis, which can readily include the continuum functions.This ensures the adequate representation in TDDOT of continuum wave functions, into which the DOs evolve in the course of time.The linear form of the TDDOT EOMs suggests a stable propagation of the DO evolution matrix elements C nj (t).
Partitioning (36) of the generic function C N (t) into ionization ''channels'' makes the application of TDDOT to real-time simulation of the electron dynamics of ionization by a strong rapidly changing electric field especially promising.The present formalism covers the single ionization events of traditional photoelectron spectroscopy.The total probabilities P 1 (t) and P 0 (t) (eqn ( 41)) of single ionization and no ionization at a particular time t as well as the probabilities of individual primary ionizations and satellites are evaluated from the densities |d j (x 1 ,t)| 2 of the propagated DOs associated with the corresponding bound cation states.The TDDOT formalism also covers multiple ionizations by very short and strong laser pulses.The total number (43) of electrons ejected in all multiple ionizations and the individual probabilities P n (t) (eqn ( 44)) of n-tuple ionization are evaluated from the densities of the DOs associated with the relevant unbound cation states.This resolves the second of the above mentioned TDDFT problems.
Real-time simulation of the paradigmatic non-sequential double ionization of the He atom offers an efficient simple test for the present TDDOT formalism.In this two-electron case all wave functions employed in TDDOT reduce to well-known orbitals.These are the the orbitals of the discrete spectrum of the He + ion as well as the unbound Coulomb waves in the +1/r field of He + and those in the +2/r field of the bare helium nucleus.Due to the above mentioned resolution of the memory problem, TDDOT is expected to reproduce a notorious memory effect, the ''knee structure'' of P 2 (t) as a function of the intensity of an applied field.

kt
times a ''free'' electron in a plane wave (or orthogonalized plane wave, or Coulomb wave) c k (x 1 ,t) with wavevector k (energy |k| 2 /2), where H ˆNÀ1 is the time-independent (N À 1)-particle Hamiltonian with the time-independent potentials v ext (x l ), l = 2. ..N and eigenfunctions C NÀ1 k .The figure brackets in (15) contain all operators acting on the first electron, while Dv ext (t) is the time-dependent change of the external potential Dv ext (x,t) = v ext (x,t) À v ext (x).(16) Inserting the expansion (11) and the partitioning (15) in (9), multiplying by ffiffiffiffi N p C ðNÀ1ÞÃ j e iE ðNÀ1Þ j t , and integrating over the coordinates x 2 ,. ..,xN , one obtains the following equations

C
nk ðtÞ d mn D jk ðtÞ þ W jk mn h i ) :

N 4 1
ðtÞ X N n¼2 P n ðtÞn ¼ N e ðtÞ h iÀ P 1 ðtÞ ¼ In this paper time-dependent Dyson orbital theory (TDDOT) has been proposed for real-time simulation of the photoionization process, i.e. the full dynamics of electron ejection from an atom or molecule.To this end a time-dependent extension of the Dyson orbital (DO) concept is introduced.The time-dependent DOs are defined by projection of the generic time-dependent wavefunction C N (t) onto the the stationary cation states C (NÀ1) j

.
yield orthogonality to high precision.The present formalism then naturally covers the traditional photoelectron experiment which measures single ionization events.We may ask at any time t for the probability to detect an emitted electron in some free electron state c k (a plane wave or, more accurately, an orthogonalized Coulomb wave in the +1 field of the cation left behind), and the cation in some state C NÀ1 hc k (x 1 ,t)|d j (x 1 ,t)i.Of course, depending on the frequency (energy) of the ionizing radiation and the energy of the ion state, d j (x 1 ,t) will develop into free electron states of specific energy, and only corresponding c k will give a finite matrix element.We have noted that the total probability of finding the system in cation state C NÀ1 j x 2 ; . ..xN ð Þ e ÀiE NÀ1 j tregardless of the state of the ejected electron is obtained by summing (integrating) over all free electron states, giving Ð dx 1