A local heat transfer analysis of lava cooling in the atmosphere: application to thermal diffusion-dominated lava flows

Journal of Volcanology and Geothermal Research 81 Ž1998. 215–243

A local heat transfer analysis of lava cooling in the atmosphere: application to thermal diffusion-dominated lava flows Augusto Neri

)

Consiglio Nazionale delle Ricerche, Gruppo Nazionale per la Vulcanologia, Centro di Studio per la Geologia Strutturale e Dinamica dell’Appennino Via S. Maria 53, I-56126 Pisa, Italy Received 20 May 1997; accepted 7 November 1997

Abstract The local cooling process of thermal diffusion-dominated lava flows in the atmosphere was studied by a transient, one-dimensional heat transfer model taking into account the most relevant processes governing its behavior. Thermal diffusion-dominated lava flows include any type of flow in which the conductive–diffusive contribution in the energy equation largely overcomes the convective terms. This type of condition is supposed to be satisfied, during more or less extended periods of time, for a wide range of lava flows characterized by very low flow-rates, such as slabby and toothpaste pahoehoe, spongy pahoehoe, flow at the transition pahoehoe-aa, and flows from ephemeral vents. The analysis can be useful for the understanding of the effect of crust formation on the thermal insulation of the lava interior and, if integrated with adequate flow models, for the explanation of local features and morphologies of lava flows. The study is particularly aimed at a better knowledge of the complex non-linear heat transfer mechanisms that control lava cooling in the atmosphere and at the estimation of the most important parameters affecting the global heat transfer coefficient during the solidification process. The three fundamental heat transfer mechanisms with the atmosphere, that is radiation, natural convection, and forced convection by the wind, were modeled, whereas conduction and heat generation due to crystallization were considered within the lava. The magma was represented as a vesiculated binary melt with a given liquidus and solidus temperature and with the possible presence of a eutectic. The effects of different morphological features of the surface were investigated through a simplified description of their geometry. Model results allow both study of the formation in time of the crust and the thermal mushy layer underlying it, and a description of the behavior of the temperature distribution inside the lava as well as radiative and convective fluxes to the atmosphere. The analysis, performed by using parameters typical of Etnean lavas, particularly focuses on the non-intuitive relations between superficial cooling effects and inner temperature distribution as a function of the major variables involved in the cooling process. Results integrate recent modelings and measurements of the cooling process of Hawaiian pahoehoe flow lobes by Hon et al. Ž1994. and Keszthelyi and Denlinger Ž1996. and highlight the critical role played by surface morphology, lava thermal properties, and crystallization dynamics. Furthermore, the reported description of the various heat fluxes between lava and atmosphere can be extended to any other type of lava flows in which atmospheric cooling is involved. q 1998 Elsevier Science B.V. All rights reserved. Keywords: volcano; crust; lava

)

Corresponding author. E-mail: [email protected]. Also at: Department of Chemical and Environmental Engineering, Illinois Institute of Technology, Chicago, IL, USA. 0377-0273r98r$19.00 q 1998 Elsevier Science B.V. All rights reserved. PII S 0 3 7 7 - 0 2 7 3 Ž 9 8 . 0 0 0 1 0 - 9

216

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

1. Introduction Lava flows involve, during their emplacement, very complex mass, momentum, and heat transfer phenomena. Relevant examples are the fractional crystallization of the chemical components forming the magma, bubble nucleation and growth, substratum erosion, formation of superficial crust and levees at the flow boundaries, generation of superimposing flows, secondary boccas, and lava tubes. The high-temperature difference between lava and the surrounding environment, in addition to the non-homogeneous nature of magma, ensures that heat transfer plays a critical role in the evolution of lava flows. On the one hand, the high-temperature lava requires the consideration of different sets of heat transfer mechanisms in different regions of the flow: examples are the simultaneous radiative and convective heat transfers between lava and atmosphere, and the conductive plus convective heat transfer contributions inside the flow. On the other hand, the multiphase and multicomponent nature of magma accounts for large variations in its chemical, physical, and rheological properties with temperature and is therefore directly responsible for the wide range of behavior-types exhibited by lava flows. As a consequence of all these considerations, most of the thermal processes we can observe on lava fields are strongly coupled with mass and momentum transport processes as well as with specific constitutive equations of the magma. Therefore, a complete description of lava flow emplacement can be achieved only by a proper integration of all these processes. Despite such a complex picture, lava flow fields exhibit quite a large number of special cases for which simple assumptions can be useful to understand the physics of the process. One of the most common cases is represented by the cooling in the atmosphere of thermal diffusion-dominated lava flows. These flows are characterized by a negligible contribution of the flow convection to the cooling process and are indeed largely controlled by heat conduction along the normal direction to the surface. Even though such a condition is almost never satisfied in a strict sense, in several cases it can be applied, for a specific period of time, to better explain the observed processes. Typical examples are several types of pahoehoe flows, such as the slabby, toothpaste, spongy, or the smooth and ropy surface pahoehoe observed, for instance, in several Hawaiian and Etnean lava fields ŽRowland and Walker, 1987; Walker, 1989; Wilmoth and Walker, 1993; Kilburn, 1989, 1993; Hon et al., 1994.. At Hawaii, pahoehoe lava flows are extremely common and cover a large part of their surface. Pahoehoe flow-units are commonly 0.1 to 5 m thick, with an average thickness on the order of a meter, and have a width ranging from a few tens of centimeters to a few meters ŽWalker, 1989; Wilmoth and Walker, 1993; Hon et al., 1994.. In most cases, pahoehoe flow sheets are formed by the merging of smaller flow tongues or lobes moving at millimeters to centimeters per second ŽKilburn, 1993; Keszthelyi and Denlinger, 1996.. At Mount Etna, even though most lavas are classified as aa individual flow units, they commonly contain secondary pahoehoe flows. They usually occur at the head of a flow and their length may extend downstream for up to tens of meters ŽChester et al., 1985; Kilburn, 1989.. Also at Etna, several types of pahoehoe surfaces have been recognized and their nature has been related to the deforming stress of the flow motion on the cooling crust ŽKilburn, 1989, 1993; Calvari et al., 1994.. Both field and experimental evidence highlights that it is mainly low discharge rates that characterize such behavior ŽPinkerton and Sparks, 1976; Kilburn, 1989; Fink and Griffith, 1990; Griffith and Fink, 1992a.. An additional application of the model is represented by those flows that exhibit an upper plug-flow region across which no velocity gradient is observed ŽCrisp and Baloga, 1990; Dragoni and Tallarico, 1994; Polacci and Papale, 1997.. Under all of these conditions, superficial effects play a major role and lava flows appear to be crust-dominated. Understanding the formation of crust and of the underlying mushy layer is therefore critical for diffusion-dominated flows in which the heat flux in the normal direction to the surface is the controlling effect. In such a context, the present work is aimed at the investigation, from a heat transfer point of view, of the local cooling of a steady portion of lava suddenly exposed to the atmosphere and at the quantification of the heat transfer process during the first few hours of cooling. Most lava flow models proposed until now have focused on the influence of different physical formulations and on assumptions regarding the thermofluid-dynamics of a generic lava flow ŽBaloga, 1987; Crisp and Baloga, 1990; Dragoni, 1989.. The present analysis

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

217

