Bruno Scaillet
email: bscaille@cnrs-orleans.fr

Alan Whittington

Caroline Martel

Michel Pichavant

François Holtz

Francois Holtz

Phase equilibrium constraints on the viscosity of silicic magmas with implications for mafic-silicic mixing processes

come

Phase equilibrium constraints on the viscosity of silicic

magmas with implications for mafic-silicic mixing processes

Abstract. Isobaric crystallization paths obtained from phase equilibrium experiments show that, whereas in rhyolitic compositions melt fraction trends are distinctly eutectic-like, dacitic and more mafic compositions have their crystallinities linearly correlated with temperature. As a consequence, the viscosities of the latter continuously increase on cooling, whereas for the former they remain constant or even decrease during 80% of the crystallisation interval, which opens new perspectives for the fluid dynamical modeling of felsic magma chambers. Given the typical dike widths observed for basaltic magmas, results of analogue modelling predict that injection of mafic magmas into crystallizing intermediate to silicic plutons under preeruption conditions cannot yield homogeneous composition. Homogenization can occur, however, if injection takes place in the early stages of magmatic evolution (i.e. at near liquidus conditions) but only in magmas of dacitic or more mafic composition. More generally, the potential for efficient mixing between silicic and mafic magmas sharing large interfaces at upper crustal levels is greater for dry basalts than for wet ones. At the other extreme, small mafic enclaves found in many granitoids behave essentially as rigid objects during a substantial part of the crystallization interval of the host magmas which implies that finite strain analyses carried out on such markers can give only a minimum estimate of the total amount of strain experienced by the host pluton. Mafic enclaves carried by granitic magmas behave as passive markers only at near solidus conditions, typically when the host granitic magma shows near-solid behavior. Thus they cannot be used as fossil indicators of direction of magmatic flow. key words: granite/rhyolite/dacite/andesite/basalts/mixing/viscosity/phase equilibria Viscosity is one fundamental parameter that controls the physical and chemical behavior of magmatic rocks. Much effort has been invested in the past 30 years to evaluate the magnitude of this parameter in terrestrial magmas as well as its major controlling factors. Since the pioneering study of [START_REF] Bottinga | The viscosity of magmatic liquids: a model for calculations[END_REF], soon modified by [START_REF] Shaw | Viscosities of magmatic silicate liquid: an empirical method of prediction[END_REF] to incorporate the effect of water, there have been several works aimed at constraining the viscosity of various silicate liquids under a wide range of P-T-melt H 2 O content conditions (see [START_REF] Lange | The effect of H 2 O, CO 2 and F on the density and viscosity of silicate melts[END_REF]). However, it has been mainly in the last five years that interest in the viscous behavior of magmas, particularly of silicic compositions, has been renewed. This stems partly from the recognised important role that viscosity plays on eruptive phenomena. It is also due to recent technological improvements made both at low and high pressures that have opened new perspectives in the field of laboratory measurement of viscous properties of silicate liquids. In particular, the development of a high perfomance furnace operated in internally heated vessels (Roux et al., 1996) has allowed viscosity determinations using the falling sphere method with greater accuracy than was possible in previous studies (e.g., [START_REF] Scaillet | Viscosity of Himalayan leucogranites: Implications for mechanisms of granitic magma ascent[END_REF][START_REF] Schulze | The influence of water on the viscosity of a haplogranite liquid[END_REF][START_REF] Holtz | The effects of silica and water on the viscosity of hydrous quartzofeldspathic melts[END_REF]. Similarly, the application of the parallel plate viscosimeter tool to hydrous glass cylinders synthesised at high pressures has allowed viscosity determinations in the high viscosity range (e.g., [START_REF] Richet | Water and the viscosity of andesite melts[END_REF][START_REF] Dingwell | The effect of water on the viscosity of a haplogranitic melt under P-T-X conditions relevant to silicic volcanism[END_REF], i.e., in a range not accessible by the more conventional, but time consuming, falling sphere method. Combination of both approaches enables the characterisation of the viscous behavior of supraliquidus liquids over very large temperature ranges, which in turn allows derivation of thermodynamic properties, at least for some, well behaved, silicate liquids (e.g., diopside; Richet, 1984). These technological improvements have also benefited from recent analytical developments that enable full chemical characterisation of the run products before and after the viscosity determination. In particular, the melt H 2 O content can be now routinely analysed accurately, using for instance the Karl Fischer titration technique, and with a high spatial resolution through microinfra red spectroscopy analysis (FTIR).

From a more geological perspective, it is well established that the viscosity of magmas depends mainly on temperature, melt water content, melt composition and melt fraction. All of these parameters can vary widely and there is not as yet empirical or thermodynamic model able to predict the controlling factors of the viscosity for all classes of natural magmas. Direct measurement of viscous properties of natural magmas (e.g., [START_REF] Pinkerton | Rheological properties of basaltic lavas at sub-liquidus temperatures: laboratory and field measurements on lavas from Mount Etna[END_REF] are difficult, hazardous and intrinsically restricted to extruded magmas. For plutonic magmas, a laboratory approach is required and we must rely on the knowledge gained from a restricted number of case studies and, provided that these are representative, attempt to find the general rules or behavior of the rheological properties of magmas. As stated above, laboratory investigations have mainly focused on silicate liquids and only a limited number of studies has been devoted to the characterisation of real magmas, i.e. of crystal+liquid+fluid mixtures (e.g., [START_REF] Van Der Molen | Experimental deformation of partially melted granite[END_REF][START_REF] Arzi | Critical phenomena in the rheology of partially melted rocks[END_REF]Rutter & Newman, 1995[START_REF] Rushmer | An experimental deformation study of partially molten amphibolite: Application to low-melt fraction segregation[END_REF]. These studies, however, were mainly conducted in the low melt fraction range, with the aim of understanding the rheological behavior of the partially molten middle to lower crust. Those studies carried out at high melt fractions were aimed at characterizing volcanological, or shallow, regimes [START_REF] Stevenson | Viscosity of microlite-bearing rhyolitic obsidians: an experimental study[END_REF].

In contrast, the rheological behavior of magmas during crystallisation at depth remains a field largely unexplored. Although the first aim of experimental phase equilibria is not to constrain the rheology of magmas, this approach is one of the most powerful that can be used to evaluate all the critical factors mentioned above. It can be used on both volcanic and plutonic rocks, and for the latter it is the only one that allows the retrieval of magmatic evolution conditions. Its main drawback, however, is that this is a static or near-equilibrium approach by definition, and as such it does not incorporate any potential effect of deformation upon crystallisation. Nevertheless, combining the information gained from phase equilibrium studies on natural rocks together with empirical rheological laws obtained from the study of model or analogue systems is a useful first step for understanding the rheological properties of magmas under relevant P-T-melt H 2 O conditions. In this paper we take this approach, with particular emphasis on the relative rheological behavior of mafic versus silicic magmas. We build upon previous work that has mainly focused on silicic magma rheology, also using a phase equilibrium approach [START_REF] Scaillet | Rheological properties of granitic magmas in their crystallization range[END_REF]1998a). We start this report by first reviewing our current knowledge about the viscosity of hydrous silicate melts as well as the effects of both bubbles and crystals on bulk viscosity. We then use recent phase equilibrium studies carried out on silicic to intermediate magmas (granite/rhyolite to andesite) as proxies for constraining the viscosity patterns of silicic to intermediate cooling magmas. These patterns are then compared against those of basaltic magma crystallizing either under dry or wet conditions to explore the most favourable conditions leading to thorough mafic-silicic magma mixing. Finally, we discuss the rheological behavior of mafic enclaves found in granitoids in which they are used as strain gauges of the finite strain that affected the host plutonic body.

Viscosity models and their application to magmas

From the rheological standpoint any magma can be usefully considered as a three phase system made of a mixture of melt, bubbles of exsolved volatiles and crystals in varying proportions. The rheological behavior of such a polyphase system will result from the interplay of the dynamical response of each of these phases with respect to any applied stress.

Decisive progress in the determination of the viscosity of hydrous silicate melts has been made in recent years and several empirical equations are now available (see Whittington et al., 1999a). These may be valid over a rather large compositional range (e.g., Hess & Dingwell, 1996;[START_REF] Baker | Granitic melt viscosity and dike formation[END_REF][START_REF] Holtz | The effects of silica and water on the viscosity of hydrous quartzofeldspathic melts[END_REF] or apply to a more specific (yet widespread) group of granitic magmas [START_REF] Scaillet | Viscosity of Himalayan leucogranites: Implications for mechanisms of granitic magma ascent[END_REF]Whittington et al. 1999a, for peraluminous granites). These empirical formulations supersede the classical model of [START_REF] Shaw | Viscosities of magmatic silicate liquid: an empirical method of prediction[END_REF], especially in the low melt water content range. These new data have shown that hydrous silicic melts are distinctly non-Arrhenian over large temperature ranges (i.e. > 1000°C) although notable exceptions do exist (ie hydrous peraluminous silicic magmas, see Whittington et al. 1999a).

However, silicic magmas sensu lato display fairly smaller melting/crystallization temperature intervals, typically 100-200°C, over which the temperature dependence of the melt viscosity can be considered as Arrhenian. The melt viscosity is also strongly dependent on the strain rate: at low strain rates, silicic melts display a newtonian behavior, while at high strain rates, non-newtonian behavior is observed, and shear viscosities decrease with increasing strain rates (Webb & Dingwell, 1990). The low strain rates prevailing in plutonic environments imply that the melt in magmas at depths obeys predominantly a Newtonian behavior. However, magmas transported through dykes, be it from the source region or from the plumbing system of volcanoes, may well experience high strain rates at which the onset of non-newtonian behavior (Webb & Dingwell, 1990) is observed. Applications of the above empirical models to such cases will likely give a maximum value for the melt viscosity. The restored fluid content of most volcanic rocks falls in the range 1-5 wt.% [START_REF] Wallace | Quantification of preeruptive exsolved gas contents in silicic magmas[END_REF]Scaillet et al., 1998b). Available experimental data indicate that the role of bubbles in magma rheology is critically dependent on the strain rate (Stein & Spera, 1993;[START_REF] Bagdassarov | Rheological investigations of vesicular rhyolite[END_REF]. Assuming that volcanic rocks are the extruded equivalents of plutonic magmas, the low fluid content of the former suggests that the rheological role of bubbles in magmas at depth is probably minor (Scaillet et al., 1998a). Only during eruption, or emplacement at very shallow levels, is the bubble fraction of a magma an important rheological factor (see [START_REF] Dingwell | The effect of water on the viscosity of a haplogranitic melt under P-T-X conditions relevant to silicic volcanism[END_REF].

The influence of crystals on magma viscosity is still a subject of active research. Experiments carried out on model systems have shown that the Eintein-Roscoe equation satisfactorily reproduces the data for crystal contents < 30% in volume [START_REF] Lejeune | Rheology of crystal-bearing silicate melts: An experimental study at high viscosities[END_REF]. Beyond a crystallinity of 30%, magma viscosities are generally observed to increase rapidly up to reaching solid behavior at around 60% crystals. Phase equilibrium studies carried out on silicic magmas have conclusively shown that they arrive in a nearly molten state at their emplacement level, their load of crystals barely exceeding 10% (e.g., [START_REF] Clemens | Crystallization and origin of some peraluminous (S-type) granitic magmas[END_REF][START_REF] Scaillet | Experimental crystallization of leucogranite magmas[END_REF][START_REF] Dall'agnol | An experimental study of a lower proterozoic A-type granite from the Eastern Amazonian Craton, Brazil[END_REF][START_REF] Scaillet | Phase relations of peralkaline silicic magmas and petrogenetic implications[END_REF]. This indicates that during extraction and ascent, the rheology of silicic magmas is mostly controlled by the melt properties. A recent compilation of silicic magma viscosities has shown that melt viscosities of felsic magmas tightly cluster at 10 4.5 Pa s, irrespective of the temperature, melt water contents and the plutonic or volcanic nature of the rock (Scaillet et al., 1998a). The main reason for such a constant viscosity is that temperature and melt H 2 O content are interrelated parameters, inversely correlated to each other (Fig. 1). Peralkaline, peraluminous and metaluminous magmas plot on the same general trend which indicates that the eruptability of silicic to intermediate magmas is not compositionally controlled. Similarly, the lack of difference between extrusive and intrusive rocks shows that the fate of silicic magmas, whether volcanic or plutonic, is not due to fundamental differences in viscous properties of the magma. In contrast, host rock rheological properties may exert a critical role on magma stoping upon ascent (e.g, [START_REF] Hogan | The A-type Mount Scott granite sheets: Importance of crustals magmas traps[END_REF]. Nevertheless, the fact that there are virtually no volcanic rocks with crystal contents higher than 50% [START_REF] Marsh | On the crystallinity, probability of occurence and rheology of lava and magma[END_REF]Scaillet et al., 1998a), corresponding to magma viscosities of 10 6 Pa s (Fig. 1), indicates that this is a fundamental rheological barrier, beyond which intermediate to silicic magmas hardly move upward.

In the following, unless otherwise stated, the [START_REF] Shaw | Viscosities of magmatic silicate liquid: an empirical method of prediction[END_REF] model has been used to calculate melt viscosities since it is still the only one available that can account for melt compositions more mafic than rhyolites (i.e., dacite, andesite and basalt), unlike all recently published models. In addition, despite the fact that several studies have shown that the [START_REF] Shaw | Viscosities of magmatic silicate liquid: an empirical method of prediction[END_REF] model may yield inaccurate results for some T-H 2 O conditions (e.g., [START_REF] Dingwell | The effect of water on the viscosity of a haplogranitic melt under P-T-X conditions relevant to silicic volcanism[END_REF], recent experimental data from hydrous acidic, intermediate and basic compositions suggests that this empirical model is quite accurate at high water contents and magmatic temperatures, although it is unreliable at viscosities higher than its intended upper limit of 10 6 Pa.s. For example, the viscosity of Crater Lake andesite containing 4 wt.% H 2 O at 1350˚C was measured as between 4 and 5 Pa.s by [START_REF] Kushiro | Viscosities of basalt and andesite melts at high pressures[END_REF], and is reproduced exactly by the Shaw model. Viscosities of synthetic phonolite and trachyte melts containing 5 wt.% H 2 O at 777 and 877 °C, extrapolated from lower temperature measurements (Whittington et al., 1999b) are reproduced to within 0.2 log units by the Shaw model. Similarly, the viscosity measurements performed on a hydrous peraluminous leucogranite melt by [START_REF] Scaillet | Viscosity of Himalayan leucogranites: Implications for mechanisms of granitic magma ascent[END_REF] are reproduced to within 0.2 log units by the Shaw model, for melt H 2 O contents between 3 and 6.6 wt.% and temperatures between 800 and 1100°C.

For the sake of simplicity, the Einstein-Roscoe equation has been used to calculate magma viscosities over the full crystallisation range, although this equation is strictly valid only at melt fractions in excess of 70 % (e.g. [START_REF] Lejeune | Rheology of crystal-bearing silicate melts: An experimental study at high viscosities[END_REF]. However, use of more appropriate equations (see [START_REF] Dingwell | Magma rheology, in: Short Handbook on Experiments at High Pressure and Applications to the Earth's Mantle[END_REF] introduces only minor changes in the calculated viscosity, that differ by less than an order of magnitude from those obtained with the Einstein-Roscoe approximation. Given that mafic and felsic magma viscosities may differ by several orders of magnitude, as illustrated below, such relatively minor deviations are unlikely to significantly alter any conclusions derived here.

Magma viscosity during crystallization

As phase equilibrium studies have shown, nearly all the crystallization of granitic magmas of crustal origin occurs isobarically, at the level of emplacement (see [START_REF] Wall | Models of granitoid evolution and source compositions[END_REF][START_REF] Clemens | Granitic magma transport by fracture propagation[END_REF]Scaillet et al., 1998a). Thus, phase equilibrium experiments are particularly well suited for constraining the rheological patterns of cooling silicic magmas.

Crystallization paths can be obtained firstly by recognizing that water has a nearly perfect incompatible behavior and secondly by assuming a closed system (e.g., [START_REF] Scaillet | Rheological properties of granitic magmas in their crystallization range[END_REF]. In all cases studied (see [START_REF] Scaillet | Rheological properties of granitic magmas in their crystallization range[END_REF], melt fraction trends, illustrated here for two compositions, a peraluminous and a metaluminous one, display a near-eutectic behavior, the magmas being still 60-80 wt.% liquid at 15 °C above their solidus (Fig. 2, see also [START_REF] Whitney | The origin of granite: The role and source of water in the evolution of granitic magmas[END_REF]). The melt and magma viscosity trends corresponding to these melt fraction trends are shown on figure 2. After a short increase, the melt viscosity decreases as crystallization proceeds. In the case of a magma beginning to crystallize with 7 wt.% H 2 O, at 750 °C and 400 MPa (peraluminous rhyolite), the viscosity decreases by nearly 1 order of magnitude during crystallisation (note that the viscosity of this specific case is calculated using the model of [START_REF] Scaillet | Viscosity of Himalayan leucogranites: Implications for mechanisms of granitic magma ascent[END_REF], calibrated on a peraluminous silcic melt, and not with the Shaw model: the decrease observed therefore not due to an inadequacy of the empirical model). In both cases, magma viscosities either remain constant or increase during cooling, but by less than a factor of 10 over more than 80 % of the crystallization interval. The reason for such behavior is simply the lowering effect of water upon viscosity that largely offsets, or at least counteracts, the opposite effect arising from cooling and increasing crystal content during crystallization. This decrease in viscosity can be even more dramatic if, in addition to water, F and B are present in the system. For instance, the addition of 1.5 wt.% F to a peraluminous liquid with 6 wt.% H 2 O decreases the melt viscosity by a factor of 4 [START_REF] Baker | The low viscosities of F+H 2 O-bearing granitic melts and implications for melt extraction and transport[END_REF]. This rheological pattern is likely to be common to most rhyolites/granites in so far as they have, by definition, a composition close to that of the eutectic and they are invariably water-rich (e.g., [START_REF] Whitney | The origin of granite: The role and source of water in the evolution of granitic magmas[END_REF]Scaillet et al., 1998a). These findings sharply contrast with the widely held belief that viscosity increases during the crystallization of rhyolitic magmas, and raise new questions about the fluid dynamical properties of felsic magma chambers. Although physical models have not yet considered the possibility of a viscosity decrease during crystallization (e.g., [START_REF] Spera | Thermal boundary layer convection in silicic magma chambers: effects of temperature dependent rheology and implications for thermogravitational chemical fractionation[END_REF][START_REF] Brandeis | On the interaction between convection and crystallization in cooling magmas chambers[END_REF][START_REF] Marsh | On the convective style and vigour in sheet-like magma chambers[END_REF] the low and constant viscosity of felsic magmas during most of their crystallization history opens the possibility that most of the fractionation mechanisms at work in basaltic magma chambers, such as crystal settling, compaction and thermal or solutal convections, can be effective in felsic magma chambers as well (see [START_REF] Petford | Porous media flow in granitoid magmas: an assessment, in Flow and creep in the solar system: Observations, modeling and theory[END_REF]. For instance the compaction length [START_REF] Mckenzie | The generation and compaction of partially molten rocks[END_REF] of magma with a melt fraction of 40 % and an interstitial melt viscosity of 10 3.3 Pa.s (corresponding in the case shown above to conditions 10 °C above the solidus for a granitic magma starting crystallizing with 7 wt.% H 2 O; Fig. 2) would exceed several hundred meters whilst the compaction time [START_REF] Mckenzie | The generation and compaction of partially molten rocks[END_REF], around 10 3 years, is well below the lifespan of felsic magma chambers several kilometers across [START_REF] Spera | Thermal evolution of plutons, a parametrized approach[END_REF]. Current physical models that postulate high and increasing viscosities for cooling silicic magmas may not, therefore, faithfully reproduce the fluid dynamics of natural felsic magma chambers. Similar crystallization paths have been calculated for a dacitic bulk composition (SiO 2 = 65 wt.%), using the experimental results obtained by [START_REF] Scaillet | The June 15, 1991 eruption of Mount Pinatubo. I. Phase equilibria and pre-eruption P-T-fO 2 -fH 2 O conditions of the dacite magma[END_REF] on the dacite of Mt Pinatubo, for which variations of melt fraction with T and melt H 2 O contents, as well as melt composition, are all well characterised from phase equilibrium experiments performed at 220 MPa. Two crystallization paths have been calculated, one for starting conditions set at 950°C with 4 wt.% H 2 O in melt, and the second for H 2 O-saturated conditions, with melt H 2 O contents of 6.9 wt.% at all temperatures (Fig. 3). Although some curvature is apparent, in both cases melt fraction trends display a much more linear relationship with temperature than those obtained with bulk granitic/rhyolitic compositions (Fig. 2). For the watersaturated case, the liquid viscosity rises steadily with lowering temperatures, having increased by nearly 2 orders of magnitude at midway between the liquidus and solidus (Fig. 3). As a consequence, magma viscosity increases throughout crystallisation, by 4 orders of magnitude at 50% of the crystallisation interval. Starting with under-saturated conditions (4 wt.% H 2 O at 950°C) has little effect. The main difference is that, after a short increase, the melt viscosity increases slowly over 80 % of the crystallisation interval, but the magma viscosity rises rapidly, by 5 orders of magnitude at 70 % of cooling (Fig. 3). Crystallisation paths calculated for an andesite composition using the phase equilibrium experiments on the Mt Pelée andesite obtained by [START_REF] Martel | Effects of fO 2 and H 2 O on andesite phase relations between 2 and 4 kbar[END_REF] are shown in Fig. 4 for two different initial melt H 2 O contents. Although not constrained in the near solidus range, the high temperature part of the crystallisation interval of andesitic magmas shows a linear correlation between melt fraction and temperature. As for the dacite, both the melt and magma viscosities continuously increase during cooling, by 1 and 3 orders of magnitude, respectively, during the first 40 % of the crystallisation interval (Fig. 4). Thus, at upper crustal levels, the rheological pattern of a dacitic/andesitic magmas does not resemble that of rhyolites, in that their viscosities increase during cooling, irrespective of the initial melt H 2 O content, although the increase remains moderate during the first half of the crystallization interval. Increasing pressure will increase H 2 O solubility, and it cannot be excluded that, beyond a certain depth, the effect of water will dominate that of an increase in both melt silica content and crystals, as seen for rhyolitic compositions [START_REF] Scaillet | Rheological properties of granitic magmas in their crystallization range[END_REF]. However, the viscosity reduction with increasing water content is less marked for more basic compositions than for rhyolitic compositions, at least at low temperatures near the glass transition (Whittington et al., 1999b;1999c). For the purposes of fluid dynamical modelling of shallow dacitic/andesitic magma chambers, the present data show that it is safe to assume that viscosity increases continuously during cooling. The experimental data available for more mafic compositions than silicic andesite does not allow us at present to carry out this type of modelling, but it is clear that, since basaltic magmas obey the same rule (i.e., linear decrease in crystal content with temperature), we can predict that such intermediate compositions will follow a rheological pattern also characterized by a continuous increase in viscosity during crystallization.

Implications for basalt-rhyolite mixing

Field evidence clearly shows that basaltic and felsic magmas often come into contact in the crust (e.g., [START_REF] Eichelberger | Origin of andesite and dacite: evidence of mixing at Glass Mountain in California and at other circum-Pacific volcanoes[END_REF][START_REF] Bacon | Magmatic inclusions in silicic and intermediate volcanic rocks[END_REF][START_REF] Barbarin | Genesis and evolution of mafic microgranular enclaves through various types of interaction between coexisting felsic and mafic magmas[END_REF][START_REF] Wiebe | Silicic magma chambers as traps for basaltic magmas: the Cadillac Mountain intrusive complex, Mount Desert Island, Maine[END_REF]. The most common situations are due either to basaltic magmas that pond at the base of the crust and subsequently melt lower continental lithologies (e.g., [START_REF] Huppert | The generation of granitic magmas by intrusion of basalt into continental crust[END_REF][START_REF] Atherton | Generation of sodium-rich magmas from newly underplated basaltic crust[END_REF] or to the injection of basaltic magmas into upper silicic magma chambers (e.g., Coulon et al., 1984;[START_REF] Nixon | Petrology of the younger andesites and dacites of Iztaccihuatl volcano, Mexico: II. Chemical composition of basaltic magma influx[END_REF][START_REF] Michael | Intrusion of basaltic magma into a crystallizing granitic magma chamber: the Cordillera del Paine pluton in southern Chile[END_REF][START_REF] Pallister | A basalt trigger for the 1991 eruptions of Pinatubo volcano?[END_REF][START_REF] Feeley | Compositional and dynamic controls on mafic-silicic magma interactions at continental arc volcanoes: evidence from Cordon El Guadal, Tatara-San Pedro complex, Chile[END_REF]. In both cases, the extent of mixing between the two end-members will depend, among other factors, on the viscosity contrast (e.g., [START_REF] Sparks | Thermal and mechanical constraints on mixing between mafic and silicic magmas[END_REF][START_REF] Campbell | The influence of viscosity on fountains in magma chambers[END_REF][START_REF] Frost | Field, chemical, and physical constraints on mafic-felsic interactions in the Lamarck granodiorite, Sierra Nevada, California[END_REF]Fernandez & Gasquet, 1994;[START_REF] Snyder | Replenishment of magma chambers: comparison of fluid-mechanic experiments with field[END_REF][START_REF] Jellinek | Mixing and compositional stratification produced by natural convection. 2. Applications to the differentiation of basaltic and silicic magma chambers and komatiite lava flows[END_REF]. Intermediate to silicic magmas have a characteristic magma viscosity of 10 4.5 Pa.s during emplacement or at pre-eruption conditions (Scaillet et al., 1998a). Trends displayed on figure 2 show that this viscosity value is maintained during a substantial part of the magmatic evolution of rhyolitic magmas. We can use this value to evaluate the effects of basalt injection into crystallizing intermediate to silicic magma chambers.

Analogue experiments carried out by [START_REF] Campbell | The influence of viscosity on fountains in magma chambers[END_REF] suggest that the effectiveness of mixing in turbulent fountains at the base of magma chambers can be characterized by the Reynolds number:

Re = w d / η (2)
where w is the velocity of the incoming fluid (basalt), d the width of the input pipe and η the viscosity of the ambient magma (rhyolite).When Re is >70, mixing is predicted to be very efficient, whereas for Re < 7 no mixing occurs at all. Adopting the set of values used by [START_REF] Campbell | The influence of viscosity on fountains in magma chambers[END_REF] for basaltic magma but taking a viscosity of 10 4.5 Pa.s obtained for silicic rocks (instead of 10 7 Pa.s), shows that this mixing criterion is met once the dike width exceeds 10 m (Re = 10). This is at the high end of the range measured for basaltic magmas, whose average dike width is close to 1 m (e.g, [START_REF] Wada | On the relationship between dike width and magma viscosity[END_REF]. Taking the average magma viscosity obtained for the pre-eruption conditions of volcanic rocks increases the minimum dike width for mixing up to 100 m. High intrusion rates may induce high strain rates at the interface between the two magmas, conditions that may promote the onset of non-Newtonian behavior in silicate melts (Webb & Dingwell, 1990;Dingwell & Webb, 1990) which is characterized by a decrease of melt viscosity as the strain rate increases (Webb & Dingwell, 1990). This decrease may be shared by both magmas, and the experiments do not indicate that acid compositions depart more significantly than mafic ones from non-newtonian behavior in such a way that the viscosity of felsic magma would be decreased more rapidly than the mafic one. Thus, the intrusion of a basalt magma into an upper solidifying silicic magma chamber will hardly result in a compositionally homogeneous intermediate magma. However, injection of basaltic magmas at early stages of the crystallization of a dacitic or andesitic magma may lead to mixing, since the near-liquidus viscosities of the latter are lower by 1-2 log units relative to those at pre-eruptive conditions (Fig. 1), decreasing the critical dike width by a factor of 10 to 100.

The previous analysis considers the case of mixing arising from vigorous injection of basaltic magma into crystallizing silicic ones, where the energy for mixing comes in part from the rate of intrusion of the incoming magma. There are many other instances where basalt and rhyolite are brought together in a more quiescent way, sharing large physical interfaces. Given that heat diffusion is many orders of magnitude faster than chemical diffusion in silicate melts (e.g., [START_REF] Watson | Diffusion in magmas at depth in the Earth: the effect of pressure and dissolved H 2 O[END_REF], the first effect of basalt juxtaposition is to force the silicic magma back up its liquid line of descent, without significant contamination (e.g., [START_REF] Sparks | Magma mixing: A mechanism for triggering acid explosive eruption[END_REF][START_REF] Sparks | Thermal and mechanical constraints on mixing between mafic and silicic magmas[END_REF]. In this case the viscosity contrast across the interface will be one critical factor in determining whether or not wholesale homogeneization can occur [START_REF] Sparks | Thermal and mechanical constraints on mixing between mafic and silicic magmas[END_REF].

Following [START_REF] Sparks | Thermal and mechanical constraints on mixing between mafic and silicic magmas[END_REF], [START_REF] Frost | Field, chemical, and physical constraints on mafic-felsic interactions in the Lamarck granodiorite, Sierra Nevada, California[END_REF] extensively modelled such a situation although, owing to the lack of pertinent experimental data, they approximated the crystallinity as an error function of T. The rheological patterns derived from recent experimental studies and summarized above allow a critical evaluation of such an approximation. [START_REF] Marsh | On the crystallinity, probability of occurence and rheology of lava and magma[END_REF] has shown that it holds well for most mafic compositions though in some cases a linear fit to the data would also give satisfactory results (e.g., the basalt of Makaopuhi, [START_REF] Marsh | On the crystallinity, probability of occurence and rheology of lava and magma[END_REF]). The experimental data show that this approximation is correct for compositions as felsic as dacites, though in detail small differences can be noted (Fig. 3,4). The modelling of [START_REF] Frost | Field, chemical, and physical constraints on mafic-felsic interactions in the Lamarck granodiorite, Sierra Nevada, California[END_REF] of the interactions between mafic and silicic magmas, the latter having 60 wt.% SiO 2 , probably closely reproduces the essential physics of the process. For more silicic rocks, however, this error function badly fails in reproducing the overall melt fraction trend (Fig. 2), irrespective of the constant used (i.e. the standard deviation parameter) for modifying the general shape. The error function greatly under-estimates the amount of crystallization that takes place from mid-way in the crystallisation interval down to sub-solidus conditions. For instance, after 80 % of cooling, the error function predicts a melt fraction of 10 % vol, whereas the experimental melt fractions vary from 60% for a hornblende-biotite rhyolite crystallizing at 300 MPa to nearly 80 % for a peraluminous leucogranite at 400 MPa (Fig. 2). This implies that for cases of injection at subliquidus conditions of the silicic end-member, the amount of mafic magma needed to prevent its being quenched by the silicic one is less than that calculated by [START_REF] Frost | Field, chemical, and physical constraints on mafic-felsic interactions in the Lamarck granodiorite, Sierra Nevada, California[END_REF] for a granitic composition, because there is less heat taken up to melt pre-existing crystals in the silicic magma.

Schematically, the masses obtained by [START_REF] Frost | Field, chemical, and physical constraints on mafic-felsic interactions in the Lamarck granodiorite, Sierra Nevada, California[END_REF] for interacting mafic and granitic magmas should be scaled down by a factor of 5. Although the interactions between both end-members can be modelled, as done by [START_REF] Sparks | Thermal and mechanical constraints on mixing between mafic and silicic magmas[END_REF] and [START_REF] Frost | Field, chemical, and physical constraints on mafic-felsic interactions in the Lamarck granodiorite, Sierra Nevada, California[END_REF], there is no particular set of initial conditions (especially the crystallinity of both end members at the time of intrusion, and their relative masses) that is of wide applicability. An evaluation of the potential for mixing between mafic and silicic magmas can, however, be done simply by considering the evolution of the viscosity contrast during crystallization. In doing so we take advantage of the fact that rhyolitic magma viscosities remain constant over a large part of the crystallization interval, which implies that, unless intrusion occurs at near-solidus conditions, variations in the viscosity contrast arise predominantly from the viscosity increase of the mafic magma.

Two limiting cases can be envisaged. One is the mixing of a silicic magma with a dry basaltic one whose rheology during crystallization is well known (e.g. [START_REF] Shaw | Rheology of basalt in the melting range[END_REF]. The evolution of magma viscosity ratio during the crystallization of dry basalt and a granite that started to crystallize with 7 wt.% H 2 O shows that, after 40 % of cooling, both magmas display similar viscosities (Fig. 5). Upon further cooling, the basalt viscosity becomes many orders of magnitude higher than that of the rhyolite because the effects of temperature drop and of increase in crystal content cannot be counteracted by that of water. Also, the crystallinity of a basaltic magma increases in a nearly linear fashion as temperature decreases (at least down to ~ 60 vol. % of crystallinity) which implies that the crystal effect on the magma viscosity happens earlier during the cooling of basaltic magmas than of rhyolitic magmas. The net result is that mixing between silicic and dry mafic magmas should be much more effective after substantial crystallization has occurred in the mafic end-member. As shown by previous studies [START_REF] Frost | Field, chemical, and physical constraints on mafic-felsic interactions in the Lamarck granodiorite, Sierra Nevada, California[END_REF], the cross-over in viscosity varies with the initial water content of granitic magmas. Decreasing initial H 2 O content down to 4.5 wt.% delays the attainment of the viscous neutral point to 50 % of cooling (Fig. 5). Assuming that the basalt becomes effectively rigid when its crystallinity is 60 % (i.e. convection ceases, see [START_REF] Marsh | On the crystallinity, probability of occurence and rheology of lava and magma[END_REF]), and given that at this point the granitic magma is largely molten (less than 30 % crystals), both magmas will be still able to convect at the viscosity cross-over. This result confirms previous findings obtained from thermomechanical modeling that, from the strict point of view of the viscosity contrast, the intrusion of a dry basalt into a granitic magma may ultimately produce a mixed/mingled rock of intermediate composition (i.e. andesitic or granodioritic), provided that the mass of the intruded basalt prevents its quenching.

The second case consists of mixing a wet basalt with a silicic magma. This case may find natural applications in active subduction zones where a growing body of evidence indicates that some arc basalts may contain a significant amount of dissolved water, up to 8 wt.% (e.g., Sisson & Grove, 1993a,b;[START_REF] Sisson | H 2 O in basalt and basaltic andesite glass inclusions from four subduction-related volcanoes[END_REF]Sobolev & Chaussidon, 1995;[START_REF] Pichavant | P-T-fO 2 -H 2 O Melt conditions, structure and dynamics of a zoned magma chamber: Mt Pelée (Martinique, Lesser Antilles Arc)[END_REF]. Obviously, the presence of water implies that such basalts will undergo a minor magma viscosity increase during cooling relative to their dry counterpart, and reach the point at which basaltic and granitic viscosities become equal later in their crystallization history. To illustrate this point, we use the experimental results that Sisson and Grove (1993b) obtained on their composition 87S35a, assuming 6 wt% H 2 O dissolved in melt at all temperatures, to calculate melt viscosities and use the modal amount of phases present to calculate magma viscosities. In the case of mixing with a normally wet-silicic magma (4.5 wt.% initial H 2 O), the calculation shows that the viscosities of the silicic magmas are always higher than that of the basalt for crystallinities lower than 60%. If the silicic magma contains an initial H 2 O content of 7 wt.%, both magmas will reach similar viscosities, but only very late in their crystallization interval, when the basaltic magma is very close to having a solid-like rheological behavior. The latter example, however, uses the rheological pattern of a granitic magma crystallizing at 400 MPa, reaching melt H 2 O contents in excess of 9 wt.% at saturation [START_REF] Scaillet | Experimental crystallization of leucogranite magmas[END_REF], and thus comparatively low magma viscosities. Silicic magma chambers of arc-volcanoes lie most often in the 200-300 MPa pressure range (e.g. [START_REF] Rutherford | The May 18, 1980, eruption of Mount St Helens.1. Melt composition and experimental phase equilibria[END_REF][START_REF] Rutherford | Pre-eruption pressure-temperature conditions and volatiles in the 1991 Mount Pinatubo magma[END_REF][START_REF] Luhr | Experimental phase relations of water-and sulfur-saturated arc magmas and the 1982 eruptions of El Chichòn volcano[END_REF][START_REF] Martel | Magma storage and control of eruption regime in silicic volcanoes: experimental evidence from Mt Pelée[END_REF][START_REF] Scaillet | The June 15, 1991 eruption of Mount Pinatubo. I. Phase equilibria and pre-eruption P-T-fO 2 -fH 2 O conditions of the dacite magma[END_REF], in which the water solubility is 7-8 wt. %. In addition, most of the erupted rocks record pre-eruptive melt H 2 O contents in the range of 3-5 wt.% (see Scaillet et al., 1998a). All these facts suggest that the conditions appropriate for thorough mafic-silicic magma mixing in the case of H 2 O-rich basalts are probably never met under shallow crustal levels. The possibility of a causal link between basalt intrusion and explosive eruption of silicic magmas [START_REF] Sparks | Magma mixing: A mechanism for triggering acid explosive eruption[END_REF] also points to a rather limited time of coexistence between the magmas, on the order of only a few days as petrological evidence suggests (e.g., [START_REF] Nixon | Petrology of the younger andesites and dacites of Iztaccihuatl volcano, Mexico: II. Chemical composition of basaltic magma influx[END_REF][START_REF] Pallister | A basalt trigger for the 1991 eruptions of Pinatubo volcano?[END_REF]. Thus, mixing between wet basalts and rhyolitic magmas under shallow crustal conditions has little chance of producing a homogeneous intermediate magma, as testified by the widespread occurrence of poorly mixed maficsilicic end-members in arc magmas (e.g., [START_REF] Eichelberger | Origin of andesite and dacite: evidence of mixing at Glass Mountain in California and at other circum-Pacific volcanoes[END_REF][START_REF] Nixon | Petrology of the younger andesites and dacites of Iztaccihuatl volcano, Mexico: II. Chemical composition of basaltic magma influx[END_REF][START_REF] Pallister | A basalt trigger for the 1991 eruptions of Pinatubo volcano?[END_REF].

Implications for the use of mafic enclaves in strain analyses

The intrusion of small dikes of more mafic composition into a large reservoir of silicic magma may eventually, upon disruption, give rise to isolated mafic enclaves (e.g., [START_REF] Bacon | Magmatic inclusions in silicic and intermediate volcanic rocks[END_REF]. Such enclaves are found in many granitoids and are commonly used either for the quantification of finite strain (e.g., [START_REF] John | Emplacement-related deformation of granitoid magmas, southern Adamello Massif, Italy[END_REF] or for evaluating strain rates [START_REF] Williams | Microgranitic enclave shapes and magmatic strain histories: Constraints from drop deformation theory[END_REF]. [START_REF] Williams | Microgranitic enclave shapes and magmatic strain histories: Constraints from drop deformation theory[END_REF] used enclave shapes to calculate magmatic flow velocities, by considering mainly the case where the enclave viscosity is lower than that of the host magma, a situation that implies overheated silicic magma (temperatures significantly higher than 1000°C). This arises when the mass of injected mafic magma is at least comparable to that of the intruded/underplated one [START_REF] Sparks | Thermal and mechanical constraints on mixing between mafic and silicic magmas[END_REF][START_REF] Frost | Field, chemical, and physical constraints on mafic-felsic interactions in the Lamarck granodiorite, Sierra Nevada, California[END_REF] such as in the case of basalt underplating of the lower continental crust [START_REF] Huppert | The generation of granitic magmas by intrusion of basalt into continental crust[END_REF]). Here we concentrate on the case where the incorporation of small blobs of mafic magmas does not affect the thermal regime of the host magma, extending the analysis done by [START_REF] Williams | Microgranitic enclave shapes and magmatic strain histories: Constraints from drop deformation theory[END_REF] to the low temperature range of the crystallisation interval of granitoid magmas. In such a case, beside transient effects related to the local superheating of the silicic magma adjacent to the mafic magma, the equilibration temperature of the resultant mixture will not differ greatly from the pre-intrusion temperature of the resident magma. The viscosity contrasts obtained under such circumstances have been calculated for the intrusion of wet basalts into andesitic, dacitic, and rhyolitic magmas, dacitic into granitic magmas, and andesitic into dacitic magmas, using the experimental results presented above to constrain the rheological patterns of each magma. The case for dry basalt is not shown since its solidus temperature, about 1000°C, is above the liquidus temperatures of most silicic to intermediate magmas which implies that small dikes of dry basalt injected will very rapidly behave as rigid objects with respect to any applied stress.

The case of mixing between wet basalt and dacite is shown in figure 6 without normalizing both temperature and viscosity to help understand the basic features common to all cases dealt with here. Two viscosity cross-overs can be observed: the first is attained when the basalt crystallinity reaches 60 %, at around 915°C, which induces a dramatic rise in the enclave viscosity.

The second cross-over is reached when the host silicic intermediate is 60 % crystallised, which occurs at about 750°C, at which point the viscosity rapidly rises and becomes equal to that of the basaltic enclave since both magmas basically behave as solids. Obviously, the temperature range between these two viscosity cross-overs is controlled by the liquidus-solidus temperature and the melt fraction-temperature trend specific to each magma, both factors depending on the phase relationships of the magma. The variation of these factors with the chemistry of the mixed magmas is illustrated in the following section.

For basaltic enclaves in either rhyolitic or dacitic magmas (Fig. 7), the rheological patterns are all similar, being characterized by small viscosity contrasts at near-liquidus conditions of the granitoids, followed at intermediate temperatures by a large domain, in some cases more than 50 % of the crystallization interval, during which mafic enclaves are many orders of magnitude more viscous than their surrounding matrix (Fig. 7). As the melt fraction of the host magma decreases, however, the viscosity contrast drops again and becomes nil as soon as the silicic magma crosses the transition from suspension-like to rigid-like behavior, thought to occur at around 30-40 % of melt by volume. In rhyolitic magmas such a transition will occur close to solidus conditions, in keeping with their near-eutectic behavior (around 680°C, Fig. 2), whereas in dacitic bulk compositions, this rheological change will occur at early stages of the crystallization interval (around 750°C at 220 MPa, Fig. 3). These calculations indicate that enclaves of mafic compositions found in granitoids will have a viscosity several orders of magnitude higher than their host magmas during an extended part of the crystallization interval.

For injection of dacitic dikes into crystallizing granitic plutons the same general pattern is observed (Fig. 8), with the viscosity of the two magmas differing by less than 1 order of magnitude during the first 30-40% of the crystallization interval, while the intermediate domain of high viscosity for the mafic enclaves is significantly reduced, as might be expected from the greater similarity in crystallization temperatures and compositions between resident and injected magmas. The case of andesite inclusions in a dacite magma is shown on figure 8, where it is apparent that the rheological behavior of enclave and host is mainly controlled by the H 2 O content of host magma. Andesite enclaves in an H 2 O-saturated dacite will display a large viscosity contrast during a longer part of the crystallisation interval of their host than when included in a H 2 O-undersaturated dacite magma. This is a direct consequence of the lower crystallisation temperature of the H 2 Osaturated dacite relative to the andesite. Although the H 2 O content of the inclusion exerts a minor influence of the viscosity contrast, figure 8 shows that, at near liquidus conditions, the viscosity contrast is smaller when the melt water content of the andesite is the same as in the dacite. Finally, injection of basalt in andesite (Fig. 7 show that the potential for mixing between a mafic and a silicic magma is optimum at near liquidus conditions of the intruded silicic end member and lasts for about 50 % of the crystallisation interval in the most favourable case.

Repeated injection of small amounts of mafic magma in a near-liquidus felsic magma chamber may thus progressively drive the composition of the ambient magma toward a more mafic bulk composition. In contrast injection of mafic dykes at later stages of the evolution of the felsic magma chamber are less likely to induce severe chemical changes or homogeneous mixtures.

Decreasing the compositional gap between the two end-members produces, on average, a substantial diminution of the domain where the viscosity contrast is very large as seen for instance for the dacite into rhyolite or basalt into andesite cases (Fig. 7,8). Differences in the rheological contrast are not simply correlated, however, to a specific compositional variable such as the bulk SiO 2 , as illustrated by the andesite/dacite situation compared to the basalt/andesite one. In the first case, the end member magmas differ only by 5 wt% in SiO 2 while there is a 10 wt% difference in SiO 2 between the silicic andesite and basalt. Yet, the viscosity difference during cooling is greater between andesite and dacite than between andesite and basalt (Fig. 7,8). This difference reflects variable liquidus-solidus temperature ranges as well as variations of cristallinity with temperature of the magmas involved in the process. More specifically, the phase equilibrium data show that, at 200 MPa, whilst the hydrous liquidi of both andesitic and basaltic magmas differ only by 10°C, the difference between liquidus temperatures of dacite and andesite reaches 55°C. For magmas having their crystallinity linearly correlated with temperature (i.e. for any magma more mafic than a dacite), the amount of crystallisation for each decrement of temperature will vary with the crystallisation interval, being larger when the crystallisation interval is smaller. Thus, a dacite with an initial melt H 2 O content of 4 wt% will have a crystallisation interval larger by 120°C than when at H 2 O-saturation. The results show clearly that the melt water content of the host exerts a critical role on the viscosity contrast as illustrated here for the andesite into dacite situation (Fig. 8). Altogether, these examples illustrate that the rheological behavior is hardly predictable from the consideration of the anhydrous composition only. Water appears to be as important as the bulk anhydrous composition of the magmas involved in the mixing/mingling process. This emphasizes the importance of a correct determination of intensive variables for a realistic evaluation of the rheological properties of natural magmas.

Considering now the possibility of using mafic enclaves as strain gauges, the previous analysis shows that mafic enclaves may experience rather complicated magmatic strain histories, reflecting changes from less to more viscous behavior during the crystallization of their host magma.

Overall, however, and in particular for rhyolitic magmas hosting enclaves of basaltic composition, enclaves are likely to behave more as rigid objects, relative to their host medium, than as passive markers, especially in the interval of viscous behavior of the granitoid. Experimental evidence suggests that once the viscosity ratio exceeds 4 (that is 0.6 in log units), enclaves are virtually unaffected by the deformation (see for more details [START_REF] Williams | Microgranitic enclave shapes and magmatic strain histories: Constraints from drop deformation theory[END_REF]) and thus cannot record it through changes of shape. This indicates that the use of mafic enclaves as finite strain markers will give only a minimum estimate of the total amount of strain that affected a plutonic intrusion. [START_REF] Williams | Microgranitic enclave shapes and magmatic strain histories: Constraints from drop deformation theory[END_REF] concluded that mafic enclaves will mainly record magmatic strain affecting the host magma at near-or superliquidus conditions, but inspection of figures 7 and 8 shows clearly that they can also record strain during the waning stages of granitoid crystallization. In fact, the state of passive marker for mafic enclaves is attained and maintained for sustained periods only very late in the crystallization history of granitoid magmas, basically when the latter reach solid behavior (Fig. 7,8). Viscosity cross-overs occurring at early stages are necessarily transient and are unlikely to last for long periods, especially in rhyolitic systems where the latent heat release upon crystallization, and thus the time of cooling, is prominent at near-solidus conditions. Indeed, experimental melt fraction trends obtained in rhyolitic systems (Fig. 2) show that the major part of the crystallization occurs at sub-solidus conditions, which implies that fabrics preserved in granitic plutons have little chance of providing clues about emplacement mechanisms, or about any strain event occuring at early stages [START_REF] Scaillet | Rheological properties of granitic magmas in their crystallization range[END_REF]. A similar conclusion can be applied to strain trajectories derived from the study of mafic enclaves when they parallel magmatic foliation patterns defined by mineral orientations. They cannot be interpreted solely as remnants of magmatic directions of flow nor in terms of effects related to emplacement mechanisms, since their shape must reflect to a large extent deformation episodes that took place when their host matrix had a rigid behavior.

Conclusions

We have attempted to place some constraints on the rheological properties of natural magmas using recent phase equilibrium data. The strength of the present approach is that full consideration is given to the P-T-melt H 2 O content conditions of evolution of the magma. These constraints are by no means definitive, however, for several reasons. First, the chemistries of magmas have been considered constant during interaction. This is clearly an over-simplification, since numerous field evidence show that chemical interaction takes place between intruded and host magmas (e.g., [START_REF] Debon | Comparative major element chemistry in various "microgranular enclave-plutonic host" pairs[END_REF][START_REF] Barbarin | Genesis and evolution of mafic microgranular enclaves through various types of interaction between coexisting felsic and mafic magmas[END_REF], decreasing the compositional gap, and thus the rheological contrast, between coeval mafic and felsic magmas. Second, in all the cases handled here, the melt viscosity has been calculated using the Shaw model [START_REF] Shaw | Viscosities of magmatic silicate liquid: an empirical method of prediction[END_REF]. As stated in the introduction, this model has been shown to be in error under some conditions for silicic magmas and it cannot be excluded that its application to mafic compositions may also lead to erroneous results. Although considerable progress has been made in the field of viscous behavior of silicic magmas, there have been few recent attempts to investigate the viscosity of hydrous mafic melts, and no direct measurements at magmatic temperatures. We consider that such data are critical for a rigorous evaluation of mixing processes such as those illustrated above and suggest that future laboratory studies take this need into account. Second, the magma viscosity is here calculated using the melt fraction obtained at equilibrium, yet a deforming magma is by definition out of equilibrium.

Available data concerning the effect of deformation upon crystallisation are restricted to dry basaltic composition. These data demonstrate that deformation induces severe changes in the crystallisation sequence of mafic magmas [START_REF] Kouchi | Effects of stirring on crystallisation kinetics of basalt: texture and element partitioning[END_REF] and it cannot be precluded that the phase relations of felsic magmas are sensitive to deformation as well. The effect of deformation on crystallising silicate liquids is clearly one aspect that should also be rigorously explored in the near future if we wish to correctly characterise the rheological properties of natural magmas seen in a dynamic perspective. [START_REF] Marsh | On the crystallinity, probability of occurence and rheology of lava and magma[END_REF]. Note the constancy or even the decrease of melt viscosity as cooling proceeds. Melt viscosities are calculated using the empirical Figure 3. Experimental melt fractions, melt and magma viscosity trends for a dacitic magma versus normalized temperature, calculated using the experimental results of [START_REF] Scaillet | The June 15, 1991 eruption of Mount Pinatubo. I. Phase equilibria and pre-eruption P-T-fO 2 -fH 2 O conditions of the dacite magma[END_REF]. Numbers on the trends correspond to the initial melt H 2 O content . Also shown is the melt fraction trend obtained with the error function using a standard deviation parameter of 1.5 [START_REF] Marsh | On the crystallinity, probability of occurence and rheology of lava and magma[END_REF]. Note the continuous increase of magma viscosities as cooling proceeds. Liquidus and solidus temperatures are 1040°C-680°C (4 wt.% H 2 O initial) and 920°C-680°C (6.9 wt.% H 2 O initial). Melt viscosities are calculated using the model of [START_REF] Shaw | Viscosities of magmatic silicate liquid: an empirical method of prediction[END_REF]. Magma viscosities are calculated using the Einstein-Roscoe equation (e.g., [START_REF] Lejeune | Rheology of crystal-bearing silicate melts: An experimental study at high viscosities[END_REF].

Liquidus and solidus temperatures are 1040°C-680°C (4 wt.% H 2 O initial) and 920°C-680°C (6.9 wt.% H 2 O initial). an andesitic magma versus normalized temperature, calculated using the experimental results of [START_REF] Martel | Effects of fO 2 and H 2 O on andesite phase relations between 2 and 4 kbar[END_REF]. Numbers on the trends correspond to the initial melt H 2 O content. Also shown is the melt fraction trend obtained with the error function using a standard deviation parameter of 1.5 [START_REF] Marsh | On the crystallinity, probability of occurence and rheology of lava and magma[END_REF]. Note the continuous increase of magma viscosities as cooling proceeds. Liquidus and solidus temperatures are 1000°C-680°C (4 wt.% H 2 O initial) and 975°C-680°C (6.9 wt.% H 2 O initial). Melt viscosities are calculated using the model of [START_REF] Shaw | Viscosities of magmatic silicate liquid: an empirical method of prediction[END_REF]. Magma viscosities are calculated using the Einstein-Roscoe equation (e.g., [START_REF] Lejeune | Rheology of crystal-bearing silicate melts: An experimental study at high viscosities[END_REF].

Liquidus and solidus temperatures are 1000°C-680°C (4 wt.% H 2 O initial) and 975°C-680°C (6.9 wt.% H 2 O initial). [START_REF] Shaw | Rheology of basalt in the melting range[END_REF], whereas the viscosity of wet basalt has been calculated from the experimental results of Sisson and Grove (1993) using the model of [START_REF] Shaw | Viscosities of magmatic silicate liquid: an empirical method of prediction[END_REF]. Rhyolitic magma viscosities are from figure 2. Numbers reported on the lines correspond to the initial melt H 2 O content of the rhyolitic magma. (1999) and for the basalt from Sisson and Grove (1993). In both cases the solidus is fixed at 680°C. Melt viscosities are calculated using the model of [START_REF] Shaw | Viscosities of magmatic silicate liquid: an empirical method of prediction[END_REF]. Magma viscosities are calculated using the Einstein-Roscoe equation (e.g., [START_REF] Lejeune | Rheology of crystal-bearing silicate melts: An experimental study at high viscosities[END_REF].

  ), is characterized by a substantial shrinking of the domain with large viscosity contrast. No fundamental difference in the rheological pattern compared to the previous examples is seen when the H 2 O content of the andesite varies, although the peak value in the viscosity contrast is the smallest of all cases considered here. Relative to the dacite/rhyolite and andesite/dacite mixing cases, however, the intermediate domain with a high viscosity contrast is displaced toward nearliquidus conditions. The various examples displayed in Figs.7-8 invariably

Figure 1 .Figure 2 .

 12 Figure 1. Temperature against the melt H 2 O content, melt and magma viscosities of silicic to intermediate magmas at pre-eruption conditions or during emplacement in the upper crust (modified fromScaillet et al., 1998a)

Figure 4 .

 4 Figure 4. Experimental melt fractions, melt and magma viscosity trends for

Figure 5 .

 5 Figure 5. Comparison of the evolution of viscosities during the crystallization of basaltic and rhyolitic magmas. The viscosity of the dry basaltic magma is from Shaw (1969), whereas the viscosity of wet basalt has

Figure 6 .

 6 Figure 6. Comparison of the viscosity evolution of dacitic and wet basaltic magmas during cooling. Data for the dacite are from Scaillet and Evans

Figure 7 .

 7 Figure 7. Evolution of the viscosity ratio of a thermally equilibrated basaltic enclave and a host rhyolitic, dacitic or andesitic magma. The temperature of the mafic enclave is that of the silicic host. Numbers reported on the lines correspond to the initial melt H 2 O content of the host magma.The viscosity in the solid behavior field (that is for melt fractions lower than 30 vol. %) is taken to 10 20 Pa.s. Liquidus and solidus temperatures of the silicic magmas

Acknowledgements

This work was supported by the IDYL and PNRN CNRS research programs and by a TMR-CEE project. The informal review of Ray Macdonald greatly improved the quality of the manuscript. Official reviews done by Don Baker and an anonymous reviewer and the careful editorial handling of Wes

Stephens are gratefully acknowledged.

are as in figures 2-4. Numbers reported on the lines correspond to the initial melt H 2 O content of the host magma. Open symbols are for 4 wt% H 2 O in the andesite enclave, filled ones for 6.9 wt% H 2 O. The viscosity in the solid behavior field (that is for melt fractions lower than 30 vol. %) is taken to 10 20 Pa.s. Liquidus and solidus temperatures of the felsic magmas are as in figures 2-4.


V Boyer
email: vboyer@laas.fr

D El Baz
email: elbaz@laas.fr

M Elkihel
email: elkihel@laas.fr

Solving knapsack problems on GPU

Keywords: Knapsack problems, Dense dynamic programming, Parallel computing, GPU computing

In this article, we propose a parallel implementation of the dynamic programming method for the knapsack problem on NVIDIA GPU. A GTX 260 (192 cores, 1.4GHz) was used for computational tests and processing times obtained with the parallel code are compared to the sequential one on a CPU with an Intel Xeon 3.0GHz. The results show a speedup up to 26 and permit one to solve large size problems within a reasonable processing time. Furthermore, in order to limit the communication between the CPU and the GPU, a compression technique is presented which decreases significantly the memory occupancy.

Introduction

Since a few years, graphics card manufacturers have developed tools to use their products for high performance computing. Graphics Processing Units, or GPUs, are high-performance many-core processors. NVIDIA GPUs are SIMT (Single Instruction, Multiple Thread) architectures which is akin to SIMD (Single Instruction, Multiple Data) architecture [START_REF]Cuda 2.0 programming guide[END_REF]. Severals parallel dynamic programming method have been proposed (see, for example, [START_REF] Baz | Load balancing methods and parallel dynamic programming algorithm using dominance technique applied to the 0-1 knapsack problem[END_REF], [START_REF] Lou | A parallel two-list algorithm for the knapsack problem[END_REF] and [START_REF] Gerash | A survey of parallel algorithms for onedimensional integer knapsack problems[END_REF]). Implementations on SIMD machine were performed on a 4K processor ICL DAP [START_REF] Kindervater | An introduction to parallelism in combinatorial optimization[END_REF] and 16K Connection Machine CM-2 systems (see [START_REF] Lin | Processor-efficient algorithms for the knapsack problem[END_REF], [START_REF] Ulm | Dynamic programming implementations on SIMD machines -0/1 knapsack problem[END_REF]) and a 4K MasPar MP-1 machine [START_REF] Ulm | Dynamic programming implementations on SIMD machines -0/1 knapsack problem[END_REF]. To the best of our knowl-edge, no parallel implementation of dynamic programming for combinatorial optimization problems has been done on a GPU. Using GPU architectures for solving combinatorial optimization problems is a great challenge in order to reduce the computing time needed to solve these NP-hard problems. In this paper, we propose a parallel implementation of the dynamic programming algorithm on a NVIDIA GPU for exactly solving the KP. Furthermore, a data compression is presented. This compression permits one to limit significantly, the communication between the CPU and the GPU and the memory needed to build the solution vector. The paper is structured as follows. Section 2 deals with the knapsack problem and its solution via dynamic programming. Section 3 focuses on GPU computing and the parallel dynamic programming method. In section 4 we propose a compression method which permits one to reduce significantly the memory occupancy and communication between CPU and GPU. Section 5 deals with computational experiences. Conclusion and future work are presented in section 6.

The knapsack problem

The knapsack problem, KP, is a NP-hard combinatorial optimization problem. It is one of the most studied discrete programming problem as it is among the simplest prototypes of integer linear programming problems and it arises in several sub-problems of many more complex problems (see, for example, [START_REF] Boyer | Heuristics for the 0-1 multidimensional knapsack problem[END_REF], [START_REF] Boyer | A dynamic programming method with list for the knapsack sharing problem[END_REF], [START_REF] Crowder | Solving large-scale zero-one linear programming problems[END_REF], [START_REF] Dietrich | More coefficient reduction for knapsacklike constraints in 0-1 programs with variable upper bounds[END_REF], [START_REF] Dietrich | New procedures for preprocessing 0-1 models with knapsack-like constraints and conjunctive and/or disjunctive variable upper bounds[END_REF], [START_REF] Kellerer | Knapsack Problems[END_REF] and [START_REF] Martello | Knapsack Problems -Algorithms and Computer Implementations[END_REF]).

Problem formulation

Given a set of n items i, with profit p i ∈ N * + and weight w i ∈ N * + , and a knapsack with the capacity C ∈ N * + , KP can be defined as the following linear integer programming problem:

(KP )              max n i=1 p i .x i , s.t. n i=1 w i .x i ≤ C, x i ∈ {0, 1}, i ∈ {1, ..., n}. (2.1)
To avoid any trivial solution, we assume that:

     ∀i ∈ {1, ..., n}, w i ≤ C, n i=1 w i > C.

Dynamic programming

Bellman's dynamic programming [START_REF] Bellman | Dynamic programming[END_REF], presented in 1957, was the first exact algorithm to solve KP. It consists in computing at each step k ∈ {1, ..., n}, the values of f k (ĉ), ĉ ∈ {0, ..., C}, using the classical recursion:

f k (ĉ) = f k-1 (ĉ), f or ĉ = 0, ..., w k -1, max {f k-1 (ĉ), f k-1 (ĉ -w k ) + p k } , f or ĉ = w k , ..., C. (2.2)
with, f 0 (ĉ) = 0, ĉ ∈ {0, ..., C}.

The algorithm, presented in this section, is based on the the Bellman's recursion (2.2). A state corresponds to a feasible solution associated with the f k (ĉ) value. Toth [START_REF] Toth | Dynamic programming algorithm for the zero-one knapsack problem[END_REF] has proposed an efficient recursive procedure in order to compute the states of a stage and used the following rule to eliminate states:

Proposition 1 [START_REF] Toth | Dynamic programming algorithm for the zero-one knapsack problem[END_REF] If a state defined at k -th stage with total weight ĉ satisfies:

ĉ < C - n i=k+1 w i ,
then the state will never lead to an optimal solution and can be eliminated.

The dynamic programming method is described in algorithm 1. The matrice M stores all the decisions and is used to build a solution vector corresponding to the optimal value by doing a backtracking. The time and space complexities are O(n.C).

Algorithm 1 (Dynamic programming) for ĉ ∈ {0, ..., C}, f (ĉ) := 0, for i ∈ {1, ..., n} and ĉ ∈ {1, ..., C}, M i,ĉ = 0,

sumW := n i=1 w i , for k from 1 to n do sumW := sumW -w k , c = max{C -sumW, w k }, for ĉ from C to c do if f (ĉ) < f (ĉ -w k ) + p k then f (ĉ) := f (ĉ -w k ) + p k , M i,ĉ := 1, end if, end for, end for. return f (C) (the optimal bound of the KP)
The high memory requirement are frequently cited as the main drawback of dynamic programming. However, this method has a pseudo-polynomial time complexity and is insensitive to the kind of instances, i.e. with correlated datas or not. In order to reduce the memory occupancy, the entries of the matrix M have been stored in integers of 32 bits. This permits one to divide by 32 the number of lines of the matrix and the memory needed. However, the memory occupancy is still important and an efficient compression method will be presented in section 4.

GPU computing

GPUs are highly parallel, multithreaded, many-core architectures. NVIDIA introduced, in 2006, CUDA, a software development kit that enables users to solve many complex computational problems on their GPU cards. This tool has been used in order to implement our parallel dynamic programming code.

NVIDIA GPU architecture

On CUDA compatible NVIDIA cards, the threads are separated in blocks and these blocks are distributed on the multiprocessors. A multiprocessor processes one block at a time. When the threads of a block terminate, a new block is launched on the idle multiprocessor. The multiprocessor executes threads in groups of 32 parallel threads called warps. Threads composing a warp start together at the same program address, they are nevertheless free to branch and execute independently.

Another important aspect in GPU computing is the memory hierarchy. Indeed, threads have access to data from multiple memory spaces. Each thread has a private local memory. Each thread block has a shared memory visible to all threads of the block which has the same lifetime as the block. Finally, all threads have access to a global memory. Furthermore, there is two other read-only memory spaces accessible by all threads which are cache memories:

• the constant memory, for constant data used by the process,

• the texture memory space, designed for graphics applications.

In order to have a maximum bandwidth for the global memory, we have to insure a coalesced memory accesses. Indeed, the global memory access by all threads is done in one or two transactions if:

• threads access:

either 32-bit words, resulting in one 64-byte memory transaction, or 64-bit words, resulting in one 128-byte memory transaction, or 128-bit words, resulting in two 128-byte memory transactions;

• all 16 words lie in the same segment of size equal to the memory transaction size (or twice the memory transaction size when accessing 128-bit words);

• threads access the words in sequence (the kth thread in the half-warp accesses the kth word).

Otherwise, a separate memory transaction is issued for each thread. For further details on the NVIDIA cards architecture and how to optimize the code, see [START_REF]Cuda 2.0 programming guide[END_REF].

Parallel dynamic programming

Parallel implementation of the dynamic programming method is optimized for GPU NVIDIA architectures. The activity that consumes processing time is the loop that processes the values of f (ĉ), ĉ ∈ {0, ..., C}. This step has been parallelized on the GPU: one thread processes one value of f . Many efforts have been made in order to limit the communication between the CPU and the GPU and ensure coalesced memory access in order to significantly reduce the processing time. The procedures implemented on the CPU and the GPU, respectively, are described in algorithms 2 and 3, respectively.

• input f are the values of f processed at the previous step,

• output f are the output values of f ,

• output m are the output values of the decisions stored as integers of 32 bits and

• c min is the minimum value of ĉ.

This function creates C -c min + 1 threads for the GPU and groups them into blocks of 512 threads (the maximum size of a block of one dimension), i.e. (C -c min + 1)/512 blocks. All threads carry out on the GPU the procedure described in algorithm 3.

Algorithm 3 (Thread processing on GPU) blocks id: the ID of the belonging block, thread id: the ID of the thread within the belonging block, k: the step number of the dynamic programming (k ∈ {1, ..., n}), i := (k + 31)%32: the rest of the division of k + 31 by 32,

ĉ:=blocks id * 512 + thread id, if ĉ < c min or ĉ > C then STOP end if, if input f (ĉ) < input f (ĉ -w k ) + p k then output f (ĉ) := output f (ĉ -w k ) + p k , output m ĉ := output m ĉ + 2 i , else output f (ĉ) := output f (ĉ), end if
In the algorithm 3, threads have to access the values of input f (ĉ -w k ), this results in un-coalesced memory accesses as described in section 3.1. In order to reduce the memory latency, the texture memory is used to access the data stored in input f . We used the texture memory since this type of memory can be allocated dynamicaly contrarily to the constant memory. output f and output m are stored in the global memory.

Reducing memory occupancy

The analysis of the values stored in the matrix M h shows that the right columns are often filled with 1 and that the left columns are filled with 0. As these bits values are grouped in integers of 32 bits, in practice it corresponds to the value 2 32 -1 for the right columns (and 0 for the left columns). Thus, the communication between the CPU and the GPU occurs every 32 iterations in order to retrieve all the decisions stored in m d into the matrix M h (see algorithm 2). This step is time consuming and we have tried to further reduce the amount of data transfered to the CPU. A simple way to reduce the vector m d is then to compress it as follows:

f or ĉ ∈ {0, ..., rc -lc}, m d c ĉ = m d ĉ+lc with lc = min{ĉ ∈ {1, ..., C} | m d ĉ-1 = 0 and m d ĉ = 0}, rc = max{ĉ ∈ {1, ..., C} | m d ĉ-1 = 2 32 -1 and m d ĉ = 2 32 -1}.
Thus, we know that:

• if ĉ < lc, then m d ĉ = 0 and • if ĉ ≥ rc, then m d ĉ = 2 32 -1.
Then we have to retrieve only the values of m d ĉ for ĉ ∈ {lc, ..., rc -1} and we process lc and rc directly on the GPU via the algorithm 4.

Algorithm 4 (Thread compression on GPU) blocks id: the ID of the belonging block, thread id: the ID of the thread within the belonging block, m d: the input vector, lc: shared variable initiate with the value C, rc: shared variable initiate with the value 0,

ĉ:=blocks id * 512 + thread id, if ĉ ≤ 0 or ĉ > C then STOP end if, if m d ĉ-1 = 0 and m d ĉ = 0 then lc := min{ĉ, lc}, end if, if m d ĉ-1 = 2 32 -1 and m d ĉ = 2 32 -1 then rc := max{ĉ, rc}, end if.
This compression method decreases the amount of data transfered from the GPU to the CPU and permits one also to decrease significantly the memory occupancy needed to store all the decisions made throughout the dynamic programming recursion. Computational experiences shows that the efficiency of the compression depends on the sorting of the variables of the KP and, in average, the best results have been obtained with the following sorting:

p 1 w 1 ≥ p 2 w 2 ≥ ... ≥ p n w n .

Computational experiences

Computational tests have been carried for randomly generated correlated problems, i.e. problems such that:

• w i , i ∈ {1, ..., n}, is randomly draw in [1, 1000],
• p i = w i + 50, i ∈ {1, ..., n},

• C = 1 2 . n i=1 w i .
For each instance, the average results displayed have been obtained with 10 problems.

A NVIDIA GTX 260 graphic card (192 cores, 1.4GHz) has been used and the parallel computational time is compared with the sequential one obtained on a CPU with an Intel Xeon 3.0GHz. Results on the memory occupancy are also presented.

Memory occupancy

In this section, we display the results obtained with the compression method presented in section 4. Table 1 shows the factor of compression computed as follow:

comp f actor = size of M c + 2. n/32 size of M ,
where M is the matrix of decision and M c the corresponding compressed matrix. 2. n/32 corresponds to the values of lc and rc needed for each lines. which leads to a very small memory occupancy as compared with the original dynamic programming algorithm. Furthermore, the factor of compression decreases with the size of the knapsack. This method of compression reduces significantly the memory occupancy of the dynamic programming algorithm and is robust when the number of variables increases. This permits one to solve larger problems that could not be solved otherwise, like problems with 100000 variables. Time spent for the compression step is presented in the next subsection, in order to be compared to the overall processing time.

Processing time

Table 2 presents the average processing time to solve KP obtained with the sequential and parallel algorithms. It also shows the corresponding average time spent during the compression step. Table 3 provides the resulting speedup.

We can see that the processing time cost of the compression step is relatively small as compared with the overall one. These results include the compression step and the transfer of data to the CPU. Thus, this approach is very efficient both in terms of memory occupancy and processing time.

The comparison of the parallel implementation with the sequential one shows that the resulting speedup increases with the size of the problem and meets a level around 26. Our parallel implementation of the dynamic programming reduces significantly processing time and shows that solving hard knapsack problems is possible on GPU. The parallel implementation of the dynamic programming algorithm on GPU combined with our compression method permits one to solve large size problems within a small processing time and a small memory occupancy.

Conclusion

In this article we have proposed a parallel implementation of the dynamic programing algorithm for the knapsack problem on NVIDIA GPU with CUDA. This algorithm has been combined with data compression techniques. Computational experiences have shown that large size problems can be solved within small processing time and memory occupancy. The proposed approaches, i.e. implementation on GPU and data compression, seems to be robust as the results are not deteriorated when the number of variables increases. The observed speedup appears to be stable (around 26) for instances with more than 40000 variables. The reduction of the size on the matrix increases with the number of variables, resulting in a more efficient compression and the overhead does not exceed 3% of the overall one.

The proposed parallel algorithm to solve knapsack problems on GPU shows the relevance of using this type of architecture for combinatorial optimization. Further computational experiences are foreseen, in particular with a NVIDIA Tesla cards and hybrid supercomputers which are dedicated to hight performance computing.

  n

Table 1

 1

	shows that in the worst case the size of the compressed data (size of M c+
	2. n/32 ) corresponds to only 0.3% of the size of the initial the matrix M ,

Table 1 :

 1 Factor of data compression.

Table 2 :

 2 Processing time (s.).

		t. //		t. seq.	t. comp. //	t. comp. seq.
	10000	3.06		58.95	0.11	1.08
	20000	11.97		226.66	0.31	4.16
	30000	26.57		536.14	0.71	9.27
	40000	47.43		1225.52	1.23	18.66
	50000	73.55		1912.43	1.85	25.66
	60000	105.93	2752.81	3.25	38.14
	70000	143.98	3739.74	3.61	50.15
	80000	183.15	4771.55	4.69	64.09
	90000	238.57	6184.28	5.95	82.56
	100000	289.21	>7200	7.16	-
	t. //: average parallel time
	t. seq: average sequential time
	t. comp. //: average parallel time for compression
	t. comp. seq.: average sequential time for compression
	n		speedup	n	speedup
	10000	18.90		60000	25.98
	20000	19.26		70000	25.97
	30000	20.17		80000	26.05
	40000	25.83		90000	25.92
	50000	26.00		100000	-

Table 3 :

 3 Speedup ( t. seq.

t. // ).

In algorithm 2, the launching of the threads on GPU is done via the following function:

Compute f and m on device(input f,output f,output m,c min)

where: