Brigitte Bidégaray-Fesquet

Kole Keita
email: kole.keita@imag.fr

Brigitte Bidegaray-Fesquet

A nonlinear Bloch model for Coulomb interaction in quantum dots

Keywords: Maxwell-Bloch model, quantum dot, Coulomb interaction, Cauchy problem, Liouville model, positiveness properties

In this paper, we first derive a Coulomb Hamiltonian for electronelectron interaction in quantum dots in the Heisenberg picture. Then we use this Hamiltonian to enhance a Bloch model, which happens to be nonlinear in the density matrix. The coupling with Maxwell equations in case of interaction with an electromagnetic field is also considered from the Cauchy problem point of view. The study is completed by numerical results and a discussion about the advisability of neglecting intra-band coherences, as is done in part of the literature.

Introduction

Bloch model is a very common model to describe the time evolution of a system of electrons in different contexts such as gases of electrons, glasses, or crystals. The very classical case of gases and glasses involves isotropic media. The electrons are supposed to be localized and non interacting. Their behavior is averaged at the mesoscopic scale. This leads to relatively simple models where matter energy levels are quantized and labelled by integers. The case of crystals [START_REF] Besse | A Maxwell-Bloch model with discrete symmetries for wave propagation in nonlinear crystals: an application to KDP[END_REF] also involves integer indexed levels, but symmetries and directions in matter have to be taken into account.

Bloch model has also been extended to the description of quantum wells [HK96, KR92, HK09], and quantum dots [START_REF] Gehrig | Mesoscopic spatiotemporal theory for quantum-dot lasers[END_REF][START_REF] Bidégaray-Fesquet | Positiveness and Pauli exception principle in raw Bloch equations for quantum boxes[END_REF]. In these models, matter is described by the state of two species of particles (electrons and holes, or equivalently conduction and valence electrons). In quantum wells, energy levels 1 are indexed by vectors, which correspond to displacements over the underlying lattice. In contrast, the confinement of electrons in quantum dots leads to integer indexed levels like in gases, which often leads to consider quantum dots as pseudo-atoms, but this is a very raw vision. In particular, among other differences, electrons at the same mesoscopic location do interact directly via Coulomb interaction.

[BF10] is a preliminary paper that derives basic Bloch equations for two species of electrons (conduction and valence) only taking into account the free electron Hamiltonian and the interaction with a laser electric field. The aim of the present paper is to include properly Coulomb interaction in this model. Beyond the sole derivation of the model, we want in particular to study its mathematical properties. In the continuation of [START_REF] Bidégaray | Contributions à l'électromagnétisme dans le domaine temporel. Modélisation classique et quantique en optique non linéaire[END_REF][START_REF] Bidégaray | Introducing physical relaxation terms in Bloch equations[END_REF], we want to show that a certain number of properties are preserved through the time evolution, such as Hermicity and positiveness of the density matrix.

Outline

The outline of this paper is as follows. We devote Section 1.2 to the description of the basic Bloch model which does not include Coulomb interaction but only the free energies of the electrons and the action of an electromagnetic field. In Section 2, we derive the Coulomb Hamiltonian in terms of the conduction and valence operators. The associated Heisenberg equation is derived in Section 3, but it ends up with an open system of equations. The system is closed using the Wick theorem, and the final Bloch equation has a Liouville form, but is nonlinear. This nonlinearity does not allow to use previous literature directly and we prove anew the Hermicity, positiveness and boundedness results in Section 4. In Section 5, we show the impact of the Coulomb contribution in numerical results and also compare our model with a vanishing intra-band coherence model defined in [START_REF] Gehrig | Mesoscopic spatiotemporal theory for quantum-dot lasers[END_REF].

The basic Bloch model

Let us first recall the main results obtained in [START_REF] Bidégaray-Fesquet | Positiveness and Pauli exception principle in raw Bloch equations for quantum boxes[END_REF] and fix the notations.

Commutators and Heisenberg equation

Let A and B be two operators, we define their commutator by [A, B] = AB-BA and their skew-commutator by {A, B} = AB+BA. For an operator A, we define the associated observable A = Tr(S 0 A) by averaging with respect to the initial state density S 0 of the system. If the system is described by a Hamiltonian H, the time-evolution for this observable is given by the Heisenberg equation

i ∂ t A = [A, H] , (1)
where is the reduced Planck constant. When the observable is the density matrix, the Heisenberg equation of motion is called the Bloch equation.

Operators for quantum dots

A quantum dot is defined as a collection of conduction and valence electrons.

There is of course no conduction in quantum dots since the electrons are confined in every direction, but this terminology is useful to distinguish between the valence electrons -which are in fact the absence of holes in the valence band -and the free, but confined, electrons. For each species, energy levels are quantized and indexed by a set of integers, I c and I v , for conduction and valence electrons, respectively. For i ∈ I c , we define the creation and annihilation operators c † i and c i . Likewise, for valence electrons, we define the creation and annihilation operators v † i and v i . Property 1.

{c i , c † j } = δ i,j , {c i , c j } = {c † i , c † j } = 0, {v i , v † j } = δ i,j , {v i , v j } = {v † i , v † j } = 0
, where δ i,j denotes the Kronecker symbol. This implies in particular that c i c i = c † i c † i = 0, which means that it is impossible to create twice or annihilate twice the same electron. This is the Pauli exclusion rule: electrons are fermions. Of course any conduction operator commutes with any valence operator.

Observables for quantum dots

The observable we are interested in is the density matrix. It includes a conduction density matrix, which elements are the ρ c ij = c † j c i . This matrix is Hermitian and positive semi-definite. Its diagonal terms ρ c ii = c † i c i are also called populations and give the probability to find an electron in state i. The off-diagonal terms, ρ c ij , i = j are called (intra-band) coherences. Of course we also define a valence density matrix, which elements are the ρ v ij = v † j v i . Besides we are interested in inter-band coherences defined by ρ cv ij = v † j c i . The entries of these matrices are the variables of the Bloch equation. They are cast in a single density matrix

ρ = ρ c ρ cv ρ vc ρ v
where ρ vc = ρ cv * , which ensures that ρ is Hermitian and positive semi-definite.

Free electron Hamiltonians and interaction with a laser

The basic Bloch equation for quantum dot is derived in [START_REF] Bidégaray-Fesquet | Positiveness and Pauli exception principle in raw Bloch equations for quantum boxes[END_REF]. It takes into account two types of Hamiltonians in the Heisenberg equation for the density matrix, namely free electron Hamiltonians and interaction Hamiltonians with a laser field. The free electron Hamiltonians read

