Jean-Michel André
email: jean-michel.andre1@upmc.fr

Philippe Jonnard

Time-dependent Bragg diffraction bymultilayer gratings

Keywords: time-dependent diffraction, short-pulse propagation, grating, multilayer

The time-dependent Bragg diffraction by multilayer gratings working by reflection or by transmission is investigated. The study is performed by generalizing the time-dependent coupled-wave theory previously developed for one-dimensional photonic crystal [André J-M and Jonnard P, J. Opt. 17, 085609 (2015)] and also by extending the Takagi-Taupin approach of the dynamical theory of diffraction. The indicial response is calculated. It presents a time-delay with a transient time that is a function of the extinction length for reflection geometry and of the extinction length combined with the thickness of the grating for transmission geometry.

Introduction

The time-dependentBragg diffraction of a one-dimensional photonic crystal (1D-PC) modelled by a periodic stack of bilayers has been studied in our recent paper [START_REF] André | Time-dependent Bragg diffraction and shortpulse reflection by one-dimensional photonic crystalsJ[END_REF]. In the present work, we extend our analysis to multilayer gratings that can be regarded as a specific case of two-dimensional photonic crystals. This kind of multilayer gratings has acquired a considerable importance in optics of short wavelength radiation (from UV to hard x-rays) [START_REF] Rife | Multilayer-coated blazed grating performance in the soft x-ray region[END_REF][START_REF] Sammar | Narrow bandpass multilayer x-ray monochromator[END_REF][START_REF] Vinogradov | Gratings and zone plates based on sliced multilayer structures[END_REF][START_REF] Yan | Hard x-ray nanofocusing by multilayer Laue lenses[END_REF]. The advent of short (femtosecond) and even ultra-short (attosecond) sources in these spectral domains has naturally led to a growing interest in the study of the temporal response of these optics. In fact, most of the works relative to the temporal diffraction have concerned the "real" crystals usually implemented in x-ray diffractive optics. Chukovskii and Förster have found analytical solutions of the time-dependent Bragg diffraction by thick crystal [START_REF] Chukhovskii | Time-dependent X-ray Bragg diffraction[END_REF] using a Tagaki-Taupin approach; Graeff has found similar solutions but for Laue geometry [START_REF] Graeff | Short X-ray pulses in a Laue-case crystal[END_REF]. More recently the problem has been treated in reflection geometry and transmission geometry by Lindberg and Shvyd'ko,also within the Tagaki-Taupin approach, with application to self-seeding free-electron lasers [START_REF] Lindberg | Time dependence of Bragg forward scattering and self-seeding of hard x-ray free-electron lasers[END_REF][START_REF] Lindberg | Spatiotemporal response of crystals in x-ray Bragg diffraction[END_REF]. Fourier analysis, which is not very adequate to treat transient phenomena, has been used by many authors to investigate spatiotemporal response of crystals or multilayer optics [START_REF] Ksenzov | Time-space transformation of femtosecond free-electron laser pulses by periodical multilayers[END_REF][START_REF] Bushuev | Diffraction of X-ray free-electron laser femtosecond pulses on single crystals in the Bragg and Laue geometry[END_REF].

In this paper, we implement the time-dependentcoupled-wave theory (CWT) in the twowave approximation,which leads to a system of coupled partial differential equations (PDEs); note that this CWT has been recently unified in [START_REF] Yang | Unified analytical theory of single-order soft x-ray multilayer gratings[END_REF] to treat the diffraction by multilayer gratings in the steady-state case.For the time-dependent case,we solve this system using the matrix methoddeveloped in [START_REF] André | Time-dependent Bragg diffraction and shortpulse reflection by one-dimensional photonic crystalsJ[END_REF]but we also use the Tagaki-Taupin approach (in quasi-specular conditions) to give a physical insight to the problemand to point out the key physical parameters of the problem.We consider multilayer gratings with various shapes (lamellar, sliced, blazed,…)diffracting in reflectionor transmission geometry. Werestrict our analysis to s-polarization (transverse electric) case and planar (non-conical) diffraction since the other situations, p-polarization and conical diffraction,are more complicated to handle because the paraxial approximation cannot be directly applied.

Time-dependent coupled-wave theory

We consider a diffraction multilayer grating as shown in Figure 1 struck by an incident plane wave under a glancing angle𝜃 0 with a wave-vector 𝒌 = 𝑐𝑜𝑠 𝜃 0 𝒙 + 𝑠𝑖𝑛 𝜃 0 𝒛 forming together with the normal to the grating surface,the incidence plane (x, z). The diffraction is considered as planar (not conical). Let us recall that in conical mounting the incident wave-vector is not perpendicular to the grating grooves and the wave-vectors of the different diffraction orders lie on a conical surface; in standard planar mounting that interests us in thiswork, the wave-vectors of the different diffracted orders remain in the incident plane. The considered geometry of the grating covers a large set of gratingshapes, from the simple (Γ = 1, 𝜙 = 0) or sliced multilayer mirror(Γ = 1, 𝜙 ≠ 0) to the lamellar multilayer grating ( Γ < 1 , 𝜙 = 0) or blazed multilayer grating (Γ < 1 , 𝜙 ≠ 0), were Γ is the ratio of the pitch to the grating period and 𝜙 defines the orientation of the multilayer with respect to the surface. Since the propagation takes place within the (𝑥, 𝑧) plane then, in the framework of the paraxial approximation, the wave time-dependent propagation equation is given in Gaussian cgs units by

𝜕 2 𝑬 𝑹 ∥ , 𝑡 𝜕𝑥 2 + 𝜕 2 𝑬 𝑹 ∥ , 𝑡 𝜕𝑧 2 - 1 𝑐 2 𝜕 2 𝑬 𝑹 ∥ , 𝑡 𝜕𝑡 2 = 4𝜋 𝑐 2 𝜕 2 𝑷 𝑹 ∥ , 𝑡 𝜕𝑡 2 (1)
where𝑬 𝑹 ∥ , 𝑡 stands for the electric field and 𝑷 𝑹 ∥ , 𝑡 for the electric polarization vector, c being the speed of light in vacuum and𝑹 ∥ = (𝒙, 𝒛) the position vector in the (𝑥, 𝑧) plane. Let us recall that the paraxial approximation neglects the term 𝛁(𝛁. 𝑬) in the Maxwell equations leading to Eq.( 1); this approximation is justified as long as the electric field remains transverse.

Since weconsider only the s-polarization (transverse electric) case, the electric field vector𝑬 𝑹 ∥ , 𝑡 is along the y axis. It is assumed that the electric field of the optical pulse is formed of a quickly varying carrier with frequency  by an envelope 𝐸 0 𝑹 ∥ , 𝑡 and we write it as follows taking into account the s-polarization geometry

𝑬 𝑹 ∥ , 𝑡 = 𝐸 0 𝑹 ∥ , 𝑡 𝑒 -𝑖 𝜔 𝑡 𝒚 (2)
We assume that the polarization is essentially electronic and follows instantly the change of the electric field, and that the media have a linear response. Hence we write the polarization 𝑷 as

𝑷 𝑹 ∥ , 𝑡 = 𝜒 𝑹 ∥ 𝑬 𝑹 ∥ , 𝑡 = 𝜒 𝑹 ∥ 𝐸 0 𝑹 ∥ , 𝑡 𝑒 -𝑖 𝜔 𝑡 𝒚 (3)
where 𝜒 𝑹 ∥ , 𝑧 is the dielectric susceptibility assumed in our model to be time independent.

