Resumen

Los lahares son flujos masivos peligrosos que han causado hasta 23,000 víctimas en el pasado reciente. El volcán Cotopaxi tiene una historia documentada de lahares destructivos por su casquete glaciar, período de erupciones y áreas densamente pobladas. Anteriores modelos mecánicos de lahares a menudo no capturan la trayectoria, velocidad y presión en la base. Se usa el modelo de flujo de escombros RAMMS con un algoritmo de arrastre para calcular el lahar de 1877 y evaluar el riesgo futuro de lahares. El modelo calibrado se aplicó a futuros escenarios de erupción y se espera un crecimiento potencial de su volumen del 50 al 400%. Se desarrolla un enfoque de Voellmy-Salm para predecir futuros lahares.

Fecha	2020
Formato	Artículo
Autores	Frimberger T, et al
Editora	Sociedad Británica de Geomorfología
Colección	Earth Surf. Process. Landforms. 2021;46:680–700
Fuente	https://onlinelibrary.wiley.com/doi/full/10.1002/esp.5056

Contenido

Modelado de futuros lahares controlado por diferentes escenarios de erupciones volcánicas del Cotopaxi (Ecuador) y calibrado con su lahar masivamente destructivo de 1877

Resumen

Los lahares se encuentran entre los procesos de flujo masivo más peligrosos del planeta y han causado hasta 23 000 víctimas en eventos únicos del pasado reciente. El volcán Cotopaxi, 60 km al sureste de Quito, tiene una historia bien documentada de lahares destructivos masivos y es un punto crítico (hotspot) para futuros lahares debido a (i) su casquete glaciar de ~10 km2, (ii) su período de retorno de 117 a 147 años de erupciones (sub)plinianas, y (iii) las zonas potencialmente inundables densamente pobladas (300 000 habitantes). Los modelos mecánicos anteriores de lahares a menudo (i) no capturan la trayectoria inicial empinada del lahar, (ii) no reproducen múltiples rutas de flujo, incluidas la bifurcación y confluencia, y (iii) no generan parámetros clave apropiados como la velocidad del flujo y la presión en la base como una medida de la capacidad de erosión. Recalculamos el bien documentado lahar de 1877 utilizando el modelo de flujo de escombros RAMMS con un algoritmo de arrastre implementado que cubre toda la trayectoria del lahar desde el edificio del volcán hasta una extensión de ~70 km desde la fuente. Para evaluar la sensibilidad y restringir el rango de entrada del modelo, exploramos sistemáticamente los valores de los parámetros de entrada, especialmente los coeficientes de fricción \mu y \xi de Voellmy-Salm. La selección objetiva de las combinaciones de parámetros más probables permite una representación realista y robusta del riesgo de lahares. Los registros históricos detallados de la altura del flujo, velocidad del flujo, descarga máxima, tiempo de viaje y límites de inundación coinciden más con una fricción \mu muy baja de tipo Coulomb (0.0025–0.005) y una fricción turbulenta \xi alta (1000–1400 m/s²). Finalmente, aplicamos el modelo calibrado a futuros escenarios de erupción (Índice de explosividad volcánica VEI 2–3, 3–4, >4) del Cotopaxi y los lahares escalados acordes. Por primera vez anticipamos un crecimiento de volumen potencial de 50 a 400% debido a la erosividad del lahar en los flancos escarpados del volcán y desarrollamos un enfoque genérico Voellmy-Salm en diferentes escalas de lahares de gran magnitud y mostramos cómo se puede utilizar el modelo para anticipar futuros lahares sineruptivos.

Palabras clave

Erosión de flujo de escombros, lahar, calibración del modelo, modelo numérico, modelado predictivo.

1. Introducción

El término indonesio ‘lahar’ se refiere a una mezcla de sedimentos y agua que fluye rápidamente y que se origina en un volcán (Smith y Fritz, 1989), ya sea directamente relacionada con la actividad eruptiva (es decir, lahar primario/sineruptivo) o durante los períodos posteriores a la erupción (es decir, lahar secundario). Estos flujos masivos se encuentran entre los peligros naturales volcánicos más desastrosos (Auker et al, 2013). El potencial destructivo de los lahares surge de su inicio repentino y características de flujo masivo (que se caracterizan por una alta capacidad de transporte de sedimentos), velocidades que superan las decenas de kilómetros por hora, descargas de hasta aproximadamente 50 000 m³/s y distancias de escape de decenas a incluso cientos de kilómetros (Mothes et al, 2004; Thouret et al, 2020). Se han documentado alrededor de 56 000 muertes relacionadas con lahares desde 1500 dC (Brown et al, 2017), incluidas alrededor de 23 000 muertes en el Nevado del Ruíz en Colombia en 1985 (Pierson et al, 1990). A medida que crece rápidamente la densidad de población y el desarrollo de infraestructura cerca de los volcanes activos, un número cada vez mayor de personas está potencialmente expuesta a peligros volcánicos (Auker et al, 2013; Small y Naumann, 2001). En las últimas décadas se han realizado esfuerzos para reducir el riesgo para la vida humana mediante la implementación de sistemas efectivos de monitoreo y alerta temprana y mediante la evaluación de áreas probables de inundación de lahar mediante el modelado y la reconstrucción de la historia eruptiva. El volcán activo Cotopaxi, coronado por un glaciar, se encuentra entre los volcanes más peligrosos del mundo, porque puede producir lahares desencadenados por erupciones muy grandes, que pueden inundar valles densamente habitados cerca del volcán. En este estudio, nuestro objetivo es contribuir a la evaluación de peligros de los lahares sineruptivos probando la aplicabilidad de un modelo de flujo de escombros basado en Voellmy que incluye el arrastre y analizando su capacidad para reproducir las características clave de los lahares.

1.1. Características de los lahares

El comportamiento del flujo de los lahares está fuertemente influenciado por la concentración de sedimentos, tipo de sedimento y distribución del tamaño de las partículas (Pierson et al, 1987; Vallance y Iverson, 2015) y puede variar desde flujos hiperconcentrados, que transportan entre 5-10 a 20–60 vol% de sólidos y muestran densidades aparentes de 1300–1800 kg/m³ (Beverage y Culbertson, 1964; Pierson, 2005), hasta flujos de detritos con concentraciones de sólidos que exceden 20–60 vol% y densidades aparentes que varían entre 1800 y 2300 kg/m³ (Costa, 1984; Iverson, 1997; McArdell et al, 2007). Varios estudios de campo han revelado transiciones de fase de lahares reversibles durante la propagación aguas abajo (Pierson y Scott, 1985; Scott, 1988) a medida que los lahares incorporan (aumentan) o depositan (disminuyen) sedimentos y pierden (infiltran) o arrastran (diluyen) agua. Especialmente los lahares no cohesivos (<3–5% de partículas del tamaño de la arcilla) son propensos a tales transiciones de fase, típicamente evolucionando desde flujos de escombros proximales a flujos hiperconcentrados o flujos de corrientes fangosas en tramos distales (Pierson y Scott, 1985).

La recopilación directa de datos cuantitativos sobre eventos de lahares es intrínsecamente difícil y especialmente rara para flujos sineruptivos a gran escala. La comprensión de las principales características físicas y reológicas de los lahares se beneficia de las observaciones de estudios experimentales (Iverson et al, 2010) y flujos de detritos naturales (McArdell et al, 2007; McCoy et al, 2012). Las características clave incluyen una progresión pulsante (Doyle et al, 2011; Iverson, 1997; McArdell et al, 2007), una presión de fluido intersticial de larga duración superior a los valores hidrostáticos que contribuyen a largas distancias de agotamiento (Iverson, 1997; Lube et al, 2009; McArdell et al, 2007) y la erosión y arrastre topográficamente controlados que conducen a un aumento de varias veces del volumen de descarga y flujo después del inicio (Berger et al, 2010; Berti et al, 1999; Manville, 2004; McCoy et al, 2012; Pierson et al, 1990; Scott et al, 2005).

1.2. Enfoques de modelado de lahares

La incorporación de la amplia gama de comportamientos de flujo y los procesos dinámicos característicos de los lahares en los análisis de descentramiento plantea una tarea desafiante, por lo que cualquier modelo teórico o herramienta práctica de simulación necesariamente se basa en simplificaciones físicas y reológicas de diversos grados. Se han propuesto diferentes sistemas informáticos que simulan lahares, incluido el modelo empírico ampliamente utilizado LaharZ (Iverson et al, 1998), que predice los límites de inundación (Delaite et al, 2005; Huggel et al, 2008; Lee et al, 2015; Muñoz-Salinas et al, 2009; Pistolesi et al, 2014; Schneider et al, 2008; Worni et al, 2011), el modelo de autómata celular LLUNPIY (Lupiano et al, 2020; Machado et al. al., 2015) y una variedad de modelos numéricos físicos que aplican diferentes supuestos reológicos. Entre estos, los modelos 1D son hoy en día reemplazados en gran medida por modelos continuos multiescala con ecuaciones promediadas en profundidad para la conservación de masa y momento. Entre las leyes reológicas, se aplican con frecuencia modelos basados en dinámica de fluidos que propagan lahares como flujos newtonianos. Dichos modelos (ej., el código dinámico de enrutamiento de inundaciones Delft3D) pudo reproducir bien los lahares de ruptura de lagos diluidos del Monte Ruapehu, Nueva Zelanda (Carrivick et al, 2008). Sin embargo, modelar los lahares como fluidos newtonianos no considera que los lahares difieren sustancialmente en la disipación de energía en comparación con los flujos de inundación newtonianos típicos. Por lo tanto, los modelos fluidodinámicos pueden extenderse mediante el coeficiente de Manning o Chézy para aproximarse a la resistencia del flujo interno. El software de simulación hidráulica Flo-2D (O’Brien et al, 1993), diseñado para inundaciones, flujos de lodo y flujos de escombros, utiliza el modelo reológico viscoplástico de Bingham y abarca términos de resistencia al flujo para tensiones de cedencia, viscosidad, turbulencia y dispersión. Flo-2D es uno de los modelos reológicos monofásicos más utilizados en el análisis de riesgo de lahares (Caballero y Capra, 2014; Worni et al, 2011), mientras que otros enfoques reológicos, como el modelo de Voellmy, se han aplicado pocas veces, aunque se ha formulado (y demostrado su éxito) para simular varios tipos de flujos cargados de sedimentos. Modelos tridimensionales más sofisticados basados en la teoría de la mezcla de Coulomb (Iverson y Denlinger, 2001) describen un flujo de dos fases y las interacciones entre fase sólida y fluida. En el código de simulación Titian2D  se incorpora una formulación simplificada promediada en profundidad de la teoría de mezcla (Pitman y Le, 2005) y se ha aplicado en algunos estudios para lahares cargados de sedimentos (Procter et al, 2010; Williams et al, 2008).

Con el aumento de la complejidad del modelo y su precisión física, es necesario definir un número creciente de parámetros de entrada, lo que requiere el muestreo y las pruebas de materiales adecuados (Iverson y George, 2016). Los modelos complejos también son computacionalmente costosos, lo que restringe el modelado de descentramiento a una pequeña sección específica de la trayectoria del flujo (Procter et al, 2010). Cuando se modelan lahares de gran magnitud a lo largo de largas distancias, es necesario lograr un equilibrio entre la solidez y su plausibilidad física: los modelos robustos y computacionalmente manejables son fundamentales para su aplicación en sistemas de drenaje a gran escala, al mismo tiempo que capturan el flujo clave. características como alturas de etapa, velocidad de flujo, presión de cizallamiento basal, tiempos de viaje, confluencia, bifurcación y aceleración. Por lo tanto, en este estudio se elige el software de simulación numérica RAMMS para flujo de escombros, ya que proporciona una buena combinación de simplificación y plausibilidad física y permite rastrear lahares cargados de sedimentos en el Cotopaxi a lo largo de toda su trayectoria en grandes cuencas.

1.3. Erosión en modelos de lahares

El arrastre de material erosionado en una masa que fluye puede tener un impacto significativo en la propagación del flujo (Iverson, 1997; McDougall y Hungr, 2005) y por lo tanto, debe considerarse en el modelado de descentramiento y evaluación de peligros (Frank et al, 2015; Iverson y Ouyang, 2015). Aunque no existe consenso sobre si el arrastre del material del sustrato reduce o mejora la movilidad del flujo (Iverson et al, 2011; Mangeney, 2011; Pudasaini y Fischer, 2016), se ha observado que un aumento en el volumen del flujo se correlaciona positivamente con la descarga máxima. área potencial de inundación, distancia de escape, altura y velocidad del flujo (Rickenmann, 1999). La cantidad de modelos numéricos capaces de simular el arrastre de flujos de masa está en aumento actualmente (Iverson y Ouyang, 2015) y el rendimiento de las teorías subyacentes se prueba y valida cada vez más frente a eventos reales de flujo másico (Dietrich y Krautblatter, 2019; Frank et al. , 2017; Hungr y McDougall, 2009; Hussin et al, 2012). La simulación de procesos de arrastre requiere la entrada de tasas de crecimiento específicas del usuario o de tasas de erosión basadas en procesos en función de la velocidad (Fagents y Baloga, 2006) o esfuerzo cortante (Frank et al, 2015; Iverson, 2012). Aunque los lahares desarrollan una fuerza erosiva considerable relacionada con la altura del flujo, densidad del flujo y esfuerzo cortante, los modelos de lahar rara vez han abordado los procesos de erosión o sus cambios volumétricos (Carrivick, 2007; Carrivick et al, 2010; Fagents y Baloga, 2006; Lupiano et al, 2020), en gran parte debido a las deficiencias de los modelos utilizados. Estos modelos se ven obligados a iniciar el cálculo ya con el volumen de flujo total hipotético, ya sea en el edificio del volcán o en la topografía canalizada aguas abajo. Este enfoque no solo puede conducir a una predicción excesiva de la descarga del flujo y la altura de la etapa, sino que también significa que la evolución dinámica de los lahares sineruptivos, caracterizados por la aceleración inicial del flujo, el rápido crecimiento volumétrico en los flancos escarpados del volcán y una mayor erosión y deposición durante la propagación aguas abajo, no se puede rastrear en el proceso de modelado. Simulaciones anteriores para lahares de Cotopaxi utilizando modelos de enrutamiento de inundaciones 1D o 2D (Aguilera et al, 2004; Barberi et al, 1992; Toapaxi et al, 2019; Vera et al, 2019) demuestran el componente fluido-mecánico aplicable en los canales distantes pero no son aplicables al escurriemiento inicial en los flancos escarpados del volcán. Simultáneamente, Ordóñez et al. (2013) y Pistolesi et al. (2014) estimaron las áreas de inundación con LaharZ (Ordóñez et al, 2013; Pistolesi et al, 2014) pero no pudieron reproducir la formación inicial del lahar porque el modelo ignora los procesos de abultamiento dinámico en el flanco del volcán. Lupiano et al (2018) incluyeron procesos de erosión y derretimiento de glaciares de lahares sineruptivos de Cotopaxi con el uso de LLUNPIY, pero allí el modelo no está calibrado y validado a un evento histórico usando restricciones reproducibles.

