Piotr Graczyk
email: graczyk@univ-angers.fr

Hideyuki Ishi
email: hideyuki@math.nagoya-u.ac.jp

Riesz measures and Wishart laws associated to quadratic maps

Keywords: 2000 Mathematics Subject Classification: 43A35, 62H05, 15B48 convex cones, homogeneous cones, Riesz measures, Wishart laws

We introduce a natural definition of Riesz measures and Wishart laws associated to an Ω-positive (virtual) quadratic map, where Ω ⊂ R n is a regular open convex cone. We give a general formula for moments of the Wishart laws. Moreover, if the quadratic map has an equivariance property under the action of a linear group acting on the cone Ω transitively, then the associated Riesz measure and Wishart law are described explicitly by making use of theory of relatively invariant distributions on homogeneous cones.

§1. Introduction.

Riesz measures and distributions on convex cones form one of fundamental tools of harmonic analysis and of the theory of the wave equation, cf. [START_REF] Faraut | Analysis on symmetric cones[END_REF] in the case of symmetric cones and [START_REF] Gindikin | Analysis in homogeneous domains[END_REF], [START_REF] Ishi | Positive Riesz distributions on homogeneous cones[END_REF] for homogeneous cones. Moreover, exponential families generated by Riesz measures are composed of Wishart laws and are of great significance in random matrix theory and in statistics.

Wishart laws are probability distributions on symmetric or Hermitian matrices with very important applications in multivariate statistics. Their role in statistics is due to two reasons:

-they are probability distributions of the maximum likelihood estimator(MLE) of the covariance matrix in a multivariate normal sample ( [START_REF] Muirhead | Aspects of multivariate statistical theory[END_REF], [START_REF] Andersson | Wishart distributions on homogeneous cones[END_REF]).

-in Bayesian statistics, Wishart laws form a Diaconis-Ylvisaker family ( [START_REF] Diaconis | Conjugate priors for exponential families[END_REF]) of prior distributions for the covariance parameter in a covariance selection model ( [START_REF] Letac | Wishart distributions for decomposable graphs[END_REF]).

On the other hand, recent developments in random matrix theory of chiral Gaussian ensembles containing Wishart laws, are intense and motivated by applications in mathematical physics, cf. [START_REF] Katori | Symmetry of matrix-valued stochastic processes and noncolliding diffusion particle systems[END_REF] and references therein.

These numerous modern applications of Wishart laws make it necessary to develop the theory of Wishart laws and Riesz measures on more general cones than in the classical case of the symmetric cones of real symmetric or complex Hermitian matrices. For example, in an r-dimensional Gaussian model X, if the marginal variables X i and X j are known to be conditionally independent given all the other variables, the statistical analysis of the covariance matrix of X must be done on the cone P of positive definite symmetric matrices Y with Y ij = Y ji = 0 and on its dual cone Q ( [START_REF] Letac | Wishart distributions for decomposable graphs[END_REF]). The cones P and Q are usually no longer symmetric. This led to some important papers in recent statistical and probabilistic literature about Wishart laws on more general cones: homogeneous cones ( [START_REF] Andersson | Wishart distributions on homogeneous cones[END_REF]) or cones related to graphical models ( [START_REF] Letac | Wishart distributions for decomposable graphs[END_REF]). In these papers, Wishart laws are introduced via their density functions (see Section 3.8).

In our paper we construct and study Riesz measures and Wishart laws on regular convex cones via quadratic maps. For a regular open convex cone Ω ⊂ R n and an Ω-positive quadratic map q : R m → R n , the Riesz measure associated to q is defined as the image of the Lebesgue measure dx on R m by q. Wishart laws studied in this paper are obtained from R m -valued normal random vectors X as the law of Y := q(X)/2. This is a classical natural approach to Riesz measures ( [START_REF] Faraut | Analysis on symmetric cones[END_REF]) and Wishart laws ( [START_REF] Muirhead | Aspects of multivariate statistical theory[END_REF], [START_REF] Faraut | Analysis on symmetric cones[END_REF]) and we propose to extend it to a much more general setting.

In Section 2 of the paper we explain the details of the quadratic construction of Riesz measures on regular convex cones and next we define the corresponding Wishart laws. We compute their Laplace transforms, what is the starting point to get formulas for their expectation, covariance and higher moments (Theorems 2.5 and 2.6). More general Riesz and Wishart distributions associated to virtual quadratic maps are introduced in Section 2.3. Moments formulas are generalized (Theorem 2.9). Group equivariance of Wishart laws is studied at the end of the section.

Section 3 of the article is thoroughly devoted to the case when Ω is a homogeneous cone and the quadratic map q is homogeneous. A crucial role in the analysis of these maps and of related Riesz measures and Wishart laws is played by a matrix realization of any homogeneous cone, coming from [START_REF]On symplectic representations of normal j-algebras and their application to Xu's realizations of Siegel domains[END_REF] and explained in Section 3.2.

It allows, among others, to define basic and standard quadratic maps in Sections 3.3 and 3.4. They play a role of generators for homogeneous quadratic maps q needed to construct all Riesz measures and Wishart laws on Ω. Next we apply the results of [START_REF] Ishi | Positive Riesz distributions on homogeneous cones[END_REF] on Gindikin-Riesz distributions on Ω and on the orbit decomposition of Ω, the closure of Ω. We explain the relation between Riesz measures related to homogeneous quadratic maps and the Gindikin-Riesz distributions on Ω (Theorem 3.9). In Section 3.7, we prove the Bartlett decomposition for the Wishart laws on homogeneous cones (Theorems 3.11 and 3.13).

Families of Wishart laws that we construct and study in Section 3 comprise Wishart distributions studied in papers [START_REF] Andersson | Wishart distributions on homogeneous cones[END_REF] and [START_REF] Letac | Wishart distributions for decomposable graphs[END_REF] (homogeneous case) and are significantly bigger: we describe all singular Wishart laws and many more absolutely continuous Wishart laws than in papers [START_REF] Andersson | Wishart distributions on homogeneous cones[END_REF] and [START_REF] Letac | Wishart distributions for decomposable graphs[END_REF]. For the symmetric cone case, our Wishart laws cover the ones studied in [START_REF] Hassairi | Riesz exponential families on symmetric cones[END_REF] as well. All the results of Section 2 apply to them, in particular the formulas for the moments. Acknowledgement. We thank Professors Gerard Letac and Yoshihiko Konno for discussions on the topic of the article. §2. Riesz measure and Wishart law on a convex cone.

2.1. Regular cones and quadratic maps. In this paper, an open convex cone Ω ⊂ R n is always assumed to be regular, that is, Ω ∩ (-Ω) = {0}, where Ω denotes the closure of Ω. Then the dual cone Ω

* := η ∈ (R n ) * ; y, η > 0 (∀y ∈ Ω \ {0})
is a regular open convex cone again in the dual vector space (R n ) * , and we have (Ω * ) * = Ω. An R n -valued quadratic map q : R m → R n is said to be Ω-positive if (i) q(x) ∈ Ω for all x ∈ R m , and (ii) q(x) = 0 implies x = 0. These (i) and (ii) are restated in a single condition q(x) ∈ Ω \ {0} (∀x ∈ R m \ {0}). For the quadratic map q, we define the associated linear map φ

= φ q : (R n ) * → Sym(m, R) in such a way that t xφ(η)x = q(x), η (η ∈ (R n ) * , x ∈ R m ).
Then the Ω-positivity of q is equivalent to the following property of φ:

(2.1) η ∈ Ω * ⇒ φ(η) is positive definite. Example 1.
Let Ω be the open convex cone in R 3 defined by

Ω :=    t 1   0 0 1   + t 2   1 0 1   + t 3   1 1 1   + t 4   0 1 1   ; t 1 , t 2 , t 3 , t 4 > 0    =      y 1 y 2 y 3   ∈ R 3 ; y 1 > 0, y 2 > 0, -y 1 + y 3 > 0, -y 2 + y 3 > 0    . (2.2) If we identify (R 3 ) * with R 3 by y, η := y 1 η 1 + y 2 η 2 + y 3 η 3 (y, η ∈ R 3 ), we have Ω * =    t 1   1 0 0   + t 2   0 1 0   + t 3   -1 0 1   + t 4   0 -1 1   ; t 1 , t 2 , t 3 , t 4 > 0    =      η 1 η 2 η 3   ∈ R 3 ; η 3 > 0, η 1 + η 3 > 0, η 1 + η 2 + η 3 > 0, η 2 + η 3 > 0    ,
see [START_REF]The gradient maps associated to certain non-homogeneous cones[END_REF]. Let q : R 4 → R 3 be the quadratic map given by

q(x) := (x 1 ) 2   0 0 1   + (x 2 ) 2   1 0 1   + (x 3 ) 2   1 1 1   + (x 4 ) 2   0 1 1   (x ∈ R 4 ).
Clearly, this q is Ω-positive. By a simple calculation, we have

φ(η) =     η 3 0 0 0 0 η 1 + η 3 0 0 0 0 η 1 + η 2 + η 3 0 0 0 0 η 2 + η 3     (η ∈ (R 3 ) * ).
Example 2. Let Π r be the set of positive definite real symmetric matrices of size r. Then Π r is a regular open convex cone in the vector space Sym(r, R) of real symmetric matrices. If we identify the space Sym(r, R) with its dual vector space by the inner product y, η := tr (yη) (y, η ∈ Sym(r, R)), then the dual cone Π * r coincides with Π r . We define q r,s : Mat(r, s; R) → Sym(r, R) by q r,s (x) = x t x (x ∈ Mat(r, s; R)).

Then q r,s is Π r -positive. We denote the (i, j) component of x ∈ Mat(r, s; R) by

x r(j-1)+i , so that Mat(r, s; R) is identified with R rs . Then we have for η ∈ Sym(r, R)

φ(η) =     η η . . . η     ∈ Sym(rs, R),
where η is put s times. In this case, the map φ : Sym(r, R) → Sym(rs, R) is a Jordan algebra representation, and q is exactly the quadratic map associated to the representation ([4, Chapter IV, Section 4]).

Example 3. Let Z be a subspace of Sym(r, R), and put P := Z ∩ Π r . Then P is a regular open convex cone in Z. Let Q ⊂ Z * be the dual cone of P. We shall construct a Q-positive quadratic map q Z : R r → Z * whose associated linear map Then the quadratic map q Z : R r → Z * is given by q Z (x) := π Z * (x t x) (x ∈ R r ). In fact, for x ∈ R r \ {0} and y ∈ P we have (2.3) y, q Z (x) = tr (yx t x) = t xyx > 0 because y is positive definite. Therefore we get q Z (x) ∈ Q \ {0}, so that q Z is Qpositive. Keeping the natural isomorphism (Z

< i 2 < • • • < i k ≤ r,
and define

(2.6)

R I := { x ∈ R r ; x i = 0 if i / ∈ I } .
We denote by q I Z the restriction of q Z to the space R I ⊂ R r . Clearly q I Z : R I → Z * is Q-positive. The associated linear map φ q I Z : Z → Sym(k, R) gives a submatrix of elements y ∈ Z, that is, φ q I Z (y) = (y iαi β ), which we denote by y I .

Riesz measures and Wishart laws associated to quadratic maps.

For a regular open convex cone Ω ⊂ R n and an Ω-positive quadratic map q : R m → R n , let µ q be the image of the Lebesgue measure dx on R m by q. Namely, the measure µ q on R n is defined in such a way that (2.7)

R n f (y) µ q (dy) = R m f (q(x)) dx
for a measurable function f on R n . The Ω-positivity of q implies that the support of µ q is contained in the closure Ω of the cone Ω. By analogy to [4, Proposition VII.2.4], we call µ q the Riesz measure associated to q.

Lemma 2.1. Let φ : (R n ) * → Sym(m, R) be the linear map associated to q.

Then, for η ∈ Ω * , the Laplace transform L µq (-η) := R n e -η,y µ q (dy) of µ q equals π m/2 (det φ(η)) -1/2 .

Proof. By definition, we have L µq (-η) = R m e -t xφ(η)x dx. Since φ(η) is positive definite, the assertion follows from a formula of the Gaussian integral.

Definition 2.2. The members of the exponential family {γ q,θ } θ∈-Ω * generated by µ q are called the Wishart laws on Ω associated to q. Namely, (2.8) γ q,θ (dy) := e y,θ L µq (θ)

µ q (dy) (y ∈ R n ).
Remark 2.3. By (2.7), (2.8) and Lemma 2.1, we have for a measurable function

f on R n R n f (y)γ q,θ (dy) = π -m/2 (det φ(-θ)) 1/2 R m f (q(x))e -t xφ(-θ)x dx.
Putting Σ := φ(-θ) -1 and replacing the variable x by x/ √ 2, we rewrite the righthand side as (2π

) -m/2 (det Σ) -1/2 R m f (q(x)/2)e -t xΣ -1 x/2 dx.
Therefore, if X is an R m -valued random variable with the normal law N(0, φ(-θ) -1 ), then γ q,θ is nothing else but the law of Y := q(X)/2. In particular, the classical Wishart law as defined in [20, Definition 3.1.3] coincides with our γ q,θ in Example 2.

Proposition 2.4. Let Y be an R n -valued random variable with the Wishart law γ q,θ .

Then the Laplace transform L γ q,θ (η) = E(e Y,η ) of γ q,θ is given by

L γ q,θ (η) = det(I m + φ(-θ) -1 φ(-η)) -1/2 for η ∈ -θ -Ω * .
Proof. By definition, we have L γ q,θ (η) = L µq (θ) -1 L µq (η + θ). Thus the formula follows from Lemma 2.1 and the observation that (det φ(-θ)

) -1 det φ(-η -θ) = det φ(-θ) -1 (φ(-θ) + φ(-η)) = det(I m + φ(-θ) -1 φ(-η)).
We shall consider the mean and the covariance of the Wishart law γ q,θ . First we fix a notation for the directional derivative of a function, and recall some basic formulas. For a (vector-valued) function f on a domain U ⊂ R N , a point a ∈ U and a vector v ∈ R N , we denote by D v f (a) the directional derivative of f at a given by D v f (a) := ( d dh ) h=0 f (a + hv). Now we consider some functions and their derivatives on the domain -Π m in the vector space Sym(m, R). First we set f

(x) = log det(-x) ∈ R for x ∈ -Π m . Then we have (2.9) D v f (a) = -tr (-a) -1 v (a ∈ -Π m , v ∈ Sym(m, R)).
Next we observe the case

f (x) = (-x) -1 ∈ Sym(m, R) for x ∈ -Π m . Then (2.10) D v f (a) = (-a) -1 v (-a) -1 (a ∈ -Π m , v ∈ Sym(m, R)).
Theorem 2.5. Let Y be an R n -valued random variable with the Wishart law γ q,θ .

(i) For η ∈ (R n ) * , one has E( Y, η ) = tr φ(-θ) -1 φ(η)/2. (ii) For η, η ′ ∈ (R n ) * , one has E(( Y, η -M)( Y, η ′ -M ′ )) = tr φ(-θ) -1 φ(η)φ(-θ) -1 φ(η ′ )/2,
where M := E( Y, η ) and M ′ := E( Y, η ′ ).

Proof. Since {γ q,θ } θ∈-Ω * is the exponential family generated by µ q , it is well known ( [START_REF] Leonov | On a method of semi-invariants[END_REF]) that the mean E( Y, η ) is given by the derivative

D η log L µ (θ), while the co- variance E(( Y, η -M)( Y, η ′ -M ′ )) equals D η D η ′ log L µ (θ).
Thus the formulas follow from Lemma 2.1, (2.9) and (2.10) because φ(θ) ∈ -Π m .

Let us discuss higher moments of the Wishart law γ q,θ . The computation is reduced to the derivations of the Laplace transform

L µq (θ) = π m/2 (det φ(-θ)) -1/2 . Namely we have for η 1 , . . . , η N ∈ (R n ) * E( Y, η 1 Y, η 2 . . . Y, η N ) = D η 1 D η 2 . . . D η N L µq (θ) L µq (θ) .
On the other hand, if f (x) = det(-x) -p for x ∈ -Π m , where p is a real constant, the derivative Lemma 5]. We review the formula briefly. For an element π of the symmetric group S N , we write C(π) for the set of cycles of π. For y ∈ Π m , we denote by r π (y) the multilinear form of v 1 , . . . , v N given by

