J Bouvier

A Chalabaev

C Charbonnel

E Falize

C Michaut

S Bouquet

C Cavet

ANALYTICAL SOLUTIONS OF SPECIFIC CLASSES OF ASTROPHYSICAL RADIATING SHOCKS

In this paper we study specific classes of radiating shocks which are widely spread in astrophysical environments. We present more general solutions of their structure and proceed to the analytical determination of physical quantities.

Introduction

Radiating shocks play a crucial role in astrophysical environments (Mignone 2005) as well as in laboratory plasmas generated by powerful facilities [START_REF] Drake | High-Energy-Density Physics[END_REF]. In astrophysics, we can find them in the head of stellar and galactic jets (Blondin & Cioffi 1989), in the first phase of supernovae explosions, in the late phase of supernova remnants (Chevalier & Blondin 1995), in accreting systems such as the magnetospherical accretion of T Tauri (Gunther et al. 2007) or magnetic cataclysmic variables (Cropper 1990). These shocks are at the basis of several models which allow the interpretation of astronomical observations. For instance, they may explain the recent observations of X spectra (Gunther et al. 2007) of classical T Tauri: TW Hya, BP Tau and V4056Sgr. Thus, it is very important to understand the structure and stability of these shocks. The structure of radiating shocks depends on the nature of pre-shock and post-shock media [START_REF] Drake | High-Energy-Density Physics[END_REF]. In this paper, we will consider shocks for which radiative effects can be modelled by entropy losses. This is the more common shock that it is encountered in interstellar phenomena. We will derive a solution by assuming a power law model for the cooling function, which generalizes the five specific solutions we find in the literature (Chevalier & Imamura 1982). Firstly, we will recall the theoretical model and discuss the approximation. Secondly, we will present the analytical solution for the one-and two-cooling processes problems. Finally, we proceed to the evaluation of fundamental physical quantities that characterize this type of shock and compare them with the results we already have.

Approximations and theoretical post-shock medium modelling

We consider a plane-parallel (∂/∂y = ∂/∂z = 0) collisional shock (i.e., t ii << t dyn where t ii is the characteristic time of collision between ions and t dyn is the dynamical time) with a post-shock medium which can be defined by a single temperature model (i.e., t ei << t cool where t ei is the characteristic time of energy exchanged between electrons and ions and t cool is the radiative cooling time). Moreover, we suppose that the local gravitational field does not modify the shock structure (i.e., g * x s << v2 s where g * , x s and v s are respectively the gravitational field, the thickness of the cooling layer and the shock velocity). Finally, we assume that the shock is stationnary (i.e., t cool << t dyn which means that cooling effects are faster than dynamical ones). Thus, the equations that give the evolution of post-shock medium are [START_REF] Kylafis | [END_REF]:

d dx [ρv] = 0 d dx [ρv 2 + P ] = 0 v dP dx -γ P ρ dρ dx = -(γ -1)Λ(ρ, P, x) (2.1)
where x, ρ, v, P , γ, Λ are respectively the spatial coordinate, the density, the velocity, the pressure, the polytropic index and the cooling function. Although these equations are theoretically consistent, we must SF2A 2007 specify a equation of state in order to connect the microscopic phenomena to the cooling function. We consider:

P = ǫ[Z]ρ α T β where T, ǫ[Z]
, α and β are respectively the temperature, a function of ionization Z and two free exponants where we must impose γ(1 -β) = (αβ) in order to make sure that the entropy is preserved in the post-shock region. On the shock front we have the Rankine-Hugoniot relations which are satisfied: ρ

[x - s ] ρ 0 = (γ + 1)M 2 (γ -1)M 2 + 2 , v[x - s ] v s = (γ -1)M 2 + 2 (γ + 1)M 2 , P [x - s ] ρ 0 v 2 s = 2γM 2 -(γ -1) γ(γ + 1)M 2 (2.2)
where ρ 0 and M are respectively the pre-shock density and the Mach number (see the figure 1). It appears that a natural rescaling is the normalisation of spatial coordinate to the thickness of the cooling zone (ξ = x/x s ), the density in the post-shock region to the pre-shock one, the velocity to the shock velocity and the pressure to the ram pressure. In order to reduce the number of equations ( 2

[ξ] = ρ 0 v 2 s (a 1 η[ξ] + a 2 )
where a 1 and a 2 are two constants. This transformation is consistent with the Rankine-Hugoniot relations provided η

[1] = [2 + (γ -1)M 2 ]/[(γ + 1)M 2 ] and P [ξ] = ρ 0 v 2 s (1 + 1/γM 2 -η[ξ]
). The determination of the intermediate function is possible through the last equation of (2.1) which is given by:

(γ + 1)η -γ 1 + 1 γM 2 dη dξ = -(γ -1) x s ρ 0 v 3 s Λ(ρ, P, ξ) (2.3)
It is this equation that we examine now.

3 Analytical solutions to one-and two-power laws cooling functions problem First, we suppose that the cooling function writes: Λ[ρ, P, x] = Λ 0 ρ ǫ P ζ (x + X 0 ) θ where Λ 0 , ǫ, ζ, θ and X 0 are five constants. This form generalizes the optically thin case where θ = 0. The arbitrarity of ǫ and ζ is motivated by the fact that Λ ∝ κ P σT 4 where κ P is the Planck opacity which can be modelled by a power law at high temperature. Introducing this form in equation (2.3) and using the integral form of Gauss hypergeometric function (that we note F 21 [a, b; c; x]), we obtain an implicit form of the compression profile:

cste + (γ + 1)η ǫ+2 1 + 1 γM 2 -ζ Γ[ǫ + 2] Γ[ǫ + 3] F 21 ζ, ǫ + 2; ǫ + 3; η 1 + 1/[γM 2 ] - (3.1) γη ǫ+1 1 + 1 γM 2 1-ζ Γ[ǫ + 1] Γ[ǫ + 2] F 21 ζ, ǫ + 1; ǫ + 2; η 1 + 1/[γM 2 ] = κ 0 ln[ξ + X 0 /x s ]; θ = -1 κ0 θ+1 [ξ + X 0 /x s ] θ+1 ; θ = -1 (3.2) where κ 0 = (γ -1)x θ+1 s Λ 0 ρ ǫ+ζ-1 0 v 2ζ-3
s and cste is defined by a specific value of the implicit function. This general solution is an extension of the five specific solutions known in the literature (Chevalier & Imamura 1982). Now we consider the case with two optically thin cooling processes. The cooling function writes: Λ(ρ, P ) = Λ 0,1 ρ ǫ1 P ζ1 + Λ 0,2 ρ ǫ2 P ζ2 where we suppose, without loss in generality, that t cool,1 < t cool,2 . Thus, the development of inverse cooling function gives:

1 Λ(ρ, P ) = 1 Λ 0,1 ∞ n=0 (-1) n Λ 0,2 Λ 0,1 n ρ [ǫ2-ǫ1]n-ǫ1 P [ζ2-ζ1]n-ζ1 (3.3)
which corresponds to a superposition of an infinity of processes. Introducing this expression into (2.3) and using once more the integral representation of Gauss hypergeometric function leads to an implicit solution:

cste

+ ∞ n=0 [-1] n Λ 0,n (γ + 1)η ǫn+2 1 + 1 γM 2 -ζn Γ[ǫ n + 2] Γ[ǫ n + 3] F 21 ζ n , ǫ n + 2; ǫ n + 3; η 1 + 1/[γM 2 ] - γη ǫn+1 1 + 1 γM 2 1-ζn Γ[ǫ n + 1] Γ[ǫ n + 2] F 21 ζ n , ǫ n + 1; ǫ n + 2; η 1 + 1/[γM 2 ] = x s ρ 0 v 3 s ξ (3.4)
where

ǫ n = ǫ 1 -[ǫ 2 -ǫ 1 ]n, ζ n = ζ 1 -[ζ 2 -ζ 1 ]n and Λ 0,n = ([Λ 0,2 /Λ 0,1 ] n /Λ 0,1 )ρ -ǫn 0 [ρ 0 v 2 s ] -ζn
We can see this solution as the one-process solution with a correction due to the presence of a second process.

Analytical evaluation of fundamental physical quantities

We will evaluate two fundamental quantities which are the thickness (which corresponds to x s ) of the post-shock medium and the accreted column density (that we note Ξ). From the previous results, it is easy to show that:

x s = ρ 0 v 3 s Λ 0 (γ -1)ρ ǫ 0 [ρ 0 v 2 s ] ζ γ 2 + (γ -1)M 2 (γ + 1)M 2 ǫ+1 1 + 1 γM 2 1-ζ Γ[ǫ + 1] Γ[ǫ + 2] × F 21 ζ, ǫ + 1; ǫ + 2; 2 + (γ -1)M 2 (γ + 1)M 2 (1 + 1/[γM 2 ]) -(γ + 1) 2 + (γ -1)M 2 (γ + 1)M 2 ǫ+2 1 + 1 γM 2 -ζ Γ[ǫ + 2] Γ[ǫ + 3] × (4.1) F 21 ζ, ǫ + 2; ǫ + 3; 2 + (γ -1)M 2 (γ + 1)M 2 (1 + 1/[γM 2 ]) = ∆[γ, M, ǫ, ζ] × v s t cool
We can apply this result to magnetic cataclysmic variables in which v s = 2GM wd /R wd where G, M wd and R wd are respectively the gravitational constant, the white dwarf mass and the white dwarf radius. Furthermore, the accretion rate ṁ is given by ṁ = ρ 0 v s . The introduction of these relations in equation (4.1) with assumptions of strong shock (M → ∞) and γ = 5/3 provides a simple relation between x s , ṁ, M wd and R wd :

x s = 7.60 × 10 6 cm ṁ 4g.cm -2 .s -1 -1 M wd 0.5M ⊙ 3/2 R wd 10 9 cm -3/2 (4.2)
This relation is consistent with Wu et al. (1994) results. Now we can evaluate the accreted column. By Laming (2004). Our results are presented in figure 2 and numerical results in table 1. These results generalize Laming ones for any cooling function. We find the trivial result which is that, for a given radiating process, smaller the polytropic index is, larger accreted matter is.

Conclusion

In this paper we present the general analytical solution to the problem of astrophysical radiating shock in the case where we have one-or two-cooling processes. These results are very important to validate the input density profile which is at the basis of all stability studies. Our solutions are totally consistent with previous evaluation of several physical quantities. Moreover it makes possible to predict more general results. This study can be considered as the counterpart of studies relative to optically thick radiative shocks (Bouquet et al. 2000).

Fig. 1 .

 1 Fig. 1. Schematic representation of radiating shock with the definition of the different physical quantities and notations.

  .1) we introduce an intermediate function η[ξ] defined by v[ξ] = -v s η[ξ]. The first and second equations of (2.1) imply ρ[ξ] = ρ 0 /η[ξ] and P

Fig. 2 .

 2 Fig. 2. Representation of the accreted column density versus Mach number for bremsstrahlung cooling for several values of polytropic index. We observe the saturation which is comptatible with values of table 1.

Table 1 .

 1 Satured values (limit M → ∞) of the Ξ-parameter for different cooling processes

	(ǫ, ζ) / γ	5/3	7/5	4/3	1.1	1.06
	( -2.35, 2.50) 1.991 2.894	3.406 11.387 19.022
	( 1.50, 0.50) 6.967 10.296 11.961 35.288 57.509
	( 2.50, -0.50) 5.869 8.651 10.047 29.629 48.292
	( 2.00, 0.00) 6.286 9.273 10.769 31.756 51.752
	( 1.00, 1.00) 8.302 12.313 14.316 42.327 68.994
	( 3.00, -1.00) 5.585 8.232	9.560 28.208 45.983
	definition, we have Ξ =

xs 0 ρ(x).dx / xs 0 ρ 0 .dx . This estimation is achieved for the former five analytical solutions by

Département de Physique Théorique et Appliquée, CEA/DIF, BP 12 91680 Bruyères-le-Châtel, France

LUTH, Observatoire de Paris, CNRS, Université Paris Diderot ; 5 Place Jules Janssen, 92190 Meudon, France c Société Francaise d'Astronomie et d'Astrophysique (SF2A) 2007


Raphaël Côte

Luis Vega

Scaling-sharp dispersive estimates for the Korteweg-de Vries group

We prove weighted estimates on the linear KdV group, which are scaling sharp. This kind of estimates are in the spirit of that used to prove small data scattering for the generalized KdV equations.