1.4. Objetivos de este estudio

Este estudio tiene como objetivo evaluar sistemáticamente los peligros de los lahares sineruptivos masivos en el volcán Cotopaxi en el sistema de drenaje norte y sur densamente habitados, mediante el análisis retrospectivo del destructivo evento de lahar de 1877, restringiendo objetivamente las entradas del modelo que mejor se ajustan y emplear los resultados para posibles lahares futuros. Por lo tanto, (i) probamos la aplicabilidad de un modelo de Voellmy-Salm para lahares sineruptivos masivos, (ii) desarrollamos un modelo genérico que simula el lahar desde la iniciación en los flancos escarpados del volcán hasta los tramos distales, y (iii) incluimos la erosividad y arrastre en las laderas empinadas de un edificio volcánico y a lo largo de la trayectoria del modelo genérico.

2. Marco teórico

2.1. Entorno geológico e historia volcánica

El Cotopaxi es un gran estratovolcán ubicado en los Andes ecuatorianos de América del Sur, a ~60 km al sureste de Quito, capital de Ecuador (Figura 1). El edificio volcánico se encuentra en el cinturón metamórfico correspondiente a la Cordillera Real, pero también está intercalado parcialmente con formaciones volcánicas de la vecina Caldera Chalupas (Hall & Mothes, 2008). El casquete glaciar del Cotopaxi se extiende desde su cumbre (5897 msnm) por todos los flancos, alcanzando los ~4900 msnm en el flanco este y ~5300 msnm en el flanco occidental. Al igual que otros glaciares tropicales de América del Sur, el casquete glaciar del Cotopaxi ha retrocedido rápidamente por varias décadas (Cáceres, 2017; Vuille et al, 2018).

La evolución geológica de Cotopaxi comenzó con un antiguo complejo llamado Cotopaxi I entre 560 000 y 420 000 años AP, seguido de un largo período de inactividad (Hall & Mothes, 2008). El volcán retomó su actividad alrededor de 13 000 años AP con la formación del actual edificio Cotopaxi II, que se caracteriza por erupciones andesíticas y riolíticas. Alrededor de 4500 años AP, tuvo lugar la última erupción riolítica y produjo el colapso parcial del flanco norte del volcán, con la formación posterior de una avalancha de escombros y un gigantesco lahar cohesivo (Mothes et al, 1998).

Los últimos 2100 años del Cotopaxi II han sido testigos de erupciones andesíticas repetitivas caracterizadas por el deposito de flujos de lava, capas de tefra en la región y lahares sineruptivos a gran escala (Pistolesi et al, 2011, 2013). Se pueden estimar períodos de recurrencia promedio entre 117 y 147 años para las erupciones del Cotopaxi (Barberi et al, 1995), aunque se ha demostrado que la ocurrencia de erupciones muestra un agrupamiento significativo (Pistolesi et al, 2011).

Los lahares primarios formados durante erupciones explosivas sin duda representan el fenómeno más peligroso en el Cotopaxi (Mothes et al, 2016; Pistolesi et al, 2013; Sierra et al, 2019). Estos flujos se forman durante erupciones explosivas, cuando las corrientes de densidad piroclástica se mueven sobre el glaciar del volcán y lo derriten, produciendo instantáneamente grandes volúmenes de agua que se mezclan rápidamente con las rocas volcánicas (Pistolesi et al, 2013). Estos flujos transitan por los tres principales sistemas de drenaje que descienden del volcán (Figura 1): río Pita al norte, río Cutuchi al sur y río Tambo al este. Los tres drenajes se caracterizan por la alternancia de profundos cañones de pendiente pronunciada que canalizan los flujos con amplias áreas planas donde se produce un deposito significativo (Mothes et al, 2004). El río Pita y el río Cutuchi transitan por las zonas densamente pobladas de Quito y Latacunga, respectivamente (Mothes et al, 2016), mientras que el río Tambo llega a la cuenca amazónica escasamente poblada (Sierra et al, 2019).

La actividad volcánica más reciente en Cotopaxi se registró en abril de 2015 e identificada por su sistema de monitoreo (Hidalgo et al, 2018). Los registros anómalos aumentaron de manera constante hasta mediados de agosto, cuando comenzó la actividad hidromagmática, lo que produjo una emisión significativa de cenizas desde el cráter de la cumbre (Bernard et al, 2016; Gaunt et al, 2016). Esta actividad disminuyó progresivamente, menguando en diciembre de 2015 (Hidalgo et al, 2018) sin producir lahares sineruptivos.

2.2. El desastre del lahar de 1877

El último lahar sineruptivo del Cotopaxi ocurrió el 26 de junio de 1877. Según informes de Sodiro (1877) y Wolf (1878), el Cotopaxi había estrado en erupción desde principios de 1877, pero en junio su actividad aumentó significativamente. Ese día, el volcán produjo una fase explosiva sostenida con corrientes de densidad piroclástica generalizadas que descendieron por todos los flancos del edificio. Estimaciones de estos autores sugieren que la erupción explosiva duró al menos 15 min y que el frente de lahar tardó alrededor de 1 hora en llegar a Latacunga (al sur) y al Valle de Chillos (cerca de Quito al norte).

En el drenaje de Pita, el lahar de 1877 experimentó un desbordamiento significativo, desviando ~20% de su volumen hacia el Río Santa Clara (Mothes et al, 2004), que transita por áreas densamente pobladas del Valle de Chillos. Ocurrió un trasvase de flujo en el sitio llamado ‘La Caldera’, donde el Río Pita forma una graganta profunda y atraviesa una curva pronunciada, mientras que la cabecera del Río Santa Clara se encuentra justo arriba. Cálculos de Mothes et al (2004) asumen que ese efecto de trasvase se produjo por un desbordamiento resultante de la superelevación del flujo en un vértice forzado, sugiriendo que la tasa de descarga fue >40 000 m3/s, pero también asumen que el caudal se pudo represar parcialmente por el cañón muy angosto inmediatamente aguas abajo.

El depósito de lahares de 1877 aún se puede distinguir en el paisaje de los flancos proximales inferiores norte y nororiental del Cotopaxi: contiene un tipo específico de proyectil volcánico escoriáceo con la característica única de fracturarse y disgregarse fácilmente cuando se golpea ligeramente con un martillo (Mothes et al, 2004) y que no se repitió en los proyectiles volcánicos de los depósitos de lahares del siglo XVIII. Esto permite un mapeo detallado de sus depósitos, especialmente en las zonas próximales (<20 km desde la fuente), así como la medición de secciones transversales precisas a partir de las cuales se pueden obtener parámetros físicos como altura, velocidad y tasa de descarga (Aguilera et al, 2004; Mothes et al, 2004). Los volúmenes totales del lahar se estiman en alrededor de 60–75 × 10⁶ m³ (0,060 – 0,075 km³) en el drenaje norte y alrededor de 80–100 × 10⁶ m³ (0,080 – 0,100 km³) en el drenaje sur, respectivamente (Mothes et al, 2004).

El lahar de 1877 estuvo influenciado por la morfología de los canales transitados (Mothes et al, 2004). En los drenajes norte y sur, los flujos transitan inicialmente por profundos cañones del edificio volcánico con supuesto comportamiento erosivo. Una primera zona de deposición ocurre en regiones planas al pie del volcán. Posteriormente, se alternan zonas de erosión y/o aumento de energía con regiones planas. Las características granulométricas de los depósitos de 1877 varían entre flujos de escombros y flujos hiperconcentrados, con muestras proximales (<40 km desde la fuente) que muestran distribuciones de flujo de escombros granulares, mientras que las distales (hasta 70 km desde la fuente) muestran características más bien hiperconcentradas (Mothes et al, 2004; Pistolesi et al, 2013). Los depósitos del lahar de 1877 son granulométricamente similares a otros depósitos de las erupciones del Cotopaxi en el siglo XVIII y carecen de cantidades significativas de partículas de arcilla de tamaño de grano muy fino (Mothes et al, 2004) como es característico de los flujos no cohesivos.

3. Métodos

3.1. Modelo RAMMS

Aplicamos el modelo dinámico basado en la física RAMMS, que se usa ampliamente en el modelado de descentramiento de diferentes tipos de movimientos masivos, incluidas avalanchas de nieve, avalanchas de rocas, flujos de escombros e inundación por desborde violento de lago glaciar o “GLOFS” (Frey et al, 2018; Hussin et al, 2012; Schneider et al. al, 2010). El modelo de flujo de escombros RAMMS está formulado como un sistema monofásico que propaga la mezcla sólido-líquido como un flujo masivo incompresible con una densidad media constante (Bartelt et al, 2017; Christen et al, 2010). Requiere una calibración cuidadosa de los dos parámetros de fricción que rigen la reología de Voellmy-Salm y brinda información sobre parámetros clave para la evaluación de peligros, como la distancia de descentramiento, altura y velocidad del flujo, presión de impacto, erosión y deposición. El núcleo del programa es un esquema de solución numérica de segundo orden para resolver ecuaciones de movimiento promediadas en profundidad en un terreno tridimensional (x, y, z) por tiempo t (Christen et al, 2010). El balance de masa está dado por

\dot{Q}(x,y,t)=\partial_{t}H+\partial_{x}(HU_{x})+\partial_{y}(HU_{y})

(1)

donde \dot{Q}(x,y,t) describe el término de fuente de producción de masa, H(x, y, t) la altura de flujo y U(x, y, t) la velocidad de profundidad promediada. El equilibrio de momento de profundidad promediado en los ejes x y y está dado por

\partial_{t}(HU_{x})+\partial_{x}\left (HU^{2}_{x}+g_{z}k_{a/p}\frac{H^2}{2} \right )+\partial_{y}\left ( HU_{x}U_{y} \right )=S_{gx}-S_{fx}

(2)

y

\partial_{t}(HU_{y})+\partial_{y}\left (HU^{2}_{y}+g_{z}k_{a/p}\frac{H^2}{2} \right )+\partial_{x}\left ( HU_{x}U_{y} \right )=S_{gy}-S_{fy}

(3)

respectivamente, donde g = (g_x, g_y, g_z) describe el vector gravitacional, k_a/p el coeficiente de presión terrestre (normalmente ajustado a 1), S_g y S_f denotan la aceleración gravitacional y desaceleración por fricción paralelas a la pendiente, respectivamente.

La resistencia friccional S_f se describe utilizando un enfoque de Voellmy-Salm (Salm, 1993; Voellmy, 1955), que incorpora una fricción seca de Coulomb \mu proporcional a la tensión normal y una velocidad de arrastre al cuadrado o fricción viscoso-turbulenta \xi (similar a Chézy fricción) (Bartelt et al, 1999; Christen et al, 2010; Fischer et al, 2012):

S_{f}=\mu N+\frac{\rho gU^{2}}{\xi } con N=\rho H\left [g\cos \left ( \varphi \right )+\kappa U^{2}\right ]

(4)

donde \rho es la densidad de la masa, \varphi el ángulo de inclinación de la superficie y \kappa la curvatura del terreno en la dirección del flujo. La fricción de Coulomb \mu representa la resistencia de la fase sólida, es independiente de la velocidad y gobierna el comportamiento del flujo cuando se mueve lentamente, mientras que la fricción turbulenta \xi representa la fase fluida, escala con la velocidad del flujo y domina a velocidades más altas en las zonas de aceleración. (Bartelt et al, 2017; Christen et al, 2010).

Un módulo de arrastre implementado en RAMMS predice la profundidad de la erosión en términos de disminución neta en la elevación del lecho del canal (Frank et al, 2017). Simula el aumento en el volumen del flujo debido al arrastre (Frank et al, 2017), pero el arrastre de sedimentos no está acoplado a cambios topográficos en el modelo digital de elevación (DEM) subyacente o cambios en el comportamiento del flujo. Las observaciones de campo en la cuenca excepcionalmente activa de Illgraben en Suiza sobre las profundidades de erosión (Schürch et al, 2011) y tasas de erosión (Berger et al, 2011) de los flujos de escombros naturales sirvieron como base para el desarrollo del algoritmo de arrastre. (Frank et al, 2015). La profundidad de erosión potencial em por celda de cuadrícula se describe como una función del esfuerzo cortante calculado que actúa sobre el lecho del canal, dado por