especially addresses the study of diffusion-dominated lava flows and tries to highlight important aspects of the heat transfer process which, to date, have received less attention. Among the most relevant influences are the effects of the superficial morphology of the crust on the heat flux from the lava surface to the atmosphere, a more complete evaluation of the radiative and convective heat transfer coefficients between crust and atmosphere, and the description of the formation in time of the crust and thermal mushy layer resulting from the evaluation of the phase change effects over a range of temperature. Further interesting information is derived from the study of the influence of other physical properties of the magma, such as those related to the crystallization process and bubble content of the lava. The model was applied to the study of the initial cooling of lavas with thermal and physical properties representative of Etnean lavas. From several points of view, the analysis completes the recently published results by Hon et al. Ž1994. and Keszthelyi and Denlinger Ž1996., regarding the initial cooling of pahoehoe lobes at Hawaii, and preliminary results by Neri Ž1994, 1995. aimed at a better definition of heat transfer coefficient for lava flows at Etna ŽDobran, 1995.. It is worthwhile mentioning again that model results are strictly applicable to diffusion-dominated flows in which convective terms are negligible compared to conductive terms. Despite such a restriction, the present analysis has a more general significance and can be used to obtain useful information for the modeling of a large number of lava flows in which atmospheric cooling plays a major role ŽStasiuk et al., 1993; Dobran, 1995; Macedonio, 1995.. Heat transfer information, appropriately combined with a description of lava rheology and flow field, should lead to the definition of more complete lava flow models and to the identification of the main emplacement mechanisms controlling lava flows. Furthermore, parametric studies performed by models of this kind can quickly supply the system response to a variable perturbation and help in carrying out field and laboratory measurements as well as in the calibration of heat transfer equipment. 2. Description of the heat transfer model The initial emplacement mechanism of most pahoehoe flows discussed above mainly consists of a series of tongues or lobes that abruptly break the surrounding crust and expose fresh lava to the atmosphere ŽKilburn, 1993; Keszthelyi and Denlinger, 1996.. Such a mechanism, as well as direct temperature measurements performed by thermocouples and infrared pyrometers at Etna ŽCalvari et al., 1994. and Hawaii ŽHon et al., 1994., gives an approximation of the initial state of the lava with a homogeneous still portion of magma at a temperature above the solidus which immediately starts cooling down from its upper surface. Fig. 1 schematically illustrates the investigated system with its main features during the two cooling periods considered in the analysis. It consists of a semi-infinite region of still lava, initially at constant temperature T l,` . The magma was modeled as a vesiculated binary melt with given liquidus and solidus temperatures, Tliq and Tsol , respectively, and with crystal fraction dependence on temperature, fsŽT ., in order to study the effect of phase change over a range of temperature and, therefore, of crystallization dynamics. Crystal and bubble contents were considered in the estimation of the mean thermal properties of the magma, as will be discussed below. At time t s 0 the lava surface starts to exchange heat by radiation and convection with the surrounding atmosphere, supposed at constant temperature Tatm . Cooling thus starts and proceeds by conduction toward the lava interior while heat is released inside the lava due to crystallization. The effect of the flow on the solidification process is neglected and therefore the heat of crystallization is the only internal heat source. The cooling is associated with the growth of a superficial mushy region, with thickness ´ m1 , and with the downward propagation of an isothermal plane delimiting the portion of magma not affected by cooling. This type of cooling lasts until time t s t 1 , when the surface temperature Tc reaches the solidus temperature, and defines Period 1 of cooling. At this time Period 2 begins; it is characterized by the presence and growth of a solid layer, i.e., the crust, of thickness ´s2 , and an underlying mushy region previously formed. During the second period, the mushy zone, of thickness ´ m2 , is limited by two isothermal fronts at the initial and solidus temperature, respectively, and grows moving downward. As time passes, the solid crust grows as well, while the heat flux at

218

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

Fig. 1. Schematic representation of the solidification process during the two investigated periods. During Period 1, just a thermal mushy layer is present at the lava surface while the superficial temperature is still higher than the solidus temperature. During Period 2, superficial temperature drops below the solidus temperature and a solid crust begins to form on the surface. The box shows the two different superficial morphologies investigated in the study.

the upper surface decreases due to the decreasing surface temperature. In the box, the figure also shows the two simplified superficial morphologies considered for the crust in order to estimate such an effect on the cooling process. Other important effects due to superficial cracks ŽCrisp and Baloga, 1990. and radiative heat transfer through the vesicles ŽKeszthelyi, 1994. have been ignored in the present analysis due to their minor importance under the specific hypotheses made above and the relatively low vesicularity of Etnean lavas. In addition, such an approximation appears to be reasonable also in the light of the high uncertainty associated with lava conductivity as discussed in the following. After setting the initial and boundary conditions for the lava portion, the model is able to determine the timewise behavior of crust and mushy layer thickness, the temperature profile within them, the radiative and convective heat fluxes to the atmosphere, and the associated heat transfer coefficients as a function of physical properties and parameters of the system. Furthermore, with minor changes, the model can implement different physical properties for the crust and mushy layer in order to account for their specific features. In particular, the addition of a volume conservation equation will allow the description of the growth or decrease of lava volume during the cooling process. However, in order to simplify the results, in this study, and also because of the lack of information, the same properties have been assumed for the two layers. Before presenting model equations, boundary conditions at the crust surface and constitutive equations for the lava will be discussed.

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

219

2.1. Heat transfer mechanisms at the cooling surface Lava cooling in the atmosphere is an excellent example of the simultaneous action of radiation and natural, or forced, convection. Even though temperature level and lava emissivity are such that radiation is certainly the main mechanism during the very first minutes of cooling, natural and forced Žwind. convection by the surrounding air become very important after just a few minutes. A detailed description of the coupling between these heat transfer mechanisms and their transient and non-linear nature is beyond the objectives of this work and would require the carrying out of specific laboratory and field experiments as well as numerical simulations by physical models. Therefore, the focus is on the main parameters controlling the cooling process, indicating major effects, assumption limits, and critical points that need further explanation. Under the basic hypothesis that the physical properties of the atmosphere surrounding the lava surface are not affected by the cooling process, one may treat radiation and convective heat transfers as independent, simultaneous mechanisms. This enables one to define the total heat flux from the surface as: qtot s qrad q qconv and in terms of the radiant and convective heat transfer coefficients, h rad and h conv : qtot s h rad Ž Tc y Tmr . q h conv Ž Tc y Tatm . where the subscript conv pertains to the natural, or forced, convective term and Tmr represents the mean radiant temperature of the surroundings as will be discussed below. From this expression, we can obtain the global heat transfer coefficient h tot : h tot s qtotr Ž Tc y Tme . with Tme representing the mean equivalent temperature of the surroundings defined by: Tme s Ž h rad Tmr q h conv Tatm . r Ž h rad q h conv . We can now focus our attention on the three heat transfer mechanisms working at the lava surface and on their quantitative estimation. 2.1.1. Radiation Under the assumption that the radiative flux density leaving the surroundings, also called radiosity, is not influenced by the flux leaving the lava surface, the net radiant heat flux from the surface can be expressed by ŽHottel and Sarofim, 1967; Edwards, 1983.: y y qrad s qqy qys eeff s Tc4 y a surr qsurr y a sun qsun

Ž 1.

where qq and qy represent the radiosity and irradiation respectively, eeff the effective emissivity of the lava surface, s the Stefan–Boltzmann constant, and a surr and a sun the absorption coefficients of the two supposed irradiations, i.e., surrounding and solar irradiations. The surroundings irradiation, assuming the lava surface can exchange radiation just with the atmosphere, in its turn, can be expressed as: y q 4 s qatm s eatm s Tatm qsurr

where the atmospheric emissivity is determined by the Brunt equation ŽEdwards, 1983.:

eatm s 0.55 q 1.8

(Ž P

H 2 O rPatm

.

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

220

in which Patm is the total atmospheric pressure equal to 1 atm and P H 2 O is the partial pressure of H 2 O. From this equation we can calculate, for instance, an emissivity of 0.84 for air at 60% relative humidity. The mean radiant temperature Tmr of the surroundings can then be determined by comparison of its definition, i.e.: 4 qrad s eeff s Ž Tc4 y Tmr .

to Eq. Ž1.. It results in: 4

Tmr s

(

y a surr qy surr q a sun qsun

seeff

4

s

(

y qy surr q qsun

s

where the validity of Kirchhoff’s law, eeff s a surr s a sun , was assumed for the latter equality. From this equation, assuming an air temperature of 300 K, air emissivity of 0.84, and solar irradiation of 500 Wrm2 Žas an average value for diurnal irradiation., we obtain a mean radiant temperature of 353.5 K. In spite of the relevant difference between the atmospheric temperature Tatm and the mean radiant temperature Tmr , no major effects were observed in the heat transfer process during the first 2 h of cooling in comparison with results obtained considering a Tmr equal to the atmospheric temperature Tatm . Some differences were observed only when the superficial temperature of lava reached a temperature below 500 K. In this case, the effects of natural and forced convections mainly control the cooling process, as will be discussed below. As a consequence, in the present analysis, the atmosphere was considered as a black body, i.e., eatm s 1, and the influence of the solar irradiation was neglected. In this analysis, attention was given to the estimation of the effective emissivity, eeff , of the lava surface. Radiative properties of a surface indeed depend very much on the ratio of the root-mean-square roughness to the radiation wavelength ŽHottel and Sarofim, 1967.. For low values of this ratio, the slope of the surface irregularities is small, multiple reflections are negligible, and the surface behaves as a partial specular reflector. For high values of the ratio, the surface is covered with deep cavities and radiations perform several reflections before leaving the surface. In this case, the radiative effect is estimated by reducing the rough surface to a smooth plane and adopting an effective emissivity eeff as a function of surface morphology and real emissivity of the material. It is evident that, for lava surfaces, the very short wavelength of infrared radiation Žon the order of microns., and the quite large value of surface irregularities Žranging from a few millimeters to several centimeters., require replacement of the real emissivity of the lava with an effective emissivity of the rough surface. At present, there are detailed studies on the morphology of several types of lava flows surfaces by various authors ŽKilburn, 1989, 1993, 1996; Walker, 1989; Kilburn and Lopes, 1991; Polacci and Papale, 1997.. A large variety of superficial morphologies has been reported and associated to different local features of the flow and cooling process. In several cases, their superficial patterns were described as a continuous rough sequence of cavities of various Žregular and irregular. shapes whose dimensions range from a few millimeters to several centimeters ŽKilburn, 1989, 1996.. However, at present, neither standard superficial morphologies have been defined, nor is it possible to exactly reproduce the complexity of a real superficial morphology. As a consequence, for the purpose of giving an idea of such an effect on the radiative properties of the surface, some results from studies on cavities with simple geometries ŽSparrow and Jonsson, 1962; Sparrow and Lin, 1963. have been adopted. Such a choice appeared to be the most reasonable in order to quantify first-order effects due to morphology and also in light of the reported description of real lava surfaces. The box reported in Fig. 1 shows the two cavity shapes investigated. Cavity length in the direction normal to the sheet was considered sufficiently great so that end effects were neglected. The first one represents a V-groove cavity with opening angle q , whereas the second is a rectangular-groove cavity with aspect ratio lra. The quantification of the cavity effect depends on the shape of the cavity as well as on the reflecting properties of the surface and the type of the surrounding medium. The nature of the magma and the above mentioned

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

221

approximation of the atmosphere to a plane black surface suggest considering the lava surface as a diffusive reflector Žin other words, the surface reflects according to Lambert’s cosine law. and the energy outcome from the cavity as diffuse Ži.e., uniformly distributed over the opening of the cavity.. The effect of the atmospheric gases inside the cavity due to the small dimensions of the system was neglected. Fig. 2 shows the values of the effective emissivity as a function of the material emissivity and geometry of the two investigated cavities. In particular, Fig. 2a illustrates the relation between effective emissivity and the opening angle of a V-groove cavity, and Fig. 2b illustrates the same relation as a function of the aspect ratio lra of a rectangular-groove cavity. It should be noted that, for average emissivity values typical of lava ranging from 0.6 to 0.95 ŽDragoni, 1989; Crisp et al., 1990; Kahle et al., 1988., the presence of superficial cavities with an opening angle of 308 Žor aspect ratio of 1. causes an increase of the effective emissivity of 30% and 5%, respectively. Therefore, surface roughness plays a greater role the lower its emissivity. In Section 3, the effects of morphology on the radiative and convective fluxes, and therefore on the cooling process, will be shown.

Fig. 2. Ža. Effect of a V-groove cavity on the effective emissivity of a surface. The figure gives the values of the effective emissivity eeff as a function of the opening angle of cavity, q , and of the emissivity of the plane surface wmodified after the work of Sparrow and Lin Ž1963.x. Žb. Effect of rectangular-groove cavity on the effective emissivity of a surface. In this case, the effect of the emissivity and aspect ratio, l r a, of the cavity is illustrated wmodified after the work of Sparrow and Jonsson Ž1962.x.

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

222

2.1.2. Natural conÕection An accurate prediction of heat transfer by natural convection from a complex high-temperature surface requires both the solution of a three-dimensional heat transfer model and the carrying out of ad hoc laboratory experiments. Up to date, some analyses of this kind have been performed, and a few laboratory experiments using air are available ŽRaithby and Hollands, 1984.. Therefore, the study investigated just the effects of simple extended surfaces on the convective heat transfer coefficient, in order to give a first order estimate of such an effect. In particular, the above presented geometries, i.e., the V- and rectangular-groove cavities, were studied. In the following, attention is restricted to isothermal fins, assuming a negligible effect of temperature drop between the top and the base of the cavity. Properties of the air in contact with the hot surface are evaluated at the mean temperature ŽTc y Tatm .r2. For an extended horizontal surface constituted by a sequence of V-groove cavities, exchanging heat with air, the convective heat transfer coefficient, h nat , has been obtained by Al-Arabi and El-Rafaee Ž1978. in the form of semi-empirical correlations of the kind: Nu s C1 RaC 2

Ž 2.

where Nu and Ra represent the Nusselt and Rayleigh numbers defined by: Nu s

h nat H k atm

Ra s GrPr s

2 H 3ratm g b Ž Tc y Tatm . catm matm 2 matm

k atm

s

H 3 g b Ž Tc y Tatm .

natm a atm

with a clear meaning of symbols Žrefer to nomenclature for explanation.. The correlations provided by the mentioned study are valid for 1.8 = 10 4 - Ra - 1.4 = 10 7 and depend on the flow regime Žstreamline or turbulent according to the Ra. and on the geometry of the cavity. These are described by the following expressions of the coefficients; for 1.8 = 10 4 - Ra - Rac : C1 s

0.46 sin Ž qr2 .

y 0.32 ;

C2 s 0.148sin Ž qr2 . q 0.187

Ž 3.

C2 s 1r3

Ž 4.

and for Rac - Ra - 1.4 = 10 7: C1 s

0.054 sin Ž qr2 .

q 0.09;

The critical Rayleigh number, Rac , that divides the streamline regime from the turbulent regime is given by: Rac s Ž 15.8 y 14.0sin Ž qr2 . . = 10 5

Ž 5.

In our case, considering a surface roughness on the order of centimeters, the characteristic Ra numbers range from 10 4 to 10 5 and, therefore, streamline conditions of the flow should be involved. 2 Under these conditions, for a smooth surface, the group Ž k atm rH .C1wŽ H 3ratm g b catm .rŽ matm k atm .xC 2 remains almost constant over a wide range of temperature involved in the cooling process Žwithin about 15% for 4008C - Tc - 11008C.. The convective heat coefficient can then be expressed in the form: h nat s lŽ Tc y Tatm .

d

Ž 6.

and using the appropriate values of the air physical properties by: h nat f 1.0 Ž Tc y Tatm .

0.33

Ž Wr Ž m2 K. . .

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

223

Such an expression gives values of about 8–10 WrŽm2 K. during the whole cooling process, which are fully consistent with previously reported estimates ŽJakob, 1957; Turcotte and Schubert, 1982; Knudsen et al., 1984; Head and Wilson, 1986.. For a corrugated surface, things change substantially. In this case, for a streamline regime of the flow, the natural convection heat transfer coefficient depends on the characteristic length. For instance, for a V-groove cavity with opening angles of 908 and 308 and characteristic length H of the order of few centimeters, it results, respectively: h nat f 3.5 Ž Tc y Tatm .

0.29

Ž Wr Ž m2 K. .

and: h nat f 10.0 Ž Tc y Tatm .

0.22

Ž Wr Ž m2 K . .

Ž 7.

It should be noted that for turbulent flow conditions Žexpected with a surface roughness of the order of decimeters. the convective heat transfer coefficient remains independent from the characteristic length H. Similar expressions have also been determined for horizontal surfaces with rectangular-groove cavities ŽRaithby and Hollands, 1984.. In this case, the reported expressions appear more complex even though the effects are comparable. We refer to the cited work for a description of this case. In any case, it should be noted that a corrugated surface has a strong effect on the natural convection heat transfer coefficient, with values of up to three to four times greater than the values typical of a flat surface. 2.1.3. Forced conÕection As previously pointed out by Head and Wilson Ž1986. and Keszthelyi and Denlinger Ž1996., forced convection by the wind is able to significantly affect the cooling of diffusion-dominated lava flows. Such an effect is even more important for submarine lavas where the high-conductive water flow represents an effective heat sink ŽGriffith and Fink, 1992b.. As in the case of natural convection, the heat transfer rates for forced convection are generally established through a combination of theoretical and experimental analysis. For practical use, theoretical analysis is predominantly based on boundary layer theory ŽSchlichting, 1968; Coulson and Richardon, 1977; Rubesin et al., 1984.. In the following analysis, an estimate of the forced convection heat transfer coefficient has been obtained in terms of the average Stanton number, St, which for a flat surface can be expressed by ŽRubesin et al., 1984.:

H Sts

A

qforc d A

ratm u atm c p ,atm Ž Tc y Tatm . A

Ž 8.

where the subscript atm refers now to the undisturbed air flow and u is its parallel-to-surface velocity. The St number, expressing the ratio between the heat lost by forced convection and the maximum amount of heat removable from the surface by the air stream, can be used to yield the average Nusselt number, Nu, by: Nus

h forc H k

sStRe H Pr .

Ž 9.

For a flat surface at a constant temperature, the average St number can be related, for both laminar and turbulent boundary layers, with the average skin friction coefficient c f through the useful correlation ŽColburn, 1933; Schlichting, 1968.: Stf

cf 2

Pry2r3

Ž 10 .

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

224

where c f accounts for the surface roughness of the plane. Schlichting Ž1968. showed that such an effect gives an important contribution to the total heat exchange and that the ratio k srX, where k s is the typical roughness dimension and X is the characteristic length of the plane, is particularly critical. Therefore, for forced convection, we will refer to a generic roughness of the surface Žusually produced by sand protuberances. abstracting from the V- and rectangular-groove cavity shapes so far considered. For Reynolds numbers greater than 10 5, as can be expected for a wind speed above 5 mrs, characteristic lengths of about 1 m, and atmospheric air at standard conditions, the skin friction coefficient reaches a constant value and the flow regime is described as ‘fully rough’. Under these conditions Schlichting gives the following correlation for the average skin friction coefficient:

ž

c f s 1.89 q 1.62log

X ks

y2 .5

/

Ž 11 .

which assuming a length scale on the order of 1 m and a roughness of the order of 1 cm, i.e., Hrk s on the order of 0.01, gives an average skin fraction of 1.4 = 10y2 . As an example, the heat transfer coefficient associated with a wind of 5 mrs and with the above computed average skin friction coefficient is about 60 WrŽm2 K. vs. 20 WrŽm2 K. associated with a smooth surface. Such a value compares reasonably well with the experimental value of about 70 " 30 WrŽm2 K. associated by Keszthelyi and Denlinger Ž1996. with a wind of 3–4 mrs, even though no information on surface roughness is reported by the authors. Furthermore, Eq. Ž9. states that the heat transfer coefficient is proportional to the wind speed and is independent from the surface temperature. Therefore, in the developed model, the heat flux from the lava surface was described by a boundary condition of the kind: 4 qtot s seeff Ž Tc4 y Tatm . q lŽ Tc y Tatm .

d

Ž 12 .

where the coefficients l and d vary according to the considered convective heat transfer mechanisms as reported above. 2.2. ConstitutiÕe equations In order to describe the main thermal processes affecting the cooling of lava, appropriate constitutive equations have to be defined. In this study, the magma was modeled as a binary melt characterized by mean thermal properties, taking into account the presence of bubbles, crystals, and phase change over a range of temperature. Such representation is certainly a simplification of the actual nature of the magma, but we will see how it is still able to represent the main first-order effects of the cooling process. Particularly critical is the representation of the crystallization dynamics. Brandeis et al. Ž1984., Cashman Ž1993., and Keszthelyi and Denlinger Ž1996. report a clear description of the several factors affecting this process. Main effects are due to cooling rate and magma composition that strongly control non-equilibrium phenomena such as undercooling or critical cooling rate for the onset of crystal nucleation and growth. Unfortunately, most of these effects are scarcely known and constrained. In spite of this, for our purposes, the modeling of magma as a binary melt supplies useful insights into the evaluation of the main effects. Such description requires the definition of liquidus and solidus temperatures and the specification of a relationship between crystal content of the mixture and its temperature. In this way, the crystallization process occurs over a range of temperatures and therefore, it accounts for a distributed release of heat. The definition of a range of temperature over which the phase change occurs also allows the description of the thermal mushy layer below the crust. The study of this layer appears particularly important because its extension can significantly affect the lava rheological properties and therefore, its emplacement mechanism.

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

225

As the solid phase is continuously formed in the mushy region as cooling propagates, the crystallization heat was considered as a volumetric generation term of the form: G Ž x ,t . s r L

d fs

Ž 13 .

dt

where r is the density, L the crystallization heat per unit mass, and fs the crystal fraction in the mushy layer. Therefore, in addition to the latent heat value, the relationship between crystal fraction and temperature is particularly important in order to take into account the real state diagram of the cooling magma. In this study, due to the lack of detailed expressions, the crystal fraction is assumed to vary linearly with temperature and is expressed by the simple relation:

ž

fs s feu 1 y

Tm y Tsol T liq y Tsol

/

Ž 14 .

where Tm ' Tm Ž x,t . is the temperature in the mushy zone and feu is the solid fraction at the eutectic temperature in the case in which the binary melt presents an eutectic state. It is evident that the crystal fraction has a value of zero at the liquidus front and of feu at the solidus front. For a binary melt outside the eutectic range feu is equal to unity. Eq. Ž14. can be easily substituted with other relationships of the crystal fraction as a function of temperature in order to account for different equilibrium or non-equilibrium conditions during the cooling process. For instance, our expression can be compared with Eq. Ž10. by Keszthelyi and Denlinger Ž1996. where the crystal variation in time is expressed by using a chain rule with respect to temperature. In our case, the derivative of crystal content with respect to temperature is expressed in term of solidus and liquidus temperatures, whereby solidus temperature was supposed to be representative of the glass transition temperature. The effect of the critical cooling rate on the onset of crystallization has been neglected in our description due to its weakly constrained value and its limited influence on the cooling process, as already shown by the above authors. Finally, the possible presence of an eutectic state allows us to investigate the net effect of a distributed or localized release of crystallization heat. The thermal properties of the crust and mushy layers have been estimated taking into account the important effect of bubbles ŽZimbelman, 1986; Aubele et al., 1988; Keszthelyi, 1994., whereas mean properties have been assumed for the melt and crystals. Particularly relevant is the effect of bubbles in the estimation of the mean conductivity and density of the two layers. In regards to the thermal conductivity, Horai Ž1991. reported a wide number of formulae that have been put forward in its estimation. Here, a well-known and tested expression based on the work of Hashin and Shtrikman Ž1962. was adopted. It consists in an arithmetic mean of the maximum and minimum conductivities estimated for a homogeneous and isotropic material. The maximum and minimum conductivities are:

w k U s k cont q

k L s kg q

Ž k g y k cont .

y1

q Ž 1 y w . Ž 3k cont .

Ž1yw .

Ž k cont y k g .

y1

q Ž w . Ž 3k g .

y1

y1

Ž 15 .

Ž 16 .

where w indicates the porosity of the layer and the subscripts g and cont are representative of the gas Žassumed water vapor. and continuum phases Žmelt plus crystals., respectively. The radiative contribution to the effective thermal conductivity of the crust and mushy layer that occurs through the vesicles has been neglected due to its limited effect in the adopted low-porosity conditions ŽKeszthelyi, 1994..

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

226

Table 1 Physical properties and parameters employed in the application to Etnean lavas ŽKilburn, 1993; Murase and McBirney, 1973. Property

Value

Liquidus temperature, Tliq ŽK. Solidus temperature, Tsol ŽK. Lava initial temperature, T l,` ŽK. Atmosphere temperature, Tatm ŽK. Melt-Crystal density, rcont Žkgrm3 . Melt-Crystal conductivity, k cont ŽWrŽmK.. Melt-Crystal specific heat, c cont ŽJrŽKg K.. Crystallization heat, L ŽJrKg. Vesicle volumetric fraction Žporosity., w f Eutectic fraction, feu Emissivity, e Dimensional constant, l ŽWrŽm 2 K dq1 .. Dimensionless constant, d