The susceptibility in the grating can be expanded in Fourier series [START_REF] Yang | Unified analytical theory of single-order soft x-ray multilayer gratings[END_REF] 𝜒 𝑹 ∥ = 𝑢 𝑚 𝑈 𝑚 ,𝑛 𝑒

with 𝜖the dielectric constant given in the Gauss unit system, by

𝜖 𝑹 ∥ = 1 + 4𝜋 𝜒 𝑹 ∥ (6)
and

𝑘 = 𝜔 𝑐 (7)
The envelope of the diffracted electric field can be represented by the Rayleigh expansion [START_REF] Petit | Electromagnetic Theory of GratingsTopics in Current Physics[END_REF] 𝐸 0 𝑹 ∥ , 𝑡 = ℇ 𝑝 𝑧, 𝑡 𝑒𝑥𝑝 proved to be very efficient in this kind of optical problem [START_REF] Yeh | Optical Waves in Layered Media[END_REF].In this context, the field can be written as the superposition of two waves propagating in opposite directions along the z-axis, so that we write, using the CWT ℇ 𝑝 𝑧, 𝑡 = 𝐹 𝑝 𝑧, 𝑡 𝑒 + 𝑖 𝜅 𝑝 𝑧 + 𝐵 𝑝 𝑧, 𝑡 𝑒 -𝑖 𝜅 𝑝 𝑧 ; 𝜅 𝑝 = 𝑘 2 𝜖 -𝑞 𝑝 2 (11a) together with the following requirementto ensure the uniqueness [START_REF] Yang | Unified analytical theory of single-order soft x-ray multilayer gratings[END_REF] 𝐹 𝑝 ′ 𝑧, 𝑡 𝑒 + 𝑖 𝜅 𝑝 𝑧 + 𝐵 𝑝 ′ 𝑧, 𝑡 𝑒 -𝑖 𝜅 𝑝 𝑧 = 0

(11b)
The term 𝜖 is the average value of the dielectric constant (term zero of the Fourier series)that means that 𝜅 𝑝 does not depend on the z variable in agreement with the essence of the CWT.

SubstitutingEqs. [START_REF] Bushuev | Diffraction of X-ray free-electron laser femtosecond pulses on single crystals in the Bragg and Laue geometry[END_REF] in Eqs. [START_REF] Ksenzov | Time-space transformation of femtosecond free-electron laser pulses by periodical multilayers[END_REF]

Only the zeroth order term 𝜖 in the Fourier expansion of 𝜖 𝑹 ∥ has been kept in the timedependent term of Eqs. [START_REF] Yang | Unified analytical theory of single-order soft x-ray multilayer gratings[END_REF]; so one can say that our calculation is a first-order (in terms of 𝜖) perturbative time-dependent model.

Two-wave theory

Matrix approach

Let us assume a strong couplingoccurringbetween the incident wave with amplitude𝐹 0 𝑧, 𝑡 and the p th wave diffracted by the grating with amplitude 𝐵 𝑝 𝑧, 𝑡 . In this condition it is possible to keep in the system given by Eqs.( 12) only the two terms 𝐹 0 𝑧, 𝑡 and 𝐵 𝑝 𝑧, 𝑡 .This is a situation that corresponds to the so-called two-wave theory (TWT) in the dynamical theory of diffraction. The conditions of validity of this approximation are discussed in several papers [START_REF] Yang | Unified analytical theory of single-order soft x-ray multilayer gratings[END_REF][START_REF] Francescangeli | On the limits of validity of the two-wave approximation in the dynamical theory of electromagnetic scattering by periodic dielectric media[END_REF] ; as mentioned in [START_REF] Yang | Unified analytical theory of single-order soft x-ray multilayer gratings[END_REF] a validity condition for this regime is that the angular width of the diffracted peak is small compared to the distance in terms of glancing angle between the neighbouring diffraction peaks, a condition generally required for a spectroscopic application of the gratings.

One finally gets a system of PDEs with terms that do not depend on the variable z,provided than one introduces the following auxiliary amplitude terms

𝑓 0 𝑧, 𝑡 = 𝐹 0 𝑧, 𝑡 exp [𝑖 (𝜅 0 + 𝜅 𝑝 ) - 𝑗 𝐺 𝑧 2 𝑧 ] (13a)
and

𝑏 𝑝 𝑧, 𝑡 = 𝐵 𝑝 𝑧, 𝑡 exp [-𝑖 (𝜅 0 + 𝜅 𝑝 ) - 𝑗 𝐺 𝑧 2 𝑧 ] (13b)
Moreover, if the system is in the vicinity of the jthBragg resonance 𝑗 𝐺 𝑧 ≈ 𝜅 0 + 𝜅 𝑝 [START_REF] Yeh | Optical Waves in Layered Media[END_REF] and 𝑓 0 𝑧, 𝑡 ≈ 𝐹 0 𝑧, 𝑡 , 𝑏 𝑝 𝑧, 𝑡 ≈ 𝐵 𝑝 𝑧, 𝑡 . Eq.( 14) gives the generalized Bragg condition.Combining Eqs. [START_REF] Yang | Unified analytical theory of single-order soft x-ray multilayer gratings[END_REF][START_REF] Petit | Electromagnetic Theory of GratingsTopics in Current Physics[END_REF],there results after some algebra the following system of time-dependent coupled PDEs satisfied by the varying amplitudes 𝑓 0 𝑧, 𝑡 and 𝑏 𝑝 𝑧, 𝑡 in the domain of the p th diffraction; using the matrix formalism, this system reads

𝜕 𝜕𝑧 Σ B 𝑧, 𝑡 = ℬ 𝐵𝑝 𝜕 𝜕𝑡 Σ B 𝑧, 𝑡 + 𝑖 ℳ 𝐵𝑝 Σ B 𝑧, 𝑡 (15)
where Σ B 𝑧, 𝑡 is the following column amplitude vector

Σ B 𝑧, 𝑡 = 𝑓 0 𝑧, 𝑡 𝑏 𝑝 𝑧, 𝑡 (16)
The subscript B stands for the Bragg-case reflection geometry. In Eq.( 15),ℳ 𝐵𝑝 is the propagation matrix in space

ℳ 𝐵𝑝 = -𝛼 𝐵 𝐾 𝐵 + 𝐾 𝐵 - 𝛼 𝐵 (17a)
𝛼 𝐵 = 𝜅 - 𝑗 𝐺 𝑧 2 - 𝑘 2 2𝜅 0 Γ 𝜒 ; 𝜅 = (𝜅 0 + 𝜅 𝑝 ) 2 (17b) 𝛼 𝐵 = 𝜅 - 𝑗 𝐺 𝑧 2 - 𝑘 2 2 𝜅 𝑝 Γ 𝜒 (17c)
𝐾 𝐵 + = - 𝑘 2 2 𝜅 0 𝑢 𝑗 𝑈 𝑗 ,-𝑝 ; 𝐾 𝐵 -= + 𝑘 2 2 𝜅 𝑝 𝑢 -𝑗 𝑈 -𝑗 ,𝑝 (17d)
ℬ 𝐵𝑝 is the propagation matrix in time

ℬ 𝐵𝑝 = - 𝑘 𝑐𝜅 0 0 0 + 𝑘 𝑐 𝜅 𝑝 (18)
We first consider the time-independent case

𝜕 𝜕𝑧 Σ B 𝑧 = 𝑖 ℳ 𝐵𝑝 Σ B 𝑧 (19)
The solution can be obtained by substituting

