Romain Bachelard

Nicola Piovella

Philippe W Courteille

Cooperative scattering and radiation pressure force in dense atomic clouds

I. INTRODUCTION

The fact that ensembles of point-like scatterers respond collectively to incident radiation is well-known since the seminal paper of Dicke [1]. The collective phenomenon, termed superradiance, has been the topic of a huge amount of theoretical and experimental work. However, the question about the exact nature of the state generated by the radiation travelling through the ensemble has only been raised very recently. Scully and coworkers pointed out [2], that the dipole moments of scatterers distributed along the incident beam's optical axis are excited in phase with beam's propagating phase front. The resulting 'timed Dicke' state emits light predominantly into forward direction, provided the ensemble size is large compared to the radiation wavelength [3].

This simple picture, however, only holds at low optical densities. At high optical densities and spherically symmetric ensembles, the refraction index of the scattering medium delays the propagation of the pump beam and distorts its phase front. This distortion can have a huge impact on any macroscopic observable of the system, such as the angular distribution of the scattered radiation, the phase-front of the transmitted beam or the force acting on the center-of-mass of the ensemble.

The aim of this paper is to calculate the correct state generated in an ensemble of two-level systems, e.g. an atomic cloud, by interaction with a weak laser beam and how this state cooperatively scatters the incident light. Cooperative scattering by many atoms has been studied extensively in the past, both classically and quantummechanically. Classically, scattering at extended objects is described as Mie scattering, showing resonances induced by the boundary conditions that the target imposes to the incident light field [START_REF] Van De Hulst | Light Scattering by Small Particles[END_REF]. Quantum-mechanically, Dicke [1] has shown that, when two-level atoms are confined inside a volume much smaller than a radiation wavelength, the emission can be superradiant or subradiant. How the classical and quantum pictures are linked has been demonstrated at the example of a sample of weakly excited atoms [START_REF] Prasad | [END_REF]6]. More precisely, when a single atom out of N is excited, the Dicke symmetric state of maximum cooperation radiates superradiantly, i.e. at a decay rate proportional to N . Cooperative effects related to the superradiant and directional emission by an extended ensemble of atoms in a timed symmetric Dicke state have been observed in the radiation pressure force acting on a large cloud of atoms driven by a resonant radiation field. Depending on the detuning of the incident radiation frequency from atomic resonance, the radiation pressure force may be either drastically reduced due to both increased forward scattering and a reduced scattering cross section [3,7], or even enhanced if the cooperative Mie scattering dominates over superradiance [8].

In this paper we revisit the scattering by N atoms driven by a constant uniform radiation field and emitting radiation into free space. Our description assumes several approximations: 1) weak excitation of the atomic ensemble (one atom out of N is excited); 2) Markov or 'rapid transit' approximation [START_REF] Prasad | [END_REF] (photon time of flight through the cloud much shorter than atomic decay time); 3) atoms frozen (zero temperature) and motionless; 4) neglected dipole-dipole interactions and collisions; 5) no resonant interactions, i.e. no resonance fluorescence, Van der Waals interactions etc. Approximation 3) excludes several cooperative effects related to atomic recoil motion (e.g. collective atom recoil lasing [9,10] or matter wave superradiance [11,12]), and thus neglects stimulated scattering processes along preferential directions, as for instance end-fire modes in Bose-Einstein condensates [13] or optical cavity modes [14]. Neglecting atomic interactions (approximation 4) is justified assuming atomic distances much larger than an optical wavelength. The yet interesting opposite regime (atoms closer than an optical wavelength) requires the solution of the atomic equations with an exponential interaction kernel (see Eq.( 3)) and will be discussed in a future publication.

We here determine the stationary solution for a spherical Gaussian distribution, as well as relevant macroscopic quantities such as medium polarization, scattered radiation intensity and radiation pressure force. Our solution is based on a solution of the eigenvalue problem in the smooth density approximation. For large size samples, a continuous spectrum limit allows to obtain explicit analytical expressions for such quantities, expressed as a function of experimentally controllable parameters such as frequency and power of the driving field, optical thickness and size of the atomic sample.

The paper is organized as follow: in Sec.II the scattering problem is expressed in terms of multi-particle coupled differential equations. In Sec.III the continuous density approximation is introduced and the solution for the atomic field is found in terms of the discrete eigenvalues of the interaction operator. Then, some macroscopic quantities of importance are calculated in Sec.IV. A continuous spectrum approximation for large atomic clouds performed in Sec.V appears to be particularly suitable to evaluate macroscopic quantities for arbitrarily large value of the atomic density. Sec. VI revises the radiation pressure force obtained assuming a symmetric timed Dicke atomic state [3] and, comparing this state to the exact solution, we find that phase shifts induced by the atomic cloud's refractive index are at the origin of important corrections for the expectation values of macroscopic quantities. Finally, numerical results and conclusions are presented in Sec.VII.

II. COOPERATIVE SCATTERING PROBLEM

A system of two-level (g and e) atoms with resonant frequency ω a and position r j , driven by an uniform laser beam with electric field amplitude E 0 , frequency ω 0 and wave vector k 0 = (ω 0 /c)ê z , is described by the interaction Hamiltonian:

H = Ω 0 2 N j=1 σj e i(∆0t-k0•rj ) + h.c. + N j=1 k g k σj e -iωat + σ † j e iωat â † k e i(ω k t-k•rj ) + âk e -i(ω k t-k•rj ) (1)
where Ω 0 = dE 0 / is the Rabi frequency of the incident laser field, σj is the lowering operator for atom j, âk is the photon annihilation operator and

g k = (d 2 ω k /2 ǫ 0 V ph ) 1/2
is the single-photon Rabi frequency, where d is the electricdipole transition matrix element and V ph is the photon volume. A special case is when a single photon is present in the mode k, as was extensively investigated in Refs. [2,15,16]. The total system (atoms+photons) is assumed to be in a state of the form [6]:

|Ψ = α(t)
where ∆ 0 = ω 0 -ω a . The first term corresponds to the initial ground state without photons, the sum in the second term are the states where a single atom has been excited by the classical field. The third term corresponds to the atom returned to the ground state having emitted a photon in the mode k, whereas the last one corresponds to the presence of two excited atoms and one virtual photon with 'negative' energy. It is due to the counter-rotating terms in the Hamiltonian (1) and disappears when the Rotating Wave Approximation (RWA) is made. In the linear regime (i.e. α ≈ 1) and in the Markov approximation (valid if the decay time is larger than the photon time-of-flight through the atomic cloud), the problem reduces to the following differential equation [START_REF] Scully | [END_REF],

βj = i∆ 0 - Γ 2 β j -i Ω 0 2 e ik0•rj + i Γ 2 m =j exp(ik 0 |r j -r m |) k 0 |r j -r m | β m , (3)
where Γ = V ph g 2 k k 2 0 /πc is the single-atom spontaneous decay rate. The kernel in the last term of Eq.(3) has a real component, -(Γ/2) m =j [sin(ρ jm )/ρ jm ] (where ρ jm = k 0 |r jr m |), describing the collective atomic decay, and an imaginary component, i(Γ/2) m =j [cos(ρ jm )/ρ jm ], describing the collective Lamb shift due to short range interaction between atoms, induced by the electromagnetic field [START_REF] Scully | [END_REF][20][21]. The latter becomes significant when the number of atoms in a cubic optical wavelength, nλ 3 , is larger than unity, in which case the contribution from the virtual processes described by the counter-rotating terms in the Hamiltonian becomes relevant. Hence, for sufficiently dilute system, such that N ≪ σ 3 where σ = k 0 σ R and σ R is the cloud size, the collective phase shift arising from the imaginary part of the kernel in Eq.( 3) can be disregarded [22] and the scattering problem reduces to

βj = i∆ 0 - Γ 2 β j -i Ω 0 2 e ik0•rj - Γ 2 m =j sin(k 0 |r j -r m |) k 0 |r j -r m | β m . (4)
with initial condition β j (0) = 0, for j = 1, . . . , N . Notice that Eq.(3) deduced in the quantum mechanical description may be also obtained classically when the two-level atoms are treated as weakly excited classical harmonic oscillators [START_REF] Prasad | [END_REF]6]. For this reason the solution of Eq.(3) (or of the approximated version (4)) has a wider interest for the general problem of collective radiation scattering.

As for the radiation field operator âk , it evolves according the Heisenberg equation

dâ k dt = 1 i [â k , Ĥ] = -ig k e i(ω k -ωa)t N m=1 σm e -ik•rm , (5)
where the fast oscillating term proportional to exp[i(ω k + ω a )t] has been neglected.

III. CONTINUOUS DENSITY APPROXIMATION

In light scattering experiments, disorder (or granularity) plays a role when the number of atoms projected onto a cross section perpendicular to the incident beam is small enough so that a light mode focused down to the diffraction limit (that is ∼ λ 2 ) would be able to resolve and count the atoms. In other words, the stochastic fluctuations induced by the random positions of the atoms can be neglected when the total number of atoms N is larger than the number of modes ∼ σ 2 that fit into the cloud's cross section, i.e. when the optical density is b 0 = 3N/σ 2 ≫ 1. Under this hypothesis, the particles can be described by a smooth density n(r), and their probability to be excited by a field β(r, t), so that Eq.( 4) turns into

∂ ∂t β(r, t) = i∆ 0 - Γ 2 β(r, t) -i Ω 0 2 e ik0•r - Γ 2 dr ′ n(r ′ ) sin(k 0 |r -r ′ |) k 0 |r -r ′ | β(r ′ , t). (6)
In what follows we will consider only spherically symmetric distributions n(r). Because of the linearity of Eq.( 6), it is convenient to introduce an eigenbasis of the coupling operator. The functions j n (r)Y nm (θ, φ), with j n s the spherical Bessel functions and Y nm (θ, φ)s the spherical harmonics, appear as a natural choice considering the following identity

sin(k 0 |r -r ′ |) k 0 |r -r ′ | = 4π ∞ n=0 n m=-n j n (k 0 r)Y * nm (θ, φ)Y nm (θ ′ , φ ′ )j n (k 0 r ′ ). (7)
In particular, the choice of the spherical harmonics guarantees the orthogonality of the basis, since

2π 0 dφ π 0 dθ sin θY * nm (θ, φ)Y n ′ m ′ (θ, φ) = δ nn ′ δ mm ′ . (8)
Therefore, assuming the following decomposition for the field,

β(r, t) = ∞ n=0 n m=-n α nm (t)j n (k 0 r)Y nm (θ, φ), (9)
the projection of Eq.( 6) along the eigenmodes leads to

αnm -i∆ 0 α nm + Γ 2 (1 + λ n )α n,m j n (k 0 r) = - i 2 Ω nm , (10)
where λ n are the eigenvalues associated to modes n

λ n = 4π ∞ 0 r 2 n(r)j 2 n (k 0 r)dr (11)
whereas Ω nm corresponds to the projection of the incident wave on mode (n, m)