\tau =\rho gH \sin \varphi

(5)

y está relacionada con un esfuerzo cortante crítico \tau_{c} τc y un factor de proporcionalidad \frac{d_{z}}{d_{\tau}} dz dτ que controla la tasa de erosión vertical, dada por

e_{m}=\left\{\begin{matrix}0&para&\tau<\tau_{c}\\\frac{dz}{d\tau}\left(\tau-\tau_{c}\right)&para&\tau\geq\tau_{c}\end{matrix}\right.

(6)

Con base en las tasas de erosión promedio registradas en Illgraben (Berger et al, 2011) se utiliza una tasa de erosión específica \frac{dz}{dt} (Frank et al, 2015):

\frac{dz}{dt}=-0.025 para e_{t}\leq e_{m}

(7)

El arrastre de sedimentos inicia una vez que se excede el esfuerzo cortante crítico \tau_{c} y finaliza cuando la profundidad de erosión real e_{t} alcanza la profundidad de erosión potencial máxima e_{m}. Las variaciones específicas de los parámetros de erosión se describen en la sección de calibración del modelo más abajo y mayores detalles del modelo RAMMS y el modelo de arrastre implementado están disponibles en Christen et al (2010) y Frank et al (2015, 2017).

3.2. Configuración del modelo

Nuestro objetivo es configurar el modelo de lahares RAMMS de manera genérica iniciando el flujo en el flanco del volcán cerca del cráter, permitiendo la acumulación de sedimentos a lo largo del camino y propagando el lahar hasta tramos distales. La configuración está diseñada con el propósito de calcular retroactivamente el evento del lahar de 1877 y se puede aplicar al modelado directo de posibles lahares futuros.

3.2.1. Premisas del modelo: DEM y límites del dominio del modelo

Los lahares de Cotopaxi generalmente no muestran límites de depósito coherentes, pero pueden trasladarse hasta 325 km, como se informó para el lahar de 1877 (Mothes et al, 2004), distalmente en una fase acuosa. Por lo tanto, la región de modelado está restringida por los centros urbanos más grandes alrededor del volcán: al norte, la región modelo se extiende hasta aproximadamente 70 km aguas abajo del volcán, cubriendo la aglomeración urbana del Valle de Los Chillos; al sur termina al sur de la ciudad de Salcedo unos 65 km aguas abajo del Cotopaxi (ver Figura 2). Dada la gran extensión, una resolución espacial de 10 m se considera suficiente para capturar patrones de flujo clave. Todas las simulaciones numéricas realizadas en este estudio se basan en un DEM de 10 m, que volvimos a muestrear a partir de un DEM original de 4 m de resolución desarrollado por el proyecto SIGTIERRAS en 2010 (Ministerio de Agricultura, Ganadería, Acuacultura y Pesca). No se dispone de información detallada sobre la morfología del valle antes del paso del lahar de 1877. Sin embargo, podemos suponer que la topografía actual de los principales canales excavados parcialmente profundos no difiere significativamente ya que la evidencia histórica y estratigráfica sugiere que confinaron todos los lahares históricos conocidos desde el siglo XVIII (Hall & Mothes, 2008).

3.2.2. Premisas modelo: Erosión y arrastre

Somos conscientes de las altas incertidumbres involucradas con el modelado de la erosividad de los lahares, ya que no existe información cuantitativa sobre el potencial de erosión de los lahares históricos en Cotopaxi. Además, al ignorar el volumen y la reducción dinámicos, sería necesario comenzar el modelo con el volumen total del lahar, lo que puede llevar a una sobreestimación de la descarga máxima y la altura del flujo (Frank et al, 2015). Por lo tanto, nuestro enfoque para simular el arrastre de sedimentos a lo largo de la trayectoria (i) se enfoca en modelar el balance de masa general y reproducir los volúmenes de flujo total promedio informados por Mothes et al. (2004) y (ii) se basa en una caracterización y mapeo sencillos de la erosionabilidad del sustrato. Usamos ortofotografías tomadas en 2011, 2014 y 2016, proporcionadas como WMS por IGM Ecuador, mapas topográficos y geológicos (Hall y Mothes, 2008; Vezzoli et al, 2017) para una caracterización directa de la erosionabilidad del sustrato relacionada con la litología (lecho rocoso/suelo) y consolidación, que verificamos en puntos críticos en el campo.

Se pueden definir tres clases de erosión en función de la erosionabilidad del sustrato. Asignamos (1) la configuración de arrastre predeterminada de RAMMS (valores promedio derivados de mediciones en la cuenca de Illgraben) a lo largo de la mayoría de los canales afectados por el lahar en los valles, donde prevalecen los depósitos de numerosos lahares históricos y flujos de escoria. Además, discriminamos entre (2) áreas con erosionabilidad comparativamente más alta ubicadas en el flanco del volcán, donde están presentes abundantes sedimentos piroclásticos sueltos o depósitos glaciares expuestos por el retroceso reciente de los glaciares; y (3) áreas con erosionabilidad limitada caracterizadas por lecho rocoso cubierto de regolito o profundos canales flanqueados por lecho rocoso.

Para cada clase, especificamos los parámetros de arrastre RAMMS que controlan el grado de arrastre de sedimentos en el flujo (Tabla 1) a partir de un conjunto de valores dados (Bartelt et al, 2017; Frank et al, 2015, 2017). En áreas con “erosionabilidad limitada”, la profundidad máxima de erosión está restringida a 0,5 m modulada por el espesor supuesto del regolito, mientras que las áreas con “erosionabilidad media” y “alta” se consideran de suministro ilimitado.

3.2.3. Estimación del volumen inicial e hidrograma de entrada

Las estimaciones del volumen y la descarga máxima de los flujos másicos son información de referencia crítica para el modelado. Para retrocalcular el lahar de 1877, reproducimos el promedio de los volúmenes totales informados por Mothes et al (2004), es decir, ~68 × 10⁶ m³ (0,068 km³) para el drenaje norte y ~90 × 10⁶ m³ (0,090 km³) para el drenaje sur, respectivamente). Los volúmenes de flujo deseados se logran mediante una aproximación gradual que equilibra la suma de los volúmenes de ingreso y erosionados. Probamos iterativamente volúmenes de ingreso entre 25 × 10⁶ y 58 × 10⁶ m³ (0,025 y 0,058 km³) y analizamos el crecimiento volumétrico dinámico resultante del lahar. Dentro de esta serie de modelos, los parámetros de fricción se mantienen constantes en \mu = 0,01 y \xi = 1800 m/s² (tomados del núcleo de la matriz de sensibilidad), pero esperamos que el cambio especial de \mu pueda afectar los volúmenes de erosión modelados, según lo informado por Frank et al (2015). En comparación con los volúmenes totales de los flujos simulados, esto puede tener solo un efecto menor.

Para una aproximación basta del hidrograma de descarga utilizamos un hidrograma de descarga de tres puntos idealizado. Los parámetros de forma se basan en las siguientes observaciones y relación empírica:

1. Descarga máxima. Mizuyama et al (1992) propusieron una fórmula empírica para estimar la descarga máxima Qp de los flujos de escombros granulares en función del volumen V, dada por

Q_{p}=0.135V^{0.78}

(8)

Pierson (1998) probó la aplicabilidad de esta correlación para lahares de gran magnitud y encontró que sirve como una aproximación adecuada para las descargas máximas cercanas a la fuente. La fórmula se aplica a los volúmenes de lahar iniciales comprobados que vayan desde 25 × 10⁶ a 58 × 10⁶ m³.

2. Momento de máxima descarga. Lahares bien estudiados en los volcanes cubiertos de glaciares del Monte St. Helens y Nevado del Ruiz mostraron una liberación de flujo relativamente instantánea (Pierson, 1995; Pierson & Scott, 1985; Pierson et al, 1990). Por lo tanto, el tiempo de descarga máxima se establece en 120 s después del inicio.

3. Duración. Podemos considerar de manera conservadora los relatos de testigos oculares sobre la duración del flujo piroclástico como un indicador de la duración de la liberación del lahar. Wolf (1878) informó que los flujos piroclásticos duraron de 15 a 20 minutos y provocaron el derretimiento inmediato de la nieve y el hielo en todo el cráter. La aplicación de la correlación volumen-descarga descrita anteriormente por Mizuyama et al (1992) a un hidrograma de forma triangular arroja duraciones de liberación de alrededor de 12 min, lo que concuerda en gran medida con las observaciones de testigos presenciales.

Se coloca un hidrograma amplio de liberación en el límite aguas arriba de cada uno de los dos sistemas fluviales, cerca de la ubicación actual de las bocas de los glaciares. Esto permite que el flujo inyectado se drene completamente en un sistema y todos los canales principales del sistema de drenaje respectivo se ven afectados simultáneamente.

3.3. Calibración de modelos y análisis de sensibilidad

Calculamos retrospectivamente el bien documentado evento del lahar de 1877 en el sistema de drenaje del norte de Cotopaxi a fin de evaluar el desempeño del modelo y derivar valores plausibles de los parámetros del modelo para futuros escenarios de lahar. Toda la información disponible sobre el evento de referencia se basa en informes de testigos presenciales de Sodiro (1877) y Wolf (1878), mapas oficiales de peligro proporcionados por IG-EPN de Ecuador (Mothes et al, 2016) e investigaciones sedimentológicas de campo de Mothes et al (2004). En total, se seleccionan 14 criterios de calibración (Figura 3) que abarcan diferentes tipos de parámetros reológicos y una distribución espacial representativa a lo largo de la trayectoria del lahar. Diferenciamos entre dos categorías de criterios de calibración sobre la base de la fiabilidad de los datos:

• Los ‘criterios primarios’ incluyen aquellos que debe tener cualquier modelo realista de lahar que retrocalcule el lahar de 1877: (i) el lahar fluyó sobre la pared del canal norte en La Caldera (correspondiente a una altura de flujo > 42 m dada por el DEM utilizado), (ii) en consecuencia inundó el pueblo de Sangolquí; y (iii) también llegó a Tumbaco aguas abajo en un tiempo razonable.
• Los ‘criterios adicionales’ comprenden datos reológicos menos confiables, como observaciones de testigos presenciales sobre los tiempos de llegada a lugares aguas abajo y estimaciones sobre velocidad del flujo, descarga máxima, alturas de las etapas y límites de inundación derivados de investigaciones de campo más de 100 años después del evento. El mapa de riesgo oficial del IG-EPN (Mothes et al, 2016) delimita áreas potenciales de inundación de lahares provocados por erupciones moderadas a fuertes, con base en el mapeo de campo de depósitos del lahar de 1877 y otros lahares del siglo XVIII.

Los criterios de calibración actúan como restricciones de simulación y se utilizan para determinar las combinaciones de parámetros de entrada que mejor se ajustan. Debido a que los parámetros de entrada pueden causar un cambio significativo en el comportamiento del modelo, analizamos sistemáticamente la sensibilidad de los dos coeficientes de fricción de Voellmy-Salm \mu y \xi, la densidad del flujo y la velocidad inicial para modelar los resultados de la distancia de descentramiento, descarga máxima, tiempo de viaje, velocidad de flujo y altura de flujo (representadas por los criterios de calibración).

Dada la ubicación de los hidrogramas de liberación al menos 1000 m pendiente abajo del cráter, lo más probable es que el lahar ya haya comenzado a evolucionar debido al derretimiento del hielo del glaciar y la mezcla con productos volcánicos recién expulsados ​​al “llegar” a la ubicación del hidrograma de liberación. Parece probable que el flujo inicial ya haya acelerado a cierta velocidad, que se aproxima explorando la velocidad inicial entre 10 y 20 m/s en intervalos de 5 m/s. La densidad del flujo modelado varía en incrementos de 200 kg/m³ de 1400 a 2200 kg/m³, lo que representa densidades conocidas de flujos hiperconcentrados y flujos de escombros. Según los coeficientes de fricción de Voellmy calibrados para flujos de escombros y GLOF en estudios previos (Frey et al, 2018; McDougall, 2006; Revellino et al, 2004; Schneider et al, 2014; Schraml et al, 2015), el parámetro espacio puede abarcar un amplio rango de 0,01 a 0,2 para valores de \mu y de 200 a 1500 m/s² para valores de \xi, con una tendencia a \mu más bajos y \xi más altos para flujos de gran magnitud. Aquí, variamos el coeficiente \mu de 0,0025 a 0,05 con intervalos de 0,0025 entre 0,0025 y 0,01 y de 0,005 entre 0,01 y 0,05. El valor de \xi varía entre 600 y 3000 m/s² con incrementos de 400 m/s².

3.4. Evaluación del rendimiento del modelo

Para evaluar objetivamente la calidad de las predicciones del modelo con RAMMS, comparamos por separado los resultados del modelo (m_i) con los valores objetivo (t_i) para las medidas de calibración (i) de la velocidad del flujo, la altura del flujo, el tiempo de viaje y la descarga máxima. El rendimiento de una simulación (\phi) utilizando un conjunto dado de parámetros de entrada (y) se cuantifica como desviación porcentual, con

\phi_{i}(y)=\frac{m_{i}(y)-t_{i}}{t_{i}}\cdot100

(9)

La normalización de las diferencias entre los valores modelados y objetivos facilita una comparación sencilla entre las diferentes medidas reológicas de los criterios de calibración.

La similitud entre áreas simulada y mapeada de impacto del lahar (Figura 3) se expresa como superposición porcentual, dividiendo el número de celdas de cuadrícula modeladas correctamente como inundadas (verdaderos positivos) por el número total de celdas modeladas (verdaderos positivos + falsos negativos). No parece significativo incluir falsos positivos en este caso porque el área mapeada de impacto del lahar abarca los límites de inundación del lahar de 1877 pero también de varios otros lahares históricos. El cálculo se aplica a la extensión compartida de ambas capas (es decir, despreciando la zona proximal alrededor del volcán) y considera áreas de impacto simuladas inundadas en más de 0,05 m. Para armonizar los resultados del modelo, los hallazgos se representan más como desviación del modelo (100 – porcentaje de superposición).

