Claudia Papandrea

Livio Battezzati
email: livio.battezzati@unito.it

A STUDY OF THE α↔γ TRANSFORMATION IN PURE IRON: RATE VARIATIONS REVEALED BY MEANS OF THERMAL ANALYSIS

Keywords: phase transformations, DTA user supplied): Iron, Interface migration Iron, massive phase transformation, interface migration, differential thermal analysis

HAL is a multi-disciplinary

During the last decade there has been continuous interest in experimental and theoretical studies on the α↔γ transformation in Fe and Fe alloys due to its relevance in steel processing [1][2][3][4][5][6][7][8][9][10][11]. Papers have appeared in relation to the kinetics of such transformation [1,3,4], the topology of grains [2], the mobility of the interface between the two phases, dilatation and enthalpic effects.

Dilatometry is often employed for the monitoring of the transformation on cooling, the kinetics of which is modelled according to a nucleation and growth mechanism [3,5]. An amount of undercooling is always reported for the start of the transformation. Assuming instantaneous nucleation at the undercooling temperature, the growth of ferrite was analysed for a series of Fe alloys containing up to 2 at% of Co, Cu, Mn, Cr, Al, and the molar mobility of the interface, M, was obtained as 2. 4•exp(-140000/RT) mol•m•(J•s) -1 when the impurity nitrogen content was below 40 ppm. It decreased about an order of magnitude for an order of magnitude increase in nitrogen content [4]. The mobility is related to the interface velocity, v, by

m G M v ∆ = • (1)
where ∆G m is the driving force per mole. Since the driving force is of the order of 10-100 J•mol -1 at the given undercoolings [1,6], the interface velocity, v, is estimated in the range from 10 -4 to 10 -6 ms -1 . Direct observation of the interface movement in a low C and low N steel provided values of v from 1•10 -4 ms -1 to 9•10 -4 ms -1 [7].

Differential thermal analysis (DTA) has been employed to obtain the transformed fraction of ferrite as a function of heating rate in Fe-1at%Mn and Fe-2at%Mn steels. From the subsequent kinetic analysis, interface velocities from 0.2 to 0.9 • 10 -6 ms -1 were calculated [1]. Further analyses 3 with a Fe-1.89at%Mn steel showed that the transformation occurs in two stages the velocities of which fall in the range just reported [8].

By performing both isochronal dilatometric and differential thermal analyses, the γ→α transformation was found to occur either normally or abnormally in Fe-Co, Fe-Mn and pure Fe according to the grain size in the sample [6,9]. With grains of the order of tens of micrometers (alloys) or 273-372 µm (pure Fe), a single step transformation was detected (normal behaviour), whereas up to three steps were detected when the grains were 439 µm (abnormal behaviour). The steps have been related to variations of the formation rate of ferrite and have been detected during cooling. As for the reverse α→γ transformation, no detailed analysis was apparently performed with either technique. It was, anyway, shown that the abnormal behaviour could be revealed by means of DTA isothermal analyses at temperatures very close to the equilibrium one. From the kinetic modelling of the dilatometric results, the interface velocities for the normal transformation were derived as 2-3•10 -7 m•s -1 for Fe-2.26 at%Mn, 4-8•10 -6 m•s -1 for Fe-1.79 at%Co and 3-6•10 -6 m•s -1 for pure Fe [6,10]. The interface velocities for the abnormal transformation were in excess of 2•10 -4 ms -1 [6]. In these studies the γ→α transformation in substitutional alloys was considered as a partitionless transformation in view of the slow diffusion of substitutional elements in Fe [9].

The topic of interface mobility has been recently reviewed by comparing all data available in the literature and concluding that the scatter of data is still very large and more work is needed to elucidate the effect of grain size, impurity content, number of transformation steps [11].

A further point of discussion refers to the technique employed to follow the transformation.

Dilatometry appears more straightforward in that the length change is immediately detected by a transducer although the accuracy of the temperature control of the sample can be limited due to its size and position in the measuring device (e. g. it is within 10 K when the instrument is operated at the heating/cooling rate of 20 Kmin -1 [5]). DTA, instead, provides much better control of the sample temperature whereas the heat effect due to the transformation can be smeared due to instrumental thermal lag [12]. The interest in the topic arose from the need of a careful instrumental calibration and interpretation of DTA results obtained with various industrial alloys as a ground support for further studies to be performed in microgravity [13]. This paper reports an analysis of the α↔γ transformation in nominally pure Fe by DTA with attention to the occurrence of the transformation on heating which was not studied either by dilatometry and thermal analysis in previous works. This is preceded and accompanied by a discussion on the calibration of the instrument and on the interpretation of the signal given by first order transformations by comparing the α↔γ transformation of Fe with the melting of pure elements and the α↔β transformation of Ti.

Experimental

Iron samples for thermal analysis were cut from a rod of 4 mm in diameter at various length in order to have sample masses from 250 mg to 2000 mg. The purity was 99.98% (Sigma Aldrich, product no. 266213) and the metal impurity content was detailed by the producer as 2.5 ppm Ti and 0.7 ppm Cu. Rolled foils were also used, made of a laboratory refined iron, with sample masses from 130 mg to 440 mg. The metal impurity content, determined by neutron activation, was (all values in ppm): Cu 3, Mo 2, Co 1, Mn 1, Ni 1, Cr 0.4, traces of W, Ag, Zn were detected.

Thermal analysis was performed with the high temperature differential scanning calorimeter of Setaram. The cell is made of alumina and the sample is contained in an alumina pan with some alumina powder to prevent it from sticking to the crucible walls. The sensor for the sample and reference cell is a thermopile and the instrumental output is a temperature difference. Therefore, the instrument will be referred to as a differential thermal analyser (DTA) The cell is evacuated and purged several times before measuring under flowing Helium filtered to remove residual oxygen and water vapour. The standard thermal cycle was as follows. The specimen was heated from room temperature to 893 K at the rate of 15 Kmin -1 and kept at this temperature for 30 min. Then it was heated at the rate of 10 Kmin -1 to 1293 K and kept at this temperature for 30 min. At last it was cooled at the rate of -10 Kmin -1 to 893 K and kept at this temperature for 30 min. For the various samples iron, one or more cycles reported before, have been performed. Various experiments were performed also at the rates of 0.1, 0.5, 2, 5 Kmin -1 . The DTA cell was calibrated for the temperature and the heat flow by recording the melting temperature and heat of fusion of various pure metal samples (Al, Ag, Au, Cu, Ni) at the rates of 2, 5, and 10 Kmin -1 . The frequency of sampling of the data points was one point every 1 to 3 s according to the heating or cooling rate for a total of 4800 point per file.

The microstructure of the samples was observed before and after thermal cycles by optical microscopy. The ferrite grain boundaries were revealed by etching with 1.5 vol.% Nital solution.

The average grain size was determined with the intercept method.

Results and Discussion

The shape of the signal and the temperature calibration of DTA: Curie point of Fe versus melting of Au

The molar heat capacity of seven samples of pure Fe has been determined on both heating and cooling at the rate of 10 Kmin -1 in the temperature range from 950 K to 1100 K. The results obtained after calibrating the temperature of the instrument with the melting point of various metals on heating, are reported in Figs. 1 and2, respectively. On heating the specific heat increases up to the ferro-to para-magnetic transition, T c , and then drops steadily (Fig. 1). The present data are in good agreement with those reported in [12]. The figures also contain the C p of Fe as computed as a function of temperature from assessed empirical equations: the one proposed by the Scientific Group Thermodata Europe (SGTE) (thin full curve) and widely employed for establishing lattice stability in phase diagram assessments, that suggested in [12] (dotted line) which uses parameters accounting for all the data available in the literature for the magnetic contribution to the specific heat of Fe in addition to those from SGTE, and that obtained in the most recent modelling of thermodynamic properties of Fe by careful consideration of both ferro-and para-magnetic contributions to the specific heat [15]. There is general agreement of data, apart from the temperature region of the transition where experiments provide a more pronounced cusp with respect to models.