D v 1 D v 2 . . . D v N f (a) (a ∈ -Π m , v 1 , . . . , v N ∈ Sym(m, R)) is given in [8,
r π (y)(v 1 , . . . , v N ) := c∈C(π) tr j∈c yv j .
Then we have

D v 1 D v 2 . . . D v N f (a) = f (a) • π∈S N p ♯C(π) r π ((-a) -1 )(v 1 , . . . , v N ),
which is deduced from (2.9) and (2.10) by induction (see [START_REF] Graczyk | The complex Wishart distribution and the symmetric group[END_REF]). Making use of the formula, we obtain Theorem 2.6. Let Y be an R n -valued random variable with the Wishart law γ q,θ .

For

η 1 , η 2 , . . . , η N ∈ (R n ) * , one has E( Y, η 1 Y, η 2 . . . Y, η N ) = π∈S N 1 2 ♯C(π) r π (φ(-θ) -1 )(φ(η 1 ), φ(η 2 ), . . . , φ(η N )) = π∈S N 1 2 ♯C(π) c∈C(π) tr j∈c φ(-θ) -1 φ(η j ) .
2.3. Wishart laws associated to virtual quadratic maps. We shall consider virtual quadratic maps, that is a 'formal linear combination' of quadratic maps, and the associated Wishart laws. First we introduce the notion of direct sum of quadratic maps. Let q i : R m i → R n (i = 1, . . . , s) be Ω-positive quadratic maps.

Then the direct sum q

= q 1 ⊕ q 2 ⊕ • • • ⊕ q s is an R n -valued quadratic map on R m 1 ⊕ R m 2 ⊕ • • • ⊕ R ms given by q(x) := q 1 (x 1 ) + q 2 (x 2 ) + • • • + q s (x s ) x = s k=1 x i , x i ∈ R m i .
It is easy to see that q is also Ω-positive. If q 1 = q 2 = • • • = q s , then the direct sum q is denoted by q ⊕s 1 . The linear map φ : (R n ) * → Sym(m, R) (m := s i=1 m i ) associated to the direct sum q = ⊕ q i is given by

(2.11) φ(η) =     φ 1 (η) φ 2 (η) . . . φ s (η)     (η ∈ (R n ) * ).
Conversely, if a symmetric matrix φ(η) is expressed by φ 1 (η), . . . , φ s (η) as above for all η ∈ (R n ) * , then the corresponding quadratic map q is the direct sum of q 1 , . . . , q s . In Example 1, the quadratic map q : R 4 → R 3 is the direct sum of 4 quadratic maps [START_REF] Faraut | Analysis on symmetric cones[END_REF], where

q i : R ∋ x → x 2 v i ∈ R 3 (i = 1, . . . ,
v 1 :=   0 0 1   , v 2 :=   1 0 1   , v 3 :=   1 1 1   , q 4 :=   0 1 1   .
In Example 2, we see that q r,s : Mat(r, s; R) → Sym(r, R) is naturally identified with q ⊕s r,1 . In Example 3 with Z given by (2.4), we have q Z = q {1} Z ⊕ q {2,3} Z , while we do not have such a decomposition for the case (2.5).

Let q i : R m i → R n (i = 1, 2) be Ω-positive quadratic maps, and q the direct sum q 1 ⊕ q 2 . Then it is easy to see that the measure µ q equals the convolution µ q 1 * µ q 2 . Thus, for θ ∈ -Ω * we have L µq (θ) = L µq 1 (θ)L µq 2 (θ) and γ q,θ = γ q 1 ,θ * γ q 2 ,θ . In general, if we set q = q ⊕s 1 1

⊕ q ⊕s 2 2 ⊕ • • • ⊕ q ⊕st t
for Ω-positive quadratic maps q i : R m i → R n (i = 1, 2, . . . , t) and positive integers s 1 , s 2 , . . . , s t , then we have (2.12)

L µq (θ) = t i=1 L µ i (θ) s i (θ ∈ -Ω * ).
Now we remark that, even though s i 's are not positive integers, there may exist a positive measure µ q on Ω for which the relation (2.12) holds. In general, for real numbers s 1 , . . . , s r , we shall call a formal sum q = q ⊕s 1 1

⊕q ⊕s 2 2 ⊕• • •⊕q ⊕st t a virtual Ω-
positive quadratic map. The measure µ q satisfying (2.12) is called the Riesz measure associated to q. By the injectivity of the Laplace transform, the associated Riesz measure µ q is unique if it exists. In this case, the Wishart laws γ q,θ (θ ∈ -Ω * ) are defined again as members of the exponential family generated by µ q .

Proposition 2.7. Let q i : R m i → R n (i = 1, . . . , t) be Ω-positive quadratic maps.

Assume that there exists a measure µ q associated to the virtual quadratic map q = q ⊕s 1 1 ⊕ • • • ⊕ q ⊕st t for certain s 1 , . . . , s t ∈ R. Let Y be an R n -valued random variable with the Wishart law γ q,θ . Then the Laplace transform L γ q,θ (η) = E(e Y,η ) of the law γ q,θ is given by

L γ q,θ (η)= t i=1 L γ q i ,θ (η) s i = t i=1 det(I m i + φ i (-θ) -1 φ i (-η)) -s i /2
for η ∈ -θ -Ω * .

Proof. We have L γ q,θ (η) = L µq (θ) -1 L µq (η + θ) by definition, and the right-hand side equals

t i=1 L µ i (θ) -1 L µ i (η + θ) s i by (2.12). Since L γ q i , θ (η) = L µ i (θ) -1 L µ i (η + θ)
, we obtain the first equality. The second equality follows from Proposition 2.4.

Since we see immediately from (2.12) that

log L µq (θ) = t i=1 s i log L µ i (θ),
the virtual versions of Theorem 2.5 is given as follows:

Proposition 2.8. Under the same assumption of Proposition 2.7, one has

(i) E( Y, η ) = t i=1 s i tr φ i (-θ) -1 φ i (η)/2 for η ∈ (R n ) * , (ii) E(( Y, η -M)( Y, η ′ -M ′ )) = t i=1 s i tr φ i (-θ) -1 φ i (η)φ i (-θ) -1 φ i (η ′ )/2 for η, η ′ ∈ (R n ) * , where M := E( Y, η ) and M ′ := E( Y, η ′ ).
As for higher moments, we generalize the formula in Theorem 2.6 as follows:

Theorem 2.9. Under the same assumption of Proposition 2.7, one has

E( Y, η 1 Y, η 2 . . . Y, η N ) = π∈S N 1 2 ♯C(π) c∈C(π) t i=1 s i tr j∈c φ i (-θ) -1 φ i (η j ) . (2.13) for η 1 , η 2 , . . . , η N ∈ (R n ) * .
Theorem 2.9 easily follows from Theorem 2.6 when s 1 , . . . , s t are positive integers, that is, q is a true quadratic map. Indeed, we can consider the associated linear map φ = φ q in this case, and we have

tr j∈c φ(-θ) -1 φ(η j ) = t i=1 s i tr j∈c φ i (-θ) -1 φ i (η j )
by virtue of (2.11). To prove (2.13) for general case, it is enough to verify that the quantity E( Y, η 1 Y, η 2 . . . Y, η N ) is a polynomial of s 1 , . . . , s t . For this purpose, we make some calculations involving the semi-invariants or the cummulants (cf. [START_REF] Leonov | On a method of semi-invariants[END_REF]).

For η ∈ (R n ) * , θ ∈ -Ω * and a positive integer k, we denote by S k (θ; η) the k-th derivative of log L µq (θ) in the direction η:

S k (θ; η) = (D η ) k log L µq (θ) = d k dh k h=0 log L µq (θ + hη).
Similarly to Proposition 2.8, we get

(2.14) S k (θ; η) = (k -1)! 2 t i=1 s i tr φ i (-θ) -1 φ i (η) k .
On the other hand, since the function L µq (θ) is analytic, we have

log L µq (θ + hη) = log L µq (θ) + ∞ k=1 h k k! S k (θ; η).
Taking the exponential, we have (2.15)

L µq (θ + hη) = L µq (θ) exp ∞ k=1 h k k! S k (θ; η) = L µq (θ) ∞ ℓ=0 1 ℓ! ∞ k=1 h k k! S k (θ; η) ℓ .
Let Y be an R n -valued random variable with the Wishart distribution γ q,θ . Then we have

E( Y, η N ) = (D η ) N L µq (θ) L µq (θ) , so that (2.16) L µq (θ + hη) = L µq (θ) ∞ N =0 h N N! E( Y, η N ).
Comparing the coefficients of h N in (2.15) and (2.16), we obtain

E( Y, η N ) = N ℓ=1 1 ℓ! k 1 +k 2 +•••+k ℓ =N N! k 1 !k 2 ! • • • k ℓ ! S k 1 (θ; η)S k 2 (θ; η) . . . S k ℓ (θ; η).
Substituting (2.14) to the equality above, we obtain a formula for the higher moments.

Proposition 2.10. One has

E( Y, η N ) = N ℓ=1 1 ℓ! k 1 +k 2 +•••+k ℓ =N N! k 1 k 2 • • • k ℓ ℓ j=1 t i=1 s i 2 tr φ i (-θ) -1 φ i (η) k j .
(2.17)

The moments

E( Y, η 1 Y, η 2 • • • Y, η N ) (η 1 , . . . , η N ∈ (R n ) *
) can be obtained now by polarization, so that it is a polynomial of s 1 , . . . , s t . Hence the proof of Theorem 2.9 is completed. We note that (2.13) with

η 1 = • • • = η N = η becomes (2.17
), which means that Proposition 2.10 together with the polarization process also yields Theorem 2.9.

2.4. Group equivariance of the Wishart laws. Let G(Ω) be the linear automorphism group { g ∈ GL(n, R) ; gΩ = Ω } of Ω. For an Ω-positive quadratic map q : R m → R n and g ∈ G(Ω), the quadratic map g • q : R m → R n is again Ω-positive.

It is easy to see that the Riesz measure µ g•q is the image of µ q by g, that is,

(2.18) µ g•q (A) = µ q (g -1 A)
for a measurable set A ⊂ R n . Let us discuss the Wishart laws γ g•q,θ for θ ∈ -Ω * .

For η ∈ (R n ) * , we denote by g * η the linear form η • g ∈ (R n ) * . If η ∈ Ω * , then g * η ∈ Ω * because y, g * η = gy, η > 0 for y ∈ Ω \ {0}. We observe

(2.19) L µq (g * θ) = R m e q(x),g * θ dx = R m e g•q(x),θ dx = L µg•q (θ).
Therefore, denoting by 1 A the characteristic function of a measurable set A ⊂ R n , we have

(2.20)

γ g•q,θ (A) = 1 L µg•q (θ) R m 1 A (g • q(x))e g•q(x),θ dx = 1 L µq (g * θ) R m 1 g -1 A (q(x))e q(x),g * θ dx = γ q,g * θ (g -1 A).
We restate (2.20) as follows.

Lemma 2.11. Let g be an element of G(Ω). If a random variable Y obeys the Wishart law γ q,θ , the law of gY is γ g•q, (g -1 ) * θ .

Let q i : R m i → R n (i = 1, . . . , t) be Ω-positive quadratic maps, and q the virtual quadratic map q ⊕s 1 1 ⊕ • • • ⊕ q ⊕st t with s 1 , . . . , s t ∈ R. Then we define g • q to be the virtual quadratic map (g

• q 1 ) ⊕s 1 ⊕ • • • ⊕ (g • q t ) ⊕st .
Proposition 2.12. If the Riesz measure µ q exists, then the Riesz measure µ g•q exists and equals the image of µ q by g. Moreover γ g•q, (g -1 ) * θ is the image of γ q,θ by g.

Proof. Let µ ′ be the image of µ q by g. For θ ∈ -Ω * , we have

(2.21) L µ ′ (θ) = R n e y,θ µ ′ (dy) =
R n e gy,θ µ q (dy) = R n e y,g * θ µ q (dy) = L µq (g * θ).

By (2.12) and (2.19), the last term equals

t i=1 L µ i (g * θ) s i = t i=1 L µg•q i (θ) s i . Thus we get L µ ′ (θ) = t i=1 L µg•q i (θ) s i which means µ ′ = µ g•q by (2.
12). Let γ ′ be the image of γ q,θ by g. Similarly to (2.21), we have L γ ′ (η) = L γ q,θ (g * η) for η ∈ -(g -1 ) * θ -Ω * , while L γ g•q i ,(g -1 ) * θ (η) = L γ q i ,θ (g * η) by Lemma 2.11. On the other hand, we see from Proposition 2.7 that

L γ q,θ (g * η) = t i=1 L γ q i ,θ (g * η) s i = t i=1 L γ g•q i ,(g -1 ) * θ (η) s i = L γ g•q,(g -1 ) * θ (η).
Thus we get L γ ′ (η) = L γ g•q,(g -1 ) * θ (η), so that γ ′ = γ g•q,(g -1 ) * θ by the injectivity of the Laplace transform.

Let Aut(Ω, q) be the set of pairs (g 1 , g 2 ) ∈ G(Ω)×GL(m, R) for which g 1 •q = q•g 2 . Then Aut(Ω, q) forms a Lie subgroup of GL(n, R) × GL(m, R), and we have a group homomorphism pr 1 : Aut(Ω, q) ∋ (g 1 , g 2 ) → g 1 ∈ G(Ω).

The condition (g 1 , g 2 ) ∈ Aut(Ω, q) is also equivalent to

(2.22) φ q (g * 1 η) = t g 2 φ q (η)g 2 (η ∈ (R n ) * ).
Then we obtain

(2.23) det φ q (g * 1 η) = C det φ q (η) (η ∈ (R n ) * ) with C = (det g 2 ) 2
, which means that det φ q (η) is a relatively invariant polynomial on (R n ) * under the contragredient action of pr 1 (Aut(Ω, q)). The following proposition describes a transformation rule of the family of the Wishart laws {γ q,θ } θ∈-Ω * under the group pr 1 (Aut(Ω, q)).

Proposition 2.13. For a measurable set A ⊂ R n and (g 1 , g 2 ) ∈ Aut(Ω, q), one has

(i) µ q (g -1 1 A) = µ g 1 •q (A) = | det g 2 | -1 µ q (A), (ii) γ q,g * 1 θ (A) = γ q,θ (g 1 A).
Proof. (i) Because of (2.18), we only have to show the second equality. By definition, we have

µ g 1 •q (A) = R m 1 A (g 1 • q(x)) dx = R m 1 A (q • g 2 (x)) dx Putting x ′ = g 2 x, the last term equals | det g 2 | -1 R m 1 A (q(x ′ )) dx ′ = | det g 2 | -1 µ q (A), whence (i) follows. (ii) By (2.20), we get for y ∈ R n γ q,g * 1 θ (dy) = γ g 1 •q,θ (g 1 dy) = e g 1 y,θ L µg 1 •q (θ) µ g 1 •q (g 1 dy), Since µ g 1 •q = | det g 2 | -1 µ q by (i), the last term equals e g 1 y,θ | det g 2 | -1 L µq (θ) | det g 2 | -1 µ q (g 1 dy) = γ q,θ (g 1 dy).
Hence (ii) is verified. §3. Homogeneous Case.

3.1. Homogeneous quadratic map. An Ω-positive map q : R m → R n is said to be homogeneous if, for any y, y ′ ∈ Ω, there exists (g 1 , g 2 ) ∈ Aut(Ω, q) for which g 1 y = y ′ . In other words, q is homogeneous if pr 1 (Aut(q, Ω)) acts on Ω transitively.

In this case, Ω is clearly a homogeneous cone, that is, a linear group on R n acts on the cone Ω transitively. Then the dual cone Ω * ⊂ (R n ) * is also a homogeneous cone on which the group pr 1 (Aut(q, Ω)) acts transitively by the contragredient action ( [START_REF] Vinberg | Homogeneous cones[END_REF]). We see from (2.1) and (2.22) that the quadratic map q is homogeneous if and only if the associated linear map φ q : (R n ) * → Sym(m, R) is a representation of the dual cone Ω * in the sense of Rothaus [START_REF] Rothaus | The construction of homogeneous convex cones[END_REF] (see also [START_REF]Representation of clans and homogeneous cones[END_REF]).

A typical example of a homogeneous cone is Π r ⊂ Sym(r, R). For A ∈ GL(r, R),

we denote by ρ(A) the linear map on Sym(r, R) defined by ρ(A)y := Ay t A (y ∈ Sym(r, R)). Then the group ρ(GL(r, R)) acts on Π r transitively. Moreover, the linear automorphism group G(Π r ) equals ρ(GL(r, R)). We see that the quadratic map q r,s : Mat(r, s; R) → Sym(r, R) in Example 2 is homogeneous. Indeed, we have a surjective homomorphism GL(r, R) × O(s) ∋ (A, B) → (ρ(A), τ r,s (A, B)) ∈ Aut(Π r , q r,s ), where τ r,s (A, B) is a linear map on Mat(r, s; R) given by τ r,s (A, B)x := AxB -1 (x ∈ Mat(r, s; R)). Thus we have pr 1 (Aut(q r,s , Π r )) = ρ(GL(r, R)), which acts on Π r transitively. For a subset I ⊂ {1, . . . , r}, we denote by q I the restriction of q r,1 : R r → Sym(r, R) to the space R I ⊂ R r defined in (2.6). The map q I coincides with q I Z in Example 3 with Z = Sym(r, R). Let us observe that q I is homogeneous in general. Let P I be the linear group consisting of A ∈ GL(r, R) for which AR I = R I .

For example, if r = 3, we have

P {1} =    A =   a 11 a 12 a 13 0 a 22 a 23 0 a 32 a 33   ; A ∈ GL(3, R)    , (3.1) P {2,3} =    A =   a 11 0 0 a 21 a 22 a 23 a 31 a 32 a 33   ; A ∈ GL(3, R)    . (3.2)
Since we have a homomorphism P I ∋ A → (ρ(A), A) ∈ Aut(Π r , q I ), it is enough to show that ρ(P I ) acts on Π r transitively. Put k := ♯I and take a permutation matrix w 0 ∈ S r ⊂ GL(r, R) sending R {r-k+1,...,r} onto R I . Then we have P I = w 0 P {r-k+1,...,r} w -1 0 , and

P {r-k+1,...,r} = A 1 0 A 2 A 3 ; A 1 ∈ GL(k, R), A 2 ∈ Mat(r -k + 1, k; R) A 3 ∈ GL(r -k + 1, R) .
Since P {r-k+1,...,r} contains the group of lower triangular matrices, ρ(P {r-k+1,...,r} ) acts on Π r transitively. Therefore ρ(P I ) = ρ(w 0 )ρ(P {r-k+1,...,r} )ρ(w 0 ) -1 also acts on Π r transitively, so that q I is homogeneous.

Coming back to the examples (3.1) and (3.2), we note that q = q {1} ⊕ q {2,3} is not homogeneous as Π 3 -positive quadratic map, while both q {1} and q {2,3} are. Indeed, the image of q generates the space Z ⊂ Sym(3, R) in (2.4). Thus, if (g 1 , g 2 ) ∈ Aut(Π 3 , q), then g 1 must preserve both Z and Π 3 . Let us take y ∈ Π 3 \ Z. Then g 1 does not send I 3 ∈ Π 3 to y because I 3 ∈ Z. Thus the action of pr 1 (Aut(q, Π 3 )) on Π 3 is not transitive.

On the other hand, if we regard q as a map from R {1} ⊕ R {2,3} to Z, then q is a homogeneous P-positive quadratic map, where P := Z ∩ Π r . In fact, since (ρ(A), A) ∈ Aut(q {1} , P) ∩ Aut(q {2,3} , P) for A ∈ P {1} ∩ P {2,3} , we have ρ(A) ∈ pr 1 (Aut(q, P)). Therefore pr 1 (Aut(q, P)) contains a group ρ(P {1} ∩ P {2,3} ) which acts on P transitively.

In Example 1, the quadratic map q : R 4 → R 3 is not homogeneous because Ω ⊂ R 4 in (2.2) is not a homogeneous cone ( [START_REF]The gradient maps associated to certain non-homogeneous cones[END_REF]).

Matrix realization of homogeneous cones.

In this section, we shall discuss a homogeneous cone realized as P V = Z V ∩ Π N with Z V ⊂ Sym(N, R) constructed from an appropriate system V = {V lk } of vector spaces in a specific way explained below, following [14, section 3.1]. The investigation of such cones is fundamental because all homogeneous cones are linearly equivalent to some P V due to [START_REF]On symplectic representations of normal j-algebras and their application to Xu's realizations of Siegel domains[END_REF]Theorem D].

Let us take a partition N = n 1 + n 2 + • • • + n r of a positive integer N, and consider a system of vector spaces V lk ⊂ Mat(n l , n k ; R) (1 ≤ k < l ≤ r) satisfying the following three conditions:

(V1) A ∈ V lk , B ∈ V kj ⇒ AB ∈ V lj (1 ≤ j < k < l ≤ r), (V2) A ∈ V lj , B ∈ V kj ⇒ A t B ∈ V lk (1 ≤ j < k < l ≤ r), (V3) A ∈ V lk ⇒ A t A ∈ RI n l (1 ≤ k < l ≤ r).
Let Z V be the subspace of Sym(N, R) defined by

Z V :=        y =     Y 11 t Y 21 • • • t Y r1 Y 21 Y 22 t Y r2 . . . . . . Y r1 Y r2 • • • Y rr     ; Y kk = y kk I n k , y kk ∈ R (k = 1, . . . , r) Y lk ∈ V lk (1 ≤ k < l ≤ r)        .
We set P V := Z V ∩ Π N . Then P V is a regular open convex cone in the vector space Z V . Let H N be the group of real lower triangular matrices with positive diagonals, and H V a Lie subgroup of H N defined by

H V :=        T =     T 11 T 21 T 22 . . . . . . T r1 T r2 • • • T rr     ; T kk = t kk I n k , t kk > 0 (k = 1, . . . , r) T lk ∈ V lk (1 ≤ k < l ≤ r)        . If T ∈ H V and y ∈ Z V , then ρ(T )y = T y t T belongs to Z V thanks to (V1)-(V3).
Moreover ρ(H V ) acts on the cone P V ⊂ Z V simply transitively (cf. [14, Proposition

2.1]).

Keeping (V3) in mind, we define an inner product on the vector space V lk (1 ≤ k < l ≤ r) by the equality

(3.3) A t A = (A|A)I n l = A 2 I n l (A ∈ V lk ).
For y, y ′ ∈ Z V , we set

(3.4) y, y ′ := r k=1 y kk y ′ kk + 2 1≤k<l≤r (Y lk |Y ′ lk ),
where y kk and Y lk (respectively y ′ kk and Y ′ lk ) denote the components of y (respectively y ′ ). Note that the inner product is not equal to tr yy ′ unless n

1 = • • • = n r = 1.
By this coupling, we identify the dual space Z * V with Z V . Let us observe that I N belongs to the dual cone

P * V of P V , that is, 0 < y, I N = y 11 + • • • + y rr (y ∈ P V \ {0}).
Indeed, since each y kk is a diagonal entry of the non-negative matrix y ∈ P V \ {0}, we have y kk ≥ 0. Suppose r k=1 y kk = 0. Then y 11 = • • • = y rr = 0 and tr y = r k=1 n k y kk = 0. This together with the non-negativity of y implies y = 0, which is a contradiction.

For T ∈ H V , define ρ * (T ) ∈ GL(Z V ) by y, ρ * (T )η = ρ(T )y, η (y, η ∈ Z V ). By [START_REF] Vinberg | The theory of convex homogeneous cones[END_REF], we have

P * V = { ρ * (T )I N ; T ∈ H V }. Moreover, the map H V ∋ T → ρ * (T )I N ∈ P * V is a diffeomorphism. For σ = (σ 1 , . . . , σ r ) ∈ C r , we define the one-dimensional representation χ σ : H V → C × by χ σ (T ) := (t 11 ) 2σ 1 . . . (t rr ) 2σr (T ∈ H V ). Note
that any one-dimensional representation χ of H V is of the form χ σ , so that χ is determined by the values on the subgroup A V ⊂ H V consisting of diagonal matrices.

Let us give some examples. When n 1 = n 2 = • • • = n r = 1 and V lk = R for all 1 ≤ k < l ≤ r, the conditions (V1)-(V3) are clearly satisfied, and we have Z V = Sym(r, R) and P V = Π r . For the case r = 3,

n 1 = n 2 = n 3 = 1, V 21 = {0} and V 31 = V 32 = R, the space Z V equals Z in (2.5). Let us set r = 3, n 1 = 2, n 2 = n 3 = 1, V 21 = v 0 ; v ∈ R , V 31 = 0 v ; v ∈ R ,
and V 32 = {0}. Then we have

(3.5) Z V =           
P V = { y ∈ Z V ; y is positive definite } = y ∈ Z V ; y 11 > 0, y 11 y 22 -(y 21 ) 2 > 0, y 11 y 33 -(y 31 ) 2 > 0 ,
which is exactly the Vinberg cone [START_REF] Vinberg | The theory of convex homogeneous cones[END_REF].

Set r = 2, n 1 = m ≥ 1, n 2 = 1 and V 21 = Mat(1, m; R). Then Z V =            y 11 v 1 . . . . . . y 11 v m v 1 • • • v m y 22     ; y 11 , y 22 , v 1 , . . . , v m ∈ R        , P V = y ∈ Z V ; y 11 > 0, y 11 y 22 -(v 1 ) 2 -• • • -(v m ) 2 > 0 , (3.7)
so that we obtain the Lorentz cone of dimension m + 2.

Basic quadratic maps. Let W i

V (i = 1, . . . , r) be the subspace of Mat(N, n i ; R) consisting of matrices x of the form

x =     0 n 1 +•••+n i-1 ,n i X ii . . . X ri     X ii = x ii I n i , x ii ∈ R X li ∈ V li (l = i + 1, . . . ,

r) .

For example, when Z V is the one in (3.5), we have

W 1 V =            x 11 0 0 x 11 x 21 0 0 x 31     ; x 11 , x 21 , x 31 ∈ R        , W 2 V =            0 0 x 22 0     ; x 22 ∈ R        , W 3 V =            0 0 0 x 33     ; x 33 ∈ R       
, while for the case that Z V is the space Z in (2.5), we have

W 1 V =      x 11 0 x 31   ; x 11 , x 31 ∈ R    , W 2 V =      0 x 22 x 32   ; x 22 , x 32 ∈ R    , W 3 V =      0 0 x 33   ; x 33 ∈ R    .
For T ∈ H V and x ∈ W i V , we see from (V1) that T x ∈ W i V , which defines a representation τ i : H V → GL(W i V ). Since W i V ⊂ Mat(N, n i ; R), we can consider the restriction q i V of the Π N -positive quadratic map q N,n i : Mat(N, n i ; R) → Sym(N, R) in Example 2 to the space W i V . Thanks to (V2) and (V3), we have

q i V (x) = x t x ∈ Z V for x ∈ W i V . Then the quadratic map q i V : W i V ∋ x → x t x ∈ Z V is P V -positive. On
can deduce (iii) recursively.

(iv) It is known ([23, Chapter 3, Section 3] and [6, Section 1]) that a homogeneous cone is described as the subset of the ambient vector space consisting of points at which all relatively invariant (appropriately normalized) functions are positive. Thus the assertion follows from (iii).

Example 4. Let P V be the Vinberg cone, that is, Z V is as in (3.5). Then we have

φ 1 V (η) =   η 11 η 21 η 31 η 21 η 22 0 η 31 0 η 33   , φ 2 V (η) = η 22 , φ 3 V (η) = η 33 for η ∈ Z V . If η = ρ * (T )I 4 for T ∈ H V , we have (t 11 ) 2 (t 22 ) 2 (t 33 ) 2 = det φ 1 V (η), (t 22 ) 2 = η 22 , (t 33 ) 2 = η 33 , so that (t 11 ) 2 = det φ 1 V (η) η 22 η 33 = η 11 - (η 21 ) 2 η 22 - (η 31 ) 2 η 33 .
On the other hand, we have by Proposition 3.1 (iv)

P * V = η ∈ Z V ; η 11 η 22 η 33 -η 33 (η 21 ) 2 -η 22 (η 31 ) 2 > 0, η 22 > 0, η 33 > 0 .
Therefore, if Z is the space in (2.5), the linear isomorphism gives a bijection from P * V onto P = Z ∩ Π 3 . The adjoint map ι * : Z * → Z * V ≡ Z V gives the matrix realization of the Vinberg cone Q as the homogeneous cone P V .

ι :Z V ∋    
3.4. Standard quadratic maps and H-orbits in P V . We call the maps q ε V (ε ∈ {0, 1} r , ε = (0, . . . , 0)) standard quadratic maps. They are of particular importance among the virtual quadratic maps q s

V = (q 1 V ) ⊕s 1 ⊕ • • • ⊕ (q r V ) ⊕sr (s = (s 1 , . . . , s r ) ∈ R r ). Let us denote by I(ε) the set { 1 ≤ i ≤ r ; ε i = 1 }. We identify q ε V with a direct sum ⊕ i∈I(ε) q i V on the space W ε V := ⊕ i∈I(ε) W i V . Recalling (3.9), we have the isomorphism W ε V ≃ ⊕ i∈I(ε) (R ⊕ ⊕ l>i V li ), which enables us to describe y = q ε V (x) ∈ Z V (x ∈ W ε
V ) as the matrix composed of the blocks

(3.10) Y lk = i≤k, i∈I(ε) X li t X ki (1 ≤ k ≤ l ≤ r),
where X ii := x ii I n i for i ∈ I(ε).

q ε V (x) = T x t T x .
For example, if r = 3 and ε = (1, 0, 1), then an element x of

W ε V = W 1 V ⊕ W 3 V is of the form x =   X 11 X 21 X 31   ⊕   0 0 X 33   X 11 = x 11 I n 1 , X 33 = x 33 I n 3 x 11 , x 33 ∈ R, X 21 ∈ V 21 , X 31 ∈ V 31 ,
and we have

T x =   X 11 X 21 0 X 31 0 X 33   ∈ Mat(N, R).
For ε ∈ {0, 1} r , let E ε be the element of V given by

E ε :=   ε 1 I n 1 . . . ε r I nr   ,
and

O ε the H V -orbit ρ(H V )E ε ⊂ Z V through E ε .
In particular, the orbit O (0,...,0)

is the origin {0}, while O (1,...,1) = ρ(H V )I N equals the cone P V . It is shown in [START_REF] Ishi | Positive Riesz distributions on homogeneous cones[END_REF]Theorem 3.5] that the H V -orbit decomposition of the closure P V is given as

P V = ε∈{0,1} r O ε . Proposition 3.2. If ε = (0, . . . , 0), the image of the quadratic map q ε V equals the closure O ε of the orbit O ε . Proof. Let W ε,+ V be the subset { x ∈ W ε V ; x ii > 0 (i ∈ I(ε)) } of W ε V . For y = ρ(T )E ε = T E ε t T ∈ O ε with T ∈ H V , we take a unique x ∈ W ε,+
V for which T x equals T E ε . Then we have y = q ε V (x) by (3.12). Conversely, for any x ∈ W ε,+ V , we put T

:= T x + I -E ε ∈ H V so that q ε V (x) = T E ε t T ∈ O ε . Therefore we obtain O ε = q ε V (W ε,+ V ). On the other hand, putting W ε, * V := { x ∈ W ε V ; x ii = 0 (i ∈ I(ε)) }, we see easily that q ε V (W ε, * V ) = q ε V (W ε,+ V ). Since W ε, * V is an open dense subset of W ε V , the orbit O ε is dense in the image of the quadratic map q ε V , which is necessarily Proof. (i) If x ′ = τ ε (T )x ∈ W ε,+ V with x ∈ W ε,+ V , we have x ′ ii = t ii x ii for i ∈ I(ε). Thus (3.16) i∈I(ε) (x ′ ii ) 2u i -1 = χ u-ε/2 (T ) i∈I(ε) (x ii ) 2u i -1 .
On the other hand, we observe that

dx ′ = | det τ ε (T )|dx = i∈I(ε) det τ i (T ) dx,
and the last term equals i∈I(ε) χ m(i)/2 (T ) dx by Proposition 3.1 (i). Since

i∈I(ε) m(i) = r i=1 ε i m(i) = ε + p(ε),
we have dx ′ = χ ε/2+p(ε)/2 (T )dx, which together with (3.16) implies (i).

(ii) By definition, we have

x 2 = i∈I(ε) (x ii ) 2 + l>i X li 2 .
Thus the left-hand side of (3.15) equals

i∈I(ε) 2 Γ(u i ) ∞ 0 e -(x ii ) 2 (x ii ) 2u i -1 dx ii l>i V li e -X li 2 dX li π n li /2 .
Therefore we obtain (3.15) from li e -X li 2 dX li = π n li /2 and ∞ 0 e -(x ii ) 2 (x ii ) 2u i -1 dx ii = Γ(u i )/2. Remark 3.5. Lemma 3.4 tells us that e -x 2 M ε u (dx) is a probability measure on W ε,+

V . Actually, we see from the proof that if X u is an W ε,+ V -valued random variable with the law e -x 2 M ε u (dx), then its components are independent and satisfy √ 2X u li ∼ N(0, I n li ), and ( √ 2X u ii ) 2 ∼ χ 2 (2u i ), where χ 2 (u) denotes the chi-square law with the density 2 -u Γ(u/2) -1 e -t/2 t u-1 (t > 0). We shall see later that any Wishart law associated to a homogeneous quadratic map is the image of this measure by an appropriate quadratic map.

3.5. Gindikin-Riesz distributions. For σ = (σ 1 , . . . , σ r ) ∈ C r , we denote (σ r , . . . , σ 1 ) by σ * . Let ∆ * σ be the function on the cone P * V given by

(3.17) ∆ * σ (ρ * (T )I N ) = χ σ * (T ) (T ∈ H V ).
By Proposition 3.1 (ii) and (iii), ∆ * σ (η) can be expressed as a product of powers of the polynomials det φ i V (η). Putting

E η :=   η 1 I n 1 . . . η r I nr   ∈ Z V for η = (η 1 , . . . , η r ) ∈ R r >0 , we have (3.18) ∆ * σ (E η ) = (η 1 ) σr (η 2 ) σ r-1 . . . (η r ) σ 1 . For each σ ∈ C r , Gindikin ([6], [7]) constructed a tempered distribution R σ ∈ S ′ (Z V ) whose Laplace transform L Rσ (θ) = R σ (e y,θ
) is given by

(3.19) L Rσ (θ) = ∆ * -σ * (-θ) (θ ∈ -P * V ).
We call R σ the Gindikin-Riesz distribution on the homogeneous cone P V . The support of R σ is contained in P V , and R σ is relatively invariant under the action of ρ(H V ), that is,

(3.20) R σ (f • ρ(T )) = χ -σ (T )R σ (f ) for T ∈ H V and f ∈ S(Z V ).
Proposition 3.6. For non-zero ε ∈ {0, 1} r and u ∈ R + (ε), put σ := u + p(ε)/2. Then R σ is the image of M ε u by the standard quadratic map q ε V .

Proof. By (3.19), it is sufficient to show that

W ε,+ V e q ε V (x),θ M ε u (dx) = ∆ * -σ * (-θ)
for θ ∈ -P * V . Take T ∈ H V for which θ = -ρ * (T )I N . Then the left-hand side is

W ε,+ V e -q ε V (x),ρ * (T )I N M ε u (dx) = W ε,+ V e -q ε V (τε(T )x),I N M ε u (dx)
by (3.13), and it is equal to

χ -(u+p(ε)/2) (T ) W ε,+ V e -q ε V (x),I N M ε u (dx)
by Lemma 3.4 (i). Since q ε V (x), I N = x 2 by (3.11), we see from Lemma 3.4 (ii) that

W ε,+ V e q ε V (x),θ M ε u (dx) = χ -(u+p(ε)/2) (T ) = ∆ * -σ * (-θ). q : R 4 → Z V by q(x) := ι • q x 1 + ix 2 x 3 + ix 4
= ι (x 1 ) 2 + (x 2 ) 2 (x 1 x 3 + x 2 x 4 ) -i(x 1 x 4 -x 2 x 3 ) (x 1 x 3 + x 2 x 4 ) + i(x 1 x 4 -x 2 x 3 ) (x 3 ) 2 + (x 4 ) 2 =   (x 1 ) 2 + (x 2 ) 2 0 x 1 x 3 + x 2 x 4 0 (x 1 ) 2 + (x 2 ) 2 x 1 x 4 -x 2 x 3 x 1 x 3 + x 2 x 4 x 1 x 4 -x 2 x 3 (x 3 ) 2 + (x 4 ) 2   .

Then we have

φ q (η) =     η 1 0 η 3 η 4 0 η 1 -η 4 η 3 η 3 -η 4 η 2 0 η 4 η 3 0 η 2    
for η ∈ Z V . It is easily checked that the map φ q • ι : Herm(2, C) → Sym(4, R) is a Jordan algebra representation. For η ∈ Ω * we have

(3.26) L µq (-η) = π 2 (det φ q (η)) -1/2 = π 2 (η 1 η 2 -(η 3 ) 2 -(η 4 ) 2 ) -1
by Lemma 2.1. On the other hand, the basic quadratic maps q i V : W i V → Z V (i = 1, 2) are given by Thus we obtain

(3.27) L µ q 1 V (-η) = π 3/2 (η 2 ) -1/2 (η 1 η 2 -(η 3 ) 2 -(η 4 ) 2 ) -1/2 , L µ q 2 V (-η) = π 1/2 (η 2 ) -1/2
for η ∈ P * V . Comparing (3.26) and (3.27), we see that µ q = µ (q 1 V ) ⊕2 ⊕(q 2 V ) ⊕(-2) , whereas the quadratic map q is by no means equal to the virtual quadratic map (q 1 V ) ⊕2 ⊕ (q 2 V ) ⊕(-2) . We see also from (3.26) and (3.27) that µ q⊕(q 2 V ) ⊕2 = µ (q 1 V ) ⊕2 , whereas the two (true) quadratic maps q ⊕ (q 2 V ) ⊕2 and (q 1 V ) ⊕2 do not coincide even up to linear transforms g 0 ∈ G(Ω), as the domains of these maps are different.

Therefore two different quadratic maps may correspond to the same Riesz measure.

3.8. Density function for the non-singular case.

Since the orbit O ε = ρ(H)E ε is contained in the boundary ∂P V of the homogeneous cone P V unless ε = (1, . . . , 1), the Gindikin-Riesz distribution R σ is a singular measure for σ ∈ Ξ(ε) with ε = (1, . . . , 1) thanks to Proposition 3.6. On the other hand, if σ ∈ Ξ(1, . . . , 1), that is, σ i > p i /2 (i = 1, . . . , r),

where p i := p i (1, . . . , 1) = l>i n li , then the Gindikin-Riesz distribution is an absolutely continuous measure with respect to the Lebesgue measure, and the density function is given explicitly in [START_REF] Gindikin | Analysis in homogeneous domains[END_REF] as follows.

Noting that the group H V acts on P V simply transitively, we define the function ∆ σ : P → C × for σ = (σ 1 , . . . , σ r ) ∈ C r by ∆ σ (ρ(T )I N ) := χ σ (T ) (T ∈ H V ).

For y = ρ(T )I N = T t T ∈ P V , we can express ∆ σ (y) as a product of powers of principal minors of y (cf. [4, p. 122]). Define d = (d 1 , . . . , d r ) ∈ Z r /2 by d k := 1 + ( l>k n lk + i<k n ki )/2. Then ∆ -d (y)dy gives a G(P V )-invariant measure on P ([6, Proposition 2.2]). Take σ ∈ Ξ(1, . . . , 1). We see from [6, Theorem 2.1] that the integral Γ P V (σ) := P V e -y,I N ∆ σ-d (y) dy converges and equals π (dim Z V -r)/2 r i=1 Γ(σ i -p i /2). By [6, Proposition 2.3], we see that (3.28) R σ (dy) = ∆ σ-d (y) Γ P V (s) dy (y ∈ P V ).

Owing to (3.24) and (3.28), we conclude the following proposition.

Proposition 3.14. Let q s V = (q 1 V ) ⊕s 1 ⊕ • • • ⊕ (q r V ) ⊕sr be the virtual quadratic map such that σ = r i=1 s i m(i)/2 belongs to Ξ(1, . . . , 1), that is, σ i > p i /2 (i = 1, . . . , r). Then one has (3.29) γ q s V ,θ (dy) = e y,θ ∆ * σ * (-θ)∆ σ-d (y) Γ P V (s) dy (y ∈ P V ).

Note that the formula (3.29) served as a definition of a Wishart law in [START_REF] Andersson | Wishart distributions on homogeneous cones[END_REF].

Example 6. Let Z V be the space defined in (3.5). If y = ρ(T )I N = T t T ∈ P V with T ∈ H V , then we see easily that

y

  Throughout the paper and from its very beginning all our concepts are illustrated on important examples, including a non-homogeneous cone (Example 1 Section 2), symmetric cones of positive definite real symmetric matrices Π r (Example 2 Section 2), 4-dimensional Lorentz cone(Example 5 Section 3), and non-symmetric but homogeneous Vinberg cone and its dual(Example 3 Section 2).

φ

  Z : Z → Sym(r, R) equals the inclusion map. Let us define the surjective linear map π Z * : Sym(r, R) → Z * by y, π Z * (S) := tr yS (y ∈ Z, S ∈ Sym(r, R)).

y 11 0 y 21 0 0 y 11 0 y 31 y 21 0 y 22 0 0 y 31 0 y 33

 33 11 , y 22 , y 33 , y 21 , y 31 ∈ R

η 11 0 η 21 0 0 η 11 0 η 31 η 21 0 η 22 33 0 η 31 0 η 22 η 21 η 31 η 21 η 11 

 2211

  For the case l = k, we have Y kk = y kk I n k and (3.11)y kk = i≤k, i∈I(ε) X ki 2 thanks to (3.3), where we putX ii := |x ii |. For each x ∈ W ε V , let T x ∈ Mat(N,R) be a lower triangular matrix whose (k, i)-block component is X ki for k > i with i ∈ I(ε), and other components are zero. Then we have(3.12)

x 1 x 4 x 1 2 1 η 3 η 4 η 3 η 2 0 η 4 0 η 2  , φ 2 V

 1222 x 3 x 1 x 4 (x 3 ) 2 + (x 4 ) (η) = η 2 .

  = {i 1 , i 2 , . . . , i k } be a subset of {1, . . . , r} with 1 ≤ i 1

	Let I
	≃ Z in mind, we see from (2.3) that φ Z (y) = y (y ∈ Z). Soon later, we shall consider the cases    y 11 0   0 
	(2.4)	Z :=			0 y 22 y 32 0 y 32 y 33	 ∈ Sym(3, R) ; y 11 , y 22 , y 32 , y 33 ∈ R	
	and
			 	 y 11 0 y 31		 
	(2.5)	Z :=		 y 31 y 32 y 33 0 y 22 y 32	 ∈ Sym(3, R) ; y 11 , y 21 , y 31 , y 22 , y 33 ∈ R	
	as concrete examples. Actually, in the latter case (2.5), the cones Q and P are called the Vinberg cone and the dual Vinberg cone respectively, which are the lowest dimen-
	sional non-symmetric homogeneous cones ([23]). We shall see another realization of
	the Vinberg cone Q in (3.6) and the last paragraph of Section 3.3.

* ) *

  11 = (t 11 ) 2 , y 11 y 22 -(y 21 ) 2 = (t 11 ) 2 (t 22 ) 2 , y 11 y 33 -(y 31 ) 2 = (t 11 ) 2 (t 33 ) 2 ,

