Emplacing a Cooling-Limited Rhyolite Lava Flow: Similarities with Basaltic Lava Flows

Accurate forecasts of lava flow length rely on estimates of eruption and magma properties and, potentially more challengingly, an understanding of the relative influence of characteristics such as the apparent viscosity, the yield strength of the flow core, or the strength of the surface crust. Consequently, even the most straightforward models of lava advance involve sufficient parameters that constraints can be relatively easily fitted within the uncertainties involved, at the expense of gaining insight. Here, for the first time, we incorporate morphological observations from during and after flow field evolution to improve model constraints and reduce uncertainties. After demonstrating the approach on a basaltic lava flow (Mt. Etna, 2001), we apply it to the 2011-12 Cordon Caulle rhyolite flow, where unprecedented observations and syn-emplacement satellite imagery of an advancing silica-rich lava flow have indicated an important crustal influence on flow emplacement. Our results show that an initial phase of viscosity-controlled advance at Cordon Caulle was followed by later crustal control, accompanied by formation of flow surface folds and large-scale crustal fractures. Where the lava was unconstrained by topography, the cooled crust ultimately halted advance of the main flow and led to the formation of breakouts from the flow front and margins, influencing the footprint of the lava, its advance rate, and the duration of flow advance. Highly similar behaviour occurred in the 2001 Etna basaltic lava flow. The processes controlling the advance of crystal-poor rhyolite and basaltic lava flow therefore appear similar, indicating common controlling mechanisms that transcend profound rheological and compositional differences.