The purpose of this short note is to give a simple proof of two dispersive estimates which are heavily used in the proof of small data scattering for the generalized Kortewegde Vries equations [START_REF] Hayashi | Large time asymptotics of solutions to the generalized Korteweg-De Vries equation[END_REF].

The proof of these estimates can be easily extended to other dispersive equations. Denote U (t) the linear Korteweg-de Vries group, i.e v = U (t)φ is the solution to

v t + v xxx = 0, v(t = 0) = φ, i.e. U (t)φ = e itξ 3 φ or (U (t)φ)(x) = 1 t 1/3 Ai x -y t 1/3 φ(y)dy,
where Ai is the Airy function

Ai(z) = 1 π ∞ 0 cos ξ 3 3 + ξz dξ. Theorem 1. Let φ, ψ ∈ L 2 , such that xφ, xψ ∈ L 2 . Then U (t)φ 2 L ∞ ≤ 2 Ai 2 L ∞ t -2/3 φ L 2 xφ L 2 , (1)
U (t)φU (-t)ψ x L ∞ ≤ Ct -1 ( φ L 2 xψ L 2 + ψ L 2 xφ L 2 ). (2)
Furthermore, the constant 2 Ai 2 L ∞ in the first estimate is optimal.

Remark 1. These estimates are often used with φ replaced by U (-t)φ : denoting J(t) = U (t)xU (-t), they take the form

φ 2 L ∞ ≤ Ct -2/3 φ L 2 J(t)φ L 2 , (3)
φψ x L ∞ ≤ Ct -1 ( φ L 2 J(t)ψ L 2 + ψ L 2 J(t)φ L 2 ). ( 4
)
Proof. Due to a scaling argument (and representation in term of the Airy function), we are reduced to show that

U (1)φ L ∞ ≤ C φ L 2 xφ L 2 ,
and similarly for the second inequality. Hence we consider (U (1)φ)(x) = Ai(x -y)φ(y)dy, and we recall that the Airy function satisfies

| Ai(x)| ≤ C(1 + |x|) -1/4 and | Ai ′ (x)| ≤ C(1 + |x|) 1/4 . Then |U (1)φ| 2 (x) = Ai(x -y)φ(y) Ai(x -z) φ(z) y -z y -z dydz = Ai(x -y)yφ(y) Ai(x -z) y -z φ(z)dz Hz →y (Ai(x-z)φ(z))(y) dy -Ai(x -z)z φ(z) Ai(x -y) y -z φ(y)dy -Hy →z (Ai(x-y)φ(y))(z) dz = 2ℜ Ai(x -y)yφ(y)H z →y (Ai(x -z)φ(z))(y)dy,
where H denotes the Hilbert transform (and with the slight abuse of notation 1 x for vp 1

x ). As H : L 2 → L 2 is isometric and hence continuous (with norm 1), and Ai ∈ L ∞ , we get

|U (1)φ| 2 (x) ≤ 2 Ai(x -y)yφ(y) L 2 (dy) H(Ai(x -•)φ)(y) L 2 (dy) (5) ≤ 2 Ai 2 L ∞ yφ L 2 φ L 2 . (6)
This is the first inequality. Let us now prove that the constant is sharp. First consider the minimizers in the following Cauchy-Schwarz inequality :

yψ(y)H(ψ)(y)dy ≤ yψ(y) L 2 (dy) ψ L 2 . (7)
There is equality if yψ(y) = λH(ψ)(y) for some λ ∈ . Then a Fourier Transform shows that ∂ ξ ψ(ξ) = λ sgn ξ ψ(ξ), hence ψ(ξ) = C exp(-λ|x|), or equivalenty, one has equality in (7) as soon as

ψ(y) = C 1 + (y/λ) 2 for some λ, C ∈ Ê. (8)
(Notice that all the functions involved lie in L 2 )

We now go back to (6). Let x 0 ∈ Ê where | Ai | reaches its maximum. Now as

Ai(x 0 ) = 0, let ε > 0 such that for all y ∈ [-ε, ε], | Ai(x 0 -y)| ≥ | Ai(x 0 )|/2
, and consider the sequence of functions

φ n (x) = √ n 1 + (ny) 2 ½ |y|≤ε Ai(x 0 -y) . Denote ψ n (y) = ½ |y|≤nε 1+y 2 . As Ai(x 0 -y)φ n (ny) = √ nψ n (ny), |U (1)φ n | 2 (x 0 ) = 2 y √ nψ n (ny)H z→y ( √ nψ n (nz))(y)dy = 2 n yψ n (ny)H(ψ n )(ny)dy.
One easily sees that ψ n (y) → 1 1+|y| 2 in L 2 and yψ n (y) → y 1+|y| 2 in L 2 , and hence, in view of (8), as H is homogeneous of degree 0 and L 2 isometric, we have

|U (1)φ n | 2 (x 0 ) ∼ 2 n y 1 + |y| 2 H( 1 1 + | • | 2 )(y)dy ∼ 2 n y 1 + |y| 2 L 2 1 1 + |y| 2 L 2 ∼ 2 y √ nψ n (ny) L 2 √ nψ n (ny) L 2 ∼ 2 y Ai(x 0 -y)φ n (y) L 2 Ai(x 0 -y)φ n (y) L 2 .
As φ n concentrates at point 0, we deduce

|U (1)φ n | 2 (x 0 ) ∼ 2| Ai(x 0 )| 2 yφ n (y) L 2 φ n (y) L 2 as n → ∞, (9)
which proves that the sharp constant in the first inequality is 2 Ai 2 L ∞ .

For the second inequality (estimate of the derivative), we have as for the first inequality :

(U (1)φU (1) ψx )(x) = Ai(x -y)φ(y) Ai ′ (x -z) ψ(z) y -z y -z dydz = Ai ′ (x -y)yφ(y) Ai(x -z)φ(z) y -z dz dy - z Ai ′ (x -y)φ(y) Ai(x -z)z ψ(z) z -y dz dy = Ai ′ (x -y)yφ(y)H z →y (Ai(x -z)φ(z))(y)dy + Ai ′ (x -y)φ(y)H z→y (Ai(x -z)z ψ(z))(y)dy. Denote ω x (y) = 1 √ 1+|x-y| : ω -1
x ∈ A 2 (with the notation of [START_REF] Stein | Singular integrals and differentiability properties of functions[END_REF]), so that there exists C not depending on x such that ∀v,

|Hv| 2 ω -1 x ≤ C |v| 2 ω -1 x .
Recall the well-know asymptotic

| Ai ′ (x)| ≤ C(1 + |x| 1/4 ). Then Ai ′ (x -y)yφ(y)H(Ai(x -•) ψ)(y)dy ≤ | Ai ′ (x -y)yφ(y)| 2 ω x dy 1/2 |H(Ai(x -•) ψ)(y)| 2 ω -1 x (y)dy 1/2 ≤ C yφ(x) L 2 | Ai(x -y) ψ(y)| 2 ω -1 y (y)dy 1/2 ≤ C yφ L 2 ψ L 2 .
In the same way,

Ai ′ (x -y)φ(y)H z→y (Ai(x -z)z φ(z))(y)dy ≤ | Ai ′ (x -y)φ(y)| 2 ω x (y) 1/2 |H z→y (Ai(x -z)zψ(z))(y)| 2 ω -1 x (y)dy 1/2 ≤ C φ L 2 | Ai(x -y)yψ(y)| 2 ω -1 x (y)dy 1/2 ≤ C φ L 2 yψ L 2 .
So that :

U (1)φU (1) ψx L ∞ ≤ 2C( φ L 2 xψ L 2 + ψ L 2 xφ L 2 ).
Up to scaling and replacing ψ by ψ, this is the second inequality.

Remark 2. This proof (especially (5)) is reminiscent of that in [START_REF] Lieb | Analysis[END_REF] (see also [START_REF] Duoandikoetxea | Some weighted Gagliardo-Nirenberg inequalities and applications[END_REF])

φ 2 L ∞ ≤ φ L 2 φ ′ L 2 ,
where the constant is sharp and minimizers are Ce -λ|x| . This has application to the Schrödinger group U(t) (i.e U(t)φ = e itξ 2 φ). We have the following Schrödinger version 4t ) :

of estimate (3) (notice that U(t)xU(-t) = e i|x| 2 4t it 2 ∂ x e -i|x| 2
ψ 2 L ∞ = e i|x| 2 4t ψ 2 L ∞ ≤ e i |x| 2 4t ψ L 2 e i |x| 2 4t ∂ x e i |x| 2 4t ψ L 2 ≤ 2 t ψ L 2 U(t)xU(-t)ψ L 2 .
From Theorem 1, we can easily obtain the optimal decay in a scaling sharp Besov like space. Let ϕ ∈ D(Ê) be non-negative with support in ]-2, 2[ and such that ϕ equals 1 in a neighbourhood of [-1.5, 1.5]. Denote ψ(x) = ϕ(2x) -ϕ(x) and ψ j (x) = ψ(x/2 j ). Finally introduce

φ Nt = ∞ j=-∞ 2 j/2 ψ j U (-t)φ L 2 .
Corollary 1. We have :

φ L ∞ ≤ Ct -1/3 φ Nt .
Proof. Notice that |xψ j (x)| ≤ 2 j+1 ψ j (x). As φ = j U (t)ψ j U (-t)φ, we have :

φ L ∞ ≤ j U (t)ψ j U (-t)φ L ∞ ≤ Ct -1/

≤

  Ct -1/3 φ Nt .

	Raphaël Côte
	Centre de Mathématiques Laurent Schwartz, École polytechnique
	91128 Palaiseau Cedex, France
	cote@math.polytechnique.fr
	Luis Vega
	Departamento de Matemáticas, Universidad del País Vasco
	Aptdo. 644, 48080 Bilbao, España
	luis.vega@ehu.es
			1/2
			L 2
	≤ Ct -1/3	1/2 L 2 xψ j U (-t)φ	1/2 L 2
	≤ Ct -1/3		1/2 L 2

3 j U (t)ψ j U (-t)φ 1/2 L 2 U (t)xU (-t)U (t)ψ j U (-t)φ j U (t)ψ j U (-t)φ j U (t)ψ j U (-t)φ 1/2 L 2 2 j/2 U (t)ψ j U (-t)φ


D Paget

T Amand

J. -P Korb

Light-induced nuclear quadrupolar relaxation in semiconductors

Keywords: 76.60. -k, 72.25. Fe, 78.55. Cr

Light excitation of a semiconductor, known to dynamically-polarize the nuclear spins by hyperfine contact interaction with the photoelectrons, also generates an intrinsic nuclear depolarization mechanism. This novel relaxation process arises from the modulation of the nuclear quadrupolar Hamiltonian by photoelectron trapping and recombination at nearby localized states. For nuclei near shallow donors, the usual diffusion radius is replaced by a smaller, quadrupolar, radius. If the light excitation conditions correspond to partial donor occupancy by photoelectrons, the nuclear field experienced by electrons trapped at shallow donors can be decreased by more than one order of magnitude.

I Introduction

In a semiconductor, the possibility to enhance the nuclear polarization by the hyperfine contact interaction with spin-polarized electrons generated by circularly-polarized light excitation is of interest both for fundamental reasons and, among others, for applications to: i) quantum computing, 1 ii) transfer of nuclear magnetization to biological systems, as an alternative to adsorption of polarized xenon, 2,3 iii) understanding of the fractional quantum Hall effect. 4 Further potential applications of the optical increase of NMR sensitivity include extension to nuclei of single spin investigations using magnetic resonance force microscopy at surfaces. 5 After the demonstration of optical nuclear polarization in silicon, 6 a number of recent investigations of the optically-enhanced bulk nuclear magnetization have been undertaken using standard NMR in Si, 7 GaAs, [8][9][10][11][12][13] InP, 14 CdTe. 15 Some of the results [11][12][13] were used to verify the predictions of a general theory for nuclear relaxation in solids according to which the presence of paramagnetic impurities, or localized centers, is crucial for relaxation of the nuclear spin system. [16][17][18] Nuclei close to the centers are relaxed by the hyperfine interaction with the spin-polarized photoelectrons trapped at these impurities, while the bulk nuclear spin system is relaxed by spin diffusion from the latter minority nuclei. A diffusion radius is defined corresponding to the distance from the impurity separating the two types of relaxation processes. 19 Optical detection of NMR, from the depolarization at resonance of the luminescence, was first reported for GaAlAs in 1974, 20 and subsequently applied to several III-V semiconductors, [21][22][23][24] as well as 2D systems 8,25 and quantum dots. 26,27 For bulk materials, this technique was shown to only detect nuclei near the sites of electronic localization, which verifies the existence of a diffusion radius. 28 The ratio of the nuclear hyperfine field acting on 14/04/2008 the electrons and of the optically measured electronic spin polarization is consistently smaller than its calculated value. The corresponding reduction of the nuclear field is found to be of 0.1 for GaAs, 29 0.02 for GaSb, 22 and of several percent for InP. 23 Such decrease is likely to significantly reduce the optical enhancement of the nuclear polarization.