On cooling, the specific heat increases on approaching the Curie point but the cusp is not as sharp as on heating. The maximum of the curve occurs at lower temperature with respect to the assessed T c which apparently coincides with the onset of the C p rise. These observations were already made in [12] and led to the de-smearing procedure described in that paper. With this procedure, the experimental curves were shifted and corrected so that they coincided on both heating and cooling.

However, a procedure for the use of the Curie temperature as a calibration point on cooling was not given. This is desirable since the solidification of metals mostly occurs with some undercooling and, therefore, there is no suggested substance for such calibration in the temperature range considered here, but only for lower temperatures [16]. If it is evident that the C p cusp at T c on heating is well suited for calibrating the cell by running it at various rates, the smeared curve on cooling does not provide an immediate indication. In fact, only after de-smearing the C p curve presents a discontinuity at the Curie point on cooling [12]. For as-recorded traces in the absence of temperature calibration, the choice is left between the maximum of the curve and the onset of its rise (Fig. 2).

In this work the Curie point of Fe samples and the melting and solidification point of Au samples were obtained as a function of heating and cooling rates in a set of experiments performed before the temperature of the instrument was calibrated. As expected [START_REF] Speyer | Thermal Analysis of Materials[END_REF], transformations points scaled with the heating rate. The melting and solidification points of Au were chosen among those available for various metals (Al, Ag, Au, Cu, Ni) because the transformation on cooling occurred with very limited undercooling, often practically undetectable with our experimental set up. 3 and4, respectively. They are discussed here, in the context of calibration but the conclusions will be useful for the general interpretation of the DTA signal given by thermodynamic first order transformations.

The DTA trace for the melting of Au (Fig. 3a) has a slightly smeared onset which is usually attributed to inevitable occasional impurities and continues with a linear trend. When melting is over, the trace decays exponentially towards the baseline * . The transformation is finished at the minimum of the curve and the transformed fraction is represented by this part of the trace. The full peak should be considered only for integration to get the enthalpy of transformation. The linear part marks the departure of the reference cell temperature which follow the imposed temperature program from that of the sample cell which remains at the melting temperature of the element. This is better illustrated by the derivative of the trace (Fig. 3b) which departs from zero at the melting onset and displays a plateau after a transient time. The end of melting is signalled by the rapid variation in the derivative. The transient time depends on heating rate and the amount of heat to be supplied. On cooling the peak is reversed (Fig. 4a) with a linear part corresponding to the progress of solidification and an exponential return to the baseline.

The difference in temperature, ∆T, between the sample and the reference cell during melting can be expressed as [20]