H c 0 = k∈I c c k c † k c k , H v 0 = k∈I v v k v † k v k ,
for conduction and valence electrons respectively. The coefficients c k and v k are the free electron energies associated with each electron level. The electron levels are described by the wave functions ψ c k and ψ v k , solutions to the free electron Schrödinger equation subject to the boundary conditions of the quantum dot (see [START_REF] Haug | Quantum Theory of the Optical and Electronic Properties of Semiconductors[END_REF]). The interaction with a laser characterized by its time-dependent electric field E(t) is described by the Hamiltonians

H Lc = 1 2 (k,l)∈(I c ) 2 (E(t) • M c kl c † k c l + E * (t) • M c * kl c † l c k ), H Lv = 1 2 (k,l)∈(I v ) 2 (E(t) • M v kl v † k v l + E * (t) • M v * kl v † l v k ), H Lcv = (k,l)∈I c ×I v (E(t) • M cv kl c † k v l + E * (t) • M cv * kl v † l c k ),
where the dipolar moment matrices are matrices with entries in C 3 and may be expressed in terms of the wave functions associated to each level:

M c kl = dr ψ c * l (r) er ψ c k (r), M v kl = dr ψ v * l (r) er ψ v k (r), M cv kl = dr ψ v * l (r) er ψ c k (r),
where e is the unsigned charge of the electron. Injecting these Hamiltonians in the Heisenberg equation, the basic Bloch equations can be cast in Liouville form

i ∂ t ρ = [V 0 (t), ρ]. (2)
In equation (2), V 0 (t) = V F +V E (E(t)) is a sum of a constant term V F stemming from the free energies collected in diagonal matrices

E c 0 = diag({ c i } i∈I c
) and

E v 0 = diag({ v i } i∈I v )
, and a time dependent term due to the interaction with the electric field:

V F = E c 0 0 0 E v 0 and V E (E(t)) = E(t) • M c E(t) • M cv E * (t) • M cv * E(t) • M v .
The scalar product of the electrical field and the dipolar moment matrix, is a scalar product in C 3 which yields a matrix with entries in C, with the same dimension as ρ, as is necessary to give a meaning to the right-hand side of (2).

Mathematical properties of the Liouville equation

Equation (2) clearly preserves the Hermitian structure of ρ. Its exact solution is

ρ(t) = exp - i t 0 V 0 (τ ) dτ ρ(0) exp i t 0 V 0 (τ ) dτ . (3)
This expression allows to prove a certain number of properties (see [START_REF] Bidégaray-Fesquet | Positiveness and Pauli exception principle in raw Bloch equations for quantum boxes[END_REF]). Let d = card(I c ) + card(I v ) be the total number of levels, given a positive semi-definite initial data ρ(0) ∈ M d (C) and a continuous electric field E(t),

• equation ( 2) is globally well-posed, i.e. there exists a unique solution ρ ∈ C 1 (R + ; M d (C)) that exists for all times t ≥ 0 and which depends continuously on the data (parameters, initial data);

• for all time ρ(t) is a positive semi-definite matrix;

• its trace is conserved through the time evolution.

Second quantification Coulomb Hamiltonian

Coulomb interaction can be introduced using field operators. We denote by ψ † c (r) and ψ † v (r) the creation field-operators of respectively a conduction electron and a valence electron at the space location r, and ψc (r) and ψv (r) the corresponding annihilation field-operators. We consider that there are N relevant electrons in the quantum dot, and can write the Coulomb Hamiltonians as

H c-c = 1 2 N i,j=1 dr i dr j ψ † c (r i ) ψ † c (r j )V c (r i , r j ) ψc (r j ) ψc (r i ), (4a)
H v-v = 1 2 N i,j=1 dr i dr j ψv (r i ) ψv (r j )V v (r i , r j ) ψ † v (r j ) ψ † v (r i ), (4b)
H c-v = N i,j=1 dr i dr j ψ † c (r i ) ψv (r j )V c-v (r i , r j ) ψ † v (r j ) ψc (r i ), (4c)
where V c , V v and V c-v are the conduction-conduction, valence-valence and conduction-valence Coulomb potentials. The Coulomb potentials have the form

V c (r, r ) = V v (r, r ) = k C |r -r | and V c-v (r, r ) = - k C |r -r | ,
where k C is Coulomb's constant. The difference of treatment of conduction and valence electrons stems from the fact that Coulomb interaction describes the interaction of electrons and holes (see e.g. [START_REF] Haug | Quantum Theory of the Optical and Electronic Properties of Semiconductors[END_REF]) and that the presence of an electron in the valence band is indeed the absence of the corresponding hole, and vice-versa, which inverts the role of creation and annihilation operators. The total Coulomb Hamiltonian is

H C = H c-c + H v-v + H c-v .
We want to derive Bloch-type equations including the Coulomb interaction. Bloch equations have the advantage not to depend explicitly on the exact form of the field-operators. To this aim, we write

ψc (r) = α∈I c ψ c α (r)c α , ψ † c (r) = α∈I c ψ c * α (r)c † α , (5a) ψv
(r) = α∈I v ψ v α (r)v α , ψ † v (r) = α∈I v ψ v * α (r)v † α . (5b)
They are weighted by the conduction and valence electron wave functions ψ c α (r) and ψ v α (r), which are the same as those who occurred in the expression of the dipolar moment matrices. In the sequel, to avoid unnecessary written complexity, we will often omit to specify which set the indices belong to.

Inserting decompositions (5) in Hamiltonians (4) we obtain

H c-c = α1,α2, α 1 ,α 2 R c α1α2α 1 α 2 c † α1 c † α2 c α 2 c α 1 , (6a)
H v-v = α1,α2, α 1 ,α 2 R v α1α2α 1 α 2 v α 1 v α 2 v † α2 v † α1 , (6b)
H c-v = - α1,α2, α 1 ,α 2 R c-v α1α2α 1 α 2 c † α1 v α 2 v † α2 c α 1 , (6c)
where

R c α1α2α 1 α 2 = N 2 2 drdr ψ c * α1 (r)ψ c * α2 (r )V c (r, r )ψ c α 2 (r )ψ c α 1 (r), (7a)
R v α 1 α 2 α1α2 = N 2 2 drdr ψ v α 1 (r)ψ v α 2 (r )V v (r, r )ψ v * α2 (r )ψ v * α1 (r), (7b)
R c-v α1α2α 1 α 2 = -N 2 drdr ψ c * α1 (r)ψ v α 2 (r )V c-v (r, r )ψ v * α2 (r )ψ c α 1 (r). ( 7c
)
The symmetries in the integrands of (7) induce the following properties.

Property 2. Since V (r, r ) is an even function of r -r , variables r and r play the same role and

R c α1α2α 1 α 2 = R c α2α1α 2 α 1 , R v α1α2α 1 α 2 = R v α2α1α 2 α 1 . Since V (r, r ) is a real valued function R c α1α2α 1 α 2 = R c * α 1 α 2 α1α2 , R v α1α2α 1 α 2 = R v * α 1 α 2 α1α2 , R c-v α1α2α 1 α 2 = R c-v α 1 α 2 α1α2 * .
The Pauli exclusion principle (skew-commutation, Property 1) also induces that some terms in H c-c and H v-v are necessarily zero.