This research was partially supported by the grant ANR-09-BLAN-0084-01.

the other hand, we observe (3.8)

which yields the group homomorphism H V ∋ T → (ρ(T ), τ i (T )) ∈ Aut(P V , q i V ). It follows that the quadratic map q i V is homogeneous. We call q 1 V , . . . , q r V basic quadratic maps for P V . Recalling the inner product on V lk given by (3.3), we define an inner product on the space W i

V via the natural isomorphism (3.9)

Taking an orthonormal basis of W i V with respect to the inner product, we identify i) , where m(i) := dim W i V . Then we consider the linear map φ i V : Z V ≡ Z * V → Sym(m(i), R) associated to the quadratic map q i V . Note that, if we write n li for dim V li , we have m(i) = 1 + l>i n li . Proposition 3.1. (i) One has

where m(i) = (0, • • • , 0, 1, n i+1,i , . . . , n ri ) ∈ Z r .

(ii) For η = ρ(T ) * I N ∈ P * V with T ∈ H V , one has

In particular, (t rr ) 2 = det φ r V (η). (iii) For 1 ≤ i < r, there exist integers c i+1,i , . . . , c ri such that

(iv) One has P * V = { η ∈ Z V ; det φ i V (η) > 0 (i = 1, . . . , r) }.

Proof. (i) Since H V ∋ T → det τ i (T ) ∈ C × is a one-dimensional representation, it is sufficient to check the equality for diagonal matrices T ∈ A V . In this case, the isomorphism (3.9) gives the eigenspace decomposition of τ i (T ), where R and V li correspond to the eigenvalues t ii and t ll respectively. Therefore we have det