R C t t T i Φ ∆ + - Φ = ∆ ) ( ( 2
)
where t is the time, t i the onset time of the peak, ∆C the difference in heat capacity between the phases in the sample (here liquid and crystal). Φ the heating rate, and R the heat resistance of the instrument. After melting is complete at t p , ∆T is given by

F o r P e e r R e v i e w O n l y [ ] [ ] τ / ) ( exp ) ( p i p t t R C t t T - - Φ ∆ + - Φ = ∆ (3)
where τ is the time constant for temperature equilibration between the sample and reference cells.

∆C is computed from known specific heats and the sample mass. The quantity R was determined by means of heat flow calibration of the instrument and ∆T computed for t = t p . Taking the relevant times, t p and t i , from the DTA trace, the time constant, τ, was determined as 27.5 s.

The derivative of the trace (Fig. 4b) displays again a plateau after a transient time which is shorter than that obtained on heating. Since heat is being evolved by the sample, the temperatures of the sample and reference cells diverge. The transient period measures the instrumental time lag for reaching the steady state of the process. It amounts to about 12 s. The shape of the derivative curve is also indicative of the fact that almost no undercooling has occurred. In fact, in the presence of undercooling the heat release would be immediate and the derivative curve would present no plateau but would contain mostly an abrupt step and a long return to the baseline (insert in Fig. 4b).

Having firmly established the parameters which determine the shape of the signal, we consider the Curie point of Fe samples and the melting and solidification point of Au as a function of heating and cooling rates as obtained in the absence of temperature calibration of the instrument.

The deviation of the transformation temperature from the assessed value, ∆T, is reported in Fig. 5.

The points are all aligned on the same straight line if the onset temperature of the specific heat rise due to the Curie transition is taken on cooling. The figure also contains the points corresponding to the maximum of the heat effect for the same transition on cooling. They are clearly misaligned with respect to those obtained on heating. It is concluded that the Curie point can be used for temperature calibration of the DTA cell and that on cooling it must be taken as the onset of the specific heat rise.

On a time base, the distance between the actual maximum of the DSC trace on cooling and T c is a measure of the time for temperature equilibration between the sample and reference cells after a specific heat jump. The corresponding temperature difference is

F o r P e e r R e v i e w O n l y [ ] { } τ / ) ( exp 1 c t t R C T - - - Φ ∆ = ∆ (4)
where t c is the time at which the Curie point is reached. With the time constant determined above, the rise of the signal as well as the decay of the curve after T c on heating are well reproduced.

The analysis of the instrumental parameters performed in this section, will be essential in assessing the results which will be reported in the following.

The α↔γ transformation of Fe by DTA

Various samples of pure Fe were cut from the same bar. One of them was used for metallography in order to determine the grain size which resulted 41 ± 4 µm. A second one was heated in the DTA cell up to 1173 K at 10 Kmin -1 , i. e. a few degrees below the transformation point and immediately cooled to room temperature. The grain size was now 91 ± 15 µm. This is taken as the measure of the grain size of ferrite undergoing transformation to austenite. The DTA traces for the α→γ and γ→α transformations are reported in Figs. 6a and7a, respectively. The traces display an asymmetric peak as in the melting transitions described above with no further discernible feature at the resolution of the drawing. If, however, the derivative is taken through the points (Figs. 6b and7b), it becomes clear that there is no plateau, contrary to Figs. 4b and5b; instead, there are steps in the curve after the initial transient and the return to the baseline takes longer. The γ→α transformation is undercooled about 10 K. The derivative of its trace has the same steps as that obtained on heating, the only difference being the higher peak in the curve due to the faster release of heat on undercooling.

The results reported in Figs. 6 and7 were reproduced with all samples taken from the same bar using heating and cooling rates ranging from 0.5 Kmin -1 to 10 Kmin -1 . A further example is given in 8 for the heating rate of 2 Kmin -1 . Here, the composite nature of the peak is clearly seen both in the DSC trace and in its derivative. Only when the DSC experiment was performed at the heating rate of 0.1 Kmin -1 the peak had a symmetric shape. Its duration spanned approximately 1300 s. On cooling at the same rate, the transformation occurred at the undercooling of 8 K and lasted about 550 s. Its DTA signal had a composite shape as the one shown in Fig. 8. It is apparent that the transformation is not finished at the minimum or maximum of the curves and that its rate varies as a function of transformed fraction. The enthalpy uptake or release for the α↔γ transformation lasts about 30 s for the heating rate of 10 Kmin -1 . This is of the same order as the time constant, τ, and definitely longer than the time needed for detection of a heat effect by the instrument. It is, therefore, concluded that the stepped signal found at various rates is due to a phenomenon occurring in the sample although smeared by instrumental effects.

After a heating and cooling cycle at 10 Kmin -1 , the grains size of ferrite was checked and resulted 276 ± 19 µm. A few samples were cycled three times through the transformation temperature range and their grain size was monitored after each cycle. The grain size ranged from 91 µm to 1100 µm (before the third cycle). In all cases the α↔γ transformation occurred with a composite peak which was more clearly resolved on heating than on cooling. A series of cycles was performed with analogous outcome using a sample of iron sliced from the same bar which was rolled in all directions to halve its thickness and reduce the grain size.

For further series of measurements a sample of high purity iron from another source (see Experimental Part) was used which confirmed once again the above results.

The α↔β transformation of Ti by DTA

In the previous paragraphs it has been shown that the derivative of the DTA trace differs for two transformations of thermodynamic first order, i. e. the α↔γ transformation of Fe and melting The former occurs fully in the solid state whereas the latter involves a liquid phase, moreover the enthalpy of transformation differs of an order of magnitude. In order to prove that the steps found for the α↔γ transformation are not due to instrumental response to a transition involving only solid phases, the same experiment has been performed at various rates with samples of pure Ti which undergoes an α↔β transformation in the solid state at 1155 K. Fig. 9 reports the DTA trace obtained on heating at 10 K•min -1 and its derivative. The DTA trace has apparently a wide linear portion as indicated by the length of the plateau in the derivative and the transformation ends at the minimum of the curve as shown by the abrupt change in sign of the derivative. The time lag of the instrument is close to that already given above.

The result obtained for the α↔β transformation of Ti confirms that the peculiar DTA features found for the α↔γ transformation of Fe are specific to this transformation.

Sample microstructure

The microstructure of the ferrite phase of all samples was observed after each DTA run and the grain size was determined. In all cases ragged grains or sub-grain boundaries were found as shown in Fig. 10 for the average grain sizes of 276 µm (a) and 1100 µm (b). Next to grains grown to large dimension, there are often fine grains apparently recrystallized in a later stage. Similar microstructures were reported previously for pure Fe, Fe-Co and Fe-Mn alloys and low carbon steel [6,7,9]. After one heating and cooling cycle the grain size of ferrite increased from 91 µm to 276 µm and to higher values for more cycles. This should be mainly due to the permanence of the sample at high temperature after the α↔γ transformation since it has been reported that the grain size decreases after the reverse transformation [9]. The ferrite grain size is thought to be representative of the average austenite grain size [6] consistent with direct observation through the temperature range of transformation [7].

The rate of the α↔γ transformation of Fe

The results illustrated by Figs 6 to 8 show that the rate of the α↔γ transformation in nominally pure Fe varies as a function of transformed fraction as already observed by others [6].

The present DTA data do not allow a direct determination of such rate. It has been verified that the time over which the transformation occurs scales linearly with the inverse square root of the heating rate employed for the experiments as expected from consideration of heat resistance models for the DTA instrument [20]. Therefore, the length of the peak is related to the need of either supplying heat to or removing heat from the sample cell. If this occurred at constant rate the peak and its derivative would appear as those due to fusion and solidification of a pure element (Figs. 3 and4).

Since the α↔γ transformation peaks are composite in shape, it is deduced that that the rate of heat uptake or emission by the sample is not constant. Only on average the transformation time should be comparable with the time expressing the length of the peak, e. g. 30 s for the heating rate of 10 Kmin -1 and 90 s for the heating rate of 2 Kmin -1 . The distance over which the transformation occurs must be of the order of the grain size in the samples reported above. Therefore, an average range of transformation rates is estimated by dividing the grain size span (276-1100 µm) for the time span (30-90 s), i. e. the time window in which the anomalies in the transformation rate can be detected.

The range is in between 4•10 -5 ms -1 and 3•10 -6 ms -1 . These values are representative of the rate of movement of the interface. If the interface would move faster, as expected at the beginning of transformation, the DTA signal would not be affected by any variation being determined by the rate of heat transfer. When the rate of transformation decreases, it becomes comparable to that of heat transfer. Actually, at the heating/cooling rate of 10 Kmin -1 the steps in the derivative of the DTA signal are better appreciated on heating than on cooling since in the latter case the transformation occurs in the undercooling regime and is, therefore, faster. At the minimum heating rate used in this work (0.1 Kmin -1 ), the step was not seen on heating since the heat transfer was slow with respect to 13 the progress of transformation; on cooling the step in the peak appeared again. The peak started at 1180 K in the undercooling regime and lasted about 230 s during which heat was spontaneously released from metastable austenite and the transformation occurred independently from the heat transfer. This finding is in agreement with the one reported in [6] where the γ→α transformation could be followed by DTA in the isothermal mode at a temperature of 1180.7 K. Here the composite peak spanned the time of 240 s.

Our results, however, show that the change in rate of the α↔γ transformation of Fe occurs in every sample whereas in [6] it was clearly detected in the dilatometric curve only for an Fe sample of large grain size (439 µm), detected in the derivative of the curve for a sample of 372 µm and not found for grains sizes of 288 µm and 273 µm. These results refer to cooling curves at rates from -5 K/min to -15 K/min. The curves on heating were not commented. The reason for the discrepancy cannot be indicated with certainty. We note that the sample of larger grain size of [6] undercooled less than the others. Accounting for the findings reported above, this can be crucial in detecting the change of transformation rate at the given cooling rates. Moreover, the temperature control within the sample used for dilatometry (rod of 5 mm diameter and 10 mm in length) was reported to be accurate within 10 K [5] whereas the temperature control is more accurate for the small sample used in thermal analysis.

The variation in transformation rate can be due to various causes. It was attributed in [6] to a variation in nucleation rate after the transformation was halted by the local development of strain in the austenite due to lattice misfit with respect of the growing ferrite. One further possibility can be the formation of the ragged microstructure which is observed frequently in all samples. The bursts of nucleation rate referred to above were suggested to cause local grain refinement and grain branching and to induce the presence of sub-grain boundaries [6]. It should also been admitted that, although high purity samples were employed in all these analyses, the inevitable impurities they contain, including interstitials, can be segregated locally at grain boundaries either as solute atoms or even ultra fine precipitates causing pinning of the boundaries during their movement. The average values of the transformation rate estimated in this work are one order of magnitude lower than the average values reported from direct observation of the frames taken in high temperature microscopy of the transformation in a extra-low carbon steel [7] but fall within the range of those collected in the literature for Fe and Fe-alloys [11]. However, the outcome of the present work is that there is a range of transformation rates and that an absolute value can hardly be determined. Therefore, the mobility of the interface should be checked in every case by careful experimentation. The range of values available to date can provide only an indicative guideline.

Conclusions

In this work the α↔γ transformation in nominally high purity Fe was analysed by DTA on heating and cooling. In all cases it occurred with a stepped peak which has been shown to be due to changes in the transformation rate occurring in the course of the transformation itself. In order to reveal this effect a detailed analysis of the instrumental output has been performed determining the time lag for the detection of the signal and for reaching the steady state of transformation by making use also of the derivative of the peak.

The change in rate occurred in all samples considered in this work irrespective of the grain size which was varied from 91 µm to 1100 µm. The range of rates of movement of the interface between the two phases is estimated in between 4•10 -5 ms -1 and 3•10 -6 ms -1 confirming values previously reported.

Possible reasons for this effect are variation in the rate induced by strain built up ahead of the moving interface, local nucleation bursts, pinning by any inevitable impurity. The microstructural evidence can be found in the formation of ragged grain boundaries. calibrated and operated at the heating rate of 10 Kmin -1 in the range from 950 K to 1100 K (bold full curve and error bars) and as computed according to Ref. [12] (dotted curve), Ref. [14] (thin full curve), and Ref. [15] (dashed curve). Fig. 2. The specific heat of Fe as obtained in the present work by high temperature DTA, fully calibrated and operated at the cooling rate of -10 Kmin -1 in the range from 1100 K to 950 K (bold full curve and error bars) and as computed according to Ref. [12] (dotted curve), Ref. [14] (thin full curve), and Ref. [15] (dashed curve). Cp on Heating Cp From Ref [14] Cp From Ref [12] Cp From Ref [15] Cp / J mol

Acknowledgement

  Examples of DTA traces for melting and solidification of Au are shown, together with their

  Fig.8for the heating rate of 2 Kmin -1 . Here, the composite nature of the peak is clearly seen both in

Fig. 1 .

 1 Figure Captions

Fig. 3 .

 3 Fig. 3. a) The DTA trace obtained for the melting of a sample of pure Au used for calibration of the