Σ B 𝑧 = 𝐴 𝐵 𝑒 𝑖 𝜓 𝑧 (20)
andasshown in Appendix 1, Some calculations show that 𝑆 𝐵𝑝 𝑧 reduces to the matrix 𝑆 𝑧 given by Eq.( 17)ofRef. [START_REF] André | Time-dependent Bragg diffraction and shortpulse reflection by one-dimensional photonic crystalsJ[END_REF]for the case of a 1D-PC. From Eqs. [START_REF] Sammar | Diffraction and scattering by lamellar amplitude multilayer gratings in the X-UV region[END_REF][START_REF] Afanas | Dynamical theory of X-ray diffraction in crystals with defects[END_REF](23)(24), which form the basis of the time-independent two-wave CW analysis, it is possible to calculate the reflection and transmission diffraction efficiencies of a grating at a given diffraction orderand the electric field distribution within the grating.

Σ B 𝑧 = 𝑆 𝐵𝑝 𝑧 Σ B 0 (21)
As for the time-dependent case, one searches the solution by analogy with the timeindependent case in the following form:

Σ B 𝑧, 𝑡 = 𝐴(𝑡) 𝐵(𝑡) 𝑒 𝑖 𝜓 𝑧 (25)
This approach can be regarded as a kind ofLagrange's method of variationof constants:

Finally, by following a way similar to the one presented for the time-independent case, we obtain

Σ B 𝑧, 𝑡 = 𝑅 B𝑝 𝑧, 𝑡 Σ B 0,0 (29)
where

𝑅 B𝑝 𝑧, 𝑡 = exp -𝐺 B𝑝 𝑡 𝑆 B𝑝 𝑧 (30)
The propagation of the electric field in time and space (field distribution within the grating, transmission, reflection, …) can be deduced from Eqs.(28-30). In the following, we consider the indicial responsein terms of reflection and transmissionunder Heaviside unit-step input. Some numerical examples will be given in Section 4.

Indicial response in the two-wave approximation

The reflection and transmission coefficientsare derived from the initial and boundary conditions: at z = 0, a Heaviside unit-step Θ(𝑡) is applied, so that 𝑓 0 0, 𝑡 = Θ(𝑡), and at z = L there is no incoming wave, so that𝑏 𝑝 𝐿, 𝑡 = 0, which gives from Eq.( 29)

𝑓 0 𝐿, 𝑡 𝑏 𝑝 𝐿, 𝑡 = 0 = 𝑅 B𝑝11 𝐿, 𝑡 𝑅 B𝑝12 𝐿, 𝑡 𝑅 B𝑝21 𝐿, 𝑡 𝑅 B𝑝22 𝐿, 𝑡 Θ(0 + ) 𝑏 𝑝 (0, 0 + ) (31)
where 𝑅 B𝑝𝑖𝑗 stand for the coefficients of the matrix 𝑅 B𝑝 𝑧, 𝑡 .Then, the calculation can becarried out as in Ref. [START_REF] André | Time-dependent Bragg diffraction and shortpulse reflection by one-dimensional photonic crystalsJ[END_REF]to give

𝑏 𝑝 0, 𝑡 = 𝑅 B𝑝21 0, 𝑡 -𝑅 B𝑝22 0, 𝑡 𝑅 B𝑝21 𝐿, 𝑡 𝑅 B𝑝22 𝐿, 𝑡 (32)
Consequently using the definition of the reflection coefficient, one finds the indicial response 𝑅 𝛩 in terms of reflectivity at the time t after switching on abruptly a constant intensity source at t = 0, or in other words, when applying a Heaviside unit-step inputΘ;

thus, it resultsfor the indicial response in terms of reflectance

𝑅 𝛩 t ≡ 𝑏 𝑝 0, 𝑡 0 t 2 (33)
From Eq.(31) one also finds that

𝑓 0 𝐿, 𝑡 = 𝐷𝑒𝑡[𝑅 B𝑝 𝐿, 𝑡 ] 𝑅 B𝑝22 𝐿, 𝑡 (34)
The indicial response in terms of transmittance is

𝑇 𝛩 t ≡ 𝑓 0 𝐿, 𝑡 0 t 2 (35)

3.3Tagaki-Taupin approach

We now consider the particular case where 𝛼 𝐵 ≈ 𝛼 𝐵 ; this situation occurs when the grating is used in specular conditionp = 0(in this case 𝛼 𝐵 and 𝛼 𝐵 are strickly equal) or when the period of the grating is large with respect to the wavelength of the diffracted radiation so that 𝜅 𝑝 ≈ 𝜅 0 and 𝜃 𝑝 ≈ 𝜃 0 . In this condition, it is useful to follow an approach adopted by Tagaki and Taupin and now usually implemented in the dynamical theory of x-ray diffraction [START_REF] Authier | Dynamical Theory of X-Ray Diffraction[END_REF],since it leads to a formulation giving aphysical insight into the problem, as we show hereafter.As shown in Appendix 2, the problem can be reduced to the following hyperbolic second order PDE ℒ 𝑣,𝑤 𝑓 0 𝑣, 𝑤 ; 𝑏 𝑝 𝑣, 𝑤 = 0

(36)
where ℒ 𝑣,𝑤 is the differential operator defined by

ℒ 𝑣,𝑤 = 𝜕 2 𝜕𝑣𝜕𝑤 + 𝜋 2 Λ 2 (37)
with Λ 2 , the quantity related to the coupling constants 𝐾 𝐵 + , 𝐾 𝐵 -by

Λ 2 = - 𝜋 2 𝐾 𝐵 + 𝐾 𝐵 - (38)
𝑣, 𝑤are the following characteristic coordinates

𝑣 = 1 2 𝑐 𝑡 𝑠𝑖𝑛𝜃 -𝑧 (39a) 𝑤 = 1 2 𝑐 𝑡 𝑠𝑖𝑛𝜃 + 𝑧 (39b)
and 𝑓 0 𝑣, 𝑤 , 𝑏 𝑝 𝑣, 𝑤 the reduced field amplitudes defined according to

𝑓 0 𝑣, 𝑤 = 𝑒𝑥𝑝 -𝑖 𝒶 𝑐 𝑡 𝑠𝑖𝑛𝜃 𝑓 0 𝑣, 𝑤 (40a)
𝑏 𝑝 𝑣, 𝑤 = 𝑒𝑥𝑝 -𝑖 𝒶 𝑐 𝑡 𝑠𝑖𝑛𝜃 𝑏 𝑝 𝑣, 𝑤

where 𝒶 ≡ 𝛼 𝐵 . Provided that the media are not absorbing, using Eqs. [START_REF] Courant | Partial Differential Equations[END_REF], it follows that

Λ = 2 𝜋 𝑠𝑖𝑛𝜃 𝑘 𝑢 𝑗 𝑈 𝑗 ,-𝑝 (41)
One recognizes that the quantity Λ is the extinction length of the dynamical theory of diffraction [START_REF] Authier | Dynamical Theory of X-Ray Diffraction[END_REF].The PDE given by Eq.(36) can be solved for given boundary conditions by implementing different methods; the most common one is the Riemann's method that we summarize hereafter.For the sake of consistency, we present in Appendix 3,a brief mathematical development of the Riemann's method applied to our problem; the rigorous mathematical foundations of the method can be found in Ref. [START_REF] Courant | Partial Differential Equations[END_REF]and some details of the calculationare available in Ref. [START_REF] Lindberg | Time dependence of Bragg forward scattering and self-seeding of hard x-ray free-electron lasers[END_REF]. The application of Riemann's methodrequires an integration contourin the characteristic coordinate plane shown in Figure 2. As shown in Appendix 3,the diffracted field in reflection geometry can be written as It is interesting to note that the Tagaki-Taupin approach shows that the indicial response in terms of reflection is conditioned by the extinction length and that it presents a transient period whose duration is given by a characteristic transient time𝑡 𝑐 approximately equal to 2 units of reduced time as shown in figure 3. One is led to think that even where the specular approximation is not valid, the indicial response of the grating still presents a transient period of the order of 𝑡 𝑐 .

