text:
G Faivre

S Bottin-Rousseau
email: sabine.bottin-rousseau@insp.jussieu.fr

S Akamatsu
email: silvere.akamatsu@insp.jussieu.fr

The trajectory of subboundary grooves during directional solidification of dilute alloys

Keywords: Solidification, Morphological instability, Grain boundaries, Spatially extended systems Solidification, Instabilité morphologique, Joints de grains , systèmes étendus spatialement

Due to their strong surface energy anisotropy, subboundaries generally do not intersect solid-melt interfaces at right angles. As a consequence, subboundary surface grooves move laterally during solidification, and thereby interact, in alloys, with the solute concentration gradients created in the liquid. We discuss the consequences of this e↵ect during thin directional solidification at speeds (V ) lower than the cellular-instability threshold Vc of the system. We show that the lateral drift of the subboundary grooves slows down (or equivalently the tilt angle of the subboundaries relative to the growth direction decreases) as V increases and vanishes as V approaches Vc.

Résumé

La trajectoire des sillons de sous-joints de grains en solidification directionnelle d'alliages dilués. Du fait de leur forte anisotropie capillaire, les sous-joints de grains ne coupent généralement pas les interfaces solide-liquide à angle droit. En consquence, les sillons de surface des sous-joints se déplacent latéralement en cours de solidification, et ainsi, dans les alliages, entrent en interaction avec les gradients de concentration de soluté créés dans le liquide. Nous discutons les conséquences de cet e↵et en solidification directionnelle d 0 échantillons minces à des vitesses de solidification (V ) inférieures au seuil d 0 instabilité cellulaire Vc du système. Nous montrons que la dérive latérale des sillons des sous-joints se ralentit (ou, ce qui revient au même, l'angle d'inclinaison des sous-joints par rapport à la direction de croissance diminue) lorsque V augmente et s'annule lorsque V atteint Vc.

Introduction

In directional solidification, i.e., when solidification is performed at an imposed speed V under a fixed unidirectional thermal gradient, the solid-melt interface of a non-faceted dilute binary alloy undergoes a morphological transition at a threshold value Vc of V : the solid-melt interface is planar (it follows an isotherm) at V < Vc, and exhibits a cellular (periodically modulated) shape at V > Vc. The physical factors at play are the solute concentration gradient ahead of the interface, the applied thermal gradient G and the solid-liquid surface energy sl . The destabilizing e↵ect exerted on the interface by the solute concentration gradient increases as V increases and exceeds the stabilizing e↵ects of the other two factors when V exceeds Vc [START_REF] Tiller | The e↵ect of growth conditions upon the solidification of a binary alloy[END_REF][START_REF] Mullins | Stability of a planar interface during solidification of a dilute binary alloy[END_REF]. Coriell and Sekerka studied the influence, on this morphological transition, of isotropic grain boundaries intersecting the solid-melt interface using linearized equations for the solidification. They showed that, after a transient, grain boundary surface grooves take on a steady profile in the form of a damped sinusoid, the amplitude of which increases as V increases and diverges as V approaches Vc [START_REF] Coriell | Morphological stability near a grain boundary groove in a solid-liquid interface during solidification of a binary alloy[END_REF].

At the melting points of alloys, high-angle grain boundaries (GBs) generally are wetted by a thin liquid film and have little surface energy anisotropy. By contrast, low-angle grain boundaries, or subboundaries (SBs) remain dry and keep, near the solid-melt interface, the high surface energy anisotropy that they have at lower temperature [START_REF] Bottin-Rousseau | Formation of grain subboundaries during directional solidification below the cellular-bifurcation threshold[END_REF]. Therefore, when local equilibrium is reached, SBs do not intersect the solid-melt interfaces at right angle. We study here some consequences of this fact on the behavior of SBs during the solidification of dilute alloys. Let us first consider a thin (two-dimensional) directional-solidification system at rest (V = 0). We take the x axis parallel to the isotherms, the z axis parallel to the growth direction, the origin of x at the abscissa of the surface-SB junction J and the origin of z at the isotherm that corresponds to the unperturbed planar solid-liquid interface (Fig. 1). We call the tilt angle of the SB with respect to z, the surface energy of the SB and sl the surface energy of the solid-liquid interface. We assume the orientation dependence of to be strong and that of sl to be negligible. We denote the entities related to the right-hand (x 0) side and the left-hand (x  0) side of the system by the indices "R" and "L", respectively. The condition of local equilibrium at J is given by the Young-Herring equation

sl (tL + tR) + = 0 . (1)
In this equation, tL and tR are the unit tangent vectors to the R and L flanks of the groove, respectively, and is the surface tension force or " vector" of the SB defined by

= t + d d n , (2)
where t and n are the unit vectors tangent and normal to the SB, respectively [START_REF] Ho↵mann | A vector thermodynamics for anisotropic surfaces[END_REF]. At rest, the solute concentration in the liquid is uniform. The condition that the surface is continuous in J imposes that the shape of the groove has mirror symmetry with respect to the z axis. Thus takes a value a (where "a" stands for "at rest") such that is parallel to z. This angle generally is a large one owing to the strong surface energy anisotropy of SBs. Incidentally, it should be stressed that the profile of the surface groove in this "symmetrical configuration" is given by the standard formulae for a meniscus [START_REF] Landau | Statistical physics[END_REF]: the anisotropy of the SB does not come into play, except in so far as it alters the Young angle of the surface at J. During solidification of a pure substance, would remain equal to a independently of V . In this note, we show that, during directional solidification of a dilute alloy, deviates from a, breaking the symmetry of the groove. Most importantly, tends to zero as V approaches Vc irrespective of the strength of the anisotropy of the SB.

t R t L melt solid ! a " t x z J SB Figure 1
. Symmetrical surface groove of a subboundary (SB) with a strong surface energy anisotropy at rest in a unidirectional thermal gradient. For the symbols, see the text.

The cellular instability

The most basic model of thin directional solidification neglects the di↵usion of solute in the solid, the rejection of latent heat by the growing solid and the kinetic di↵erence in thermodynamic potentials between liquid and solid at the interface. The temperature distribution is assumed to be ⇥ = ⇥ liq + Gz (⇥ liq : liquidus temperature). The "thermal length" l th = ⇥/D ( ⇥: solid-liquid thermal gap; D: di↵usion coe cient of the solute in the liquid) is taken as a unit of length and l 2 th /D as a unit of time. V and sl are represented by the dimensionless control parameters µ = l th V /D and do = l 1 th ⇥ 1 S 1 v sl ( Sv: entropy of fusion per unit volume), respectively. The unknowns are the profile of the solid-melt interface, denoted by z = ⇣(x), and a reduced concentration field u(x, z) defined such that u(x, z = 0) = 1 and u(x, z = +1) = 0 in the planar steady state of the system. The equations for the solidification in the absence of grain boundary are three in number. They include the di↵usion equation at z > ⇣ (@xx + @zz + µ@z) u = @tu ,

and two equations at z = ⇣, namely, the mass conservation equation

ru.n = [K + (1 K) u] (µ+ @t⇣)z.n , (4)
and the Gibbs-Thomson equation

⇣ = 1 u + do , ( 5
)
where  is the curvature of the interface and K is the partition coe cient. These equations admit a planar stationary solution denoted by a superscript "0", which reads 8 <

:

⇣ 0 = 0 , u 0 = exp ( µz) . (6)
The linear stability analysis of this solution was famously performed by Mullins and Sekerka [START_REF] Mullins | Stability of a planar interface during solidification of a dilute binary alloy[END_REF] and was studied in greater detail by several authors afterwards [START_REF] Wollkind | A nonlinear stability analysis of the freezing of a dilute binary alloy[END_REF][START_REF] Caroli | On the emergence of one-dimensional front instabilities in directional solidification and fusion of binary mixtures[END_REF]. We give an outline of this analysis. A perturbed solution of the form ⇣ 0 + Re ( ⇣), u 0 + Re ( u), where

:

⇣ = Z exp (ikx + !t) , u = U exp ( qz) exp (ikx + !t) , (7)
is plugged into the equations for the solidification. The constant Z is the amplitude of the surface deformation and U is the amplitude of the perturbation to the concentration field. As the system under consideration is infinite in the x and z directions, the wavevector k must be real, entailing that ! and q are real, and q must be positive. We may assume k to be positive without loss of generality. Eq. ( 3) leads to

! = q 2 µq + k 2 . ( 8
)
After linearization with respect to the perturbation, Eqs. ( 4) and ( 5) lead to a homogeneous linear system for Z and U , which reads 8 < :

`! + Kµ 2 ´Z [q (1 K) µ] U = 0 `µ 1 dok 2 ´Z U = 0 . (9)
The characteristic equation for System ( 9) is

! = [q (1 K) µ] `µ 1 dok 2 ´ Kµ 2 . ( 10
)
The concentration-to-deformation amplitude ratio is given by

M = U/Z = [q (1 K) µ] 1 `! + Kµ 2 ´. ( 11
)
The problem is now solved, at least implicitly. An elementary calculation leads from Eqs. ( 8) and ( 10) to an algebraic equations for q as a function of k (more precisely, k 2 ) and µ, and then to the function !(k, µ). The planar solution is stable against those modes, for which ! < 0 and unstable against those, for which ! > 0. The neutral (! = 0) modes, denoted by a subscript "n", are solutions of the system of equations 8 < :