Fig. 4 .

 4 Fig. 4. a) The DTA trace obtained for the solidification of the sample of pure Au used for

Fig. 5 .

 5 Fig. 5. The difference between the recorded and assessed values of the Curie temperature of Fe and

Fig. 6 .

 6 Fig. 6. a) The DTA trace obtained for the α γ transformation of a sample of pure Fe at the heating rate of 10 Kmin -1 . b) The time derivative of the trace reported in a).

Fig. 9 .

 9 Fig. 9. a) The DTA trace obtained for the α β transformation of a sample of pure Ti at the heating rate of 10 Kmin -1 . b) The time derivative of the trace reported in a).

Fig. 10 .

 10 Fig. 10. a) The microstructure of a sample of pure Fe after one heating and cooling run through the

  Figure 3a

  Figure 10b

http://mc.manuscriptcentral.com/pm-pml Philosophical Magazine & Philosophical Magazine Letters

This work has been performed in the framework of the ESA-Thermolab project "High-Precision Thermophysical Property Data of Liquid Metals for Modelling of Industrial Solidification Processes".


Paolo Cermelli
email: paolo.cermelli@unito.it

Antonio J Di Scala

Constant-angle surfaces in liquid crystals

Keywords: liquid crystals, disclinations (user supplied): interfaces in nematics, disclination cores, layers in smectics

We discuss some properties of surfaces in R 3 whose unit normal vector has constant angle with an assigned direction field. The constant angle condition can be rewritten as an Hamilton-Jacobi equation correlating the surface and the direction field. We focus on examples motivated by the physics of interfaces in liquid crystals and of layered fluids, and examine some properties of the constant-angle surfaces when the direction field is singular along a line (disclination) or at a point (hedgehog). We also show how our results may be used to study the shape of disclination cores.

Introduction

We study here surfaces in R 3 whose unit normal vector forms a constant angle with an assigned direction field. The main idea is to show that the constant angle condition is equivalent to an Hamilton-Jacobi equation correlating the surface and the direction field. We discuss examples motivated by the physics of interfaces in liquid crystals and of layered fluids, and study the properties of constant-angle surfaces when the direction field is singular along a line or at a point. Precisely, we show how constant-angle surfaces may be used to describe interfaces occurring in special equilibrium configurations of nematic and smectic C liquid crystals, and to determine the shape of disclination cores in nematics (for an extensive study of boundaries of nematic drops or tactoids in nematics, in which the constant angle condition is relaxed, see Virga [START_REF] Virga | Drops of nematic liquid crystals[END_REF], [START_REF] Virga | Variational theories for liquid crystals[END_REF], Prinsen and van der Schoot [START_REF] Prinsen | Shape and director-field transformations of tactoids[END_REF], [START_REF] Prinsen | Continuous director-field transformations of nematic tactoids[END_REF], [START_REF] Prinsen | Parity breaking in nematic tactoids[END_REF] and Calderer, Liu and Shen [START_REF] Calderer | Axisymmetric configurations of bipolar liquid crystal droplets[END_REF]).

The Hamilton-Jacobi approach to the determination of constant angle surfaces has been also applied extensively to the shape-from-shading problem (cf. e.g., Lions, Rouy and Tourin [START_REF] Lions | Shape-from-shading, viscosity solutions and edges[END_REF]), and specific numerical methods are now available for the computation of viscosity solutions of Hamilton-Jacobi equations (cf. e.g., Osher and Sethian [START_REF] Osher | Fronts propagating with curvature-dependent speed: algorithms based on Hamilton-Jacobi formulations[END_REF] and Sethian [START_REF] Sethian | Level set methods and fast marching methods[END_REF]).

Precisely, a direction field on a domain Ω ⊂ R 3 can be locally described by a vector field m such that |m| = 1. However, since the vectors m and -m define the same direction, it is customary to assign a direction field through a map (which we continue to denote by the same symbol) m : Ω → P 2 (R), with P 2 (R) the real projective plane, i.e., the set of all lines in R 3 passing through the origin or, equivalently, the set of equivalence classes of unit vectors of the form {m, -m}.

Consider a piecewise smooth surface S ⊂ Ω whose unit normal n has a constant angle α with m:

|n • m| = cos α, (1)
for fixed α ∈ [0, with A α = A α (x) = cos 2 α1m(x) ⊗ m(x). Writing H(x, p) = p • A α (x)p for the Hamiltonian, it is well known that local solutions of (2) may be obtained by the method of characteristics, which are the solutions of the canonical system (see Arnold [START_REF] Arnold | Lectures on partial differential equations[END_REF])

ẋ = ∂H ∂p , ṗ = - ∂H ∂x . ( 3
)
When α = 0 the method of characteristics fails, since the right-hand side of (3) vanishes identically for H(x, p) = 0. Indeed, when α = 0, A α is positive semidefinite, and the equation H(x, p) = 0 leads to an overdetermined system that must be studied by Frobenius theory. Thus, solutions of (1) for α = 0 exist if and only if m

•curl m = 0 in Ω.
In what follows we restrict to α ∈ (0, π/2) when using the Hamilton-Jacobi approach, with the exception of Section 3.1 where the integrability condition is automatically satisfied and complex-variable methods can be applied.

In this work we study [START_REF] Virga | Variational theories for liquid crystals[END_REF] for some special choices of the direction field m. Namely, we characterize (i) cylindrical and conical surfaces arising when m is singular along a straight line and orthogonal to that line (line defects in nematics), and (ii) conical surfaces corresponding to a constant orientation field m = m 0 (layered structures in smectics C).

We now discuss briefly two physical examples which lead to equations of the form (1).

Applications to nematics

A nematic liquid crystal is an ordered fluid whose constituents are macromolecules which locally tend to align parallel to each other (cf., e.g., De Gennes and Prost [START_REF] De Gennes | The physics of liquid crystals[END_REF]). A nematic liquid crystal can be described by a direction field m that measures the average local orientation of the macromolecules. In the simplest approximation, the field m at equilibrium is a stationary point of the (one-constant) Frank energy functional

Ω κ 2 |∇m| 2 dv, (4)
with κ > 0, whose Euler equation is the harmonic mapping equation (see Eells and Sampson [12] where (5) is introduced and studied in a general setting)

m + |∇m| 2 m = 0. (5)
Often, the direction field describing the orientational order in a nematic fluid is not defined on the whole Ω, the most common singularities being straight disclinations and hedgehogs. Concisely, a (straight) disclination is a line in R 3 such that m is not defined on but is smooth in Ω \ { }. Restricting to situations such that m is a plane vector field orthogonal to , a disclination is characterized by its integer strength k, i.e., the degree of the mapping m : γ → P 1 (R), with γ a loop in Ω \ { } encircling once. Roughly, k is the net number of π-turns of m when γ is traversed once (cf. Section 3).

A hedgehog is a point defect in R 3 such that m is not defined at some point x 0 , but is smooth in Ω \ {x 0 }. The simplest example of hedgehog, and the one we shall study here, is m(x) = x/|x| for x = 0.

When singularities are present, it is usually assumed that, at equilibrium, m satisfies (5) away from the singularity.

Equations of the form (1) arise in the study of free surfaces of nematic liquid crystals and interfaces between nematic and isotropic phases. Indeed, denoting by n : ∂Ω → S 2 the outward unit normal to ∂Ω, which we view as the interface between the nematic and isotropic phase, it is often assumed that the strong anchoring boundary condition holds at ∂Ω with α ∈ (0, π/2), which is an equation of the form [START_REF] Virga | Drops of nematic liquid crystals[END_REF].

|n • m| = cos α,
In Sections 2 and 3 we restrict to two special cases: first, assuming that m is a constant solution of (5), we show that surfaces satisfying (1) are flat and ruled, and can be explicitly described. Second, assuming that m has a line singularity with arbitrary integer strength, we study cylindrical surfaces which satisfy (1), and discuss how their topology depends on the strength.

However, since in general surfaces solutions of (1) are not bounded, they cannot be used to describe the boundaries of nematic drops (tactoids). Actually, (1) must be relaxed by requiring that the director field at the surface ∂Ω satisfies weak anchoring conditions, i.e., that m and ∂Ω minimize the total energy, given by the sum of the Frank energy (4) and a surface energy term of the form (cf. Rapini and Papoular [START_REF] Rapini | Distorsion d'une lamelle nématique sous champ magnétique. Conditions d'ancrage aux parois[END_REF])

∂Ω (ψ 0 + ψ(n • m)) da ( 6
)
where ψ 0 > 0 is a fixed constant and ψ is the anchoring energy, which is assumed to have absolute minima at n • m = ± cos α.

The energy-minimizing approach to determine the optimal shape of tactoids has been extensively developed by Virga [START_REF] Virga | Drops of nematic liquid crystals[END_REF], [START_REF] Virga | Variational theories for liquid crystals[END_REF] and, more recently, by Prinsen and van der Schoot [START_REF] Prinsen | Shape and director-field transformations of tactoids[END_REF], [START_REF] Prinsen | Continuous director-field transformations of nematic tactoids[END_REF], [START_REF] Prinsen | Parity breaking in nematic tactoids[END_REF] and Calderer, Liu and Shen [START_REF] Calderer | Axisymmetric configurations of bipolar liquid crystal droplets[END_REF].

In this paper we do not discuss the optimal shape of nematic drops, but we make use of the energyminimizing approach to study the shape of the core of a straight disclination. In an oversimplified model, the core of a defect may be viewed as an isotropic inclusion in a nematic environment: in this case, since the inclusion has a volume which is negligible with respect to the environment, scaling arguments show that we may assume that the nematic field is a given singular solution of (5) outside the defect core, whose boundary is a minimizer of the surface energy (6) (cf. the related work by Rosso and Brunelli [START_REF] Rosso | Forces on nematic disclinations with optimal core[END_REF], who determine the optimal shape of a disclination core neglecting surface energy).

Applications to smectics C

Smectic C liquid crystals are ordered fluids characterized by a layered structure (cf. De Gennes and Prost [START_REF] De Gennes | The physics of liquid crystals[END_REF]). When density fluctuations are neglected, the description of equilibrium configurations of smectics requires, in addition to the orientation field m, a scalar order parameter f : Ω → R whose level sets define the layers. In smectics C, the director field m is constrained to a constant angle α with the layer normal n = ∇f /|∇f |, i.e.,

|m • n| = cos α, with α ∈ (0, π/2).
The energy of a smectic C liquid crystal may be written as the sum of a nematic energy term, which in the simplest approximation is the Frank energy (4), and a smectic term of the form (cf. [START_REF] De Gennes | The physics of liquid crystals[END_REF], [START_REF] Biscari | Telephone-cord instabilities in thin smectic capillaries[END_REF])

Ω A|∇f -qm| 2 + B|m • ∇f -q| 2 dv, ( 7
)
where A, B, q are positive constants.

An energetically-based approach to the determination of layered structures in smectics, which allows the director field to vary, is developed by Kidney, McKay and Stewart [START_REF] Kidney | Energy considerations for parabolic cyclides in SmA liquid crystals[END_REF], McKay and Leslie [START_REF] Mckay | A continuum theory for smectic liquid crystals allowing layer dilation and compression[END_REF], and McKay [START_REF] Mckay | Spherical, cylindrical and toroidal layering of smectic C liquid crystals[END_REF]. We restrict here to special configurations such that:

• the orientation field m is assigned; • the smectic order parameter is an absolute minimizer of the smectic energy density [START_REF] Lions | Shape-from-shading, viscosity solutions and edges[END_REF] subject to the constraint [START_REF] Virga | Drops of nematic liquid crystals[END_REF].

Letting P = 1 -m ⊗ m, then |∇f | 2 = |P ∇f | 2 + (m • ∇f ) 2 and
|P ∇f | 2 = tan 2 α(m • ∇f ) 2
, so that the integrand of ( 7) may be rewritten as

(A + B + A tan 2 α)(m • ∇f ) 2 -2q(A + B)(m • ∇f ) + q 2 (A + B),
and minimization with respect to x = (m • ∇f ) yields the system of PDEs

|P ∇f | = γ tan α, m • ∇f = γ, (8)
with γ = (A+B)q (A+B)+A tan 2 α . We show below that, contrary to (2), ( 8) has no solution for a large class of orientation fields, and discuss two cases when the solution does exist, i.e., for a uniform orientation field, and a line disclination with strength 2.

Constant-angle conical surfaces in a uniform orientation field

For a uniform orientation field m = m 0 on R 3 , (1) reduces to

|n • m 0 | = cos α, ( 9
)
and we assume that α ∈ (0, π/2).

Let S be a solution of ( 9): notice first that smooth portions of S ⊂ S are flat, since the Gaussian curvature K identically vanishes on these portions. By consequence, S is also ruled.

To see this, let S ⊂ S be smooth, and notice that ( 9) implies that the image of S by the Gauss map n : S → S 2 is a circle. By consequence, the Gaussian curvature vanishes, since

K = lim area(S )→0 area(n(S )) area(S ) = 0,
and S is flat.

To solve [START_REF] Sethian | Level set methods and fast marching methods[END_REF], choose Cartesian coordinates (x, y, z), with associated basis (i, j, k) such that k = m 0 . Assume that S is smooth and is locally the graph of a function g

:= D ⊂ R 2 → R, i.e., S = {(x, y, z) / z = g(x, y)}. Since n = (-∇g + k)/ 1 + |∇g| 2 , with ∇g = (∂ x g)i + (∂ y g)j, (9) reduces to the eikonal equation |∇g| = tan α. ( 10
)
Let C 0 = S ∩ {z = 0} be the trace of S on the (x, y)-plane. Then S is a piece of a conical surface generated by C 0 . To see this, let x 0 (s) be a parametrization of C 0 with arc parameter s, and let τ = τ (s) and ν = ν(s) be the unit tangent and normal vectors to C 0 . Writing x = xi + yj, let also t = t(x) := dist(x, C 0 ) be the distance of a point x from C 0 (which is well defined in a neighborhood of the curve if C 0 is smooth), and notice that the solution of [START_REF] Arnold | Lectures on partial differential equations[END_REF] with boundary condition g(x, y) = 0 on C 0 is

g(x, y) = t(x, y) tan α. ( 11
)
In fact, consider the local change of coordinates ((s, t) are known as Fermi coordinates)

x(s, t) = x 0 (s) + tν(s), Since ∇g•ν = ∂ t ĝ, and ∇g•τ = ∂ s ĝ/(1-ct), with c the curvature of C 0 , it follows, by taking ĝ(s, t) = t tan α, that ∇g = tan α ν, so that (10) is satisfied.

As is well known, ( 11) is a local solution of [START_REF] Arnold | Lectures on partial differential equations[END_REF], since the distance is in general a multi-valued function. However, the method of characteristics actually yields a surface S in R 3 which is a global solution of the original equation ( 9), given by the global parametrization

x 0 (s) + tν(s) + (t tan α) m 0 , ( 12
)
in terms of the Fermi coordinates (s, t) as t varies in R. The surface S parametrized by ( 12) may have self-intersections, for instance at points (s, t) corresponding to the caustics of [START_REF] Arnold | Lectures on partial differential equations[END_REF]. In Figure 1 a portion of a conical surface with C 0 an ellipse is shown: notice the sharp edge where the generatrices meetcorresponding to the caustics of the eikonal inside the ellipse.

Cylindrical surfaces in the presence of a straight disclination

Let Ω = D × R be a cylindrical domain, with D ⊂ R 2 , and consider an orientation field m : Ω → P 2 (R). Assume that m is planar and independent of the vertical coordinate, so that

m = m(x), m • k = 0, (13)
with x = xi + yj ∈ D, and (i, j, k) an orthonormal basis with k orthogonal to the plane containing D.

Under the above assumptions the problem is 2-dimensional and we may view m :

D → P 1 (R).
Since m is planar, we may write

m(x) = cos ϕ(x) i + sin ϕ(x) j.
Notice that, for nematic orientation fields, m satisfies (5) if and only if ϕ is harmonic:

∆ϕ = 0, (14)
with ∆ the laplacian in R We first look for surfaces S ⊂ Ω such that

|n • m| = 1, (15)
with n a choice of unit normal for S. Since m•curl m = 0, Frobenius theorem is trivially satisfied and local solutions to (15) exist. Moreover, the same argument of Section 2 shows that the corresponding surfaces are flat and ruled, actually cylinders S = C × R with vertical axis, with C ⊂ D a curve in the plane orthogonal to k.

Noting that the plane curve C has unit normal n, it follows that C is still defined by [START_REF] Biscari | Telephone-cord instabilities in thin smectic capillaries[END_REF], which we henceforth view as an equation in the unknown curve C.

Equation ( 15) defines the family of (orthogonal) integral curves of the direction field m. Describing C as a level set of a function f : D ⊂ R 2 → R:

C = {x ∈ D : f (x) = c}, with c ∈ R, (15) is equivalent to ∇f × m = 0. ( 16
)
Assume first that D is simply connected and m is a smooth solution of (5) in Ω, so that no defect is present and ϕ is harmonic in D. Denote by w the harmonic conjugate of ϕ, solution of the Cauchy-Riemann equations

∂ϕ ∂x + ∂w ∂y = 0, ∂ϕ ∂y - ∂w ∂x = 0,
and define u = e -w m = e -w (cos ϕi + sin ϕj), [START_REF] Mckay | A continuum theory for smectic liquid crystals allowing layer dilation and compression[END_REF] so that u is a non-zero vector field on D. The function e -w is an integrating factor for m, and it is easy to see that Div u = 0 and Curl u = 0, where Curl u = u 2 xu 1 y for u = u 1 (x, y)i + u 2 (x, y)j. By consequence, ∆u = 0 in D. In particular, there exists f : D → R such that ∇f = u, [START_REF] Mckay | Spherical, cylindrical and toroidal layering of smectic C liquid crystals[END_REF] so that n = ∇f /|∇f | = m and {f (x) = c} is an orthogonal integral curve of the direction field m, i.e., is a solution of [START_REF] Biscari | Telephone-cord instabilities in thin smectic capillaries[END_REF]. Equivalently, identify R 2 to C and, for z = x + iy, define

m(z) = e -iϕ(z) and u(z) = e -(w(z)+iϕ(z)) , ( 19
)
where w(z) and ϕ(z) are w(x, y) and ϕ(x, y) expressed in terms of the complex variable z. By construction u is holomorphic in D.

Consider now a holomorphic function in

D F = f + ig, ( 20
)
with f and g real valued. Since F is holomorphic, Comparing ( 18) and ( 19) we conclude that f , as a function on D ⊂ R 2 , satisfies ( 15) and ( 16) if and only if F as defined by [START_REF] Bonnans | Numerical optimization[END_REF] satisfies

F = 2 ∂f ∂z = ∂f ∂x -i ∂f ∂y .
F (z) = λ(z)u(z),
for some real function λ(z) on D. However, since F and u are holomorphic, λ is a real constant, so that (15) and ( 16) reduce (modulo a multiplicative real constant) to

F (z) = u(z). ( 21
)
Therefore, the plane curves solutions of ( 15) are the level sets of the real part f of any primitive F of u, i.e.,

C = {z ∈ C : Re(F (z)) = c}. ( 22
)
Going back to three dimensions, assume that the director field m is singular along a line parallel to the axis of the cylinder Ω: precisely, assume that is a straight disclination of strength k. Assuming that intersects D at x = 0, we may henceforth restrict to the plane direction field m : D \ {0} → P 1 (R), which is singular at 0. Recall now that the strength k ∈ Z of the disclination at 0 is defined as the net number of π-rotations of m along any simple closed loop encircling 0, given by the relation

k 2 = 1 2π γ (-m 2 ∇m 1 + m 1 ∇m 2 ) • dx = 1 2π γ dϕ.
Fix now k ∈ Z: any solution of ( 14) with given strength k may be written in the form

ϕ( , ϑ) = k 2 ϑ + g( , ϑ), (23)
with g : D \{0} → R harmonic (and single-valued) and ( , ϑ) polar coordinates centered at {0}. Therefore, by means of a conformal transformation 1 , we may restrict to director fields that have the form

m = cos kϑ 2 i + sin kϑ 2 j. ( 24
)
Switching to complex-variables, and letting ϕ = kϑ/2 in ( 19) we obtain

u(z) = z -k 2 . ( 25
)
Now recall that the orthogonal integral curves of m are given by ( 22), where F is a primitive of u,

F (z) = z -k 2 .
which yields

F (z) = 2 2-k z 2-k 2 , k = 2, log z, k = 2.
f (z) = Re z 2-k 2 = 1 2 z 2-k 2 + z 2-k 2 , k = 2, log |z|, k = 2.
To avoid discussing the sign of c, consider the curves C defined by |f (z)| = c > 0. The topology of C depends on k as follows:

• For k ≤ 1, k = 0, the curve |f (z)| = c is not bounded, does not contain the defect z = 0 and has 2k connected components with 2k asymptotic directions (cf. Figure 2). In polar coordinates C has the explicit representation

= c cos 2-k 2 ϑ 2 2-k . ( 26
)
• For k = 2, the field m is radial and the curves |f (z)| = c are circles centered at 0.

• For k ≥ 3 the curve C contains the defect at z = 0. C is bounded, and has a clover-like shape with k -2 leafs. In polar coordinates C has the explicit representation

= c cos k-2 2 ϑ 2 k-2 . ( 27
)
Figure 2. From left to right: (i) a section of a cylindrical constant-angle surface for α = 0, with a disclination of strength k = 1 located along the axis of the cylinder (directed out of the page); (ii) same as before with k = -1;

(iii) a section of a cylindrical domain K with boundary a constant angle surface: a bump for k = 1 and (iv) a meniscus for k = -1.

Surfaces with α arbitrary

The above results may readily be extended to surfaces S ⊂ R 3 such that

|m • n| = cos α, ( 28
)
with α ∈ (0, π/2). Again, solutions of ( 28 Consider now, as in Section 3.1, a straight disclination with strength k parallel to the axis of the cylinder, and the corresponding singular plane direction field m as in (24). In this case, (25) takes the form

u(z) = e ±iα z -k 2 ,
so that

F (z) = 2 2-k e ±iα z 2-k 2 , k = 2, e ±iα log z, k = 2.
and it follows that the curve C has the form f (z) = c, with

f (z) = 2-k 2 cos 2-k 2 θ ± α , k = 2, log ± θ tan α, k = 2.
The curves C defined by |f (z)| = c > 0, are now such that

• For k = 2, the curve C is a rigid rotation (of angle 2 2-k α) of the curves defined by ( 26) and ( 27). • For k = 2,the curves |f (z)| = c are logarithmic spirals of the form = ce ±θ tan α .

Disclination pairs

The above argument may be extended to any number N of parallel dislocations with strength k i , intersecting the horizontal plane at z i , i = 1, . . . , N (complex representation z i ). Modulo a conformal transformation, at equilibrium the nematic angle may be written in the form

ϕ = ϕ 0 + 1 2 N i=1 k i ϑ i ,
with ϑ i (x) = arctan y-yi x-xi a polar angle centered at z i and ϕ 0 a constant angle. Proceeding as in the above section, we find in the complex representation

u(z) = e -iϕ0 N i=1 (z -z i ) -k i 2 .
To be specific, consider two parallel disclinations of strength k 1 = k 2 = 1, located at z 1 = 0 and z 2 = -ai, with a > 0 (complex representation z 1 = 0 and z 2 = -a), and assume that the nematic angle is given by

ϕ = 1 2 ϑ + 1 2 θ, ( 29
)
where θ(x) = arctan y x+a . Then

u(z) = z -1 2 (z + a) -1 2 , ( 30
)
and a primitive F of u is easily calculated: From left to right: (i) a section of cylindrical orthogonal surface surrounding a pair of parallel disclinations with strength 1/2 and mutual separation a = 2; (ii) a section of a capillary tube with shape a constant-angle surface: for orthogonal (homeotropic) anchoring at the walls an equilibrium configuration has two disclinations of strength k = 1.

F (z) = log z + a 2 + z(z + a) .
It follows that the orthogonal curves of the direction field with angle (29) are given by

z + a 2 + z(z + a) = c, (31)
with c a real positive constant. The explicit form of this curve in polar coordinates follows from a straightforward calculation, and will be omitted here. Figure 3 shows the plot of one such curve surrounding a pair of disclinations with separation a = 2.

Disclination cores in nematics

Consider a singular orientation field m as in (24). Domains whose boundaries are constant-angle surfaces correspond to natural equilibrium interfaces between nematic and isotropic phases, since the orientation field inside the domain is by construction a solution of (5), and strong anchoring conditions are automatically satisfied at the interface. Specifically, consider a cylindrical domain K = D × R ⊂ D × R, with C = ∂D ⊂ C a portion of a constant-angle surface. From the discussion in Section 3.1 it follows that

• For k ≤ 1, K is not bounded. This class of domains corresponds to curved bounding surfaces of domains containing a disclination (k = 1), or to menisci connecting films containing a disclination (k = -1)(Figure 2). When k = 0 we obtain a film in which the orientation field is aligned. • For k = 2, K is either a circular cylindrical filament centered at the disclination, or its complement, which corresponds to a cylindrical isotropic core about a straight disclination. • For a pair of disclinations of strength k = 1, the domain K is a filament with elongated section containing the disclination lines. This configuration may be realized for instance in a capillary tube with elongated section and orthogonal anchoring at the walls (Figure 3).

A second issue that can be discussed in this contest is the equilibrium shape of disclination cores. As mentioned above, in a crude approximation we may view the core as an isotropic inclusion in a nematic environment. For very small sectional areas of the core, the nematic field can be assumed to be a fixed equilibrium solution of (5) outside the core, and the shape of the core itself is a minimizer of the interfacial energy [START_REF] Calderer | Axisymmetric configurations of bipolar liquid crystal droplets[END_REF].

Precisely, assume that the nematic angle ϕ is assigned, a fixed solution of ( 14) of the form (24). Let D be the planar section of the disclination core, write t for a choice of unit tangent to ∂D , and consider the surface energy functional (per unit height of the cylinder) with ω, ψ 0 > 0 and s the arc-parameter on the curve ∂D . Clearly, orthogonal anchoring is preferred at the core-boundary. Consider the following problem: find the plane region D , containing the origin, that minimizes (32) under the constraint that area(D ) = 1. Noting that m = ∇f /|∇f |, (32) may be written as

∂D ψ 0 + ω(m • t) 2 ds, (32)
∂D ψ 0 + ω |∇f | 2 ∂f ∂s 2 ds.
We only discuss the cases k = ±1. In polar coordinates centered at the origin (the trace of the disclination), let = R(ϑ) be the parametric representation of the surface f = const., i.e., ∂D , and notice that we may restrict attention to curves that are symmetric with respect to the x-axis. Then the minimization of (32) is equivalent to the minimization of

π 0 AR 2 -2CRR + BR 2 √ R 2 + R 2 dϑ, ( 33
)
on the set of all smooth functions R = R(ϑ) ≥ 0 such that

π 0 R 2 dϑ = 1, (34)
and where

A = ψ 0 + ω cos 2 ϑ 2 , B = ψ 0 + ω sin 2 ϑ 2 , C = sign(k)ω cos ϑ 2 sin ϑ 2 .
We have performed a numerical minimization of (33) as described in the Appendix. The results are shown in Figure 4: planar sections of disclination cores have elongated oval shapes, whose aspect ratio depends on ω. As ω → 0 the planar section of the core tends to be a disk, while as ω increases the core tends to flatten along the x-axis (Figure 5). We can use the same procedure to study the shape of the cores of a pair of disclinations located at z = 0 and z = -ai with strength k = 1: the energy functional of the surface surrounding the disclination at z = 0 has still the form (33), but now

A = ψ 0 + ω 2 1 + a cos ϑ + (cos 2 ϑ -sin 2 ϑ)R(ϑ) a 2 + R 2 (ϑ) + 2aR(ϑ) cos ϑ , B = ψ 0 + ω 2 1 - a cos ϑ + (cos 2 ϑ -sin 2 ϑ)R(ϑ) a 2 + R 2 (ϑ) + 2aR(ϑ) cos ϑ , C = ω 2 a sin ϑ + 2R(ϑ) cos ϑ sin ϑ a 2 + R 2 (ϑ) + 2aR(ϑ) cos ϑ .
The results of the numerical minimization of (33) are shown in Figure 6: again sections of disclination cores have elongated oval shapes along the line joining the defecs.

Applications to smectics C

For layered structures such as smectics C, the constraint (8) 2 poses severe restrictions to the existence of solutions of the Hamilton-Jacobi equation (8) 1 . Also, while for nematics the basic unknown is the surface S and f is just an auxiliary function of which S is a level set, in smectics the function f is an order parameter that should be uniquely determined by [START_REF] Osher | Fronts propagating with curvature-dependent speed: algorithms based on Hamilton-Jacobi formulations[END_REF] and suitable boundary conditions. Indeed, fix a direction field m as before and consider the local coordinate system (s, ξ) on R 3 , with ξ a coordinate system on a surface S 0 transversal to m, and s the arc parameter along the integral curves of m, measured from the intersection of the corresponding integral curve and the initial surface S 0 . In other words, s = s(x) and ξ = ξ(x) are defined so that ξ(x) is the point of intersection of the integral curve of m passing through x and the surface S 0 , and s is the length of the arc of this curve between ξ and x.

Writing f = f (ξ, s), it follows that (8) 2 is equivalent to

∂ f ∂s = γ, so that there exists g = g(ξ) such that f (ξ, s) = γ(s + g(ξ)), (35)
and using the fact that m • ∇s = 1, (8) 1 can be rewritten as

|∇ ξ g| 2 + ∇s • ∇ ξ g + |∇s| 2 -1 = tan 2 α. ( 36
)
We discuss here a few special cases of (36). Consider first a constant orientation field m = m 0 in a cylindrical domain D × R, with coordinates (x, y, z) with z along the axis of the cylinder. By (35) we may write the solution of (8) in the form

f (x, y, z) = γ(z -g(x, y)), ( 37
)
and a direct calculation shows that (36) reduces to the two-dimensional eikonal equation [START_REF] Arnold | Lectures on partial differential equations[END_REF] for g. Assume now, to fix ideas, that the boundary condition g(x, y) = 0 on ∂D (38) holds: arguing as in Section 2, we obtain a local solution [START_REF] De Gennes | The physics of liquid crystals[END_REF] in the form f (x, y, z) = γ(zt(x, y) tan α), with t(x, y) the distance function with respect to ∂D.

As discussed in section 2, the distance function t(x, y) is in general multivalued, with branching at the caustics, and a selection rule is necessary, for instance by requiring that g be a viscosity solution of the two-dimensional eikonal equation (36) (cf. e.g., Lions [START_REF] Lions | Generalized Solutions of Hamilton-Jacobi Equations[END_REF] and Lions, Rouy and Tourin [START_REF] Lions | Shape-from-shading, viscosity solutions and edges[END_REF]). Figure 1 shows the viscosity solution corresponding to a smectic layer in a cylindrical tube with elliptic cross-section. The smectic layers are equidistant portions of conical surfaces obtained by translating the surface in Figure 1 along the vertical axis.

As a second example, consider a singular direction field m corresponding to a disclination line with strength k = 2, as in Section 3. Specifically, let m = e r , with (r, ϑ, z) cylindrical coordinates with z along the disclination. Then s = r, ξ = (ϑ, z) and ∇ ξ g = (∂ z g)e z + 1 r (∂ ϑ g)e ϑ , so that (36) yields the eikonal equation on the cylinder

(∂ z g) 2 + 1 r 2 (∂ ϑ g) 2 = tan 2 α. ( 39
)
Since g is independent of r, ∂ ϑ g ≡ 0 and As before, g is independent of R, but now this implies that the above equation has no solutions for α = 0. Therefore, there is no smectic C configuration for which the orientation field has a hedgehog defect and the layers are constant-angle surfaces.

L (k) = ∇ 2 R F (R (k) ) + λ (k) ∇ 2 R c(R (k)
). The computation of the Hessian of F is straightforward but cumbersome and its explicit expression will not given here.

The genetic algorithm proceeds as follows: first, restrict to the class of polynomial radial functions on [0, π] of the form R(ϑ) = N i=1 a i ϑ i-1 , so that on this class the surface energy functional (33) reduces to a function F (a 1 , . . . , a N ) of the N coefficients a i . This function can be evaluated numerically, for given a i , by approximating the integral in (33) by a Newton-Cotes method. The resulting finite-dimensional minimization problem is then solved numerically in terms of the a i , by the continuous genetic algorithm discussed in [START_REF] Haupt | Practical genetic algorithms[END_REF], by taking into account the constant-area constraint by renormalizing the coefficients at every generation.

We have checked the results obtained by the Newton method with the results obtained by the genetic algorithm, and they turned out to be consistent, in that the computed minimum surface energies differ by 0.5% after 10 iterations of the Newton method and 5000 generations of the genetic algorithm (for ω = 0.7 and ψ 0 = 1).

π 2 ]

 2 . Viewing S as a level set of a Lipschitz continuous function f : Ω → R, S = {x ∈ Ω / f (x) = const.}, and letting n = ∇f /|∇f |, (1) is equivalent to the quadratic Hamilton-Jacobi equation ∇f • A α ∇f = 0

  , t) := g(x(s, t)).