𝑏 𝑝 𝑣, 𝑣 = -𝑖 𝐾 𝐵 -Λ 𝜋 𝐽 1 2 𝜋 Λ 𝑣 -𝑣 ′ 𝑣 -𝑣 ′ 𝑓 0 𝑣 ′ , 𝑣 ′ 𝑄 𝑃 𝑑𝑣 ′ (
As also shown in Appendix 3,the diffracted field in transmission geometry is approximately given by :

𝑓 0 𝑣, 𝑤 𝑇 ≈ - 2 𝜋 2 𝐿 Λ 2 𝑓 0 𝑣 ′ , 𝑣 ′ 𝐽 1 2 𝜋 Λ 𝑣 + 𝐿 -𝑣 ′ 𝑣 -𝑣 ′ 2 𝜋 Λ 𝑣 + 𝐿 -𝑣 ′ 𝑣 -𝑣 ′ 𝑄 𝑃 𝑑𝑣 ′ (48)
To determine the impulse response in terms of transmission, we insert in Eq.( 48) the expression of the incident pulse given by Eq.( 44) as for the reflection caseand we performthe integration as for the reflection geometry; it resultsthat the impulse response in terms of transmission 𝑓 0 𝑧 = 𝐿, 𝑇 δ is Strejc [START_REF] Strejc | [END_REF][START_REF] Larminat | Analysis and Control of Linear Systems[END_REF]. The temporal dependence of the reflection or transmission of a short pulse can, then, be determined using Eq.( 51) and then performing an inverse Laplace transform. Details can be found in Ref. [START_REF] André | Time-dependent Bragg diffraction and shortpulse reflection by one-dimensional photonic crystalsJ[END_REF].

𝑓 0 𝑧 = 𝐿, 𝑇 δ ≈ 𝑔 𝑇 𝑇 = 𝜋 2 𝐿 𝑠𝑖𝑛𝜃 Λ 2 𝐽 1 𝜉 𝑇 𝜉 𝑇 Θ 𝑇 ; 𝜉 𝑇 = 𝜋 Λ 𝑐 𝑇 2 𝐿 𝑠𝑖𝑛𝜃 + 𝑐 𝑇 𝑠𝑖𝑛𝜃 2

Numerical simulations

We present a numerical example obtained with the time-dependent matrix approach in the continuation of the one presented in Ref. [START_REF] André | Time-dependent Bragg diffraction and shortpulse reflection by one-dimensional photonic crystalsJ[END_REF]. We consider a lamellar unslanted (𝜙 = 0) grating formed from multilayer bars. The multilayer structure consists of a periodic stack of N= 20Fe/C bilayers; the period d is equal to 5.0 nm and the ratio is equal to 0.5, that is to say the thicknesses of the Fe and C layers are the same. Figure 4 displays this structure and the parameters.

Figure 4: Scheme of the unslanted lamellar grating considered in the numerical applications. D is the period of the grating andD the width of the multilayer bars.

The energy of the incident radiation is 8 keV. In all calculations, we use for the optical indices the values tabulated in the CXRO database [START_REF]CXRO X-Ray Interactions With Matter[END_REF]. One assumes that the grating period D equals 1 m.Figures5 and 6show the steady state reflectance and transmittancerespectively, for the diffraction order p = -1 at different values of the  parameter: 0.1, 0.2, 0.4 and 0.6; they arecalculated by the TWT presented in section 3.1Let us emphasize that in the Bragg domain, the results of the TWTare in very good agreement with results given by a rigorous coupled wave analysis RCWA (see for instance Eq.( 12) of Ref. [START_REF] Yang | Unified analytical theory of single-order soft x-ray multilayer gratings[END_REF])or a modal theory (MT) [START_REF] Sammar | Diffraction and scattering by lamellar amplitude multilayer gratings in the X-UV region[END_REF] as illustrated by

Conclusion

We have generalized the time-dependent coupled-wave theory initially developed for one-dimensional photonic crystals [START_REF] André | Time-dependent Bragg diffraction and shortpulse reflection by one-dimensional photonic crystalsJ[END_REF] to the spatiotemporal diffraction by multilayer gratings. The results obtained with the matrix formalism are in agreement with the Tagaki-Taupin theory originally developed in the framework of the dynamical theory of crystal diffraction that we have extended to grating diffraction. From Tagaki-Taupin theory, it appears that in reflection geometry, the key quantity is the extinction length while in transmission geometry, the thickness of the grating is also an important parameter. Although a very interesting means for obtaining a physical insight of the problem, Tagaki-Taupin approach is not relevant for predicting accurately the nonspecular case behaviour; instead a matrix approach is very efficient.

This theoretical work gives a useful tool to predict the temporal response of optics implemented with short sources such as free electron lasers, high harmonic generation sources, …Numericalapplications of the present theory to gratings other than the lamellar one will be presented in a forthcoming paper. Our time-dependent matrix approach should be applicable beyond the two-wave approximation that is in the framework of the rigorous coupled-wave analysis; this work is underway. Finally this theory has to be extended to p-polarizationand conical diffraction that requires going beyond the paraxial approximation.

(A3.1)

To treat the Bragg-case geometry, we consider the contour 𝒞 formed by the triangle (PQR) formed in the plane of characteristic coordinates 𝑣, 𝑤as shown in Figure 2

Figure 1 :

 1 Figure 1: Sketch of the grating of period D; each grating pitch has a width equal to Γ𝐷 and is formed by a periodic stack of N bilayers with thickness d; the bilayer is made up of a material a of dielectric susceptibilityχ 𝑎 and material bof dielectric susceptibilityχ 𝑏 , with layer thickness da = d and db = (1-)d respectively. The incoming radiation strikes the

  𝐴 → 𝐴 𝑡 , 𝐵 → 𝐵(𝑡).Inserting Eq.(25) in Eq.(15) gives after derivation with respect to space 𝜕 𝜕𝑡 𝐴(𝑡) 𝐵(𝑡) = -𝐺 B𝑝 𝐴(𝑡) 𝐵(𝑡) (26) with 𝐺 B𝑝 = 𝑖 ℬ B𝑝 -1 𝜓 𝐼 -ℳ B𝑝 (27) Integration of Eq.(25)with respect to time gives 𝐴(𝑡) 𝐵(𝑡) = exp -𝐺 B𝑝 𝑡 𝐴(0) 𝐵(0)

Figure 2 :

 2 Figure 2: Grating geometry in the characteristic coordinate reference frame (v,w); the grating front surfaceis given by w = v(line PQ) while the rear surface is given by w = v + L (line RT).

  42) where𝑓 0 𝑣 ′ , 𝑣 ′ corresponds to the incoming wave at the front surface (z =0). 𝐽 𝑛 is a Bessel function of the first kind.From Eq.(42) one can deduce the impulse and indicial response in terms of reflection coefficient.The impulseincident reduced field amplitude 𝑓 0 𝑣, 𝑣 δ reads 𝑓 0 𝑣, 𝑣 δ = 𝑒𝑥𝑝 + 𝛿stands for the Dirac peak. Inserting Eq.(43) in Eq.(42) and performing the integration gives for the reduced diffracted field𝑏 𝑝 𝑣, 𝑣 δ under the incidence of a Dirac pulse 𝑏 𝑝 𝑣, 𝑣 δ = -𝑖 𝑠𝑖𝑛𝜃 𝐾 𝐵 -that is, for the diffracted field 𝑏 𝑝 𝑧 = 0, 𝑇 δ = 𝑖 𝑠𝑖𝑛𝜃 𝐾 𝐵 introduced the time delay T measured with respect to the diffracted wave plane 𝑇 = 𝑐 𝑡 -𝑥 𝑐𝑜𝑠 𝜃 𝑐 (46) The quantity𝑔 𝑅 𝑇 = 𝑏 𝑝 𝑧 = 0, 𝑇 δ is the temporal Green function. For time coherent radiation with time-dependent causal distribution Ξ (normalized to unity), the indicial response 𝑅 𝛩 𝑡 in terms of reflection coefficient is given by 𝑅 𝛩 𝑡 = 𝑔 𝑅 (𝑇) 47) allows one to draw a "universal" curve for the indicial response 𝑅 𝛩 𝑡 in terms of reduced time 𝑡 = 𝑡 𝜋 sin 𝜃 𝑐 Λ , see Figure 3.