ri , whence we closed. Indeed, introducing the projective imbedding ι V of a vector space V by

In the proof of Proposition 3.2, we see that the quadratic map q ε V gives a surjective map from W ε,+ V onto O ε . The map is also one-to-one thanks to [START_REF] Ishi | Positive Riesz distributions on homogeneous cones[END_REF]Lemma 3.3 (ii)], while the map

Actually, a large part of the content of this section is presented in language of normal j-algebra in [START_REF] Ishi | Positive Riesz distributions on homogeneous cones[END_REF]Sections 3 and 4].

We define the representation τ ε :

Then we have by (3.8)

given by

(ii) One has

We set

If ε = (0, . . . , 0) and σ ∈ Ξ(ε), then R σ is a positive measure on the orbit O ε by Proposition 3.6. For the case ε = (0, . . . , 0), we have Ξ(0, . . . , 0) = {(0, . . . , 0)} and R (0,...,0) is the Dirac measure at the origin {0}. It is proven in [START_REF] Ishi | Positive Riesz distributions on homogeneous cones[END_REF] that they exhaust all the cases that R σ is a positive measure.

The parameter set Ξ is also described as

3.6. Riesz measures and Gindikin-Riesz distributions. Let us investigate a relation of the Riesz measures µ q associated to homogeneous P V -positive quadratic maps q and the Gindikin-Riesz distributions on P V .