Figure 1 .

 1 Figure 1. Side and top views of a portion of a conical constant-angle surface with vertical uniform orientation field m = m 0 . The plane curve generating the cone is an ellipse.

2 . 6 3. 1

 261 Surfaces orthogonal to m: α = 0

  curve C has the form f (z) = c, with

  ) are cylinders S = C × R with vertical axis, with C ⊂ D a plane curve solution of (28) in D.Since (28) is equivalent to | m • n| = 1, where m = cos(ϕ ± α)i + sin(ϕ ± α)j, all considerations in Section 3.1 hold, with u(z) = e -(w(z)+i(ϕ(z)±α)) .

10 Figure 3 .

 103 Figure 3.From left to right: (i) a section of cylindrical orthogonal surface surrounding a pair of parallel disclinations with strength 1/2 and mutual separation a = 2; (ii) a section of a capillary tube with shape a constant-angle surface: for orthogonal (homeotropic) anchoring at the walls an equilibrium configuration has two disclinations of strength k = 1.

0. 5 1Figure 4 .

 54 Figure 4. Left: section of a cylindrical surface minimizing the interfacial energy (33) for k = 1. Right: core of a disclination of strength k = 1 parallel to the axis of the cylinder. The surface-energy constants are ψ 0 = 1 and ω = 0.7. Also, area(D ) = 1.