Figure 3 :

 3 Figure 3: Indicial response in terms of reflectance versus reduced time for a grating working in Bragg geometry close to the specular condition. The response is normalized to unity.

3 . 4

 34 radiation with time-dependent causal distribution Ξ (normalized to unity), the indicial response 𝑇 𝛩 𝑡 in terms of transmission coefficient is given by 𝑇 𝛩 𝑡 = 𝑔 𝑇 (𝑇) From indicial response to time-dependent reflection and transmission of a short pulse Let 𝐸 𝑡 be the temporal envelope shape of any incident pulse and S(t) the envelopeof the time-dependent reflected or transmitted pulse. Then the Laplace transform 𝑆 𝑠 of S(t) is related to the Laplace transform 𝐸 𝑠 of 𝐸 𝑡 by means of the convolution theorem 𝑆 𝑠 = Ζ 𝛿 𝑠 𝐸 𝑠 (51) whereΖ 𝛿 𝑠 is the transfer function in terms of reflectance or transmittance that is the Laplace transform of the impulse response Ζ 𝛿 𝑡 . This impulse response corresponds to the instantaneous reflectance or transmittanceobtained when the system is struck by a Dirac pulse at t=0 and should not be confused with the indicial response Ζ Θ 𝑡 . Nevertheless Ζ 𝛿 𝑡 and consequently Ζ 𝛿 𝑠 can bedetermined from Ζ Θ 𝑡 by different methods; we have applied a basic but efficient graphic method introduced by

Figure 7 .

 7 In this figure, as an example of the accuracy of the TWT in the Bragg domain, the steadystate reflectance for the unslanted grating with a value of 0.2is presented: the data D GD from TWT are given by the solid line and the data from RCWA and MT (which are identical) by the dots; the calculations in the RCWA or MT are done with 15 Fourier terms(or modes). Similar agreements are observed for the other reflectances and transmittances displayed in figures 5 and 6.

Figure 5 :

 5 Figure 5: Steady-state reflectance versus glancing angle0of an unslanted lamellar multilayer grating shown in figure 4 for different values of the parameter . The photon energy is 8 keV and the parameters of the grating are given in the text: =0.6 green thick solid line, =0.4 black dotted line, =0.2 blue dashed line, =0.1 red thicksolid line. The diffraction order p = -1.

Figure 7 :

 7 Figure 7: Steady-state reflectance versus glancing angle 0of the unslanted lamellar multilayer grating shown in figure 4 for =0.2 calculated from TWT (solid line) and rigorous theories RCWA and MT (dots). The diffraction order p = -1.

Figure 8 Figure 9

 89 Figure 8 showsthe indicial response of the peak reflectance𝑅 𝛩 𝐵𝑟𝑎𝑔𝑔 for the grating diffraction order p = -1, that is the reflectance at the generalized Bragg angleBragg ,for different values of the  parameter: 0.1, 0.2, 0.4 and 0.6. It appears that the response presents a transient period in agreement with the Tagaki-Taupin approach; the transient time is approximately equal to two units of reduced time and depends on the value of the parameter: the smaller the  the  (Eq.(42)) since 𝑈 𝑗 ,-𝑝 is proportional to (Eq.(4 f))  consequently, the longer the transient period as indicated by the Tagaki-Taupin theory.Let us outline that the generalized Bragg angle varies slightly with the value of the parameter according to Eq.(14). After the transient period the peak reflectance reaches the steadystate valueRsp which also depends on the value of  Figure 9 shows the indicial response 𝑇 𝛩 𝐵𝑟𝑎𝑔𝑔 for the grating diffraction order p = -1, now in terms of transmittance at the generalized Bragg angle Bragg , for different values of the  parameter: 0.1, 0.2, 0.4 and 0.6. The indicial response in terms of transmittance presents also a transient period but contrary to the reflectance case, the transient time depends not only on the   but also on the thickness of the grating.

Figure 8 :

 8 Figure 8: Indicial response in terms of peak reflectance 𝑅 𝛩 𝐵𝑟𝑎𝑔𝑔 of an unslanted lamellar multilayer grating shown in figure4for different values of the parameter . The photon energy is 8 keV, the grating diffraction order p = -1 and the parameters of the grating are given in the text: =0.6green thick solid line,=0.4black dotted line, =0.2blue dashed line, =0.1 red thin solid line.

Figure 9 :

 9 Figure 9: Indicial response in terms of transmittance𝑇 𝛩 𝐵𝑟𝑎𝑔𝑔 of an unslanted lamellar multilayer grating shown in figure 4 for different values of the parameter  at the Bragg angle. The photon energy is 8 keV, the grating diffraction order p = -1 and the parameters of the grating are given in the text: =0.6 green thick solid line, =0.4 black dotted line, =0.2 blue dashed line, =0.1 red thin solid line.

Figures 10 and 11 displaythe 3 -

 3 Figures 10 and 11 displaythe 3-dimensionaltemporal and spectral indicial response in terms of reflectance and transmittance for the unslanted lamellar grating described in figure 4 for =0.2. The spectral dependence versus the glancing angle at a given photon energy (8keV) is given; we choose this representation in order to avoid taking into account the energy dispersion of the optical indices. The figures are zoomed on the short values of time in order to focus on the transient period. It clearly appears that the transient time is shorter for transmittance than for reflectance; as mentioned previously the transient time in transmittance is mainly dependent on the thickness Lof the grating, in agreement with the Tagaki-Taupin theory (see Eq.(48)).

Figure 10 :

 10 Figure 10: Temporal and spectral indicial response in terms of reflectance of the lamellar grating shown in figure 4 for =0.2. The photon energy is 8 keV, the grating diffraction order p = -1.

Figure 11 :

 11 Figure 11: Temporal and spectral indicial response in terms of transmittance of the lamellar grating shown in figure 4 for =0.2. The photon energy is 8 keV, the grating diffraction order p = -1.

Figure 12

 12 Figure 12 displays the reflected peak height of an incident Gaussian pulse with unit amplitude and temporal width equal to 10 fs for the grating considered in figure 4. The calculations are performed according to the method given in section 3.4. The influence of the  parameter is clearly evident, particularly concerning the stretching of the reflected pulse. The calculations are carried out according to the method given in section 3.4.

Figure 12 :

 12 Figure 12: Time dependence of the reflected pulse height for an incident pulse of width equal to 10 fs, of an unslanted lamellar multilayer grating shown in figure 4 for different values of the parameter . The photon energy is 8 keV, the grating diffraction order p = -1 and the parameters of the grating are given in the text: =0.6 green thick solid line, =0.4 black dotted line, =0.2 blue dashed line, =0.1 red thin solid line.