Ω nm = Ω 0 2π 0 dφ π 0 dθ sin θ Y * nm (θ, φ) e ik0r cos θ = 2Ω 0 δ m0 π(2n + 1)i n j n (k 0 r). (12)
Assuming the cloud is initially unexcited, i.e. α mn (0) = 0, only spherically symmetric components with m = 0 will grow so that, defining α n (t) ≡ α n0 (t), Eq.( 10) reduces to

αn -i∆ 0 - Γ 2 (1 + λ n ) α n = -i n+1 π(2n + 1)Ω 0 (13)
Eq.( 13) straightforwardly integrates and, inserted in Eq.( 9), leads to the following expression for the excitation field

β(r, θ, t) = Ω 0 Γ ∞ n=0 i n (2n + 1)j n (k 0 r)P n (cos θ) 2δ + i(1 + λ n ) 1 -e i∆0t e -(Γ/2)(1+λn)t , (14)
where the scaled detuning δ = ∆ 0 /Γ was introduced. Hence, each mode n relaxes toward the steady-state with a characteristic time τ n = 1/Γ(1 + λ n ): the first modes relax very quickly since λ n is proportional to N , yet for the highest modes, τ n ∼ Γ -1 , even if their macroscopic contribution is usually small. Eventually, for times much longer than the single-atom decay time Γ -1 , the field tends toward a stationary state fully characterized by the spectrum

β s (r, θ) = Ω 0 Γ ∞ n=0 i n (2n + 1) 2δ + i(1 + λ n ) j n (k 0 r)P n (cos θ). (15)
Notice that the set of eigenvalues ( 11) is complete since from the identity n≥0 (2n + 1)

j 2 n (z) = 1, it follows that ∞ n=0 (2n + 1)λ n = 4π ∞ 0 r 2 n(r)dr = N, (16)
which corresponds to the trace of the coupling operator.

IV. MACROSCOPIC QUANTITIES

The description of the field β s (r, θ) in terms of spectrum also provides expressions for any macroscopic quantities, the most relevant of which are calculated here below. These formulae will be specialized to Gaussian clouds in the subsequent section.

A. Average amplitude and probability of excitation

The average 'phased' probability of the timed Dicke state [START_REF] Scully | [END_REF] and the excitation probability are respectively

β s e -ik0•r = 2π N π 0 dθ sin θ ∞ 0 drr 2 n(r)β s (r, θ)e -ik0r cos θ ( 17
)
|β s | 2 = 2π N π 0 dθ sin θ ∞ 0 drr 2 n(r)|β s (r, θ)| 2 . ( 18
)
Inserting Eq.( 15) and using the identities

1 -1 dxP n (x)e iαx = 2i n j n (α) , 1 -1 dx P m (x)P n (x) = 2 2n + 1 δ mn , (19)
we obtain

β s e -ik0•r = Ω 0 ΓN ∞ n=0 (2n + 1)λ n 2δ + i(1 + λ n ) ( 20
)
|β s | 2 = Ω 2 0 Γ 2 N ∞ n=0 (2n + 1)λ n 4δ 2 + (1 + λ n ) 2 . ( 21
)
B. Scattered field

The electric field radiated by the excited atoms reads, in the smooth density limit (see Appendix A),

E S (r, θ, φ, t) = - dk 2 0 4πǫ 0 r e ik0(r-ct) dr ′ n(r ′ )β(r ′ , t)e -iks•r ′ . ( 22
)
Using the stationary solution (15) for a spherically symmetric distribution and using the integral

π 0 dθ ′ sin θ ′ P n (cos θ ′ )J 0 (k 0 r ′ sin θ sin θ ′ )e -ik0r ′ cos θ cos θ ′ = 2i -n j n (k 0 r ′ )P n (cos θ) (23)
we obtain

E S (r, θ) = - E 0 4πk 0 r e ik0(r-ct) ∞ n=0 (2n + 1)λ n 2δ + i(1 + λ n ) P n (cos θ), (24)
where we used the relation Γ = d 2 k 3 0 /(2π ǫ 0 ). In the forward direction (θ = 0) the scattered field is proportional to the 'phased' probability amplitude (20). The expression (24) provides the angular distribution of the scattered radiation field. We can also calculate the scattered intensity as

I S (r, θ, φ) = cǫ 0 Ê † S ÊS = cǫ 0 dk 2 0 4πǫ 0 r 2   N j=1 |β j | 2 + j =m β * m β j e -iks•(rj -rm)   (25)
Passing to the continuous limit and using Eqs.( 21), ( 22) and (24), we obtain:

I S (r, θ) = I 0 (4πk 0 r) 2   ∞ n=0 (2n + 1)λ n 4δ 2 + (1 + λ n ) 2 + ∞ n=0 (2n + 1)λ n 2δ + i(1 + λ n ) P n (cos θ) 2   (26)
where

I 0 = cǫ 0 E 2 0 .
The scattered intensity is the sum of the incoherent contribution, proportional to N (since λ n ∝ N , see Eq.( 11)) and isotropic, and the superradiant contribution, proportional to N 2 and directed (for extended clouds) mainly in the forward direction. Integrating over the solid angle, the total scattered power is

P S = 2πr 2 π 0 dθ sin θI S (r, θ) = I 0 4πk 2 0 ∞ n=0 (2n + 1)λ n (1 + λ n ) 4δ 2 + (1 + λ n ) 2 (27)
where we used the second of the identities [START_REF] Scully | [END_REF].

C. Radiation pressure force

The radiation force operator acting on the jth atom is calculated from Eq.(1) as Fj = -∇ rj Ĥ = Faj + Fej where [3]

Faj = i k 0 Ω 0 2 σj e i(∆0t-k0•rj ) -h.c. ( 28
) Fej = i k kg k â † k σj e i(ω k -ωa)t-ik•rj -σ † j âk e -i(ω k -ωa)t+ik•rj (29)
are the forces due respectively to the absorption and emission processes. In Eq.(29) we have neglected the counter-rotating terms proportional to exp[±i(ω k + ω a )t]. We are here interested in the average radiation force F = (1/N ) j Fj , which stands for the force acting on the center-of-mass of the atomic cloud along the direction of the incident field, k 0 = k 0 êz . This average force is relatively easy to measure by time-of flight techniques in cold atomic clouds released, for instance, from magneto-optical traps (MOTs) and has recently revealed cooperative effects in the scattering by extended atomic samples [7,8]. It may provide a convenient measurement (aside from the scattered radiation) of the effects that cooperative scattering imprints on the atoms. The average absorption force along the z-axis, resulting from the recoil received upon absorption of a photon from the incident laser, is

Fa = = i 2N k 0 Ω 0 N j=1 σj e i∆0t-ik0•rj -h.c. . (30)
The second contribution Fe = (1/N ) j Fej results from the emission of a photon into any direction k. Inserting âk from Eq.( 5) into Eq.( 29) and approximating the sum over the modes k by an integral, it is possible to obtain, in a way similar as done for the radiation field operator ÊS of Eq.(A5), the following expression for the average emission force along the z-axis [3]:

Fe = - k 0 Γ 8πN 2π 0 dφ π 0 dθ sin θ cos θ N j,m=1 e -ik•(rj-rm) σ † m σj + h.c. . (31)
Evaluating their expectation values on the state (2) (neglecting virtual photon contributions), the emission force in the discrete model is

F e = - k 0 Γ 8πN 2π 0 dφ π 0 dθ sin θ cos θ N j,m=1 β j β * m e -ik•(rj -rm) + c.c. = i k 0 Γ 2N N j,m=1 (z j -z m ) |r j -r m | j 1 (k 0 |r j -r m |)(β j β * m -c.c.). ( 32
)
where we used the identity (B1) and j 1 (z) is the first order spherical Bessel function. Then, passing to the continuous distribution limit, Eqs.( 30) and (31) are approximated by

Fa = -k 0 Ω 0 Im βe -ik0•r (33) Fe = -k 0 Γ 8πN 2π 0 dφ π 0 dθ sin θ cos θ drn(r) dr ′ n(r ′ ) β(r)β * (r ′ )e -ik•(r-r ′ ) + c.c. . (34)
The absorption stationary force is readily obtained from Eq.( 20),

Fa = k 0 Ω 2 0 ΓN ∞ n=0 (2n + 1)λ n (1 + λ n ) 4δ 2 + (1 + λ n ) 2 (35)
whereas a longer calculation, reported in Appendix B, yields the emission stationary force,

F e = -k 0 2Ω 2 0 ΓN ∞ n=0 (n + 1)λ n λ n+1 [4δ 2 + (1 + λ n )(1 + λ n+1 )] [4δ 2 + (1 + λ n ) 2 ][4δ 2 + (1 + λ n+1 ) 2 ] . (36)
As expected, the absorption force pushes the atomic cloud in the direction of the driving field, whereas the emission force is oppositely directed and is proportional to N . Both forces depend on N also through the collective decay rate Γ(1 + λ n ). We observe that the absorption force (35) is linked to the scattered power (27) by the relation

N Fa = 4π P S c , (37)
i.e. the absorption force is proportional to the scattered power per atom, P S /N . Fig. 1 compares the (a) absorption, (b) emission and (c) total forces vs. N for δ = 10 and a Gaussian density profile with σ = 8 (where σ = k 0 σ R ), calculated using the series (35) and (36) (circles) and N -body simulations (see Eq.4) (dots). The eigenvalues λ n for the Gaussian density profile are given by Eq.(38). The forces are reported as a ratio between the cooperative force and the single-atom force, F 1 = k 0 ΓΩ 2 0 /(4∆ 2 0 + Γ 2 ). The error bars correspond to the standard deviation of the observables over 8 realizations. Fig. 1d shows the total force vs. N for a different choice of parameters, σ = 5 and δ = 200, for which the force exhibits a maximum as a function of N [8]. We observe a good agreement between the analytical solution (35) and the N -body simulations at large N , when effects due to the discreteness of the system become negligible. 35), (36) and their sum, the circles to the analytical expressions (46), (36) and their sum, the dots to N -body simulations (see Eq.4) and the dash-dotted lines to the STD state. The error bars correspond to the standard deviation of the observables over 8 realizations.

We recall that our expressions have been obtained in the continuous density approximation, i.e. assuming a sample with large optical thickness, i.e. N ≫ σ 2 , but sufficiently dilute to neglect the collective Lamb shift, i.e. with a small number of atoms in a cubic optical wavelength volume, N ≪ σ 3 . These two conditions imply that the cloud size should be much larger than the optical wavelength, σ ≫ 1. Nevertheless, our results remain valid also in the limit n(r) → 0 where, in particular, the standard radiation pressure force, F 1 , is recovered. However, at low densities, fluctuations will strongly affect individual measurements of the radiation pressure force and in this regime our expressions represent only the average expectation values.

V. LARGE GAUSSIAN CLOUDS

Large clouds (σ ≫ 1) behaves fundamentally differently from small clouds (σ 1), as can be deduced from their spectrum. Let us for example consider the case of Gaussian clouds, with density n(r) = (N/(2π) 3/2 σ 3 R ) exp(-r 2 /2σ 2 R ). The spectrum then reads