Figure 5 .

 5 Figure 5. Sections of disclination cores for different values of ω, with area(D ) = 1.

0. 5 1Figure 6 .

 56 Figure 6. Left : (i) section of a cylindrical surface surrounding one of the two disclinations (the rightmost one) of a pair with strength k = 1 and separation a = 4, and minimizing the surface energy. Right: rescaled view of both disclination cores surrounded by the nematic field. The surface-energy constants are ψ 0 = 1 and ω = 0.7, and area(D ) = 1.

g 14 eikonal equation on the sphere 1 R 2 (∂ ϑ g) 2 + 1 R 2

 141212 = γs + zγ tan α + const. (40)Hence, the layers are families of equidistant half-cones with axis the disclination line.Consider finally a hedgehog defect described by the orientation field m = e R , with (R, ϑ, ψ) spherical coordinates in R 3 . Here s = R, ξ = (ϑ, ψ) and ∇ ξ g = 1 R (∂ ϑ g)e ϑ + 1 R sin ϑ (∂ ψ g)e ψ , so that (36) yields the sin 2 ϑ (∂ ψ g) 2 = tan 2 α.(41)

  k) = c(R (k) ), A (k) = ∇ R c(R (k) ), G (k) = ∇ R F (R (k) )), and

Page 24 of 27 http

	November 9, 2006 10:14 Philosophical Magazine cermelli-discala-revised
	Philosophical Magazine & Philosophical Magazine Letters
	13
	F o r
	P
	e
	e r
	R
	e
	v i e
	w
	O n l
	y