5 )(A3. 6 )(A3. 8 ) 1 (A3. 9 ) 15 )Figure 2 ;(A3. 16 )

 5681915216 Figure2; the segment𝑅𝑇 corresponds mathematically to the condition v = w-Land physically to the rear surface of the grating. In a way similar to the reflection case, contour integration along ℬ leads to

(A3. 21 )

 21 along the rear surface 𝑤 = 𝑣 + 𝐿and

  Implementing these conditions Eqs.(A3.17-A3.19) in Eq.(A3.16) gives 𝐹 𝑣 ′ , 𝑤 ′ 𝑇 = 𝑅 𝑣, 𝑤; 𝑣 ′ , 𝑤 ′ 𝐹 𝑣 ′ , 𝑤 ′ 𝑄 -𝑅 𝑣, 𝑤; 𝑣 ′ , 𝑤 ′ 𝜕𝐹 𝑣 ′ , 𝑤 ′ 𝜕𝑣 ′Afanas'ev and Kohn[START_REF] Sammar | Diffraction and scattering by lamellar amplitude multilayer gratings in the X-UV region[END_REF]have given the Riemann 𝑅 𝑣, 𝑤; 𝑣 ′ , 𝑤 ′ function that fulfills the transmission conditions given by Eqs.(A3.17-A3.[START_REF] Larminat | Analysis and Control of Linear Systems[END_REF]

				′	𝑇 = 1
							(A3.19)
						𝑄
							𝑑𝑣 ′
					𝑃
	-	𝑄 𝑃	𝐹 𝑣 ′ , 𝑤 ′ 𝜕𝑅 𝑣, 𝑤; 𝑣 ′ , 𝑤 ′ 𝜕𝑣 ′	𝑤 ′ = 𝑣 ′ 𝑑𝑣 ′
							(A3.20)
	𝑅 𝑣, 𝑤; 𝑣 ′ , 𝑤 ′ = 𝐽 0	2 𝜋 Λ	𝑤 -𝑤 ′ 𝑣 -𝑣 ′ +	𝑣 -𝑣 ′ 𝑤 -𝑤 ′ 𝐽 2	2 𝜋 Λ	𝑤 -𝑤 ′ 𝑣 -𝑣 ′

APPENDIX 1

Inserting Eq. [START_REF]CXRO X-Ray Interactions With Matter[END_REF] in Eq. [START_REF] Francescangeli | On the limits of validity of the two-wave approximation in the dynamical theory of electromagnetic scattering by periodic dielectric media[END_REF] To the two eigenvalues 𝜓 ± are associated two eigenvectors

The solutions for Σ B 𝑧 can be derived using the eigenmatrix

that is

Since the general solution is a linear combination of the eigensolutions Σ B 𝑧 = 𝐶 + 𝑒 + 𝑖 𝜓 + 𝑧 𝑉 + + 𝐶 -𝑒 𝑖 𝜓 -𝑧 𝑉 -= 𝑃 𝑒 + 𝑖 𝜓 + 𝑧 0 0 𝑒 𝑖 𝜓 -𝑧 𝐶 + 𝐶 - (A1.6)

At z = 0, one has

Putting Eq.(A1.7) in Eq (A1.5), it follows that

where 𝑆 𝐵𝑝 𝑧 is obtained by the product 𝑃 𝑒 + 𝑖 𝜓 + 𝑧 0 0 𝑒 𝑖 𝜓 -𝑧 𝑃 -1 .

APPENDIX 2

With 𝛼 𝐵 ≈ 𝛼 𝐵 ≡ 𝒶 and 𝜃 𝑝 ≈ 𝜃 0 ≡ 𝜃, the system given by Eq.( 15) can be written


Professor, Ana Maria Cruz

Eric Piatyszek

Junlei Yu

Akihiko Hokugo

Carine El Hajj

Michel Lesbats

Professor, Alicja Tardy

The ResTO-TerRiN Project: Contribution to the Systemic Modeling of Technical and Organizational Issues of a Territory Exposed to Natech Risk

This paper presents a brief summary of work being carried out under the Franco-Japanese research project ResTO TerRiN. The main goal of the project is to produce relevant knowledge and effective methods and tools to improve the resilience of a territory against Natech accidents (chemical accidents triggered by natural hazards) especially those due to flood / tsunami. The work is based on a posteriori (in Japan and France) and a priori analysis of the industrial as well as the local governments' emergency management to the Natech accidents during the Great East Japan earthquake and tsunami (GEJET) and during severe flooding events in France. Survey questionnaires in Japan and France are being applied to chemical facilities, to government agencies, and to citizens in France and Japan The data collected are used to model the impact of the natural hazard events on the facilities and the safety barriers, as well as the community and overall social impacts. These results are then used to understand societal and territorial resilience to these complex disasters and propose a Natech resilience model. In this paper we present the preliminary results of the ongoing work in Japan.

INTRODUCTION

Chemical accidents triggered by natural hazards (known as Natechs) are relatively rare. In fact they represent between 2 -7 % of all chemical accidents reported in databases in Europe and the United States, and show an increasing trend in the last 30 years [START_REF] Sengul | Analysis of hazardous material releases due to natural hazards in the United States[END_REF][START_REF] Santella | Empirical estimation of the conditional probability of Natech events within the United States[END_REF]). Growing urban populations and industrialization in areas subject to natural hazards, coupled with environmental change due to climate change and other factors will contribute to an increase in the number of Natech accidents in the years to come [START_REF] Cruz | Vulnerability of the Oil and Gas Sector to Climate Change and Extreme Weather Events. Special Issue on "Climate Change, Extremes, and Energy Systems[END_REF].

Natural hazards are generally taken into account in the design of chemical plants by insuring that buildings and structures are built to codes. Nonetheless, they are mostly absent in process hazard assessment and chemical accident prevention programs. Thus, the safety and mitigation measures, and emergency response plans to prevent and respond to chemical accidents concurrent with natural disasters may not be effective under the natural disaster conditions. Furthermore, local government and civil protection authorities are often not aware of or not familiar with the dangers posed by certain industrial activities lacking knowledge and/ or training about chemicals and their potential health effects [START_REF] Cruz | Emerging Issues for Natech Disaster Risk Management in Europe[END_REF].

Natech risk management is relatively complex. These technological accidents occur while the "system which should provide resilience" are already under "stress" due to the concurrent natural hazard event. In the case of a Natech the industrial plant is both impacted by the natural hazard event, and may become the source of additional "hazards" due to the release/ potential release of hazardous materials which may have direct and indirect impacts offsite to nearby residents and communities.

This project studies these aspects by considering two levels of analysis: the scale of the industrial site, called microscopic scale, and the scale of the territory or society, called macroscopic scale. These two spatial scales look at the problem from the perspective of "vulnerability" and via a comprehensive analysis of the notion of territorial or societal resilience to Natech accidents. With regard to the microscopic scale, we aim to improve the understanding of the physical and organizational vulnerability of industrial facilities exposed to natural hazards (particularly tsunamis and floods), to identify and define prevention and protection measures, to ensure that plant owners/operators and emergency organizations are able to cope with Natech events. The idea is to provide to industrial risk managers tools and guides for improved Natech risk management during both prevention and response. With regard to the macroscopic scale, the objective of this project is to build a Natech resilience model through the identification of Natech resilience criteria. In this paper we present the preliminary results of research work being carried out in Japan.

METHODOLOGY