1473 1223 1373 300 2650 2.0 1200 0%3.5=10 5 0%0.3 0.0%1.0 0.6%0.9 1%10 1.22%1.33

The mean density has been equivalently calculated, on the basis of lava vesicularity and mean values of crystals and melt densities, by the relation:

r s wrg q Ž 1 y w . rcont f Ž 1 y w . rcont Finally, the effects of bubbles on other thermal properties, such as specific and crystallization heat, have been neglected due to their minor importance. It is important to outline that almost all the described properties show variations with temperature ŽMurase and McBirney, 1973; Touloukian et al., 1989.. As far as lava density and specific heat are concerned, variations are limited to few percentages of the average value that can therefore represent a good approximation for the entire cooling process. In regards to latent heat and thermal conductivity, values are more sensitive to magma composition and have a lower accuracy. Sometimes they also show contradictory behavior: for instance, increasing and decreasing conductivities with temperature are shown by Birch and Clark Ž1940., Murase and McBirney Ž1973., Peck et al. Ž1977. and Touloukian et al. Ž1989. for similar basaltic materials. As a consequence, for our purpose, mean constant values of these two properties were assumed, although the model allows implementation of temperature-dependent relationships. Physical properties employed in the model are reported in Table 1 and can be considered representative of Etna’s lavas. These values are not referred to a specific magma composition but represent mean values reported in the literature for Etnean lavas or similar magmas ŽKilburn, 1993; Williams and McBirney, 1979; Hon et al., 1994.. Such a choice, even though forced by the lack of data, was considered adequate for the described objectives of the work. Anyway, for some properties and parameters, the effect of their variations within appropriate ranges were studied due to the great uncertainty associated with them. 2.3. Model equations and solution procedure The two cooling periods described above were studied separately because of the different nature and number of the required heat transfer equations. We will now illustrate the implemented equations and the solutions obtained for the two periods. The non-linearities present in the physical process, such as for instance the strongly non-linear boundary condition Ž12. or the possible variations of physical properties with temperature, prohibit exact analytical solutions, and allow only approximated analytical and numerical solutions. In this work, the solution was obtained by an original application of the integral heat balance method ŽGoodman, ¨ 1964; Ozisik, 1993. and numerically solving the reduced system of ordinary differential equations, in terms of

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