Para cada conjunto de parámetros de entrada (y), se suman los valores de los análisis de rendimiento de restricción para dar un único valor de rendimiento de simulación (\phi_{total}), con

\phi_{total}(y)=\frac{\sum^{n}_{i=1}\left|\left(\frac{m_{i}(y)-t_{i}}{t_{i}}\cdot100\right)\right|}{n}

(10)

Este número resalta las desviaciones mínimas del modelo dentro del espacio de parámetros de entrada analizados y por lo tanto, nos permite seleccionar las combinaciones de parámetros más representativas para el evento retrocalculado.

3.5. Simulación de futuros lahares potenciales

Con base en investigaciones de la historia eruptiva del Cotopaxi (Hall & Mothes, 2008; Mothes et al, 1998, 2004), Andrade et al (2005) propusieron cuatro posibles escenarios futuros de erupción (Figura 9a), que pueden vincularse con el tamaño de la erupción expresada en el Índice de Explosividad Volcánica (VEI por Volcanic Explosivity Index) (Newhall & Self, 1982): la erupción esperada más débil (VEI 1–2, N. del T.: Escenario 1) se caracteriza por actividad estromboliana, emisiones de tefra de pequeñas a moderadas y eyección de proyectiles balísticos, mientras que durante las erupciones de VEI 2–3 (Escenario 2), VEI 3–4 (Escenario 3) y VEI > 4 (Escenario 4) con actividad vulcaniana, subpliniana y pliniana, respectivamente, la erupción probablemente esté acompañada de flujos piroclásticos. Debido a que los flujos piroclásticos son agentes necesarios para derretir cantidades significativas de nieve y hielo requeridas para la generación de lahares sineruptivos, el escenario hipotético de lahares 1 se excluye aquí del modelo de descentramiento. A modo de comparación, durante la última erupción en Cotopaxi en 2015, clasificada como VEI 1–2 (Escenario 1) (Mothes et al, 2017), no se produjo ningún lahar sineruptivo, sino varios lahares secundarios de pequeña escala dentro de la zona proximal. Es probable que un aumento en el VEI corresponda a mayores volúmenes y dispersión de tefra emitida y al aumento de la energía y el impacto espacial de los flujos piroclásticos (Andrade et al, 2005). Un escenario VEI 3-4 (Escenario 3), equivalente a la erupción de 1877, puede caracterizarse por una actividad de ‘rebosamiento’ que produce flujos de escoria, mientras que durante erupciones plinianas altamente explosivas VEI > 4 (Escenario 4), un colapso de la columna de erupción podría afectar a todo el glaciar de la cumbre (Pistolesi et al, 2013). Por lo tanto, estudios previos han concluido que un mayor grado de interacciones volcán-glaciar con el aumento de VEI da como resultado lahares potenciales más voluminosos y destructivos (Mothes et al, 2004; Ordóñez et al, 2013; Pistolesi et al, 2013).

Un enfoque común, pero simplificado, para estimar los volúmenes de lahar del Cotopaxi (Barberi et al, 1992; Mothes et al, 2004; Ordóñez et al, 2013; Toapaxi et al, 2019; Vera et al, 2019) se basa en tres controles: los dos controles dinámicos incluyen (i) espesores variables de hielo y nieve derretidos uniformemente desde la superficie del glaciar durante una erupción como indicador del grado de interacción volcán-glaciar; y (ii) el tamaño de la cumbre del glaciar. Un tercer control estático para la estimación del volumen total del lahar es (iii) la suposición de un crecimiento del triple del volumen del lahar debido a la incorporación de sólidos. La escala del espesor del hielo derretido en cada escenario se ha realizado con un factor de aproximadamente 2, asumiendo un derretimiento uniforme de 2, 4–5 y 8–10 m para una erupción tipo VEI 2–3, VEI 3–4 y VEI > 4, respectivamente. Con respecto a la extensión de los glaciares del Cotopaxi, Cáceres (2017) informó que la superficie ha disminuido en promedio alrededor de 52% entre los años 1977 y 2016. Con el cálculo de los volúmenes de lahar basándose en el espesor de fusión y la extensión del glaciar, un retroceso de los glaciares de alrededor del 50% se traduciría directamente en una reducción del 50% en los volúmenes iniciales de posibles lahares futuros en comparación con la situación durante la erupción de 1877 (Toapaxi et al, 2019; Vera et al, 2019).

Adaptamos algunos de estos supuestos simplificados a la cuantificación de futuros escenarios de lahares en este estudio utilizando el siguiente enfoque:

1. El modelo de lahar calibrado de la erupción de 1877 con descarga de liberación estimada y volumen de entrada (que define el hidrograma de entrada) proporciona la base para el modelado de escenarios. El evento modelado fue equivalente a una erupción VEI 3–4 (N. del T.: Escenario 3).
2. Debido a una reducción del 50% en el área de la superficie del glaciar, se espera que una futura erupción VEI 3–4 produzca lahares con descargas de liberación de la mitad de las del lahar de 1877 (suponiendo que las áreas glaciares totales en 1877 y 1977 fueran de la misma extensión). Esto representa claramente una suposición simplificada; el grado en que el retroceso de los glaciares afecta las futuras magnitudes de los lahares también puede necesitar una mirada más cercana a la dinámica de las erupciones y los mecanismos desencadenantes (Pistolesi et al, 2013, 2014) como se analiza más adelante.
3. Se utiliza un factor de escala de 2 para expresar el grado creciente de interacciones volcán-glaciar con el aumento de VEI. Partiendo de un Escenario 3 (VEI 3–4), reducimos a la mitad las descargas de liberación para futuros Escenario 2 (eventos VEI 2–3) y las duplicamos para futuros Escenarios 4 (VEI > 4). Los volúmenes de entrada correspondientes se estiman utilizando la relación volumen-descarga de Mizuyama et al (1992) (Ecuación 8).
4. El lahar inicial puede cambiar dinámicamente su volumen debido al arrastre y deposición de material. La erosión del flujo está controlada por la erosionabilidad del sustrato.

4. Resultados

4.1. Estimación del volumen de entrada, volumen de arrastre e hidrograma de entrada para el lahar de 1877

El análisis de los volúmenes de entrada del modelo para el cálculo retroactivo del lahar de 1877 se inició en el drenaje norte y sur con un volumen inicial estimado de 25 × 10⁶ m³, pero como los volúmenes de flujo total modelados eran demasiado bajos, los volúmenes de entrada se elevaron a 34 × 10⁶ m³ y sucesivamente aumentado en incrementos de 3 × 10⁶ m³. Los volúmenes de flujo total objetivo (~68 × 10⁶ m³ en el drenaje norte y ~90 × 10⁶ m³ en el drenaje sur) están bien aproximados con un volumen de entrada de 46 × 10⁶ m³ en el drenaje norte y 49 × 10⁶ m³ en el drenaje sur, respectivamente (Figura 4). Al aplicar estos valores en la ecuación volumen-descarga dada por Mizuyama et al (1992) (Ecuación 8), obtenemos descargas máximas de 128 000 m/s³ en el drenaje norte y 134 500 m/s³ en el drenaje sur. De esta forma, se estiman todos los parámetros de forma que definen el hidrograma de liberación de tres puntos.

Dadas las configuraciones de arrastre (Figura 2 y Tabla 1), un flujo con un volumen inicial de 46 × 10⁶ m³ erosiona y arrastra alrededor de 23 × 10⁶ m³ a lo largo de su desagüe hacia el norte, mientras que en el drenaje sur, un flujo inicial de 49 × 10⁶ m³ crece por 40 × 10⁶ m³ debido al arrastre.

Esto corresponde a un crecimiento de volumen de 50 y 80% en el sistema de drenaje norte y sur, respectivamente. En ambos sistemas fluviales podemos observar que la cantidad de material erosionado varía con el volumen inicial, lo que revela una tendencia general hacia un menor crecimiento del volumen atribuido a la erosión a medida que aumenta el volumen de entrada. Esta tendencia es más pronunciada en el drenaje sur, donde la erosión y el arrastre dan como resultado un crecimiento del volumen que va desde el 135% para un volumen inicial de 25 × 10⁶ m³ hasta el 70% para un volumen inicial de 58 × 10⁶ m³.

4.2. Calibración de modelos y análisis de sensibilidad

Dentro de una serie de ejecuciones piloto, encontramos que la densidad del flujo y la velocidad inicial solo afectan levemente la propagación del flujo. Las modificaciones de la densidad de flujo dan como resultado resultados del modelo para tiempos de llegada (Sangolqui), alturas de flujo (La Caldera) y caudales máximos (Proaño) que difieren en 5, 2 y 15%, respectivamente. Un cambio en la velocidad inicial del flujo hace que los tiempos de llegada (Sangolqui), las alturas de flujo (La Caldera) y las descargas máximas (Proaño) difieren en 1, 0.4 y 4%, respectivamente. Con base en estas observaciones mantuvimos la densidad del flujo en un valor constante de 2000 kg/m³ y comenzamos todas las simulaciones con una velocidad de 10 m/s. Por lo tanto, el análisis de sensibilidad se limitó a los parámetros de fricción gobernantes \mu y \xi de la reología de Voellmy-Salm.

Primero, analizamos si los resultados del modelo cumplen con los “criterios de calibración primarios” requeridos. En simulaciones usando \mu > 0.02, el lahar no experimenta el desbordamiento observado en La Caldera hacia la cabecera del Santa Clara y al utilizar \mu > 0.015, el lahar no llega a Tumbaco en un tiempo de cálculo razonable de 4 h. En consecuencia, pudimos restringir en gran medida el espacio de parámetros para el análisis de sensibilidad a coeficientes \mu que oscilan entre 0,0025 y 0,02 y coeficientes \xi entre 600 y 3000 m/s², lo que da como resultado un total de 42 combinaciones \mu/\xi diferentes. Aunque las ejecuciones del modelo con \mu = 0.02 no cumplen con los requisitos de todas las “restricciones del modelo principal”, este valor aún se incluye para obtener una comprensión más amplia del comportamiento del modelo.

En segundo lugar, utilizamos los “criterios de calibración adicionales” para analizar la sensibilidad del modelo y encontrar las combinaciones de parámetros más adecuadas. Los resultados de las 42 ejecuciones del modelo para cada criterio se resumen en las Figuras 6, 7 y 8 y se presentan de cuatro maneras diferentes, mostrando (a) la sensibilidad del resultado del modelo en función de \mu y \xi, (b) una estadística frecuencia de agrupación de análisis de los resultados del modelo, (c) una matriz visual de los resultados del modelo, y (d) una matriz visual que representa el rendimiento de una simulación (\phi) en términos de desviación del modelo del valor objetivo. Nos referimos a la Figura 3 para la localización de las restricciones de calibración posteriores a lo largo del drenaje norte. Los resultados del modelo para la velocidad muestran que esta medida es sensible a los coeficientes \mu y \xi, pero la sensibilidad al coeficiente de fricción \mu aumenta en las ubicaciones aguas abajo de Proaño y Hacienda San Rafael. En general, las simulaciones dieron velocidades que se comparan favorablemente con las tres velocidades objetivo. Los resultados del modelo para la altura del flujo demuestran que esta medida de calibración es muy sensible a las variaciones de \mu y en gran medida insensible a las variaciones de \xi. Sin embargo, todos los resultados del modelo en las dos ubicaciones distales, Hacienda Valencia y San Rafael (Megamaxi) claramente sobrestiman los valores objetivo, lo que puede atribuirse a grandes incertidumbres en el valor objetivo mismo. En San Rafael (Megamaxi), las desviaciones son del orden del 25%, mientras que en la Hacienda Valencia, las alturas de flujo simuladas superan los valores objetivo en un 130–170%. Como se indica en la distribución de frecuencias (Figura 7, columna b), la mayoría de las ejecuciones del modelo (60%) predicen alturas máximas de flujo en la Hacienda Valencia entre 35 y 37 m en comparación con el valor objetivo de 15 m. Usando los tiempos de viaje como restricciones de calibración, podemos observar una alta sensibilidad de esta medida a \mu y una influencia menor de \xi. Los resultados muestran además que los tiempos de viaje simulados coinciden en gran medida con los informes históricos de testigos oculares. Los resultados de la simulación para las medidas de descarga máxima reaccionan de manera sensible a las variaciones en \mu, pero son en gran parte insensibles a las variaciones en \xi. Destacamos de manera particular que la simulación del desbordamiento en la cabecera del río Santa Clara es una restricción crucial para la evaluación del desempeño del modelo. Observamos que todas las ejecuciones del modelo reprodujeron un desbordamiento en ‘La Caldera’, aunque aproximadamente el 90% de ellas subestiman la descarga máxima prevista. Curiosamente, podemos observar particularmente bien en las animaciones de nuestro modelo que la bifurcación del flujo en La Caldera puede no ser el resultado de un avance real del frente de flujo contra la pared del canal norte casi vertical. El lahar es capaz de pasar la curva pronunciada en La Caldera, pero a medida que el cañón Pita aguas abajo se estrecha significativamente, el flujo se ralentiza, lo que conduce a un relleno progresivo aguas arriba del cañón hasta superar la altura de la pared del canal norte después de unos 90 s. Los resultados del modelo que representan áreas de impacto simuladas (altura máxima de flujo > 0,05 m) se corresponden bien con las áreas de impacto mapeadas, con una buena similitud general entre las ejecuciones del modelo que oscila entre el 79 y el 86%.