The project consortium involves a multi-disciplinary team of researchers from Japan and France. In both countries data has been collected through interviews, field visits, survey questionnairesin case study areas: in the Tohoku area (Japan) and in several industrial areas in France, particularly the Presqu'île d'Ambès (Gironde). Furthermore, we review reports, documents and the scientific literature. This peninsula, in the north of Bordeaux, is located between two rivers (the Garonne and the Dordogne) not far from the Atlantic coast. Therefore, this territory can be inundated by both river floods and storm surges. There are several small and large facilities (11 SEVESO sites and ICPE plants) along the Garonne banks. The overall population of the peninsula, which is exposed to both technological and natural risks, is about 30.000. Thus, this territory is appropriate for the analysis of the flood natech events, although there have been no natech accidents to date (several SEVESO sites were inundated by the 1999 storm, and some technological accidents have taken place, but never a technological accident has been triggered by flood events).

In Japan, the main objective is to collect and analyze data and extract lessons learned from Natech events triggered by the Great East Japan earthquake and tsunami disaster of March 2011. Through field visits and person-to person interviews with industrial and government officials it has been possible to construct event trees and fault trees of the most important accidents triggered by the earthquake and tsunami. Furthermore, through mail surveys of industry and local residents we are assessing Natech disaster risk reduction practices, risk perception and disaster preparedness for these types of events. In this paper we present preliminary results of two mail surveys. One to industrial facilities in affected areas in Miyagi Prefecture, and the other to residents living near a refinery in Sendai that suffered fires after the tsunami.

An industrial mail survey was sent out through the Japan Industrial and Medical Gases Association (JIMGA) in Sendai to member companies. JIMGA had the contact information of these companies and was able to contact their managers directly even if the facilities were not in operation. We present here a summary of the responses of 21 companies. The industrial survey asked questions regarding damages and losses, performance of safety and mitigation measures, and the opertional readiness/ effectiveness of the emergency response. The resident/ household survey asked questions regarding risk perception of Natech hazards, disaster preparedness and the evacuation for the tsunami as well as the Natech accidents.

The household survey was led by research team members at Kobe University to residents living within an area of 2km in radius from the source of a Natech accident that occurred at the JX Refinery following the tsunami. The purpose of the survey was to understand residents' attitude and protective actions such as evacuation due to the Natech. The questionnaire also inquired about damages to homes or injuries to residents, and the evacuation order issuing process.

In France, a priori analysis of the industrials as well as the local governments' (e.g., city and prefecture governments, fire division/ departments, industrial emergency responders) involved in the management of the Natech disasters during severe flooding events or storm surges is led. By the way of field visits, survey questionnaires, and interviews, we identify and map all players/actors and their responsibility, mechanisms, public policies and tools which are used in planning for natural and technological emergencies. The aim is to represent the organizational and management frameworks to identify criteria which could characterize the global natech resilience of the territory. Moreover, we are capitalizing the knowledge of the industrials to improve our understanding of how they manage natural disasters, technological accidents and eventually natech events. These data are used to construct various models (systemic model, functional representation, statistical models and risk models). Subsequently, these models are used to elaborate the guidelines/tools, identify the resilience criteria. In this paper, we present the first elements of a model built to assess the flood-induced Natech vulnerability of an industrial site.

PRELIMINARY RESULTS

Survey of Industrial Facilities

The interviews, field visits and mail survey to industrial facilities have provided important information concerning damage and losses, chemical accidents and spills triggered by the earthquake and tsunami, the vulnerability of the facilities, the accident consequences, the performance of safety and mitigation measures, and the overall emergency response to contain them. The interviews and visits to affected companies have provided detailed information about the accident mechanisms and cascading events, and has allowed the construction of flow diagrams, event trees and fault trees. Figure 1 presents an event tree of the impact of the earthquake and tsunami at a refinery and other nearby facilities at an industrial park in Sendai. It is important to note that in Figure 1, the cascade of events, direct and indirect, that contributed to the triggering of releases and/ or spills, and escalation into a large fire. Some of the steps in the chain (e.g., oil spill spread by flood waters) of cascading events would not normally occur during day-to-day plant operation, and therefore may be overlooked if the potential for natural hazard loads are not considered in the risk assessment of chemical process plants. Detailed fault trees have been constructed to identify common cause failures and key line of events.

The mail surveys provided additional important information about impacted facilities. Table 1 presents a summary of the characteristics of the industrial facilities surveyed thus far. The majority of facilities are nationally owned, small and medium size, chemical companies with 50 or less employees. 40% of the companies were subject to strong ground shaking (Shindo scale: 6 plus), 30% were subject to ground shaking intensity of 5 or higher, and the other 30 % did not know. Six companies said their plants were subject to tsunami inundation heights of 3.0 meters or higher, with one of these reporting 9.5 m. See Table 2.

Fig. 1 Event tree showing the impact of the Great East Japan earthquake and tsunami at a refinery and an industrial park in Sendai.

All of the facilities reported some kind of damage. 76% of the facilities said they were directly damaged by the earthquake while 52% said the facilities were directly damaged (e.g., direct damage due to strong ground shaking; direct impact by tsunami wave and damage due to hydrodynamic and hydrostatic loads) by the tsunami. Indirect damages (e.g., collapse of nearby building on equipment; impact by floating debris) from the earthquake and the tsunami were also reported. Table 3 shows the statistics concerning overall damage at the facilities surveyed. It is important to note that the region experienced several large aftershocks causing additional direct and indirect damages at more than half of the facilities.

Furthermore, all of the facilities reported that they had to partially or completely shut down after the disaster. Complete shutdown of operations occurred at 18 facilities (86%) and 3 reported partial shutdown. Table 4 shows the causes of total or partial shutdown. 81% of facilities indicated 'loss of electricity' as one of the contributing factors for the total or partial shutdown, followed by 76% of facilities that indicated that direct / indirect damage contributed to the shutdown.

Table 1. Industrial facility characteristics

Despite the fact that 17 facilities suffered total or partial damage, only five reported the release of hazardous materials. Two of the releases occurred due to the earthquake, one due to the earthquake and tsunami, and two due to the tsunami. Two of the releases occurred sometime after the earthquake, while the tsunami-triggered releases occurred during the tsunami. Earthquake damage of storage tanks, pipelines and other equipment occurred due to structural failure of support structures, overturning of equipment and buckling, among other causes. Tsunami damage occurred due to inundation and water intrusion, debris impact and floating-off of equipment. Table 5 summarizes these findings.

!

Table 6 summarizes damages to emergency response systems (including buildings, equipment and resources) by the earthquake and the tsunami. Table 5. Failure modes leading to chemical releases Damage to emergency response systems included direct damage to fire stations and fire trucks onsite, inundation of low lying emergency response equipment such as pumps and motors, floating off of emergency water pipelines, loss of utilities at emergency response/ services buildings, damage to IT equipment and damage to communication systems. Most facilities completely evacuated the plant sites before the arrival of the tsunami. Two facilities indicated that they had staff who stayed onsite. During the person-to-person interviews, several staff indicated that in the future, facilities in flood prone/ tsunami prone areas need to maintain emergency boats to move around once the danger has passed to rescue stranded workers or residents, or perform emergency response actions. Table 6. Number of facilities reporting damage to emergency response systems

Preliminary Findings of the Household Survey

The mail survey was sent to a random sample of 1,632 households, and 484 households responded, for a response rate of 29.4 %. Figure 2 presents the distribution of (a) sampled households and (b) returned questionnaires [START_REF] Yu | Households' Risk Perception and Behavioral Responses to Natech Accidents[END_REF].