227

surface temperature and crust and mushy layer thicknesses, by a Runge–Kutta method. Such a method represents a straightforward and widely applicable technique to study transient, non-linear, heat transfer problems involving a phase change. Solutions obtained by this method, when tested with problems showing an exact solution, demonstrate a good accuracy Žwithin few percentages. for many practical and engineering applications. In the following, the solution procedure will be outlined, leaving in Appendix A the description of method details. Period 1. As discussed above, this first period is characterized by the formation of a mushy layer at the top of the investigated semi-infinite region and by a surface temperature comprised between the initial temperature Tl,` and the solidus temperature of the magma. In other terms, this period describes the very first cooling moments when no crust has formed yet and a mushy layer starts growing. Assuming constant thermal properties, Fourier’s conduction equation relative to the mushy layer can be expressed by:

E 2 Tm 1 Ex

2

1 E Tm 1

G q

s

for 0 F x F ´s

Et

am

km

Ž 17 .

where the temperature has been referred to the solidus temperature, k m and a m represent the thermal conductivity and diffusivity of the mushy zone, and the subscript m1 refers to the mushy region during Period 1. The corresponding boundary conditions are: x s 0, 0 - t F t 1 ;

km

E Tm 1

E Tm 1

x s ´m , 0 - t F t1 ;

4 s qtot s seeff Ž Tc4 y Tatm . q lŽ Tc y Tatm .

Ex Ex

d

s0

Ž 17a. Ž 17b.

and the initial condition is:

´m Ž t s 0 . s 0 or equivalently: Tc Ž t s 0 . s Tl ,` s Tm1 Ž x ,0 . In the previous equations Tc Ž t . corresponds to the superficial temperature, i.e., Tm1Ž0,t ., while Eq. Ž17a. defines the total heat flux from the lava surface toward the atmosphere as defined above and: DTl ,` s T l ,` y Tsol

According to the integral heat balance method the integration of Eq. Ž17. over the thermal mushy layer from x s 0 to x s ´m results in:

E Tm 1 Ex

y

E Tm 1

´m

1 s

am

Ex

´m

H0

1 km

0

E Tm 1 Et

´m

q

H0

G Ž x ,t . d x

dx

Ž 18 .

where, using Eqs. Ž13. and Ž14., the crystallization heat is expressed by: G Ž x ,t . s y

r Lfeu E Tm DT

Et

.

Ž 19 .

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

228

To generalize the solution, the following dimensionless variables are defined: xs s

x D

´m s

;

´m D

´s s

;

´s D

; Ts

T y Tsol DT l ,`

;

ts

t as D

2

;

ks

km ks

;

cs

cm cs

;

as

am as

;

Ls

L) cs DTl ,`

where D is the characteristic length and L) s LŽ1 y feu .. Assuming a quadratic profile of the temperature, after some manipulations Žsee Appendix A., we obtain the final ordinary differential equation in Tc and the associated expressions for the thickness of the mushy layer, ´ m 1 , and the temperature profile inside it: dTc

2 qtot

sy

dt

´m s y

Tm 1 s

2 Z Ž Tc y 1 .

Y Ž Tc y 1 .

3

qtot

Ž 20 . y2

2 Ž Tc y 1 . qtot

Ž Tc y 1 . ´ m2

x2y

2 Ž Tc y 1 .

´m

x q Ž Tc y TsrDTl ,` .

Ž 21 .

where: qtot s

ž

qtot D l ,`Tk m

Ž 1ras .

Ž 22 .

/

D Z

Zs

,

; with Z s

1 q

am

r Lfeu

Ž 23 .

k m DT

and the variable Y can be expressed by: 4seeff Tc3 q dl Ž Tc y Tatm .

Y