λ n = N 2 π ∞ 0 ρ 2 e -ρ 2 /2 j 2 n (σρ)dρ = N π 2σ 2 e -σ 2 I n+1/2 (σ 2 ), (38)
where σ = k 0 σ R is the scaled size of the cloud, and I n (x) the n-th modified Bessel function. Recently the spectrum λ n for an uniform spherical cloud has been calculated by A.A. Svidzinsky et al. [6,16], also for the exponential kernel of Eq.( 3). However, a Gaussian distribution is certainly more realistic for experiments with cold dense atomic ensembles. Generally, the spectrum of small clouds (σ ≤ 1) is composed of a few significant eigenmodes, whereas for σ large, all the eigenmodes for n < σ are significant and the spectrum can be treated as a continuum. In particular, in the latter case the λ n can be approximated, for n < σ, by λ n ∼ (N/2σ 2 ) exp[-(n + 1/2) 2 /2σ 2 ] (see e.g. [START_REF] Abramowitz | Handbook of Mathematical Functions[END_REF]). Switching to a continuous treatment of the spectrum, we define η = n + 1/2 and get

λ η = N 2σ 2 e -η 2 /(2σ 2 ) . ( 39
) ∞ n=0 (2n + 1) → 2 ∞ 0 ηdη. ( 40
)
Remark that using these definitions, the completeness condition ( 16) is still preserved. The continuous spectrum limit allows for the evaluation of the sums in Eqs.( 20) and ( 21) as continuous integrals

β s e -ik0•r = 2Ω 0 ΓN ∞ 0 λ η ηdη 2δ + i(1 + λ η ) = Ω 0 Γ 6 b 0 b0/6 0 dx 2δ + i(1 + x) ( 41
)
|β s | 2 = 2Ω 2 0 Γ 2 N ∞ 0 λ η ηdη 4δ 2 + (1 + λ η ) 2 = Ω 0 Γ 2 6 b 0 b0/6 0 dx 4δ 2 + (1 + x) 2 (42)
where we have set x = (b 0 /6) exp(-η 2 /2σ 2 ), with b 0 = 3N/σ 2 the optical thickness. The above expressions integrate as

β s e -ik0•r = Ω 0 Γ 3 b 0 2 arctan b 0 3 δ 1 + 4δ 2 + b 0 /6 -i ln 1 + b 0 3 1 + b 0 /12 1 + 4δ 2 ( 43
)
|β s | 2 = Ω 0 Γ 2 3 δb 0 arctan δb 0 /3 1 + 4δ 2 + b 0 /6 . ( 44
)
These formulae highlight the prominent role of the parameters b 0 and δ in the high-density limit. In a similar way, we calculate the total scattered power as:

P S = I 0 N 4πk 2 0 1 0 dx 1 + (b 0 /6)x 4δ 2 + [1 + (b 0 /6)x] 2 = I 0 σ 2 R 4π ln 1 + b 0 3 1 + b 0 /12 1 + 4δ 2 (45)
As expected, for small optical thickness (b 0 ≪ 1) the scattered power is incoherent, P S ≈ [I 0 N/(4πk 2 0 )]/(4δ 2 + 1). However, for large optical thickness it shows a logarithmic dependence on N . The superradiant character of the radiation is visible only observing the scattered intensity in the forward direction (see the second term of Eq.( 26)), but not in the total scattered power.

The absorption force is deduced from Eqs.( 33) and (43) as

Fa = k 0 3Ω 2 0 b 0 Γ ln 1 + b 0 3 1 + b 0 /12 1 + 4δ 2 . ( 46
)
The emission force can be written, in the continuous spectrum approximation in the integral form

Fe = -k 0 Ω 2 0 b 0 6Γ e -1/(4σ 2 ) ∞ 1/σ dy ye -y 2 4δ 2 + 1 + (b 0 /6)e -(y+1/2σ) 2 /2 1 + (b 0 /6)e -(y-1/2σ) 2 /2 4δ 2 + 1 + (b 0 /6)e -(y+1/2σ) 2 /2 2 4δ 2 + 1 + (b 0 /6)e -(y-1/2σ) 2 /2 2 , ( 47
)
where we have set y = (η + 1)/σ. As can be observed in Fig. 1, the continuous-spectrum approximation gives excellent results compared to the full series ( 35) and (36). In the limit σ → ∞ and finite b 0 , Eq.( 47) would lead to

Fe ≈ -Fa + k 0 Γ |β s | 2 , ( 48
)
which has a transparent interpretation: for a very large cloud, the atoms scatter radiation to forward direction and the recoil received by the atoms upon emission cancels out with the recoil received by absorbing a photon from the driving field. The net force remaining after the substraction is the non collective contribution to the emission. The net force is equal to the photon momentum k 0 times the emission rate Γ |β s | 2 . This emission rate depends indirectly on N and σ through the enhanced superradiant decay, Γ(N/4σ 2 ), which decreases the emission rate when the optical thickness increases. A more accurate expression of the emission force valid for large but finite cloud size would require the exact evaluation of the integral in Eq.(47). 44), the circles to series expressions ( 20) and ( 21), the dots to N -body simulations (see Eq.4) and the dash-dotted lines to the STD state. The error bars correspond to the standard deviation of the observables over 8 realizations. Note that for the STD state, σ β = 0; for too small N , the approximations (43) and (44) provide inconsistent results, that is negative σ β . Simulations realized for σ = 10 and δ = 10.

VI. SYMMETRIC TIMED DICKE STATE

A particular ansatz used by Scully and coworkers [2,6,16] is the symmetric timed Dicke (STD) state, given by

β(r, t) = β T D (t)e ik0•r . (49)
After integration over space of Eq.( 6), one obtains the following evolution equation

dβ T D (t) dt = i∆ 0 - Γ 2 (1 + N s ∞ ) β T D (t) -i Ω 0 2 , (50)
where s ∞ is the integrated structure factor of the cloud defined as

s ∞ = 1 N 2 n(r)dr dr ′ n(r ′ ) sin(k 0 |r -r ′ |) k 0 |r -r ′ | e -ik0•(r-r ′ ) = 1 4π 2π 0 dφ π 0 dθ sin θ| e i(k0-k)•r | 2 . ( 51
)
This ansatz is of particular interest since it allows to evidence e.g. the superradiant nature of the decay of such state when the pump is turned off [START_REF] Scully | [END_REF]. As for its steady-state, it reads [3,7]

β T D = Ω 0 Γ 1 2δ + i(1 + N s ∞ ) , (52)
and for a large cloud with Gaussian distribution, since s ∞ ≈ 1/4σ 2 , we get

β T D = Ω 0 Γ 1 2δ + i(1 + b 0 /12) . (53)
Thus, the STD solution (53) approximates the exact result (44) only for b 0 ≪ 12,

|β T D | 2 ≈ Ω 0 Γ 2 1 4δ 2 + (1 + b 0 /6) . (54)
Fig. 2 shows the average excitation probability |β s | 2 (left) and its variance,

σ 2 β = |β s | 2 -| β s e -ik0
.r | 2 vs. N for σ = 10 and δ = 10. The plain curves refer to the analytical expressions (43) and (44), the circles to series expressions (20) and (21), the dots to N -body simulations (see Eq.4) and the dash-dotted lines to the STD state (49). The error bars correspond to the standard deviation of the observables over different realizations. Note that for the STD state, σ β = 0; for too small N , the approximations (43) and (44) provide inconsistent results, that is negative σ β . We observe excellent agreement between the series and the analytical solutions, and a consistent reduction of the excitation probability decrease vs N with respect to the TDS prediction (dashed line in fig. 2,left). Also the fluctuations obtained from the N -body simulations converge for large values of N toward the results obtained in the continuous density approximation (fig. 2, left), showing the presence of a shot-noise contribution for small N .

The radiation force for the STD ansatz (49) is [3,7,8]

F = Fa + Fe = k 0 -Ω 0 Im(β T D ) -Γ|β T D | 2 N f ∞ (55)
where

f ∞ = 1 4π 2π 0 dφ π 0 dθ sin θ cos θ| e i(k0-k)•r | 2 . ( 56
)
Since for a spherically Gaussian distribution N (s ∞ -f ∞ ) ≈ N/(8σ 4 ) = b 0 /(24σ 2 ), from Eq.( 52) the stationary radiation force is [3,7,8]

F = k 0 Γ Ω 0 Γ 2 1 + N (s ∞ -f ∞ ) 4δ 2 + (1 + N s ∞ ) 2 (57)
≈ k 0 Γ Ω 0 Γ 2 1 + b 0 /24σ 2 4δ 2 + (1 + b 0 /12) 2 . ( 58
)
As can be observed in Fig. 1, the STD state yields a good agreement with the full-spectrum approach and the N -body simulations only for small values of the optical thickness. Finally, for that state (49) the scattered radiation electric field, Eq.(A5), and intensity, (25), become

E S (r, θ, φ, t) = dk 2 0 4πǫ 0 r e ik0(r-ct) β T D (t) e i(k0-k)•r (59)
and

I S (r, θ, φ) = cǫ 0 dk 2 0 4πǫ 0 r 2 |β T D | 2 N + N 2 | e i(k0-k)•r | 2 . (60)
In particular, for a Gaussian distribution, exp

[i(k 0 -k) • r] = exp[-σ 2 (1 -cos θ)/2].
Fig. 3 shows the phase of the excitation amplitude β s (left) and the excitation probability, |β s | 2 , (center) in the (x, z) plane (y = 0), calculated from the exact solution (15), for N = 10 4 , σ = 10 and δ = 10. Fig. 3 (right) shows the contribution to the electric field radiation in the same plane, weighted by the local atomic density. The simplification of the STD, as compared to the exact solution (15), resides in the assumption that all atoms are equally excited and oscillate in phase. According to the exact calculation, the atomic dipoles appear to be in phase only in the core of the cloud (see Fig. 3, left), but this phase profile has strong distortion away from it. This phenomenon is all the more important as the atoms are much more excited in the peripheral region than in the core (see Fig. 3, center). In particular, even when this excitation probability is weighted by the particle density, two areas at the cloud entrance and exit contribute significantly to the radiation electric field (24) (see peaked structures in Fig. 3, right). For a timed Dicke state, both the phase profile and average excitation remain flat throughout the cloud. From a macroscopic point of view, the STD state neglects phase shifts imprinted into the pump beam by the cloud's reflective index. This can be seen in Fig. 4(a), which compares the phase of the STD state, β T D exp(ik 0 • r), (linear curve, no phase shift) and β s (r) (additional phase-shift) along the optical axis across the cloud. Fig. 4(b) shows the pump beam phase shift after transmission through the atomic cloud as a function of atom number. This phase shift is at the origin of the deviation between the radiation pressure forces calculated for the STD state and the exact solution. The pump beam phase shift leads to a reduction of the absorption and the emission forces. This can be understood as destructive interference of forward radiation emitted from different atoms, located at the same plane z = z 0 but different x or y.

As for the emitted wave, it is concentrated in the forward direction (see Fig. 5, left), and there is no backscattering (not shown here). The wavefront phase does not exhibit significant distortion in the central region of the radiated beam (see Fig. 5, right).