q 2 n µqn + k 2 n = 0 , [qn (1 K) µ] `µ 1 dok 2 n ´ Kµ 2 = 0 . (12)
The calculations show that neutral modes are represented by a closed curve in the (µ, k) plane. Planar solutions are stable against any perturbation outside that "neutral curve", and unstable against some set of perturbations inside it. Let (µc, kc) be the lower critical point of the neutral curve. µc is the cellular-instability threshold of the system, and kc the wavevector of the critical mode. When do << 1, as is generally the case in the experiments,

kc ⇡ (K/2) 1/3 d 1/3 o and µc ⇡ 1 + 3 (K/2) 2/3 d 1/3
o . The speed at which µ = 1 is called the "constitutional supercooling" speed. It is the speed above which the solute concentration gradient overcomes the stabilizing e↵ects of G, but not yet those of sl .

Isotropic grain boundary surface grooves

We wish to calculate ⇣(x) and u(x, z) below the cellular threshold for a stationary solid-melt interface intersected by a grain boundary at x = 0 (Fig. 2). As the emergence of the grain boundary generates a slope discontinuity of the interface, ⇣(x) and u(x, z) must be sought in the form of piecewise-defined (R or L) functions. Let us first consider an isotropic GB. This problem was dealt with by Coriell and Sekerka [START_REF] Coriell | Morphological stability near a grain boundary groove in a solid-liquid interface during solidification of a binary alloy[END_REF], as already mentioned, but we use here a more straightforward approach consisting of replacing the equations for the solidification by the linearized Eqs. ( 8) to ( 12) from the outset. Naturally, the exact solution could be found numerically using, for instance, a boundary-integral method [START_REF] Langer | Studies in the theory of interfacial stability[END_REF]. However, our approach is likely to yield a good approximation of the exact solution in the limit of small p, i.e. for shallow surface grooves. Moreover, it can be used as a method for constructing trial functions for the numerical calculations. When µ < µc, the perturbation due to the GB grove must vanish at large x. Therefore, the wavevector in (7) must be complex, as well as the amplitudes. To be more explicit, the R-side of the interface profile should be written as

Re ( ⇣R) = exp ( kix) [Zr cos (krx) Zi sin (kix)] (13)
with ki > 0, where the subscripts "r" and "i" designate real and imaginary parts, respectively. On the L-side, we have ⇣L(x) = ⇣R( x), or equivalently, kL = kR. The calculation of ⇣(x) and u(x, z) at fixed µ can be performed in two stages. Firstly, qn and kn are calculated using Eqs. [START_REF] Akamatsu | A theory of thin lamellar eutectic growth with anisotropic interphase boundaries[END_REF]. Note that a solution with ki > 0 always exists since only k 2 appears in the equations. Secondly, the amplitudes Zn and Un are calculated using Eq. ( 11) and the boundary conditions at x = 0. Figure 3 shows kn(µ) calculated with values of K and do corresponding to the model transparent alloy CBr4-C2Cl6 at a concentration of 0.15mol%C2Cl6 under a thermal gradient of 10 Kmm 1 [START_REF] Mergy | Quantitative determination of the physical parameters relevant to the thin-film directional solidification of the CBr4-Cr2Clr6 eutectic alloy[END_REF].

As could be expected, kn equals the (real) critical wavevector kc at µ = µc, and tends to ka = id 1/2 o (i.e. the purely imaginary wavevector of surface grooves at rest) as µ ! 0. The evolution of Z and M with µ also deserves attention (right panel of Fig. 3). Note, in particular, that Re(Mn) and Im(Mn) increase very slowly as µ varies from 0 to a value of about 1. Above this value, Re(Mn) increases almost linearly while Im(Mn) goes through a maximum to fall back to zero as µ ! µc. This crossover is linked to the crossing of the constitutional-supercooling threshold. The boundary conditions at x = 0+ include the Young-Herring condition, which reads

kiZr krZi = p , (14)
or, in a compact form, Im (knZn) = p , (15) where p (which is independent of µ and positive) is derived from Eq. ( 1) with = z. Other equations are provided by the continuity conditions. In fact, as the continuity of ⇣ and u is already assured by symmetry, it only remains to meet the condition of continuity on @u/@x, which reads Im (knMnZn) = 0 .

(16)

The solution of the system formed by Eqs. ( 15) and ( 16) is

Zn = p M ⇤ n knMni , (17)
where M ⇤ n is the complex conjugate of Mn. By expanding Eqs. [START_REF] Akamatsu | A theory of thin lamellar eutectic growth with anisotropic interphase boundaries[END_REF] in powers of µ = µ µc, it can be shown that ki ⇠ | µ| 1/2 and Mi ⇠ | µ| 1/2 while Mn tends to a real value as µ ! 0. Thus Eq. ( 17) predicts that Zr diverges like | µ| 1/2 as µ ! µc. Simultaneously, kiZr and krZi both tend to constant values, allowing Eq. ( 14) to be fulfilled until µ = µc. This means that, if the existence of a GB groove greatly a↵ects the nature of the transient leading to the cellular instability, it leaves unchanged the value of the instability threshold, as previously noted by Coriell and Sekerka [START_REF] Coriell | Morphological stability near a grain boundary groove in a solid-liquid interface during solidification of a binary alloy[END_REF].

Anisotropic subboundary surface grooves

When the grain boundary that intersects the interface has a strong surface energy anisotropy, the mirror symmetry between the two sides of the surface groove is broken, and the groove travels along the surface (right panel of Fig. 2). We must therefore look for solutions to the equations for the solidification, which travel as a whole at a constant velocity V T , where T = tan , along the x axis. By making the change of variable x ! x µT t in Eq. ( 7), it can be seen that travelling modes (designated by a subscript "t" if useful) correspond to ! = ikµT . Plotting this relation into Eqs. ( 8) and ( 10) leads to 8 < :

q 2 t µqt + k 2 t = iktµT , [qt (1 K) µ] `µ 1 dok 2 t ´ Kµ 2 = iktµT, (18)
which replaces the neutral-mode equations [START_REF] Akamatsu | A theory of thin lamellar eutectic growth with anisotropic interphase boundaries[END_REF]. The rhs terms break the mirror symmetry and entail that ktL 6 = ktR and qtL 6 = qtR. The calculations show that a R L pair of solutions with kRi > 0 and kLi < 0 exists at any value of µ below µc for su ciently low values of T . Once kD and kL are known, the amplitudes can be determined using the Young-Herring equation and the continuity conditions at x = 0, as explained above. However, Ut is related to Zt by Ut = MtZt, where Mt = [qt (1 K) µ] 1 `iµktT + Kµ 2 ´. There are thus only four independent unknowns, for instance, the real and imaginary parts of ZR and ZL. On the other hand, the Young-Herring equation now yields two independent equations (for the R-and L-flanks of the groove, respectively), which, added to three continuity conditions, bring to five the number of conditions at x = 0. Therefore, the problem has no solution except perhaps for a particular value of T , which, if it exists, depends on µ. This is the central conclusion of this report. At this point, it should be noted that, as qR 6 = qL, the continuity between uR and uL can only be ensured at z = ⇣, but not at z > ⇣ along the x = 0 line. This imperfection is inherent in the linearized theory. While its quantitative consequences are unknown, its is clear that it cannot a↵ect the fact that the lateral drift velocity is uniquely determined at given control parameters, which is a general property symmetry-broken one-dimensional out-of-equilibrium patterns [START_REF] Coullet | Parity-breaking transitions of modulated patterns in hydrodynamics systems[END_REF][START_REF] Akamatsu | A theory of thin lamellar eutectic growth with anisotropic interphase boundaries[END_REF][START_REF] Akamatsu | Lamellar eutectic growth with anisotropic interphase boundaries: Experimental study using the rotating directional solidification methods[END_REF].

We shall now proceed with the calculation of T as a function of µ. The five above-mentioned conditions at x = 0 are:

Im (kRZR) = pR , (19)
Im (kLZL) = pL , (20)
Re (ZL) = Re (ZR) , (21)
Im (MRZR) = Im (MLZL) , (22)
Im (kRMRZR) = Im (kLMLZL) . (23)
As mentioned, they represent a system of five linear equations for the four unknowns ZRr, ZRi, ZLr and ZLi.

Elementary but lengthy calculations show that the condition for this system to have a solution can be written in the form

Fsym (pR + pL) + Fanti (pR pL) = 0 , (24)
where

Fsym = `|kR| 2 + |kL| 2 2kDikLi ´MDiMLi `|MR| 2 + |MR| 2 2MRrMLr ´kDrkLr , (25)
and

Fanti = `|kR| 2 |kL| 2 ´MDiMLi + `|MR| 2 |MR| 2 ´kDrkLr + 2 (kDrMDrkLiMLi kLrMLrkRiMRi) . (26)
The functions Fsym and Fanti are symmetrical and antisymmetrical with respect to permutations of R and L, respectively, hence their names. The quantities kR, kL, MR and ML are functions of T and µ via System (18) while pR and pL are functions of T via Eq. ( 1). Therefore Eq. ( 24) can be viewed as an equation for T at given µ. This equation can be shown to have a single solution provided that the SB has a moderate, non-singular (without facet or unstable orientation range) anisotropy. A calculated T (µ) curve is shown in Figure 4. The anisotropy function taken for the calculation was / sl = 1 0.5 ⇥ 10 2 cos (4 ⇡/4) . With this anisotropy, pa ⇡ 0.05 and Ta ⇡ a ⇡ 1.4 ⇥ 10 2 rad. The values of K and do were the same as in Fig. 3. It can be seen that, as µ increases from 0 to µc, T first remains close to Ta over a large range of µ, and then falls rapidly to zero. Most of the decrease occurs between µ = 1 and µ = µc. The right panel of Fig. 4 illustrates the amplification and shift of u linked to the lateral propagation of the SB groove.

Knowing that T ! 0 as µ ! µc, we can study in more detail the behavior of the T (µ) function in this limit. By expanding Eqs. ( 18), ( 25) and (26) in powers of T it can be shown that Fsym ⇡ 4f anti T , where

f (0) sym = M 2 ni |kn| 2 , (27)
and

f (1) anti = knrRe `m0 k ⇤ n M ⇤ n ´+ MniIm `0 k ⇤ n M ⇤ n ´. ( 28
)
where m 0 and  0 depend on µ = µ µc, but not on T . At small T , Eq. ( 1) can be written in the form pR(T ) ⇡ pR0 + p 0 R T and pL(T ) ⇡ |pL0| + p 0 L T . Plugging these expressions into Eq. ( 24), we obtain

T ⇡ f (0) sym (pR0 |pL0|) f (1) anti (pR0 |pL0|) + f (0) sym (p 0 R + p 0 L ) . ( 29
anti ⇠ | µ| 1/2 . Thus, we can finally rewrite Eq. ( 29) in the form

T ⇡ A pR0 |pL0| pR0 + |pL0| (µc µ) 3/2 , ( 30
)
where A only depends on K and do. In the limit of small do, A reads

A ⇡ 2 4/3 3 1/2 K 1/2 d 1/6 o ( 31
)
and generally is of the order of magnitude of unity. It should be noted that the ratio (pR0 |pL0|) / (pR0 + |pL0|) may be termed "strength of anisotropy" of the SB near = 0. Eq. 30 indicates that ! 0 as µ ! µc does not depend on the strength of anisotropy of the SB. Finally, it should be reminded that the steady state of a SB-free interface in the limit of small µ is di↵erent from its state at rest. Therefore the finding (displayed in Fig. 4) that, in the presence of a SB, T ! Ta as µ ! 0 is not quite obvious and requires justification. We have seen that Mnr ! 0, Mni ! 0, knr ! 0 and kni ! d 1/2 o as µ ! 0. The fact that T does not tend to zero imposes that mRi ! T and mLi ! T . An elementary calculations shows that, as a consequence, Fsym contains a term that does not tend to zero as µ ! 0 (namely, the term 4k 2 ni mDimLi ! 4d 1 o T 2 ) , while Fanti has none. Therefore, according to Eq. ( 24), the SB must approach an orientation such that pR + pL = 0 as µ ! 0. In other words, this proves that, at low V , the tilt angle of the SB is close to the symmetrical configuration corresponding to T = Ta.

Conclusion

The theory presented in this study predicts that, during the directional solidification of a dilute alloy, the tilt angle of a SB should vary from a finite value (corresponding to the symmetrical configuration of its surface grooves and therefore depending on the strength of the anisotropy of the SB) to zero as the solidification speed varies from 0 to the cellular-instability threshold of the system. More specifically, it has been established that most of the decrease in the tilt angle takes place in the interval between the constitutional-supercooling velocity and Vc, indicating that the process is mostly linked to the overcoming of the surface energy barrier by the concentration gradients ahead of the interface. The experimental data about this phenomenon are not very numerous although the question of the interaction between solute and SBs can be traced back to the origins of modern solidification science [START_REF] Atwater | The influence of impurities on the macromosaic structures of tin and lead[END_REF]. It seems to be established that the tilt angle of the SBs is essentially V -independent at low V and that SBs run essentially parallel to the growth direction above the cellular-instability threshold [START_REF] Bottin-Rousseau | Formation of grain subboundaries during directional solidification below the cellular-bifurcation threshold[END_REF]. However, little is known about the transition between these two extreme regimes from the experimental standpoint. We conclude with a remark about SBs in the cellular regime. Several authors have noted that, on crossing the cellular-instability threshold, the SBs become attached to inter-cellular grooves and then run parallel to z, as mentioned. We wish to point out that those inter-cellular grooves, which contain a SB generally exhibit asymmetric comma-shaped profiles that contrast with the symmetrical profiles of perfect inter-cellular grooves. As explained in Fig. 5, this means that the angles between the SBs and their vectors are large, which clearly supports the basic presuppositions (strong capillary anisotropy of the SBs, local equilibrium at the junctions of the SBs with the solid-melt interface) of this work.

Figure 2 .

 2 Figure 2. Steady-state surface grooves during directional solidification. Left panel: isotropic grain boundary (GB). Right panel: anisotropic subboundary (SB). The groove is drifting towards the left. The deformation of the solid-melt interface has been strongly exaggerated.

Figure 3 .

 3 Figure 3. Subcritical (µ  µc) and supercritical (µ µc) neutral modes calculated with K = 0.75 and do=6.5⇥10 3 . Left panel: wavevector kn . Hatched area: instability range of the planar steady solution. Right panel: deformation amplitude Zn (top graph) and amplitude ratio Mn (bottom graph).

) 6 !Figure 4 .

 64 Figure 4. Left panel: Tilt angle of an anisotropic SB with respect to the direction of growth as a function of the solidification speed. For the values of the material parameters assumed in the calculation, see the text. Right panel: Profile of the surface (⇣) and perturbation to the solute distribution along the x axis ( u) at µ = 0.15. aniso: : same SB as in the left panel. iso: isotropic SB, everything else equal. The scales of the x and ⇣ axes are in a ratio of 2, 000 : 1 By expanding Eqs. (27) and (28) in powers of µ it can be shown that f (0) sym ⇠ | µ| and f

Figure 5 .

 5 Figure 5. Directional-solidification front in a shallow-cell regime (V slightly higher than Vc) with a SB attached to an intercellular groove. Sketch from Fig.4of Ref.[START_REF] Bottin-Rousseau | Formation of grain subboundaries during directional solidification below the cellular-bifurcation threshold[END_REF].

  Figure


domain: phys nlin chim

text:
Mikaël Perrut

Sabine Bottin-Rousseau

Gabriel Faivre

Silvère Akamatsu

Dynamical instabilities of rodlike eutectic growth patterns: A real-time study

Keywords: Eutectic solidification, directional solidification, solidification microstructures

We present a real-time experimental study of the rodlike growth patterns formed during directional solidification in a non-faceted transparent eutectic alloy, succinonitrile-(d)camphor. Slightly convex isotherms were used to slowly increase the pattern spacing from an appropriate starting value to the threshold spacing for the rod elimination or rod splitting instabilities allowing a quantitative determination of these thresholds as a function of the solidification rate V. We show that the threshold spacing for rod splitting obeys the general ⇠ V 1/2 scaling law of eutectic growth whereas the threshold spacing for rod elimination deviates from this law at low V, exhibiting the same overstability e↵ect as that previously reported for lamellar eutectic patterns. We demonstrate that topological defects (walls between hexagon domains) play an important role in rod splitting processes. We also describe a spatially incoherent mode of oscillation that we observed in disordered rodlike patterns.

Introduction

Nonfaceted eutectic growth is one of the oldest known experimental example of 2D pattern-forming out-of-equilibrium system [1,2,3]. It shares the same ⇤ Corresponding author Email address: akamatsu@insp.jussieu.fr () translational and rotational symmetry as, for instance, Rayleigh Benard convection experiments and, like the latter, is multistable [START_REF] Manneville | Dissipatives structures and weak turbulence[END_REF]. At given control parameters, it can settle in space-periodic stationary patterns either lamellar (periodic in one direction) or rodlike (consisting of 2D periodic arrays). The spacing and orientation of these patterns are not sharply selected but are boundary condition-and history-dependent. Various "defective" patterns, ranging from space-periodic patterns traversed by a few topological defects to fully disordered patterns, can also be more or less stationary. The selection processes between these di↵erent types of growth patterns as a function of boundary conditions and solidification history have not yet been fully elucidated. The investigation of such processes mostly relies on real-time observations of well-controlled systems over su ciently long period of times. Because eutectic patterns are by their nature immersed in a solidifying melt, and, moreover, usually have spacing values ( ) in the micrometer range, their in-situ study has long been limited to quasi 1D systems (leading, incidentally, to an almost complete elucidation of the dynamics of such systems [START_REF] Karma | [END_REF]6]). It is only a few years ago that an experimental method of observing insitu 2D eutectic patterns in transparent alloys was proposed [7]. This method was applied to lamellar eutectic patterns, revealing that the stability of such patterns is limited by a zigzag bifurcation followed by lamella disruption and branching on the large-side [8], in agreement with the results of phase-field numerical studies [9]. It was also found that a necessary condition for lamellar patterns to form was the use of a solidification setup with a slight misalignment of the direction of growth and the axis of the thermal gradient. Without this additional ingredient, called "thermal bias", the eutectic growth patterns remained indefinitely labyrinthic [10].

The same method of in-situ observation is currently being used to study the dynamics of rodlike eutectic patterns. The experiments are performed in the transparent alloy succinonitrile-(d)camphor (SCN-DC), which is known to be a nonfaceted rodlike eutectic [11,12]. The question of the long-term growth regime of this system was addressed in a previous report, which pointed out the crucial role played by "imperfections" of the experimental setup in the selection of this regime [13]. More specifically, it was shown that a slight forward curvature of the isotherms drove the rodlike pattern towards an operating point close to its large-stability limit, which corresponds to a rod splitting instability. During this study, we used the same mild external forcing to study quantitatively the threshold spacings for both the rod splitting and the rod elimination instabilities, the latter corresponding to the small-stability limit of the growth patterns. Although curved isotherms are a common feature in bulk directional solidification, it is only relatively recently that their e↵ects on growth patterns have begun to be studied in detail [14,15,16]. Before presenting our findings about the instability threshold spacings, we shall therefore sum up our previous conclusions about the e↵ects of convex isotherms on the dynamics of rodlike eutectic patterns.

Experimental Methods

The eutectic point of the SCN-DC alloy is at T E = 38.3 o C and C E =13.9 mol%DC [11]. The eutectic plateau connects an almost pure-SCN body centered cubic phase and an almost pure-DC hexagonal phase. Thereafter, these phases are called SCN and DC, respectively. The alloys we used were within less than ±0.1mol% of C E . A solid at this concentration has a DC volume fraction of about 0.23 and grows from the melt by forming a rodlike eutectic pattern, with rods of the DC phase embedded in a SCN matrix [12]. Directional-solidification experiments were conducted as follows. The alloy is put in cartridges made of two 300 mm-thick glass plates separated by plastic spacers delimiting an empty space of cross section 6 ⇥ 0.4mm 2 . A device acting as grain selector is placed near one end of the cartridges allowing one to grow single, or at least large, eutectic grains. A filled cartridge is placed between a cold oven and a hot oven adjusted to give a thermal gradient of 8 ± 1 Kmm 1 in the region of interest. Solidification is carried out by pulling the sample towards the cold end of the gradient at a constant rate V. Let z be the direction opposite to that of the pulling, y the direction perpendicular to the sample plane and x the direction perpendicular to y and z. In theory, z coincides with the axis of the thermal gradient and the isotherms are planes perpendicular to z. In practice, the isotherms are generally slightly curved and tilted with respect to z.

During the pulling, the average position of the solid-liquid interface remains locked onto an isotherm close to T E and is observed with an optical microscope. We use a dark-field method based on the fact that the three (SCN, DC, liquid) phases of interest have di↵erent refractive indices. This method can be briefly described as follows (For details, see [7]). The solid-liquid interface is observed from the outside of the crucible through the liquid under oblique view. The angle between the y axis and the direction of observation is of about 40 o . The solid is illuminated obliquely and slightly from below. The light rays that travel along a DC rod are refracted by the curved interface ("DC cap") between this rod and the liquid to form a caustic at a short distance above the interface. The images collected essentially consist of arrays of bright spots, each of which represents the caustic associated to a DC cap (Fig. 1). The as-received images are processed for subtracting background and stretched in the y direction in order to correct for the contraction due to the inclination of the direction of observation. The resolution limit of this method is of about 3µm. It is well known that the order of magnitude of the spacing of a eutectic pattern is given by the Jackson-Hunt scaling law

⇠ JH = p K JH V 1/2
, where K JH is a material constant [3]. In near-eutectic SCN-DC, K JH = 10.2 ± 1.5 µm 3 s 1 [12]. During this study, the V-range used was of 0.018 to 0.25µms 1 , corresponding roughly to a -range of 5 to 30µm.

Following common practice, we shall present the measurements of threshold spacings in terms of the dimensionless ("reduced") spacing ⇤ = / JH . When illustrating local instability processes, we shall make use of Voronoi diagrams with a gray scale ranging from white to black as the number of sides of the Voronoi cells ranges from 4 to 8. With this representation, hexagon domains appear as medium-grey areas, while domain walls mostly consist of strings of penta (lightgrey )-hepta (dark-grey) pairs of cells (Fig. 2). We define the local spacing (t) of a particular DC cap in a rodlike pattern as the average of the distances between this DC cap and its nearest neighbors at time t, as determined by the Voronoi construction. A number of other local parameters, such as the distances and angular orientations of the di↵erent nearest neighbors, are of course needed for a complete geometrical description of a rod splitting or rod elimination event in an irregular environment. Empirically, however, we found that a single parameter, namely, (t), was su cient for a study of the dynamics of these events under the conditions of our experiments [17]. We therefore limit ourselves to this parameter thereafter. We define the mean spacing ¯ (t) as the average of the local spacings of all the DC rods observed at time t. The experimental uncertainty on each individual measurement of (t) was of about ±1.5µm. For an essentially hexagonal (more exactly, triangular) arrangement of the rods, ¯ ⇡

h A 1 N(t) sin ⇡/3 i 1/2
, where A is the area of the observation window and N(t) is the number of rods in this window. Since N(t) > 300 in our experiments, the uncertainty on ¯ was negligible.

Directional solidification with convex isotherms

In our experiments, the T E -isotherms were flat in the x direction, but slightly convex (bulging into the liquid) in the y direction due to the di↵erences in thermal conductivities of the various materials of which the samples were composed [13,17]. The corresponding radius of curvature R was measured to range between 2.5 and 5mm in the di↵erent samples used. Due to the fact that DC rods (i.e. the trajectories of the DC caps) remain approximately perpendicular to the envelope of the growth front during solidification, the growth patterns, whatever their intrinsic dynamics, were continually stretched in the y direction at a rate V/R. Consequently, when the pattern was quiescent (i.e. as long as no rod creation or rod elimination events took place, except of course at the sample walls), the mean spacing increased over time according to

¯ (t) = ¯ o e t/⌧ , (1)
where ⌧ = 2R/V and ¯ o is a constant. As R was much larger than all the characteristic distances of the growth patterns, ⌧ was long compared to the characteristic times of the processes of interest. Let s (V) and c (V) be the threshold spacings for the rod splitting and rod elimination instabilities, respectively. These quantities are not uniquely defined when the rodlike pattern is disordered (see below) but we may neglect this fact for the moment. Let us consider a long-duration solidification run at constant V starting from ¯ o < c . Such a run can be divided into four stages (Fig. 3a). During Stage I, rod elimination is active and ¯ increases rapidly. Observations showed that, during this stage, the growth pattern is very disordered and the width of the -distribution is large. Rod elimination stops completely when the lower tail of this distribution is above c .

During Stage II, the pattern is quiescent and ¯ rises exponentially, as explained above. Simultaneously, the degree of order of the pattern increases. After a suciently long solidification time, the growth pattern becomes organized in hexagon domains close to either of two particular orientations, namely, the H-orientation, which has a nearest-neighbor direction parallel to the sample walls, and the Vorientation, which has a nearest-neighbor direction perpendicular to the sample walls (Fig. 2). The formation process of such structures will not be considered here. We simply note that they are certainly linked to wall and thermal-bias effects, as are also the structures perpendicular to samples walls in lamellar eutectics [10]. Near the end of Stage II, the hexagon domains were about 10 in size and were separated from each other by almost ideally thin walls. The hexagons were then substantially elongated (by about 20%) in the y direction. This is the highest degree of order that was observed during this study.

Stage III corresponds to the gradual increase of the rod splitting frequency from zero to the value (R 1 V per rod) at which it balances out the stretching e↵ect due to convex isotherms. Finally, Stage IV corresponds to a statistically steady, indefinitely long, growth regime, in which the rod splitting frequency and the mean spacing remain constant, or, more exactly, fluctuate about constant values. In fact, during Stage IV, the whole histogram of local spacings averaged on relatively short time scales and, in particular, its mean value < ⇤ > are constant and history-independent (Fig. 4). Experiments performed at di↵erent solidification rates yielded < ⇤ >= 1.03 ± 0.04, indicating that < ⇤ > obeys the ⇠ V 1/2 scaling law. Finally, we note that, during Stages III and IV, the growth pattern retains a structure in V-and H-domains in spite of the recurring rod-splitting events. However, the size of the domains is smaller and the domain walls are more di↵use than in Stage II. Examples are given below.

The rod splitting instability

The solidification rate programs we used during this study are sketched in Fig. 3b. After a long transient in a quiescent state at constant V, we observed rod splitting events to occur either inside hexagon domains or in domain walls. Voronoi analyses of intra-domain rod splitting processes are shown in Figure 5. The most noteworthy facts are as follows. We noted no important di↵erence in rod splitting processes between V-and H-domains. In both cases, the rod splitting events occurred in quick succession at neighboring sites, which suggests the existence of a cooperative mechanism (cascade) with a propagation rate on the same order of magnitude as V. The value of ⇤ when the splitting occurred was of 1.15±0.02 and the deformation of the hexagons in the y direction was of about 25%. The direction of splitting (i.e. the direction of the vector between the two daughter branches resulting from the splitting of a rod) was that of the most elongated nearest-neighbor direction. After the passage of a rod splitting cascade, the hexagonal ordering of the domains was restored at a smaller spacing (⇡1.11). To sum up, the intra-domain rod splitting process is the closest approximation we could observe to an instability in a perfect hexagonal rodlike pattern. The observations suggest that a hexagonal rodlike pattern, when subjected to unidirectional deformation, becomes unstable against rod splitting for a deformation slightly less than 25%, or, in terms of local spacing, at ⇤ hd s (hd: hexagon domain) slightly smaller than 1.15.

Rod splitting processes located in domain walls are shown in Figure 6. The most characteristic feature of such processes is that they only a↵ected rods with particularly large ⇤-values. These rods were always located inside domain walls and were mostly, but not exclusively, represented by heptacoordinated cells in the Voronoi constructions. In the examples under consideration, the average spacing inside the hexagon domains was of about 1.05, and was thus clearly smaller than ⇤ hd s , whereas the local spacing of the rods that were about to split was substantially larger than 1.1 (see below). It should also be noted that the sequencing of splitting events did not suggest that a cooperative mechanism was at play. In other words, the di↵erent rod splitting events seemed to occur independently from each other, supporting the view of a strictly localized mode of splitting. This mode of splitting was largely predominant in Stage IV. Broadly speaking, its localized character explains the persistence of a multidomain structure with a small grain-size and di↵use walls during this stage.

Figure 7a shows the time evolution of the local spacing of a rod during a splitting process. The measured data clearly deviated from a ⇤(t) = ⇤ o e t/⌧ law after a time point at which ⇤ ⇡ 1.14. Thus ⇤ s was smaller or equal to that value. No theory of morphological instabilities in disordered hexagonal patterns being available, we fitted various elementary models of instability processes to the data points in order to obtain a lower bound for ⇤ s . In brief, we concluded that ⇤ s = 1.2 ± 0.1.

During the transients leading up to the rod splitting events, we observed a gradual ovalization of the optical images of the rods, which was undoubtedly due to an elongation of the DC rods that were about to split (Figs. 6b and7b). This rod elongation process coincided in time with the transitory deviation of the local spacing from a e t/⌧ law, and stopped immediately after the splitting had occurred and is therefore clearly part of the instability process itself. It has no obvious link with the stationary rod elongations that exist in stable non-equilateral hexagonal patterns [18].

The rod elimination instability

The initial condition ¯ o < c was implemented by applying a downward Vjump (See Section 3). Rod elimination started at a fast rate almost immediately after the V-jump and led to fully disordered growth patterns through a process that could not be resolved with our method of observation. Observable isolated rod elimination events continued to occur in this disordered state for some time (Fig. 8). These events were not accompanied by any detectable rearrangement of the pattern beyond the nearest neighbors of the eliminated rod. In particular, we did not observe any analogue of the phase-di↵usion dynamics that is, in theory, the precursor to rod elimination in regular hexagon patterns [19]. Figure 9 shows a ⇤(t)-curve measured after a downward V-jump together with the best-fit exponential curve to the data points of Stage II. The spacing value ⇤c (0.78 ± 0.03, in this example) at the time when the data points stop to deviate significantly from the exponential curve gives us a reasonable (over)estimate of the average local value of the rod elimination threshold spacing in disordered states.

Figure 10 displays ⇤c -values measured at di↵erent values of V. It is striking that ⇤c , contrary to ⇤ s , does not follow the ⇠ V 1/2 scaling law, but decreases rapidly as V decreases. A similar e↵ect was recently demonstrated, both experimentally and numerically, for lamellar eutectic patterns [20,21] and was called "overstability" in reference to the fact that ⇤ c = 1 within the usual normal-growth approximation [22]. A complete theory of the overstability e↵ect is still lacking, but it was argued that ⇤ c should obey the following semi-empirical formula:

1 1 ⇤ 2 c + B G⇤ c K r V = 0 , ( 2
)
where G is the thermal gradient, K r is a material constant and B is an adjustable constant on the order of 0.1. In near-eutectic SCN-DC, K r = 0.028Ksµm 2 [12]. Least-squares fitting of Eq. ( 2) onto experimental and numerical data for various lamellar eutectic alloys yielded values of B ranging from 0.067 to 0.15 [21,12,23]. As shown in Fig. 10, Eq. ( 2) with B = 0.075 ± 0.03 (including the uncertainties on the values of G and K r ) fits the experimental data for the rod elimination threshold reasonably well. This value of B falls in the range of values found for lamellar patterns, indicating that the overstability e↵ect has the same physical origin in disordered rodlike patterns as in regular lamellar eutectic patterns. In other words, the overstability e↵ect seems to be essentially insensitive to the geometrical particulars of the patterns. It is therefore reasonable to argue that the results shown in Fig. 10 should apply, at least semi-quantitatively, to regular rodlike patterns as well.

Oscillations of disordered rodlike patterns

Two kinds of spatially coherent modes of oscillation have recently been reported in hexagonal patterns. One was a "sublattice" (or period-tripling) oscillation, found in hexagonal cellular solidification patterns of dilute alloys [24,25]; the other was a "row" oscillation, occurring in a hexagonal liquid column array [26]. A variety of other coherent modes of oscillation are in theory possible in hexagonal patterns, but have not yet been observed experimental or numerically [27]. We looked for, but did not find any spatially coherent mode of oscillation in rodlike patterns during this study. There are several possible reasons for this, among which an insu cient size of the hexagon domains in our system. Surprisingly, however, we observed a spatially incoherent mode of oscillation which lasted for the duration of Stage IV in two independent experiments performed at V = 0.018µms 1 , which was the lowest V-value we used. The two experiments gave similar results. We limit ourselves here to a brief characterization of this mode of oscillation.

A snapshot and a local space-time diagram of a disordered oscillatory rodlike pattern are shown in Fig. 11. It can be seen that the oscillation manifested itself mostly through an essentially periodic, although somewhat irregular, variation in brightness of the image of each DC rod. A closer analysis showed that this was accompanied by a small-amplitude periodic displacement of the images, indicating that the variation in brightness was in fact mostly due to a variation in the angle with z of the rods. A rod-by-rod Fourier analysis of the oscillation revealed that the time period T varied slightly from rod to rod. It ranged from 19 to 23.5min, and was thus close to the value of < ¯ > /V (⇡ 22min). This is most probably meaningful given that the time period of the 1 O oscillations of lamellar eutectic patterns obeys a T ⇠ /V scaling law (See Fig. 30 of Ref. [6]). There was no obvious spatial correlation between the di↵erences in period, phase or amplitude of the oscillations of the di↵erent rods. There has been no previous report of phenomena related to these spatially incoherent oscillations, as far as we know, except perhaps for the transient oscillation localized on penta-hepta pairs that was noted by Pirat and Gil in imperfect hexagonal liquid column array [26].

Conclusion

The main conclusions of this study are summed up in Fig. 12. First, the threshold spacing for rod splitting obeys a ⇠ V 1/2 scaling law and is clearly above the Jackson-Hunt scaling length (which corresponds to ⇤ = 1). Second, the threshold spacing for rod elimination exhibits a large deviation from the ⇠ V 1/2 scaling law at low V (overstability) in conformity with Eq. ( 2). The best-fit value of the empirical constant appearing in this equation (B = 0.075) is comparable to those previously found in lamellar eutectic patterns. Third, the steady-state value of the mean spacing during the last stage of directional solidification with convex isotherms is significantly smaller than s . This is consistent with the fact that, in this stage, rod splitting is localized in domain walls and, more precisely, only a↵ects those rods, which have a much larger local spacing than average. Fourth, fully disordered rodlike patterns exhibit, under conditions which remain to be clarified, a spatially incoherent mode of oscillation, whose time period T probably obeys a T ⇠ /V scaling law.

It should be kept in mind that these results were established in multidomain rodlike eutectic patterns, for which no theory of morphological instabilities is available yet. As mentioned above, there are reasons to believe that the thresh-old spacings for rod splitting and rod elimination in perfect hexagonal rodlike patterns are not very di↵erent from those presented here, but this is not certain. It now seems possible, and would indeed be very interesting, to study this question, and, more generally, explore the dynamics of multidomain rodlike patterns, by phase-field simulations [18]. Future experimental studies should undertake to grow single-domain rodlike eutectic patterns using thicker samples and appropriate solidification rate programs [13,28,29].

List of figure captions

Figure 1 :

 1 Figure 1: Top view of a rodlike eutectic pattern during directional solidification of a near-eutectic SCN-DC alloy. V = 0.035µms 1 . The bright spots are the optical images of the DC rods. The upper and lower edges of the field of view correspond to the inner faces of the crucible. Vertical dimension of the snapshot: 400 µm. Insert: Schematic 3D drawing.

Figure 2 :

 2 Figure 2: Schematic Voronoi representation of a multidomain rodlike eutectic pattern. H and V refer to the orientations of the hexagon domains. (See the text).

Figure 3 :

 3 Figure 3: (a) Definition of the successive stages of a long-duration solidification run at constant V. ¯ : mean spacing. c : rod elimination threshold. s : rod splitting threshold. Dash line: ¯ (t) = ¯ o e t/⌧ . (b) Solidification rate programs for the study the rod splitting (1) and rod elimination (2) instabilities. Hatched area: range without rod splitting or elimination.

Figure 4 :

 4 Figure 4: Histogram of local spacings averaged over a short period of time (⇡ 0.02⌧) during the last stage of a long-duration solidification run[13,17]. < ⇤ >: mean value. ⇤ s : rod splitting threshold (this study).

Figure 5 :

 5 Figure5: Intra-domain rod splitting processes. From left to right: Voronoi diagrams of two micrographs taken immediately before and after the process, respectively. The rod splitting events have been labeled alphabetically in chronological order. (a) V-domain. V = 0.11µms1 . Time interval between snapshots: 2min. Horizontal dimension of a snapshot: 55µm. (b) H-domain. V = 0.053µms 1 . Time interval between snapshots: 98min. Horizontal dimension of a snapshot: 140µm.

Figure 6 :

 6 Figure 6: Rod splitting process located in domain walls. V = 0.035µms 1 . (a) Same representation as in Fig. 5. Time interval between snapshots: 66min. Horizontal dimension of a snapshot: 150µm. (b) Successive snapshots of a single rod splitting process. From left to right: t = 0, 14, 21, 28min. Horizontal dimension of a snapshot: 65µm.

Figure 7 :

 7 Figure 7: (a) Time evolution of the local spacing of a DC rod during a rod splitting process. V = 0.035µms 1 . Symbols: measured values. Continuous line: ⇤ = ⇤ o e t/⌧ (⇤ o = 1.096; ⌧ = 2.4 ⇥ 10 3 min). ⇤ s : Largest possible value for the rod splitting threshold. (b) Time evolution of the ovalization (aspect ratio minus 1) of the optical image of the DC rod. The nonzero (⇡ 0.5) baseline is due a residual misalignment of the optical system.

Figure 8 :

 8 Figure8: Successive snapshots of an isolated rod elimination event following a V-jump from 0.042 to 0.028µms 1 . From left to right: t = 0, 50, 100min. Horizontal dimension of a snapshot: 50µm.

Figure 9 :

 9 Figure 9: Mean spacing as a function of time after a V-jump from 0.071 to 0.042 µms 1 applied at t = 0. Symbols: measured values. Continuous line: ⇤ = ⇤ o e t/⌧ (⇤ o = 0.768, ⌧ = 3.3 ⇥ 10 3 min). ⇤c : rod elimination threshold.

Figure 10 :

 10 Figure 10: Rod elimination threshold spacing as a function of solidification rate. Symbols: measured values. Continuous line: Eq. (2) with B = 0.075.

Figure 11 :

 11 Figure 11: (a) Snapshot of a disordered oscillating rodlike pattern. V = 0.018µms 1 . Horizontal dimension: 440µ. (b) Time evolution of the intensity distribution along a 180 µm-long, 1pixelthick area. The rods that are apparently created (eliminated) are in fact entering (leaving) the selected area.

Figure 12 :

 12 Figure 12: Morphology diagram of multidomain rodlike eutectic patterns in near-eutectic SCN-DC alloys. Hatched area: range without rod creation or elimination. Square symbols: oscillatory patterns. < ⇤ >: mean spacing at long times in the presence of convex isotherms.

Acknowledgments

We warmly thank V. T. Witusiewicz, L. Sturz, and S. Rex from ACCESS (Aachen, Germany) for kindly providing purified succinonitrile-(d)camphor alloy and for useful advices. We are glad to acknowledge many fruitful discussions with M. S ¸erefoglu and M. Plapp. This work was supported by the Centre National d' Etudes Spatiales, France.


Sourabh Dharurkar
email: sourabh1289@gmail.com

Ganesh Fulari
email: 2ganesh.fulari33@gmail.com

Ram Agarwal

Suraj Gawhale

S Pattanaik

Real Time Hand Gesture Recognition For Content Based Image Retrieval

•

Navigating and manipulating in virtual environments.

•

Can be used to develop Medical based application.

• Vision based application.

• Gaming application.

VIII. CONCLUSION

This paper describes the design and implementation of a hand gesture recognition system using just one color camera. The developed system can obtain a high recognition rate of hand gestures. Our current research mainly focuses on the single hand static gestures for the simplicity. However, we are far from building a generalpurpose gesture recognition system. The dynamic gestures and two-handed gestures should also be explored in the future, since they are more expressive and allow more natural interaction. The system implementation shown above adopted by the detection of hand or highlighting hand posture in the input image using calibration, motion detection, and resulted with segmented hand posture is defined with results. As in case of recognition techniques video processing needs speed and higher intensity in providing experience of virtual and real time experience. So we have studied and described how efficiently and effectively we can retrieve images in database by using real time hand gestures.

IX. ACKNOWLEDGMENT

We would like to thank Prof. S.Pattanaik (RSCOE, computer Department), Prof. V.M.Barkade (RSCOE, Computer Department) for judging and assessing the design of proposed system and her valuable guidance.

X. REFERNCES


Laurent Théry
email: laurent.thery@sophia.inria.fr

A Formalisation of the Generalised Towers of Hanoi

This notes explains how the optimal algorithm for the generalised towers of Hanoi has been formalised in the Coq proof assistant using the SSReflect extension.

Introduction

The famous problem of the towers of Hanoi was proposed by the french mathematician Édouard Lucas. It is composed of three pegs and some disks of different size. Here is a drawing of the initial configuration for 5 disks 1 : Initially, all the disks are pilled-up in decreasing order of size on the first peg. The goal is to move them all to another peg. There are two rules. First, only one disk can be moved at a time. Second, a larger disk can never be put on top of a smaller one. 1 We use macros designed by Martin Hofmann and Berteun Damman for our drawings.

1

The towers of Hanoi one of the classical example that illustrates all the power of recursion. If we know how to solve the problem for n disks, then the problem for n + 1 disks can be solved in 3 steps. Let us suppose we want to transfer all the disks to the last peg. The first step uses recursion and moves the n-top disks to the intermediate peg.

The second step moves the largest disk to its destination The last step uses recursion and moves the n disks on the intermediate peg to their destination. This simple recursive algorithm is also optiomal: it produces the minimal numbers of moves. In particular, if we look at each recursion depth, the key idea is that largest disk always moves once from its current peg to its destination.

The generalised version of the towers of Hanoi considers an arbitrary initial configuration and an arbitrary final configuration. These two configurations must be valid : there is no larger disk on top of a smaller disk. The problem is to find an algorithm that generates the minimal nunber of moves that connects the two configurations. Here, the naive recursive algorithm is still applicable to solve the problem but does not lead to an optimal algorithm. This can be illustrated by 3 disks when trying to go from the initial configuration:

to the final position If the naive recursive approach, that tries to move the largest disk only once, leads to a 7-move long solution as depicted in Figure 1 at page 12. The optimal solution requires instead to move the largest disk twice and is 5-move long as depicted in Figure 2 at page 2. In the following, we explain how the generalised towers of Hanoi has been formalised in the Coq proof assistant and how an algorithm that solves this problem has been proved correct.

Pegs

The set of natural numbers strictly smaller than three I 3 is used to represent the three pegs.

Definition peg := I 3 .

An operation that is frequently used in the algorithm is, having two arbitrary peg p 1 and p 2 to get the third one, This is done by the function opeg using some arithmetic:

Definition opeg p 1 p 2 : peg := inord (3 -(p 1 + p 2 )).
where inord is the function that injects a natural number into the type I n .

Disks

A parameter n is used. A disk is then an element of I n .

Definition disk := I n .

The comparison of the respective size of two disks is simply performed by comparing their natural number.

Configurations

Disks are ordered from the largest to the smallest on a peg. This means that a configuration just needs to record which disk is on which peg. It is then defined as a finitite function from disks to pegs.

Definition configuration := {ffun disk → peg}.
Note that in this encoding, we do not have invalid configurations.

A configuration is called perfect if all its disks are on a single peg. It is the case for the initial and final configurations in the standard towers of Hanoi. The constant function that always returns p is then the perfect configuration where all the disks are on the peg p.

Definition perfect p := [ffun d ⇒ p].
We also need two functions to build configurations. The first one builds a new configuration from a configuration c by putting all the disks of size strictly smaller than m on the peg p: this new configuration is perfect at depth m.

Definition mk perfect m p c := [ffun d ⇒ if d < m then p else c d].
The second function performs a single change on the configuration c : the disk d is moved to the peg p.

Definition setd d p c 2 := [ffun d 1 ⇒ if d 1 = d then p else c d 1 ].

Moves

A move is defined as a relation between configuration. A parameter m is introduced. It represents a bound on the size of the disk that has been moved. Its purpose is to let us perform proof by induction: the induction is then performed on the depth of the moves m rather than on the number of disks n.

Definition move m : rel configuration := [rel c 1 c 2 | [∃d 1 : I n , [&& d 1 < m, c 1 d 1 = c 2 d 1 , [∀d 2 , d 1 = d 2 ⇒ c 1 d 2 = c 2 d 2 ], [∀d 2 , c 1 d 1 = c 1 d 2 ⇒ d 1 ≤ d 2 ] & [∀d 2 , c 2 d 1 = c 2 d 2 ⇒ d 1 ≤ d 2 ]]]].
The definition simply states that there is a disk d 1 htat fulfills five conditions. Its size is smaller than m. It is has moved. It is the unique disk that has moved. No disk is in on top of

d 1 in c 1 . No disk is on top of d 1 in c 2 .
Trivial facts first need to be derived. For example, the move relation is cumulative and symmetrical.

Lemma moveW m 1 m 2 : m 1 ≤ m 2 → subrel move m1 move m2 . Lemma move sym m c 1 c 2 : move m c 1 c 2 = move m c 2 c 1 .
A first interesting is then that before and after a move, all the disks smaller than the disk that has moved are all pilled up on the same peg:

Lemma move perfectl m d c 1 c 2 : move m c 1 c 2 → c 1 d = c 2 d → c 1 = mk perfect d (opeg (c 1 d) (c 2 d)) c 1 . Lemma move perfectr m d c 1 c 2 : move m c 1 c 2 → c 1 d = c 2 d → c 2 = mk perfect d (opeg (c 1 d) (c 2 d)) c 2 .
An important corollary is that if a disk d moves twice on different pegs then the two configurations after the first move and before the second one are perfect configurations at depth d.

Lemma move twice d c 1 c 2 c 3 c 4 : c 1 = c 4 → c 1 d = c 2 d → c 3 d = c 4 d → move d.+1 c 1 c 2 → connect move d c 2 c 3 → move d.+1 c 3 c 4 → c 2 = mk perfect d (opeg (c 1 d) (c 2 d)) c 2 ∧ c 3 = mk perfect d (c 1 d) c 2 .
This is a key lemma that is used to get a direct lower bound (2 d -1) on the number of moves that are needed for connecting c 2 are c 3 . We will explain this later.

In order to be able to decompose paths under the move m.+1 relation, two inversion lemmas are needed. The first one checks if the disks m moves. If it is the case, it singles out its first move. The decomposition is presented as an inductive predicate in order to get the decomposition by the direct application case: pathSP tactic on a path. The second inversion lemma considers a path at depth d + 1 for which the disk d may have moved but at the end it remains on the same peg. It builds the "restricted" path for which the disk does not move.

Lemma pathS restrictE d c cs :

path move d.+1 c cs → last c cs d = c d → {cs 1 | [∧ path move d c cs 1 , last c cs 1 = last c cs & size cs 1 ≤ size cs ?= iff cs 1 = cs ]}.
The number of moves of the restricted path gets strictly smaller only if there is a move of the disk d in the path cs.

Naive algorithm

Our definition of the naive algorithm works at depth m, starts with a configuration c and tries to move all the disks of size less than m to the peg p.

The strategy works as follows. If the disk m is already on the peg p, there is nothing to do at depth m: the algorithm is recursively called at depth m -1. Otherwise, the disk m needs to be moved. The algorithm is first called at depth m -1 to move all the disks of size smaller than m -1 to the intermediate peg p 1 . Then, the disk of size m is moved to the peg p and finally the algorithm is called a second time to move the disk of size smaller than m -1 to the peg p. Formally, this gives2 : This gives us directly that it computes the actual size of the naive algorithm.

Fixpoint rpeg path rec m c p := if m is m 1 . + 1 then if c m = p
Lemma size rpeg path rec pr m c p : size (rpeg path rec m c p) = size rpeg path rec m c p.

As a matter of fact 2 m -1 is the maximum a naive solution can get

Lemma size rpeg path rec pr m c p : size rpeg path rec m c p ≤ 2 m -1.

With these results, it is possible to prove the optimality of the naive algorithm for the special case where the initial configuration is perfect. Note that we also prove that the optimal solution is unique. The proof works by a double induction: one induction on the depth m and another strong induction on the size of cs. So, when proving the step case at depth m, the property is known to hold for all paths at depth m -1 and for paths at depth m which size is strictly smaller than cs. We then simply do a discussion on the number of moves the disk m does in the path cs:

-If it does not move, the inductive hypothesis for m -1 gives directly the result.

-If it moves once, the path cs mimics the strategy of the naive algorithm at depth m. So, combining the two applications of the inductive hypothesis for m -1 (one before and one after the move) gives the result.

-If it moves more than once, there are two possibilities. Either the disk visits a peg more than once (this is always the case if the disk moves more than two times). In this case, the lemma pathS restrictE gives us a stricly smaller path on which we can apply the inductive hypothesis on the size. Either the disk has moved twice on different pegs. The lemma move twice tells us that if we consider the path the first move and before the second, it connects two perfect configurations at depth m -1. So the inductive hypothesis for m -1 and the lemma size rpeg path rec 2p tell us that its size is greater than 2 m-1 -1. The same holds for the path after the second move and the final configuration. Altogether, this gives a size for cs (not considering the path before the first move of the disk m) that is larger than 1 + (2 m-1 -1) + 1 + (2 m-1 -1) = 2 m . This is strictly more than the bound 2 m -1 for the naive algorithm given by the lemma size rpeg path rec pr.

This ends the proof.

Optimal algorithm

In order to define the optimal algorithm, we first define the symmetric of the naive algorithm that goes from a perfect configuration to any configuration by simply reversing the path. This algorithm is clearly correct and optimal. The optimal algorithm is defined recursively in order to find the first disk that has to be moved. When this disk is found, it simply chooses the best solution between moving it directly to where it has to go (going from c 1 to c 3 then c 2 ) and moving it twice (going from c 1 to c 3 then c 4 and fimally c 2 ) using the intermediate peg p. The computation of these two solutions can use the naive algorithm since one of the two configurations that are connected is perfect, so we know it is optimal.

Inductive

  pathS spec (d : I n ) c : seq configuration → bool → Type := pathS specW : ∀cs, path move d c cs → pathS spec d c cs true | pathS spec move : ∀c 1 cs 1 cs 2 , path move d c cs 1 → c d = c 1 d → move d.+1 (last c cs 1 ) c 1 → path move d.+1 c 1 cs 2 → pathS spec d c (cs 1 ++ c 1 :: cs 2 ) true | pathS spec false : ∀cs, pathS spec d c cs false. Lemma pathSP d c cs : pathS spec d c cs (path move d.+1 c cs).

  then rpeg path rec m 1 c p else let p 1 := opeg (c m) p in let c 1 := setd m p (mk perfect m 1 p 1 c) in rpeg path rec m 1 c p 1 ++ c 1 :: rpeg path rec m 1 c 1 p else [::] Definition rpeg path c p := rpeg path rec n c p. Lemma size rpeg path rec 2p m p 1 p 2 c (c 1 := mk perfect m p 1 c) : size rpeg path rec m c 1 p 2 = (2 m -1)(p 1 = p 2 ).

Lemma

  rpeg path rec min m c 1 p cs (c 2 := mk perfect m p c 1 ) : path (move m) c 1 cs → last c 1 cs = c 2 → size rpeg path rec m c 1 p ≤ size cs ?= iff (cs = rpeg path rec m c 1 p).

Definition

  lpeg path rec m p c := rev (belast c (rpeg path rec m c p)). Definition lpeg path p c := lpeg path rec n p c.

The formalisationIn this section, we present the different elements our formalisation, starting with pegs, disks, configurations and moves, then we describe the naive recursive algorithm and finally the optimal one.

This is a simplified version. The real code is less readable because m is an ordinal.

Note that c 1 is configuration after the disk m has been moved to the peg p: all the disks smaller than m are on the intermediate peg p 1 .

The first basic property that needs to be proved is that this algorithm is correct: what is build is a path that goes from the configuration c to the perfect configuration on peg p Lemma rpeg path correct c p (cs := rpeg path c p) :

path

This directly gives the fact that any configuration is connected to any perfect configuration

Lemma move connect rpeg c p : connect (move n) c (perfect p).

Since the relation is symmetric, this gives that any two configurations are connected.

There is always a solution to the generalized tower of Hanoi.

If we are only interested by the size of the solution, it is possible to give an algorithm that computes the size of the connection given by the naive algorithm.

Note that in this version, there is only one recursive call. The justification for this comes from the fact that this algorithm returns 2 m -1 when called on a perfect configuration c that is different from p.

It is then easy to derive that this algorithm is correct. The proof for optimality is similar to the one for the naive algorithm : we simply show that the largest disk cannot move three times in the optimal solution.

Lemma hanoi path correct c 1 c 2 (cs := hanoi path c 1 c 2 ) :

path

Conclusion

We have presented a formalisation of the generalised towers of Hanoi. The formalisation clearly benefits from the SSReflect library. In particular, finite function have been a convenient tool to encode configuration. Most of the proofs are elementary. Without surprise, the difficult part is to get the optimality results. We had to device two dedicated inversion principles in order to mechanise the case distinctions that were needed. The complete proof is available at http://www-sop.inria.fr/marelle/Laurent.Thery/Hanoi.


Abdul-Hamid Soubra

Richard Kastner

A Benmansour

Passive earth pressures in the presence of hydraulic gradients

Keywords: earth pressure, ®lters, limit state design/analysis, pore pressures, seepage, sheet piles and cofferdams

de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Passive earth pressures in the presence of hydraulic gradients A. H. SOUBRA, Ã R. KA STNER { and A. BENMANSOUR Ã

The paper describes a variational approach applied to the limit equilibrium method for calculating the effective passive pressures of a cohesionless soil, taking into consideration the seepage ¯ow. It is shown that in the general case of non-homogeneous and non-isotropic hydraulic properties of the soil medium, the shape of the slip surface which veri®es the three limiting equilibrium equations of the soil mass at failure is a log-spiral. It is also shown that the passive earth pressure calculation is independent of the normal stress distribution along this surface. The variational limit equilibrium method is equivalent to the upper bound method in limit analysis for a rotational log-spiral mechanism. Numerical results of the coef®cients of passive earth pressures in the presence of seepage ¯ow are presented and discussed.

INTRODUCTION

The design of deep sheeted excavations is often dominated by the ¯ow of water around the sheet piles. The seepage ¯ow in¯uences the stability of the excavation where bulk heave or piping may occur. While the piping takes place at the excava tion level, the heaving is more catastrophic and its risk is usually evaluated by considering a rectangu lar failure mechanism adjacent to the wall [START_REF] Terzaghi | Theoretical soil mechanics[END_REF]. The vertical force equilibrium of this soil mass is then considered by neglecting the vertical frictional forces along the vertical faces of this mechanism.

Based on laboratory model tests, [START_REF] Kastner | Excavations profondes en site urbain: proble Ámes lie Âs a Á la mise hors d'eau. Dimensionnement des soute Ánements butonne Âs[END_REF] has shown that the failure of the sheet piling structures in the presence of seepage ¯ow is not only due to the heaving phenomenon but may also occur due to the reduction of the passive earth pressures in front of the wall. Our aim in this paper is to propose an outline for the calculation of the effective passive pressures, taking into ac count the seepage forces.

Looking for a simple model capable of correctly describing the behaviour of soil in the passive state in the presence of seepage ¯ow, we opted for the limit equilibrium method. This approach is based on a priori hypotheses concerning the shape of the slip surface (kinematic function) and the normal stress distribution (static function) along this sur face. The variational approach applied to the limit equilibrium method has been employed in order to avoid the restrictions of such a priori hypotheses.

OVERVIEW OF PREVIOUS VARIATIONAL ANALYSIS

The variational limit equilibrium method has been used by several authors in geotechnical en gineering. [START_REF] Kopacsy | Three-dimensional stress distribution and slip surface in earth work at rupture[END_REF] applied this approach to the three dimensional slope stability problem; however, no explicit solution is offered. Several investigators, for example [START_REF] Dorfman | Variational method of investigating the stability of slopes[END_REF], [START_REF] Garber | Variational method for investigating the stability of slopes[END_REF], [START_REF] Revilla | The calculus of variations applied to stability of slopes[END_REF] and [START_REF] Ly | Determination of active earth pressure by a variational method[END_REF], have used the calculus of variations to avoid introducing an assumption concerning the shape of the slip surface but they made an assump tion concerning the normal stress distribution along this surface. A more interesting analysis consists of ®nding the two unknown functions without any a priori assumptions. Thus, Baker & Garber (1977, 1978) applied the variational approach to the two dimensional slope stability problem, and Garber & Baker (1977) and [START_REF] Castillo | One application of the calculus of variations to bearing capacity of foundations[END_REF] ap plied this approach to the problem of the bearing capacity of a strip footing. Then, [START_REF] Garber | Extreme value problems of limiting equilibrium[END_REF] treated the problems of slope stability, bear ing capacity and earth pressure distribution in a uni®ed manner as a single problem. Later, this method was used by Leshchinsky et al. (1985), [START_REF] Ugai | Three-dimensional stability analysis of vertical cohesive slopes[END_REF] and [START_REF] Leshchinsky | Three dimensional slope stability: end effects[END_REF] to study the three dimensional slope stability problem, and by [START_REF] Leshchinsky | Stability of membrane reinforced slopes[END_REF], who applied it to the reinforced slope stability problem. Finally, [START_REF] Leshchinsky | Pseudostatic seismic stability of slopes: design charts[END_REF] applied the variational limit equilibrium method to the seismic stability of slopes, and presented interesting results in the form of design charts. It is to be noted here that the variational approach has been the subject of some controversy. Particularly, the work by [START_REF] Castillo | A critical analysis of some variational methods in slope stability analysis[END_REF] shows that the functional has no minimum. Notice, however, that experience indicates that the slip surface determined by the variational analysis reasonably duplicates reality [START_REF] Leshchinsky | Pseudostatic seismic stability of slopes: design charts[END_REF]. Consequently, the solutions given by the variational limit equilibrium method are very interesting. We present in the following the application of this method to the passive earth pressure problem, taking into account the seepage forces.

VARIATIONAL APPROACH OF THE PASSIVE EARTH PRESSURE PROBLEM

Figure 1 shows a double walled cofferdam sub jected to a seepage ¯ow where H is the total head loss and u(x, y) represents the distribution of the pore water pressures in the soil medium.

The assumptions made in the analysis can be summarized as follows:

(a) The soil is cohesionless. It is homogeneous and isotropic with respect to its angle of internal friction ö. (b) The soil medium is non homogeneous and non isotropic with respect to the hydraulic properties. It is composed of n permeable layers overlying impermeable rock. Each layer is characterized by its coef®cients of per meability K hi and K vi . (c) The breadth B 0 (Fig. 1) is large enough so that there is no interaction of the two failure mechanisms which develop in front of the two walls of the cofferdam.

(d) The resultant P P of the effective passive pressures is assumed to act at the bottom third of the penetration depth (Fig. 2(a)). This force can be expressed as follows:

P P K P ã9 f 2 2 (1)
where K p is the coef®cient of passive earth pres sure in the presence of seepage ¯ow, ã9 is the submerged unit weight of the soil, and f is the penetration depth. The variational approach is brie¯y presented in this paper. For more details, refer to [START_REF] Soubra | Application de la me Âthode variationnelle au proble Áme de de Âtermination des pressions passives des terres[END_REF]. Fig. 2(a) illustrates the formulated problem and shows the notation used. A potential slip sur face y(x) is subjected to a total normal stress ó (x). Both functions y(x) and ó (x) are assumed to be continuous. Using Coulomb's failure criterion ô(x) [ó (x) u(x)]tan ö ó 9(x)tan ö, the global limiting equilibrium equations for the soil mass (Fig. 2(b)) can be written as

P p cos ä U 1 x1 x0 [ó 9(tan ö y) u y] dx (2a) P p sin ä x1 x0 [ó 9(1 tan öX y) u ã sat ( f y)] dx (2b) P p cos äX 1 U 1 X 2 x1 x0 [ó 9(1 tan öX y)x ó 9(tan ö y) y uxuy y ã sat ( f y)x]dx (2c)
Notice that one of the endpoints x 0 is null 2b) and (2c), one rea lizes that the effective passive force P p is a func tional of two functions, y(x) and ó 9(x). The mathematical problem of the passive earth pressure is to ®nd those functions which give the minimum value of the passive force functional and simulta neously satisfy all three equations of limiting equi librium (equations (2a), (2b) and (2c)).

H f B 0 /2 K h1 , K v1 K h2 , K v2 K hi , K vi K hn , K vn
Using equation (2b) to de®ne P p , while consider ing the other two equilibrium equations (equations (2a) and (2c)) as constraints (i.e. equations that must be satis®ed), the passive earth pressure problem is a variational isoparametric one with a variable end point. This variational problem is equivalent to the minimization of an auxiliary functional G:

G L 0 ë 1 L 1 ë 2 L 2 (3)
where L 0 , L 1 and L 2 are given as follows:

L 0 ó 9(1 tan ö y) u ã sat ( f y) L 1 ó 9(tan ö y) u y L 2 ó 9(1 tan ö y)x ó 9(tan ö y) y ux uy y ã sat ( f y)x (4)
ë 1 and ë 2 are the Lagrange undetermined multi pliers. Finally, the two extremal functions y(x) and ó 9(x) must satisfy the following conditions:

(a) The system of Euler's differential equations for the functional G:

σ(x ) P p (b) δ U 1 x o X 1 X 2 X x 1 x u σ ′ τ W sat y
A, we have x A y A 0. At the variable end point B, a variational condition must be satis ®ed. This condition is called the `transversality condition' and it can be written as follows:

G y dG d y ó 9 dG d ó 9 x x1 äx 1 dG d ó 9 x x1 äó 9 1 dG d y x x1 ä y 1 0 (6)
where ä is a variational operator.

The ®rst Euler equation

Combining equations ( 3) and (5a) yields a dif ferential equation. The solution of this equation in a polar coordinate system is a log spiral (Fig. 3) whose equation is given as follows:

r r 0 exp(è è 0 )tan ö (7)
Note that the log spiral function has a particular property, that the resultant of the forces (ó 9 dl) and (tan öó 9 dl) passes through the pole of the spiral. Hence, the moment equation about the pole is independent of the stress distribution ó 9(x), and may be used for the determination of the effective passive force. The two remaining equilibrium equa tions may be satis®ed by every ó 9(x) distribution that has two degrees of freedom. Thus, one has to ®nd the critical è 0 and è 1 angles which satisfy the moment equilibrium equation and give the mini mum value of the effective passive force P p . This is done by a two dimensional minimization proce dure of P p with respect to è 0 and è 1 .

The independence of the effective passive force from the normal stress distribution can also be shown, due to the special property of the present functional. This functional can be written as fol lows:

G ó 9 f (x, y, y) g(x, y, y) (8)
G is linear in ó 9 and is independent of ó 9. The ®rst Euler equation implies that f (x, y, y) 0. Substituting this equation into equation ( 8), one can see that this functional becomes independent of ó 9 as follows: G g(x, y, y). This result is a direct consequence of the shape of the slip surface. Thus, the ®rst Euler equation transforms the pas sive functional from a functional of two unknown functions to a functional of a single function. Therefore, it is possible to solve the passive earth pressure problem by simply minimizing the new functional G without specifying the normal stress distribution. Finally, it is easy to see that the moment equation of the rotational log spiral me chanism around the centre is identical to the work equation for the same mechanism in the upper bound method in limit analysis. Thus, solving the passive earth pressure problem by writing the moment equation around the centre of the log spiral will give an upper bound solution of the exact solution for an associated ¯ow rule Coulomb material.