Ys

s

Ž k m rD .

dy1

Ž k m rD .

where: d qtot dt

sY

dTc

Ž 24 .

dt

Period 2. This period is characterized by the presence of a solid layer, the crust, overlying the thermal mushy layer formed during the first period and therefore by values of superficial temperature below the solidus temperature. Fourier’s conduction equations for the solid and mushy layers are respectively expressed by:

E 2 Ts 2 Ex

2

E 2 Tm 2 Ex

2

s

1 E Ts 2

as E t G

q

s km

for 0 F x F ´s 1 E Tm 2

am

Et

for ´s F x F Ž ´s q ´ m .

Ž 25 . Ž 26 .

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

229

The boundary conditions for the solid region are:

E Ts 2

4 s qtot s seeff Ž Tc4 y Tatm . q lŽ Tc y Tatm . Ex E Ts 2 E Tm 2 d ´s ks s km q r L Ž 1 y feu . Ex Ex dt

x s 0, t ) t 1 ;

ks

x s ´s , t ) t 1 ;

d

Ž 25.a . Ž 25.b .

while those for the mushy region are: x s ´s , t ) t 1 ;

ks

E Ts 2

s km

Ex

E Tm 2 Ex

E Tm 2

x s Ž ´s q ´ m . , t ) t 1 ;

Ex

q r L Ž 1 y feu .

d ´s

Ž 26.a .

dt

s0

Ž 26.b .

with initial conditions: ´sŽ t s t 1 . s 0 and Tc Ž t s t 1 . s Tsol . The application of the integral heat balance method to the two layers, and the use of boundary conditions and Eq. Ž19. lead to the following equations: k m E Tm 2

y

ks

Ex

y

k s E Ts 2

´s

E Ts 2 Ex

q

km E x

q

r L Ž 1 y feu . d ´s ks

0

km

´s

as d t

dt

r L Ž 1 y feu . d ´s

1 dQs 2

y

r Lfeu

dQm 2

k m DT

dt

1 y

q

am

dt

s0