Property 3. If α 1 = α 2 or α 1 = α 2 , c † α1 c † α2 c α 2 c α 1 = 0, v α 1 v α 2 v † α2 v † α1 = 0.
The order of the operators (which has to be read from the right to the left) in the Coulomb Hamiltonians (4) has a meaning: in order that two particles interact via Coulomb interaction they have to pre-exist at locations r i and r j . Then they are annihilated while interacting and recreated at the same locations. Definition 1. A product of operators will be said to be in the normal order, if the annihilation operators are on the right and the creation operators on the left.

H c-c already follows a normal ordered form and we can keep it untouched:

H c-c = α1,α2, α 1 ,α 2 R c α1α2α 1 α 2 c † α1 c † α2 c α 2 c α 1 . (8a)
Hamiltonians H v-v and H c-v given by equations (6b) and (6c) do not follow normal ordered forms.

To express H v-v we need to compute a normal ordered form of

v α 1 v α 2 v † α2 v † α1 : v α 1 v α 2 v † α2 v † α1 = δ α1α 1 δ α2α 2 -δ α1α 1 v † α2 v α 2 -δ α 2 α1 δ α 1 α2 + δ α 2 α1 v † α2 v α 1 +δ α 1 α2 v † α1 v α 2 -δ α2α 2 v † α1 v α 1 + v † α1 v † α2 v α 2 v α 1 . Thanks to Property 2, R v α1α2α 1 α 2 = R v α2α1α 2 α 1 , therefore -δ α1α 1 v † α2 v α 2 and -δ α2α 2 v † α1 v α 1 lead to the same contribution. The same argument can be applied to δ α 2 α1 v † α2 v α 1 and δ α 1 α2 v † α1 v α 2 . Hence H v-v = 2 α,α ,β (R v βαα β -R v βαβα )v † α v α + α1,α2, α 1 ,α 2 R v α1α2α 1 α 2 v † α 1 v † α 2 v α2 v α1 . (8b)
In the definition of H v-v we have dropped the δ α1α 1 δ α2α 2 and -δ α 2 α1 δ α 1 α2 terms which would lead to zero contributions in the Heisenberg equation.

In the same way

c † α1 v α 2 v † α2 c α 1 = δ α2α 2 c † α1 c α 1 -c † α1 v † α2 v α 2 c α 1 , hence H c-v = - α,α ,β R c-v αβα β c † α c α + α1,α2, α 1 ,α 2 R c-v α1α2α 1 α 2 c † α1 v † α2 v α 2 c α 1 . ( 8c
)
3 Formulation of the Heisenberg equation

We now write Heisenberg equation (1) where A are operators c † j c i , v † j v i or v † j c i and H are the Hamiltonians defined by Equation (8).

Wick theorem

For the free Hamiltonians and the laser interactions, the computation of the commutators led to expressions in terms of the two-operator densities, and therefore to a closed set of equations [START_REF] Bidégaray-Fesquet | Positiveness and Pauli exception principle in raw Bloch equations for quantum boxes[END_REF]. This will not be the case any more here, since the commutators stemming from Coulomb Hamiltonians will give rise to four-operator densities. To go further we should a priori have evolution equations for these observables via the Heisenberg equation, computing commutators with the already defined Coulomb Hamiltonians. This would lead inevitably to six-operator densities, and so on. To avoid this endless procedure, we have to close the system at some point. This is the goal of the Wick theorem [START_REF] Wick | The evaluation of the collision matrix[END_REF] which amounts in our case to write the four-operator densities as sums of products of two-operator densities following e.g. the rule

c † α1 c † α2 c α 2 c α 1 (WT) = c † α1 c α 1 c † α2 c α 2 -c † α1 c α 2 c † α2 c α 1 ,
where the symbol

(WT)
= means "is approximated through Wick theorem".

Commutators involving c † j c i

In order to derive the Heisenberg equation, we have to compute

[c † j c i , H c-c ] and [c † j c i , H c-v ] (since [c † j c i , H v-v ] is clearly zero).

Commutator with H c-c

According to Equation (8a)

[c † j c i , H c-c ] = α1,α2, α 1 ,α 2 R c α1α2α 1 α 2 [c † j c i , c † α1 c † α2 c α 2 c α 1 ].
Remark 1. We already know many situations where [c † j c i , c † α1 c † α2 c α 2 c α 1 ] is necessarily zero:

• if none of the indices α 1 , α 2 , α 1 , α 2 is equal either to i or j,

• if α 1 = α 2 or α 1 = α 2 (see Property 3).
We compute separately each commutator:

[c † j c i , c † α1 c † α2 c α 2 c α 1 ] = δ iα1 c † j c † α2 c α 2 c α 1 + δ iα2 c † j c † α1 c α 1 c α 2 -δ jα 2 c † α2 c † α1 c α 1 c i -δ jα 1 c † α1 c † α2 c α 2 c i .
Using Property 2, we see that

[c † j c i , H c-c ] = α1,α2, α 1 ,α 2 R c α1α2α 1 α 2 2δ iα1 c † j c † α2 c α 2 c α 1 -2δ jα 1 c † α1 c † α2 c α 2 c i = 2 α, α 1 ,α 2 R c iαα 1 α 2 c † j c † α c α 2 c α 1 -2 α1,α2, α R c α1α2jα c † α1 c † α2 c α c i .
We now apply Wick theorem which leads to

[c † j c i , H c-c ] (WT) = 2 α, α 1 ,α 2 R c iαα 1 α 2 c † j c α 1 c † α c α 2 -c † j c α 2 c † α c α 1 -2 α1,α2, α R c α1α2jα c † α1 c i c † α2 c α -c † α1 c α c † α2 c i = 2 k, α,α (R c iαkα -R c iαα k ) c † α c α c † j c k -2 k, α,α R c kαjα -R c αkjα c † α c α c † k c i .
Defining matrix Λ c (ρ) as

Λ c ik (ρ) = 2 (α,α )∈(I c ) 2 (R c iαkα -R c iαα k ) ρ c α α , (9)
we can cast the result as

[c † j c i , H c-c ] (WT) = [Λ c (ρ), ρ c ] ij .

Commutator with H c-v

The same sort of computation as in Section 3.2.1 is performed on Equation (8c) to compute [c † j c i , H c-v ]. Using the matrices ζ v (ρ), γ c-v (ρ), and η c-v ik where

ζ v ik (ρ) = (α,α )∈(I v ) 2 R c-v iαkα ρ v α α , (10)
γ c-v ik (ρ) = - (α,α )∈I v ×I c R c-v iαα k ρ cv α α , (11)
η c-v ik = - β∈I v R c-v iβkβ , (12)
we obtain that