The second Euler equation and the transversality condition

Since the aim of this study is the determination of the critical effective passive force, the results obtained so far are enough to solve the problem. Indeed, it has been shown by Baker & Garber (1977) that the second Euler equation and the transversality condition give the normal stress dis tribution. In this work, this is also the case, but since the ó 9(x) distribution is not necessary to assess the effective passive force, we will not express these equations.

Existence of a minimum of the functional

The second variation of the passive earth func tional shows that this functional is degenerated. Thus, we cannot say mathematically that there is a minimum [START_REF] Castillo | A critical analysis of some variational methods in slope stability analysis[END_REF]. This dif® culty concerning the existence of a minimum can be overcome due to the equivalence between the variational limit equilibrium method and the upper bound method in limit analysis for a rotational mechanism as has been shown before. This equiva lence has been established in a more general way by [START_REF] Castillo | Variational methods and upper bound theorem[END_REF] and Leshchinsky et al. (1985).

CALCULATION SCHEME OF THE PASSIVE EARTH PRESSURES

As mentioned before, the solution of the passive earth pressure problem by a variational limit equi librium method is equivalent to saying that the equation of moment equilibrium must be satis®ed for the soil mass bounded by the log spiral and the ground surface. As shown in Fig. (2b), the forces acting on the collapse mechanism are (a) the saturated weight of the soil mass between the log spiral surface and the ground surface, (b) the effective passive force P p which is inclined at ä to the normal of the sheet pile, (c) the pore water pressures along the penetration depth and the log spiral surface, and (d) the effective normal and tangential stress distributions along y(x). The moment equilibrium equation can be written as follows:

W sat (r 0 cos è 0 X ) P p sin ä(r 0 cos è 0 )

P p cos ä r 0 sin è 0 f 3 M 1 M 2 (9)
where M 1 and M 2 represent the moments of the force U 1 and of the pore water pressures along y(x) respectively. They are given as follows:

M 1 f 0 u( r 0 sin è 0 y) dy r 0 sin è 0 f 0 u d y f 0 uy d y (10) M 2 è1 è0 (u sin ö) r dè cos ö r r 0 2 tan ö exp ( 2è 0 tan ö) 3 è1 è0 u exp(2è tan ö)dè (11)
From equation ( 9), one can easily see that P p is a function of the two parameters è 0 and è 1 which completely describe the failure mechanism. The most critical K p value is obtained by a minimiza tion procedure of the K p coef®cient given by equa tion (1), that is [K P (2P P )a(ã9 f 2 )] with respect to the two parameters mentioned above. A compu ter program has been developed for assessing the minimal K p values and the corresponding critical slip surfaces.

NUMERICAL RESULTS OF THE PASSIVE EARTH PRESSURES

Case of no seepage ¯ow Table 1 compares the passive earth pressure coef®cient K p for ö 408 and äaö 1a2 ob tained from the present analysis with that of other authors in the case of no seepage ¯ow. The com parison of the present upper bound solution with the upper and lower bound solutions given respec tively by [START_REF] Chen | Limit analysis solutions of earth pressure problems[END_REF] and [START_REF] Lysmer | Limit analysis of plane problems in soil mechanics[END_REF] shows that the difference with the lower bound solution is smaller than 3%, which means that the present solution is very close to the exact solution for an associated ¯ow rule Coulomb material. On the other hand, the currently accepted values given by [START_REF] Sokolovski | Statics of soil media[END_REF] and [START_REF] Caquot | Tables for the calculation of passive pressure, active pressure and bearing capacity of foundations[END_REF] lie in the range between the best upper and lower bound solutions given by the limit analysis method. Finally, the comparison between the present result and the ones given by the limit equilibrium methods (Coulomb, 1776;[START_REF] Shields | Passive pressure coef®cients by method of slices[END_REF] shows that the traditional limit equilibrium method may greatly overestimate or underestimate the passive earth pressure coef® cients due to the a priori assumptions concerning the shape of the slip surface and the normal stress distribution along this surface.

Case of seepage ¯ow

From equation ( 9), the determination of the effective passive force P P requires the determina tion of the terms M 1 and M 2 . The hydraulic head distribution j(x, y, z) in the soil medium is gov erned by the following equation:

d dx K x dj dx d d y K y dj d y d dz K z dj dz 0 (12)
In some simple cases, such as the case of the single sheet pile driven into a semi in®nite homo geneous soil medium, the pore water pressure distribution is given analytically (Soubra & Kast ner, 1992). For more complex geometry or for a multi layered soil medium with different coef® cients of permeability, the pore water pressure dis tribution cannot be known analytically and it requires a numerical resolution of equation ( 12).

For the double walled cofferdam considered in this paper, the numerical method used for the determination of the potential ®eld in the soil medium is the well known ®nite difference meth od where the differential equation (equation 12) is approximated by a ®nite difference equation. The boundary conditions used are shown in Fig. 4. The medium is discretized by a rectangular mesh (Fig. 5). The ®nite difference equations written at the different nodes form a system of linear equations whose unknowns are the values of the hydraulic head at the nodes. This system is solved by the Gauss Seidel method, using over relaxation in order to accelerate the rate of convergence.

The determination of the moments M 1 and M 2 has been made by numerical integration using the Gaussian quadrature method where the pore water pressures along both the sheet pile and the log spiral surface are determined by numerical interpolation. Let us pass now to the presentation of some numerical results obtained from the com puter program. Notice that all subsequent results concern the case when ã sat aã w 2 and B 0 a2 10 m (Fig. 1).

Case of a homogeneous and isotropic soil medium

Case of a single layer of in®nite depth. [START_REF] Soubra | In¯uence of seepage ¯ow on the passive earth pressures[END_REF] published the results of the passive earth pressure coef®cients in the presence of seepage ¯ow in the case of a single sheet pile wall driven into a homogeneous and isotropic semi in®nite soil medium where the hydraulic head can be known analytically. Note that the same results have also been obtained by the present analysis using the ®nite difference method for the determi nation of the hydraulic heads.

Figure 6 shows the variation of the passive earth pressure coef®cient with Ha f for ö 308 and for four values of äaö (äaö 0, 1a3, 1a2, 2a3). For a zero K p value, the corresponding Ha f value is the same for different äaö values ( Ha f 2 X 78). This means that the angle of friction at the soil structure interface has no effect on the Ha f value causing failure by heaving. This fact can Fig. 5. Finite difference mesh 4. Boundary conditions for the seepage ¯ow of the double-walled cofferdam be explained as follows: when the effective passive force vanishes, there is no interaction at the soil structure interface and we have the tradi tional heaving phenomenon. For the same case, Terzaghi's approach gives a value of Ha f at fail ure equal to 2´82. Furthermore, the piping phe nomenon which appears for the critical hydraulic gradient at the point D (Fig. 4) occurs for a value of Ha f equal to 3´14. It should be mentioned that the numerical results have shown that in the case of a homogeneous and isotropic semi in®nite soil medium, the failure by heaving will occur before the piping phenomenon as long as ö is smaller than 458.

∂ ∂ν ϭ 0 ϕ ∂ ∂ν ϭ 0 ϕ h H A f F G E v Z B D ϭ Z ϭ H ϩ Z ϭ H ϩ Z ∂ ∂h ϭ 0 C ϕ ϕ ∂ ∂ν ϭ 0 ϕ ϕ ϕ Fig.
Figure 7 shows some charts of the variation of the passive earth pressure coef®cient as a function of Ha f for different values of ö (ö 208, 258, 308, 358, 408) and äaö (äaö 0, 1a3, 2a3). From these ®gures, the reduction of this coef®cient is quasi linear for the Ha f values varying from 0 to 2´5.

Finally, it should be mentioned that the K p value increases with ã sat . This increase is to be expected, since the soil weight has the favourable effect of increasing the stability of the soil mass in front of the sheet pile. For example, when ö 358, äaö 2a3 and Ha f 2, the increase of the passive earth pressure coef®cient attains 22% when ã sat aã w increases from 2 to 2´2. Case of a single layer of ®nite depth D. The numerical results obtained from the present pro gram have shown that the passive earth pressure coef®cient increases with the Da f decrease. The minimum relative depth Da f necessary to obtain the results of the semi in®nite case must be greater than or equal to 6.

Figure 8 shows some charts of the variation of the passive earth pressure coef®cient as a function of Ha f for different values of ö (ö 208, 258, 308, 358, 408) and äaö (äaö 0, 1a3, 2a3) and for Da f 2.

Case of a homogeneous and non isotropic soil medium. Figure 9 shows the variation of the passive earth pressure coef®cient with K h aK v when ö 358, Ha f 2 and äaö 0, 1a3, 2a3 and 1. There is a large decrease of the passive earth pressure coef®cient up to K h aK v 100. Beyond this limit, the passive earth pressure coef®cient tends to an asymptote. This can be explained by the fact that the equipotential lines in front of the sheet pile become quasi horizontal beyond a certain value of K h aK v and, thus, the potential ®eld does not change any more in the zone concerned with the failure mechanism.

Figure 10 shows some charts of the variation of the passive earth pressure coef®cient as a function of Ha f for different values of ö (ö 208, 258,

δ/φ ϭ 2/3 δ/φ ϭ 1/3 δ/φ ϭ 0 φ ϭ 35Њ, H/f ϭ 2 log (K h /K v )
Fig. 9. K p versus K h aK v for ö 358 and Haf 2 in the case of a homogeneous and non-isotropic semi-in®nite medium 308,358,408) and äaö (äaö 0, 1a3, 2a3) and for Da f 2 when the permeability ratio K h aK v 50.

Case of an isotropic two layered soil medium. In this section, we consider the frequent case of a cofferdam driven into a two layered soil medium.

Case where the bottom of the sheet pile lies in the upper layer (case A). Figure 11 shows the case of an isotropic two layered soil medium where the permeability coef®cients of the upper and lower layers are respectively K 1 and K 2 .

Figure 12 shows the variation of the passive earth pressure coef®cient with K 1 aK 2 when ö 358, äaö 2a3 and Ha f 2.

For the case of a single layer (K 1 aK 2 1), the passive earth pressure coef®cient is equal to 3´81. When K 1 aK 2 . 10, one obtains the passive earth pressure coef®cient corresponding to the case of a single layer of limited depth (ö 358, äaö 2a3, Ha f 2 and Da f 2), since the lower layer can be considered as an impermeable substra tum. Notice, however, that for cases when the lower layer has a greater permeability coef®cient than the upper layer (K 1 aK 2 , 1), most of the head loss is concentrated in the upper layer, resulting in greater pore water pressures. Consequently, the passive earth pressure coef®cient decreases with the K 1 aK 2 decrease.

Case where the bottom of the sheet pile lies in the lower layer (case B). Figure 13 shows the case of an isotropic two layered soil medium. The unique difference from the previous case is that the bottom of the sheet pile wall lies in the lower layer.

Figure 14 shows the variation of the passive earth pressure coef®cient with K 1 aK 2 when ö 358, äaö 2a3 and Ha f 2. As in the previous section (case A), the passive earth pressure coef®cient in the case of an isotro pic single layer (K 1 aK 2 1) is equal to 3´81. For K 1 aK 2 . 100, the head loss takes place solely in the lower layer, and the upper layer can be consid ered as a ®lter increasing the global stability of the soil in front of the sheet pile. Consequently, the increase of the passive earth pressure coef®cient with the K 1 aK 2 increase is a logical phenomenon. On the other hand, for cases when the lower layer has a greater permeability coef®cient than the upper layer (K 1 aK 2 , 1), most of the head loss takes place in the upper layer, resulting in a sig ni®cant reduction of its resultant body force (buoy ant weight seepage force). Consequently, the passive earth pressure coef®cient decreases with the K 1 aK 2 decrease. However, it should be noted that this calculation scheme is only valid as long as the vertical hydraulic gradient in the upper layer is smaller than the critical gradient i c ã9aã w , as otherwise failure by heaving of the upper layer will occur due to the fact that the saturated weight of this layer is equal to the resultant of the pore water pressures on the base of this layer. This limitation is shown in Fig. 14 by a vertical dotted line.

CONCLUSION

The variational limit equilibrium method was applied to the passive earth pressure problem, tak ing into account the seepage ¯ow due to dewater ing. It showed that the failure mechanism, in the general case of non homogeneous and non isotro pic hydraulic properties of the soil medium, is a log spiral. This method is equivalent to the upper bound method in limit analysis for a rotational log spiral mechanism.

In the case of no seepage ¯ow, the present upper bound solution is the smallest one existing in the literature and is very close to the currently ac cepted results of [START_REF] Caquot | Tables for the calculation of passive pressure, active pressure and bearing capacity of foundations[END_REF].

In the case of seepage ¯ow, the present mechan ism allowed us to determine the reduction of the passive earth pressures. For the limiting case corre sponding to zero passive pressures, the present mechanism describes the traditional heaving phe nomenon in front of the sheet pile. The numerical results have shown that: shown that: (i) for great values of K 1 aK 2 , when the bottom of the sheet pile lies in the upper layer, one obtains the case of a single layer of limited depth, and when the bottom of the sheet pile lies in the lower layer, the upper layer can be considered as a ®lter; and (ii) for small values of K 1 aK 2 , most of the head loss takes place in the upper layer, thus resulting in a signi®cant reduction of the effective passive pressures.

Finally, one can see that the effect of soil aniso tropy and non homogeneity is signi®cant for the reduction of the passive earth pressures. Thus, the assessment of the reduction of the effective passive pressures taking into account these parameters is of great interest in the practice of geotechnical engi neering.

NOTATION B 0 breadth of the double-walled cofferdam D layer depth dl elementary surface along the slip surface f penetration depth H total hydraulic head loss h, v coordinate system whose origin is at point F i c critical hydraulic gradient K 1 , K 2 isotropic permeability coef®cients of layers 1 and 2 K h , K v horizontal and vertical permeability coef®cients K p coef®cient of passive earth pressure K x , K y , K z permeability coef®cients in the principal directions x, y and z M 1 moment of the force U 1 M 2 moment of the pore water pressures acting on the slip surface P p effective passive force r 0 , r 1 initial and ®nal radius of the log-spiral slip surface U 1 resultant of pore water pressures at the soil structure interface u pore water pressure W sat saturated weight of the soil mass OAB X distance between the wall and the line of action of the force W X 1 distance between the bottom of the wall and the normal component of the effective passive force X 2 distance between the bottom of the wall and the line of action of force U 1 y(x) equation of the slip surface in the (x, y) coordinate system y d y dx ä friction angle at the soil structure interface ö angle of internal friction of the soil j hydraulic head ã9 submerged unit weight of the soil ã sat saturated unit weight of the soil ã w unit weight of water ë relaxation factor ë 1 , ë 2 Lagrange's undetermined multipliers è 0 , è 1 angles de®ning the log-spiral slip surface ó (x) total normal stress distribution along the slip surface ó 9(x) effective normal stress distribution along the slip surface ó 9 dó 9 dx ô(x) tangential stress distribution along the slip surface

Fig. 1 .

 1 Fig. 1. Double-walled cofferdam in a multi-layered soil medium

  The constraint equations. (c) The boundary conditions: at the ®xed endpoint

Fig. 2 .

 2 Fig. 2. (a) Slip surface for passive earth pressure analysis. (b) Free body diagram

  Fig. 3. Log-spiral slip surface for passive earth pressure analysis

Fig. 6 .Fig

 6 Fig. 6. K p versus Haf for ö 308 and äaö 0, 1a3, 1a2 and 2a3 in the case of a homogeneous isotropic semi-in®nite medium

Fig. 8 .

 8 Fig. 8. K p versus Haf for different values of ö and ä in the case of a homogeneous isotropic single layer of ®nite depth (Daf 2)

Fig. 12 .

 12 Fig. 10. K p versus Haf for different values of ö and ä in the case of a homogeneous and non-isotropic single layer (Daf 2)

  (a) The heaving of a soil mass in front of the sheet pile occurs before the piping phenomenon in the case of a homogeneous and isotropic semi in®nite soil medium as long as ö is smaller than 458. (b) The effective passive pressures increase with the decrease of the layer depth in the case of a single layer problem.

Fig. 14 .

 14 Fig. 13. Case of an isotropic two-layered soil medium: case B

Table 1 .

 1 Passive earth pressure coef®cient K p as given by different authors for ö 408, äaö 1a2

	Authors	K p