Ž 27 . y DTl ,`

ž

d ´s dt

q

d ´m dt

/

s0

Ž 28 .

where the auxiliary variables Qs2 and Qm2 are defined in Appendix A. By assuming again a quadratic temperature profile in each layer, we obtain a final system of ordinary differential equations in ´s and Tc which, in terms of dimensionless variables, can be expressed in the form:

° ž N qN P / ~ 1

2

1

d ´s dt

q N2 P2 y Q2

ž

d ´s

/

dTc dt

s Q1

Ž 29 .

dTc

¢ž M q M P / d t q Ž M P . d t s M 1

2

1

2

2

3

where the expressions of symbols are reported in Appendix A. This system has been solved by a Runge–Kutta method using the above reported initial conditions and a variable time step to guarantee a relative error of less than 10y6 . Similarly, the thickness of the mushy layer and the temperature distributions inside the two layers result in: k 2 qtot ´s q 2 kTc

´m s

A) s

y

2 qtot ´s k 2

2

y

2Tc2

´s

q

L)Tc

´s

B)

Ž 30 .

q qtot k L) y 2Tc

2

2

Ts 2 s A s 2 Ž ´s y x . q Bs 2 Ž ´s y x . q Cs 2 ; Tm 2 s A m 2 Ž x y ´s . q Bm 2 Ž x y ´s . q Cm 2

Ž 31 .

where the constants A, B and C of the two distributions are independent of x, and can be determined by the prescribed boundary conditions. In particular they result in: As 2 s

1

´ s2

Tc q Ž b 1 q b 2 . ;

Bs2 s y

1

´s

Ž b1 q b2 . ;

C s2 s 0

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

230

with:

b1 s

k ´s

´m

y L);

(

b 2 s b 12 y 2 L)Tc

Ž 32 .

for the crust and: Am 2 s y

1

´m2

;

Bm 2 s

2

´m

;

Cm 2 s 0.

Ž 33 .

for the mushy layer.

3. Results and discussion In this section, some results of the described model are presented mainly to establish the relative importance of different physical properties and parameters of the system and to highlight their effect on the cooling process. In particular, the effects of variations of the lava’s superficial morphology, vesicle content, and crystallization dynamics on the cooling system were investigated. Types of behavior are reported over the first 2 h of cooling in order to cover a wide range of cooling conditions and to investigate both the initial and long-term features of the heat transfer process. Fig. 3 illustrates, in dimensionless form, the temporal evolutions of crust and mushy layer thicknesses ŽFig. 3a. and temperature profiles inside lava ŽFig. 3b.. The figure was obtained assuming the physical properties reported in the caption and a characteristic length, D, equal to 0.1 m. A first interesting result shown in Fig. 3a is that the thermal mushy layer initially grows more rapidly than the crust, whereas at longer times, the crust becomes thicker. From the figure reported in the box, we can observe that crust formation on the lava surface is very rapid and therefore, the associated Period 1 is quite short. If we transform the variables into dimensional form, we can estimate the time of crust formation at about 10 s from the start of cooling. We can estimate the time at which the crust becomes thicker than the underlying mushy layer at 15 min. In the same way, we can estimate that the crust and the mushy layer thicknesses have reached, after 2 h of cooling, about 7.5 and 5.5 cm, respectively. Fig. 3b shows the very high temperature gradient Žup to 100 Krmm. existing on the lava surface during the first minutes of cooling. Another interesting and non-intuitive result consists of the progressively slower decrease of superficial temperature in time and how, after the first few minutes, the highest temperature variations in time are observed within the lava and not on its surface. The net effect of lava surface morphology on the cooling process is illustrated in Figs. 4 and 5. These figures compare results obtained for a perfectly smooth surface and a surface characterized by the presence of V-groove cavities on it. The characteristic length of cavities was considered to be on the order of some millimeters and the opening angle of the cavity was assumed to be 308. Hypothesizing a lava emissivity of 0.8 and the described geometry of the cavity, Fig. 2a gives an effective emissivity of the surface of about 0.9. In the same way, the effect of surface roughness on the natural convection heat transfer coefficient can be determined from Eqs. Ž3. – Ž5.. For the described geometries, the values of l and d reported in Eq. Ž7. were computed. The other physical properties are reported in the caption. Fig. 4a shows the growth of crust and mushy layer as a function of time for the smooth and rough surfaces. From the figure, it can be noted that surface conditions significantly change the relative proportion of crust and mushy layer but only slightly influence the distance reached by the cooling front inside the lava Ži.e., ´ m q ´s .. In more detail, the rough surface causes a more rapid cooling, an increase of the crust thickness Žalmost 1 cm in 2 h. and an equivalent decrease of the mushy layer thickness. More evident is the effect of surface roughness on the temperature profiles within the lava even after only a few minutes of cooling Žsee Fig. 4b.. In fact, the more effective radiative and convective heat transfers strongly influence the superficial temperature with variations of

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

231

Fig. 3. Ža. Timewise distribution of the thickness of the solid Žcrust. and thermal mushy layer during the first 2 h of cooling. The box represents an enlargement of the first minutes of cooling. The plot is represented in dimensionless form using a characteristic length Ds 0.1 m. Model parameters correspond to w s 0.1, eeff s 0.8, l s1.0, d s1.33, Ls 3.5=10 5 Jrkg, and feu s1.0. Žb. Temperature profile inside the lava at different dimensionless times correspond to 60, 300, 900, 3600, and 7200 s.

about 1008C in the first minutes of cooling. Such a variation increases at longer times up to about 1308C after 2 h Žsee also Fig. 5b.. However, as already shown by Fig. 4a, such a large variation of superficial temperature only slightly affects the transmission speed of the cooling front inside the lava. In fact, this effect was expected due to the dependence of the cooling depth on thermal diffusivity, which, in this case, is the same for both conditions. Therefore, a rough superficial morphology of the lava strongly affects the surface temperature and the most superficial layers of the lava, but no influence is expected on the thickness of the portion affected by cooling. An increase of the crust thickness, with respect to the thermal mushy layer, is also observed. Similar, and even greater effects, were found for lava with lower emissivity and a different geometry of the cavities. Fig. 5 addresses the effect of natural and forced convection on the total heat transfer from the lava surface to the atmosphere. Both superficial morphology and wind effects have been investigated and quantified from the study of convective fluxes and lava superficial temperature. In detail, comparisons were made between a smooth surface and two rough surfaces characterized by different opening angles of their V-groove cavities. A third comparison was made considering a rough surface affected by a 5 mrs wind. In the first two cases, radiation and natural convection are acting on the surface, whereas, in the third case, radiation and forced convection due to the wind define the total heat transfer. Fig. 5a illustrates the behavior of the ratio between the convective and

232

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

Fig. 4. Ža. Timewise distribution of the thickness of the crust and thermal mushy layer during the first 2 h of cooling for a smooth and rough surface. The smooth surface is characterized by e s eeff s 0.8, l s1.0, and d s1.33. The rough surface is composed by V-groove cavities with opening angle q s 308 and is characterized by eeff s 0.9, l s10, and d s1.22. Other model parameters are w s 0.1, Ls 3.5=10 5 Jrkg, and feu s1.0. Žb. Temperature profiles inside the lava at 60, 300, 900, 3600, and 7200 s for the smooth and rough surfaces.

the total flux from the surface as a function of time for the four cases. For the smooth surface, the convective flux increases rapidly in the first minutes of cooling, reaching a value of 10 of the total flux after about 3 min and over 20% after 2 h. A large increase of the natural convection flux is observed for rough surfaces. The supposed superficial V-groove cavities with q s 908 and q s 308 increase the above defined ratio up to 30% and 40% after 5 min and 40% and 60% after 2 h, respectively. The presence of a 5 mrs wind acting on the rough surface with a cavity opening angle of 308 further increases the convective contribution to the total flux with ratios ranging from 0.4 up to 0.8. Similar effects can be observed in Fig. 5b, in which the superficial temperature is plotted vs. time. The plot clearly shows, as already discussed above, the drop in superficial temperature due to the significant effect of superficial roughness on the convective heat transfer coefficient. A further cooling of the surface is observed when the 5 mrs wind acts on the surface. In this case, the temperature drop is about 50–808C with respect to the case of simple natural convection acting on the rough surface. The effect of forced convection on the growth of the crust and mushy layer is similar to that reported in Fig. 4 for the superficial morphology effect and is not shown here. The above reported estimates for natural convection acting on rough surfaces are significantly greater than those previously applied in the literature ŽHead and Wilson, 1986; Oppenheimer, 1991; Keszthelyi and Denlinger,

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

233

Fig. 5. Ža. Timewise distribution of the ratio between the convective and total flux for smooth and rough surfaces characterized by different opening angles of the V-groove cavities and atmospheric conditions. For the rough surface, affected by the 5 mrs wind, the opening angle of the cavities is 308. Superficial features and lava properties are reported in the text. Žb. Timewise distribution of the superficial temperature for the four investigated cases.

1996., whereas, for forced convection, the suggested correlation appears to confirm the importance of this heat transfer mechanism ŽKeszthelyi and Denlinger, 1996. and able to relate it to important variable as surface roughness and wind speed. It is also important to highlight that the presence of superficial roughness makes the estimation of natural and forced convection heat transfer coefficients extremely critical from the very first minutes of cooling. The estimation of these effects is even more relevant over a longer period of time when these mechanisms control the cooling process of a wide range of lava flows. As already pointed out by several authors ŽZimbelman, 1986; Jones, 1993; Keszthelyi, 1994., vesicularity is certainly among the most important physical properties controlling lava cooling, due to the strong effect of bubbles on both thermal conductivity and density. Fig. 6 illustrates the effect of different lava vesicularities, chosen in a range commonly observed in Etnean lavas ŽKilburn, 1993; Polacci and Papale, 1997., on the growth of the crust and mushy layer, and on the temperature distribution inside them. As we can see from Fig. 6a, a variation of vesicularity from 0 to 0.3 brings about a considerable decrease in the growing rate of both the crust and the mushy layer. Lower superficial temperatures and larger temperature gradients are also associated with higher vesicularities. These effects are qualitatively consistent with experimental data on Hawaiian pahoehoe flows by Jones Ž1993. and Keszthelyi and Denlinger Ž1996., who reported a clear sensitivity of surface temperature on lava porosity.

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

234

Fig. 6. Ža. Timewise distribution of crust and mushy layer thickness during the first 2 h of cooling for vesicle volumetric fractions equal to 0 and 0.3. The surface is characterized by eeff s 0.9, l s10, and d s1.22. Other model parameters are Ls 3.5=10 5 Jrkg, and feu s1.0. Žb. Temperature profiles inside the lava at 60, 300, 900, 3600, and 7200 s for the two values of vesicularity.

This effect can be clearly observed in Fig. 6b, which shows the temperature profiles at different times. In this case, such effects are the result of two combined factors that occur at high vesicularity: Ža. a significant decrease in surface temperature, and Žb. a slower propagation of the cooling front inside the lava. These two factors can be easily related to the thermal inertia and thermal diffusivity of the lava, respectively expressed by:

'

I s kc r ; a s

k cr

.

Ž 34 .

Thermal inertia is a measure of the surface temperature variations induced by changes in the thermal boundary condition of the system ŽZimbelman, 1986; Keszthelyi, 1994.. In our case, the supposed change in vesicularity causes a change in the thermal inertia of the system from about 2500 JrŽm2 K s 0.5 . to less than 1400 JrŽm2 K s 0.5 ., with a drop in the superficial temperature of about 1508C only after a few minutes of cooling. On the other hand, the lower penetration of the cooling front is mainly related to the thermal diffusivity of the system. In our case, the stronger decrease in thermal conductivity with respect to density causes a decrease of thermal diffusivity Žfrom 6.3 = 10y7 to 3.9 = 10y7 m2rs. for more vesiculated lavas and, therefore, a slowing down of their cooling rate. Anyway, whereas thermal inertia always decreases by increasing vesicularity, the diffusivity trend is strongly related to the relative effect of vesicularity on thermal conductivity and density. For instance, Wilmoth and Walker Ž1993. observed a reversed behavior of crust thickness vs. porosity probably due to the

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

235

much higher range of porosity investigated and to the non-negligible effect of radiation inside the lava in that case. In the next paragraph, we will see how thermal diffusivity has to be corrected in order to account for crystallization dynamics and to better describe the penetration depth of the cooling front. Such simple considerations highlight again how important the knowledge of the thermal properties of lava is and how different and non-intuitive behavior can be simply explained in terms of a consistent framework. Another important factor that plays a major role in the cooling of lava is the crystallization process. In the next two figures we will examine both the effects of the amount of crystallization heat and its release route associated, for instance, with the presence of a eutectic. As far as crystallization heat is concerned, as already stated above, its knowledge is weakly constrained but mainly seems to depend on lava composition. In this context, the net effect of this important variable is shown by assuming a value in the range usually assumed for Etnean lavas ŽKilburn, 1993; Dragoni, 1989.. In detail, Fig. 7a illustrates the timewise growth of the crust and thermal mushy layer with values of 0 and 3.5 = 10 5 Jrkg for the crystallization heat. Fig. 7b shows the corresponding temperature profiles within the lava at different times. From the figure we can note that crystallization heat strongly limits the thickness of the crust and mushy layer, largely reducing the transmission of the cooling front inside the bulk. Such behavior can be easily explained in terms of an effective diffusivity

Fig. 7. Ža. Timewise distribution of the crust and mushy layer thickness during the first 2 h of cooling for crystallization heats equal to 0 and 3.5=10 5 Jrkg. The surface is characterized by eeff s 0.9, l s10, and d s1.22. Other model parameters are w s 0.1, and feu s1.0. Žb. Temperature profiles inside the lava at 60, 300, 900, 3600, and 7200 s for the two values of crystallization heat.

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

236

Fig. 8. Ža. Timewise distribution of crust and mushy layer thickness during the first 2 h of cooling for crystallization heats distributed over a range of temperature or localized at the solidus temperature Žpresence of an eutectic.. The surface is characterized by eeff s 0.9, l s 10, and d s 1.22. Other model parameters are w s 0.1, and L s 3.5 = 10 5 Jrkg. Žb. Temperature profiles inside the lava at 60, 300, 900, 3600, and 7200 s for the two kinds of heat release.

acting in the mushy layer. In fact, in this layer, the above reported diffusivity of Eq. Ž34., has to be corrected, due to the release of crystallization heat, according to the expression:

a eff s

1

rL

1 q

a

ž / kDT

Ž 35 .

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

237

where the denominator coincides with the variable Z reported in Eq. Ž23. if we assume feu s 1. In our application, the thermal diffusivity is reduced from about 5 = 10y7 m2rs, corresponding to a null latent heat, down to 2.3 = 10y7 m2rs. As already mentioned above, such a variation strongly affects the thickness of the two investigated layers whereas only a small influence is observed on the superficial temperature. The effects of a different release route of the crystallization heat is investigated in Fig. 8, where the presence of a eutectic has been supposed. In such a case, the heat release is all localized at the solidus temperature. In other words, we no longer consider a phase change over a range of temperatures, but only a phase change at a fixed temperature. This case can be directly analyzed by setting feu s 0 instead of 1, as so far assumed. Fig. 8a illustrates again the timewise behavior of the crust and mushy layer thicknesses, whereas Fig. 8b shows the temperature profiles. The model results show that the localized release of heat reduces the crust thickness and simultaneously increases the thickness of the mushy layer by about the same amount, thus determining no major changes in the propagation of the cooling front within the lava. Greater effects are observed in the temperature distribution that now shows a break-in-slope in correspondence with the solidus temperature, whereas only a limited change is observed in the surface temperature. At the crust–mushy layer interface, the temperature difference between the two conditions can range from about 40–508C, after a few minutes of cooling, to almost 1008C after 2 h.

4. Summary and conclusions The study performed allowed investigation of the main variables and parameters affecting the local cooling process of a thermal diffusion-dominated lava flow. The analysis particularly focused on the effect of surface morphology on the radiative and convective heat fluxes leaving the lava surface, taking into account the several non-linearities of the process. The solution method of the heat transfer equations allowed description of the evolution in time of the crust and thermal mushy layer as well as the behavior of their temperature profiles. The relative importance of radiative and convective contributions to the total heat flux was estimated during the first 2 h of cooling and a description of the associated heat transfer coefficients was given trying to better define such critical effects. Results suggest relevant effects of the superficial morphology. Surface features can strongly affect the global heat transfer between lava and atmosphere as well as the local cooling history of the flow. In particular, natural and forced convections appear to play a very relevant role from the very first minutes of cooling. The transient one-dimensional heat transfer model also allowed investigation of the influence of the main physical and thermal properties of the lava on the cooling process. The effect of thermal conductivity, vesicularity, and crystallization dynamics appears particularly relevant. With these variables it is possible to estimate the thermal diffusivity and inertia of the lava which control, respectively, the propagation velocity of the cooling front within the bulk and the cooling rate of its surface. Finally, it is important to outline that, in spite of the simple analysis performed, the study presented can already give an idea of the complexity of heat transfer processes occurring during lava cooling and can be useful in planning future investigations and research. Further work should focus on a more complete description and modeling of specific and local phenomena observed in the lava field, such as for instance the growth or decrease in volume of the investigated lava flow parcel and the linking of the heat transfer mechanisms described here with the fluid-dynamics of the flow. In this way, new important and non-intuitive effects could be discovered and a more realistic definition of the heat transfer effects would be achieved. For Etnean lavas, it is also particularly urgent to make an accurate definition of their thermal properties as a function of composition and temperature. Greater attention needs to be focused on magma emissivity, thermal conductivity, and crystallization dynamics. Laboratory measurements of these variables are strongly recommended and would allow a significant reduction in the number of variables needed in the model applications. On the other hand, field measurements appear to be fundamental in the description of real morphological features and in the validation of model results.

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

238

5. Nomenclature

a A A) b B B) c C C1,2 D

Geometrical parameter Integration constant of temperature profile Variable defined by Eq. Ž30. Geometrical parameter Integration constant of temperature profile Defined by Eq. Ž30. Specific heat Integration constant of temperature profile Constants of Eq. Ž2. Characteristic length

Tsol Tliq x X Y Z

fs g G Gr h H I k ks kU, kL

Crystal fraction Gravity acceleration Volumetric generation term Grashoff number Heat transfer coefficient Characteristic length Thermal inertia Thermal conductivity Roughness characteristic dimension Thermal conductivities reported in Eqs. Ž15. and Ž16. Geometrical parameter Crystallization latent heat Defined by Eq. Ž43. Defined by Eqs. Ž41. and Ž42. Nusselt number Pressure Defined by Eqs. Ž44. and Ž45. Prandtl number Net heat flux Radiosity Irradiation Defined by Eqs. Ž46. and Ž47. Rayleigh number Critical Rayleigh number Stanton number Time Temperature Surface temperature Initial temperature of lava Mean radiant temperature Mean equivalent temperature

b1, b 2 d e ´ w l m r s q Q

l L M1,2,3 N1,2,3 Nu P P1,2 Pr q qq qy Q1,2 Ra Rac St t T Tc Tl,` Tmr Tme