The identification of the relaxation mechanisms responsible for this loss of nuclear polarization remains an open problem. In the absence of light excitation, the hyperfine coupling with the unpolarized holes, 30 or the quadrupolar interaction modulated by lattice phonons 31 are negligible at low temperature. The total hyperfine field of nuclei near shallow donors is decreased because of the competition between spin-lattice relaxation and spin diffusion, but only by a factor 3. 28 Another possibility is that the averaging caused by spin exchange between trapped electrons and free electrons reduces the effective nuclear field measured experimentally. 21 Interestingly, in addition to the dynamic nuclear polarization, light excitation also creates an intrinsic leakage mechanism for the same nuclei as the ones which are dynamically-polarized. The nuclei close to shallow donors experience a very strong electric field from the ionized donor. Since the latter field is modulated by trapping and recombination of photoelectrons, there results a significant nuclear depolarization.

The present work is devoted to an evaluation of the efficiency of such light-induced nuclear relaxation for the case of nuclei near shallow donors. In Sec. II, the characteristic time of the quadrupolar-induced evolution of the nuclear spin temperature is calculated using the semi-classical rate equation for the nuclear spin density matrix. [START_REF] Abragam | The Principles of Nuclear Magnetism[END_REF] Quantitative estimates of the nuclear magnetization as a function of distance to the donor and of the nuclear field experienced by electrons trapped at shallow donors are performed in Sec. III, using the known magnitudes of quadrupolar [START_REF] Gill | [END_REF][34][35][36] and hyperfine couplings. 29 Provided the light power density is such that shallow donors are partially occupied, the light-induced quadrupolar relaxation is found to induce a decrease of the nuclear hyperfine field by as much as one order of 14/04/2008 magnitude. The corresponding effect in quantum dots and the resulting dependence of the nuclear field as a function of temperature and light excitation power will be discussed elsewhere. 37 II Light-induced quadrupolar nuclear relaxation time and nuclear polarization value

In the absence of a trapped photoelectron, the electric field experienced by nuclei near a shallow donor is given by

2 0 1 . 4 ) ( r e r E off πεε = ( 1
)
where e is the electronic charge, ε is the static dielectric constant and r is the distance from the donor. Photoelectron trapping and recombination induces a modulation of the electric field between Eq. ( 1) and such that

) (r E on [ ) ( 1 ) ( ) ( r s r E r E off on - = ] (2)
where the expression for , found using Gauss's theorem and the shape of the electronic wavefunction, is

) (r s * 0 / 2 2 * 0 2 * 0 2 2 1 1 ) ( a r e a r a r r s - ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ + + - = (3)
Here * 0 a is the electronic Bohr radius. In GaAs, one has and is of the order of 10

3 . 0 ) ( * 0 ≈ a s ) ( * 0 a E off 6 V/m. The modulation amplitude ( ) ( ) r E r E on off -
induced by photoelectron trapping and recombination is very large. Unlike the usual quadrupolar relaxation, the corresponding relaxation process does not rely on phonons for modulation and can be relevant at low temperature. The present section is devoted to the calculation of the corresponding relaxation time and of the resulting decrease of the nuclear polarization. ( )

[ ∑ = + + = 2 1 . k Qk Qk Q Q A A r F H ] (4)
Taking the normal z to the surface as the quantization axis, the spin operators are given by

Qk A ( ) [ ]

Zeeman spin levels separated by energies given by

Qk A (k=1,2) (6) B k k γ ω h h =
where γ is the nuclear gyromagnetic ratio. It is convenient to write 38,39 (7) where e is the electronic charge, Q is the quadrupolar moment of the bare nucleus of spin I.

( ) ( ) ( ) r E b r E I I Q eR r F Q Q γ h = - = ) 1 2 ( 4 14
The factor R 14 , which includes the electrostatic antishielding, is in the present frame of 14/04/2008 coordinates Oxyz the value of the only nonzero components of the third rank tensor relating the electric field gradient to the electric field. [START_REF] Gill | [END_REF][34][35][36] The quantity is the ratio of a magnetic to an electric field. It is calculated in Appendix A for different compounds and is given in Table I. The Hamiltonian

H [ ] 1 14 ) 1 2 ( 4 - - = I I Q eR b Q γ h
Q can be rewritten as the sum of a static and of a modulated part

( ) [ ] ( ) ( ∑ + + + = Qk Qk Q Q A A r F t h H . 1 0 ) (8)
where is given by

) ( 0 r F Q [ ] [ ] ( ) r E b r s r F r s r F off Q t off Q t Q Γ - = Γ - = ) ( 1 ) ( . ) ( 1 ) ( 0 γ h ( 9
)
and Γ t is the fraction of the time during which the electron is present at the donor site. The function describes temporal fluctuations due to the trapping and recombination of an electron at the localized site. This function has a time average equal to zero and varies randomly between and -

( ) t h [ ] 1 ) ( 1 ) ( - Γ - Γ t t r s r s ( )[ ] 1 ) ( 1 1 ) ( - Γ - Γ - t t r s r s
. Its correlation function, as found in Appendix B, is given by

( ) ( ) [ ] cQ e r s r s t h t h g t t t τ τ τ τ / 2 2 . ) ( 1 ) ( 1 ) ( ) ( - Γ - Γ - Γ = > - < = (10)
The latter result expresses the fact that the interaction is not modulated for s = 0 or = 0 or

Γ t Γ t =1.
The correlation time τ cQ for the quadrupolar interaction is the sum of two independent

contributions c r cQ τ τ τ / 1 / 1 / 1 + = ( 11
)
where is the recombination time of the electron trapped at the donor and is the lifetime of the ionized donor due to capture of a free electron. Following a semi-classical treatment, the quadrupolar-induced evolution of the nuclear spin density matrix for the nuclear spin system, in the interaction representation and within the secular approximation, is given by

* σ 32 [ ] [ ] [ ] [ ( k Q k Qk Qk Q Q J A A r F H i dt d ω σ σ σ σ ∑ - - - = + 0 * 2 2 0 * 0 * , , ) ( , h h ] ) (12)
where is the total static Hamiltonian and J ω , taken for k ω defined by Eq. ( 6), is given by 40

( ) ( ) ∫ ∞ + ∞ - - = τ τ ω τ ω d g e J k i k Q (13)
Here, we assume that the existence of interactions between nuclei results in the establishment of a temperature among the nuclear spin system. With the latter hypothesis, justified in Sec. IIID, the nuclear spin density matrix is, in the high temperature limit, of the form 29

(

) [ ] ( )

1 / 1 Tr H H H Z SS Q IS + + + - ≈ β σ ( 14
)
where β = 1/k B T n , k B is the Boltzmann constant and T n is the temperature of the nuclear spin system. It is then found that the nuclear mean spin lies along the direction of the magnetic field independently on the relative magnitudes of Zeeman and quadrupolar interactions.

Since the latter operator commutes with the static Hamiltonian, the density matrix in the interaction representation is and also the first term of Eq. ( 12) vanishes. An equation for evolution of the inverse nuclear spin temperature β is obtained, after multiplication of Eq. ( 12) by I σ σ σ = * z , taking the trace, and using Eq. (14). Assuming that ( ) ( )

L k cQ k cQ t t t Q Q r s r s r F t β β τ ω τ θ β - ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ + Γ - Γ - Γ - = ∂ ∂ ∑ 2 2 k 2 2 0 2 1 ) ( K 2 ) ( 1 1 ) ( ) ( 1 h (15)
Here β L = 1/k B T L, T L being the temperature of the lattice. The numerical, angle-dependent,

quantity K k (θ), defined by [ ] [ ] { } [ ] 2 z z Qk Qk z k I Tr / I , A , A I Tr ) ( K + = θ (16)
is calculated in Appendix C. Its value is as expected zero for I = ½ and is given by

[ ] ( ) ( ) r E r E I I K off off 2 2 1 , 3 ) 1 ( 4 5 2 ) ( θ θ ⊥ - + = (17) [ ] ( ) ( ) ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ + - + = + r E r E I I K K off off 2 2 // 2 1 , 3 1 . 3 ) 1 ( 4 5 2 ) ( ) ( θ θ θ ( 18
)
where we recall that the parallel and perpendicular components of the electric field, defined with respect to the normal z to the surface, are equal to E off cosθ and E off sinθ, respectively.

The quadrupolar relaxation rate is finally given by

( ) ( ) ( ) ( ) [ ] [ ] ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ + + + - - = 2 2 2 2 2 2 1 1 2 1 1 2 1 2 1 1 cQ cQ cQ cQ on off Q t t Q ) ( K ) ( K . r E r E b , r T τ ω τ θ τ ω τ θ γ Γ Γ θ (19)
Its value is proportional to the square of the amplitude of the modulated electric field and further depends on the angle θ which defines the direction of the electric field.

C Steady-state nuclear mean spin

The rate equation for the evolution of the nuclear mean spin along z, neglecting the thermodynamic nuclear and electronic polarizations in the applied magnetic field as well as nuclear spin lattice relaxation processes other than the hyperfine and quadrupolar ones, is given by to the donor by spin-lattice relaxation, is too short to allow this polarization to be transferred to the bulk nuclei by spin diffusion. As shown in Sec. IIIC, in the latter case, spin diffusion only marginally modifies the results of the present section so that this term will not be considered here. The steady-state value of the nuclear mean spin under the sole effect of spinlattice relaxation is given by + . The quantity

[ ] > < ∆ + > < - ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ + - > < - = > < ) , ( ) , ( ) , ( 1 ) 1 ( 3 4 ) , ( ) ( 1 )
[ ] [ ] z z z S I I r f r f S I I r p r I ) 1 ( ) , ( 1 ) , ( 3 4 ) 1 ( 3 4 ) , ( ) , ( + +
= + = θ θ θ θ (21)
) , r ( f θ , given by ) ( ) , ( ) , ( 1 1 r T r T r f H Q θ θ = (22)
is equal to ) , r ( p θ in the extreme case where the quadrupolar relaxation is much more efficient than the hyperfine one. The relaxation time T 1H (r) is given by 28

( ) [ ] 2 2 2 1 1 2 1 cH H cH * t H ) r ( b r T e τ ω τ γ Γ + = ( 23
)
Here, is the instant electronic hyperfine field acting on the nuclei. The time ) ( * r b e cH τ is the correlation time of the hyperfine interaction. The energy , corresponding to the flip-flop of an electronic and a nuclear spin, is given by

H ω h ) B B ( n e H ± ≈ γ ω h h (24)
14/04/2008 where B n is the nuclear hyperfine field acting on the electrons, which is added or subtracted to B depending on the sign of the electronic spin. The latter energy, which depends on the electronic gyromagnetic ratio γ e , is larger than Assuming that , , and are small with respect to unity, which sets an upper limit to the magnetic field value, the quantity f(r, θ) is finally given by

2 2 cH H τ ω 2 2 1 cQ τ ω 2 2 2 cQ τ ω ( ) [ ] [ ] 1 2 1 1 - ∑ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ - - ≈ k k on off Q * t cQ cH ) ( K . ) r ( E ) r ( E b ) r ( b . ) , r ( f e θ Γ τ τ θ (25)
Note that, since the spatial dependence of the electric fields and does not appear explicitly, Eq. ( 19) and Eq. ( 25) are valid for any localized electronic state. For nuclei near a donor one has Using Eq. ( 2), Eq. ( 17), and Eq. ( 18), it is possible to separate ( )

) (r E off ) (r E on ( ) 2 * 0 * 0 / . ) ( ) ( r a a E b r E b off off Q Q = (26) (27)
θ , r f
into the product of a radial dependence ) (r ϕ , of an angular one, and of a numerical coefficient f 0 which is a measure of the relative strengths of hyperfine to quadrupolar relaxations:

θ ϕ θ 2 0 3 1 cos ) r ( f ) , r ( f + ≈ (28) ( ) 4 * 0 2 1 / 4 ) ( ) ( * 0 ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ = - - a r r s e r a r ϕ (29) [ ] ( ) 2 * 0 * 0 * 1 1 0 ) ( ) ( 1 3 ) 1 ( 4 2 5 ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ Γ - - + = - - a E b a b I I f off cQ cH t Q e τ τ ( 30
)
The implications of the latter equations are discussed in the following section. A key parameter for the value of f 0 is the rate of occupation of the donors by photoelectrons, which depends on the light excitation power. Indeed:

t Γ a)
The correlation times cH τ , and cQ τ depend on the free electron density n f. The time cQ τ can be written using Eq. ( 11)