VII. CONCLUSIONS

In summary, we studied collective scattering from a dense and large atomic cloud with gaussian density profile in terms of the eigenvalues of the interaction operator. This enabled us to calculate the state generated by the interaction with a laser beam. We found that this state considerably deviates from the 'timed Dicke' state at high optical densities. In order to characterize this state, we calculated the phase front of the excited atomic dipole moments and the force due to cooperative scattering by the atomic cloud.

Our approach consisted in expanding the solution on spherical harmonics, which form a complete orthogonal basis of the Hermitian interaction matrix. Under the assumption that the system is sufficiently dilute to neglect short-range effects due to dipole-dipole interactions, the continuous spectrum approximation allows to transform infinite series into solvable integrals, and eventually to derive analytical expressions for the most relevant observables, such as the scattered intensity and radiation pressure force. These analytical expressions show good agreement with the numerical solution of the N -body problem, and highlight the dependence on the accessible experimental parameters such as optical thickness, atomic cloud size and laser frequency.

The analytical solution appears particularly useful for studying the thermodynamic limit when N → ∞ and V → ∞ with N/V fixed, until collisions or non resonant interactions come into play. The thermodynamic limit is hardly accessible to N -body simulations, since the latter are highly CPU consuming.

In contrast, the eigenvalue approach opens the possibility to study the fascinating link between microscopic and macroscopic domains of light scattering, and in particular between single-atom scattering and Mie scattering for extended continuous samples [START_REF] Prasad | [END_REF].

For large optical thickness, the refraction index of the cloud acts back on the driving field and shifts its phase. For this reason, the solution for large optical thickness shows appreciable deviations from that obtained assuming a symmetric timed Dicke state for the atomic sample, since the latter corresponds to a completely degenerate eigenvalue spectrum. The exact solution (15) takes into account the induced phase shift and reduces to the symmetric timed Dicke state in the limit of relatively small optical thickness.

An important further development of the present study should be to understand how the dipole-dipole interactions contribute to the observed cooperative effects, completing in this way the cooperative scattering description for highly compressed and dense atomic clouds. In the stationary case,

α n = Ω 0 Γ i n 4π(2n + 1) 2δ + i(1 + λ n )
and the stationary force is ) [1] R.H. Dicke, Phys. Rev. 93, 99 (1954).

F e = -k 0 2Ω 2 0 ΓN ∞ n=0 (n + 1)λ n λ n+1 [4δ 2 + (1 + λ n )(1 + λ n+1 )] [4δ 2 + (1 + λ n ) 2 ][4δ 2 + (1 + λ n+1 ) 2 ] . ( B10

  FIG. 1: (a) Absorption, (b) emission and (c) total forces vs N for δ = 10 and a Gaussian cloud with σ = 8. (d) Total force for σ = 5 and δ = 200. Forces are relative to the single-atom force F1. The plain curves refer to the series (35), (36) and their sum, the circles to the analytical expressions (46), (36) and their sum, the dots to N -body simulations (see Eq.4) and the dash-dotted lines to the STD state. The error bars correspond to the standard deviation of the observables over 8 realizations.

2 FIG. 2 :

 22 FIG. 2: Average excitation probability |βs| 2 (left) and variance σ 2 β = |βs| 2 -| βse -ik 0 .r | 2 (right) vs N . The plain curves refer to the analytical expressions (43) and (44), the circles to series expressions(20) and(21), the dots to N -body simulations (see Eq.4) and the dash-dotted lines to the STD state. The error bars correspond to the standard deviation of the observables over 8 realizations. Note that for the STD state, σ β = 0; for too small N , the approximations (43) and (44) provide inconsistent results, that is negative σ β . Simulations realized for σ = 10 and δ = 10.

FIG. 3 :FIG. 4 :

 34 FIG. 3: Phase of the atomic radiation ang(βs) (left), Excitation probability |βs| 2 of the atoms in the (x, z) plane (y = 0) (center) and contribution to the radiation of the electric field in the same plane, which corresponds to the level of excitation of the atoms weighted by their local density (right). Simulations realized using the analytic expressions (15), for N = 10 4 , σ = 10 and δ = 10.

FIG. 5 :

 5 FIG. 5: Transverse intensity profile (left) of the wave emitted by the cloud, at z = 3σ and its wavefront (right). Simulations realized using the analytic expression (24), for N = 10 4 , σ = 10 and δ = 10.

dφ π 0 ××dθ 0 dθ sin 3

 003 dθ sin θ cos θe -ik•(r-r ′ ) = 4πi z -z ′ |rr ′ | j 1 (k 0 |rr ′ |). k 0 |rr ′ |) = -k 0 (z -z ′ ) |rr ′ | j 1 (k 0 |rr ′ |),where j 0 (x) = sin(x)/x, Eq.(34) can be written asFe = -i k 0 Γ 2N drn(r) β(r) ∂ ∂(k 0 z) dr ′ n(r ′ )j 0 (k 0 |rr ′ |)β * (r ′ ) -c.c. . (B2)Using the expansions (7) and (9) and Eqs.(8) and(11), we obtaindr ′ n(r ′ )j 0 (k 0 |rr ′ |)β * (r ′ ) = λ n j n (k 0 r)Y * nm (θ, φ). j n (k 0 r)Y * nm (θ, φ) -c.c. ∂j n (k 0 r) ∂(k 0 r) cos θ Y pq (θ, φ)Y * nm (θ, φ) + j n (k 0 r) k 0 r sin 2 θ Y pq (θ, φ) ∂Y * nm (θ, φ) ∂(cos θ) i (α pq α * nm -c.c.) . (B5)Assuming as beforeα nm = α n δ m,0 ∂j n (k 0 r) ∂(k 0 r) cos θ P p (cos θ)P n (cos θ) + j n (k 0 r) k 0 r sin 2 θ P p (cos θ) ∂P n (cos θ) ∂(cos θ) i (α p α * n -c.c.) . sin θ cos θ P p (cos θ)P n (cos θ) = 1 -1 dxxP n (x)P p (x)and π θP p (cos θ) ∂P n (cos θ) 1)xP p (x) = (p + 1)P p+1 (x) + pP p-1 (x),(x 2 -1) dP n (x) dx = n [xP n (x) -P n-1 (x)] n (x)P p (x) = 2 (2p + 1)(2n + 1) {(p + 1)δ n,p+1 + p δ p,n+1 } and 1 -1 dx(x 2 -1)P p (x) d dx P n (x) = 2n (2n + 1)(2p + 1) {p δ p,n+1 -(n + 1)δ n,p+1 } .By substituting these expressions in Eq.(B6) we obtainFe = -k 0 Γ 2N ∞ p=0 ∞ n=0 i(α p α * n -c.c.) λ n (2n + 1)(2p + 1) ∞ 0 drr 2 n(r)j p (k 0 r) × × (p + 1) ∂j p+1 (k 0 r) ∂(k 0 r) + (p + 2) j p+1 (k 0 r) k 0 r δ n,p+1 + p ∂j n (k 0 r) ∂(k 0 r) p α * n -c.c.) λ n λ p (2n + 1)(2p + 1) {(p + 1)δ n,p+1 -pδ p,n+1 } . (B8)Eliminating one of the two sums, we obtain:Fe = -k 0 Γ 4πN ∞ n=0(n + 1)λ n λ n+1 (2n + 1)(2n + 3) i(α n α * n+1 -c.c.). (B9)

  |g 1 . . . g N |0 k + e -i∆0t

		N
		β j (t)|g 1 . . . e j . . . g N |0 k +	γ k (t)|g 1 . . . g N |1 k
		j=1	k
	N
	+	ǫ m<n,k (t)|g 1 . . . e m . . . e n . . . g N |1 k ,
	m,n=1

VIII. ACKNOWLEDGEMENTS

We acknowledge fruitful discussions with R. Kaiser and T. Bienaimé. We are also grateful to F. Staniscia for support for the numerical simulations. The work has been supported by the COSCALI network within the IRSES program of the European Commission under grant number PIRSES-GA-2010-268717.

Appendix A: Derivation of the scattered electric field

The scattered radiation is provided by the positive-frequency part of the electric field operator

where E k = ( ω k /2ǫ 0 V ph ) 1/2 . Integrating Eq.( 5) with a k (0) = 0, inserting it in Eq.(A1) and approximating the sum over the modes k by an integral, we obtain

Introducing spherical coordinates, dk = dkk 2 dφdθ sin θ, and integrating the angular part, Eq.(A2) becomes

Assuming the radiation spectrum centered around k ≈ k 0 , we approximate kE k g k ≈ k 0 E k0 g k0 . Then, extending the lower limit of integration of k to -∞, we obtain for t < |r -

Neglecting the radiation retard in the limit t ≫ σ R /c, where σ R is the cloud size, and approximating in the far field limit |rr m | ≈ r -n • r m where n = r/r, Eq.(A4) becomes ÊS (r, θ, φ, t) ≈ -dk 2 0 4πǫ 0 r e ik0(r-ct)+i∆0t N m=1 σm (t)e -iks•rm , (A5)

where k s = k 0 (sin θ cos φ, sin θ sin φ, cos θ). When applied on the state (2), neglecting virtual transitions, it yields ÊS |Ψ = E S |g 1 . . . g N and where E S , in the continuous density approximation will be given by Eq. (22).

Appendix B: Derivation of Eq.(36)

The angular integral in Eq.( 34) is 2π 0


Michel Perrin
email: michel.perrin@ensmp.fr

Laura S Mastella
email: laura.mastella@gmail.com

Olivier Morel
email: o.morel@brgm.fr

Alexandre Lorenzatti
email: alexandre.lorenzatti@inf.ufrgs.br

Geological Time Formalization: an improved formal model for describing time successions and their correlation

Keywords: Geological time scales, Geological dating, Ontologies, Codification, Stratigraphic correlation

Geological time description largely rests on an event based chronology based on the stratigraphical model. It uses a hierarchy of chronologically ordered geochronological units and boundaries. In order to be easily dealt with within large databases used by complex engineering systems, the geological time chronology must be formalized. Stratigraphical time successions should accordingly be described by using adequate semantic tools (ontologies) complemented by a set of logical rules. At present, geological time formalization mainly rests on the GeoSciML model. This model is fit for describing individual geological time scales but does not provide all the necessary tools for comparing various time successions and for operating full stratigraphic correlations. For complementing the GeoSciML model, we define two ontologies for geological time description and for geological dating. They extend the GeoSciML model, so that it becomes possible to fully use the Allen rules for operating time correlations between any couple of time scales or stratigraphic successions. We additionally propose a codification resting on the defined ontologies, which allows operating all age identification and correlation by means of simple computation rules.

Introduction

Geological time description and rock dating