Proposition 3.8. For i = 1, . . . , r, one has µ q i V = π m(i)/2 R m(i) . Proof. ¿From Lemma 2.1 and Proposition 3.1 (ii), we have

Assume that there exists the Riesz measure µ q s V associated to a virtual quadratic map q s

. By (2.12) and (3.21), we have for θ

We put

Then the equality above can be rewritten as

where |σ| := σ 1 +• • •+σ r . Thus µ q s V equals π |σ| R σ , so that σ belongs to Ξ(ε) for some ε ∈ {0, 1} r owing to Theorem 3.7. The converse argument is also valid. Therefore we obtain Theorem 3.9. For a virtual quadratic map q s

⊕sr , there exists the associated Riesz measure µ q s V if and only if σ := r i=1 s i m(i)/2 belongs to Ξ. In this case µ q s V = π |σ| R σ , and there exist ε ∈ {0, 1} r and u ∈ R + (ε) for which

Let q : R m → Z V be any homogeneous P V -positive quadratic map. As is noted in Section 3.1, the group pr 1 (Aut(P V , q)) acts on the cone P V transitively. Assume first that pr 1 (Aut(P V , q)) contains ρ(H V ). Then the polynomial det φ q (η) is relatively invariant under the action of ρ(H V ) by (2.23). Namely, for each T ∈ H V , there exists

, where C := det φ q (I N ). By (3.18), we have

which gives a practical way to determine m i . Indeed, we see from this formula that m i are non-negative integers. Comparing the degrees of both sides, we obtain m = m 1 + • • • + m r . Similarly to Proposition 3.8, we have

Let us consider the virtual quadratic map q ⊕s . The associated Riesz measure exists if and only if sm/2 ∈ Ξ, and in this case

As for the general case, we have the following result.

Proposition 3.10. Let q : R m → Z V be a homogeneous P V -positive quadratic map. Then there exist g 0 ∈ G(P V ), m ∈ Z r , and C > 0 for which

The Riesz measure µ q ⊕s associated to the virtual quadratic map q ⊕s exists if and only if sm/2 ∈ Ξ. In this case, µ q ⊕s equals the image of C -s/2 π sm/2 R sm/2 by g 0 .

Proof. We note that pr 1 (Aut(P V , q)) acts on the cone P V transitively, and that the identity component of pr 1 (Aut(P V , q)) equals the identity component of an algebraic group (cf. [START_REF]Representation of clans and homogeneous cones[END_REF]Theorem 2]). It follows that an Iwasawa subgroup (maximal connected split solvable subgroup) H of pr 1 (Aut(P V , q)) acts on P V simply transitively ([23, Chapter 1]). Since H is also an Iwasawa subgroup of G(P V ), it is conjugate to another Iwasawa subgroup ρ(H V ) ⊂ G(P V ). Namely, there exists g 0 ∈ G(P V ) for which g -1 0 Hg 0 = ρ(H V ). Let q ′ be the P-positive quadratic map g -1 0 • q : R m → Z V . We have

0 pr 1 (Aut(P V , q))g 0 ⊃g -1 0 Hg 0 = ρ(H V ). Thus we can apply the argument preceding Proposition 3.10 for q ′ , so that we have