INTRODUCTION
Lava flows present a significant hazard in the immediate vicinity of volcanic complexes and can, in some cases, also present a hazard to more distal communities. However, despite the sophisticated numerical models that have been developed for basaltic lava flows (Crisci et al., 1986;Hidaka et al., 2005;Vicari et al., 2007;Hérault et al., 2011), a full understanding of the factors controlling the rate and extent of lava flow advance, which is essential for adequate hazard forecasting over a broad range of lava geochemistries, remains elusive due to the complexity of lava flow rheology, internal architecture, and interactions with topography. The frequency of basaltic eruptions has provided many examples for modeling low-silica lavas (Crisci et al., 1986;Hidaka et al., 2005;Vicari et al., 2007;Hérault et al., 2011), but equivalent studies of high-silica flows are relatively rare and poorly constrained, and thus reflect weaknesses in our universal understanding of lava emplacement processes. Here, we use observations of the 2011-2012 Cordón Caulle rhyolite lava flow as an unparalleled opportunity to study a high-silica flow, and we present the first modeling study to our knowledge that examines the advance of a high viscosity and crystal-poor lava. Our results significantly improve constraints on the rheological controls of high-viscosity flows.
Lava flow advance is thought to be initially controlled by the effusion rate at the vent and topographic slope (Pinkerton and Wilson, 1994). Volume-limited flows are those lavas that cease to advance when fresh material is no longer supplied to the flow front, due to cessation of effusion (Walker, 1971). If, instead, a flow continues to advance, it will eventually reach its cooling-limited length, at which point cooling has formed a crust sufficiently strong to inhibit flow advance (Cashman et al., 1998;Blake and Bruno, 2000;Hon, 2003). Pressure builds within the flow as new material is supplied to the flow front, causing inflation and eventually breaking the crust to form a secondary flow or breakout (Calvari and Pinkerton, 1998;Blake and Bruno, 2000). This is a common process in basaltic eruptions (Walker, 1971;Kilburn and Guest, 1993;Hon et al., 1994;Calvari and Pinkerton, 1998;Applegarth et al., 2010) but has been observed rarely in high silica lava flows (Harris et al., 2004), partly due to the infrequency with which such eruptions occur.
The 2011-2012 eruption of Cordón Caulle in southern Chile (low-silica rhyolite, 69.7 wt% SiO 2 ) permitted the first detailed scientific observations of the emplacement of a crystal-poor, silica-rich lava flow (Silva Parejas et al., 2012;Castro et al., 2013;Schipper et al., 2013;Tuffen et al., 2013;Bertin et al., 2015). Prior to this eruption the only detailed observations of rhyolitic lava emplacement came from the 2008-2009 Chaitén eruption in Chile, in which a composite rhyolite dome was emplaced without the development of a laterally extensive flow field (Bernstein et al., 2013;Pallister et al., 2013). Initial observations of the flow at Cordón Caulle indicated that over the duration of the ∼10 month eruption the flow front stalled, in some places against topographic barriers, and breakouts formed at the flow front and margins, creating a complex flow field . With flow fronts typically 30 m thick, models of conductive cooling suggest that the flow interior may have retained sufficient heat to permit advance for 36-48 months after emplacement (Farquharson et al., 2015), and breakouts were observed to advance for up to 10 months after the eruption ceased Farquharson et al., 2015). The formation of breakouts in both basalts and rhyolites suggests that both compositional endmembers can generate cooling-limited lava flows and, hence, flow front advance may be controlled in the latter phases of emplacement by a cooled crust.
Previous studies of rhyolitic lava flows have focused mostly on higher-silica (72-74 wt% SiO 2 ) Holocene flows in the US such as Big and Little Glass Mountain at Medicine Lake (Fink, 1980b(Fink, , 1983Grove and Donnelly-Nolan, 1986;Manley and Fink, 1987) and Big Obsidian Flow in Newberry Caldera (Castro et al., 2002). Detailed studies of lava flow morphology led to conceptual models in which slowly-spreading crystalpoor lava moved over talus derived from the flow front (Fink, 1980a(Fink, , 1983 before stopping when supply ceased (i.e., volumelimited). This style of advance has been observed in active dacitic flows at Santiaguito, Guatemala (Harris and Flynn, 2002). However, the characteristics of the surface facies of Holocene lava flows additionally inspired models that included a cooled surface crust to account for features such as ogives, formed by crustal folding as the flow front slowed (Fink, 1980a;Castro and Cashman, 1999), crease structures Lescinsky and Merle, 2005), and large gas cavities (Castro et al., 2002). Such crusts have been inferred to significantly influence lava flow growth (Griffiths and Fink, 1993;Fink and Griffiths, 1998), and would also act to thermally insulate the flow, facilitating prolonged mobility and advance (Manley, 1992), as later observed at Cordón Caulle .
Straightforward analytical models to forecast lava flow length have been based on a number of different lava rheologies (Huppert, 1982;Lister, 1992;Kerr and Lyman, 2007;Castruccio et al., 2013). In these models, lava flow advance is assumed to be controlled by either an apparent Newtonian viscosity, a non-Newtonian flow core, or a surface crust of defined strength (Figure 1). Control may switch between rheology types as the flow advances, as suggested by a transition between best-fit models during emplacement (Castruccio et al., 2013). However, for most studied examples of silicic lava emplacement, relatively few parameter values or constraints are available (e.g., syn-eruptive flow lengths and widths only) and models are thus vulnerable to correlated parameter uncertainties which weaken confidence in the resulting interpretations of rheological control. Nevertheless, previous modeling studies have proposed that advancing crystal-rich high viscosity lava flows are yield strength controlled (Castruccio et al., 2013).
Here we examine model robustness by using field data in addition to morphological observations from sequential satellite images of evolving flow fields, to help constrain the time of transition between rheological models. We test the application of these additional observations on a basaltic lava from the 2001 Mt Etna eruption and then extend their application to the 2011-2012 Cordón Caulle rhyolitic flow field. Our overall aims are to (1) assess the consilience of such numerical models with morphological evidence for rheological change and (2) apply lava flow models to the Puyehue-Cordón Caulle lava flow field to explore the role of a thickening crust on the advance of a crystal-poor, high viscosity rhyolite flow. We show that whereas both the Etna and Cordón Caulle lava flows were initially under viscous control, their subsequent development was predominantly controlled by the evolving surface crust. FIGURE 1 | The three flow regimes that can control lava flow behavior: viscosity dominated, crust yield strength dominated, and core yield strength dominated (after Castruccio et al., 2013). The control on flow advance may change throughout the duration of the lava flow emplacement.

FLOW LENGTH MODELS
The models used in this study are given in Kerr and Lyman (2007) and Castruccio et al. (2013), where their derivations are presented in detail. If a lava flow is treated as a Newtonian viscous fluid emplaced into an environment of negligible density, its length L can be given by where C vs is a constant defined as (3/2) 2/3 = 1.31, g is gravitational acceleration, ρ is lava density, β is the ground slope, q is the volume erupted per unit width of the flow front, t is the elapsed time, and η is the lava viscosity (Table 1; Huppert, 1982;Lister, 1992;Kerr and Lyman, 2007). As a lava flow advances, its surface cools to form a crust of thickness δ, which can be assumed to thicken diffusively as δ ∼(kt) 1/2 (Griffiths and Fink, 1993), where k is the thermal diffusivity. If the strength of the crust controls the flow motion, the lava flow length is given by where C c is an unknown constant of the order 1, and σ c is the crustal yield strength (Griffiths and Fink, 1993;Lyman et al., 2005;Kerr and Lyman, 2007). Equation (2) does not include ground slope because, as a cosβ multiplier term (Griffiths and Fink, 1993), for the low slope angles examined here, it can be assumed cosβ ≈ 1 (Lyman and Kerr, 2006). Although useful for first order estimates, Equations (1) and (2), termed here the fixed parameter approach, do not enable down-flow variations in parameters (e.g., slope) to be accounted for. Thus, Castruccio et al. (2013) derived equivalent models, but using a summative approach to enable parameter variation during emplacement. In this case, the length of an advancing flow dominated by a Newtonian viscosity is given by where C 1 is a constant in the viscous regime defined as (3/2) 2/3 = 1.31 (in line with Equation 1), V i is the volume added at each time step, β i is the terrain slope, t i is the time interval between each step, and η i and W i represent the lava viscosity and flow width at a given time step (Castruccio et al., 2013). In the case of flow advance limited by a crustal yield strength, the evolving flow length can be given by where t is the total time elapsed, C 3 is a constant of order 1 (Castruccio et al., 2013), and σ ci is the crustal yield strength at a given time interval (Castruccio et al., 2013). A third case represents flow behavior as a bulk non-Newtonian fluid rather than Newtonian. In this case a yield strength must be exceeded before the flow can advance (Castruccio et al., 2013), and flow length can be given by where C 2 is a constant of order 1, and σ yi is the bulk yield strength at a given time interval (Castruccio et al., 2013). These three models are collectively termed the flexible parameter approach throughout this paper. Here we attempt to improve upon the application of these models by incorporating observations of lava flow advance and morphology to better constrain the timing of key rheological changes during lava flow emplacement.

BASALTIC LAVA FROM 2001 MT ETNA ERUPTION
The 2001 eruption of Mt Etna produced several lava flows from numerous vents (Ingv, 2001;Behncke and Neri, 2003;Coltelli et al., 2007;Applegarth, 2008). The eruption of the largest lava flow began from an altitude of ∼2,100 m near La Montagnola (Figure 2) in the early hours of 18 July 2001 (Behncke and Neri, 2003;Applegarth, 2008). The flow that issued from this vent was the largest of the eruption, reaching over 6 km long and with a volume of 0.02 km 3 . The lava flow was active for 23 days but reached its maximum length after 7 days, after which activity was dominated by breakout development farther up-flow (Behncke and Neri, 2003;Coltelli et al., 2007). As described in Applegarth (2008), the lava initially advanced as a thin (up to 2 m) sheet at ∼150 m per hour before becoming concentrated in a channel at its eastern margin. By the evening of the 19 July the flow had attained a length of 3.4 km and was advancing at a rate of ∼45 m per hour. On the morning of the 20th the flow front had reached an altitude of 1,350 m, and was 100 m wide and ∼5 m high with an incandescent interior visible beneath a cover of black scoria. The flow continued to move down the volcano flanks at a rate of ∼50 m per hour, eventually reaching an area of slightly flatter topography at 1,055 m on 23 July (Behncke and Neri, 2003;Applegarth, 2008). The flow reached its maximum length on 25 July, at which point the flow front was 200-300 m wide and 10-20 m high (Applegarth, 2008). The eruption continued and numerous small breakouts formed at the most distal flow margins and larger breakouts formed further up-flow until the eruption ceased on 9 August 2001 (Behncke and Neri, 2003;Coltelli et al., 2007;Applegarth, 2008). Between 26 and 27 July, a lava flow that had erupted from higher elevations merged with the lower 2001 lava flow (Lanzafame et al., 2003;Coltelli et al., 2007). However, because this was after the main flow had attained its final length , it did not contribute to the overall flow advance, but likely contributed to the later breakout emplacement.

Field Evidence for a Late Crustal control
Field observations of the lava flow were made in April 2015 and focused on morphology and late stage flow features. The main portion of the lava is a channelized 'a'ā flow with welldeveloped levees and a significant portion of drained channel (>2 km long), whereas the distal region lacks levees and has numerous breakouts from its margin. The 2001 flow deposit is almost entirely covered by loose clinker (rubble), making it difficult to observe any surface crust. However, breakouts are common and typically found at flow margins ( Figure 3A) and fronts, together with some areas of drained channel (Applegarth et al., 2010). Most breakouts are hundreds of meters long but some are substantially longer.
The surface crust is visible at the boccas (feeder points) of some breakouts, where it can be buckled upwards due to inflation ( Figure 3B). Where observed, the crust is generally ∼tens of centimeters thick, with a layer of incipiently welded clinker adhering to what would have been the top surface. However, this crustal thickness, formed in the latter stages of the Frontiers in Earth Science | www.frontiersin.org FIGURE 2 | Location of Mt Etna and the 2001 lavas. The flow examined in this work (dark gray and highlighted by the dashed box) was erupted from a vent at 2,100 m elevation. This flow can be seen to be joined by others erupted from higher elevations (light gray), but this occurred after it had attained its final length and so did not contribute to the main flow advance . eruption, may not be representative of the crust earlier in the flow evolution.

Inferred Changes in Rheological Control
As there are relatively sparse published observations of the eruption and spatio-temporal development of the lava flow field, the timing of inferred changes in rheological control is uncertain. However, a reported reduction in flow front advance rate occurred after ∼40 h on a constant slope, and at a time when the vent effusion rate ( Table 3) was increasing (Behncke and Neri, 2003;Coltelli et al., 2007). This suggests that surface cooling had generated a flow-retarding crust that had begun to exert significant influence. The presence of breakouts, formed after the main flow front stalled, indicates the presence of a cooled crust that acted to halt, or substantially slow, the flow front and enabled breakout formation from the flow margins. Flow front halting was likely assisted by a reduced slope gradient. Direct observations of the crust at boccas show the development of a flow surface crust that could have influenced the flow development. The crust is visible only at areas of inflation or breakout formation, providing further evidence for the effectiveness of this crust, as it was able to accommodate an internal pressure increase within the flow due to material continuing to be supplied to the stalled front. It is worth noting that the breakout locations may be substantially controlled by inherent weaknesses in the surface crust. Here we assume flow behavior was initially controlled by its apparent Newtonian viscosity and there was insufficient cooling and crystallinity to acquire a yield strength in the first 40 h.

Application of Flow Models to the 2001 Mt Etna Flow
Previous work based on surface velocity measurements from videos of channelized sections of the flow have estimated a Newtonian flow viscosity of 3.5 × 10 4 -5.6 × 10 5 Pa s for flow in a wide channel, where ρ is the density, g is the acceleration due to gravity, d is the flow depth, β is the slope angle, n is an empirical constant for channelized flow and U is the flow surface velocity. It is assumed that the apparent viscosity η A is equivalent to the viscosity η (Hulme, 1974;Stevenson et al., 2001;Harris et al., 2004;Farquharson et al., 2015). Here we use the central viscosity value of this range −3.0 × 10 5 Pa s (Applegarth, 2008).

Models of Flow Lengthening
To account for uncertainty in the timing of the transition from viscosity to crust-controlled flow, lava flow models are used with transition times of 12, 24, and 36 h (Tables 2, 3), along with observed flow length changes from Applegarth (2008). The modeled lengths using the fixed parameter approach (Equations 1, 2) suggest that the transition from viscosity to crustal control occurred after 12-24 h, although in the case of a transition in control after 24 h the final lava flow length was overestimated by ∼1 km (Figure 4A). The model may overestimate the advance rate during initial lava emplacement but this cannot be verified due to a lack of observed flow length data from this early stage.
Modeled flow lengths using the flexible parameter approach (Equations 3, 4) show best fit to observational data with a viscosity to crustal-control transition after 36 h (Figures 4B). The reduced advance rate in Figure 4B in the crust controlled regime, compared with the model in Figure 4A, is due to changes in the flow width and a reduction in the ground slope that is not accounted for in the fixed parameter approach. The modeled results using values from Castruccio et al. (2013) are shown in Figure 4C and discussed fully in Section Comparison to previously modeled results.

Model Sensitivities
Most of the parameters used in the models are presented or calculated as a range. Here we briefly explore the effect of varying viscosity and crustal yield strength on the models. For viscosity, the plausible range of 3.5 × 10 4 -5.6 × 10 5 Pa s (Applegarth, 2008) leads to final modeled flow lengths that vary by ∼10 km (Figures 5A,B). These substantial differences highlight the model sensitivity to viscosity and emphasize that strong constraints on the timing of rheological control transitions are vital, as well as for parameter values. However, within this range, some model results provide a good match to the observed flow lengthening. Models are similarly sensitive to the crustal yield strength value; for the fixed parameter approach (Equations 1. 2) decreasing yield strength from 10 7 to 10 6 Pa results in a final lava flow length overestimated by >10 km (Figure 5C). For the flexible parameter approach (Equations 3, 4) the final flow length increased by ∼2 km (Figure 5D), with a transition in flow control after 36 h. However, with transition in flow control after 12 or 24 h, the model results give a closer fit to the observed flow front advance, reflecting a correlation between the parameters for crustal yield strength and the inferred timing of the onset of crustal control. Results from the fixed parameter approach (Equations 1, 2) assuming the flow is controlled by its viscosity for 12, 24, or 36 h. (B) Results from flexible parameter approach (Equations 3, 4) assuming that the flow is controlled by its viscosity for 12, 24, or 36 h. (C) Results from flexible parameter approach (Equations 3, 4) using flow viscosities, crust yield strengths and transition times from Castruccio et al. (2013). Changes in flow control are inferred to occur after 3 days.

RHYOLITIC LAVA FROM THE 2011-2012 CORDÓN CAULLE ERUPTION
Puyehue-Cordón Caulle is a basaltic-rhyolitic volcanic complex in the southern volcanic zone of the Chilean Andes (Figure 6). One of the world's most frequently active rhyolitic/rhyodacitic volcanoes, it has most recently erupted in 1921-1922(Katsui and Katz, 1967), 1960(Lara et al., 2004, and 2011-2012. The VEI 5 2011-2012 eruption of crystal-poor rhyolite (69-74 wt% SiO 2 ) began on 4 June 2011, initially producing a 15 km high Plinian column . At the onset of lava effusion, from 15 June 2011, the lava flux exceeded 50 m 3 s −1 (Bertin et al., 2015). Lava effusion was first identified from radar satellite images on 19 June 2011 (Bignami et al., 2014), with lava flowing into a gently-sloping (5-9 • ) topographic depression. Much of the eruption constituted hybrid explosive-effusive activity, with simultaneous lava emission and pyroclastic venting from a common vent Schipper et al., 2013). The eruption ended on 15 March 2012, but parts of the lava flow continued to advance into January 2013 . In total ∼0.4 km 3 lava was emitted , alongside ∼1 km 3 tephra (Pistolesi et al., 2015). We acknowledge that the 2011-2012 Cordón Caulle lava flow lies at the low-silica end of the rhyolitic spectrum, but its high viscosity and crystalpoor nature does allow direct comparison with more typical higher-silica rhyolitic lava flows.

Lava Flow Advance and Morphology
Satellite images from GeoEye-1 and Earth Observing-1 ALI (EO-1 data available from the U.S. Geological Survey, https://eo1. usgs.gov/) show the evolution of the rhyolite flow during the eruption (Figure 7). The images were used to derive the position of the flow front and supplemented by flow outlines derived from COSMO-SkyMed radar data from Bignami et al. (2014).
By 26 June 2011 (d = 22 where d is the number of days from the eruption onset) the flow had split into two channels, one heading north and one heading south, hereafter referred to as the northern and southern flows. These flows advanced synchronously but the northern flow stopped against a topographic high on around 22 October 2011 (d = 140) after advancing ∼3 km. This portion of the lava flow was impeded by topography at both its margins and flow front. The southern flow advanced throughout the ∼285 day eruption, and continued to advance after eruption cessation, reaching ∼4 km in total length. The lava advanced into a shallow topographic basin across a substrate of loose volcaniclastic material cut by numerous small gullies and ridges up to ∼10 m in height, which are minor compared to the 30-40 m thickness of the lava flow.
The first 2 months of the eruption produced lava with a blocky appearance and no obvious large scale structures (Figure 7A), similar in appearance to other rhyolitic and high silica flows such as Big Glass Mountain in the US (Fink, 1980b;Fink et al., 1992;Castro et al., 2002). By 18 August 2011 the flow had developed surface folds hundreds of meters in length, with a wavelength of ∼30-40 m (Figure 7B, d = 75).
By 9 October 2011 (d = 127, Figure 7C) large-scale fractures (hundreds of meters or longer) had formed in the surface crust of the flow (Figures 7C-F). Two fracture types developed-tensional fractures generated by coupling of the crust to the spreading flow core, and shear zones that separated the stalled flow margins from the mobile channel center.
The first breakouts formed from stalled flow fronts and margins by 4 November 2011 (d = 153, Figure 7D). Thereafter, almost all advance in the northern flow, from the front or margins, was via breakout formation ( Figure 7D). Such prevalent breakout formation from this portion of the flow was likely aided by the early stalling of the northern flow against large (up to 100 m high) topographic barriers. In contrast, the southern flow continued to advance during and beyond the remainder of the eruption . Breakouts in the southern flow mostly formed from the flow margins, whilst some small breakouts formed at the flow front (Figures 7F, 8).
A secondary flow, 2 km in length, formed from one breakout on the northern flow after stalling of the main flow front (Figures 7E,F, 8). It began to form around 23 Dec 2011 (d = 202) on the eastern margin of the northern flow (Figures 7E,F, 8). The secondary flow was narrow (100-400 m across), contained within a steep gully (∼15 • slope), and advanced into a forested area on the eastern volcano flank before halting. The southern flow began to advance into topographic barriers around 23 Dec 2011 (d = 202) and this may have impeded lava flow advance. A 50-150 m high ridge spanning the entire western side of the northern and southern flows severely restricted westward spreading of the flow front.
Other flow surface features visible in the satellite imagery include light brown/yellow patches several tens of meters across, and mostly confined to the main flow channel, which represent rafts of flow top pumice. Numerous rough, dark brown patches up to 500 m in length, visible in Figures 7F, 8, are elongate and ogive-parallel far from the vent, but occur as large, more equant patches closer to the vent. These consist of poorlyconsolidated, coarse-grained pyroclastic deposits thought to derive from gravitational collapse of the tephra cone at the vent. Lava effusion ended around 15 March 2012 (d = 285) but local advance of breakouts in the southern flow continued for several months .

Field Evidence for Late Crustal Control
Field campaigns at Cordón Caulle were conducted in January 2013, 2014, and 2015. Fieldwork focused on the morphological features of the flow as well as ground-truthing features observed in satellite images. The lava flow is blocky and ∼30-40 m thick, with a surface consisting of loose debris tens of centimeters to meters across that comprise fragments of surface crust and pumiceous material. This cover of loose debris may have contributed to the thermal insulation of the flow but is unlikely to have had a substantial mechanical influence on the bulk flow advance. Direct evidence for an underlying surface crust is presented by extensional fractures meters to hundreds of meters in length (Figure 3D). The depth of downward penetration of these fractures (several meters to >10 m) indicates a plausible crustal thickness. Further evidence for vertical rheological variation is given by the presence of well-developed ogives (or surface ridges) that are 30-40 m in wavelength and several meters in amplitude.
Breakouts are widely distributed around the flow front and margins, where they display a range of morphologies. Many breakouts developed extensive surface crusts tens of centimeters thick (Figure 3C), which then fragmented as the breakout grew and inflated . Their presence indicates that although the main flow advance had halted, mobile lava remained available in regions close to the flow edges and was likely still being effused from the vent. Breakouts that initiated after vent effusion ceased were fed by lava drainage through the flow field, as observed at other flows (Dietterich et al., 2012).

Inferred Changes in Rheological Control
In its initial stages, flow advance may have been controlled either by an apparent Newtonian viscosity or a non-Newtonian (yield strength) rheology. However, in its latter phases the field and satellite observations of folds and fractures suggest that a surface crust exerted a strong influence. Surface folds (ogives) form in response to a slowing flow front and subsequent buckling of a rheologically stiff surface crust behind it (Fink, 1980a). The extensional fractures that formed after ogive formation appear similar to those generated in analog models of lava flows with a strong brittle crust (Lescinsky and Merle, 2005;Applegarth, 2008), and their orientation suggests that some fractures relate to the formation and spreading of flow front lobes. We infer that when surface fractures were produced the flow crust was significantly influencing flow advance. Furthermore, the presence of breakouts points toward sufficient crustal strength to halt overall advance, despite the availability of mobile material in the flow core Farquharson et al., 2015).

Calculation of Bulk Flow Viscosity
In order to apply Newtonian flow lengthening models to the Cordón Caulle lava, the bulk flow viscosity must first be estimated using Equation (6) (Jeffreys, 1925). We use a lava density of 2,300 kg m −3 , a flow depth of ∼30 m at the flow front (Farquharson et al., 2015), a slope of 7 • (determined from a pre-eruption Aster DEM along the length of the southern flow), and n = 2 for flow in a wide channel (Hulme, 1974;Stevenson et al., 2001;Harris et al., 2004;Farquharson et al., 2015).
In order to determine the lava flow surface velocity, surface features in the main channel of the flow (ogives, fractures in the flow crust, and pale pumice rafts) were tracked between two georeferenced satellite images taken in 9 October 2011 and 26 January 2012 (Figure 8). Surface velocities of 2.2-6.2 × 10 −5 m s −1 were derived (mean 3.8 × 10 −5 m s −1 ), giving a viscosity range of 2.0-5.6 × 10 10 Pa s (mean 3.6 × 10 10 Pa s).

Models of Flow Lengthening
Due to the northern flow being topographically constrained at its front, lava flow advance models are only relevant to the initial advance of the southern flow. We assume that the lava flow was initially controlled by either a Newtonian or non-Newtonian bulk rheology for the first 75-127 days. Henceforth, the presence of ogives and surface fractures indicate that surface crust had attained sufficient thickness and strength to exert a controlling influence on flow behavior. Model results should be disregarded for periods after the southern flow front began to interact with a topographic barrier (202 days). However, by this time the major changes in rheological control had already occurred, with the lava firmly in the crustal-control regime and breakouts forming at many points around the lava flow. The effusion rates (Table 3) are cited for the entire flow field (Bertin et al., 2015). As a best estimate, we assume that the north and south flow branches were fed approximately equally-i.e., each at half the rate given in Table 3.
Models based on the fixed parameter approach (Equations 1, 2, utilizing values in Table 2) underestimate final lava flow length by ∼1,500-1,000 m for transitions from viscosity to crustal control after 75 and 127 days respectively ( Figure 9A). The biggest disparity is in the early stages of the eruption, when the observed lava flow advance rate far exceeded model predictions, suggesting that early-erupted lava was substantially less viscous than the later-erupted lava for which we have viscosity estimates (Farquharson et al., 2015). Such increases in the lava flow viscosity with time likely relate to cooling, degassing and crystallization of the lava flow.
Models that employ the flexible parameter approach (Equations 3-5) and assume a core yield strength control underestimate the final flow length by ∼2,000 m ( Figure 9B). When an initial viscosity control and 75 day crustal control transition are used the final lava flow length is underestimated by ∼900 m. However, when the crustal control transition occurs after 127 days, the modeled final lava flow length closely fits with observed values. Not surprisingly, the flexible parameter approach provides a better fit for the early stages of the lava evolution than the fixed parameter approach, because it accounts for variations in lava flow properties, such as the flow width.

Inferring Cordón Caulle Flow Properties
The Cordón Caulle lava flow gained most of its length in the first 2-3 months of the eruption, during inferred viscosity-controlled advance. However, reliable viscosity estimates could only be made from the latter stages of the eruption when recognizable surface features could be tracked and the flow was less obscured by the eruption plume. Thus, viscosity estimates are upper bounds because they are likely to be affected by the presence of a crust, and the later lava was likely cooler, more degassed and more crystalline than during the earlier stages (Schipper et al., 2015). Flow length models offer an alternative approach to constrain plausible values of flow viscosity, crust yield strength and internal yield strength.
The optimized fixed parameter model provides a best fit to observed flow lengthening for a viscosity of 3.3 × 10 9 Pa s, a crustal yield strength of 2.5 × 10 8 Pa and a transition from viscosity to crustal control after 60 days ( Figure 10A). This lower viscosity, an order of magnitude less than estimates for later lava (Farquharson et al., 2015), could be representative of initially  effused lava that was hotter, less crystalline and more volatilerich. A viscosity value can also be estimated directly using the Einstein-Roscoe equation: where η 0 is the melt viscosity which can be determined from Giordano et al. (2008), R and Q are constants equal to 1.67 and 2.5 respectively, and Φ is the crystal packing fraction. Using an eruptive temperature of 900 • C, a melt composition from Schipper et al. (2015) with a volatile content of 0.1%, and an initial crystal content of ∼30% at the vent (Schipper et al., 2015) gives a viscosity of ∼1.5 × 10 9 Pa s. This is in close agreement to the fixed parameter model suggested here (Equations 1, 2, Figure 10A). We acknowledge that the Einstein-Roscoe equation is based on the assumption that the crystal population is isotropic. However, in the Cordón Caulle flow most crystals are rod-like microlites (Schipper et al., 2015), and such high aspect ratio crystals increase the magma viscosity more strongly than isotropic crystals (Mueller et al., 2011;Mader et al., 2013). Accordingly, this estimate could be taken as a minimum value, and viscosity also likely increased with time due to degassing and crystallization of the lava flow. Alternatively, the alignment of an initial randomly orientated microlite population could lead to shear thinning behavior during lava flow emplacement and a reduction in apparent viscosity. Despite these limitations, the value independently derived from the Einstein-Roscoe equation does support the viscosity value derived through the optimized lava flow models, Equations 1, 2.
Optimizing a flexible parameter model with Newtonian viscosity and subsequent crustal control (Equations 3, 4) gives the best fit to the observed flow lengthening for a rather higher viscosity of 1 × 10 10 Pa s, a lower crustal yield strength of 4 × 10 7 Pa and a similar flow control transition time of 50 days ( Figure 10B). A flexible parameter model that assumes initial core yield strength control provides best fit for a core yield strength of 9.8 × 10 4 Pa, a crustal yield strength of 4 × 10 7 Pa and a transition in flow control at 50 days. As these models are inherently non-unique and other combinations of flow properties and transition time in flow control can yield equally good fits to observations, it is valuable to compare the determined flow parameters to previously published results. The inferred values of crustal yield strength correspond with those proposed in other studies of silicic lava flows and domes (Griffiths and Fink, 1993;DeGroat-Nelson et al., 2001) and the lower core yield strength determined here is in line with other estimations of high silica content lava yield strength (Blake, 1990;Fink and Griffiths, 1998). The inferred viscosity from the idealized flexible parameter approach (Equations 3, 4, Figure 10B) is in line with the lower viscosities calculated by Farquharson et al. (2015) for breakouts (1.21 × 10 10 -4.01 × 10 10 Pa s) from the flow margin but is slightly lower than the inferred bulk lava viscosity inferred from satellite observations.

Model Sensitivities
The uncertainties associated with model input result in uncertainties in model outputs; for example varying slope angle by 2 • in all models leads to final flow length changes of ∼2-10%, corresponding to flow length variations up to a few hundred meters. However, variations due to uncertainties in crustal yield strength and flow viscosity are much greater. Applying the range of viscosity values derived from initial range in flow channel velocities leads to a range in the final modeled flow lengths of ∼1.5 km for the fixed parameter approach (Figure 11A) and 2 km for the flexible parameter approach (Figure 11B).
Previously published values for the crustal yield strength of silicic lavas range from 10 6 to 10 8 Pa (Griffiths and Fink, 1993;Bridges, 1997;Fink and Griffiths, 1998;DeGroat-Nelson et al., FIGURE 11 | Model results for the Cordón Caulle lava flow. The black arrow shows the time at which the flow starts to interact with a topographic barrier. (A) Range of model results when viscosity is varied using the fixed parameter approach (Equations 1, 2). The red area represents a transition in rheological control after 127 days and the gray area represents a transition after 75 days. Dashed lines show the high and low viscosity end members and the solid lines show the mean viscosity. (B) Similar to (A) but for a flexible parameter approach (Equations 3, 4). (C) Fixed parameter approach (Equations 1, 2) where the crust yield strength is reduced by an order of magnitude to 10 7 Pa. (D) Similar to (C) but for a flexible parameter approach (Equations 3-5). Flow is initially controlled by either its viscosity or core yield strength. Two core yield strength curves give similar results for a transition in control after 75 days and 127 days.
2001; Kerr and Lyman, 2007;Castruccio et al., 2013). Reducing the crustal yield strength by an order of magnitude to 10 7 Pa in the fixed parameter approach (Equations 1, 2) produces a significant increase in modeled lava flow advance rate in the crustal control regime (Figure 11C) compared to the original model (Figure 9A), and increases the final lava flow length by >1.5 km. An earlier transition in flow control (at day 75), results in a longer lava flow than when the transition occurs after 127 days. This is because at this time in the modeled scenario, the crust provided a lesser retarding force than the flow viscosity, leading to accelerated lava advance rates upon flow control transition. In reality though, the lava flow would have retained the viscous retardation and remained viscosity controlled until the crust became sufficiently strong to dominate over the viscous forces. Thus, this apparent increase in lava flow advance is an artifact of the modeling approach that considers rheological controls in isolation rather than combination. The effect on the flexible parameter models of lowering crust yield strength (Equations 3, 4) is not as pronounced (Figure 11D), but leads to an increase in the final flow length compared to the original model (Figure 9B), and provides a closer fit to the actual flow lengthening.

DISCUSSION
Field and remote sensing observations, as well as straightforward flow advance models, suggest that both the 2001 Etna basaltic flow and the 2011-2012 Cordón Caulle rhyolitic lava flow were controlled in their latter stages by a cooled crust after an initial viscous control. The cooled crust acted to retard the lava flows, ultimately halting them before breakouts formed at the flow front and along flow margins (Behncke and Neri, 2003;Coltelli et al., 2007;Applegarth et al., 2010;Tuffen et al., 2013).

Model Limitations
The models are limited by the uncertainties in parameter values and inherent in-model assumptions (aleatory and epistemic uncertainties). The models themselves are straightforward and their ability to simulate flow lengths across a broad compositional range demonstrates their generality (Kerr and Lyman, 2007;Castruccio et al., 2013). The flexible parameter approach (Equations 3-5) will be most appropriate for lava flows that undergo substantial widening, such as the Cordón Caulle rhyolite flow, and likely explains the disparity among the model results. None of the models are valid for sections of uphill slope (i.e., negative β values) and so breakdown when topographic barriers are encountered. For this reason, averaged slope angles are typically used.
Values of crustal strength are poorly reported in the literature and often rely on scaling analysis from analog experiments (Griffiths and Fink, 1993;Fink and Griffiths, 1998). The evolution of crustal strength through time is also unknown and very difficult to determine. The effusion rate is relatively well constrained for the 2001 Mt Etna eruption (Behncke and Neri, 2003;Coltelli et al., 2007). However, the effusion rate for the 2011-2012 Cordón Caulle eruption is less well constrained, with only preliminary values currently available. Due to lack of data we assume that the north and south flows received an equal share of the vent effusion rate. Because both flows were continuously fed and active throughout the eruption (either through breakout formation or main flow advance) this assumption appears reasonable.
Error in imagery-derived lava flow velocity estimates mostly relate to potential georeferencing problems, although physical effects, such as compression of the flow surface during ogive formation giving an apparent shortening of the flow's surface, will also contribute. Such errors are likely small compared with other uncertainties and the velocity values calculated encompass those estimated by Farquharson et al. (2015) for breakouts in the latter phases of flow emplacement (3.57 × 10 −5 m s −1 ). The topographic slope may also have changed during the eruption, as Castro et al. (2016) have shown that considerable uplift (>200 m) occurred around the Cordón Caulle vent in the first month of the eruption. This uplift, which is attributed to the emplacement of a shallow laccolith intrusion, may have increased the average gradient along the length of the southern flow by up to ∼3 • . However, such a magnitude of slope change has a minor effect on the model outputs when compared with other influences, as discussed in Section Model sensitivities.
The frequency of available observational data remains a substantial factor within the uncertainties of model outputs. For Cordón Caulle, the precision in the timing of the transition between viscous and crustal control (identified through observations of crustal fracturing or the formation of features such as ogives and breakouts), is limited by the infrequent satellite image acquisitions (one to two partial images per month).
Despite these uncertainties, the models provide reasonable estimates of the evolution of compositionally diverse lava flows. Nevertheless, is clear that good knowledge of the evolving flow architecture can help constrain models and identify when underpinning assumptions breakdown. Morphological considerations can thus aid model selection and improve the accuracy of estimated parameter values, as well as helping identify where the greatest remaining uncertainties lie.

Comparison to Previously Modeled Results
Due to the additional constraints that we impose, our model results for the 2001 Mt Etna flow have a poorer fit to the observed flow lengthening (Figure 4) than those shown previously by Castruccio et al. (2013), in which more parameters were allowed to vary within the fitting process. Castruccio et al. (2013) suggests that the 2001 Etna flow advance was either controlled by an initially low viscosity (3.8 × 10 6 Pa s), which increased after 3 days (2.4 × 10 8 Pa s), or by an initial low viscosity (4.5 × 10 6 Pa s) and then by the crustal yield strength (1 × 10 6 Pa, Figure 4C). We have initially constrained crustal yield strength to a greater value in our models. However, when a value closer to that of Castruccio et al. (2013) is used, the modeled fit in the flexible parameter approach is much improved (Figure 5D). This suggests that the crust of the 2001 Etna flow may have been weaker than indicated by Griffiths and Fink (1993). The models presented here, as well as flow observations, suggest that (1) the transition in the lava flow control at Mt Etna occurred much sooner, within 36 h, than the 3 days suggested by Castruccio et al. (2013), and (2) the initial lava flow viscosity was considerably lower.
For high viscosity (dacite) lava flows, Castruccio et al. (2013) suggested a core yield strength control. However, our findings suggest that the high viscosity and crystal poor Cordón Caulle lava flow was initially viscosity-controlled. This arguably reflects a fundamental difference between the emplacement of crystal-poor and crystal-rich high-silica flows. Rhyolitic lavas are typically extruded as degassed, crystal-poor melt, because their high melt viscosity makes diffusion sluggish (Eichelberger et al., 1986;Eichelberger, 1995). In contrast, other high viscosity lavas including andesites and dacites typically have a significantly higher crystal content (Sparks et al., 2000), making a yield strength control far more likely. Indeed, the dacitic lava at Santiaguito is significantly more crystal-rich, 47% by volume (Rose, 1987;Holland et al., 2011), than the initially effused Cordón Caulle rhyolite.

Development of a Yield Strength
Thus, the advance rate of some high silica lava flows is thought to be controlled by their core yield strength (Castruccio et al., 2013), which develops when the lava crystal content exceeds ∼30% volume (Saar et al., 2001;Castruccio et al., 2010). However, it is currently unclear whether the Cordón Caulle rhyolite could have developed substantial core yield strength during emplacement, as measured crystal contents in samples taken from the vent region are only ∼10-30% (Schipper et al., 2015), and generally beneath the threshold proposed by Saar et al. (2001). Nevertheless, such vent samples only represent the last-erupted material and may not typify earlier-effused material, which could have had a lower crystal content.
In contrast, samples from breakouts exhibit significantly higher crystal content of 50-60% (Schipper et al., 2015), clearly sufficient to generate a yield strength. Nevertheless, given the slow crystallization kinetics of rhyolites (Swanson et al., 1989), it is unclear whether this crystal population could have formed before lava flow advance became dominated by a crust. Swanson et al. (1989) suggest the bulk of crystallization in rhyolites occurs only after lava has halted and had time to cool. Furthermore, much of the crystal population within a rhyolite lava flow may be products of devitrification, such as spherulites or globulites (Schipper et al., 2015), which largely form after lava has halted. Further work on detailed textural studies is needed in order to reconstruct the crystallization history relative to the emplacement duration, and unravel to what effect, if any, crystal populations had on the evolution of an effective crust and broader flow processes.

Breakouts and Crustal Control
Observations of Cordón Caulle demonstrate that not only was the lava flow controlled in its latter phases by a cooled crust, but the crust was able to halt the flow leading to breakout formation at the stalled front and flow margins. Such a strong crustal control had been inferred at other high silica flows Castro et al., 2002), but the absence of breakouts at these flows is intriguing. As the first Cordón Caulle breakouts formed only after ∼150 days, the absence of breakouts at other rhyolite lava flows may reflect briefer emplacement timescales of lava flows that were volume-limited. Alternatively, a lack of breakouts may reflect the inability of other lava flows to develop a substantial surface crust, potentially indicative of higher flow advance rates than at Cordón Caulle caused either by higher effusion rates (Loewen et al., 2017) or steeper topography. In both cases, the greater lava flow velocity is likely to apply greater stress to the surface crust through developing greater velocity gradients, causing it to fail more regularly and potentially preventing development into a pervasive surface capable of halting the flow (Anderson et al., 1995;Kilburn, 2004). Similarly, pulses of higher effusion rate could result in periods of increased flow rate, leading to increased disruption of the flow surface crust. Furthermore, higher effusion rates would allow the flow to advance farther before cooling became a substantial controlling factor (Fink and Griffiths, 1998). Alternatively, surface crusts at other high silica flows could be sufficiently strong to prevent breakout formation in the first place, with the internal pressure failing to exceed the crustal yield strength.
A combination of these processes may have contributed to the lack of breakouts from the late north-east lobe of the Cordón Caulle flow (Figure 8). The north-east lobe displays pronounced levees and a possible area of drained channel, suggesting that its emplacement was volume limited. This lobe advanced down a much steeper slope (∼15 • ) than the main flow field, potentially leading to increased deformation of any surface crust by the greater stresses involved. In contrast, the bulk of the flow field advanced over very gentle topography and, in some cases, abutted topographic barriers, which could have aided crust formation through putting the top surface under compression, reducing the degree to which the surface crust was ruptured by extensional fracturing.
Given the presence of very large fractures in the Cordón Caulle flow surface, the crust clearly behaved brittly for part of its lifetime. The models used here assume that the crust acts as a deforming viscoelastic layer that controls the lava flow advance. This layer is likely bounded by a cooled, strong brittle layer that would fragment. The energy used in the fragmentation of a brittle surface layer can be accounted for (Kilburn, 2004) and future models could attempt to incorporate this additional factor in lava flow lengthening models. It is also worth noting that, in many large basaltic lava flows, the crust not only halts flow advance but enables prolonged supply to the flow front through effective insulation of the material in the lava flow core (Anderson et al., 1999;Vye-Brown et al., 2013). This is also true of the Cordón Caulle rhyolite flow, where material continued to be supplied beneath the stalled crust to breakouts at the flow front and margins, enabling continued lava flow transport long after the eruption ended Farquharson et al., 2015).