Como era de esperar, las combinaciones de \mu/\xi que mejor se ajustan para cada medida de calibración con menos desviaciones de los valores objetivo se agrupan en ubicaciones ligeramente diferentes dentro de la matriz visual (Figuras 6, 7 y 8, columna d). Sin introducir una ponderación subjetiva para criterios de calibración específicos o tipos de medidas de calibración, consideramos igualmente todos los criterios para derivar un valor final del rendimiento total del modelo (Ecuación 10). Nuestro análisis muestra que las ejecuciones del modelo con \mu = 0,005/\xi = 1400 m/s² y \mu = 0,0025/\xi = 1000 m/s² produjeron las simulaciones generales más representativas del lahar de 1877, reveladas por un valor de rendimiento cuantificado como desviación total del modelo de 23% para ambas combinaciones (Figura 8, recuadro rojo). Es importante destacar que los resultados demuestran que los parámetros de mejor ajuste caen en una amplia zona de valores bajos de desviación.

Las simulaciones de la altura del flujo, velocidad y erosión del lahar de un modelo ejecutado con \mu = 0,005 y \xi = 1400 m/s² se ilustran en la Figura 5 y el rendimiento del modelo para esta combinación de parámetros de mejor ajuste se puede resumir en las Figuras 6, 7 y 8 de la siguiente manera. (i) El área de impacto está bien reproducida (desviación del modelo del 16%) y traza las rutas de flujo principales, lo que indica que el modelo captura características clave como la bifurcación y la confluencia del flujo. En el área norte de baja pendiente entre las quebradas El Salto y Río Pita se sobreestima la expansión lateral del lahar, lo que contribuye a las desviaciones observadas del modelo. (ii) Las velocidades de flujo están muy bien reproducidas con desviaciones de 14, −1 y 3% respectivamente en Jatabamba, Proaño y Hacienda San Rafael. (iii) Mientras que las alturas de caudal simuladas en Hacienda San Rafael y San Rafael (Megamaxi) muestran desviaciones de los valores objetivo de 25 y 21%, las alturas de caudal en Hacienda Valencia están sobreestimadas en un 138%, posiblemente debido a las dificultades para delimitar los límites de caudal en tramos distales. (iv) Los tiempos de viaje a Sangolquí y Tumbaco caen muy bien en el rango de los tiempos de llegada informados con una desviación de 0 y 10%, respectivamente. (v) El caudal máximo en Proaño está ligeramente sobreestimado por el modelo con un valor del 10%. El desbordamiento al río Santa Clara en La Caldera se reproduce con una desviación del -20%.

4.3. Simulación de lahares futuros potenciales

Aplicamos el modelo de lahar calibrado con parámetros de entrada de mejor ajuste derivados previamente (\mu = 0.005 y \xi = 1400 m/s²) y estimaciones del volumen inicial y la descarga del lahar de 1877 a posibles lahares futuros durante tres diferentes escenarios de erupción del Cotopaxi. Nuestro enfoque para evaluar los futuros volúmenes de lahar considerando el efecto del retroceso de los glaciares arrojó las siguientes estimaciones hidrográficas (Figura 9): para los Escenarios 2, 3 y 4 en el drenaje norte, calculamos descargas máximas de 32 000, 64 000 y 128 000 m³/s, respectivamente y volúmenes de entrada correspondientes de 7,8, 18,9 y 46 × 10⁶ m³. En el drenaje sur, derivamos para los escenarios 2, 3 y 4 descargas máximas de 33 600, 67 300 y 134 500 m³/s, respectivamente y volúmenes de entrada correspondientes de 8,3, 20,2 y 49 × 10⁶ m³.

Los resultados de las simulaciones de escenarios se resumen en la Tabla 2 y la Figura 10; se proporcionan resultados de simulación más detallados de alturas de flujo, velocidad de flujo y erosión en la información de apoyo en línea. Observamos que los lahares simulados en el drenaje norte alcanzan volúmenes de flujo total de ~20 × 10⁶ m³ (~0,020 km³) en el Escenario 2, ~35 × 10⁶ m³ (~0,035 km³) en el Escenario 3 y ~70 × 10⁶ m³ (~0,070 km³) en el Escenario 4, lo que corresponde a un crecimiento de volumen debido a arrastre de 150, 80 y 50%, respectivamente. Los lahares potenciales en los tres escenarios de erupción transitan a través de la aglomeración urbana densamente poblada del Valle de Los Chillos (incluido el sector de San Rafael) a través del Río Pita, luego se canalizan en el Río San Pedro y llegan a la ciudad de Tumbaco y otros lugares aguas abajo (fuera de la región de modelado). Los tiempos de viaje simulados a San Rafael oscilan entre 84 min (Escenario 2) y 60 min (Escenario 4). Un desbordamiento en La Caldera y entrada al río Santa Clara ocurre en grado menor, en el Escenario 3, pero solo en el Escenario 4 la descarga de desbordamiento es lo suficientemente grande como para que el lahar pueda propagarse al pueblo de Sangolquí en 52 minutos y luego unirse con el lahar en el canal del Río Pita.

Los futuros lahares que descienden por los canales del sur alcanzan volúmenes de flujo total de ~40 × 10⁶ m³ (~0,040 km³) en el Escenario 2, ~55 × 10⁶ m³ (~0,055 km³) en el Escenario 3 y ~90 × 10⁶ m3 (~0,090 km³) en el Escenario 4, porque los volúmenes de lahar aumentan respectivamente en 400, 180 y 90% debido a la incorporación de sedimentos. El crecimiento del volumen del lahar, especialmente durante las erupciones VEI 2–3 (Escenario 2) y VEI 3–4 (Escenario 3), es notablemente mayor en el drenaje sur. Los lahares en los tres escenarios llegan a la ciudad de Salcedo y aguas abajo, pero el ancho de la inundación en el valle interandino densamente poblado varía significativamente en la sección norte de Latacunga, en la ciudad de Latacunga y en la ciudad de Salcedo. Los tiempos de llegada simulados de los lahares del sur son en gran medida insensibles a la magnitud de la liberación, ya que difieren muy poco entre los tres escenarios de erupción. Los tiempos de llegada a los pueblos y aldeas ubicados cerca del volcán, como Mulaló y Chasqui, varían menos (solo 1–2 min), lo que indica que el relieve abrupto del edificio del volcán favorece una fuerte aceleración de los lahares provocados por pequeños VEI 2–3 así como por grandes erupciones VEI> 4 (consulte la información de apoyo en línea sobre velocidades de flujo para comparar). Además, los tiempos de tránsito simulados a Latacunga coinciden con los informados por testigos presenciales del evento de 1877 (<1 h).

5. Discusión

Este estudio simuló un flujo completo desde el flanco superior del volcán Cotopaxi hasta las áreas de impacto distal utilizando el software 2D de flujo de escombros RAMMS por primera vez. Esto implica investigar la aplicabilidad de una reología de Voellmy-Salm para fines de simulación de lahar, lo que logramos mediante la calibración del modelo y la evaluación del rendimiento. Una parametrización realista del comportamiento del flujo implica necesariamente estimar la capacidad de erosión de los lahares y analizar la influencia de la erosión y el arrastre en los volúmenes de flujo.

5.1. Estimación de la capacidad de erosión

Para estimar la erosión y el arrastre por lahares, utilizamos un enfoque de modelado directo, que se basa en estimar la erosionabilidad del sustrato y es independiente de la magnitud de un evento eruptivo. Demostramos que los lahares de Cotopaxi pueden experimentar una acumulación sustancial durante la propagación aguas abajo (para obtener detalles, consulte la Tabla 2).

Nuestras simulaciones implican (Figura 11) que las adiciones volumétricas importantes al flujo ocurren casi en su totalidad en las laderas empinadas del volcán entre 0 y ~20 min desde el inicio del evento. Esto refleja cualitativamente la comprensión común de cómo evolucionan los lahares sineruptivos en los volcanes cubiertos de hielo: una transición rápida de un flujo inicialmente más diluido a un flujo hiperconcentrado o flujo de escombros debido al arrastre coincidente de material piroclástico en las regiones proximales (Major & Newhall, 1989; Pierson et al, 1990). Las razones de la concentración espacial modelada del impacto erosivo en los flancos del Cotopaxi están relacionadas con (i) el relieve escarpado (Figura 1c), que da como resultado grandes esfuerzos cortantes simulados (Ecuación 5), y (ii) el predominio de rocas piroclásticas sueltas y fácilmente erosionables (Figura 2), lo que facilita altas tasas de erosión (Cuadro 1). El impacto de la topografía en las fuerzas erosivas modeladas también es evidente en la reducción significativa del volumen erosionado a medida que los lahares entran en áreas de bajo gradiente al pie del volcán (ver Figura 11: después de 7–8 y 9–13 min en los drenajes norte y sur, respectivamente).

A diferencia de los modelos anteriores del Cotopaxi, esta investigación calcula la erosividad a lo largo de toda la trayectoria empinada inicial y canalizada más plana, mientras que los modelos anteriores solo anticipan esta última. Para comparar la erosión en tramos de bajo gradiente, estimamos que un lahar de magnitud similar al de 1877 (Figura 11 y Tabla 2, Escenario 4) crece solo ~4.5 × 10⁶ m³ (~6% del volumen total) al ingresar a la parte más plana de la trayectoria norte; y ~6 × 10⁶ m³ (~6% del volumen total) a lo largo de la trayectoria sur en el Valle Interandino. Nuestros números son muy similares a los valores de volumen teóricos calculados de modelos de lahar anteriores en Cotopaxi: Aguilera et al (2004) calcularon para el lahar norte de 1877 un crecimiento teórico relacionado con la erosión de 3 × 10⁶ m³ (<5% del volumen total asumido de 60 × 10⁶ m³), relacionando la altura del flujo y la pendiente del canal con los valores de volumen observados de los lahares del Nevado del Ruíz (Pierson et al, 1990). Usando la misma relación empírica, Pistolesi et al (2014) encontraron para el lahar del sur de 1877 un valor de volumen de 8,5 × 10⁶ m³ (7% del volumen total supuesto de 120 × 10⁶ m³).

En cuanto a la cuantificación de los volúmenes erosionados totales de los futuros lahares de Cotopaxi, se encontró que pueden erosionar y arrastrar entre ~12 × 10⁶ y ~43 × 10⁶ m³ (~0,012 y ~0,043 km³) de material de sustrato y potencialmente pueden crecer en volumen entre un 50 y un 400% en comparación con volumen inicial sobre una distancia de viaje modelada de ~70 km (Tabla 2). Sin embargo, no podemos validar estos valores en números absolutos porque no hay datos cuantitativos de erosión disponibles para los lahares de Cotopaxi. Sin embargo, podemos compararlos con los lahares del Nevado del Ruíz de 1985 debido a características similares relacionadas con (i) el volumen de flujo total (~90 × 10⁶ m³ o ~0,090 km³), (ii) el mecanismo desencadenante (erupción VEI 3), que generó flujos piroclásticos que interactuaron con el casquete glaciar, y (iii) la distancia de descentramiento (70–100 km desde la fuente). Además, (iv) los modelos recientes para lahares del Cotopaxi (Mothes et al, 2004; Ordóñez et al, 2013; Vera et al, 2019) utilizan un factor de volumen de 3 como regla general en todos los escenarios de erupción modelados, que es el media de los factores de volumen de 2 a 4 reportados para los lahares del Nevado del Ruíz (Pierson et al, 1990). Sin embargo, el crecimiento volumétrico observado en el evento del Nevado del Ruiz de 1985 puede no ser directamente transferible a los lahares del Cotopaxi por dos razones. Primero, la topografía específica del sitio alrededor de los dos volcanes y, por lo tanto, las presiones de corte basales ejercidas por el flujo son notablemente diferentes: en el Nevado del Ruíz los canales de los ríos mantienen pendientes pronunciadas en promedio del 30% a distancias de 20 a 50 km desde el volcán, mientras que en el Cotopaxi los gradientes de pendiente de este orden solo se extienden de 4 a 8 km desde el cráter y descienden después a alrededor de 1 a 10%. Los gradientes pronunciados del canal y la presencia de material erosionable aguas abajo del Nevado del Ruíz favorecieron un gran crecimiento del volumen desde inicialmente 20 × 10⁶ m³ (0,020 km³) hasta aproximadamente 70 × 10⁶ m³ (0,070 km³) debido al arrastre de sedimentos durante la propagación aguas abajo (Pierson et al, 1990). Esto es aproximadamente el doble de sedimento erosionado en comparación con nuestras simulaciones en Cotopaxi, donde un lahar con un volumen inicial de ~20 × 10⁶ m³ (Escenario 3, drenaje sur) puede erosionar y arrastrar ~36 × 10⁶ m³ (~0,036 km³). Y en segundo lugar encontramos que el crecimiento volumétrico de los lahares puede variar notablemente entre los diferentes escenarios de erupción futura: en el drenaje sur, los volúmenes iniciales de los lahares pequeños (Escenario 2) aumentan en ~400%, mientras que los lahares con volúmenes iniciales grandes (Escenario 4) crecen un ~100%. Esto indica que incluso los lahares iniciales pequeños del Cotopaxi son potencialmente capaces de arrastrar grandes volúmenes de sedimentos, favorecidos por la presencia de sedimentos piroclásticos fácilmente erosionables en los flancos del volcán con altas tasas de erosión y bajos esfuerzos de fluencia (Tabla 1). Las variaciones en el crecimiento del volumen tienen implicaciones importantes para la estimación de las futuras magnitudes de lahar. Nuestras simulaciones sugieren que la aplicación de un factor de crecimiento de volumen de 3 puede resultar en una subestimación de los volúmenes generales de lahares de magnitud inicialmente más pequeña (VEI 2–3) del Cotopaxi.