Since q ⊕s = g 0 • (q ′ ) ⊕s , we get the last statement from Proposition 2.12. Proposition 3.10 states that the Riesz measure µ q associated to a homogeneous q is equal to some Gindikin-Riesz distribution up to a linear transform on G(P V ). On the other hand, Theorem 3.9 tells us that if a Gindikin-Riesz distribution is a positive measure, then it equals a Riesz measure associated to the virtual sum q s V = (q 1 V ) ⊕s 1 ⊕• • •⊕(q r V ) ⊕sr of basic quadratic maps up to a constant multiple. For example, let us recall the homogeneous Π r -positive quadratic map q I : R I → Sym(r, R) with I ⊂ {1, . . . , r} and the permutation matrix w 0 ∈ S r ⊂ GL(r, R) in Section 3.1. Putting g 0 := ρ(w 0 ), we have q I = g 0 • q {r-k+1,...,r} (k := ♯I), while q {r-k+1,...,r} is exactly the basic quadratic map q r-k+1 V for Z V = Sym(r, R). Therefore, the Riesz measure µ (q I ) ⊕s exists if and only if (0, . . . , 0, s/2, . . . , s/2

In this case, µ (q I ) ⊕s equals the image of π s/2 R (0,...,0,s/2,...,s/2) = µ (q r-k+1 V ) ⊕s by g 0 .

Let q 1 : R m 1 → Z V and q 2 : R m 2 → Z V be two homogeneous P V -positive quadratic maps. As we have seen in Section 3.1, the direct sum q 1 ⊕ q 2 is not necessarily homogeneous. Let us assume that the group pr 1 (Aut(P V , q 1 )) ∩ pr 1 (Aut(P V , q 2 )) acts on P V transitively. In this case, we see easily that q 1 ⊕ q 2 is homogeneous. As in Proposition 3.10, we can take g 0 ∈ G(P V ) for which g 0 ρ(H V )g -1 0 ⊂ pr 1 (Aut(P V , q 1 ))∩ pr 1 (Aut(P V , q 2 )). Then we have det φ q 1 ((g -1 0

. Now we consider a virtual quadratic map q = q ⊕s 1 1 ⊕ q ⊕s 2 2 . We see that the associated Riesz measure µ q exists if and only if s 1 m ′ + s 2 m ′′ /2 ∈ Ξ, and in this case, µ q is the image

by g 0 . Obviously, the same argument is valid for general quadratic maps q = q ⊕s 1 1

3.7. Bartlett decomposition of the Wishart laws. Let q : R m → Z V be a homogeneous P V -positive quadratic map. Then the Wishart law γ q,θ (θ ∈ -P * V ) is the image of the normal law N(0, φ(-θ) -1 ) on the vector space R m by the quadratic map q/2, see Remark 2.3.

However, this description of the Wishart law does not permit us to determine its support in general. In this section, we shall give another construction of the Wishart random matrices , which is a generalization of the Bartlett decomposition ( [START_REF] Bartlett | On the theory of statistical regression[END_REF], [START_REF] Muirhead | Aspects of multivariate statistical theory[END_REF]Theorem 3.2.14]) and has the advantage of controlling the support of the underlying Wishart law. Moreover, the result is valid for virtual quadratic maps.

First we consider the virtual quadratic map q s V = (q 1 V ) ⊕s 1 ⊕ • • • ⊕ (q r V ) ⊕sr whose associated Riesz measure µ q s V exists. Then σ = r i=1 s i m(i)/2 belongs to Ξ and we have µ q s V = π |σ| R σ by Theorem 3.9. Moreover, we have θ). Therefore we obtain from (2.8) that

We remark that distributions of this type are considered in [START_REF] Hassairi | Riesz exponential families on symmetric cones[END_REF] for the case when P V is a symmetric cone. Assume that γ q s V ,θ is not the Dirac measure. Then σ = (0, . . . , 0), so that we can take a non-zero ε ∈ {0, 1} r and u ∈ R + (ε) for which σ = u + p(ε)/2. Recall the standard quadratic map q ε V :

V is identified with a lower triangular matrix T x for which q ε V (x) = T x t T x . Thus, the W ε,+ V -valued random variable X u in the following theorem can be regarded as a triangular random matrix, similarly to the Bartlett decomposition of the classical Wishart distribution. Theorem 3.11. Let σ = u + p(ε)/2 and X u be an W ε,+ V -valued random variable whose components are independent and satisfy (X u ii ) 2 ∼ χ 2 (2u i ) and X u li ∼ N(0, I n li )

Proof. For a measurable function f on Z V , we see from (3.24) and Proposition 3.6 that

, by the change of variable of x by x/ √ 2, we rewrite the last term as

Keeping the Remark 3.5 in mind, we see that the law of the random variable X u is e -x 2 /2 M ε u (dx/2). Hence (i) holds. To show (ii), it suffices to check that γ q s V ,θ is the image of γ q s V ,-I N by ρ(T ) -1 . Since γ q s V ,-I N (dy) = e -y,I N R σ (dy), we have γ q s V ,-I N (ρ(T ) dy) = e -ρ(T )y,I N R σ (ρ(T ) dy) = e y,θ χ σ (T )R σ (dy) by (3.20). Therefore (3.24) together with (3.17) leads us to the assertion (ii). Now we consider the Wishart distribution γ q ⊕s ,θ , where q is a general homogeneous P V -positive quadratic map. First we show a refinement of the first part of Proposition 3.10. Lemma 3.12. Let q : R m → Z V be a homogeneous P V -positive quadratic map, and θ an element of -P * V . Then there exist g 0 ∈ G(P V ), m ∈ Z r and C > 0 for which (3.25) det φ q ((g

Proof. It is shown in the proof of Proposition 3.10 that there exists a 0 ∈ G(P V ) for which a 0 ρ(H V )a -1 0 ⊂ pr 1 (Aut(P V , q)). Since -a * 0 θ ∈ P * V , we take T 0 ∈ H V for which ρ * (T 0 )I N = -a * 0 θ. Put g 0 := a 0 ρ(T 0 ) -1 ∈ G(P V ). Then g 0 ρ(H V )g -1 0 ⊂ pr 1 (Aut(P V , q)) and g * 0 θ = ρ * (T 0 ) -1 a * 0 θ = -I N . Similarly to Proposition 3.10, we see that g 0 together with an appropriate m and C > 0 satisfies the required properties.

Assume that there exists the Riesz measure µ q ⊕s associated to a virtual quadratic map q ⊕s , and that µ q ⊕s is not the Dirac measure. Then we can take non-zero ε ∈ {0, 1} r and u ∈ R + (ε) such that sm/2 = u + p(ε)/2 as in Theorem 3.9. Using these data together with g 0 in Lemma 3.12, we obtain the Bartlett decomposition of the Wishart distribution γ q ⊕s ,θ . Theorem 3.13. Let sm/2 = u+p(ε)/2 and X u be the W ε,+ V -valued random variable in Theorem 3.11. Then the Wishart law γ q ⊕s ,θ is the law of Y = g 0 • q ε V (X u )/2.

Proof. Put q ′ := g -1 0 • q ⊕s . As is seen in the proof of Proposition 3.10, the Riesz measure µ (q ′ ) ⊕s equals C -s/2 π sm/2 R sm/2 . Thus, similarly to the proof of Theorem 3.11 (i), we see that γ (q ′ ) ⊕s ,-I N (dy) = e -y,I N R sm/2 (dy), and that γ (q ′ ) ⊕s ,-I N is the law of Y = q ε V (X u )/2. Since q ⊕s = g 0 • (q ′ ) ⊕s , Theorem 3.13 follows from Proposition 2.12.

We have seen that Riesz measures and Wishart laws associated to a homogeneous quadratic map are obtained (up to linear transforms as in Proposition 3.10 and Theorem 3.13) as the ones associated to a virtual quadratic map q s

V ) ⊕sr , that is, a virtual sum of basic quadratic maps. However, it does not mean that every homogeneous quadratic map is equal to a direct sum of basic quadratic maps. The structure of homogeneous quadratic maps is more rich than the maps generated by basic quadratic maps. Let us study the following example. Since ρ(GL(2, C)) acts on Ω transitively, the quadratic map q is homogeneous. Keeping the natural isomorphism C 2 ≃ R 4 in mind, we define the quadratic map so that ∆ σ (y) = (t 11 ) 2σ to the virtual quadratic map (q 1 V ) ⊕s , where q 1 V : W 1 V → Z V is the basic quadratic map. Since m(1) = (1, 1, 1), we observe that sm(1)/2 ∈ Ξ if and only if s ∈ {0, 1} ∪ (1, +∞). If s = 0, the associated Wishart law is the Dirac measure. If s = 1, the associated Wishart law γ q 1 V , θ (θ ∈ -P * V ) is described as the image of the normal law N(0, φ 1 V (-θ) -1 ) on W 1 V ≡ R 3 by the quadratic map q 1 V /2, where φ 1 V (-θ) is given in Example 4 after Proposition 3.


C Krzeminski

Q Brulin

V Cuny

E Lecat

E Lampin

F Cleri

Molecular dynamics simulation of the recrystallization of amorphous Si layers: Comprehensive study of the dependence of the recrystallization velocity on the interatomic potential

The molecular dynamics method is applied to simulate the recrystallization of an amorphous / crystalline silicon interface. The atomic structure of the amorphous material is constructed with the method of Wooten, Winer and Weaire. The amorphous on crystalline stack is afterwards annealed on a wide range of temperature and time using five different interatomic potentials : Stillinger-Weber, Tersoff, EDIP, SW115 and Lenosky. The simulations are exploited to systematically extract the recrystallization velocity. A strong dependency of the results on the interatomic potential is evidenced and explained by the capability of some potentials (Tersoff and SW115) to correctly handle the amorphous structure while other potentials (Stillinger-Weber, EDIP and Lenosky) lead to the melting of the amorphous. Consequently, the interatomic potentials are classified according to their ability to simulate the solid or the liquid phase epitaxy.

I. INTRODUCTION

The pre-amorphization of silicon by high dose ion implantation is one of the identified processes for the realization of the ultra-thin junctions of the technological nodes 65 nm and below of the silicon microelectronics industry. The transformation of crystalline (c-Si) into amorphous silicon (a-Si) on a layer close to the surface reduces the channeling of the dopant consecutively implanted. Such a technological stage is therefore required to limit the penetration into the substrate of small radius atoms such as boron, the standard p dopant in microelectronics [1,2]. The pre-amorphization and the dopant implant might be followed:

• by an annealing at a temperature low enough (< 700 • C) to recrystallize the amorphous layer and activate the doping, but insufficient to cause a significative diffusion of the dopant. Reference 3 gives such an example, with a depth of the junction at 10 18 cm -3 of [16][17][18][19] nm compatible with the values [10][11][12][13][14][15][16][17][18][19] nm targeted by the ITRS roadmap [START_REF]International Technology Roadmap for Semiconductors[END_REF]. Likewise, the high level of activation allows to reach sheet resistances of [405 -785] Ω/ to compare with the range [760-830] Ω/ of the ITRS Roadmap.

• by an annealing at very high temperature (laser or spike annealing). In the case where the melting temperature of a-Si is reached, the profile after annealing has an ideal "square-like" shape as the diffusivity of the dopant is higher in the liquid than in the crystal. An example of such a profile is given in the reference 5. While ideal for a full wafer implantation, this solution however is not always compatible with the other technological stages [6]. The other alternative is to anneal just beneath the melting temperature (ex : 1300 • C) during a very short duration (ex : spike of 330 µs wide at 1250 • C) to limit the dopant diffusion [7].

In the two cases, it is critical to have a good understanding of the recrystallization process of the amorphous layer to develop such technologies. If, on the one hand, the recrystallization velocity by solid phase epitaxy (SPE) of an amorphous layer is well known from experiments in the intrinsic case [16], this is not the case in presence of impurities where accelerations, slowdown, compensation and even hysteresis are observed [16]. Moreover in this case, the final distribution of the dopant is intimately related to the velocity of recrystallization. Now, this case of doped a-Si is precisely of current interest for the formation of junctions. In the absence of an experimentally established law, is it possible to calculate at the atomic level the recyrstallisation velocity of doped a-Si ? Atomistic simulations of a-Si, of its interface with c-Si and of its recrystallization have already been published [8][9][10][11][12][13][14][15]. However in these works the atomistic simulations of the recrystallization of a-Si are not always exploited to extract the recrystallization velocity v cryst . The comparison with experiments is most often limited to the activation energy of v cryst . Moreover, the description of silicon by means of empirical interatomic potential is not unique and the question of the most adapted one(s) for this particular case is still open. The aim of this paper is to present a systematic study of the recrystallization velocity as a function of the annealing temperature for the various interatomic potentials available for silicon. The present work has been devoted to the intrinsic case, that already appears an issue by itself.

In this paper, we begin by a presentation of the method of simulation we have used to model the a-Si, its recrystallization and the extraction of the recrystallization velocity. The results are then presented for five different interatomic potentials in the second part. In the last part, the wide differences obtained as a function of the interatomic potential are analyzed and an explanation is proposed based on a possible transformation of the amorphous layer into a liquid one.

II. METHOD

This section is first devoted to a discussion about the different empirical interatomic potentials used in the molecular dynamics simulation. Next, the description of the method used to build the a-Si cluster and to its properties is presented. The generation of the a-Si/c-Si structure is also reported. Finally, the simulation of the recrystallization of a-Si on c-Si with molecular dynamics is described. The systematic exploitation of these atomistic simulations to determine the recrystallization velocity is also explained in details.

A. Molecular dynamics using interatomic potentials

In this article, the annealing of the a-Si on c-Si stack was investigated by means of molecular dynamics simulation employing empirical interatomic potentials. Interatomic potentials are an efficient way to describe the interaction between atoms using a continuous and derivable functional of the atomic positions. As they do not contain deep physical effects as in ab-initio description, their domain of validity is not universal. Most often, interatomic potential are developed on a particular application such as the study of surfaces or the description of a phase transition. Five empirical potentials developed for silicon were tested in this study : the Stillinger-Weber potential in its original parametrisation [23] and in a recent parametrisation called SW115 [26] designed for a better description of the melting of a-Si and c-Si, the third formulation of the Tersoff's potential [24], the Environment Dependent Interatomic Potential (EDIP) [25] and the Lenosky's pontential [27]. To our knowledge, these are the main potentials used in the literature to describe silicon. The 4 functional forms (Stillinger-Weber and SW115 have the same formulation) contain an attractive and a repulsive part given in the form of continuous functions of the relative distance between pairs of atoms, with a two-body and a three-body character, respectively. These functional forms are adjusted to reproduce the T = 0 equilibrium crystalline structure, cohesive energy, elastic constants and possibly other structural or energetic parameters of crystalline Si. The Stillinger-Weber potential has a three-body term that imposes the angle θ between the three neighbors to cos θ = -1/3 and therefore strictly enforces the diamond structure. The EDIP potential is inspired from the Stillinger-Weber potential, but it adapts to the coordination (4, 3, 2,...) of the atom and it is therefore more general. The formulation of the Tersoff potential is slightly different from the two previous ones, and has the form of a pair potential with a bonding strength that contains the dependence upon the local atomic environment.

Finally, the Lenosky is a modified embedded-atom method including a pairwise embedding functional and an angular term. Which potential should be more adapted to our case is not obvious a priori. Therefore simulations with these 5 different potentials were performed.

B. Amorphous cluster

Amorphous silicon in silicon devices is created by ion implantation of heavy ions. Ideally, the best approach would be to mimic the experiment and therefore to simulate the ion implantation. However, it is not currently realistic since the simulation of the implantation of an ion and of all the interactions and collisions with the atoms of the target by molecular dynamics is already computer resources consuming, and many events of that kinds should be simulated to obtain an a-Si layer of good quality. The other alternative is to construct a cluster of amorphous silicon and to glue it on crystalline silicon. The main two methods to form an amorphous structure starting from the crystal are a) the disordering of tetrahedral interatomic bonds by the method of Wooten, Winer and Weaire (WWW) [17] and its derivatives [18,19] and b) the simulation of an idealized a-Si, obtained by a melting followed by a quenching to a solid but amorphous state [20,21]. We will call this second method the "melting-quenching". Each of these two methods has its advantages and drawbacks. The main disadvantage of the WWW method is to enforce a tetravalent coordination, while it is experimentally known [22] that a-Si contains defects which may be floating or dangling bonds. The simulation of the melting-quenching on the other hand, done using a molecular dynamics simulation with an empirical potential, relies first on the ability of the potential to correctly describes a liquid at elevated temperature, and second to correctly describe the amorphous structure upon quenching, which is not unconditionally ensured for neither of these two aspects as emphasized later (see Section IV).