Geology is in many ways a historical science. Information concerning the date in the course of geological times at which some geological process took place or the age of some given geological formation is a significant part of geological information in general, essential in many cases for geologists and for geology users. Geological time is currently described in two different ways. The first one rests on a quantitative chronology based on "absolute ages" expressed in millions years (My), which are established by means of radiometric measurements [START_REF] Hardenbol | Mesozoic and Cenozoic sequence chronostratigraphic framework of European basins: Mesozoic and Cenozoic sequence stratigraphy of European basins[END_REF]. The second one, more frequently used by geologists, uses an event based chronology based on stratigraphic time scales. Event based chronology is commonly used in archaeology and in history. It relies on time scales composed of a succession of time units, each corresponding to the span of time which separates two definite events. Moreover, as shown in figure 1, a period of rank 1 is likely to be recursively divided in sub-periods of higher ranks limited by less important events.

Historically, the main events that were considered for defining stratigraphic time scales were the appearance or disappearance of significant fossil species. However, since adequate fossil occurrence data were not available everywhere, other events were also used for establishing these stratigraphic time scales, such as lithostratigraphic facies occurrences, sequential stratigraphy and geochemistry, palaeomagnetism, solar flux cycles (cf. [START_REF] Ogg | A geologic time scale[END_REF]. The universal reference for stratigraphic dating is the International Stratigraphic Scale (ISS). This standard time sale and the reference rock records on which it rests -formerly stratotypes presently replaced by GSSPs (Global Stratotype Section and Point), each representing the point in time at which a particular stage is starting [START_REF] Ogg | A geologic time scale[END_REF]-are established by the International Commission on Stratigraphy 1 of the International Union of Geological Sciences (IUGS) 2 . Since the ISS can only be used for dating successions containing adequate fossil associations, complementary time scales were defined at the regional scale in various parts of the world and were correlated with the ISS by all available means. Finally, at local scales, rock dating is currently operated by reference to local stratigraphic successions, whose units are themselves correlated with the ISS or with some regional time scale. Fig1. Example of an event based chronology in the field of history: the Pharaoh dynasties3

Compared to the time scales used in history, the peculiarity of the time scales used in geology is the fact that they are based on the stratigraphic model [START_REF] Jackson | Glossary of Geology[END_REF]. Used for describing all types of stratigraphic successions, this model is based on the "superposition principle" [START_REF] Tarbuck | Earth: an introduction to physical geology[END_REF]. It establishes a correspondence between space and time, considering that any sedimentary succession observed in the field materializes a particular slice of geological time. As pictured in figure 2, each time unit or time boundary belonging to a given time succession, corresponds to a particular geological unit or to a particular geological boundary. The model of time succession pictured in figure 2, rests on the assumption than any time limit corresponds to one atomic time instant in the course of geological times and that, consequently, each time unit begins and terminates at definite time instants. However, this is a drastic simplification of reality since the actual ages of these various limits expressed in My can be determined only with some approximation. One reason for this is the uncertainty attached to radiometric measurements. Another reason, even more important, is the fact that the geological boundaries to which these measurements are related, do not always have one definite age. This notably happens when the geological limits considered are based on lithostratigraphy with no direct reference to fossil occurrences. For these various reasons, rock dating can hardly be seen as a straightforward operation. It always depends in some way on more or less complex geological interpretations, which may generate contradictory solutions. We will further examine how some these difficulties can be taken into account.

Need of geological time formalization

Management and exchange of information within large databases and complex systems is a capital issue in earth sciences as in many other fields. Typical examples concern for instance; -the edition of geological maps, which synthesize a very large number of atomic field and rock sample observations that need to be periodically incremented, updated and attached to the various geological objects represented on the maps; -practical activities such as exploration of water or mineral resources, underground waste storage, civil engineering, development planning, which require considering a large number of geology data that must be easily accessed by various kinds of users such as geologists, engineers, decision makers; -3D or 4D earth models performed for water or hydrocarbon resource exploration or for underground storage, which result of a complex chain of operations involving data of various kinds: geophysical, geological, petrophysical [START_REF] Perrin | Knowledge-driven applications for geological modeling[END_REF], [START_REF] Mastella | Knowledge Management for Shared Earth Modelling 69 th EAGE SPE/EUROPEC Conference[END_REF].

In the age of internet, managing the information required for dealing with these problems should be preferably envisaged as a computer-assisted knowledge management activity using web resources. This supposes building research engines and specialized data bases for dealing with geological knowledge. And, as [START_REF] Richard | Geoscience concept models[END_REF] rightly points out, this further supposes that "geology as a science [should be object of] a greater degree of formalization in order to take advantage of evolving computer-aided knowledge representation and analysis systems". Knowledge management thus supposes knowledge formalization

3. Goal of the present paper

Examining the issues related to the formalization of the event based chronology currently used for dating rock successions is the exclusive subject of the present paper. Formalization supposes drastic simplification. For this reason, we will not be concerned here with the complex geological issues related to establishing time scales and more generally stratigraphic successions but only with the formal aspects of their description and of their correlation. For accommodating the uncertainty attached to geological ages, we will just consider that the ages of time boundaries may be determined only with some approximation and that, consequently, the correlation of two stratigraphic successions may eventually be object of contradictory interpretations. This formal issue is illustrated in figure 3. In the present paper, after having listed the requirements for geological time formalization, we will examine how these requirements are met by the GeoSciML model currently used. We will then propose a few improvements, which open the possibility of not only describing individual time scales but also operating correlations between stratigraphic scales and/or stratigraphic successions in general. We will finally propose a codification method that enables actual stratigraphic successions to be described and correlated.

C C 3 R 3 R C 2 ? R 2 C 1 ? ? R 1 B ? ? Q ? P A C C 3 R 3 R C 2 R 2 C 1 R 1 B Q A P C C 3 R 3 R C 2 R 2 R 1 C 1 B Q P A a) Correlation issues b) Interpretation 1 c) Interpretation 2

Issues attached to geological time formalization.

We have listed in section 1.2, some fields of activity and some issues for which geological time formalization is required. We should now examine which requirements are attached to these various issues.

Geochronologic hierarchy

For being accessed within a database, objects need to be classified and, in our view, such a classification should preferably be knowledge based. This supposes considering the semantics of the various objects and attaching them to concepts, which must themselves be hierarchically classified. In our case, considering that, within a definite time scale, any geochronologic unit of rank i corresponds to a fraction of the time span of the related unit of rank i-1, this hierarchy corresponds to a partonomy (cf. infra § 3.1). Moreover, since the various units of rank i attached to a geochronologic unit of rank i-1 are chronologically ordered with respect to each other, this partonomy is chronologically ordered.

Synonymies

A user who searches documents related to a particular geological period, for instance late Triassic, will be wanting to retrieve from the data base, not only the documents, which explicitly mention "late Triassic" but also those containing the word "Keuper". For a broader research, he/she may also be interested in retrieving documents which contain the word "Triassic" or, on the contrary, documents containing the names of some particular divisions of late Triassic such as "Carnian", "Norian", "Rhaetian" or even "Lettenkohle".

There exists, attached to a given set of geochronologic units, a dual set comprising the associated time boundaries. In order to be accessed by users within a database, these various boundaries must also be described and classified by specifying their links with the limiting geochronologic units. This supposes to take into account many synonymies since the top boundary of a given unit is equivalent to the bottom boundary of the unit, which directly overlies it and since a given limit is likely to have different names, depending on the rank of the units to which it is attached. An example is shown in figure 4, where all the listed terms correspond to one same geochronologic boundary corresponding to the beginning of the Triassic period.

Base of

Mesozoic

Chronological relationships

Since geochronologic units and boundaries are chronologically ordered, it is necessary to specify the temporal relationships, which exist between them. The various possible relationships between time intervals were examined by [START_REF] Allen | Maintaining knowledge about temporal intervals[END_REF], who proposed 13 basic relations between two time intervals. Figure 5 illustrates the application of these relations to stratigraphic units and boundaries. This set of relationships will be largely used hereafter in our ontological model and for exploiting the codification of geological scales that we will propose in section 4 of the present paper. [START_REF] Allen | Maintaining knowledge about temporal intervals[END_REF]. Seven relationships are shown here. Six additional relationships could be deduced by inverting U1 and U2 in those of the above relationships that are not symmetric (i.e. all relationships except the last one)

Geological dating

Formally, dating consists in establishing a link between a geological age (stratigraphical or absolute) and some geological object. This supposes that not only geological ages but also geological objects be described in a formal way and that the significance and formats of the established dating links be specified. Geological objects are various and may be simple (for instance: a stratigraphical surface, a single sedimentary strata, an elementary fault,) or complex (for instance a full stratigraphical column, a set of deltaic formations, a fault network). Their formal description supposes that these various objects should be considered as combinations of elementary types defined within some conceptual model [START_REF] Richard | Geoscience concept models[END_REF]. . Specifying the age of a definite object is not always a simple issue since some objects were eventually the results of complicated successions of geological events. In order to avoid all possible ambiguities, it is thus necessary that geologists specify which event(s) the age of a given object should be attached to. We should also consider that the age of a given object may correspond either to a given geological date corresponding to some event considered as "instantaneous" at geological time scale or to a geological time span that begun at a given geological date and ended at another one. Formalization should take into account the above mentioned peculiarities and allow age attribution according to two different time formats, one referring to a quantitative chronology expressed in absolute ages in My and the other to an event based chronology defined in reference to some stratigraphic time scale.

Accommodating various hypotheses

Since the age of a geological object is not an objective datum but the result of geologists' interpretation, different ages may eventually be attributed to one same geological object. At a larger scale, this is likely to generate different, possibly contradictory versions of the maps and models related to a given geological area. For comparing and evaluating these different versions, a record should be kept of all the elementary interpretations, which generated them [START_REF] Rainaud | A Short History of the Last 15 year's Quest for IT Interoperability in the Petroleum E&P Industry[END_REF]. Consequently, the formalization of geological time that is required should be able to accommodate age uncertainties and to represent many interpretations related to the ages of elementary objects and/or of more or less complex geological assemblages.

State of the art of geological time formalization

Ontologies