Greek a b

Solidus temperature Liquidus temperature Spatial coordinate Lava flow characteristic length Defined by Eq. Ž24. Defined by Eq. Ž23.

Thermal diffusivity, absorption coefficient Thermal expansion coefficient Constants of Eq. Ž32. Constant of Eq. Ž6. Emissivity Thickness of layer Vesicularity Žporosity. Dimensional constant of Eq. Ž6. Viscosity Density Stefan–Boltzman constant V-groove opening angle Auxiliary variable

Subscripts 1 Relative to period 1 2 Relative to period 2 atm Atmospheric cont Continuum phase Žmelt plus crystals. conv Convective Žnatural and forced. eff Effective eu Eutectic phase forc Forced convection Žwind. g Gas phase m Mushy region nat Natural convection rad Radiative s Solid phase sun Solar surf Superficial surr Surrounding tot Total

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

239

Acknowledgements The study was performed within the context of the Volcanic Simulation Group of the Gruppo Nazionale per la Vulcanologia, Consiglio Nazionale delle Ricerche, Italy. I wish to thank S.M. Baloga for the helpful review and F. Dobran, C. Kilburn, G. Macedonio, and P. Papale for many interesting discussions during the work. I also like to thank A. Bilancieri for drawing some figures. Funds were provided by the Commission of the European Communities, Contract no. CEC-EV5V-CT92-0190, and Gruppo Nazionale per la Vulcanologia.

Appendix A In the following, the application of the integral heat balance method to the two investigated periods of lava cooling is described to the purpose to illustrate the method and allow the reproduction of results. Period 1. In order to solve the integral Eq. Ž18. it is necessary to express its right-hand term, by the use of Leibnitz’s rule, in the form: ´m

H0

E Tm 1 Et

d xs

Qm 1 dt

d ´m

y DTl ,`

dt

.

where Qm1 is an auxiliary variable defined by:

Qm 1 s

´m

H0

Tm 1d x.

Ž 36 .

Eq. Ž18. can thus be expressed by substitution of the previous Eqs. Ž19. and Ž36. and taking into account boundary condition Ž17.b., as: Z

EQm 1 Et

y ZDT l ,`

d ´m

q

dt

dTm 1 dx

s0

Ž 37 .

0

where Z is reported in Eq. Ž23.. Assuming a quadratic temperature profile of the form: Tm 1 s A m 1 x 2 q Bm 1 x q Cm 1

Ž 38 .

and using boundary condition Ž17.a. the temperature profile reported, in dimensionless form, by Eq. Ž21. is obtained. With regards to the quadratic approximation, it has to be noted that the assumption of cubic or quadric temperature profiles slightly improves the accuracy of the results but greatly increases the computational manipulation. As a consequence, for our purpose, the assumption of a quadratic profile for each thermal layer Ži.e., crust and mushy layer. appears to be fully suitable. After determining Qm1 by Eq. Ž36., and its derivative with respect to time, Eq. Ž37. can be rewritten as: Z

Ž T y Tliq . 3 c

d ´m dt

Z q 3

´m

dTc dt

y