f c r cQ vn σ τ τ + = 1 1 ( 31
)
where v is the velocity of free electrons and σ c is the cross section for their capture at donors.

The correlation time cH

τ of the hyperfine interaction is given by

f e ex r cH vn T σ τ τ τ ≈ + + = 1 1 2 1 1 1 (32)
as obtained in Appendix B, assuming that the electronic polarization is weak with respect to unity. Here 1 T is the electronic spin-lattice relaxation time and ex τ is the characteristic time for spin exchange between trapped and free electrons. In GaAs, for above bandgap light excitation, it has been found that the latter process is dominant by several orders of magnitude, so that cH τ has a simple approximate expression, also given in Eq. (32), where σ e is the spin exchange cross section. 21 b) The rate Γ t of donor occupation is obtained by writing the rate equation for the

Using Eq. ( 28), Eq. ( 29) and Eq. ( 30), f 0 is given by

14/04/2008 [ ] ( ) ( ) t t t t off e c f a E b a b I I f Q e Γ - Γ = Γ - Γ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ - + = - 1 1 1 ) ( ) ( 3 ) 1 ( 4 2 5 00 2 * 0 * 0 * 1 0 σ σ (34)
Eq. ( 34) has a simple form in which the quantity f 00 , which is a measure of the maximum magnitude of the quadrupolar-induced loss of nuclear magnetization, is independent of experimental conditions such as excitation power. The latter dependence is concentrated in the donor occupation rate . According to Eq. ( 34), the quadrupolar-induced loss of nuclear polarization occurs when the donors are partially occupied, which can be easily characterized from the power dependence of the donor luminescence. For a density of conduction electrons much smaller than (

, one has

Γ t Γ ) 1 - v r c τ σ
t << 1 and the quadrupolar effects are small since the correlation time τ cH is large. Conversely, if Γ t = 1, the quadrupolar interaction is not modulated and cannot relax the nuclear spins.

B Order of magnitude estimates

For As 75 in GaAs, the efficiency of the quadrupolar relaxation process comes from the fact that the spin exchange cross section [σ e ∼9x10 -16 m 2 ] 21 is three orders of magnitude larger than the one for electron capture at donors [σ c =5.1 x 10 -19 m 2 ]. 41 Using Table I and Ref. (29), we obtain and we find f -3 , and f 0 ∼ 10 -2 for Γ t = 1/2. As found from Eq. ( 28), (

) ( 2 . 5 5 . 1 ) ( * 0 * 0 * a E b mT a b off Q e ≈ = 00 ∼ 2 x10
)

-1 * 0 10 2 / , ≈ π a f and ( ) -2 * 0 x10 2 0 , ≈ a f
, so that the nuclei at the Bohr radius are depolarized by the quadrupolar relaxation.

Nuclei such as In 115 in InP and Sb 121 in GaSb are believed to exhibit stronger quadrupolar effects because of their larger spin values (9/2 for In 115 and 5/2 for Sb 121 ).

However, as seen in Table I, the quantity b Q is smaller than for As 75 in GaAs. Using Table I and Eq. ( 34) and assuming that both σ c and σ e scale like the Bohr radius, so that their ratio is independent on material, we estimate that f 00 is equal to 9.4x10 -3 and 1.2x10 -2 for In 115 in InP 14/04/2008 and Sb 121 in GaSb respectively. This implies that the latter materials should also exhibit nuclear polarization losses of quadrupolar origin, although slightly smaller than for GaAs.

C Radial and angular dependences of the nuclear polarization: quadrupolar diffusion radius Shown in Fig. 1 are the radial dependences of p(r, 0) and p(r, 2 / π ) using f 0 ∼10 -2 .

Close to the donor position, one has p(r, θ) =1, as the quadrupolar relaxation is inefficient because

( 3 * 0 / 3 4 ) ( a r r s ≈
) so that the electric field is not modulated. As a function of distance, although the quadrupolar rate first increases and then decreases, f(r, θ) exhibits a monotonic, decreasing behavior. The nuclei are depolarized above a distance to the donor corresponding to f = 1. As seen in Fig. 1, this distance is smaller in the direction z of the magnetic field ( ) than in the perpendicular directions ( ). . As found using Eq. ( 21) and Eq. ( 28), this quantity is given by

[ ] ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ + + = > < ) ( 1 3 ) ( 1 3 ) ( ) ( 0 0 0 r f Arctg r f r f r p ϕ ϕ ϕ (35)
for which the radial dependence, also shown in Fig. 1, is intermediate between those of p(r, 0) and p(r, π/2). The second approximation consists in replacing <p(r)> by a step function at r

= Q ρ such that <p( Q ρ )> = 2 1 ( 36
)
The nuclear hyperfine field, defined by 29 is approximated by 14/04/2008 ( )

, dr r p e r a B B a r n n ∫ ∞ - - > < = 0 / 2 2 3 * 0 0 ) ( 4 * 0
Q n n s B B ρ 0 ≈ ( 37
)
where B n0 is the nuclear field value for a homogeneous nuclear polarization and s(r) is defined in Eq. ( 3). The latter approximation implies that the quadrupolar relaxation is inefficient for distances smaller than Q ρ and dominant for larger distances (f=0). Such approximation is usual in analyses of nuclear polarizations near shallow donors, [16][17][18] and results in defining a sphere around the donor inside which the nuclear polarization is not affected by the quadrupolar relaxation. The radius of this sphere, which will be called the quadrupolar radius, replaces the usual diffusion radius for the estimate of the nuclear hyperfine field. Shown in We now discuss the effect of spin diffusion between neighboring spins, which appears as the third term of Eq. ( 20) and has so far been neglected. A diffusion radius D ρ is defined, corresponding to the distance from the donor at which the efficiencies of direct relaxation and of spin diffusion are equal. This radius is given by Spin diffusion should however be taken account of in two extreme cases concerning the light excitation power. In the case of a very weak efficiency of the overall spin-lattice relaxation, Eq. ( 38) does not have a real solution so that spin diffusion becomes predominant at all distances from the donor. Since the efficiencies of both the hyperfine and quadrupolar relaxation processes are proportional to Γ t , such case is obtained for a weak light excitation density.[for As 75 in GaAs, we estimate that this situation corresponds to a threshold characterized by Γ t <0.15] This situation is outside the scope of the present work and is not considered here. Conversely, for a high light excitation power, if Γ t ∼ 1, the quadrupolar radius ρ Q increases because the quadrupolar spin-lattice relaxation becomes negligible with respect to the hyperfine one. When ρ Q is larger than ρ D , since the polarization of nuclei at a distance larger than ρ D is decreased by spin diffusion towards the bulk nuclei, the nuclear field is obtained by replacing ρ Q by ρ D in Eq. (37). Thus, the maximum nuclear field obtained for a negligible quadrupolar relaxation is given by B n0 s(ρ D ). As a result, the relative decrease of the nuclear field produced by the light-induced quadrupolar relaxation is s(ρ D )/s(ρ Q ) which, in the conditions of Fig. 1, is of the order of 15.

( ) ( ) 2 1 1 1 1 D D H D Q D T T ρ ρ ρ = + (

D Magnetic field effects

The present Subsection is devoted to the justification of three hypothesis made in Sec.

II B, for which the validity depends on the magnetic field value.

1) The zero magnetic field expression of the quadrupolar-induced decrease of nuclear magnetization has been used in Eq. (25). With the values of the cross sections σ c and σ e given in Sec. III. B and taking n f = 10 21 m -3 , one finds that for B= 20T. The same magnetic 43 It is concluded that Eq. ( 25) is valid up to very large magnetic field values.

= cQ τ ω r ∼1 nsec.
2) It has been assumed that the heat capacity of the Zeeman reservoir is larger than those of the quadrupolar and spin-spin ones. Such assumption is obviously not valid at very low magnetic field. The lower magnetic field limit is obtained by expressing the heat capacities of the various reservoirs using the following relation 29

( ) ( ) 2 2 2 2 3 1 1 h γ + - = + > < + > < = > < I I T k B B H H B Z n B Q L Q SS (39)
where the electronic field acting on the nuclear spins has been neglected. Here is the local field. For a magnetic field larger than , only flip-flops between nuclei of the same isotopic species are energetically allowed. In the particular case of As

L B L B
75 in GaAs, this leads to B L ∼0.03 mT. 29 The local field of quadrupolar origin Q B , equal to zero for I = ½, is given by

( )( )( ) ( )( ) ( ) [ ] 3 1 4 1 5 4 1 2 1 3 2 2 2 2 2 - + Γ - = + + > < = I I s E b I I I H Tr B t off Q Q Q h γ (40)
We conclude that the high field limit discussed in Sec. IIB is valid provided

(41) 2 2 2 Q L B B B + >>
Thus the effective local field is larger than the spin-spin local field. For a magnetic field along the z direction, assuming for simplicity Γ t = 0 and taking , we calculate

* 0 5 . 0 a r ≈ B Q ∼ 1.6 mT
which is more than one order of magnitude larger than B L . 29 3) The hypothesis made in Sec. IIB of a nuclear spin temperature gives an independent low magnetic field limit for the validity of the present model. In the absence of quadrupolar couplings, there is no doubt that there exists a spin temperature since the time T 2 of establishment of the nuclear temperature is of the order of 1/γB L ∼300 µsec. Inclusion of 14/04/2008 quadrupolar interactions does not change the latter picture provided the following condition is

fulfilled [ ] [ ] L j Q j Q i Q i Q B E E E E m m m m γ δ δ δ δ h < - - - = Ε ∆ + - ' 1 ' 1 (42)
where ∆E is the energy balance of a difference between transition energies, expressed as a function of the quadrupolar-induced shift of the level m of spin i. Here, j is the nearest neighbor of nucleus i, of the same isotopic specie. There are two distinct reasons for which the latter condition is likely not to be fulfilled:

i Q m E δ i)
For a homogeneous electric field, although the quadrupolar shifts of nuclei i and j are the same, one has for

j Q i Q ' m m E E δ ≠ δ ' m m ≠ .
This may prevent some flip-flops between neighboring nuclei and therefore induce a decrease of the local field. Following Abragam, 44 it is found that such effect leads to a decrease by only 15%, so that quadrupolar interactions weakly affect the time T 2 of establishment of a nuclear spin temperature.

ii) Near a donor, flip flops between nearest neighbors may be prohibited because of the strong spatial dependence of the electric field so that for a given quantum number m. Since the effect of distinct quadrupolar shifts of states with distinct m values has been examined in i) above, we replace Eq. ( 42) by

j Q i Q m m E E δ ≠ δ L j Q i Q mB E E m m γ δ δ h < -
, which states that, for a fixed m, the difference in quadrupolar shifts of neighboring nuclei is smaller than the Zeeman energy in the local field. The latter condition allows us to estimate the characteristic radius r Q of the zone outside which a spin temperature exists, using the following value of , obtained by second order perturbation

i Q m E δ ( ) ( ) ( ) [ ] ( ) ( ) ( ) [ ] ⎪ ⎭ ⎪ ⎬ ⎫ ⎪ ⎩ ⎪ ⎨ ⎧ - - + - - - + = ⊥ 1 2 1 2 , 1 8 1 4 , ) ( 2 ) ( 2 2 2 // 2 2 2 m I I r E r E m I I r E r E B r E mb r E b E i i i i i Q i Q i Q m θ θ γ δ h (43) 14/04/2008 where [ ] ( ) i off t i r E r s r E Γ - = ) ( 1 ) (
is the time average of the electric field. Due to the presence of the magnetic field B at the denominator of the latter equation, the radius r Q decreases with increasing magnetic field. Derivating Eq. ( 43) with respect to distance, we find for r Q an expression of the type:

(44) 5 / 1 - = B r Q η
The radius r Q is largest when the electric field is parallel to z, when Γ t <<1 so that ( )

i off i r E ) r ( E ≈
and for m=3/2 in the case of I=3/2. In the latter case, taking account of the fact that in GaAs the i-j direction is along the [110] crystal axis, and using the known interatomic spacing, we find . It is concluded that, at a distance r from the donor, the hypothesis of nuclear spin temperature is valid provided

5 / 1 3 T nm ≈ η 5 ' ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ = > r B B Q η ( 45
)
Since most of the effects discussed above occur for distances larger than r = , B' 2 / * 0 a Q is evaluated to ∼0.09 T. Although much larger than B L and B Q defined by Eq. ( 41), the latter value is smaller than magnetic fields used in most experiments.