Popularized by authors such as [START_REF] Gruber | A translation approach to portable ontology specifications: Knowledge acquisition[END_REF], [START_REF] Uschold | Ontologies: principles, methods, and applications[END_REF], [START_REF] Guarino | Formal ontology in information systems[END_REF], [START_REF] Noy | Ontology development 101: A guide to creating your first ontology[END_REF] and many others, ontologies have become in two decades a classical AI tool for formalizing technical knowledge. Presently, a widely accepted definition of an ontology in the field of computer science is that given by [START_REF] Gruber | A translation approach to portable ontology specifications: Knowledge acquisition[END_REF]: "an ontology is an explicit specification of a conceptualization". In this definition, the word "conceptualization" designates an abstract model of things that are assumed to exist in some area of interest (objects, relations) and the expression "explicit specification" intends to specify that the concepts and relationships in the abstract model are given explicit names and definitions. So, in contrast to other kinds of models, such as numerical models for instance, an ontology is a symbolic model of the types of objects attached to some domain. It represents a domain conceptualization by means of words and of their meanings. Various ontologies have been defined for earth sciences (cf. [START_REF] Sinha | Geoinformatics: Data to Knowledge[END_REF]. Ontologies defined by [START_REF] Richard | Geoscience concept models[END_REF] and more specifically by [START_REF] Cox | A formal model for the geologic time scale and global stratotype section and point, compatible with geospatial information transfer standards[END_REF] for the description of knowledge attached to geological time, were at the origin of the GeoSciML model , which is the only significant knowledge formalization presently in use in the field of geosciences.

Before entering into the details of GeoSciML model, it is necessary to recall a few important points concerning ontologies.

i. An ontology is most often built in view of a definite practical goal. For this reason, "there is no one correct way to model a domain [but] there are always viable alternatives" [START_REF] Noy | Ontology development 101: A guide to creating your first ontology[END_REF]. This is a key point for understanding the choices that were made by the various categories of geoscientists, who have already proposed solutions for geological knowledge formalization. ii. An ontology is not a software model but merely an abstract model. For using such a model in a software application, it is necessary to formalize it into a logical model with the help of an ontology language (Gomez-Perez & al., 2004). OWL 4 (McGuinness & Van Harmelen, 2004) for instance is the standard language proposed by W3C to formalize knowledge, i.e. make the knowledge from a domain processable by a computer. OWL is based on Description Logics, which is considered one of the most important knowledge representation (KR) formalism, which unify all the well known KR formalisms5 . iii. In order for end users to check its relevance, the structure of an ontology must be visualized in some way.

Since ontology languages such as OWL cannot be easily read by end users, ontologies must be visualized using graphical representations (for instance, the Protégé platform (Noy & al., 2001) for visualizing ontologies expressed in OWL). However, there presently exists no solution allowing end users to visualize ontologies in a fully satisfactory way.

A formerly developed geological time formalization: the GeoSciML model

Till now, only one significant concept model built by the Arizona Geological Survey [START_REF] Cox | A formal model for the geologic time scale and global stratotype section and point, compatible with geospatial information transfer standards[END_REF] has been proposed for formalizing Geological Time Scales. This same model is presently being reused by the working group GeoSciML, which involves geological surveys from various countries (notably Australia, Canada, France, Italy, UK, USA) with the goal of unifying and formalizing geological knowledge for the sake of geological map editors. The model, which is in free access on the GeoSciML website6 , will be described hereunder as the GeoSciML model for Geological Time Systems or more simply as the GeoSciML model. The GeoSciML model is written in the UML graphical language and consists of a set of UML diagrams visualized by means of XML representations. The part of the model describing Geological Time Systems refers to the ISO model for Temporal Reference System (ISO 19108) and is thus compatible with geospatial information transfer standards. The ISO 19108 model distinguishes 3 types of reference systems: one based on years, months and days (TM_Calendar), one based on hours, minutes and seconds (TM_Clock) and one based on named intervals (TORS = TM_OrdinalReferenceSystem). The latter is of major importance in our case since a TORS is composed of an ordered sequence of one or more component TOE (TM_OrdinalEra) elements (the term "Era" having here a generic meaning and being thus different from the geological concept "Era", which refers to geological units of rank 1). A TOE can be for instance some geological time interval (such as a geological Eon), which can be recursively decomposed into ordered TOE elements (Eras, Periods, Epochs, Ages, Chrons). This allows a hierarchical age system to be constructed.

With respect to the ISO 19108 Temporal Reference Model, some variations were introduced into the GeoSciML model represented in figure 6. The major one consists in defining the era limits as TOEBs (TimeOrdinalEraBoundary), each TOEB corresponding to one TM_Instant, which is a point in the span of Geological Time. A Geologic Timescale is thus a kind of TORS, which includes as first-order elements: 1) GeochronologicEra, which is a kind of TOE with boundaries defined as TOEBs. GeochronologicEra instances may have names such as Triassic, Oxfordian etc.

2)

TOEB having two specializations: a. GeochronologicBoundary defined with reference to some geologic evidence such as, for instance, fossil appearance/disappearance or significant change in local rock lithology. b. NumericEraBoundary corresponding to geological boundaries dated by absolute age measurements.

The correlation of the timescale elements GeochronologicBoundary and GeochronologicEra with actual geological records is established by means of the two relationships:

GeochronologicBoundary ---isManifestedBy---> ChronostratigraphicBoundary GeochronologicEra ----------isManifestedBy---> ChronostratigraphicUnit

The concept of ChronostratigraphicBoundary corresponds to an actual reference Stratigraphic_Boundary (for instance a GSSP) and that of ChronostratigraphicUnit to an actual reference Geological_Unit (for instance a stratotype). Moreover the concept of ChronostratigraphicUnit may eventually be further characterized as a LithostratigraphicUnit with reference to its lithological content.

The GeoSciML model is conceived to record the characteristics of various geological data from which the Geological Time Scale was built, i.e. the data attached to the stratotypes and to the Global Stratotype Sections and Points [START_REF] Ogg | A geologic time scale[END_REF]. Generally speaking, a geological object described in the GeoSciML model as a GeologicalFeature, can be given an age by using the relations:

GeologicalFeature ----is linked to -----> GeologicalEvent ---has---> Age

Geological dating thus goes through associating a geological object to the GeologicalEvent that generated it. This is a minor difference with respect to the concept model formerly proposed by the Arizona Geological Survey [START_REF] Richard | Geoscience concept models[END_REF], which establishes a direct link:

GeologicalUnit ---has ----> GeologicalAge

Need of complementing the existing GeoSciML ontology

The main goal of the GeoSciML model for Geological Time Systems is formalizing the knowledge attached to the International Stratigraphic Scale (or possibly to some other Geological Time Scale) and to the geological records, stratotypes and GSSPs, which define it. However, as we mentioned in section 1.2., in their everyday practice, geologists also need to operate correlations between stratigraphic successions by specifying chronological relationships or eventual synonymies between the various items attached to these successions. At present, this cannot be done in a simple way by using the GeoSciML model.

Developing new ontologies or significantly modifying those which already exist is an issue that must be carefully thought about in a field like geology, where significant efforts are being made to unify and formalize vocabulary. But, considering the needs that are not covered by the GeoSciML model and which are of importance for geomodelers like us, we have considered the possibility of proposing improved ontologies for Geological Time description and for Geological Dating. The ontologies that we will hereafter present in section 3 are fully compatible with the GeoSciML model and are in no way redundant or contradictory with this model. Furthermore, these improved ontologies open the possibility of formulating new rules for solving practical issues concerning geological dating and stratigraphic correlation by means of a codification methodology that will be presented in section 4.

Proposal for Geological Time Formalization

The ontologies that we present hereafter for geological time description and for geological dating were partly inspired by the GeoSciML model. Their interest is that they open the possibility, which does not exist in GeoSciML, of performing full mutual comparison between stratigraphic successions of any type. A detailed comparison with the GeoSciML model will be presented at the end of the present section. A preliminary version of these ontologies was given in Mastella (2010)7 .

The Geological Time Ontology

The Geological Time ontology was developed for describing the hierarchy of the geological periods of time as it appears in stratigraphical time scales and for further allowing establishing correspondence between different time scales based either on the use of fossils or on absolute ages. The two main concepts that were defined for formalizing the main elements of a geological time scale are GeochronologicUnit, which corresponds to a geological time interval, and GeochronologicBoundary, which represents a geological time boundary corresponding to an instant having no temporal duration.

Fig 7. Proposed Geological Time ontology

The Geological Time ontology is represented in figure 7. It can be described in the following way.

i. The abstract concept GeochronologicElement is the superclass of GeochronologicUnit and GeochronologicBoundary. ii. The class GeochronologicUnit is defined by reference to a StratigraphicReferenceSystem (a TORS as defined in GeoSciML). iii. Age relationships between two GeochronologicElements of any nature can be established by the relations isYoungerThan and isOlderThan. This allows specifying the order of occurrence of these elements in the course of geological times, even in the case when they are not attached to the same StratigraphicReferenceSystem. iv. GeochronologicUnit instances relate to other GeochronologicUnit instances by means of detailed interval relationships (e.g. overlaps, meets, starts, etc.). These relations are those defined by the Allen's rules as mentioned in section 2.1.3. As illustrated in figure 8, it is thus possible to precisely describe age relationships between two GeochronologicUnits eventually belonging to different StratigraphicReferenceSystems, notably in the cases when the two units meet or are included one in the other. v. GeochronologicUnit instances (the equivalents to GeochronologicEra instances) are organized in a partonomy (i.e. by a partOf relation). Units such as Eon, Era, Period, etc. are sub-concepts of GeochronologicUnit and are organized in specific partonomies: an instance of Chron is part of some Age, which is part of some Epoch. Actual GTS units, such as Triassic, Jurassic, and so on, are represented as instances of the concept GeochronologicUnit (in fact, as instances of the concepts Eon, Era, and so on). Actual boundaries between units are represented as instances of the concept GeochronologicBoundary. One can stipulate that definite GeochronologicBoundaries correspond to the base or to the top of some GeochronologicUnit by using hasBase and hasTop relationships. vi. GeochronologicInstant is a generalization of GeochronologicBoundary representing one particular instant within a GTS, which may correspond or not to a boundary between GeochronologicUnits. The age of GeochronologicInstant may be expressed by an AbsoluteAge (for 1.5 My) or as a stratigraphical age (Lower Pleistocene) or remain purely virtual.

Fig. 8. Examples of relationships between time intervals

(in reference with the interpretation 1 envisaged in figure 3b).

Unit Q begins at the same time as unit B is expressed by the relation:

starts(Q,B)
Unit R 2 begins after the beginning of unit C and terminates before the termination of unit C is expressed by the relation: during(R 2 ,C).

Equivalence between units P and A is established by the relation: isEquivalentTo(P,A).

C C 3 R 3 R C 2 R 2 C 1 R 1 B Q A P
The position of this Geological Time ontology with respect to the GeoSciML model can be appreciated in the following way.

i. The concepts of GeochronologicUnit and GeochronologicBoundary are respectively equivalent to the GeoSciML classes GeochronologicEra and GeochronologicBoundary. We prefered the name "Unit" to "Era" in order to avoid any possible confusion with Geological Era (= GeochronologicUnit of rank 2). ii. The class GeochronologicUnit is more general than the GeoSciML class Geologic Time Scale, since it includes stratigraphic successions of any type. The class StratigraphicReferenceSystem is a TM_OrdinalReferenceSystem (TORS), as defined in GeoSciML. iii. Relative age relationships between GeochronologicEra and/or GeochronologicBoundary within a given Time Scale can be established in the GeoSciML model by considering the hierarchy of the TM_OrdinalReferenceSystem. This is also possible in our proposed Geological Time ontology. However, the relations isYoungerThan and isOlderThan completed by the detailed interval relationships issued from Allen's rules open the possibility of operating detailed age comparisons between GeochronologicElements, which do not belong to the same Time Scale. Such possibility is an original addition with respect to the GeoSciML model. iv. The partonomy of our GeochronologicUnits is strictly equivalent to that of GeoSciML GeochronologicEras. The relationships between GeochronologicBoundary instances and the associated GeochronologicUnit are also strictly equivalent to the GeoSciML relationships: GeochronologicBoundary ---starts/ends---> GeochronologicEra. v. Our class GeochronologicInstant is equivalent to the TM_Instant class of GeoSciML. However, we proceeded to some simplification by ignoring the GeoSciML concept of NumericEraBoundary , which did not appear necessary in view of our goal of stratigraphic correlation and, as a matter of consequence, by also ignoring the concepts of TimeOrdinalEraBoundary (TOEB) and NumericEraBoundary. Compared to the GeoSciML model, the main peculiarity of the Geological Time ontology that we have defined is that it can be applied for geological time correlation, that is for creating a correspondence between elements of different GTSs. For depicting geologists' interpretations, correlations between units can be created by means of the time interval relations.