[c † j c i , H c-v ] (WT) = [ζ v (ρ) + η c-v , ρ c ] ij + k γ c-v ik (ρ)ρ vc kj - k ρ cv ik γ c-v * kj (ρ).

Commutators involving

v † j v i 3.3.1 Commutator with H v-v
According to Equation (8b), we have to evaluate two types of commutators to compute

[v † j v i , H v-v ].
The first commutator is clearly computed in the same way as [c † j c i , H c-c ] replacing conduction electron operators by valence ones. For the second part, we have to compute

[v † j v i , v † α v α ]. We obtain [v † j v i , H v-v ] (WT) = [Λ v (ρ) + κ v , ρ v ] ij ,
where Λ v (ρ) and κ v are defined as

Λ v ik (ρ) = 2 (α,α )∈(I v ) 2 (R v iαkα -R v iαα k ) ρ v α α , (13)
κ v ik = 2 β∈I v (R v βikβ -R v βiβk ). ( 14
)
3.3.2 Commutator with H c-v
According to Equation (8c), the commutator [v † j v i , H c-v ] only involves the first term of H c-v , which we can write as

[v † j v i , H c-v ] = α1,α2, α 1 ,α 2 R c-v α1α2α 1 α 2 c † α1 c α 1 [v † j v i , v † α2 v α 2 ].
The simplifaction of the commutators and the Wick theorem leads to

[v † j v i , H c-v ] (WT) = [ζ c (ρ), ρ v ] ij + k γ c-v * ik (ρ)ρ cv kj - k ρ vc ik γ c-v kj (ρ),
using the previously defined γ c-v (ρ) and the new matrix ζ c (ρ) defined by

ζ c ik (ρ) = (α,α )∈(I c ) 2 R c-v αiα k ρ c α α . (15)

Commutators involving v † j c i

The commutators involving v † j c i are all possible to write with already defined matrices, indeed we first compute [v † j c i , H c-c ] according to Equation (8a) and recognize

[v † j c i , H c-c ] (WT) = k Λ c ik (ρ)ρ cv kj .
Then Equation (8b) leads to

[v † j c i , H v-v ] (WT) = - k ρ cv ik Λ v kj (ρ) - k ρ cv ik κ v kj ,
and Equation (8c) to

[v † j c i , H c-v ] (WT) = k ζ v ik (ρ)ρ cv kj + k γ c-v ik (ρ)ρ v kj - k ρ c ik γ c-v kj (ρ) - k ρ cv ik ζ c kj (ρ) + k η c-v ik ρ cv kj .
evolution equation of ρ c ij , the ρ c ij term only involves E c 0 , and the other entries of ρ play a role through M c and M cv .

We can therefore easily analyze the Coulomb contributions in terms of shifts on the free electron energies and off-diagonal terms. Hence we compute

δ c i (ρ) = 2 (α,α )∈(I c ) 2 (R c iαiα -R c iαα i ) ρ c α α + (α,α )∈(I v ) 2 R c-v iαiα ρ v α α - β∈I v R c-v iβiβ , δ v i (ρ) = 2 (α,α )∈(I v ) 2 (R v iαiα -R v iαα i ) ρ v α α + (α,α )∈(I c ) 2 R c-v αiα i ρ c α α + 2 β∈I v (R v βiiβ -R v βiβi ).
We can define the energy shift matrices as δE c (ρ) = diag({δ c i (ρ)} i∈I c ) and

δE v (ρ) = diag({δ v i (ρ)} i∈I v ). Hence V (t, ρ(t)) = E(ρ(t)) + R(t, ρ(t))
, where

E(t, ρ(t)) = diag(E c (t, ρ(t)), E v (t, ρ(t))) and E c (ρ(t)) = E c 0 + δE c (ρ(t)) and E v (ρ(t)) = E v 0 + δE v (ρ(t)).

Vanishing intra-band coherence assumption

The model for quantum dots given in [START_REF] Gehrig | Mesoscopic spatiotemporal theory for quantum-dot lasers[END_REF] has clearly been derived mimicking quantum well models as in [START_REF] Kuhn | Monte Carlo simulation of ultrafast processes in photoexcited semiconductors: Coherent and incoherent dynamics[END_REF]. In quantum wells, electrons and holes can interact only if they "see each other" long enough, i.e. if they have (and are indexed by) the same wave vector. This leads morally to weakly coupled twolevel systems. Therefore, in [START_REF] Gehrig | Mesoscopic spatiotemporal theory for quantum-dot lasers[END_REF], the variables are the level populations and the inter-band coherences. In our model, this means that intra-band coherences ρ c ij and ρ v ij for i = j are not considered. We can therefore wonder what becomes of our model if we set intra-band coherences to zero. First we notice that in the evolution equation of, e.g., ρ c ij , there is, e.g., a contribution of V c ij (ρ c jj -ρ c ii ). Hence, even if intra-band coherences are initially zero, they are not always zero through the evolution with Equation (16). Setting artificially intra-band coherences to zero therefore destroys the Liouville structure of the system.

In our model, this hypothesis also changes the definition of Λ c (ρ), ζ v (ρ) and Λ v (ρ), which become

Λ c ik (ρ) = 2 α∈I c (R c iαkα -R c iααk ) ρ c αα , ζ v ik (ρ) = α∈I v R c-v iαkα ρ v αα , Λ v ik (ρ) = 2 α∈I v (R v iαkα -R v iααk ) ρ v αα .
4 Mathematical analysis

Estimates on the density matrix

In the same way as in the usual Bloch case [START_REF] Bidégaray | Introducing physical relaxation terms in Bloch equations[END_REF], the Liouville structure allows to state that the density matrix remains Hermitian through time evolution.

Its trace is conserved and it remains a positive semi-definite matrix. Hence for a given electric field E(t), the elements of the density matrix are bounded, more precisely populations are bounded

|ρ c ii (t)|, |ρ v jj (t)| ≤ Tr(ρ(t)) = Tr(ρ(0)),
for all i ∈ I c and j ∈ I v , as well as coherences

|ρ c ij (t)| ≤ ρ c ii (t)ρ c jj (t) ≤ Tr(ρ(0)) 2 , for all i, j ∈ I c ,
and likewise

|ρ v jk (t)|, |ρ c-v ij (t)| ≤ Tr(ρ(0)) 2 ,
for all i ∈ I c and j, k ∈ I v .

Coupling with the Maxwell equations

The density matrix governed by Bloch equations in Section 3 was only depending on time. We can now consider a collection of quantum dots which are scattered in space and interacting, not directly but through the interaction with an electromagnetic wave that propagates through the medium. This can be modeled by a density matrix, that now depends on time and space, which is coupled with Maxwell equations for the laser field through the expression of polarization. We therefore address the system

             µ∂ t H = -curl E, ε∂ t E = curl H -∂ t P, P = N b Tr(Mρ), ∂ t ρ = - i [V (ρ), ρ], (17)