2 Ž Tc y T l ,` .

´m

s0

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

240

where Tc is an arbitrary function of time. Such a degree of freedom can be eliminated by using boundary condition Ž17.a. and solving with respect to the superficial temperature. Expressing the total flux as:

E Tm 1

qtot s k m

Ex

s

y2 k m Ž Tc y Tliq .

´m

0

we can obtain ´ m and its derivative with respect to time and, by substitution of dimensionless variables, the final ordinary differential equation in Tc reported in Eq. Ž20.. Period 2. Similarly to the procedure shown for Period 1 the auxiliary variables Qs2 and Qm2 can be introduced:

Qs 2 s

´s

Ts 2 d x ; Qm 2 s

H0

H´ ´ q´ Ž

s

m.

Tm 2 d x.

Ž 39 .

s

Assuming a quadratic temperature distribution in the two layers and by using the prescribed boundary conditions the expressions reported in Eqs. Ž31. – Ž33. are obtained. For the temperature distribution in the solid layer, constants can be obtained from boundary condition Ž25.b. and the equation: d ´s E Ts E x sy dt E t E Ts obtained by differentiation of Ts2 Ž ´sŽ t .,t . at x s ´s . From the knowledge of the temperature profiles, it is possible now to compute the auxiliary variables of Eq. Ž39., their derivatives with respect to time, and all other derivatives appearing in Eqs. Ž27. and Ž28.. Finally, we obtain the following system of ordinary differential equations in ´s and ´ m in which the superficial temperature Tc is again an arbitrary temperature:

°N d ´ q N d ´ s

~

m

s N3 dt d ´m M1 q M2 s M3 dt dt 1

¢

2

dt d ´s

Ž 40 .

where: 1

N1 s

6

N2 s y 1

N3 s

´s q

k´

Ž b1 q b 2 . y 2Tc q ´ s m

1 k ´s2 6

´s

M2 s

L c 1

1q

1 q feu Lfeu

3 cTliq

b2

q L) ,

b1

ž / 1q

Ž 41 .

b2

y2Tc y Ž b 1 q b 2 . y

3 M1 s

´m2

b1

ž / 1q

L)

dTc

2 b2

dt

ž

1 Tliq

y1

/

2k

´m

Ž 42 . q 1,

a q1 ,

M3 s

k ´s

Ž b1 q b 2 . .

Ž 43 .

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

241

Now the superficial temperature Tc can be constrained by boundary condition Ž25.a. that in this case becomes: qtot sk s

E Ts 2 Ex

ks

s

´m

0

2 Ž Tsol y Tc . y Ž b 1 q b 2 . .

Resolving this equation for ´ m and introducing dimensionless variables we obtain the expression reported in Eq. Ž30.. Deriving ´m with respect to time, we obtain: d ´m dt

s P1

d ´s dt

q P2

dTc dt

with: P1 s

P2 s

k 2 qtot

2 qtot k2

A) q

B)

B)2

k 2 2 q k ´s Y

ž

2

/q

B)

A) B)2

y

2Tc2

´s2

q

L) Tc

Ž 44 .

´s2

qtot Y´s k 2 q

4Tc

´s

L) y

´s

y L) kY q 2YkTc q 2 kqtot

Ž 45 .

and by substitution of Eq. Ž31. in the system Ž40. we obtain the final system of differential Eq. Ž29.. Finally, it should be noted that in this form N3 has been decomposed into: N3 s Q1 q Q2 where: Q1 s

Q2 s

1

´s ´s 3

y2Tc y Ž b 1 q b 2 . y

1q

L)

dTc

2 b2

dt

2k

´m

Ž 46 .

Ž 47 .

References Al-Arabi, M., El-Rafaee, 1978. Heat transfer by natural convection from corrugated plates to air. Int. J. Heat Mass Transfer 21, 357–359. Aubele, J.C., Crumpler, L.S., Elston, W.E., 1988. Vesicle zonation and vertical structure of basalt flows. J. Volcanol. Geotherm. Res. 35, 349–374. Baloga, S., 1987. Lava flows as kinematic waves. J. Geophys. Res. 92, 9271–9279. Birch, F., Clark, H., 1940. The thermal conductivity of rocks and its dependence upon temperature and composition. Am. J. Sci. 238, 529–558. Brandeis, G., Jaupart, C., Allegre, C.J., 1984. Nucleation, crystal growth, and the thermal regime of cooling magmas. J. Geophys. Res. 89, 10161–10177. Calvari, S., Coltelli, M., Neri, M., Pompilio, M., Scribano, V., 1994. The 1991–1993 Etna eruption: chronology and lava flow-field evolution. Acta Volcanol. 4, 1–14. Cashman, K.V., 1993. Relationship between plagioclase crystallization and cooling rate in basaltic melts. Contrib. Mineral. Petrol. 113, 126–142.

242

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

Chester, D.K., Duncan, A.M., Guest, J.E., Kilburn, C.R.J., 1985. Mount Etna. Chapman & Hall, London, 300 pp. Colburn, A.P., 1933. A method for correlating forced convection heat transfer data and a comparison with fluid friction. Trans. Am. Inst. Chem. Eng. 29, 174–210. Coulson, J.M., Richardon, J.F., 1977. Chemical Engineering: I. Fluid Flow, Heat Transfer, and Mass Transfer. Pergamon, Oxford, p. 450. Crisp, J., Baloga, S., 1990. A model for lava flows with two thermal components. J. Geophys. Res. 95, 1255–1270. Crisp, J., Kahle, A.B., Abbott, E.A., 1990. Thermal infrared spectral character of Hawaiian basaltic glasses. J. Geophys. Res. 95, 657–669. Dragoni, M., 1989. A dynamical model of lava flows cooling by radiation. Bull. Volcanol. 51, 88–95. Dragoni, M., Tallarico, A., 1994. The effect of crystallization on the rheology and dynamics of lava flows. J. Volcanol. Geotherm. Res. 59, 241–252. Dobran, F., 1995. Magma and lava flow modeling and volcanic system definition aimed at hazard assessment at Etna. Commission of the European Communities, Final Report on Contract CEC-EV5V-CT92-0190, Giardini, Pisa. Edwards, D.K., 1983. Heat Transfer by Radiation. Hemisphere Publishing, New York. Fink, J.H., Griffith, R.W., 1990. Radial spreading of viscous gravity currents. J. Fluid Mech. 221, 485–501. Goodman, T.R., 1964. Application of integral method to nonlinear heat transfer. In: Irvine, T.F., Hartnett, J.P. ŽEds.., Advances in Heat Transfer, Vol. 1. Academic Press, New York, pp. 52–120. Griffith, R.W., Fink, J.H., 1992a. The morphology of lava flows in planetary enviroments: prediction from analog experiments. J. Geophys. Res. 97, 739–748. Griffith, R.W., Fink, J.H., 1992b. Solidification and morphology of submarine lavas: a dependence on extrusion rate. J. Geophys. Res. 97, 19729–19737. Hashin, Z., Shtrikman, S., 1962. A variational approach to the theory of effective magnetic permeability of multiphase material. Appl. Phys. 33, 3125–3131. Head, J.W., Wilson, L., 1986. Volcanic processes and landforms on Venus: theory, prediction, and observation. J. Geophys. Res. 91, 9407–9446. Hon, K., Kauahikaua, J., Denlinger, R., McKay, K., 1994. Emplacement and inflation of pahoehoe sheet flows: observations and measurements of active lava flows on Kilauea Volcano, Hawaii. Geol. Soc. Am. Bull. 106, 351–370. Horai, K., 1991. Thermal conductivity of Hawaiian basalt: a new interpretation of Robertson and Peck’s data. J. Geophys. Res. 96, 4125–4132. Hottel, H.C., Sarofim, A.F., 1967. Radiative Transfer. McGraw-Hill, New York, p. 520. Jakob, M., 1957. Heat Transfer. Wiley, New York, p. 411. Jones, A.C., 1993. The cooling rates of pahoehoe flows: the importance of lava porosity. Lunar. Planet. Sci. Conf. 24, 731–732. Kahle, A.B., Gillespie, A.R., Abbott, E.A., Abrams, M.J., Walker, R.E., Hoover, G., Lockwood, J.P., 1988. Relative dating of Hawaiian lava flows using multispectral thermal infrared images: a new tool for geologic mapping of young volcanic terrains. J. Geophys. Res. 93, 239–251. Keszthelyi, L., 1994. Calculated effect of vesicles on the thermal properties of cooling basaltic lava flows. J. Volcanol. Geotherm. Res. 63, 257–266. Keszthelyi, L., Denlinger, R., 1996. The initial cooling of pahoehoe flow lobes. J. Volcanol. Geotherm. Res. 58, 5–18. Kilburn, C.R.J., 1989. Surfaces of aa flow-fields on Mount Etna, Sicily: morphology, rheology, crystallization and scaling phenomena. In: Fink, J.H. ŽEd.., Lava Flows and Domes. IAVCEI Proceeding in Volcanology, Vol. 2, pp. 129–156. Kilburn, C.R.J., 1993. Lava crust, Aa flow lengthening and the pahoehoe-aa transition. In: Kilburn, C.R.J., Luongo, G. ŽEds.., Active Lava Flows: Monitoring and Modeling. UCL Press, London, pp. 263–280. Kilburn, C.R.J., 1996. Patterns and predictability in the emplacement of subaerial lava flow and flow fields. In: Tilling, R.I., Scarpa, R. ŽEds.., Monitoring and Mitigation of Volcano Hazards. Springer, Berlin, pp. 491–537. Kilburn, C.R.J., Lopes, R.M.C., 1991. General patterns of flow field growth: aa and blocky lavas. J. Geophys. Res. 96, 19721–19732. Knudsen, J.C., Bell, K.J., Holt, A.D., Hottel, H.C., Sarofim, A.F., Standiford, F.C., Stuhlborg, D., Uhl, V.W., 1984. Heat transmission. In: Perry, R.H. ŽEd.., Perry’s Chemical Engineers’ Handbook. McGraw-Hill, New York. Murase, T., McBirney, A.R., 1973. Properties of some common igneous rocks and their melts at high temperature. Geol. Soc. Am. Bull. 84, 3563–3592. Macedonio, G., 1995. A 3-dimensional model for lava flows. In: Dobran, F. ŽEd.., Magma and Lava Flow Modeling and Volcanic System Definition Aimed at Hazard Assessment at Etna. Commission of the European Communities, Final Report on Contract CEC-EV5VCT92-0190, Giardini, Pisa. Neri, A., 1994. A heat transfer analysis of lava cooling in the atmosphere. In: Barberi, F., Casale, R., Fratta, M. ŽEds.., The European Laboratory Volcanoes. Workshop Proceedings, Acicastello ŽCatania., June 18–21, p. 171. Neri, A., 1995. A heat transfer analysis of lava cooling in the atmosphere. In: Dobran, F. ŽEd.., Magma and Lava Flow Modeling and Volcanic System Definition Aimed at Hazard Assessment at Etna. Commission of the European Communities, Final Report on Contract CEC-EV5V-CT92-0190, Giardini, Pisa. Oppenheimer, C., 1991. Lava flow cooling estimated from Landsat thematic mapper infrared data: the Lonquimay eruption ŽChile, 1989.. J. Geophys. Res. 96, 21865–21878.

A. Neri r Journal of Volcanology and Geothermal Research 81 (1998) 215–243

243

¨ Ozisik, M.N., 1993. Heat Conduction. Wiley, New York, p. 692. Peck, D.L., Hamilton, M.S., Shaw, H.R., 1977. Numerical analysis of lava lake cooling models: II. Application to Alae lava lake, Hawaii. Am. J. Sci. 277, 415–437. Pinkerton, H., Sparks, R.S.J., 1976. The 1975 sub-terminal lavas, Mount Etna: a case history of the formation of a compound lava field. J. Volcanol. Geotherm. Res. 1, 167–182. Polacci, M., Papale, P., 1997. The evolution of lava flows from ephemeral vents at Mount Etna: insights from vesicle distribution and morphological studies. J. Volcanol. Geotherm. Res., in press. Raithby, G.D., Hollands, H.G., 1984. Natural convection. In: Rohsenow, W.M., Hartnett, J.P., Ganic E.N. ŽEds.., Handbook of Heat Transfer Fundamentals, 2nd edn., Chap. 6. McGraw-Hill, New York. Rowland, S.K., Walker, G.P.L., 1987. Toothpaste lava: characteristics and origin of a lava structural type transitional between pahoehoe and aa. Bull. Volcanol., Vol. 49, pp. 631,641. Rubesin, M.W., Inouye, M., Parikh, P.G., 1984. Forced convection: external flows. In: Rohsenow, W.M., Hartnett, J.P., Ganic, E.N. ŽEds.., Handbook of Heat Transfer Fundamentals, 2nd edn., Chap. 8. McGraw-Hill, New York. Schlichting, H., 1968. Boundary Layer Theory. McGraw-Hill, New-York, p. 748. Sparrow, E.M., Jonsson, V.K., 1963. Thermal radiation absorption in rectangular-groove cavities. J. Appl. Mech., pp. 237–244. Sparrow, E.M., Lin, S.H., 1963. Absorption of thermal radiation in a V-groove cavity. Int. J. Heat Mass Transfer 5, 1111–1115. Stasiuk, M.V., Jaupart, C., Sparks, R.S.J., 1993. Influence of cooling on lava-flow dynamics. Geology 21, 335–338. Touloukian, Y.S., Judd, W.R., Roy, R.F., 1989. Physical Properties of Rocks and Minerals. Hemisphere Publishing, New York, p. 548. Turcotte, D.L., Schubert, G., 1982. Geodynamics. Applications of Continuum Physics to Geological Problems. Wiley, New York, p. 450. Walker, G.P.L., 1989. Spongy pahoehoe in Hawaii: a study of vesicle distribution patterns in basalt and their significance. Bull. Volcanol. 51, 199–209. Wilmoth, R.A., Walker, G.P.L., 1993. P-type and S-type pahoehoe: a study of vesicle distribution patterns in Hawaiian lava flows. J. Volcanol. Geotherm. Res. 54, 129–142. Williams, H., McBirney, A.R., 1979. Volcanology. Freeman, Cooper, San Francisco, p. 396. Zimbelman, J.R., 1986. The role of porosity in thermal inertia variations on basaltic lavas. Icarus 68, 366–369.