La clasificación espacial del sustrato en clases de erosionabilidad se llevó a cabo de forma remota en su mayor parte, respaldada por evidencia de campo en ubicaciones seleccionadas. Esto arrojó estimaciones aproximadas razonables de las cantidades y la ocurrencia espacial del arrastre de lahares; pero el mapeo detallado basado en el campo de los tipos de sustrato puede mejorar aún más el poder predictivo. No obstante, queda abierto si los parámetros del modelo para las tasas de erosión, las profundidades de erosión potenciales y el límite elástico asignado a cada clase de erosionabilidad (Tabla 1) pueden usarse para inferir la erosión y el arrastre por lahares de gran magnitud, ya que estos parámetros de erosión se derivaron de flujos de escombros naturales en los Alpes (Frank et al, 2015). El material piroclástico suelto presente en los flancos del volcán, que controla en gran medida el crecimiento volumétrico del lahar, se caracterizó por tasas de erosión dos veces más altas que las tasas promedio observadas en la cuenca del Illgraben (Frank et al, 2017), pero posiblemente el arrastre es aún más extremo. Se necesita su aplicación a eventos de lahares con patrones de erosión espacial documentados para cuantificar las posibles desviaciones. También notamos que el modelo se limita a predecir el arrastre de material de la erosión vertical en el lecho del canal, pero no considera fuentes de volumen adicionales como el colapso del banco del canal, deslizamientos de tierra importantes (re)activados a lo largo de los canales y en las áreas de origen del volcán o material volcaniclástico recién expulsado, que puede alcanzar un espesor en zonas proximales (<10 km) de 1 a 1,5 m durante una futura erupción VEI > 4 (Andrade et al, 2005).

5.2. Rendimiento del modelo para la simulación de lahar

Calculamos retroactivamente el evento lahar de 1877 del Cotopaxi y realizamos una calibración multicriterio utilizando medidas de calibración como el área de impacto, altura del flujo, velocidad del flujo, tiempos de viaje y descarga máxima (Mothes et al, 2004; Sodiro, 1877; Wolf, 1878). Al evaluar los resultados del modelo de todas las restricciones de calibración simultáneamente, pudimos derivar un buen nivel general de ajuste, no solo para una sola zona, sino para áreas más amplias de combinaciones de parámetros de entrada friccionales, donde las desviaciones del modelo son siempre menores al 25% (Figura 8). Esto corrobora la solidez del modelo y sugiere además que el comportamiento general de los lahares se puede reproducir adecuadamente con la reología probada de Voellmy-Salm.

El proceso de calibración se ve afectado por incertidumbres en el volumen total informado del lahar norte de 1877 (60 × 10⁶ a 75 × 10⁶ m³) en nuestra estrategia para estimar el hidrograma de entrada y en la precisión del DEM utilizado con respecto al cambio morfológico desde 1877. Asimismo, la evaluación del rendimiento del modelo se basa en gran medida en la comparación de los valores objetivo de los criterios de calibración con los resultados del modelo. Sin embargo, existe cierta incertidumbre en los valores objetivo en sí, que surge de las fuentes de error en las mediciones de peralte para estimar la velocidad máxima del flujo (Pierson & Scott, 1985; Prochaska et al, 2008) y del grado de preservación de los depósitos de lahar para estimar los límites de inundación y las alturas de los escenarios. El número total de rendimiento del modelo está especialmente influenciado por un valor atípico, la restricción de altura en la Hacienda Valencia, donde los valores simulados se desvían fuertemente del objetivo en un 130-170%. En esta ubicación, las alturas de flujo simuladas en todo el espacio de parámetros \mu/\xi tienen valores entre ~35 y 41 m, lo que indica que los resultados de la simulación parecen estadísticamente sólidos (Figura 7). La fuente de esta desviación podría tener su origen en (i) la complicada determinación de los extremos del depósito y (ii) mayores desviaciones en el DEM de la topografía del canal anterior al evento a medida que aumenta la distancia desde el volcán, lo que resulta de las modificaciones antropogénicas y el extenso crecimiento de la vegetación en valles aguas abajo climáticamente favorecidos. Es importante destacar que las desviaciones en la altura del flujo también se pueden atribuir al hecho de que (iii) las alturas objetivo se estiman a partir de los depósitos de flujo final, mientras que el modelo calcula las alturas máximas de las etapas. Es probable que estas incertidumbres en las estimaciones de altura también se apliquen a los resultados en San Rafael (Megamaxi). Para tales casos, el enfoque de calibración desarrollado por Aaron et al (2019) permitiría incluir una incertidumbre específica para cada criterio, pero estas no han sido cuantificadas en Mothes et al (2004) y la estimación posterior parece inapropiada. Si excluimos la restricción de altura en la Hacienda Valencia, el valor de rendimiento total se reduciría al 13% en comparación con el 23% actual para \mu = 0,005/\xi = 1400 m/s².

En general, encontramos que los valores extremadamente bajos para la fricción de tipo Coulomb \mu y los valores relativamente altos para la fricción turbulenta \xi reproducen mejor el comportamiento de los lahares sineruptivos masivos en Cotopaxi y que \mu tiene una influencia más fuerte en los resultados del modelo que \xi (Figuras 6, 7 y 8). Como esta es la primera aplicación a lahares sineruptivos, actualmente no hay parámetros de Voellmy analizados retrospectivamente disponibles para la comparación y se necesita la aplicación del modelo a otros estudios de casos para una mayor interpretación. Los valores de \mu calibrados (0,0025 y 0,005) son al menos un orden de magnitud más pequeños que los utilizados para modelar flujos de detritos no volcánicos (Hussin et al, 2012; Schraml et al, 2015), GLOF (Frey et al, 2018; Schneider et al, 2014) o avalanchas de hielo y rocas (Schneider et al, 2010; Sosio et al, 2008). Sin embargo, Sosio et al (2011) utilizaron combinaciones similares de \mu/\xi (\mu entre 0,001 y 0,02, \xi entre 700 y 1200 m/s²) para simular avalanchas de escombros volcánicos que evolucionan a partir de la falla repentina de los flancos del volcán. Una razón indiscutible de nuestro resultado de calibración es que solo las combinaciones de \mu pequeño/\xi grande pueden reproducir las distancias extremadamente largas de los lahares del Cotopaxi. Un valor pequeño de \mu y un valor alto de \xi indican una baja resistencia al flujo (Ecuación (4)) y conducen a un mecanismo de detención lento, que a su vez es la razón por la cual el modelo no logra reproducir la deposición de material en las llanuras alrededor del volcán. Los valores grandes de \xi indican además que el término de fricción turbulenta en la Ecuación (4) es marginalmente relevante para la simulación de lahares sineruptivos en este estudio.

5.3. Simulación de futuros lahares del Cotopaxi

Predecir la magnitud de futuros lahares sineruptivos es un desafío intrínsecamente difícil porque implica estimar la eficacia de las interacciones volcán-glaciar durante una erupción y el consiguiente suministro de agua, así como el potencial de incorporación de material y el crecimiento del volumen del lahar (Pierson, 1995). Para anticipar esto último, propusimos un enfoque de modelo genérico e iniciamos el lahar en las laderas empinadas del volcán superior para proporcionar una primera estimación de la capacidad de erosión de los lahares. Como se discutió anteriormente, las simulaciones de escenarios revelaron un grado variable de crecimiento de volumen entre los tres escenarios futuros potenciales y los dos drenajes. Este enfoque de modelo puede, por tanto, respaldar la estimación de los factores de volumen de los lahares específicos de los volcanes.

En nuestro modelo, la eficacia de las interacciones volcán-glaciar se cuantifica en el hidrograma de liberación de cada escenario de erupción, definiendo el volumen de lodo producido por la mezcla rápida de productos calientes de erupción con hielo y nieve. Similar a estudios previos de lahares en Cotopaxi (Aguilera et al, 2004; Barberi et al, 1992; Mothes et al, 2004), nuestras estimaciones de volúmenes iniciales de futuros lahares están directamente relacionadas con la extensión del glaciar. Debido a que la extensión de los glaciares del Cotopaxi ha disminuido en al menos ~50 % desde la última gran erupción en 1877 (Cáceres, 2017), asumimos que esto se convertiría en una reducción del 50% en los volúmenes iniciales de futuros lahares. Sin embargo, el grado en que el retroceso de los glaciares influye en los volúmenes de lahares resultantes es un tema abierto a discusión. Pistolesi et al (2013, 2014) argumentaron que las magnitudes de lahar esperadas en Cotopaxi dependen en gran medida de la dinámica de la erupción y del tipo de mecanismo desencadenante del lahar: las erupciones del tipo de 1877 (VEI 3–4 o Escenario 3) con actividad de ebullición producen principalmente flujos de escoria que erosionan cañones profundos y lineales en el glaciar. Por lo tanto, la extensión del glaciar podría tener solo un efecto menor en los volúmenes iniciales del lahar. En consecuencia, es de suponer que futuras erupciones generarán lahares de magnitud similar a la de 1877, a pesar de una reducción del 50% en la superficie del glaciar. Por el contrario, las erupciones VEI > 4 producen oleadas o flujos piroclásticos que interactúan con todo el casquete glaciar y, por lo tanto, el retroceso del glaciar ejercería un control importante sobre los volúmenes iniciales del lahar. En pocas palabras, nuestro estudio muestra que los modelos de Voellmy-Salm pueden anticipar de manera efectiva y sólida los lahares futuros, pero sigue siendo esencial comunicar las incertidumbres involucradas en la estimación de las magnitudes de los lahares futuros y, en consecuencia, las áreas de impacto de los lahares simulados a las autoridades encargadas de la toma de decisiones y las comunidades río abajo.

6. Conclusión

En este estudio simulamos lahares masivos de larga distancia del volcán Cotopaxi, desde la iniciación en los flancos superiores escarpados del volcán hasta tramos distales a unos 70 km de la fuente. El modelo de lahar genérico en el flujo de escombros RAMMS se calculó de nuevo con base en los registros del evento destructivo de 1877 (drenaje norte de 60–75 × 10⁶ m³, drenaje sur de 80–100 × 10⁶ m³) y estimaciones sobre la capacidad de erosión. Luego, el modelo calibrado se aplicó a posibles lahares futuros controlados por escenarios de erupción basados en VEI para el volcán Cotopaxi. Los hallazgos clave se pueden resumir de la siguiente manera.

La reología probada de Voellmy-Salm reproduce de manera plausible y sólida las características de flujo clave de los lahares sineruptivos masivos. Esto está respaldado por una evaluación objetiva del rendimiento del modelo basada en 14 restricciones de calibración de las propiedades físicas del flujo del lahar de 1877, como la altura del flujo, velocidad del flujo, descarga máxima, tiempo de viaje y área de impacto. Los resultados del modelo coinciden especialmente bien con el área de impacto observada (desviación del 16 %), velocidades de flujo (desviaciones de -1 a 14 %), tiempos de viaje (desviaciones de 0 y 10 %) y caudales máximos (desviaciones de -20 a 10 %), mientras que las alturas de escenario simuladas muestran mayores desviaciones de los valores objetivo. Debido al largo recorrido del lahar, los valores extremadamente bajos de fricción de tipo Coulomb \mu (0.0025–0.005) y los valores altos de fricción turbulenta \xi (1000–1400 m/s²) reproducen mejor la propagación de lahares sineruptivos.

La erosión del lahar se estima mediante un enfoque de modelo directo, que se basa en la erosionabilidad del sustrato subyacente dada por el entorno geológico y topográfico actual. Las simulaciones demuestran que el grado de crecimiento del volumen disminuye con el volumen de liberación, oscilando entre el 50% para volúmenes iniciales grandes (46 × 10⁶ m³) y el 400% para volúmenes iniciales comparativamente más pequeños (8 × 10⁶ m³), respectivamente. A pesar de las incertidumbres involucradas en este enfoque, la predicción cualitativa de la erosividad del lahar simulado en los flancos escarpados del volcán utilizando la propuesta estimación directa de la capacidad de erosión es prometedora en comparación con los patrones de erosión espacial estudiados. Se obtendrá más información mediante la transferencia de este enfoque a eventos de lahar más recientes, donde la información espacial cuantitativa y detallada sobre la erosión por flujo esté fácilmente disponible.

Una estimación confiable de las magnitudes de los futuros lahares desencadenados por erupciones y de acuerdo con los peligros para las comunidades río abajo es un desafío y debe considerar el grado de interacciones volcán-glaciar así como el grado de incorporación de sólidos en el flujo. El enfoque del modelo genérico desarrollado en este estudio incluye la erosividad del lahar y puede ayudar a anticipar el crecimiento volumétrico y los patrones de agotamiento de los lahares sineruptivos. Sin embargo, los procesos críticos que determinan las magnitudes iniciales del lahar (es decir, la magnitud de la erupción, dinámica de la erupción, interacciones entre los productos volcánicos y el hielo/nieve, impacto del retroceso de los glaciares) no se comprenden completamente y necesitan una investigación más dedicada.

Agradecimientos

Este estudio se ha realizado como parte del proyecto RIESGOS (Grant No. 03G0876D), financiado por el Ministerio Federal Alemán de Educación e Investigación (BMBF). El proyecto involucra el análisis de situaciones de riesgo múltiple y efectos en cascada relacionados en regiones piloto seleccionadas en Chile, Ecuador y Perú. En este contexto, reconocemos especialmente la colaboración con IG-EPN Ecuador y también agradecemos a Franziska Petry y Andreas Dietrich por varias discusiones y aportes. Los autores agradecen a Marc Christen (SLF) y Brian McArdell (WSL) por las discusiones y el apoyo del modelo RAMMS.