where all the variables depend on time and space in 3 dimensions: the electric and magnetic fields E and H, the polarization P, and the density matrix ρ. In Equation ( 17), ε and µ denote the electromagnetic permittivity and permeability of the underlying medium. They both can depend on the space variable. The density of quantum boxes is given by N b . The Bloch and Maxwell equations are coupled via the polarization that involve the dipolar moment matrix M. Such models have already been written and studied mathematically and numerically in a few physical contexts. Here the specificity is the fact that the Bloch equation is nonlinear in ρ. Note that even in the case when V does not depend on ρ, the full coupled model is already nonlinear since V is affine in E.

This system can be cast in the abstract setting of [START_REF] Dumas | Cauchy problem and quasi-stationary limit for the Maxwell-Landau-Lifschitz and Maxwell-Bloch equations[END_REF]. In this paper, they introduce a general abstract setting able to treat both Maxwell-Landau-Lifschitz and classical Maxwell-Bloch equations. In this setting, the electromagnetic field is supposed to exist in all space R 3 . Matter described by the density matrix is only occupying a subdomain Ω of R 3 . The variables are gathered in one variable U = (u, v), where u = (u 1 , u 2 ) = (H, E) and v = ρ. The variable u can be viewed as 6 real variables and variable v as d 2 real variables. This variable U is supposed to be in

L 2 = L 2 (R 3 ; R 6 ) × L 2 (Ω; R d 2 ). The abstract system reads (∂ t + B)u = (κ -1 • l)F (v, u), for x ∈ R 3 , ∂ t v = F (v, u), for x ∈ Ω. ( 18
)
In this formulation κ(x) = (κ 1 (x), κ 2 (x)) = (µ, ε). It is uniformly positive as needed in [START_REF] Dumas | Cauchy problem and quasi-stationary limit for the Maxwell-Landau-Lifschitz and Maxwell-Bloch equations[END_REF]. Let H curl be the space of functions

f ∈ L 2 (R 3 ; R 3 ) with curl f ∈ L 2 (R 3 ; R 3 ). The linear differential operator B is defined on H curl ×H curl by B(u 1 , u 2 ) = (κ -1 1 curl u 2 , -κ -1 2 curl u 1 )
. The variable v is extended by v on the whole R 3 and is zero outside Ω. We can identify l 1 = 0, l 2 = -N b Tr(M•), and

F (v, u) = -i [V (ρ), ρ].
The system (17) verifies the hypotheses given in [START_REF] Dumas | Cauchy problem and quasi-stationary limit for the Maxwell-Landau-Lifschitz and Maxwell-Bloch equations[END_REF], namely,

• F is affine in u: F (v, u) = F 0 (v) + F 1 (v)u, • for j = 0, 1, F j (0) = 0, • for all R > 0 there exists C F (R) such that for all v ∈ B R (ball of radius R in R d ), |F j (v)| + |∂ v F j (v)| ≤ C F (R),
• there exists

K ≥ 0 such that for all (u, v) ∈ R 6 × R d 2 , F (v, u) • v ≤ K|u| 2 .
We have in particular used the L ∞ bounds of Section 4.1, more precisely, we look for