Geological dating

Defining geological objects

The geological dating procedure consists in assigning an "age" to some geological object. As we explained before, this first supposes that geological objects should be formally defined. For this, we use the ontology of Basic Geology defined by the e_Wok Hub consortium8 [START_REF] Aït Ameur | Semantic Hubs for Geological Projects: Workshop on Semantic Metadata Management and Applications[END_REF], which is represented in figure 9. This simplified model is strictly in accordance with the GeoSciML model.

The Basic Geology ontology is built around the concept GeologicalObject, which is constituted by GeologicalUnits and GeologicalBoundaries. Geological objects are very various (examples among many others are: a stratified sedimentary unit, a reef, a diapir, a fault network etc.) and can be simple or complex. Complex geological objects can be described as associations of a various number of atomic geological objects belonging to two categories: A given geological object is the result of some geological event (represented by the concept GeologicalEvent). A geological event may consist of a single geological process (e.g. the deposition of a sedimentary unit) or be composed of multiple geological processes (example: a metamorphic formation deformed by late tectonics). The various specializations of GeologicalProcess correspond to thecreation, destruction or transformation of geological matter. They were detailed by the e_Wok Hub consortium in a subontology entitled Geological Process. A geological event is also responsible for geological object structures (e.g. Synform Fold, Reverse Fault), which are detailed by the e_Wok Hub consortium in a sub-ontology Geological structures.

Fig 10. Geological units and boundaries

The Geological Dating Ontology

The Geological Dating ontology that we have defined is based on the data model of the Geological Time Scale presented in the GeoWhen Database9 . GeoWhen is a compilation of the main naming schemes for the geologic timescale that appear in the literature. We modified it in many points, and we imported the temporal relations already formalized in the OWL-Time ontology [START_REF] Hobbs | Time ontology in OWL, Ontology Engineering Patterns Task Force of the Semantic Web Best Practices and Deployment Working Group[END_REF]. The Time Ontology is an ontology of temporal concepts originally developed for describing the temporal content of Web pages and the temporal properties of Web services. This ontology can be easily reused for adding temporal aspects to any domain ontology.

Our Geological Dating ontology represented in figure 11 allows representing the main characteristics of geological dating and notably the two different ways in which geological ages can be expressed i.e. absolute dating and relative dating using stratigraphic scales. The ontology introduces abstract concepts, which make the link between concepts of the Geological Time ontology and of the Basic Geology ontology. The concepts imported from the Geological Time ontology were given the prefix GeoTime, and those imported from the Basic Geology ontology the prefix BasicGeo.

Fig 11. Proposed Geological Dating ontology

The concepts and relationships related to the Geological Dating ontology can be described as follows: i.

The geological objects that can be dated are: BasicGeo:GeologicalBoundary, BasicGeo:GeologicalUnit, BasicGeo:GeologicalEvent and BasicGeo:GeologicalObject.

Ii

All geological objects imported from the Basic Geology ontology have a defined stratigraphic age, which is represented by the relation hasStratigraphicAge with the Geological Time ontology concept GeoTime:GeochronologicUnit. iii. GeoIntervalObjects are abstract classes that are superclasses of all interval-like objects imported from the Basic Geology ontology. They are related to GeoTime:GeochronologicUnit through the relations meets/overlaps/finishes/during. These objects also have a duration. They can be dated by specifying a stratigraphic age, by defining an interval relation with some GeoTime:GeochronologicUnit or by stating that the object begins or ends on some GeoTime:GeochronologicUnit. iv. A GeologicalBoundary can be also dated by stating a stratigraphic age. Moreover, since it is a boundary, it may be dated by a specific instant in time, represented by the concept GeoTime:GeochronologicInstant. This instant can be characterized by an absolute age value or be considered as equivalent to a GTS boundary. It may also remain a purely virtual point in GTS, since GeoTime:GeochronologicInstant is a superclass of GeoTime:GeochronologicBoundary.

Since we were not concerned by the issue of linking timescale elements to the geological records, which define them, we have defined, in contrast with GeoSciML, geological dating in one simple way:

by associating any GeologicalBoundary to a GeochronologicInstant, -by associating GeochronologicUnit to GeologicalUnits (or more generally to GeologicalObjects) and to GeologicalEvents, This solution allows considering that actual geological objects on the one hand and items belonging to some time scale one the other, can be respectively attached to the same concepts in the following way:

GeologicalBoundary --associated with ---> GeochronologicInstant GeochronlogicBoundary GeologicalUnit --associated with ---> Geochronologic Unit GeochronologicUnit

Conclusion on the proposed formalization

Although being largely inspired by the GeoSciML model and being fully compatible with it, the two ontologies that we developed for formalizing Geological Time and Geological Dating open new possibilities. By introducing a class GeochronologicUnit which not only represents time scales but any type of stratigraphic succession and by formally expressing detailed interval relationships issued from Allen's rules, we open the possibility of comparing the ages of GeochronologicElements eventually belonging to different time successions. In contrast to the GEoSciML, our improved model thus enables the user to easily operate stratigraphic correlations. The full set of ontologies formalized in OWL language is presented and can be downloaded from the web site of the e-Wok project.10

Proposal for an Ontology Based Codification of Geological Time

Introduction

For operating the concepts and properties described by means of the two above defined ontologies, we propose to characterize the various instances of GeochronologicUnit and those corresponding to the associated GeochronologicBoundaries within one time succession by means of a codification bearing in itself all the information related to these entities. Such a codification should enable one to chronologically classify the various GeochronologicElements belonging to one time succession and to determine all the instances associated to a given GeochronologicInstant. Moreover, when applied to two time successions, this codification should enable establishing time correlations between these entities.

The codification that we propose, intends to be applicable to any chronological succession based on named intervals i.e. to any TM_OrdinalReferenceSystem (TORS) as defined in the ISO 19108 model and in the GeoSciML model. We will describe here the principle of this codification and briefly explain how it can be used for classifying GeochronologicElements and for establishing correlations between different time successions

Codification rules for GeochronologicUnits and GeochronologicBoundaries

The codification that we propose basically rests on a chronologically ordered numbering. Two sets of code numbers are defined one for time units (U_code) and one for time boundaries (B_code). Figure 12 shows how the related concepts can be inserted in the Geological Time ontology.

Let us consider a TORS comprising time units of n ranks. The U_code corresponding to some time unit U i of rank i that is a division of a parent time unit U i-1 of rank i-1 will consist in :

i. a prefix for identifying the considered TORS and for specifying that the code is a U_code. For instance, for a GeochronologicUnit belonging to the ISS, this prefix could be ISU = IS (International Scale) + U (U_code). ii a sequence of n pairs of numbers corresponding to the n ranks of the considered TORS. The i-th pair of numbers ≠ 0 comprises from left to right : a. an ID number indicating the position of U i within the division of the parent unit U i-1 b. a ND number equal to the total number of divisions of U i-1 .

In a similar way, the (i-1)th pair will consist in an ID number describing the position of unit U i-1 within the division of the parent unit U i-2 and in a ND number equal to the total number of divisions of U i-2 .

In the case when i < n, the last n-i pairs of numbers of the code are given the values 00. Moreover additional pairs 00 can be eventually added to the code for facilitating time scale correlations as it will be explained later in section 4.3. For attributing B_codes to GeochronologicBoundaries, we propose to apply two rules. Rule 1 specifies that the B_code attached to the lower time boundary related to any "ultimate" time unit of the TORS, (i.e. to any time unit corresponding to a subdivision, which is not itself divided and which is thus locally the "ultimate" in the partonomy), bears the same number as the ultimate time unit to which it is attached. This B_code will just be differentiated from the U_code to which it is related by the ending of its prefix where the letter U will be replaced by the letter B (the prefix ISU being for instance replaced by the prefix ISB in the case of the ISS).

However, the lower time boundary of some "ultimate" time unit may also be the boundary of time units of lower ranks (i.e. of time units that are not ultimate). In this case, rule 2 specifies that the B_code attached to any boundary possibly related to several time units of different ranks is the one that is defined by considering the ultimate time unit to which this boundary is attached and by applying rule 1. This allows all synonymy problems between boundaries to be easily solved.

Figure 13 illustrates a practical example of codification for a time unit divided into two subunits and seven sub-subunits.

Correlation between two time successions

The codification that we have just described applies to any time succession and, more generally, to any TORS. It can also be used for establishing correlations between two time successions.

In the logical frame that we are using for codification based on GeochronologicUnit numbering, correlation between two time scales A and B is necessarily an asymmetric procedure. This means that correlation can be established in two different equally valid ways, either by correlating B with A or A with B. Figure 14 illustrates the methodology that we propose in the case of a correlation of a time scale B with a time scale A. We suppose that the two scales were codified as described in section 5.1.2, scale A with a code having the headers AU and AB and scale B with headers BU and BB. We also suppose that the user has determined the time relationships (youngerThan/olderThan), which exist between any two time boundaries respectively belonging to scale A and scale B, these relationships being eventually the result of his/her interpretation. In order to allow the correlation, we first add to the B_codes attached to scale A a pair of numbers, which will be given zero default values. Correlation of B with A then basically consists in solving the two following cases concerning time boundaries: i. Some time boundary of scale B has the same time position as some boundary of scale A. In this case, the two boundaries are put as equivalent. In figure 14 left, we thus have the equivalences:

BB 22 00 ≡ AB 22 33 00 , BB 12 55 ≡ AB 22 23 00, BB 12 15 ≡ AB 12 00 00

ii. One or several boundaries bB belonging to scale B correspond to GeochronologicInstants positioned within a given GeochronologicUnit uA of scale A (which means that boundaries bB are younger than the bottom limit and older than the upper limit of uA). In such a case, we create within scale A additional boundaries that we put equivalent to boundaries bB. The codes of these boundaries will be deduced from the code of the upper boundary of uA by replacing the last ID zero number by numbers 1, 2… in the chronological orders of the bBs from the youngest to the oldest one. Beyond that, the number of additional boundaries must also be set in the corresponding number of subdivisions. Considering the example illustrated in figure 14, the following boundaries can be added to time scale A, as shown on the third column of the left part of the figure:

AB 12 00 11 ≡ BB 12 25, AB 22 13 12 ≡ BB 12 35, AB 22 13 22 ≡ BB 12 45

Having solved the two above issues concerning GeochronologicBoundaries, existing relationships between any two time units respectively belonging to scales A and B can then be easily determined by applying the Allen's rules.

Integration of two time successions

Considering the above presented methodology, we can go one step further by operating full integration of the two time scales into one. This is again an asymmetric operation, which can be operated either by integrating B into A or A into B with two different and equally valid results.