Crustal Thicknesses
Lava flow crustal thicknesses (δ) can be estimated at the time of transition from viscosity to crustal control if we assume that the crust grows diffusively by δ∼(kt) 1/2 (Griffiths and Fink, 1993), and with values of k given in Table 1. In the case of the 2001 Mt Etna flow, this suggests crustal thicknesses of 0.1, 0.15, and 0.2 m after 12, 24 and 36 h respectively, and 0.4 m by the time the lava flow reached its maximum length; broadly similar to thicknesses observed in the field. The Cordón Caulle flow would have attained a crustal thickness of ∼1.9 and 2.5 m thickness after 75 and 127 days, and 4 m after 285 days at the eruption end. This is also a similar scale to the depths of fractures observed in the field and the estimated crust thickness at other rhyolite flows from fractures and ogive wavelengths (Fink, 1980a).
The point at which a lava flow becomes crustal controlled could be considered as a ratio between crustal thickness and the bulk flow thickness (Lescinsky and Merle, 2005). The timescale of crustal thickening is primarily controlled by the thermal diffusivity, with higher diffusivities, e.g., associated with low vesicularity, leading to a faster thickening and hence more rapid transition to crustal control. A crust with a higher yield strength ought to become a controlling factor more rapidly than a lower yield strength crust because it would not have to gain as great a thickness to have a commensurate effect on lava flow emplacement. The timing of the onset of crustal control must therefore be sensitive to both the initial crystallinity of the lava and the rate of in-flow crystallization.
The ratio of crustal thickness to flow thickness (Ψ ) can be calculated for the moment of inferred viscosity to crustal control transition. In the case of the 2001 Etna flow, Ψ is equal to 0.02-0.04, depending on the inferred transition time. Similarly, for Cordón Caulle Ψ is equal to 0.06-0.08. We thus infer that the crust becomes a dominant influence when it accounts for 2-8% of the bulk flow thickness, but this proportion could be lower for crust with a higher yield strength. Such dimensionless numbers help to systematically define the timing of viscosity to crustal control transition. Crustal strength is clearly a key, but poorly understood, aspect of lava flow emplacement. Improved constraints are required to define and quantify the factors that influence it, which must include lava texture (i.e., crystallinity, vesicularity), composition, and cooling rate. Such an improved understanding could aid in defining the value of Ψ required to control lava flow advance rate for a broad compositional range.