v ∈ L ∞ ((0, ∞); L ∞ (Ω; R d 2 )). (19)
Besides we suppose to have at time t = 0 the conditions div(κ j u j -l j v) = 0, for j = 1, 2, (

which are indeed the physical conditions div(µH) = 0 and div(εE + P) = 0 [START_REF] Dumas | Global existence for Maxwell-Bloch systems[END_REF]. The structure of Equation (18) ensures that this condition holds for all time if it is valid at the initial time.

Cauchy problem

In this section, we state without proof the results obtained in [START_REF] Dumas | Cauchy problem and quasi-stationary limit for the Maxwell-Landau-Lifschitz and Maxwell-Bloch equations[END_REF] and that we can apply to our context. The first result addresses the existence of global finite energy solutions.

Theorem 1 (Theorem 3, [START_REF] Dumas | Cauchy problem and quasi-stationary limit for the Maxwell-Landau-Lifschitz and Maxwell-Bloch equations[END_REF]). For any initial data

U 0 = (u 0 , v 0 ) ∈ L 2 (R 3 ; R 6 )× (L 2 (Ω; R d 2 ) ∩ L ∞ (Ω; R d 2 )) satisfying (20), there exists U ∈ C([0, ∞); L 2
) which is a solution to (18)-(20), and satisfies the finite energy condition (19). Moreover, for all T > 0, there exists a constant C that only depends on T , F , l and

v 0 L ∞ , such that for all t ∈ [0, T ], U (t) L 2 ≤ C U 0 L 2 .
To have a uniqueness result, we need some regularity on ε and µ. Usually in physical contexts, ε and µ may be discontinuous across the boundary of Ω, but we do not know how to tackle with this problem. We hence assume that Ω is bounded and

κ i -1 ∈ C ∞ Ω (R 3 ), for i = 1, 2, (21)
the space of C ∞ functions on R 3 with compact support included in Ω, which means in particular that κ i is 1 outside Ω and the transition across the boundary is smooth. We also assume that the initial data for the electromagnetic wave is smooth enough.

Theorem 2 (Theorem 5, [START_REF] Dumas | Cauchy problem and quasi-stationary limit for the Maxwell-Landau-Lifschitz and Maxwell-Bloch equations[END_REF]). Under the assumptions of Theorem 1 and (21), and assuming that curl u 0i ∈ L 2 (R 3 ) for i = 1, 2, there exists only one solution to (18)-( 20) with initial data U 0 , given by Theorem 1. It satisfies

curl u i ∈ C([0, ∞); L 2 (R 3 )) for i = 1, 2.
This stems from the fact that l 1 = 0 and F does not depend on u 1 but only u 2 .

5 Numerical experiments

Self-Induced Transparency

Self-Induced Transparency (SIT) is a typical two-level phenomenon: using a light pulse which is resonant with the transition, absorption and stimulated emission are combined to obtain exact population inversion and an unchanged electric field. This phenomenon has been predicted theoretically and confirmed experimentally [AE87, MH67, GS70].

The propagating field is a pulse given by

E(t, z) = E(t, z) sin(ω 0 (t -z/v)),
where v is the velocity of the pulse, ω 0 is both the center frequency of the pulse and the transition frequency of the medium, and E(t, z) is the pulse envelope. It is shown that the envelope is not reshaped by the medium, only if it is a symmetric hyperbolic secant

E(t, z) = E 0 sech t -z/v τ , where E 0 = 2 mτ .
In this expression, τ is the pulse duration and m = M/ , where M is the dipolar moment associated to the transition. According to the Area Theorem [START_REF] Mccall | Self-induced transparency by pulsed coherent light[END_REF], the medium undergoes k exact inversions if

A = m ∞ -∞ E(t, z)dt = kπ.
The corresponding pulse is called a kπ-pulse. It is easy to compute that

A = mτ E 0 [arctan(sinh t)] ∞ -∞ = mτ E 0 π,
and hence a kπ-pulse is obtained for the amplitude E 0 = k /M τ . For our test-cases, we will use a 2π-pulse, for which the medium is inverted and goes eventually to its original state. The return to the initial state is an easy-tocheck criterion to validate numerical approaches, as has been already done in [START_REF] Bidégaray-Fesquet | Hiérarchie de modèles en optique quantique. De Maxwell-Bloch à Schrödinger non-linéaire[END_REF][START_REF] Ziolkowski | Ultrafast pulse interaction with two-level atoms[END_REF].

Adaption to the quantum dot context

In this paper we want in particular to investigate the validity of the vanishing intra-band coherence assumption. To this aim we need a minimum of three levels and we therefore adapt the SIT experiment to our framework. We absolutely do not claim that SIT has any practical application in the quantum dot context, but only choose this test-case because of the easiness to interpret the results. In Figure 1(a), the original two-level test case is represented, for which there is a single conduction level and a single valence level, separated by the energy corresponding to the field frequency. The upper plot represents the (normalized) time-evolution of the electric field. The time-evolution of the population of the initially empty conduction level is given by the lower curve. We observe that the medium undergoes two complete population inversions.

In Figure 1(b), we have two three-level test cases with a single conduction level and two valence levels. In the first place we do not take Coulomb interaction into account. In the first test case (represented by solid lines both on the plot and on the scheme) the transition between the two valence levels is also resonant with the field and this destroys the SIT phenomenon. It suffices to get both valence levels far apart enough (e.g. 2ω 0 as in the second case, represented by dashed lines) to recover SIT. We use this last configuration as basis test case for the following numerical experiments.

Numerical features

The simulations are performed using a code based on a finite difference Yee scheme and a relevant choice for the time discretization of the Bloch equation (see [START_REF] Bidégaray | Time discretizations for Maxwell-Bloch equations[END_REF]). It allows to keep the good properties of the original Yee scheme: second order and explicitness. A splitting scheme, first described in [START_REF] Bidégaray | Introducing physical relaxation terms in Bloch equations[END_REF], allows to preserve positiveness at the discrete level. It is strongly based on the exact solution given by Equation (3). It has been adapted to include also Coulomb interaction, still preserving positiveness, but at the cost of a loss of approximation order, which becomes one.

Integrating the zero intra-band coherences assumption is a priori a problem since it destroys the Liouville structure and an exact solution is no more available. It is however possible to solve a Liouville-like equation and set artificially intra-band coherences to zero. This adds a step at each time iteration but allows to preserve the general structure of the numerical code.

To determine the right envelope amplitude for numerics, we use the argument of [START_REF] Ziolkowski | Ultrafast pulse interaction with two-level atoms[END_REF]: in practice the input pulse is cut off on an interval t ∈ [-10τ, 10τ ], therefore the numerical area is

A n = mτ E 0 [arctan(sinh t)]
10 -10 = mτ E 0 (0.999942π), which slightly changes the value of E 0 .

Results

Impact of Coulomb terms

To include Coulomb interaction in the SIT test case, we have to give values to the coefficients given by Equation (7). Their exact computation is not in the scope of the present paper. We choose to take them of the same order R 0 , taking into account the symmetries described by Property 2, but not equal (which would lead to vanishing Λ c , Λ v and κ v ). The test is performed using the full Coulomb terms (no vanishing intra-band coherence assumption).

For small values of R 0 (see the evolution of ρ c 11 described for R 0 = 10 -21 in Figure 2(a), solid plot) SIT is only slightly affected. In this figure, the dashed plot corresponds to the reference case (R 0 = 0) and is the same as the dashed plot of Figure 1(b). The effect is clearer for stronger values of R 0 (e.g., R 0 = 3 × 10 -21 in Figure 2(b)). The total inversion is prevented by Coulomb interaction. Although inversion is not complete, the return to zero of ρ c 11 is observed in this test case. We notice that the first valence level (which is not supposed to take part in the SIT experiment) is slightly populated, which makes this test case not really a two-level experiment.

Impact of vanishing intra-band coherences

We first test the impact of the vanishing intra-band coherence assumption on the Coulomb-free model. We always use the same experimental setting (see Figure 3 The result is displayed in Figure 3(b). The final equilibrium state for matter is slightly changed and ρ c 11 , which is given by the solid curve, does not eventually return to zero. Inversion is not total. Now we combine both Coulomb interaction and the vanishing intra-band coherence assumption. If R 0 is low (e.g., 10 -21 ), the result is not much affected by this assumption and is essentially the same as that plotted on Figure 2(a).

In the case when R 0 = 3×10 -21 , we see in Figure 4 that the final equilibrium state is not physical (ρ v 11 < 0). This is due to the destruction of the Liouville structure. We can easily explain why this is not observed when R 0 is low. When we are close to the SIT experiment we have a typical two-level phenomenon: ρ v 11 0 during the whole experiment. For a two-level system, the positiveness of the density matrix is equivalent to

• the positiveness of each diagonal term (populations),

• the estimation of coherence by populations, here: |ρ c-v 12 | 2 ≤ ρ c 11 ρ v 22 , if the second valence level would be the only relevant one (see [START_REF] Bidégaray | Introducing physical relaxation terms in Bloch equations[END_REF]).

Setting intra-band coherences to zero within the numerical process does not affect these properties, and the iteration used in the proof of the positiveness of the density matrix applies.

But for a three-level system (and the case R 0 = 3 × 10 -21 is a true threelevel case), the positiveness of the matrix involves some more properties, which are affected by setting ρ v 12 to zero. Although trace is still conserved, the positiveness of the population is affected. We would of course have the same result with a dedicated code where intra-band coherences would simply not be computed. Besides the effect is clear enough in Figure 4 not to be attributed to simple round-off errors (there are only 600 time-steps in this computation). The conclusion is that even if intra-band coherences seem not to be very relevant for some physical applications, it is very important to include them in the mathematical description and in the numerical computation to keep the natural mathematical structure of the density matrix. Remark 2. In absence of electromagnetic field the evolution equation for the conduction electrons is reduced to i ∂ t ρ c 11 = 0. This is only due to the fact that there is only one conduction level in our test case, and does not depend on the intra-band coherence vanishing assumption or on specific values of the Coulomb coefficients. Hence, when the population of the conduction level has been set into a non-physical state (and this is due to the intra-band coherence vanishing assumption), it remains in this state for ever.

Conclusion

In this paper, Bloch-type equations have been derived considering Coulomb effects in quantum dots. We have shown analytically and numerically that Coulomb effects are not negligible in some quantum dot structures, and we have given the link between mathematical properties and physical relevancy of the Bloch model and more specifically in the treatment of intra-band coherences. Then this model has been coupled with the description of laser propagation in the quantum dot structures, leading to a Maxwell-Bloch system for which we have studied the global Cauchy problem. This system has been implemented numerically and simulations have been performed on a self-induced transparency test-case. In particular, we have tested the impact of Coulomb parameters and intra-band coherences. We have illustrated numerically that the modification of the equation structure when intra-band coherences are neglected can lead to non-physical solutions. Further work will include additional effects in the same Bloch-type framework.

Figure 1 :

 1 Figure 1: Adaption of the SIT test case to the quantum dot context. (a) Original two-level case; (b) 2 three-level test cases.

  (a) R 0 = 10 -21 (b) R 0 = 3 × 10 -21

Figure 2 :

 2 Figure 2: Impact of Coulomb terms for two interaction strengths. (a) R 0 = 10 -21 -comparison with the Coulomb-free case; (b) R 0 = 3×10 -21 -excitation of the conduction level and of the first valence level.

Figure 3 :

 3 Figure 3: Impact of vanishing intra-band coherences the Coulomb-free model.

Figure 4 :

 4 Figure 4: Impact of vanishing intra-band coherences on the full model for R 0 = 3 × 10 -21 .

Acknowledgements

LJK is partner of the LabEx PERSYVAL-Lab (ANR-11-LABX-0025-01) funded by the French program Investissement davenir. The authors wish to thank Eric Dumas for fruitful discussions.

Matrix formulation

We would like to cast the former results as

where

Identifying the coefficients computed in Sections 3.2, 3.3, and 3.4, we obtain that the result can indeed be cast as [V C (ρ), ρ] where

and the various matrices have been defined by Equations ( 9) to (15). The evolution equation for the density matrix including also the free electron Hamiltonians, the interaction with a laser field, and Coulomb interaction can therefore be cast in a Liouville form

where V (t, ρ(t)) = V 0 (t) + V C (ρ(t)) and V 0 (t) has been introduced in Equation (2).

Energy shifts

Recall (see Equation ( 2)) the basic Bloch equation

where

The symmetries in the definition of M c and M v imply that their diagonal entries are zero. On the contrary for E c 0 and E v 0 there are only diagonal entries. In the


J Rosato

H Capes

R Stamm

Divergence of the Stark collision operator at large impact parameters in plasma spectroscopy models

Keywords: numbers: 32.70.Jz, 32.60.+i, 52.25.Dg, 52.70.Kz

The divergence that occurs at large impact parameters in Stark collision operators is examined for low-density hydrogen plasmas. In a previous work [Rosato, Capes, and Stamm, Phys. Rev. E 86, 046407 (2012)], we showed that the correlations between a radiating atom and the charged particles surrounding it affect the mean evolution of the atom, resulting in a mitigation of the Stark broadening near the line center. In this work, we examine the physical mechanism underlying this mitigation with an approach inspired from the standard semi-classical impact model. Our approach accounts for the atom-perturber correlations in a simple fashion, through a cut-off at large impact parameters, and embraces the impact model in the weakly coupled plasma limit. Comparisons with numerical simulations are performed and indicate a good agreement.

Recent theoretical works carried out in the framework of magnetic fusion research have led to successive improvements on ion collision operator models, devoted to provide fast numerical routines for hydrogen line shapes, for diagnostic applications but also for radiative transfer simulations (e.g. [1]). Although the ion Stark broadening is strongly dynamic for lines with a low upper principal quantum number n (such as Ly-α, D-α), the impact approximation for ions is only marginally valid and its use can lead to significant mistakes, with an overestimate of the line width at large densities or low temperatures (typically when N e,i > ∼ 10 14 cm -3 and/or T e,i < ∼ 1 eV, e.g. [2]). The inadequacy of the impact approximation rests in several assumptions that are not fulfilled. This is sometimes referred to as the "static effects". One such an effect concerns the reduction of the atom-perturber interaction to a duration of the order of |∆ω| -1 (with ∆ω being the frequency detuning), which is effective as soon as |∆ω| exceeds the plasma frequency ω p of the perturbing species. The latter corresponds to the characteristic time scale for an individual collision given the Debye shielding. Incomplete collisions have been investigated in the past for the electron broadening by using refined models for the collision operator, either based on kinetic theory (such as the "unified theory" [3,4]) or semiempirical procedures (e.g. using the Lewis cutoff [5]). Such models reproduce the result of the impact approximation at the line center (i.e., with the frequency detuning ∆ω = 0) and yield an asymptotic behavior identical to that expected within the quasi-static approximation, which assumes a constant electric field. Another failure of the impact approximation concerns the line center in moderately coupled plasma conditions. The typical time for the atom's dipole decorrelation (estimated as the inverse of the collision operator's characteristic matrix elements) can be shorter than ω -1 p , yielding a re- * Electronic address: joel.rosato@univ-amu.fr duction of the collision duration. Recently, we have examined this issue within a kinetic theory treatment inspired from the Bogoliubov-Born-Green-Kirkwood-Yvon hierarchy [6]. The reduction of the collision duration is a feature of emitter-perturber correlations and involves a new time scale which is in competition with the collision duration. In this Brief Report, we examine the physical mechanism underlying the emitter-perturber correlations with an approach inspired from the standard semiclassical impact model [7,8].

We consider the mean evolution of an ensemble of emitters under the influence of the ion microfield. The ions are assumed statistically independent to each other, independent with respect to the electrons, and the electron perturbation is neglected. A "collision" refers to the perturbation of the atomic wavefunction due to one ion only. In the standard semi-classical impact model, it is customary to write down a finite difference equation for the atomic evolution operator and to average it along the following scheme

{U (t + ∆t, 0)} = {U (t + ∆t, t)U (t, 0)} ≃ {U j (t + ∆t, t)}{U (t, 0)}. (1)
Here, U (t 2 , t 1 ) stands for the evolution operator of the emitter from time t 1 to t 2 , the brackets {...} denote average with respect to the ions, the time interval ∆t is chosen sufficiently large so as to include one collision only, and U j refers to the evolution of the emitter under the influence of this collision, which is labeled by j (this includes the time of closest approach, the impact parameter, and the velocity). The average factorization stems from the independent ion assumption. In the standard impact model, the collisions are assumed complete during ∆t, so that the evolution operator U j (t + ∆t, t) is replaced by U j (+∞, -∞). This quantity is identical to the S-matrix used in the theory of collisions, hence, we will refer to it as S j in the following. The finite difference Eq. ( 1) for the mean evolution operator {U (t)} becomes a differential equation at the limit ∆t → 0

d{U (t)} dt = (-iL 0 -K){U (t)}, (2)
where L 0 is the atomic Liouvillian and K = ∑ j ν j (1-S j ) is the collision operator, ν j being the frequency of collision of type j. If the Stark effect is linear, the sum over collisions diverges logarithmically at large impact parameters. This result is a feature of the Coulomb interaction and is especially relevant for hydrogen lines. The usual way to prevent this divergence consists in putting an upper cut-off at the Debye length, accounting phenomenologically for the screening of the Coulomb field due to the other perturbers (which is not retained in the independent perturber assumption). There is another regularization that prevents the divergence, which stems from the fact that each collision cannot be strictly completed, even at the line center. The collisions yield a decorrelation of the atomic dipole over a finite time scale t i (estimated as the inverse line width and usually referred to as the "time of interest"), which provides a maximum time for the effective action of each collision. Although this time can be much larger than the characteristic duration of each collision (of the order of ω -1 p ), the dipole decorrelation should be retained in a realistic calculation. An illustration is given in Fig. 1. The microscopic electric field has been simulated at plasma conditions relevant to tokamak edge, during the time of interest for Ly-α, assuming the ions moving along straight lines and generating a Debye electric field (due to the perturber-perturber correlations). For simplicity no magnetic field is considered here. The figure shows the dipole autocorrelation function obtained from numerical integration of the time-dependent Schrödinger equation, averaging with 1000 realizations of the electric field. Also shown in the figure is the effective duration of a collision, only slightly shorter than the time of interest.

In order to retain this regularization we extend the semi-classical impact model, by relaxing the complete collision assumption. We replace the S-matrix by the evolution operator U j (t j , -t j ) that accounts for one collision during a finite time. This time is described as a random variable distributed according to an exponential distribution of parameter γ, which is assumed of the order of the line width and is evaluated as some relevant matrix element of the collision operator (e.g., γ ≃ ⟨200|K|200⟩ for Ly-α). With this model, the collision operator reads

K = N ∫ +∞ 0 dvf (v)v ∫ +∞ 0 db2πb × ∫ λD 0 dtγe -γt {1 -U (t, -t)} angle . (3)
Here, the j-index has been removed to clarify the notation, and the sum has been written explicitly in terms of integrals over the velocity (module) v and the impact parameter (module) b. The brackets denote an angular average, and N , f stand for the ions' density and their velocity distribution function, respectively. To illustrate the regularization we evaluate explicitly the cross section

operator σ = ∫ +∞ 0 db2πb ∫ λ D 0
dtγe -γt {1 -U (t, -t)} angle and focus on the weak collision contribution, i.e., that provided by the second order expansion of the evolution operator. The calculation is similar to that performed in the standard impact model and it leads to the following analytical expression

σ = 4πh 2 3m 2 e v 2 ⃗ r 2 Λ(x)| γbst/v γλD/v , (4)
where ⃗ r is the position operator of the atomic electron (atomic units) restricted to the Liouville subspace of the line under consideration, b st is the strong collision radius, of the order of the Weisskopf radius b W = hn 2 /m e v, and Λ is defined as an integral

Λ(x) = ∫ +∞ x dy { 1 y + Im[e iy E 1 (iy)] } . ( 5
)
This function is plotted in Fig. 2. It behaves as -ln xγ E (with γ E being the Euler constant) at small x. The presence of the first exponential integral E 1 regularizes the integral for large x. The integral can be evaluated explicitly in terms of the sine and cosine integrals [START_REF] Abramowitz | Handbook of Mathe-matical Functions with Formulas, Graphs, and Mathematical Tables[END_REF]:

Λ(x) = -Ci(x) cos x -si(x) sin x. (6)
We have compared the model to the result of a simulation excluding strong collisions, i.e., the perturbers with an impact parameter smaller than b st have not been retained. For simplicity we have assumed all perturbers with the same velocity module, evaluated at the relative thermal velocity v i = √ 4T /m i (the factor 4 accounts for the reduced mass of the emitter-perturber system). Figure 3 shows the dipole autocorrelation function for N = 10 15 cm -3 and T = 1 eV, obtained from the simulation and from the model. In the latter case, it has been evaluated as ⟨210| exp (-N σv i t)|210⟩. The γ parameter has been estimated graphically from the simulation result, assuming an exponential decrease. Setting the equality ⟨d(0)d(γ -1 )⟩ = exp(-1) yields γ ≃ 4.2×10 10 s -1 . The model provides a good estimate of the dipole decorrelation, with only a small deviation (by about 20% regarding the decay rate N ⟨210|σ|210⟩v i ) from the simulation, which falls in the typical range of uncertainties inherent to the description of the Debye screening and the estimate of strong collision radius. In contrast, the impact approximation leads to a strong deviation at large times. The faster decrease of the autocorrelation function stems from an overestimate of the Stark perturbation due to the complete collision assumption. This result is in agreement with that obtained from the kinetic theory treatment presented in [6]. In summary, we have shown that the semi-classical impact approximation can be extended to account for incomplete collisions in a rather simple fashion, using an appropriate cut-off at large impact parameters. The resulting collision operator has a structure similar to that obtained within the usual impact model, making it suitable for implementation in line shape codes (e.g. [START_REF] Ferri | [END_REF][11][12]) and methods (e.g. [13][14][15]). The γ parameter involved in the cut-off can be estimated from numerical simulations or analytically. In the latter case a characteristic matrix element of the collision operator was previously suggested as a good candidate [6]. The model is of particular interest for lines affected by linear Stark effect, given their sensitivity to the weak collisions. It can be applied to the ion broadening of hydrogen lines, in the case where ion dynamics effects are strong, but also to the electron broadening. This is of particular interest for lines with a high principal quantum number, in tokamaks but also in astrophysics. The recent Workshop on Spectral Line Shapes in Plasmas [16] provides a motivation to proceed with further investigations from the theoretical point of view.

FIG. 1 :

 1 FIG. 1: Atomic dipole autocorrelation function of hydrogenLy-α obtained from a simulation, at conditions relevant to tokamak edge plasmas. The collisions' effective duration is of the same order as the time of interest, suggesting an inadequacy of the complete collision assumption at the line center.

FIG. 2 :

 2 FIG. 2: Plot of the Λ function. It has a logarithmic behavior at small x.

FIG. 3 :

 3 FIG. 3: Atomic dipole autocorrelation function of hydrogenLy-α corresponding to perturbation by weak collisions. The complete collision assumption done in the impact approximation yields an overestimate of the Stark perturbation, which results in an overestimate of the decorrelation. The model proposed in this work is in a better agreement.

Acknowledgments

This work was carried out within the framework of the European Fusion Development Agreement and the French Research Federation for Fusion Studies. It is supported by the European Communities under the contract of Association between Euratom and CEA. The views and opinions expressed herein do not necessarily reflect those of the European Commission.