In the present work, the results are obtained with a WWW a-Si. The reason for this choice is first to construct an amorphous cluster of good quality, and next paragraph will discuss of this point. Moreover, the objective was to test which of the numerous empirical interatomic potential developed to model silicon [23][24][25][26][27] is the more adapted to keep the amorphous structure and to give the recrystallization velocity measured by solid phase epitaxy experiments [16]. Using the melting-quenching method would have required to construct for each potential its corresponding amorphous cluster, which, while more rigorous than the use of a unique WWW amorphous cluster, would have introduced different amorphous structures. Therefore all the results presented in the section III have been obtained with the same initial a-Si cluster constructed by the method presented just below. On the other side, we believe that, once the choice of the most adapted interactomic potential has been made, the better approach should on the contrary be the simulation of the melting-quenching-recrystallization with this potential.

Construction

The WWW [17] a-Si structure is obtained using an algorithm of bond switching between pairs of first neighbours. The total energy associated to this bond reordering is originally calculated using a Keating potential [28]. A Maxwell-Boltzmann probability is then determined with this energy and a temperature. In this way, an elevated value of the temperature can be used at the beginning of the amorphisation to strongly disorder the structure, and be subsequently decreased to reduce the number of defects in order to obtain a very realistic a-Si. This procedure was applied to an initial cluster of 4096 atoms of c-Si.

Properties

The current paragraph presents the properties of the WWW a-Si cluster of 4096 atoms and the comparison with published experiments.

The side of the initial cubic cluster of c-Si is of L c-Si =43.44 Å. The density before amorphisation is therefore equal to d c-Si =4096/L c-Si = 0.0500 at/ Å3 . After amorphisation, the dimensions of the cluster are slightly reduced: L a-Si =43.38 Å. The density is therefore 0.4% higher (d c-Si =0.0502 at/ Å3 ). This feature is in disagreement with experiments where a reduction of the density of approximately 1.8 % has been measured [29]. This difference might be due to the dilatation and to the strain induced by the implantation at high dose and high energy (> MeV) used in the experimental articles [31], both effects neglected in the model.

It is not obvious to simulate these macroscopic effects at the atomic level. Moreover, it is not sure that this rather slight difference is significative of a problem in the model.

The angular distribution of the bonds between first neighbours is shown in Fig. 2. The mean angle of the distribution is equal to 109.1 • and the standard deviation is equal to 10.6 • , in excellent agreement with the work of Djordevic et al. [30].

The pair correlation function g(r), i.e. the probability to find an atom between the spheres of radius r and r + dr centered on any atom is given in Fig. 3. g(r) is calculated using the radial distribution G(r) :

g(r) = G(r) 4πr 2 ρ 0 (1)
where ρ 0 is the density of c-Si. The radial distribution is the number of neighbours at the distance r :

G(r) = 1 N i,j δ(r -|R i -R j |) (2)
where |R i -R j | is the distance between atoms i and j and N the total number of atoms.

The abscissa of the maxima of the pair correlation function (see Fig. 3) corresponds to the distances between first (2.35 Å), second (3.83 Å), third (5.43 Å) and even fourth (6.62 Å) neighbours. The peaks are of course enlarged due to the disorder in a-Si compared to are clearly separated, with a gap to g(r)=0 between the two peaks. The pair correlation function is not directly comparable to experiments. Therefore its Fourier transform was calculated to extract the structure factor S(Q), given in Fig. 4. The shape, the position and the intensity of the peaks is in perfect agreement with the structure factor measured by Laaziri et al. [32,33] on a-Si obtained by ion implantation, except a slight underestimation of the intensity of the first peak. This agreement is significative of the good description of the local order in the WWW a-Si.

To conclude, although the discrepancy on the density, the WWW a-Si has excellent geometrical characteristics and the structure factor is in good agreement with experiments.

Therefore the quality of the a-Si structure is adapted to the objective of this study.

C. Construction and properties of the a-Si/c-Si structure

The resulting amorphous cluster of 4096 atoms (43.38 Å each side) was then truncated and glued along the [100] surface of the c-Si cell to form the interfacial structure [11,15,43].

The structure contains 1905 atoms and the crystal part of the interfacial structure corresponds to 50% of the global thickness. It has been observed that the crystalline thickness is an important factor since it appears in some cases that a recrystallization of the a-Si layer was obtained for the bigger structure while the disorder propagated in the c-Si part of the smaller structure (21% of the global thickness) although the two atomic layers at the bottom of the structure were frozen at their perfect diamond lattice site to eventually compensate the thin thickness of the crystalline part. Tests on bigger stack (8034 atoms) did not give significantly different recrystallization velocities than with the stack of 1905 atoms. Since the computational time is drastically increased for this stack of 8034 atoms, we dropped it and focused the work on the 1905 atoms interfacial system. Table I presents the mean coordination number and the silicion defects distribution of the initial interface.

The statistics have been estimated for all the silicon atoms located at less than 2.35 Å of the interface both for the cristal and the amorphous part. A cut-off distance of 2.7 Å has been used to calculate the first neighbours. Large amount of threefold and twofold silicons defects are present initially due to the difficulty to perfectly align an amorphous interface and a perfect crystal.

In order to improve the interface quality, prior to the main recrystallisation step, an equilibration step has been assumed. An objective is to relax the interface and to reach the temperature targeted value for the main recrystallisation simulation. Periodic boundary conditions only along the lateral directions have been used allowing the system to relax along the regrowth direction. The equilibration phase is performed during 4 ps at the targeted temperature value using a 2 fs time step for the velocity Verlet algorithm. During this equilibration, the temperature is imposed by the scaling of the velocity which ensures an efficient convergence of the temperature as shown in Fig. 5. The empirical atomic potential used in this pre-annealed step is identical to the one used for the recrystallisation simulation.

Table I presents the evolution of the mean coordination number distribution for the different interatomic potentials. First of all, it has been checked that during this equilibration step, no real recrystallisation occurs (the interface variation is less than 0.4 Å for all the interatomic potential). The equilibration step improves the quality of the interface since the fraction of fourfold silicon atoms is enhanced for almost all interatomic potentials. On the other side, the choice of the interatomic potential has a clear influence on the defects at the interface even for a very short annealing. After equilibration, both Tersoff and SW115 tend to favorize the presence of threefold defects. On the other side, fivefold silicon atoms are the most dominant defects with the SW, EDIP and Lenosky potentials.

D. Recrystallization of the amorphous layer

The recrystallization of the amorphous layer has been simulated after the equilibration step in an second step where a Nosé-Hoover thermostat [34] is applied. The acquisition begins during a number of time steps sufficient to observe the recrystallization of all or of one a part of the a-Si layer. During this phase, the positions of the atoms are saved to allow a post-treatment explained in the following part. The interval between two backups is adapted to the speed of the recrystallization : from 10 ps to 1000 ps. Once the simulation has ended, i. e. the amorphous layer is partially or totally crystallized, the positions stored during the simulation are exploited to extract the recrystallization velocity. Several methods have been used in the literature to quantify the motion of the recrystallization front. Bernstein et al.

[10] classify the atoms as "crystalline" or "amorphous" on the basis of a geometric criterion called χ and equal to the sum of the absolute values of the dot products between bond directions around each atom and the bond directions in the underlying crystal lattice. This parameter χ is equal to 4 in the crystal and 2.7 in the a-Si. The velocity of recrystallization is extracted from the number of crystal atom versus time. Marqués et al. [14] extract the velocity in atom/ps from the slope of the potential energy curves, since the recrystallization is accompanied by a decrease of this energy. Motooka et al. [11] only specify that the a-Si/c-Si interface is defined as a middle of their rough interface, without additional information on a criteria used. Albenze et al. [35] define an order parameter equal to 1 in solid, 0.7-0.8 in an amorphous state and 0.2-0.4 in liquid, from the angles between each atom and its first neighbours. Next, the sample is divided into slices, and the parameter is calculated inside each slice which is consecutively classified into one of the three states. This method helps to follow the displacement of a crystalline/amorphous interface. In our opinion, the most rigorous analysis is made by Mattoni and Colombo [15]. They calculate a structure factor along the growth direction S(z) that simply quantify the crystalline order in 2D slices in the z direction. This structure factor is equal to 1 when the plane is ordered (case of a crystal at 0 K) and around 0 when the order is lost. A detection of S(z) = 0.5 is used to easily place the a-Si/c-Si interface and further deduce the recrystallization velocity. The same criterion has been used to extract the recrystallization velocity from our simulations. Fig. 6 gives an example of structure factor for a stack of 8034 atoms. The structure factor along the growth direction is aligned with the representation of the atoms to emphasize the link between order/disorder and the value of S(z). The interface position obtained at S(z) = 0.5 is reported versus time in Fig. 7. Until ≈ 240 ps, the position increases linearly with time :

the disordered layer recrystallizes at a roughly constant velocity. After 240 ps, the position saturates, the whole layer has recrystallized. The velocity of recrystallization is estimated from the fit of the linear region (see Fig. 7). In this particular case, a value of 0.015 nm/ps was found. This procedure is applied to 50 to 100 files of 1905 atomic positions, for 4 to 6 temperatures for each of the 5 potentials. The evolution of the a-Si on c-Si was studied up to ≈ 100 ns for the slowest cases, a computer job of two months.

III. RESULTS AND COMPARISON WITH PUBLISHED SIMULATIONS

The results of our molecular dynamics simulation of the recrystallization of an a-Si on c-Si stack are presented in Fig. 8 for the Stillinger-Weber, Tersoff, EDIP, SW115 and Lenosky potential. This figure evidences a strong dependence of the velocity of recrystallization on the interatomic potential. At 1200 • C, there is for example 4 orders of magnitude separating the results with the Stillinger-Weber potential from the results with the Tersoff potential.

These differences might also be even more pronounced at low temperature, but the limit of our simulations at time being is 2×10 -6 nm/ps (one monolayer recrystallized after 100 ns)

and does not allow to get additional points for the Tersoff potential.

Very few data are available in the literature to perform a comparison with our calculations.

As for the Stillinger-Weber potential, the only point to our knowledge can be deduced from Fig. 1 of Mattoni and Colombo [15] : at 1200K, we have estimated a velocity of ≈ 2×10 -3 nm/ps. However, their calculation was made i) on a system obtained by melting-quenching with the Stillinger-Weber potential, although this procedure is known to give an amorphous layer of poor quality [36] and ii) including B atoms in the Si lattice to study their clustering upon annealing. It is known that dopant atoms modify the velocity of the recrystallization [16]. The direct comparison with our calculations would therefore not be conclusive, although the order of magnitude is roughly comparable (around 10 -2 nm/ps in our case versus the value of 2×10 -3 nm/ps in the paper of Mattoni and Colombo [15]).

Concerning the simulation of the recrystallization with the Tersoff potential, the velocities calculated by Marqués et al. [14] in atom/ps is not useful for the comparison with our simulations. Motooka et al. [11] give the velocity as a function of the inverse temperature.

Their results are reported in Fig. 9. Their simulations with an a-Si system obtained by melting-quenching give two regimes, a low temperature regime with an activation energy of 2.6 eV and a high temperature regime with an activation energy of 1.2 eV. On the contrary, experiments give an Arrhenius behavior from 470 to 1350 • C with a single activation energy of 2.68 eV [16]. We do not observe these two regimes, but a single one with an activation energy of 2.99 eV. Nevertheless, our results are in perfect agreement with the results of Motooka et al. in the low temperature regime (see Fig. 9).