://mc.manuscriptcentral.com/pm-pml

http://mc.manuscriptcentral.com/pm-pml Philosophical Magazine & Philosophical Magazine Letters

Acknowledgments P.C. was supported by the Italian M.I.U.R. research grant "Modelli matematici per la scienza dei materiali" (2005).

Appendix

The numerical minimization of the surface energy functional (33) has been performed by two different techniques, which yield compatible results.

For the single-disclination problem, we have used both a modified Newton method for constrained minimization [START_REF] Bonnans | Numerical optimization[END_REF], and a continuous genetic algorithm [START_REF] Haupt | Practical genetic algorithms[END_REF]. For the disclination-pair problem, we have used a genetic algorithm.

To apply the modified Newton method, we may proceed as follows: first, subdivide the interval (0, π) in N -1 subintervals of length h = π/(N -1), and discretize R and R on the resulting mesh of N points by centered finite differences (or, equivalently, approximate R by piecewise-linear functions). Let R = (R 1 , . . . , R N ) be the vector of discretized values of R: computing the integral in (33) and (34) by a Newton-Cotes method, we reduce to the problem of minimizing a function F (R) subject to the constraint c(R) = 0, where

Local minimizers of F may be computed by a Newton method which traces the values of the Lagrange multiplier associated with the constraint, as discussed in [START_REF] Bonnans | Numerical optimization[END_REF]. Specifically, a sequence R (k) approximating the minimizer, and a sequence λ (k) approximating the Lagrange multiplier are constructed by the iterative formulas R

, and the Newton steps (d (k) , λ Q (k) ) are found by solving the linear system