Considering the same theoretical case as in section 4.3, our integration methodology is illustrated in figure 15. It consists in considering that the extra GeochronologicBoundaries that were added to time scale A for correlating time scale B, correspond to extra divisions of the GeochronologicUnits of scale A. If the maximum rank of the scale A divisions was originally n (2 in the case considered in figure 14), it will then become n+1 (i.e. 3 in the case of figure 14). It then suffices to recompute the codification of time scale A according to the rules presented in section 4.1, taking into account these new divisions of rank n+1. The right part of figure 15 shows the result that is obtained for the integrated time scale AB corresponding to the integration of B with A. Figure 16 shows an additional example of stratigraphic correlation and integration corresponding to the example presented above in figure 3b.

X Y C C 3 R 3 R C 2 R 2 C 1 R 1 B Q A P

Use of the proposed classification

Thanks to the properties of the proposed classification, it is possible to solve with simple algorithms, all practical questions related to relationships between time units and/or boundaries such as specifying:

-the units of lower ranks into which a given time unit U is inserted and the relationships of U with these units (begins / during / finishes) -the units of higher ranks inserted in unit U and their relationships with U -the unit of the same rank immediately laying over U -the unit of the same rank immediately overlaid by U, -the time units overlaying or overlaid by a given time boundary B To take but one example, let us consider the issue of identifying the time units, which overlay a given time boundary B. This problem can be solved by first identifying the ultimate time unit U overlaying B, whose U_code number is the same as the B_code of B and then by identifying all the successive time units of lower ranks, which are begun by U. For finding those, one can apply the rule:

If, in the U_code of a time unit of rank i begun by U, the ID and ND numbers of the last pair of numbers non equal to zero are equal (i.e. have the format kk), this time unit corresponds to the oldest division of the associated time unit of rank i-1, which is itself, consequently, begun by U.

A software application using an OWL version of the ontologies that we defined as well as our proposed codification, is presently being developed at IFP. Software modules have already been realized for enabling the user to interactively build time successions and correlate them. The codes corresponding to time units and boundaries are generated in an automatic way. The correlation of a time succession B with a time succession A and the merging of A and B in one succession AB are simply operated by interactively inserting with the mouse time boundaries belonging to B into time units belonging to A. The codes of the new succession AB are also automatically generated.

These first results already demonstrate that the proposed codification is a flexible tool to easily manage time successions and for correlating and integrating them according to one or eventually to several interpretations. This is likely to help geologists to manage uncertainties attached to the time position of Geochronologic_Boundaries by rapidly checking various interpretative hypotheses.

Conclusion

At present, many activities using geology rest on large databases and depend on more or less complex engineering systems. One condition for establishing and maintaining such data bases and for efficiently managing the engineering systems, which use them, is adequate knowledge formalization. These issues are the ones which motivated us for reexamining how geologists presently formalize the specific event based chronology that is currently used for geological time description and for dating rock successions and for suggesting some possible improvements

The ISS and the various other event based time successions used in geology rest on the stratigraphic model. Their architectures and mutual relationships depend on the relative chronological order of the various boundaries, which separate time units. For many reasons, mostly geological, the ages of these boundaries can only be determined with some uncertainty. Establishing, maintaining and correlating stratigraphic time scales is thus a difficult task, which rests for a good part on geological interpretation Since formalization supposes drastic simplification, we have decided however not to consider these various geological difficulties but only the formal issues related to the description and correlation of geological time scales and of stratigraphic successions in general.

Considering these formal aspects, we have mentioned that the semantic and logical tools used for describing time scales and stratigraphic successions should enable users:

to correctly describe the hierarchy of geological time units belonging to any time succession as a partonomy of chronologically ordered elements, -to take into account the many synonymies between GeochronologicElements that exist within a given time scale and notably those related to GeochronologicBoundaries, -to operate a detailed comparison of the relative ages of two GeochronologicElements (GeochronologicUnit or GeochronologicBoundary) possibly issued from two different time successions by operating the Allen's rules, -to operate full chronostratigraphic correlation between two time successions as well as eventual merging of one of them into the other.

For fulfilling these requirements and notably those related to chronostratigraphic correlation, we have shown that there is the need of modifying and complementing the widely used formalization of geological time available within the GeoSciML model. We have thus developed two ontologies for describing Geological Time and for operating Geological Dating. These ontologies are derived for a good part from the GeoSciML model. With respect to this model, we only have introduced minimum modifications and extensions, which consisted in defining a limited number of new concepts and a full set of relationships resting on Allen's rules. The practical result is that geological time scales are now considered as being particular cases of stratigraphic successions.

Any two stratigraphic successions can then be mutually correlated one with another by establishing detailed chronological relationships between the GeologicalUnits and GeologicalBoundaries related to each of them. The developed ontologies have been expressed as UML schemas and also implemented as OWL code.

A practical tool has also been defined for allowing users to describe time successions and to easily correlate and integrate them by means of software applications. Our proposed codification is applicable to any time scale or stratigraphic succession and more generally to any TM_OrdinalReferenceSystem (TORS). Within a given time scale, the codes attributed to each time unit (U_code) and to each time boundary (B_code) bear in themselves all the age information related to these entities. GeochronologicalObjects are described in accordance with the model proposed in the two ontologies that we have defined. By means of simple computation rules operated on code numbers, any GeochronologicalObject can be easily retrieved (even when it eventually bears different names) and its age can be easily compared with that of any other GeochronologicalObject. In addition, age relationships established by the geologist between two time successions can be expressed as code number relations. These code number relations are used for expressing specific age relationships between GeochronologicalObjects respectively belonging to the two successions to be compared and for eventually merging these two successions into one.

The solutions that we propose for describing and formally correlating time scales and stratigraphic successions do not contradict the widely used GeoSciML model but complement it. The cores of the two ontologies that we have defined are in strict adequation with the GeoSciML model and the new concepts that we have introduced are mere generalizations of already existing ones. Our proposed codification itself strictly reflects the existing ontologies available for geological time description and for geological dating, so that it can be used with no major difficulty by most existing systems.

Although we paid no attention to the many geological difficulties that exist for establishing and maintaining usable stratigraphic time scales, the knowledge based tools that have developed present, we think, some interest by allowing geologists to operate easy stratigraphic correlations. Considering this particular issue, we have started studying knowledge based user interfaces that will hopefully help geologists to easily appreciate the practical results of many possibly correlation assumptions and consequently to accommodate in some way the many uncertainties, which affect geological ages.

Considering these various aspects, we hope that our present contribution will be judged acceptable by the geoscience community and useful to many.

Fig 2 .

 2 Fig 2. The stratigraphic model: (a) geological succession (b) time succession

Fig 3 .

 3 Fig 3. The stratigraphic correlation issue. Thick lines in (a) indicate correlations that can be established considering lateral continuity or other criteria such as fossil appearance or disappearance while question marks indicate boundaries whose position relative to the other column remains questionable. Boxes (b) and (c) show two possible interpretations for (a).

Fig 4 .

 4 Fig 4. Synonyms for the base of Triassic. The figure shows the various possible names of the boundary corresponding to the beginning of the Triassic in the ISS and in the European continental Triassic stratigraphic scale. Data are extracted from Callec et al. (2006).

Fig 6 .

 6 Fig 6. GeoSciML model for Geologic Time Scales (from Cox & Richard, 2005)

-

  Fig 9. Proposed ontology for Basic Geology

Fig 12 .

 12 Fig 12. Insertion of the U_code and B_code concepts in the Geological Time ontology

Unit to Unit relationships Boundary to Boundary relationships

		Upper U 2 older than Lower U 1
	U 2 older than U 1
		Upper U 2 equivalent to Lower U 1
	U 2 meets U 1
		Lower U 2 older than Lower U 1
		AND Upper U 2 younger than Lower U 1
	U 2 overlaps U 1	AND Upper U 2 older than Upper U 1
		Lower U 2 equivalent to Lower U 1
	U 2 starts U 1	AND Upper U 2 older than Upper U 1
	U 2 during U 1	Lower U 2 younger than Lower U 1 AND Upper U 2 older than Upper U 1
		Lower U 2 younger than Lower U 1
	U 2 finishes U 1	AND Upper U 2 equivalent to Upper U 1
		Lower U 2 equivalent to Lower U 1
	U 2 equals U 1	AND Upper U 2 equivalent to Upper U 1
	Fig 5. Allen's relationships

Rank n-2 Rank n-1 Rank n Fig 13.

  Principle of the proposed codification for time units (TM_OrdinalEras) and time boundaries (TimeOrdinalEraBoundaries) in the case of a TORS comprising n ranks.

			U ...23 12 14	B ...23 12 14
		U ...23 12 00	U ...23 12 24	B ...23 12 24
			U ...23 12 34	B ...23 12 34
	U ...23 00 00	B ...23 12 44	U ...23 12 44	B ...23 12 44
			U ...23 22 13	B ...23 12 44
		U ....23 22 00	U ...23 22 23	B ...23 22 13
	B ...23 22 33	B ...23 22 33	U ...23 22 33	B ...23 12 44

  Integration of time scale B into time scale A resulting into a time scale A

	Time succession A	Integrated time succession AB
	AU 12 00		ABU 12 00 00
		AB 12 00 00		ABB 12 00 00
				ABU 22 13 12
	AU 22 13	AB 12 00 11	ABU 22 13 00	ABB 22 13 12
				ABU 22 13 22
		AB 22 13 00		ABB 22 13 22
				ABU 22 23 13
	AU 22 00	AB 22 13 12	ABU 22 00 00	ABB 22 23 13
	AU 22 23	AB 22 13 22	ABU 22 23 00	ABU 22 23 23 ABB 22 23 23
				ABU 22 23 33
		AB 22 23 00		ABB 22 23 33
	AU 22 33		ABU 22 33 00
		AB 22 33 00		ABB 22 33 00
	Fig 15.

http://www.stratigraphy.org/column.php?id=Chart/Time%20Scale

http://www.iugs.org/

http://www.touregypt.net/kings.htm

http://www.w3.org/TR/owl-ref/

http://dl.kr.org/

http://www.geosciml.org/

http://pastel.archives-ouvertes.fr/pastel-00005770/fr/

e_Wok Hub : Environmental Web Ontology Knowledge Hub, http://www.inria.fr/sophia/edelweiss/projects/ewok/

http://www.stratigraphy.org/bak/geowhen/index.html

http://www-sop.inria.fr/edelweiss/projects/ewok/ontologyview/ontologies.html (in French).

Acknowledgments

We wish to thank the many persons, who have been indirectly involved in the work here presented thanks to their participation to joint research projects with the authors (notably the French ANR project e_Wok Hub) or through various informal discussions. We wish to specially thank Mara Abel (UFRGS, Brazil), Jean-François Rainaud (IFP), Beiting Zhu (IFP), Philippe Verney (IFP), Priscille Durville (INRIA), Olivier Corby (INRIA), Sandrine Grataloup (BRGM), Yamine Aït Ameur (ENSMA), Stéphane Jean (ENSMA).