Bernstein et al. [10] have published a figure of the recrystallization velocity as a function of the inverse temperature between 700 and 1100 K for the EDIP potential. They also obtained two regimes. The low temperature regime ( < 950 K) has a very low energy of activation (0.4 eV) [START_REF]Atomistic scale structure and dynamics of amorphous-crystal interfaces in silicon: simulations with empirical and quantum-mechanical approaches[END_REF]. The regime at high temperature has an activation energy of ≈ 2.0 eV. The fit of their points in this regime is given in Fig. 9 together with our results.

Our result is 3 times higher at 1100 K and the difference increases when the temperature decreases. However a direct comparison is difficult since different methods have been used.

For example, the velocity of recrystallization is calculated by identification of the "crystal"

and "amorphous" atoms, not with a structure factor. Moreover the amorphous part of the atomic structure is generated by a procedure of melting-quenching.

To our knowledge, there are no published values for the SW115 and Lenosky potentials.

To conclude, the comparison between our simulations and the literature is difficult. Many differences exist in the method to generate the amorphous layer or to estimate the interface position. However it can be noticed a relatively good agreement with the simulations of Motooka et al. [11] for the low temperature regime.

IV. DISCUSSION

Although the results of the simulations presented in Fig. 8 seems at first glance very disparate, a deeper analysis allows to regroup them under two categories. The first category is made up of the results obtained with the Tersoff and SW115 potentials. Both interatomic potentials give an Arrhenius behavior and smaller values of the velocity of recrystallization.

The second category includes the Stillinger-Weber, the EDIP and the Lenosky potential.

Using these three potentials, the velocity of recrystallization does not follow an Arrhenius law any more but presents a maximum followed by a decrease of the velocity at high temperature.

Moreover, the values of the velocities are higher in this second case. In this section, we firstly compare the results of each category with the measurements of the recrystallization velocity in the case of a SPE and LPE (Liquid Phase Epitaxy) respectively, before presenting an analysis of the order parameter and melting temperatures to evidence SPE in the first category and LPE in the second.

A. Comparison with experiments

The first category of results, obtained with the Tersoff and SW115 potentials, is characterized by an Arrhenius behavior similarly to the numerous measurements of the velocity of recrystallization of a-Si on c-Si reviewed by Olson and Roth [16]. Fig. 10 presents the experimental results of Olson and Roth, our simulations with the Tersoff and SW115 potential and their fits. The activation energies are 2.99 eV for Tersoff and 1.87 eV for SW115. Not only the activation energy obtained with the Tersoff potential is close to the value of 2.68 eV obtained by measurement, but also the absolute value of the velocity. As for the SW115, the behavior is globally similar, although the recrystallization is accelerated in comparison to the experiments, especially at low temperature (lower activation energy).

The second category includes the Stillinger-Weber, the EDIP and the Lenosky potentials.

These three potentials have in common higher values of the velocities, up to 0.02 nm/ps i.e. 20 m/s. These results are orders of magnitudes far from the law of Olson and Roth [16],

they don't seem to have anything to do with the phenomenon measured by these authors.

However, velocities of m/s are typical of the recrystallization of silicon melted by pulse laser [START_REF] Galvin | [END_REF] and of explosive recrystallization of a-Si during pulsed laser irradiation [39]. In both cases, the speed of the recrystallization is governed by the transformation of liquid silicon (l-Si) into c-Si. In the case of explosive recrystallization of a-Si, the liquid is formed at the a-Si/c-Si interface due to a local heating by the release of heat generated by the transition from a-Si to c-Si [39]. The speed of crystallisation of l-Si has been measured by Galvin et al. [START_REF] Galvin | [END_REF]. These measurements are presented in Fig. 11 together with the results obtained with the Stillinger-Weber potential. The simulations at high temperature are in excellent agreement with these measurements. Moreover, the transition-state theory of crystal growth gives the following expression for the velocity of recrystallization [40]:

v(T ) = C exp - Q kT 1 -exp - L(T 0 -T ) kT ( 3
)
where Q is the activation energy for viscous or diffusive motion in the liquid, L is the heat of fusion of the phase transformation, T 0 is the equilibrium melting temperature and C is a constant. Our simulations were fitted to this law, and Fig. 11 shows that an excellent adjustment can be obtained with T 0 = 1600 K and Q = 0.54 eV. The results obtained with the EDIP and Lenosky potentials can also be perfectly calibrated with Eq. 3, but the agreement with experiments is less good in these two cases : EDIP gives a slower recrystallization, and Lenosky is shifted towards smaller temperatures.

The conclusion is therefore that the recrystallization velocities obtained with the Stillinger-Weber, EDIP and Lenosky potentials are very different of the results obtained for a solid phase epitaxy, but on the other hand they are typical from the solidification of a liquid. The question is now on the origin of this behavior, since the initial system is indeed an a-Si / c-Si stack, and not a l-Si / c-Si stack. This problem is studied in the following section.

B. Analysis : SPE versus LPE

In order to analyze in more details the structure of the initial a-Si layer upon annealing, we have calculated the order parameter defined by Albenze et al. [35]. First the parameter called A is calculated for each atom:

A = i cos θ i + 1 3 2 (4)
where the sum is performed on the 4 nearest neighbors of the atom and θ i is the angle made by each triplet. A to classify the atom as a solid-like (A < 0.4) or a liquid-like (A > 0.4) one. Consecutively, the structure under study is cut in slides, and the percentage of solidlike or liquid-like atoms is calculated in each slide. This order parameter is close to 1 in a crystal, around 0.2-0.4 in a liquid, and around 0.7-0.8 in an amorphous material. The order parameter as a function of z, the direction perpendicular to the disordered Si/c-Si, is represented for a simulation at 1273 K with the Stillinger-Weber potential (Fig. 12) and at 1673 K with the Tersoff potential (Fig. 13). Fig. 12 clearly demonstrates that the disordered silicon for z > 30 Å, while amorphous at the origin, transforms to a liquid during the initial equilibration stage. When the simulation of the annealing proceeds, the interface moves towards the right, as the recrystallization occurs, and the disordered part remains liquid. On the contrary with the Tersoff potential, the order parameter remains in the range [0.6-0.8] during all the process of recrystallization. This is a first proof that the a-Si initially constructed is transformed into l-Si when the Stillinger-Weber potential is used, while it remains amorphous with the Tersoff potential.

One could argue that these differences of behavior are due to the ability of the interatomic potentials to correctly predict the experimental melting temperature of a-Si. Indeed, the melting temperature of a-Si, T a m , is measured in the range [1420-1530] K [39,41]. However, the Tersoff potential is known for its large overestimation of the melting temperature, for c-Si as for a-Si : Marqués et al. [42] have determined a value of that the melting temperature of T a m = 2050 K, while for the Stillinger-Weber potential, a value of T a m = 1400 K was found by Brambilla et al. [43] and 1075 K by Albenze et al. [26]. Indeed, Fig. 14 shows that the temperature used with the Stillinger-Weber potential must be decreased to [473-673] K to keep an order parameter of the amorphous kind after the initial equilibration stage. Therefore we think that not only the melting temperature enter into account but more globally the ability of the potential to handle or not an amorphous structure. The Stillinger-Weber potential clearly favors the liquid state, as already observed in the literature [36].

In order to further corroborate the transformation to a liquid state, the radial and angular distribution of the atoms have been calculated in the disordered layer after the equilibration stage. The results are respectively presented in Fig. 15 and16. After a simulation with the Tersoff and SW115 potentials, the typical characteristics of an amorphous material are preserved: distinct peak of the second neighbours with a gap between first and second neighbours, and angular distribution comparable to the initial one (WWW a-Si). With Stillinger-Weber, EDIP and Lenosky on the other hand, there is no gap between the peaks of the first and second neighbors in the radial distribution (see Fig. 15), these second peak being even much less pronounced. Moreover, the angular distribution is drastically smeared out (see Fig. 16). All such structural features are far more typical of a liquid [44,45].

Consequently, our understanding of MD simulations of the recrystallization of a-Si using empirical potentials is the following. The Tersoff and SW115 potentials are adapted to the description of the SPE. The amorphous character is preserved with these potentials and the velocity of recrystallization follows an Arrhenius law as observed experimentally. The SW115 potential gives however a recrystallization velocity two orders of magnitude faster than the experiment. The Tersoff potential gives an excellent agreement with the law of Olson and Roth [16]. Its drawback is its excessive value of the melting temperature (for a-Si as for c-Si).

We have tried to apply a simple rescaling of the temperature T rescaled = T simu ×T c m,exp /T c m,simu , T simu × 1685/2400 in the case of the Tersoff potential. This approach is not satisfying since the simulation has moved away from the experiment by nearly 4 order of magnitudes as shown in Fig. 17.

The Stillinger-Weber, EDIP and Lenosky potentials are not adapted to the simulation of the SPE of a-Si. These three potentials seem to have a tendency to transform the amorphous structure into a liquid one, even below the nominal bulk melting temperature. Nevertheless, the resort to the Stillinger-Weber potential is perfectly justified to simulate the LPE, the experimental results obtained by pulse laser annealing [START_REF] Galvin | [END_REF] being in the margin of incertitude of our simulations. The sluggish behavior of the EDIP potential, already mentioned by environments. This environment dependency is a disadvantage for this particular problem, the return to the crystal arrangement is not enough enforced. Finally, the Lenosky potential is also a good candidate for the simulation of LPE if the temperature is scaled to account for the differences in metling temperatures (see Fig. 17).

V. CONCLUSION

The recrystallization of an amorphous silicon layer on a crystalline silicon substrate has been simulated by means of molecular dynamics. The initial a-Si on c-Si stack has been constructed by the WWW method and the structural properties of the resulting amorphous were carefully analyzed to validate the atomic arrangement in the structure with respect to the available experimental results. Subsequently, the annealing of this amorphous-crystal stack was studied using five of the most widely used interatomic potentials for silicon. A systematic campaign of molecular dynamics simulations as a function of temperature in the appropriate interval to observe amorphous recrystallization was carried out to extract the velocity of recrystallization versus the temperature. A wide disparity is obtained as a function of the interatomic potential, underscoring the limitations of the different empirical potentials that are not able to treat with an equivalent and consistent degree of accuracy the three Si phases of interest in the technological process (diamond crystal, amorphous and liquid). A similar shortcoming had been already pointed out by Ding and Andersen [47] for the case of germanium. However, the results can be classified in two categories. The first category is the expected solid phase epitaxy of the a-Si, obtained with the Tersoff and SW115 potentials. The second category corresponds rather to a liquid phase epitaxy as obtained with the potentials of Stillinger-Weber, EDIP or Lenosky, due to the transformation of the starting a-Si cluster into a liquid one. The reliability of the empirical potential molecular dynamics compared to experiments is shown to be very good, provided that a suitable potential is chosen to investigate each different phenomenology. Based on the results of the present study, the best candidates turn out to be Tersoff for the SPE and Stillinger-Weber for the LPE. potential is also represented for the same condtions (velocity scaling, 1273K). The equilibration step improves the quality of the interface since the fraction of fourfold silicon atoms is obviously enhanced. On the other side, the choice of the interatomic potential has a clear influence on the defects at the interface even for a very short annealing.

FiguresFIG. 1

 1 Figures

FIG. 2 :

 2 FIG.2: Angular distribution in the a-Si cluster.

FIG. 3 :FIG. 4 :FIG. 5 :FIG. 6 :FIG. 7 :FIG. 8 :FIG. 9 :FIG. 10 :

 345678910 FIG.3: Radial distribution in the a-Si cluster.

FIG. 11 :FIG. 12 :FIG. 13 :FIG. 14 :FIG. 16 :

 1112131416 FIG.11: Simulations of the recrystallization velocities with the Stillinger-Weber potential (+), measurement of the velocity of recrystallization of melted Si[START_REF] Galvin | [END_REF] (×) and fit of the simulations according to Eq. 3.

FIG. 17 :

 17 FIG. 17: Experimental and calculated velocities of recrystallization versus the temperature. The temperature has been rescaled for the simulations following T rescaled = T simu ×T c m,exp /T c m,simu , where the experimental melting temperature of c-Si T c m,exp is equal to 1685 K and the simulated melting temperature of c-Si T c m,simu is potential dependent.

TABLE I :

 I Distribution of the mean coordination number for the initial a-Si/c-Si interface. The evolution of this distribution during the equilibration step of4 ps with the various interatomic

	Tables
		Conditions		Coordination Number Twofold Threefold Fourfold Fivefold Sixfold
	Initial Interface		3.73	10.6%	16.9%	55.7% 16.8%	0%
	Duration Temperature Potential
	4ps	1273K	Tersoff	3.92	0.6%	10.7% 84.20% 4.50%	0%
	4ps	1273K	SW115	3.72	4.3%	18.50% 77.20%	0%	0%
	4ps	1273K	SW	3.86	0%	15.4%	79%	5.6%	0%
	4ps	1273K	EDIP	4.08	0%	3.0%	85.3% 11.7%	0%
	4ps	1273K	Lenosky	4.63	0%	6.4%	49.7% 23.1% 20.8%

VI. ACKNOWLEDGEMENTS

This work is part of the common laboratory STMicroelectronics/IEMN supported by the MINEFI. The numerous computer simulations would not have been possible without the reliable park of workstations managed by Jean-Michel Droulez. The authors thank W.

Windl for the suggestion to look at the Lenosky potential and Guy Allan for discussion about the amorphous cluster.