For smaller magnetic fields, the evolution of the mean nuclear spin value > < I , calculated using Eq. ( 12) and ( )

I Tr I σ >= <
, is found to be nonexponential as a function of time so that the calculation of the steady-state nuclear magnetization becomes intricate. However, qualitatively, the decrease of the steady-state nuclear magnetization is still expressed by an equation of the same type as Eq. ( 25) with distinct numerical factors, so that the magnetization is still strongly reduced. Furthermore, because of the very weak magnetic field dependence of r Q , the conclusions of the present work are still qualitatively correct for : For a magnetic field equal to B'

' Q B B <
Q /3 , one has . As seen using Eq. ( 37), such increase is of moderate impact on the nuclear field value since the relative increase of Here we summarize the results of the preceding subsections and obtain the value of the nuclear field B n , using Eq. ( 37) and further considering the various processes on which depends B n . The shallow acceptor concentration N A is assumed to be larger than the donor one so that the donor levels are unpopulated in the absence of light excitation. As shown in Appendix D, the donor occupancy factor t Γ , on which depends the quadrupolar radius ρ Q , is related to the excitation power density P by

( ) ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ + Γ - = D A t N N P P ξ ξ 2 0 1 ( 46
)
where

D A N N k h L P ν = 0 (47) t t c v k Γ - Γ = 1 σ ξ (48)
Here ν h is the photon energy (assumed to be above bandgap) and L is the electron diffusion length. The quantity is the coefficient for bimolecular electron-hole recombination. k Shown in Fig. 3 is the specific case of As 75 in GaAs, using for illustration purposes 47 and , which corresponds to T∼40K.

N D =10 22 m -3 , N A =5x10 22 m -3 , L∼ 5µm, 46 , /s m 10 1x 3 -14 ≈ k
1 10 / - ≈ v k c σ 41
Shown in Curve a is the dependence of t Γ as a function of light excitation power, obtained using Eq. ( 46). The light excitation power range P < 0.2 P 0 corresponds to <0.15 and is not considered here because, as discussed in Sec. IIID, spin diffusion becomes dominant. Shown in Curve b of Fig. 3 is the power dependence of the t Γ 14/04/2008 quantity ) ( Q s ρ , obtained using Eq. (3), Eq. ( 34), Eq. ( 35), Eq. ( 36). Immediately apparent is the fact that, in a relatively broad power range near P 0 , corresponding to =0.5, the nuclear field is decreased with respect to the maximum value of 0.5 imposed by spin diffusion, by a factor of about 20. Note that, at high power, the light-induced quadrupolar relaxation still decreases the nuclear field because, as seen from Eq. ( 46) and Eq. ( 47), , equal to , is smaller than unity. In the present case, one finds

P t Γ t Γ ( 1 / 1 - + v k c σ ) 0 =2.4x10 6 W/m 2 , a realistic value.
The effective nuclear field, measured experimentally from its effect on the electronic polarization, is further modified by the efficient spin exchange processes between free and trapped electrons considered in Sec. IIIA and is equal to the average of the nuclear fields experienced by the two electronic species. Such effect does not modify the quadrupolar radius, but only the multiplicative factor B n0 defined by Eq. (37). Assuming, as performed throughout the present work, that the bulk nuclei are weakly-polarized, the nuclear field experienced by free electrons is very small so that B no is proportional to ( )

D t f D t N n N Γ + Γ /
where the concentration of free electrons is given by Eq. (D6). One has finally

> < = 0 0 S b B n n n α (49)
where does not depend on light power and doping. Since the common mean spin of free and trapped electrons can be measured from the luminescence polarization, the quantity

0 n b > < 0 S n α , given by ) ( Q D t f D t n s N n N ρ α Γ + Γ = (50)
is the reduced nuclear field for which we now consider the power dependence.

Shown in Curve c of Fig. 3 is the light excitation dependence of n f /N A , obtained using Eq. (D6) and assuming N D =N A /5. Shown in Curve d is the light excitation dependence of n α .

14/04/2008

For P< P 0 , the dominant mechanism for nuclear field reduction is the light-induced quadrupolar relaxation. The nuclear field reduction due to spin exchange becomes significant for P> P 0 since n f increases while Γ t is nearly constant and induces an overall decrease of the nuclear field with light excitation power. [START_REF]D. Paget, unpublished] that nuclear effects decrease when the excitation power density[END_REF] The effect of a change of doping of the quadrupolar-induced reduction of nuclear field is limited to the sole variation of the quantity P 0 . As a result, an increase of acceptor and donor doping levels simply shifts Curve b of Fig. 3 along the X axis by a similar factor without any change of shape. As seen from Eq. (50), this is still true if one includes the effect of spin exchange, provided the ratio N A /N D remains constant. Note finally that resonant excitation of donor states might enable to increase the nuclear field value: As seen from Eq. ( 23), Eq. ( 32), and Eq. ( 50), the subsequent decrease of the free electron concentration should induce an increase of the efficiency of the hyperfine relaxation and a decrease of the effect of spin exchange.

We now discuss the possibility of experimental demonstration of the light-induced quadrupolar relaxation. In order to separate the contribution to the nuclear field value of the light-induced quadrupolar relaxation from that of spin exchange with free electrons, it is crucial to analyze the dependence of the nuclear field as a function of light excitation and donor concentration. However, among the experimental works which have estimated the leakage factor f, 22,23,29 none of them has performed the latter analysis, so that experimental proof of the present mechanism is lacking. Such analysis is beyond the scope of the present paper and will be published elsewhere for the case of quantum dots. 37

IV Conclusion

We now summarize the main results of the present work: 14/04/2008 a) The effect of the light-induced quadrupolar relaxation is evaluated assuming that there exists a temperature among the nuclear spin system. The latter hypothesis implies that the external magnetic field is sufficiently large to decrease the difference between the quadrupolar shifts of neighboring nuclei so that flip flops are allowed. In the latter case, the time evolution of the nuclear spin temperature is found to be exponential, so that a relaxation time can be defined. The latter time , within numerical factors, depends on the product of the square of the modulation amplitude and of the correlation time of the modulation.

Comparison of

with the relaxation time due to the hyperfine contact interaction gives the expression for the nuclear polarization as a function of the distance to the donor under the combined effects of quadrupolar and hyperfine relaxations.

I I r V I I eQ r H Y X Y Y X X z z Z Y z z Z X z Z Z Q (A1)
where the quantization axis Z' is the magnetic field direction, the spin operators are equal to

' ± I ' Y ' X iI I ± and ( ) ( ) j i ij X X r E r V ' ' 2 ∂ ∂ ∂ = (A2)
and X' i stands for X', Y', or Z'. These directions are distinct from the xyz directions of the cubic crystal lattice, z being also the normal to the sample surface. The components of the electric field gradient tensor in the X'Y'Z' frame are obtained by using elementary rules for tensor transformation and are given by 35

⎟ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎜ ⎝ ⎛ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ ⎛ - - - - - = ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ ⎛ z y x Y X Z X Z Y Z Z Y Y X X E E E R V V V V V V 14/04/2008
ionicity of the solid, on 2 -n, where n is the infra-red optical index, and on the antishielding factor) and of the covalent contribution (which further depends on the bandgap value). [START_REF] Gill | [END_REF] The expression of the quadrupolar Hamiltonian is then obtained from Eq. (A1) and (A3). For an arbitrary orientation of the magnetic field, this expression is intricate and depends both on θ' and ϕ'. Eq. ( 4) is readily obtained.

In order to estimate b Q, it is necessary to determine R 14 . One of the first determinations was performed for GaAs, where the effect of application of an electric field along the [111] direction on the quadrupolar splitting of the NMR line was studied. [START_REF] Gill | [END_REF] Here, we take the more recent measurements of Ref. 36, which give slightly larger values, arguing that the smaller values obtained in Ref. 33 were due to sample inhomogeneities. For GaAs, InAs and GaSb, independent estimates of R 14 were obtained using the broadening of the nuclear acoustic resonance. 34 For GaAs, they differ from the latter value by about a factor of 3-4. As a result, for a nucleus α of InAs or GaSb, we have chosen to determine according to the following scaling involving Ref. (36) The modulation of the quadrupolar interaction is described by the function h(t), given by Eq. ( 9). This function is of zero average and takes two discrete values (where

α h 2 , 1 = α )
given respectively by

( ) 1 1 1 - Γ - Γ = t t s s h or ( )( ) 1 2 1 1 - Γ - Γ - - = t t s s h
, with respective probabilities and . The correlation function is written under the form The result of Eq. ( 10) is obtained after replacing by their values in Eq. (B1).

t w Γ = 1 t w Γ - =1 2 ) ( - h(t)h(t ) ( τ τ τ αβ β β α α α P h w h g Q ∑ ∑ = > ) < = (
f ( ) 1 1 2 , 1 - = Π τ , ( ) 1 2 1 , 2 - = Π τ , , , where
( ) 1 1 2 , 2 - - = Π τ ( ) 1 2 1 , 1 - - = Π
) ( ) [ ] ( ) ( ) ( ) [ ] ( ) [ ] ( ) ( ) [ ] 1 2 1 1 22 1 2 1 1 12 1 2 1 1 21 1 2 1 1 11 exp 1 exp exp 1 1 exp 1 - - - - - - - - + - Γ - + Γ = + - Γ - Γ = + - Γ - - Γ - = + - Γ + Γ - = τ τ τ τ τ τ

αβ

P

The same procedure can be applied to calculate the correlation function for the hyperfine interaction. Here three states, labeled +1, -1, or 0 are considered, depending on the absence or presence of an electron of a spin equal to +1/2 or -1/2. In addition with the recombination time τ r , the correlation time also depends on the spin-lattice relaxation time T 1 and of the characteristic time τ ex due to possible spin-exchange processes with delocalized electrons. The final expression for the correlation valid in the limit of small electronic polarizations (i. e. ), is

1 1 1 1 - - - >> + r ex T τ τ ( ) cH e g t H τ τ τ / . - Γ = (B4)
where τ cH is given by Eq. ( 32). Eq. (B4) expresses the fact that, unlike for the quadrupolar coupling, the hyperfine relaxation is inefficient in the only case where the probability of occupation of the localized state is zero. where k describes the bimolecular electron-hole recombination. further that the total concentrations of photocreated holes and electrons are equal, one obtains successively

[ ] [ ] { } [ ][ ] { z k Q k Q z z k Q k Q z I A A I Tr I A A I Tr , . , , , , , , , + +
= } (C1) One finds [ ] [ ] { } ( ) [ ] [ ] 2 2 2 2 , 2 , 2 sin , , z z z Q Q z I Tr I I I I I Tr I A A I Tr + - - + + + = θ (C2) [ ] [ ] { } ( ) [ ] [ ]
A D t f N p N n - = Γ + (D4) ( ) [ ] [ ] { } A f D t t c f N n k N k v n g + + Γ + Γ - = 1 σ (D5) ( ) ( ) t t c t A D t f k v k N N n Γ - Γ - Γ + Γ = 1 σ (D6)
The excitation power density P corresponding to a given value of g is given by

υ gLh P = (D7)
where L is the electron diffusion length and υ h is the photon energy. The latter equation assumes that the diffusion length is larger than optical absorption length, and that the surface recombination velocity is negligible. Eq. ( 46) is readily obtained using Eq. (D5), Eq.

  conjugate operators are obtained by replacing i by -i and therefore by . Here θ is the angle between z and the direction Z of the electric field ϕ is the angle between the x direction and the zZ plane. The operators induce transitions between

B

  Calculation of the nuclear relaxation time 14/04/2008

  time due to the hyperfine coupling and <S > is the mean electronic spin along the direction z of light excitation. The third term of the latter equation describes spin diffusion, due to flip-flops between neighboring spins. z 10, 16-18, 28 Here ∆ is the Laplacian operator and D is the diffusion constant. Throughout the present work, it will be considered that the duration of light excitation, although sufficient to polarize the nuclei close

  <1 expresses the reduction of nuclear mean spin with respect to its maximum value [ ]

A

  Effect of the donor rate of occupation.

  population of electrons trapped at donors, of concentration N D . The latter equation, given in Appendix D considering above bandgap light excitation, yields . 1

  For calculation of the nuclear field experienced by trapped electrons, two approximations will be made. First, we shall use for simplicity the angular average of the nuclear polarization, defined as

Fig. 2

 2 Fig.2are the variations of ρ Q and of s(ρ Q ), as a function of f 0 . For f 0 = 10 -2 , one finds ρ Q ∼

4 has 5 .

 45 38)At a distance smaller than ρ D , all the considerations of the present work are valid , since spin diffusion does not affect the nuclear polarization. On the other hand, for nuclei situated beyond ρ D , the polarization is strongly decreased because of the efficient diffusion towards the bulk nuclei. In the absence of quadrupolar relaxation, the maximum value of the reduced nuclear field is s(ρ D ). Since a value ρ D ∼ 1.As seen from Eq. (38), the quadrupolar relaxation results in a decrease of ρ D . For f 0 =10 -2 and using the value of D of Ref.(28), we calculate a modified diffusion radius ρ D value close to . Since the latter value is still larger than ρ * 0 a Q , spin diffusion is generally negligible with 14/04/2008 respect to spin-lattice relaxation for the nuclei which contribute to the nuclear field experienced by trapped electrons.42

E

  , is much smaller than the decreases found in the present work, of more than one order of magnitude. Nuclear field dependence on light excitation power and doping

  shallow donors in semiconductors, the angular-averaged effect of the quadrupolar relaxation is to replace the diffusion radius ρ D up to which the nuclei are spinpolarized by a novel, smaller, radius called the quadrupolar radius ρ Q .c) The quadrupolar-induced decrease of the nuclear field occurs in conditions of light excitation corresponding to partial donor occupancy by photoelectrons. This should induce a decrease of the nuclear field by more than one order of magnitude in GaAs and by slightly smaller factors for InP and GaSb. In addition, the effect of averaging of the nuclear field between free and trapped electrons, due to spin exchange, produces a further decrease of the nuclear field for larger light excitation powers.

  If the magnetic field B direction coincides with a [100] crystal axis z,(θ'=ϕ'=0) the only nonzero components of V ij in the xyz frame are

For

  In 115 in InP no estimate of R 14 has to our knowledge been published. However, R 14 of In 115 in InP should not differ from that of In 115 in InAs by more than a factor of 50% since the ionicities of InAs and InP are identical and since the effect of bandgap should be similar to 14/04/2008 the ratio of the R 14 of As 75 between GaAs and InAs. The final results are shown in Table I. Appendix B: Correlation functions of the quadrupolar and hyperfine interactions