Financiamiento de acceso abierto habilitado y organizado por Projekt DEAL.

Conflicto de intereses

Los autores declaran no tener ningún conflicto de intereses.

Declaración de disponibilidad de datos

Los datos que respaldan los hallazgos de este estudio están disponibles del autor correspondiente a pedido razonable.

ORCID

Theresa Frimberger: https://orcid.org/0000-0002-0248-7292

S. Daniel Andrade: https://orcid.org/0000-0002-4814-8691

Samuel Weber: https://orcid.org/0000-0003-0720-5378

Michael Krautblatter: https://orcid.org/0000-0002-2775-2742

Referencias bibliográficas

• Aaron J, McDougall S y Nolde N. (2019) Two methodologies to calibrate landslide runout models. Landslides, 16(5), 907–920. https://doi.org/10.1007/s10346-018-1116-8
• Aguilera E, Pareschi M, Rosi M y Zanchetta G. (2004) Risk from lahars in the northern valleys of Cotopaxi Volcano (Ecuador). Natural Hazards, 33(2), 161–189. https://doi.org/10.1023/B:NHAZ.0000037037.03155.23
• Andrade D, Hall M, Mothes P, Troncoso L, Eissen JP, Samaniego P, et al (2005) Los peligros volcánicos asociados con el Cotopaxi. Corporación Editora Nacional, CEN/Escuela Politécnica Nacional: Quito.
• Auker MR, Sparks RSJ, Siebert L, Crosweller HS y Ewert J. (2013) A statistical analysis of the global historical volcanic fatalities record. Journal of Applied Volcanology, 2, art. 2. https://doi.org/10.1186/2191-5040-2-2
• Barberi F, Caruso P, Macedonio G, Pareschi M y Rosi M. (1992) Reconstruction and numerical simulation of the lahar of the 1877 eruption of Cotopaxi volcano (Ecuador). Acta Vulcanologica, 2, 35–44.
• Barberi F, Coltelli M, Frullani A, Rosi M y Almeida E. (1995) Chronology and dispersal characteristics of recently (last 5000 years) erupted tephra of Cotopaxi (Ecuador): Implications for long-term eruptive forecasting. Journal of Volcanology and Geothermal Research, 69(3–4), 217–239. https://doi.org/10.1016/0377-0273(95)00017-8
• Bartelt P, Bieler C, Buehler Y, Christen M, Deubelbeiss Y, Graf C, et al (2017) RAMMS::DEBRISFLOW User Manual. ETH: Zürich.
• Bartelt P, Salm B y Gruber U. (1999) Calculating dense-snow avalanche runout using a Voellmy-fluid model with active/passive longitudinal straining. Journal of Glaciology, 45(150), 242–254. https://doi.org/10.1017/S002214300000174X
• Berger C, McArdell BW, Fritschi B y Schlunegger F. (2010) A novel method for measuring the timing of bed erosion during debris flows and floods. Water Resources Research, 46, W02502. https://doi.org/10.1029/2009WR007993
• Berger C, McArdell BW y Schlunegger F. (2011) Sediment transfer patterns at the Illgraben catchment, Switzerland: Implications for the time scales of debris flow activities. Geomorphology, 125(3), 421–432. https://doi.org/10.1016/j.geomorph.2010.10.019
• Bernard B, Battaglia J, Proaño A, Hidalgo S, Vásconez F, Hernandez S y Ruiz M. (2016) Relationship between volcanic ash fallouts and seismic tremor: Quantitative assessment of the 2015 eruptive period at Cotopaxi volcano, Ecuador. Bulletin of Volcanology, 78(11), 80. https://doi.org/10.1007/s00445-016-1077-5
• Berti M, Genevois R, Simoni A y Tecca PR. (1999). Field observations of a debris flow event in the Dolomites. Geomorphology, 29(3–4), 265–274. https://doi.org/10.1016/S0169-555X(99)00018-5
• Beverage JP y Culbertson JK. (1964). Hyperconcentrations of suspended sediment. Journal of the Hydraulics Division, ASCE, 90, 117–128.
• Brown SK, Jenkins SF, Sparks RSJ, Odbert H y Auker MR (2017) Volcanic fatalities database: Analysis of volcanic threat with distance and victim classification. Journal of Applied Volcanology, 6, 1–20.
• Caballero L y Capra L. (2014) The use of FLO2D numerical code in lahar hazard evaluation at Popocatépetl volcano: A 2001 lahar scenario. Natural Hazards and Earth System Sciences, 14(12), 3345–3355. https://doi.org/10.5194/nhess-14-3345-2014
• Cáceres B. (2017) The role of geosciences in societal development: A German–Latin American perspective. In Evolución de los glaciares del Ecuador durante los últimos 60 años y su relación con el cambio climático, p. 149.
• Carrivick JL. (2007) Modelling coupled hydraulics and sediment transport of a high-magnitude flood and associated landscape change. Annals of Glaciology, 45, 143–154. https://doi.org/10.3189/172756407782282480
• Carrivick JL, Manville V y Cronin SJ. (2008) A fluid dynamics approach to modelling the 18th March 2007 lahar at Mt. Ruapehu, New Zealand. Bulletin of Volcanology, 71, 153–169.
• Carrivick JL, Manville V, Graettinger A y Cronin SJ. (2010) Coupled fluid dynamics–sediment transport modelling of a crater lake breakout lahar: Mt. Ruapehu, New Zealand. Journal of Hydrology, 388(3–4), 399–413. https://doi.org/10.1016/j.jhydrol.2010.05.023
• Christen M, Kowalski J y Bartelt P. (2010) RAMMS: Numerical simulation of dense snow avalanches in three-dimensional terrain. Cold Regions Science and Technology, 63(1–2), 1–14. https://doi.org/10. 1016/j.coldregions.2010.04.005
• Costa JE. (1984) Physical Geomorphology of Debris Flows. Developments and Applications of Geomorphology. Springer: Berlin.
• Delaite G, Thouret J, Sheridan M, Labazuy P, Stinton A, Souriot T y van Westen C. (2005) Assessment of volcanic hazards of El Misti and in the city of Arequipa, Peru, based on GIS and simulations, with emphasis on lahars. Zeitschrift für Geomorphologie, 140, 209–231.
• Dietrich A y Krautblatter M. (2019) Deciphering controls for debris-flow erosion derived from a LiDAR-recorded extreme event and a calibrated numerical model (Roßbichelbach, Germany). Earth Surface Processes and Landforms, 44(6), 1346–1361. https://doi.org/10.1002/esp.4578
• Doyle E, Cronin S y Thouret JC. (2011) Defining conditions for bulking and debulking in lahars. Bulletin, 123, 1234–1246.
• Fagents SA y Baloga SM. (2006) Toward a model for the bulking and debulking of lahars. Journal of Geophysical Research, 111, B10201. https://doi.org/10.1029/2005JB003986
• Fischer JT, Kowalski J y Pudasaini SP. (2012) Topographic curvature effects in applied avalanche modeling. Cold Regions Science and Technology, 74, 21–30.
• Frank F, McArdell BW, Huggel C y Vieli A. (2015) The importance of entrainment and bulking on debris flow runout modeling: Examples from the Swiss Alps. Natural Hazards and Earth System Sciences, 15 (11), 2569–2583. https://doi.org/10.5194/nhess-15-2569-2015
• Frank F, McArdell BW, Oggier N, Baer P, Christen M y Vieli A. (2017) Debris-flow modeling at Meretschibach and Bondasca catchments, Switzerland: Sensitivity testing of field-data-based entrainment model. Natural Hazards and Earth System Sciences, 17(5), 801–815. https://doi.org/10.5194/nhess-17-801-2017
• Frey H, Huggel C, Chisolm RE, Baer P, McArdell B, Cochachin A y Portocarrero C. (2018) Multi-source glacial lake outburst flood hazard assessment and mapping for Huaraz, Cordillera Blanca, Peru. Frontiers in Earth Science, 6. https://doi.org/10.3389/feart.2018.00210
• Gaunt HE, Bernard B, Hidalgo S, Proaño A, Wright H, Mothes P, et al. (2016) Juvenile magma recognition and eruptive dynamics inferred from the analysis of ash time series: The 2015 reawakening of Cotopaxi volcano. Journal of Volcanology and Geothermal Research 328, 134–146. https://doi.org/10.1016/j.jvolgeores.2016.10.013
• Hall M y Mothes P. (2008) The rhyolitic–andesitic eruptive history of Cotopaxi volcano, Ecuador. Bulletin of Volcanology, 70(6), 675–702. https://doi.org/10.1007/s00445-007-0161-2
• Hidalgo S, Battaglia J, Arellano S, Sierra D, Bernard B, Parra R, et al. (2018) Evolution of the 2015 Cotopaxi eruption revealed by combined geochemical and seismic observations. Geochemistry, Geophysics, Geosystems, 19(7), 2087–2108. https://doi.org/10.1029/2018GC007514
• Huggel C, Schneider D, Miranda PJ, Delgado Granados H y Kääb A. (2008) Evaluation of ASTER and SRTM DEM data for lahar modeling: A case study on lahars from Popocatépetl volcano, Mexico. Journal of Volcanology and Geothermal Research, 170(1–2), 99–110. https://doi.org/10.1016/j.jvolgeores.2007.09.005
• Hungr O y McDougall S. (2009) Two numerical models for landslide dynamic analysis. Computers y Geosciences, 35(5), 978–992. https://doi.org/10.1016/j.cageo.2007.12.003
• Hussin HY, Quan Luna B, van Westen CJ, Christen M, Malet JP y van Asch TWJ. (2012) Parameterization of a numerical 2-D debris flow model with entrainment: A case study of the Faucon catchment, Southern French Alps. Natural Hazards and Earth System Sciences, 12 (10), 3075–3090. https://doi.org/10.5194/nhess-12-3075-2012
• Iverson RM. (1997) The physics of debris flows. Reviews of Geophysics, 35 (3), 245–296. https://doi.org/10.1029/97RG00426
• Iverson RM. (2012) Elementary theory of bed-sediment entrainment by debris flows and avalanches. Journal of Geophysical Research – Earth Surface, 117. https://doi.org/10.1029/2011JF002189
• Iverson RM. y Denlinger RP. (2001) Flow of variably fluidized granular masses across three-dimensional terrain: 1. Coulomb mixture theory. Journal of Geophysical Research – Solid Earth, 106(B1), 537–552. https://doi.org/10.1029/2000JB900329
• Iverson RM y George DL. (2016) Modelling landslide liquefaction, mobility bifurcation and the dynamics of the 2014 Oso disaster. Géotechnique, 66(3), 175–187. https://doi.org/10.1680/jgeot.15.LM.004
• Iverson RM, Logan M, LaHusen RG y Berti M. (2010) The perfect debris flow? Aggregated results from 28 large-scale experiments. Journal of Geophysical Research – Earth Surface, 115. https://doi.org/10.1029/2009JF001514
• Iverson RM y Ouyang C. (2015) Entrainment of bed material by earth-surface mass flows: Review and reformulation of depth-integrated theory. Reviews of Geophysics, 53(1), 27–58. https://doi.org/10.1002/2013RG000447
• Iverson RM, Reid ME, Logan M, LaHusen RG, Godt JW y Griswold JP. (2011) Positive feedback and momentum growth during debris-flow entrainment of wet bed sediment. Nature Geoscience, 4(2), 116–121. https://doi.org/10.1038/ngeo1040
• Iverson RM, Schilling SP y Vallance JW. (1998) Objective delineation of lahar-inundation hazard zones. Geological Society of America Bulletin, 110(8), 972–984. https://doi.org/10.1130/0016-7606(1998)110<0972:ODOLIH>2.3.CO;2
• Lee SK, Lee CW y Lee S. (2015) A comparison of the Landsat image and LAHARZ-simulated lahar inundation hazard zone by the 2010 Merapi eruption. Bulletin of Volcanology, 77(6), 46. https://doi.org/10.1007/s00445-015-0920-4
• Lube G, Cronin SJ y Procter JN. (2009) Explaining the extreme mobility of volcanic ice-slurry flows, Ruapehu volcano, New Zealand. Geology, 37(1), 15–18. https://doi.org/10.1130/G25352A.1
• Lupiano V, Chidichimo F, Machado G, Catelan P, Molina L, Calidonna CR, et al. (2020) From examination of natural events to a proposal for risk mitigation of lahars by a cellular-automata methodology: A case study for Vascún valley, Ecuador. Natural Hazards and Earth System Sciences, 20(1), 1–20. https://doi.org/10.5194/nhess-20-1-2020
• Lupiano V, Machado G y Di Gregorio S. (2018) Revisiting the 1877 cataclysmic lahars of Cotopaxi volcano by a cellular automata model and implications for future events. In Proceedings of the 2nd International Conference on Computer Science and Application Engineering, pp. 1–7.
• Machado G, Lupiano V, Avolio MV, Gullace F y Di Gregorio S. (2015) A cellular model for secondary lahars and simulation of cases in the Vascún Valley, Ecuador. Journal of Computational Science, 11, 289–299. https://doi.org/10.1016/j.jocs.2015.08.001
• Major JJ y Newhall CG. (1989) Snow and ice perturbation during historical volcanic eruptions and the formation of lahars and floods. Bulletin of Volcanology, 52(1), 1–27. https://doi.org/10.1007/BF00641384
• Mangeney A. (2011) Landslide boost from entrainment. Nature Geoscience, 4(2), 77–78. https://doi.org/10.1038/ngeo1077
• Manville V. (2004) Palaeohydraulic analysis of the 1953 Tangiwai lahar; New Zealand’s worst volcanic disaster. Acta Vulcanologica, 16(1–2), 137–152.
• McArdell BW, Bartelt P y Kowalski J. (2007) Field observations of basal forces and fluid pore pressure in a debris flow. Geophysical Research Letters, 34(7). https://doi.org/10.1029/2006GL029183
• McCoy SW, Kean JW, Coe JA, Tucker GE, Staley DM y Wasklewicz TA. (2012) Sediment entrainment by debris flows: In situ measurements from the headwaters of a steep catchment. Journal of Geophysical Research – Earth Surface, 117, F03016. https://doi.org/10.1029/2011JF002278
• McDougall S. (2006) A new continuum dynamic model for the analysis of extremely rapid landslide motion across complex 3D terrain. PhD thesis, University of British Columbia.
• McDougall S y Hungr O. (2005) Dynamic modelling of entrainment in rapid landslides. Canadian Geotechnical Journal, 42(5), 1437–1448. https://doi.org/10.1139/t05-064
• Mizuyama T, Kobashi S y Ou G. (1992) Prediction of debris flow peak discharge. Internationales symposion: Interpraevent: 99–108.
• Mothes P, Espin P, Hall M, Vasconez F, Sierra D y Andrade D. (2016) Mapa regional de amenazas volcánicas potenciales del volcán Cotopaxi, zona norte. Instituto Geofísico de la Escuela Politécnica Nacional.
• Mothes P, Hall ML, Andrade D, Yepes H, Pierson TC, Gorki Ruiz A y Samaniego P. (2004) Character, stratigraphy and magnitude of historical lahars of Cotopaxi volcano (Ecuador). Acta Vulcanologica, 16(1–2), 85–108.
• Mothes PA, Hall ML y Janda RJ. (1998) The enormous Chillos Valley Lahar: An ash-flow-generated debris flow from Cotopaxi Volcano, Ecuador. Bulletin of Volcanology, 59(4), 233–244. https://doi.org/10.1007/s004450050188
• Mothes PA, Ruiz MC, Viracucha EG, Ramón PA, Hernández S, Hidalgo S, et al. (2017) Geophysical footprints of Cotopaxi’s unrest and minor eruptions in 2015: An opportunity to test scientific and community preparedness. In , Unrest V, Gottsmann J, Neuberg J, Scheu B. (eds). Springer: Cham; 241–270.
• Muñoz-Salinas E, Castillo-Rodríguez M, Manea V, Manea M y Palacios D. (2009) Lahar flow simulations using LAHARZ program: Application for the Popocatépetl volcano, Mexico. Journal of Volcanology and Geothermal Research, 182(1–2), 13–22. https://doi.org/10.1016/j.jvolgeores.2009.01.030
• Newhall CG y Self S. (1982) The volcanic explosivity index (VEI) an estimate of explosive magnitude for historical volcanism. Journal of Geophysical Research, Oceans, 87(C2), 1231–1238. https://doi.org/10.1029/JC087iC02p01231
• O’Brien JS, Julien PY y Fullerton W. (1993) Two-dimensional water flood and mudflow simulation. Journal of Hydraulic Engineering, 119 (2), 244–261. https://doi.org/10.1061/(ASCE)0733-9429(1993)119: 2(244)
• Ordóñez J, Samaniego P, Mothes P y Schilling S. (2013) Las potenciales zonas de inundación por lahares en el volcán Cotopaxi. Publicación Instituto Geofísico-Escuela Politécnica Nacional, N. 20. Pierson, T.C. (1995) Flow characteristics of large eruption-triggered debris flows at snow-clad volcanoes: Constraints for debris-flow models. Journal of Volcanology and Geothermal Research, 66(1–4), 283–294. https://doi.org/10.1016/0377-0273(94)00070-W
• Pierson TC. (1998) An empirical method for estimating travel times for wet volcanic mass flows. Bulletin of Volcanology, 60(2), 98–109. https://doi.org/10.1007/s004450050219
• Pierson TC. (2005) Hyperconcentrated flow – transitional process between water flow and debris flow. In Debris-flow Hazards and Related Phenomena, Jakob M, Hungr O. (eds). Springer: Cham; 159–202.
• Pierson TC, Costa JE y Vancouver W. (1987) A rheologic classification of subaerial sediment-water flows. In Debris Flows/Avalanches, Vol. 7, Costa, JE., Wieczorek, GF. (eds). Geological Society of America: Boulder, CO. https://doi.org/10.1130/REG7-p1
• Pierson TC, Janda RJ, Thouret JC y Borrero CA. (1990) Perturbation and melting of snow and ice by the 13 November 1985 eruption of Nevado del Ruiz, Colombia, and consequent mobilization, flow and deposition of lahars. Journal of Volcanology and Geothermal Research, 41(1–4), 17–66. https://doi.org/10.1016/0377-0273(90)90082-Q
• Pierson TC y Scott KM. (1985) Downstream dilution of a lahar: Transition from debris flow to hyperconcentrated streamflow. Water Resources Research, 21(10), 1511–1524. https://doi.org/10.1029/WR021i010p01511
• Pistolesi M, Cioni R, Rosi M y Aguilera E. (2014) Lahar hazard assessment in the southern drainage system of Cotopaxi volcano, Ecuador: Results from multiscale lahar simulations. Geomorphology, 207, 51–63. https://doi.org/10.1016/j.geomorph.2013.10.026
• Pistolesi M, Cioni R, Rosi M, Cashman KV, Rossotti A y Aguilera E. (2013) Evidence for lahar-triggering mechanisms in complex stratigraphic sequences: The post-twelfth century eruptive activity of Cotopaxi Volcano, Ecuador. Bulletin of Volcanology, 75(3), 698. https://doi.org/10.1007/s00445-013-0698-1
• Pistolesi M, Rosi M, Cioni R, Cashman KV, Rossotti A y Aguilera E. (2011) Physical volcanology of the post–twelfth-century activity at Cotopaxi volcano, Ecuador: Behavior of an andesitic central volcano. Bulletin, 123, 1193–1215.
• Pitman EB y Le L. (2005) A two-fluid model for avalanche and debris flows. Philosophical Transactions Series A: Mathematical, Physical, and Engineering Sciences, 363(1832), 1573–1601. https://doi.org/10.1098/rsta.2005.1596
• Prochaska AB, Santi PM, Higgins JD y Cannon SH. (2008) A study of methods to estimate debris flow velocity. Landslides, 5(4), 431–444. https://doi.org/10.1007/s10346-008-0137-0
• Procter JN, Cronin SJ, Fuller IC, Sheridan M, Neall VE y Keys H. (2010) Lahar hazard assessment using Titan2D for an alluvial fan with rapidly changing geomorphology: Whangaehu River, Mt. Ruapehu. Geomorphology, 116(1–2), 162–174. https://doi.org/10.1016/j.geomorph.2009.10.016
• Pudasaini SP y Fischer JT. (2016) A mechanical erosion model for two-phase mass flows. arXiv preprint arXiv:1610.01806.
• Revellino P, Hungr O, Guadagno FM y Evans SG. (2004) Velocity and runout simulation of destructive debris flows and debris avalanches in pyroclastic deposits, Campania region, Italy. Environmental Geology, 45(3), 295–311. https://doi.org/10.1007/s00254-003-0885-z
• Rickenmann D. (1999) Empirical relationships for debris flows. Natural Hazards, 19(1), 47–77. https://doi.org/10.1023/A:1008064220727
• Salm B. (1993) Flow, flow transition and runout distances of flowing avalanches. Annals of Glaciology, 18, 221–226. https://doi.org/10.3189/S0260305500011551
• Schneider D, Bartelt P, Caplan-Auerbach J, Christen M, Huggel C y McArdell BW. (2010) Insights into rock-ice avalanche dynamics by combined analysis of seismic recordings and a numerical avalanche model. Journal of Geophysical Research, 115, F04026. https://doi.org/10.1029/2010JF001734
• Schneider D, Granados HD, Huggel C y Kääb A. (2008) Assessing lahars from ice-capped volcanoes using ASTER satellite data, the SRTM DTM and two different flow models: Case study on Iztaccíhuatl (Central Mexico). Natural Hazards and Earth System Sciences, 8 (3), 559–571. https://doi.org/10.5194/nhess-8-559-2008
• Schneider D, Huggel C, Cochachin A, Guillén S y García J. (2014) Mapping hazards from glacier lake outburst floods based on modelling of process cascades at Lake 513, Carhuaz, Peru. Advances in Geosciences, 35, 145–155. https://doi.org/10.5194/adgeo-35-145-2014
• Schraml K, Thomschitz B, McArdell BW, Graf C y Kaitna R. (2015) Modeling debris-flow runout patterns on two alpine fans with different dynamic simulation models. Natural Hazards and Earth System Sciences, 15(7), 1483–1492. https://doi.org/10.5194/nhess-15-1483-2015
• Schürch P, Densmore AL, Rosser NJ y McArdell BW. (2011) Dynamic controls on erosion and deposition on debris-flow fans. Geology, 39 (9), 827–830. https://doi.org/10.1130/G32103.1
• Scott KM. (1988) Origins, behavior, and sedimentology of lahars and lahar-runout flows in the Toutle-Cowlitz River system. Professional Paper 1447-A.
• Scott KM, Vallance JW, Kerle N, Luis Macías J, Strauch W y Devoli G. (2005) Catastrophic precipitation-triggered lahar at Casita volcano, Nicaragua: Occurrence, bulking and transformation. Earth Surface Processes and Landforms, 30(1), 59–79. https://doi.org/10.1002/esp.1127
• Sierra D, Vasconez F, Andrade SD, Almeida M y Mothes P. (2019) Historical distal lahar deposits on the remote eastern-drainage of Cotopaxi volcano, Ecuador. Journal of South American Earth Sciences, 95, 102251. https://doi.org/10.1016/j.jsames.2019.102251
• Small C y Naumann T. (2001) The global distribution of human population and recent volcanism. Global Environmental Change Part B: Environmental Hazards, 3(3–4), 93–109. https://doi.org/10.1016/S1464-2867(02)00002-5
• Smith GA y Fritz WJ. (1989) Volcanic influences on terrestrial sedimentation. Geology, 17(4), 375–376. https://doi.org/10.1130/0091-7613(1989)017<0375:VIOTS>2.3.CO;2
• Sodiro L. (1877) Relación sobre la erupción del Cotopaxi acaecida el día 26 de Junio de 1877.
• Sosio R, Crosta GB y Hungr O. (2008) Complete dynamic modeling calibration for the Thurwieser rock avalanche (Italian Central Alps). Engineering Geology, 100(1–2), 11–26. https://doi.org/10.1016/j.enggeo.2008.02.012
• Sosio R, Crosta GB y Hungr O. (2011) Numerical modeling of debris avalanche propagation from collapse of volcanic edifices. Landslides, 9, 315–334.
• Thouret JC, Antoine S, Magill C y Ollier C. (2020) Lahars and debris flows: Characteristics and impacts. Earth-Science Reviews, 201, 103003. https://doi.org/10.1016/j.earscirev.2019.103003
• Toapaxi J, Torres C, Hamad K, Vera P y Méndez J. (2019) Actualización de la modelación numérica del flujo de escombros producto de una erupción del volcán Cotopaxi–Flanco Norte. Revista Politécnica, 44(1), 7–14. https://doi.org/10.33333/rp.vol44n1.01
• Vallance JW y Iverson RM. (2015) Lahars and their deposits. In The Encyclopedia of Volcanoes. Elsevier: Amsterdam.
• Vera P, Ortega P, Casa E, Santamaría J y Hidalgo X. (2019) Modelación numérica y mapas de afectación por flujo de lahares primarios en el drenaje sur del volcán Cotopaxi. Revista Politécnica, 43(1), 61–72. https://doi.org/10.33333/rp.vol43n1.971
• Vezzoli L, Apuani T, Corazzato C y Uttini A. (2017) Geological and geotechnical characterization of the debris avalanche and pyroclastic deposits of Cotopaxi Volcano (Ecuador). Journal of Volcanology and Geothermal Research, 332, 51–70. https://doi.org/10.1016/j.jvolgeores.2017.01.004
• Voellmy A. (1955) Ueber die Zersstorungskraft von Lawinen. Schweizerische Bauzeitung, 73, 212–217.
• Vuille M, Carey M, Huggel C, Buytaert W, Rabatel A, Jacobsen D, et al. (2018) Rapid decline of snow and ice in the tropical Andes – impacts, uncertainties and challenges ahead. Earth-Science Reviews 176, 195–213. https://doi.org/10.1016/j.earscirev.2017.09.019
• Williams R, Stinton A y Sheridan M. (2008) Evaluation of the Titan2D two-phase flow model using an actual event: Case study of the 2005 Vazcún Valley Lahar. Journal of Volcanology and Geothermal Research, 177(4), 760–766. https://doi.org/10.1016/j.jvolgeores.2008.01.045
• Wolf T. (1878) Memoria sobre El Cotopaxi y su ultima erupcion, acaecida el 26 de junio de 1877.
• Worni R, Huggel C, Stoffel M y Pulgarín B. (2011) Challenges of modeling current very large lahars at Nevado del Huila Volcano, Colombia. Bulletin of Volcanology, 74, 309–324.

Información de soporte

Puede encontrar información de apoyo adicional en línea en la sección Información de apoyo al final de este artículo.

Cómo citar este artículo:

Frimberger T, Andrade SD, Weber S, Krautblatter M. Modelling future lahars controlled by different volcanic eruption scenarios at Cotopaxi (Ecuador) calibrated with the massively destructive 1877 lahar. Earth Surf. Process. Landforms. 2021;46:680–700. https://doi.org/10.1002/esp.5056

Medios