We used risk likelihood and severity to measure the risk perception levels by asking: how likely did you think a Natech would be a threat to your lives or property before the Tohoku earthquake, just after the earthquake shaking, and in the future 10 years. Furthermore, we asked to what extent residents felt that a Natech would affect their lives or property before the Tohoku earthquake, just after the earthquake shaking, when they perceived the danger of the Natech, when they received the Natech evacuation order, while staying at the evacuation shelters, and in the future 10 years. One of the more interesting results of the survey is the fact that residents evacuated more than once. In fact 65% of respondents said they evacuated at least one time; over 30% evacuated the first time because of the Natech accident not the tsunami. Almost 50 % of the total that evacuated the first time, evacuated a second time, and over 66% of those that evacuated a second time, had to evacuate a third time. Some residents had to change shelter more than four times. Many residents evacuated on their own initiative without waiting for authorities to issue an evacuation order. During person-to-person interviews with residents, some indicated that there was a lack of clear information regarding where it was safe to evacuate to, or what they had to do in order to protect themselves and their families.

The survey asked about households' risk perception. The results show that people's risk perception in terms of the perceived likelihood that a Natech would cause harm to their lives or property increased (as expected) after experiencing the Natech during the Tohoku earthquake. However, as presented in Fig 3, no change was found for risk perception in terms of the perceived severity of the impacts of a Natech accident when compared with their responses before the Tohoku earthquake and in the next 10 years. This may be because no deaths, injuries or severe environmental damage due to the Natech was reported. We found that households felt that the Natech was very serious when: a. they perceived its occurrence, b. while they were staying at the first/second evacuation shelter, and c. when they received the Natech evacuation order.

Using logistic regression models, we also analyzed how risk perceptions together with geographic and demographic variables influenced people's evacuation decisions due to the Natech danger in order to comply with official orders.

!

Results showed that wind direction (B=2044, exp(B)=11.41, p≤0.01) and people's perception of the Natech severity when they perceived its occurrence motivated (B=1.77, exp(B)=5.89, p≤0.05) more people to evacuate. And we also found that being a female (B= 1.22, exp(B)=3.40, p≤0.01) made it more likely to comply with the Natech evacuation order. At last, to our surprise, results showed that people's perception about the Natech severity when they received the Natech evacuation orders played a negative role in evacuation motivation (B=-0.95, exp(B)=0.39, p≤0.05). This may due to the fact that many evacuees were reluctant to evacuate again for a Natech because they had already evacuated due to the earthquake/tsunami. Even those who have higher severity risk perception about the Natech when they received the Natech evacuation order were less likely to evacuate again. This finding may be important to emergency mangers or local officers to make a more effective plan to protect evacuees who may need to evacuate more than once.

First Elements of a Model to Assess Flood-Induced Natech Vulnerability of an Industrial Site (El Hajj et al., 2015)

The tools used to assess the vulnerability of industrial sites versus flood-induced Natech, will use final generic accidental scenarios (model) in the form of fault trees [START_REF] Vesely | Fault Tree Handbook[END_REF].

This model is based, on the one hand, on all the data collected during the project, and on the other hand, on a systematic risk analysis methodology. The systematic approach uses the MADS model, which is the French acronym for "Analysis Method of Dysfunctional Systems" [START_REF] Périlhon | La gestion des risques-Méthodes MADS MOSAR II[END_REF], (Lesbats et al. 2014). These models are built in three steps : a. Conceptual modeling of an industrial installation: This model is divided into several sub-systems responsible for triggering major accidents in case of flooding of an industrial site. The decomposition is an important step to facilitate the risk analysis: It is used to organize the information obtained through experience feedback and to conduct a rigorous analysis. b. MADS analysis. By adopting the MADS MOSAR methodology and building on the analysis of past accidents, each sub-system of the industrial installation will be analyzed in order to understand how it would be responsible for major accidents in case of floods.

Therefore, the identification of potential hazards and their process are established. c. Accident scenarios development: In order to identify potential accidents triggered by flood events, a fault tree is associated to each sub-system.

Building on the analysis of data collected during the project, only parts of an industrial installation which are involved, directly or indirectly, in the generation of major accidents during floods were considered (sub systems). We have identified 6 sub systems (fig. 5 In a second step, an analysis of each sub-system of the industrial installation was done in order to find how it can be impacted by the natural hazard (target of hazard flow) and how it could be considered, after a flood event, as a system causing mechanical shocks (ch), overpressure (Pr), heat radiation (th), pollution (P), electrical hazards (el) and regulation and energetic hazards (En) (source of hazard). For instance, for the critical equipment atmospheric storage tank, the following questions are asked systematically: how the atmospheric storage tank can cause mechanical shocks, overpressure, heat radiation…following a direct contact with floodwaters or, indirectly after the inundation of the site. Therefore, the main hazard flows between all the sub-systems that can be observed after a flood episode of an industrial installation and all the potential hazards in an industrial site were identified. The results of this systematic analysis were presented in a table containing the different elements of the MADS hazard flow.

In a third step, based on the various tables containing the different elements of the MADS hazard flow, fault tree are built for each sub system. For example Fig. 6 presents the fault tree of the sub system Electricity utility. Consequently, the final results of our work will be the development of generic accident scenarios triggered by flood events. These scenarios will be built using all the fault trees of the various sub systems which can be impacted during a flood event. They can be used for flood Natech risk analysis, a diagnosis tool based on these scenarios to find the weak points regarding technological accidents triggered by floods and finally recommendations on how to mitigate these accidents in the form of a list of preventive and protective measures.

CONCLUSIONS

These preliminary results have highlighted that natural hazard events can directly and/ or indirectly impact industrial facilities that handle hazardous materials resulting in spills or smaller releases, which may escalate to major accidents through cascading events including due to the loss of lifeline systems and emergency response capabilities.

The survey facilities suffered damages and losses, but less than 20 % actually experienced releases. It is not clear at this point whether this is due to their location (located in a less vulnerable area), good risk management practices, or other causes. This needs to be further studied.

Further research is continuing with a focus on how to improve the effectiveness of the evacuation orders by identifying when, where, and how to issue and disseminate information about chemical accidents, and what emergency response actions residents can take (particularly in the absence of warnings or evacuation orders from local authorities). Another key point is to understand the dynamics of multiple evacuations such as identifying which groups are more likely willing to evacuate more than once and what factors are more effective in motivating them.

Fig. 2

 2 Fig. 2 Distribution of (a) sampled households and (b) number of households completing the questionnaire (Source : Yu et al. 2015).

Fig. 3

 3 Fig. 3 Mean risk likelihood rating for natech (where 0 is not likely, and 5 is very likely) (Source: Yu et al. 2015).

Fig. 4

 4 Fig. 4 Mean risk severity ratings for Natech throughout the evacuation process (where 0 is no severe impacts expected to 4 very severe impacts expected) (Source: Yu et al. 2015).

Fig. 5

 5 Fig. 5 Conceptual modeling of an industrial installation during a Natech event.

Fig. 6

 6 Fig.6 Fault tree of Electricity utility.

Table 2 .

 2 Reported tsunami inundation heights at the surveyed companies.

Table 3 .

 3 Overall damage

Table 4 .

 4 Causes of total and partial shutdown

		No. of facilities	Percent
	Multinational	1	5%
	Plant size
	Small	13	62%
	Medium	5	24%
	Large	2	10%
	No. of employees
	(0-50)	18	86%
	(51-200)	0	0%
	(201+)	3	14%
	Type of industry
	Chemical (High pressure gas )	18	86%
	Oil refinery	2	10%
	Metallurgical	1	5%
	N	21
	!

ACKNOWLEDGEMENTS

The ResTO-TerRiN project is funded by the French Ministry of Ecology, Sustainable Development and Energy (Ministre de l'Écologie, du Développement durable et de l'Énergie -MEDDE).