τ 2 ]

 2 is the lifetime of stateα .[with the definitions of Sec. IIA, one has r Using the latter values, resolution of Eq. (B2) yields

(

  Expression of K k (θ) defined by Eq. (16) Applying the relations Tr(ABC)= Tr(BCA) and A, B, C and D are spin operators, one obtains

  calculation proceeds using the following relations, where m is the quantum number of I z

Fig. 1 :Fig. 2 :Fig. 3 :Fig. 4 :

 1234 Figure captions

  The nuclear spin Hamiltonian, given by H = Z + H IS + H SS + H Q , is the sum of the Zeeman term Z, of the hyperfine Hamiltonian H IS , of the nuclear spin-spin interaction H SS and of the quadrupolar interaction H Q . The expressions for the first three terms can be found in Ref.(29). For a cubic semiconductor, the expression for the quadrupolar Hamiltonian is given in Appendix A for arbitrary magnetic field B and sample surface orientations. If the magnetic field is perpendicular to a (001) sample surface, which is the case of a wide majority of experimental situations, the quadrupolar Hamiltonian is simpler:

	14/04/2008
	A Quadrupolar Hamiltonian

  and Ref.(34)

									α 14 R
	(	.	34	)	.	( (	. .	34 36	) )

  Assuming that the fluctuation process is Markovian and stationary, the quantity

										B1)
	where	αβ P	(τ	)		is the conditional probability that
	αβ P	(τ	)	is given by 49
				dt dP αβ	=	γ ∑ Π	( ) (t P αβ β γ ,	)	(B2)
	where ( β γ , Π	)	is a numerical factor, equal for	γ ≠	β	to the probability per unit time that h(t)
	goes from the value	γ h	to the value	h	β	. The quantity -	Π	( ) β β ,	is the probability that h(t)
	goes from			β	to the other value. One has

β h h = at time τ , under the condition that at time t = 0. α h h =

Acknowledgements

We are grateful to B. Urbaszek, X. Marie and D. Petit for useful discussions and to A.

C. H. Rowe for a critical reading of the manuscript.

Appendix A : Form and magnitude of the quadrupolar Hamiltonian

We consider here the general case, described by Fig. 4, where both the electric field direction Z and the magnetic field one Z' do not coincide with a crystal axis z, taken perpendicular to the crystal surface. The quadrupolar Hamiltonian of a given nucleus at position Q H r is related to the components of the electric field gradient by 38,39 ( )

and gives the results shown in Eq. ( 17) and Eq. ( 18).

Appendix D : Calculation of the rate Γ t of donor occupation Complete calculation of Γ t requires considering the kinetics of generation and recombination for the conduction band, the valence band and the donor and acceptor levels.

Although tractable, such calculation leads to intricate results. We assume here, for simplicity and for illustration purposes, that the kinetics of generation and recombination of acceptor levels and of valence holes are similar. Such assumption is reasonable because donor-acceptor recombination, which is specific to holes trapped at shallow acceptors, is known to be less efficient than band to band or exciton recombination. 50 As a result, we consider only one hole specie, for which the total concentration p is the sum of those of valence holes and of neutral acceptors. In steady-state, the rate equations for the concentrations n f of free electrons and Γ t N D of electrons trapped at donors are, respectively ( )

Eq. ( 33) is readily obtained using Eq. (D2). The recombination time τ r of free and trapped electrons is given by The quantity b Q, which has the dimension of the ratio of a magnetic field to an electric field, is given by Eq. ( 7) and characterizes the strength of the quadrupolar relaxation. This quantity estimated in Appendix A, is given below for several nucleus/semiconductor matrix combinations (the isotopic specie under consideration is indicated in bold).

Nucleus R 14 (10


Marie Chavent

Hervé Guégan

Vanessa Kuentz

Brigitte Patouille

Jérôme Saracco

Air pollution sources apportionment in a french urban site

Keywords: air pollution data, Principal Component Analysis (PCA), Positive Matrix Factorization (PMF), rotation

Introduction

Air pollution is a complex mixture of extremely small particles and liquid droplets suspended in the air we breathe. Various sources as factory and utility smokestacks, vehicle exhaust, wood burning, mining, construction activity and agriculture, are known to generate particulate pollution, also called particulate matter (PM). High concentrations of particles have been found to present a serious danger to human health [START_REF] Pope | Lung cancer, cardiopulmonary mortality, and long-term exposure to fine particulate air pollution[END_REF][START_REF] Samet | Fine particulate air pollution and mortality in 20 US cities, 1987-1994[END_REF]. In this study, particles of special preoccupation to the protection of lung health are those known as fine particles, less than 2.5 microns in diameter and called PM2.5 in the rest of this paper. Thus, a wide preoccupation of environmental protection agencies concerns the development of PM2.5 control strategies. One of the main goals of these strategies is to improve ambient air quality. In consequence, this involves the reduction of emissions from primary sources. Therefore it is important to be able to identify these air pollution sources and apportion the contributions of these sources.

A reliable way to provide information regarding source characteristics is often obtained from a receptor modelling approach, using measurements of chemical composition data for particles at a sample site; see [START_REF] Hopke | Receptor Modeling for Air Quality Management[END_REF] for some details and useful references. Most of the multivariate receptor models are based on the analysis of the correlations between measured concentrations of chemical species, assuming that highly correlated compounds come from the same source. Principal Component Analysis (PCA) is commonly used as multivariate receptor model, and this multivariate method has been successfully applied to identify sources in several studies. However, PCA fails to quantify sources contributions. Then specific methods are needed to address this problem. One of them is Positive Matrix Factorization (PMF); see for instance [START_REF] Paatero | Positive Matrix Factorization: a non-negative factor model with optimal utilization of error estimates of data values[END_REF].

The case study presented here corresponds to the statistical part of the scientific program PRIMEQUAL (Projet de Recherche Interorganisme pour une MEilleure QUalité de l'Air à l'échelle Locale), initiated by the French ministry of Ecology and Sustainable Development and the ADEME (Agence de l'Environnement et de la Maîtrise de l'Energie, that is French Environment and Energy Management Agency), about atmospheric pollution and its impacts. In this statistical work, a methodology for determining particulate emission sources and their concentrations at the urban site of Anglet located in the South-West of France has been proposed and applied. This paper is based on a more complete work presented in [START_REF] Chavent | PCA and PMF based methodology for air pollution sources identification and apportionment[END_REF]. The following three steps process has been implemented:

(i) The air pollution data (that is PM2.5) were collected with sequential fine particle samplers on the receptor site and the chemical composition of each sample was measured with PIXE (Particle Induced X-ray Emission) method. A data matrix of chemical compounds concentrations in each sample has been obtained after several pretreatments.

(ii) To identify possible air pollution sources, we have implemented PCA approach to this data matrix and we have rotated the standardized principal components in order to obtain more interpretable results.

(iii) For the sources apportionment issue, we have applied PMF to the same data matrix and normalized the results so as to find components with physical interpretations (concentration of each source in each sample). It is interesting to note that steps (ii) and (iii) are numerically and computationally independent. Because the PMF method can be used for both identifying and quantifying the pollution sources, case studies usually don't mix PCA and PMF. But in practice, it can be difficult when using PMF to identify potential sources without some sort of profile to which to compare the final results. In this case study, we mix PCA and PMF in the sense that we check that each source quantified with PMF is clearly correlated with a single source identified with PCA. Sources difficult to identify with PMF are then clearly identified with the help of PCA and the identity of the others are validated with PCA. The rest of this case study is organized as follows. A description of the air pollution data set is given in Section 2. The sources identification step via PCA and Varimax rotation is developed in Section 3. Section 4 is devoted to the sources apportionment step via PMF. Finally, Section 5 gives conclusion and summarizes the paper.

The data set

PM2.5 samplers were collected by AIRAQ1 during December 2005 and July 2006 at the French site of Anglet located in the South-West of France (see Figure 1). In this case study, we only exhibit the results corresponding to the winter data set. This sampling site located at "Station fixe d'Anglet" (see the map given in Figure 1) was chosen because of its proximity with:

a high traffic road in red on the map, -three cities (Bayonne, Anglet, Biarritz) with a total of 170.000 inhabitants, -an industrial area with a steelworks and a refinery at North-East, -Atlantic Ocean at West. This receptor site is thus subject to different pollution origins: traffic road, urban and industrial activities, and natural dust. The knowledge of potential origins has been determinant in the choice of the site. Indeed, it allowed to check if the sources identified with the statistical methodology (without using any information about the sources) were coherent with the expected ones. The n = 61 samples of PM2.5 were collected every twelve hours: one for the day (7AM:19PM) and one for the night (19AM:7PM). The mass, the volume and the concentration C in ng/m3 of each particle sampler were measured with the PIXE method by ARCANE-CENBG2 , as well as the concentrations of p = 15 chemical elements (Al, Si, P, S , Cl, K, Ca, Ti, Mn, Fe, Ni, Cu, Zn, Br, Pb) Elements Ni and Ti that were frequently present at concentrations below the detection limits (BDL) were excluded and only 13 elements were selected. Then the few BDL data remaining in this data set were replaced by values corresponding to one-half the appropriate analytical detection limit. In nature, the elements Al, Si, S and Fe are usually found in the following oxidized forms: Al 2 O 3 , Si0 2 , S0 4 , Fe 2 0 3 . For this chemical reason, Al, Si, S and Fe were replaced by the compounds Al 2 O 3 , Si0 2 , S0 4 , Fe 2 0 3 : we added the mass of the measured element to the mass of oxygen of its oxidized form. Then, the remaining concentration, called C org , which was not measured by the previous compounds and elements was calculated for each particle sampler:

C org = C-(Al 2 O 3 + Si0 2 + P + S0 4 + Cl + K + Ca + Mn + Fe 2 0 3 + Cu + Zn + Br + Pb).
The addition of the column C org in the data matrix is a key point specific to this case study. We have indeed noticed during the data pre-treatment that the remaining concentrations not measured by PIXE could be an important part of the samples concentrations. We can see Figure 2 that the proportion of C org in the samples which have a concentration greater than 5 µg/m 3 , is at least 50% of the total concentration. C org will then be used to distinguish among the sources identified, those mostly participating to the concentration in PM2.5. The other 15 columns will be used to identify the sources. For instance, it is known that the elements Zn and Pb are in particulates emitted by industrial sources. A meteorological data set was also used to interpret and validate some results. Hourly temperatures and wind directions (in degrees) were collected during the sampling period at a meteorological station located 2.5 km far from the sampling site. Temperatures were averaged to match with the 12-hours samples. The 360° of the wind directions were split into 8 categorical wind directions (North, North-East, North-West, South….) and a wind direction data matrix of 61 rows (12-hours samples) and 8 columns (wind directions) were constructed. An element of this matrix is the percentage of hours during which the wind direction has been observed (see Table 3).

Date

N N-E E S-E S S-W W N-W 23-11-05 day 23-11-05 night 24-11-05 day 24-11-05 night

! M 17% 25% 0% 0% ! M 0% 17% 0% 0% ! M 25% 33% 8% 0% ! M 8% 25% 17% 0% ! M 8% 0% 50% 8% ! M 0% 0% 25% 9% ! M 0% 0% 0% 58% ! M 42% 0% 0% 25% ! M Table 3: Wind directions data

PCA and Varimax rotation for sources identification

In order to identify the sources of fine particulate emission, we applied PCA to the concentration matrix X. The idea was to find groups of correlated chemical compounds that are characteristics of air pollution sources. For instance, if the elements Zn and Pb are strongly correlated to the same factor, since these elements are known to have industrial origin, this factor will be associated to industrial pollution source. In PCA, we consider a (n,p) numerical data matrix X where n objects are described on

! p < n variables ! x 1 ,K, x p . We note ! x j a column of X. Let ! ˜ X = ˜ x ij ( ) n, p be the standardized data matrix with ! ˜ x ij =
x ij " x j s j where !

x j and ! s j are the sample mean and the sample standard deviation of !

x j . The basic idea underlying Factor Analysis (using correlation matrix) is that the p observed standardized variables

! ˜ x 1 ,K, ˜
x p , can be expressed, except for an error term, as linear functions of ! q < p unobserved variables or common factors ! f 1 ,K, f q . Given the observed standardized matrix ! ˜ X , Factor Analysis model can be expressed in its simplified form as:

! ˜ X = FA'+ E
, where F is the (n,q) matrix of unobserved values of the factors, whereas the (p,q) matrix A is the unknown loading matrix which provides information that relates the factors ! f 1 ,K, f q to the original variables ! x 1 ,K, x p . Several approaches were developed to estimate the model (principal factor, maximum likelihood, …) but PCA is often used in practice. With PCA estimation method, the q columns of F are the first q standardized principal components (which are mutually orthogonal and of variance equal to 1) and each element ! a jk of A is equal to the correlation between the variable j

x and the factor ! f k . In this case study, we applied PCA to the concentration matrix X where the 61 samples (in row) are described by 14 compounds (in column). We see Figure 3 that each k th column of F and k th row of A' obtained with PCA will be associated to a source. The approach is the following: we search in the k th row of the loading matrix A', the compounds strongly correlated with the k th factor. If these compounds are known to be characteritisc from a source, this source is associated to this factor. Because the n samples are chronologically ordered, the k th column of the factor score matrix F gives an idea of the evolution of the quantity of fine particulate emitted by the source associated with the k th factor. Of course, depending of the correlations in the loading matrix, it is not always possible to associate clearly a source to a factor. In this study, after several trials, we have chosen q=5. Indeed, it was not possible to associate clearly a source to each factor with the decompositions in more than 5 factors. Moreover, with q=5, we have 90,93% of total variance that is explained. where:

-data=hiverorg is the sas dataset constructed from the concentration data matrix (Table 2) -method=prin because the factor analysis method is PCA -outstat=load is the sas dataset with the loadings reported Table 4 f 1 f 2 f 3 f 4 f 5 Al 2 O 3 0,672 -0,663 0,221 -0,193 0,023 SiO 2 0,649 -0,669 0,254 -0,203 -0,058 P 0,682 -0,629 0,223 -0,238 0,006 SO 4 0,589 0,449 -0,413 -0,220 0,175 Cl -0,474 -0,281 0,192 0,681 0,344 K 0,888 -0,154 -0,209 0,004 0,284 Ca 0,638 -0,405 0,099 0,399 -0,195 Mn 0,384 0,776 0,183 0,093 -0,247 Fe 2 O 3 0,793 0,319 0,038 0,269 -0,360 Cu 0,796 0,248 -0,098 0,232 -0,360 Zn 0,352 0,663 0,589 -0,072 0,246 Br 0,746 -0,182 -0,126 0,352 0,363 Pb 0,428 0,659 0,519 -0,088 0,296 C org 0,600 0,297 -0,613 -0,025 0,222

Table 4: Correlations between the chemical compounds and the 5 factors

Since most compounds are positively correlated with the first standardized principal component, it is difficult to detect groups of correlated elements. However the graphical representation of these compounds according to their correlations with f 1 and f 2 (see Figure 4) or according to their correlations with f 2 and f 3 (see Figure 5) shows groups of compounds that seem to be correlated with each other (Zn and Pb for instance or P, Al 2 O 3 and SiO 2 ). Because it is known that Zn and Pb for instance come from fine particulates with industrial origin, we would like to see clear correlations between those two elements and a factor. In order to identify more clearly groups of correlated compounds, a Varimax rotation has been applied to the standardized principal components. The idea of rotation in Factor Analysis is the following. Let T be an orthogonal transformation matrix, ! TT'= T'T = I q . The factor analysis model can then be expressed as: where:

! ˜ X = FT ˜ F { T' A' ˜ A { + E with: - ! ˜ F = FT the
-outstat=load is the sas dataset containing the rotated loading -out=fact is the sas dataset containing the rotated factor scores

! ˜ f 1 ! ˜ f 2 ! ˜ f 3 ! ˜ f 4 ! ˜ f 5
Al 2 O 3 0,981 0,087 -0,042 0,070 -0,038 SiO 2 0,979 0,012 -0,055 0,104 -0,074 P 0,972 0,090 -0,017 0,071 -0,092 SO 4 -0,028 0,765 0,247 0,180 -0,345 Cl -0,153 -0,274 -0,136 -0,181 0,879 K 0,597 0,716 0,111 0,233 0,031 Ca 0,608 0,091 -0,113 0,560 0,272 Mn -0,279 0,119 0,604 0,582 -0,238 Fe 2 O 3 0,198 0,282 0,289 0,848 -0,112 Cu 0,213 0,359 0,161 0,816 -0,149 Zn -0,029 0,053 0,977 0,129 -0,044 Br 0,490 0,615 0,097 0,281 0,392 Pb 0,004 0,163 0,969 0,126 -0,054 C org -0,018 0,893 0,021 0,222 -0,160

Table 5: Correlations between the chemical compounds and the 5 rotated factors

The loadings after rotation are reported Table 5. These new loadings are correlations between the rotated factors and the compounds. We notice now that the five rotated factor

! ˜ f 1 , …,
! ˜ f 5 are clearly correlated with some compounds (highlighted in grey in the table). Because some of these compounds are known to be characteristic from specific pollution sources, we have il is possible here to associate five sources to the 5 rotated factors: ! ˜ f 5 is clearly correlated with Cl, which is known to have sea salt origin. In the same way, the elements Zn and Pb correlated with ! ˜ f 3 have industrial origin and the compounds Al 2 O 3 and Si0 2 correlated with

! ˜ f 3 are coming from soil dust. Factors ! ˜ f 1 , ! ˜ f 3 and ! ˜ f 5 can
then be clearly associated to those three sources (see Table 6).

Rotated factor

Source Characteristic correlated elements or compounds

! ˜ f 1 Soil dust Al 2 O 3 , Si0 2 ! ˜ f 2 Combustion S0 4 ! ˜ f 3 Industry Zn, Pb ! ˜ f 4 Vehicle Fe 2 0 3 , Cu ! ˜ f 5

Sea

Cl

Table 6: Factor-source associations

In the same way, S0 4 correlated with ! ˜ f 2 is usually linked to combustion and Fe 2 0 3 and Cu correlated with ! ˜ f 4 can be linked linked traffic road. In order to confirm these two last associations we have confronted the rotated factors ! ˜ f k with external informations such as meteorological data (temperatures and wind directions) and the periodicity night/day of the sampling. The coefficients ! ˜ f ik of the column represent here a "relative" contribution of the source k to the sample i. Figure 6 gives the evolution of the relative contribution of the source associated with ! ˜ f 4 . The night have been distinguished from the day ones, which enables to notice that the contribution of this source is stronger during the day than at night. It is then a confirmation that this source corresponds to vehicle pollution. We notice an increase in the contribution of this source at the middle of the sampling period, which corresponds to a decrease in the temperature measured on the sampling site, see Figure 8. Finally, correlations between the rotated factor scores and the eight vectors of wind directions (see Table 3) have been calculated. Figure 9 (resp. Figure 10) is a graphical representation of the correlations between

! ˜ f 3 (resp ! ˜ f 5
) and wind directions. The correlation of ! ˜ f 3 with East and North-East wind directions confirms that this source corresponds to industrial pollution. Indeed, the main industries are located at North of the sampling site (see Figure 1). In the same way, the correlation of ! ˜ f 5 with West and North-West wind directions is a confirmation that this source corresponds to Atlantic Ocean located at West of the sampling site, and then to sea pollution.

Sources apportionment via PMF

We have seen how fine particulate pollution sources can be identified by applying PCA to a concentration data matrix. But the identification is not sufficient. The danger for the health is indeed linked to the quantity of PM2.5 we breath. The problem is then not only to identify the PM2.5 sources but also to determine in which proportion these sources participate to the global dust contamination. In order to quantify unknown sources of fine particulate emission, we have approximated a receptor model by first defining a Positive Matrix Factorisation of the concentration matrix X, and then by normalizing the results to find components with physical interpretations.

The basic problem of receptor modelling is to estimate, from the data matrix X and the number q of sources, their compositions and their contributions. To address this problem, we consider the mass balance equation:

! = = q k jk ik ij b g x 1 where -ij
x is the concentration of the jth chemical species in the ith sample, ik g is the concentration in particles from source k in the sample I, is a mass fraction or, in other words, it is the percentage of species j emitted by source k when sample i was collected. Since the mass fraction jk b is independent from i, it means that the sources profiles are assumed to be constant during the sampling period.

-
In matrix form, the mass balance equation can be written:

' GB X =
where G is a (n,q) matrix of sources contributions and B is a (p,q) matrix of sources compositions (see Figure 11).

Concentration matrix X

Source k Source k Compound j  Matrix G of sources contributions Sample i  Compound j  = Source ?
Matrix B of sources profiles In this case study, 5 sources were clearly identified with PCA. We have then chosen to approximate the matrices G and B for q = 5 in the receptor model. Two steps (PMF and normalization) have been necessary to approximate G and B. Once these approximations G ˆ and B ˆcalculated, the user has an approximation of the quantities of fine particulate emitted by 5 sources in each sample. He has also an approximation of the profiles of the sources. But no name is associated to each source. As we will a supplementary step is necessary to identify clearly the sources quantified in G ˆ.

PMF step.

The matrix X is factorized in a product HC' of rank q under constraints of positivity of the coefficients. This condition is required by physical reality of non-negativity of sources compositions and contributions:

0 ! ik g and 0 ! jk b
. The PMF algorithm developed by [START_REF] Paatero | Positive Matrix Factorization: a non-negative factor model with optimal utilization of error estimates of data values[END_REF] in the context of receptor modelling minimizes

!! ! = = = " " " " # $ % % % % & ' ( n i p j ij jk q k ik ij c h x 1 2 1 1 ) subject to 0 ! ik h and 0 ! jk c
. The coefficient ij ! is a measure of uncertainty of the observation ij x . In this case study, we are dealing with variables measured on very different scales which can cause problems when approximating X globally on all the variables (for instance, minimizing an unweighed quadratic error, that is when

, 1 = ij !
will give better approximations for the columns of X corresponding to variables with large dispersion). Hence we have opted for

j ij s = !
, the sample standard deviation of the jth variable. The numerical results have been obtained with the PMF algorithm and program poposed by [START_REF] Jianhang | Technical details and programming guide for a general two wat positive matrix factorization algorithm[END_REF]. This program takes as input:

the concentration matrix

! ! = = = = q k jk k k ik q k jk ik ij c h c h x 1 1 " "
. Then, the matrix X ˆ can be written:

' ˆC H X ( ( = with k ik ik h h ! = ( and k jk jk c c ! = (
. The objective of scaling is then to define the scaling

constants q k k , , 1 , K = !
such that H ( and C ( verify the physical conditions of the matrices G and B of the mass balance equation. In this study, we consider the two following conditions:

(a) The sum of the concentrations of the sources adds up to the total concentration of the samples; that is for each sample,

! = = q k ik i g 1 "
where i ! is the total concentration in the ith sample.

(b) If the sum of the concentrations of the observed species adds up to the total concentration of the samples, then the sum of all species in each source profile is equal to unity, that is:

! b jk = 1 j=1 p " if ! x ij = " i j=1 p # .
Note that the introduction in this case study of the 15 th column C org (the concentration in the samples not measured by PIXE) calculated for each particle sampler yields

! x ij = " i j=1 p # .
From the physical constraint (b), the scaling coefficients k ! can be directly calculated: from

' Ĉ H X =
, we get

! = = p j jk k c 1 "
. Details can be in [START_REF] Chavent | PCA and PMF based methodology for air pollution sources identification and apportionment[END_REF].

In practice, we have imported in Excel the numerical values of the matrices the 5 scaling coefficients (given Table 7) by summing the 5 columns of the matrix gives the percentages reported in Table 8. We deduce possible identifications of the five pollution sources, see Table 9.

! k = 1 ! k = 2 ! k = 3 ! k = 4 ! k = 5
! k = 1 Soil dust ! k = 2 Vehicles ! k = 3 Industry ! k = 4 Combustion ! k = 5
Sea One can remark that the sources identified in Table 9 are the same than those found with PCA in Table 6. To verify the coherence of these sources identifications, we have calculated, in Table 10, the correlations between the factors obtained with PCA and the sources obtained by receptor modelling (the columns of G ˆ). We observe that the factors match well with the receptor model sources.

Sources apportionments.

From matrix G ˆ of the source contributions, we can deduce some interesting comments. First we can focus on the relative contribution of each source in each particle sampler. For instance, Figure 13 represents the relative contributions of the combustion source in the 61 particle samplers. We can notice the increase in the part of this source in the second period of sampling, corresponding to a decrease in the temperature (see Figure 8). We can also focus on the contribution of the sources to the PM2.5 dust contamination during the sampling period. Figure 14 shows the domination of the combustion source during this winter sampling period.

Conclusion

We exhibited in this case study a methodology in order to determine air pollution sources apportionment in a French urban site. The data correspond to measurements of chemical composition data for particles. Our approach is divided into two parts: first we identify the sources profiles via Principal Component Analysis followed by a Varimax rotation. Then, we evaluate the apportionment of the sources via a receptor modelling approach based on Positive Matrix Factorization as estimation method. The corresponding numerical results allow to characterize and apportion five principal sources of fine particulate emission: soil dust, vehicles, industry, combustion and sea. These results are not hardly surprising, they confirm the environmental knowledge of the phenomenon of air pollution by fine particles. What is interesting here is the fact that we do not use any prior information in order to retrieve the usual five air pollution sources. This kind of methodology can then be applied to any dataset and problems of identification and apportionment of sources. For mathematical details on the proposed methodology, the reader can refer to [START_REF] Chavent | PCA and PMF based methodology for air pollution sources identification and apportionment[END_REF].

Figure 1 :

 1 Figure 1 : The French urban site

Figure 2 :

 2 Figure 2: The proportion of C org in the sample according to the total concentration C of the sampleTable 2 shows the data set after having applied the transformations described above. The (n,p) concentration matrix

Figure 3 :

 3 Figure 3: Decomposition of the standardized concentration matrix

Figure

  Figure 4: Factor 1-2 correlation circle

  matrix of the rotated standardized principal components (rotated factors) which are still mutually orthogonal and of variance equal to 1, -! ˜ A = AT the matrix of rotated loadings which are correlations between the variables and the rotated principal components. From a practical point of view, the orthogonal transformation matrix T is defined in order to construct a matrix! ˜ A such that each variable ! x j is clearly correlated to one of the rotated factor ! ˜ f k (that is! ˜ a jk close to 1) and not to the other rotated factors (that is rotation technique is Varimax. It seeks rotated loadings that maximize the variance of the squared loadings in each column of ! ˜ A . The matrix ! ˜ A of the loadings after rotation and the matrix ! ˜ F of the rotated factor scores are obtained with the procedure FACTOR of SAS:PROC FACTOR data=hiverorg method=prin nfactors=5 outstat=load out=fact rotate=varimax; var Al2O3 SiO2 P SO4 Cl K Ca Mn Fe2O3 Cu Zn Br Pb Corg; run;

Figure 6 :

 6 Figure 6: Evolution of the rotated factor 4 associated to cars pollution Figure 7 gives the evolution of the relative contribution of the source associated with

Figure 7 :

 7 Figure 7: Evolution of the rotated factor 2 associated to heating source

Figure 9 :

 9 Figure 9: Correlations between the rotated factor 3 associated to industrial source and the 8 wind directions

  jk b is the mass fraction (percentage) of species j in source k. In the common parlance of receptor modelling, the ! b jk 's are the sources compositions (or sources profiles) and the ik g 's are the sources contributions. The product jk ik b g is then the approximation of the concentration in the sample i in particles from the jth species coming from the source k. Let ijk m be the mass, in the sample i, of species j from source k, and let ! m ik be the mass in the sample i from source k.

Figure 11 :

 11 Figure 11: Decomposition of the concentration matrix

  of the above constrained minimization problem. Obviously, this solution is not unique and other physical constraints have been used to calculate from these two matrices G ˆ and B ˆof the contribution matrix and of the profile matrix. This is the scaling step. Let ' The columns of the approximations H ˆ and C ˆ obtained in the previous step, must be scaled in order to get the approximations G ˆ and B ˆ. The scaling coefficients are defined to fulfil other physical constraints of the sources compositions and contributions.Let us first remark that

  of the PMF program and calculated:

--

  the (14,5) matrix B ˆ of the approximated compositions (profiles) of the 5 sources on the 14 compounds, by dividing the five columns of the (61,5) matrix G ˆ of the approximated concentrations of the 5 sources in the 61 samples, by multiplying the five columns of

Table 7 :

 7 The scaling coefficientsQuality of the model approximation. We can evaluate the quality of the approximation of X by ' B G : the Figure12shows a good fitting of the i ! 's by the i !

Figure

  Figure 12: Adjustment of

Figure 13 :

 13 Figure 13: Relative contribution of the source combustion to the samples

Figure 14 :

 14 Figure 14: Global sources contributions to the PM2.5 dust contamination

Table 1 : Subset of the original data table

 1 3 . Table1gives a subset of the data in their initial form. We notice on this data table that the concentrations of the 15 elements measured with PIXE are very small comparing to the total concentrations C in the samples. Obviously the sum of the 15 concentrations measured with PIXE doesn't add up to the total concentration of the samples. The elements (H, C, N, O) 4 not measured with PIXE represent almost all the remaining concentration.

		Date	C	Al Si … K Ca … Br Pb
		23-11-05 day	7264.2	92		75		…	163	35	…	7	10
		23-11-05 night	9633.0	135		90		…	211	23	…	7	77
		24-11-05 day	10952.4	175	137		…	241	69	…	8	19
		24-11-05 night	5333.3	36		31		…	94	44	…	9	7
		M	M	M		M		M	M	M	M	M	M
		24-12-05 day	20978.3	<2		<1		…	266	<1	…	7	18
	!	24-12-05 night !	18130.8 !	18 !	!	<1	!	…	307 ! !	<1 !	… !	7	19
		25-12-05 day	23297.9	37		22		…	311	12	…	7	14
		25-12-05 night	36105.3	<2		<1		…	277	<1	…	10	19

Table 2 : The final data set

 2

Table 4

 4 gives the loading matrix A' obtained with the following procedure FACTOR of SAS:

	PROC FACTOR data=hiverorg method=prin nfactors=5 outstat=load;
	var Al2O3 SiO2 P SO4 Cl K Ca Mn Fe2O3 Cu Zn Br Pb Corg;
	run;

Table 8 : Relative contributions of the sources to the chemical compounds

 8

			5	3.9	8.2	6
	S0 4	4.5	9.5	10.7	67.9	7.5
	Cl	0.0	0.0	0.0	0.0	100.0
	K	38.8	0.0	4.4	56.7	0.2
	Ca	42.0	39.6	0.0	0.0	18.4
	Mn	0.0	54.9	33.1	8.5	3.5
	Fe 2 0 3 19.0	59.2	14.4	7.4	0.0
	Cu	18.5	56.8	9.1	15.6	0.0
	Zn	9.0	0.5	87.5	0.0	3.1
	Br	19.4	12.1	5.7	33.4	29.4
	Pb	10.7	0.0	81.4	7.9	0.0
	C org	0.0	8.0	0.0	92.0	0.0

  Table8is used to identify the nature of the sources. For instance, Al 2 0 3 and SiO 2 are exclusively emitted by source 1. Because Al 2 0 3 and SiO 2 are known to have natural origin, this source is associated to the soil dust pollution source. We proceed in the same way for the other sources.

Table 9 : Receptor model sources identification

 9

Table 10 : Correlations between the sources of the receptor model and the factors of PCA Sources descriptions.

 10 The matrix B ˆ of the sources profiles is reported in Table11. We notice that, according to these profiles, C org which represents almost the total concentration in PM2.5, is only emitted by the vehicle and combustion sources. It is also possible to describe the composition of the sources using Table11: for instance, sea pollution source is made of around 75% of Chlorine and 23% of S0 4 .

		Source 1 Source 2 Source 3 Source 4 Source 5
	Factor 1	0.98	-0.18	-0.11	-0.02	-0.18
	Factor 2	0.11	0.12	0.06	0.95	-0.30
	Factor 3	-0.05	-0.09	0.98	0.02	-0.15
	Factor 4	0.12	0.96	0.10	0.11	-0.22
	Factor 5	-0.02	-0.13	-0.10	-0.27	0.88
		Soil dust Vehicles Industry Combustion Sea
	Al 2 O 3	41.6	0.0	0.0	0.0	0.0
	Si0 2	18.5	0.0	0.0	0.0	0.0
	P	6.2	0.0	0.7	0.0	0.6
	S0 4	10.1	15.3	59.6	12.2	22.6
	Cl	0.0	0.0	0.0	0.0	74.5
	K	12.9	0.0	3.6	1.5	0.1
	Ca	2.4	1.6	0.0	0.0	1.4
	Mn	0.0	0.2	0.3	0.0	0.0
	Fe 2 0 3	6.7	15.0	12.7	0.2	0.0
	Cu	0.3	0.7	0.4	0.0	0.0
	Zn	0.7	0.0	16.3	0.0	0.3
	Br	0.2	0.1	0.2	0.0	0.5
	Pb	0.3	0.0	6.2	0.0	0.0
	C org	0.0	67.1	0.0	85.9	0.0

Table 11 : The sources profiles

 11

Réseau de surveillance de la qualité de l'air en Aquitaine

Atelier Régional de Caractérisation par Analyse Nucléaire Elémentaire -Centre d'Etudes Nucléaires de Bordeaux Gradignan

Aluminum (Al), Silicon (Si), Phosphorus (P), Sulphur (S), Chlorine (Cl), Potassium (K), Calcium (Ca), Titanium (Ti) Manganese (Mn), Iron (Fe), Nickel (Ni), Copper (Cu), Zinc (Zn), Bromine (Br), Lead (Pb).

Hydrogen (H), Carbon (C), Nitrogen (N), Oxygen (O).