CONCLUSIONS
Observations of lava flow emplacement from published studies, field data and satellite imagery indicate that a strong crust may significantly influence the growth of both basalt and high viscosity, crystal poor, rhyolite lava flows. By using these observations to help constrain numerical models of lava flow length change over time, we show that both the 2001 Mt Etna basaltic lava flow and the 2011-2012 Cordón Caulle rhyolite lava flow were initially controlled by their apparent Newtonian viscosities and then by their surface crusts. In the case of the rhyolite, this change in control occurred only after several months and our results underscore that relatively straightforward lava advance models, based on bulk rheological properties, can be effective across a large compositional range of lava flows. The emplacement of high silica flows is not always dominated by their internal yield strength, as previously inferred, but crystal-poor high-silica lavas, including rhyolites, may behave in a manner similar to initially viscosity controlled basalts, albeit over different temporal and spatial scales.

AUTHOR CONTRIBUTIONS
All authors contributed to the project design and preparation of the manuscript. NM, MJ and HT conducted fieldwork at Cordón Caulle. NM and MJ conducted fieldwork at Mt Etna. NM analyzed satellite imagery and carried out the numerical modeling assisted by MJ and HT.

FUNDING
NM is supported by a NERC Envision studentship and a BUFI grant from the British Geological Survey. HT is supported by a Royal Society University Research Fellowship. CV publishes with the permission of the Executive Director of the British Geological Survey.