Métodos con base física#

Los modelos con base física evalúan la ocurrencia de movimientos en masa en términos de factor de seguridad o probabilidades de ocurrencia, a través de modelos numéricos. Sin embargo existen diferentes aproximaciones a la física del problema con simplificaciones diferentes, por lo cual es fundamental utilizar el modelo que se ajuste a las necesidades o particularidades del área de estudio. Los métodos con base físisca son recomendables para áreas pequeñas a escalas detalladas, ya que requieren una gran cantidad de información de entrada, que se obtiene de ensayos de laboratorio o mediciones de campo. Por lo que la parametrización de los modelos puede ser complicada, en especial la distribución espacial de la profundidad del suelo, la cual juega un papel fundamental. Los modelos con base física permiten tener en consideración la complejidad de los factores detonantes (sismo y/o lluvia) y las características geomecánicas del terreno.

Factor detonante lluvia#

El componente hidrológico de estos modelos explica el proceso de infiltración del agua en el suelo y el avance del frente húmedo o aumento de un nivel freático colgado. Para esto utilizan modelos de flujo estático horizontal o modelos de flujo transitorio vertical. Desde el componente geotécnico, la mayoría de modelos utilizan métodos de equilibrio límite de talud infinito. Algunos han incoporado modelos de dovelas que permiten análisis de superficies de falla circulares.

Como parte de los modelos físicos una gran variedad de modelos hidrológicos se han utilizado considerando el flujo de agua en el suelo en estado estático o transitorio, y en condiciones saturadas o parcialmente saturadas, utilizando las ecuaciones de Richards para el flujo vertical y las ecuaciones de Darcy para el flujo lateral. Entre estos se encuentra TOPOG (O’Loughlin, 1986), TAPES (Moore et al., 1988), TOPMODEL (Beven & Kirkby, 1979), GEOtop (Bertoldi & Rigon, 2004), tRIBS (Ivanov et al., 2004), SHIA (Vélez 2001), y los modelos de Green-Ampt (1911), Pradel & Raad (1993) e Iverson (2000). Por el contrario, los modelos geotécnicos utilizados son muy similares, generalmente utilizando el método talud infinito y análisis de equilibrio límite, suponiendo la superficie de falla paralela a la superficie del terreno y que la longitud de falla es mucho mayor que el espesor de la capa desplazada (Borga et al., 2002).

modelo hidrologico

Fig. 63 Dos posibles mecanismos de saturación para deslizamientos planares superficiales detonados por lluvias. Izquierda: incremento del nivel freático colgado paralelo a la ladera a partir de una superficie de contraste de permeabilidad. Derecha: avance de frente húmedo a partir de la superficie del terreno.#

Uno de los primeros modelos físicos fue propuesto por Wu y Sidle (1995) denominado dSLAM (distribuited Shallow LAndslide Model), y posteriormente modificado como IDSSM (Integrated Dynamic Slope Stability Model) por Dhakal & Sidle (2003). Corresponde a un modelo físico, dinámico y distribuido, que utiliza un análisis de estabilidad de talud infinita con un modelo hidrológico de onda cinemática propuesto por Takasao y Shiiba (1988). Este modelo considera la vegetación y asume que la infiltración es siempre mayor a la intensidad de la lluvia.

Sin embargo uno de los modelos físicos más reconocido por su sencillez y fácil aplicación es el modelo propuesto por Montgomery & Dietrich (1994) y Montgomery et al. (1998), denominado SHAllow Landslide STABility (SHALSTAB). Este modelo emplea el modelo hidrológico TOPOG para estimar la altura de la porción saturada de suelo, el cual asume que el control dominante de la distribución espacial de los movimientos está dado por la topografía. Estos autores definen para su análisis un índice topográfico, el cual es utilizado para predecir el nivel freático en función del flujo del agua en el suelo y la intensidad de la lluvia. El modelo de estabilidad propuesto utiliza el criterio de Mohr-Coulomb, en el cual por simplificación asumen la cohesión igual cero.

Otros modelos conceptualmente similares al SHALSTAB se han desarrollado. El modelo SINMAP, desarrollado por Pack et al. (1998), evalúa las condiciones de estabilidad para movimientos en masa, utilizando la misma ecuación del factor de seguridad y las ecuaciones de Darcy para flujos saturados, pero a diferencia del SHALSTAB considera la cohesión del suelo. Y el modelo LISA (por sus siglas en inglés Level I Stability Analysis), desarrollado por Hammond et al. (1992), realiza un análisis probabilístico basado en el factor de seguridad, considerando la vegetación. Los valores para cada parámetro en la ecuación son definidos por una función de distribución de probabilidades y los resultados son presentados en un histograma que muestra la distribución del factor de seguridad calculado usando el método de Monte Carlo.

El modelo GEOtop –FS, desarrollado por Simoni et al. (2008), combina un análisis de estabilidad de talud infinito con el modelo hidrológico distribuido espacialmente y de diferencias finitas, GEOtop. El modelo estima el contenido de agua en el suelo y la evolución de la presión de poros resultante simulando los procesos de infiltración, escorrentía superficial, flujo subsuperficial saturado y no saturado, y flujos en canales. Para lo cual todas las celdas de la cuenca son divididas en celdas tipo canal y tipo laderas. El factor de seguridad se calcula con un enfoque probabilístico, donde los parámetros del suelo se asignan con funciones de distribución.

El modelo CHASM (Combined Hydrological And Stability Model), desarrollado por Wilkinson et al. (2000), utiliza el método simplificado circular de Bishop (1955) con un modelo hidrológico bidimensional de diferencias finitas para simular el cambio en las presiones de poros como respuesta a eventos de lluvia individuales. En cada intervalo de tiempo de simulación, las presiones de poro, tanto negativas como positivas, se incorporan directamente en la determinación de los esfuerzos efectivos y la resistencia al corte del suelo. Este análisis implica una búsqueda numérica de la superficie de deslizamiento que reduce al mínimo el factor de seguridad.

El modelo Transient Rainfall Infiltration and Gridbased Slope-Stability (TRIGRS), desarrollado por Baum et al. (2002), utiliza el modelo hidrológico transitorio unidimensional para infiltración vertical propuesto por Iverson (2000), con un modelo de estabilidad de talud infinita asumiendo condiciones saturadas. El modelo evalúa el factor de seguridad como respuesta de las presiones de poros durante un evento de lluvia función de la profundidad y el tiempo, para lo cual divide el factor de seguridad en un componente estático (Fs0) asociados al estado estático o de largo plazo, y un componente que varía con el tiempo (FS)’ asociado al estado transitorio del evento de lluvia simulado.

Arnone et al. (2011) proponen un modelo físico distribuido, en el cual acoplan un modelo de talud infinita con el modelo hidrológico tRIBS (Triangulated Irregular Network Real-Time Integrated Basin Simulator) desarrollado por Ivanov et al. (2004). El modelo simula la intercepción de la lluvia, la evapotranspiración, la infiltración del agua en el suelo, la escorrentía superficial, y el flujo subsuperficial lateral a la ladera en condiciones saturadas y parcialmente saturadas, representando la heterogeneidad espacial de la topografía por una red irregular triangular (TIN) utilizada ampliamente por Sistemas de Información Geográfica (SIG) que reduce el tiempo computacional sin pérdidas aparentes de información. El modelo es también capaz de simular el volumen del material fallado y su posible trayecto.

Finalmente el modelo SHIA_Landslide propuesto por Aristizábal et al. (2014) acopla un modelo de talud infinita en condiciones saturadas con el modelo hidrológico distribuido de tanques a escala de cuenca que considera los flujos verticales y horizontales, denominado SHIA (Vélez, 2001). El modelo hidrológico está conformado por dos componentes: un balance de agua que simula los proceso hidrológicos dominantes en la cuenca, tales como infiltración, percolación, escorrentía superficial, flujo sub superficial, y flujo subterráneo; y un componente que simula el flujo del agua a través de la red de drenaje, lo que permite calibrar con valores observados la hidrología de la cuenca y las laderas previo a la simulación de estabilidad de las vertientes.

modelos

Fig. 64 modelos físicos.#

Modelo SHALSTAB#

El modelo SHALSTAB, desarrollado por Montgomery & Dietrich, (1994), emplea el modelo hidrológico TOPOG (O’Loughlin, 1986) en condiciones de lluvia estacionaria para construir un mapa del patrón de la humedad basado en el área aferente a cada punto, la pendiente y la transmisividad del suelo. En el modelo la humedad del suelo es calculada como:

\(W=\frac{Qa}{bTsinθ}\)

Donde \(Q\) es la lluvia en condiciones estacionarias (\(mm/d\)), a es el área de drenaje (\(m^2\)), \(b\) es la longitud de cada celda (\(m\)), \(T\) es la transmisividad del suelo en condiciones saturadas (\(m^2/\)), \(\theta\) y es la pendiente. Asumiendo que la transmisividad no varía con la profundidad, se puede entonces asumir:

\(W= \frac{h}{z_t}\)

Donde \(h\) es el espesor del suelo saturado y \(z_t\) el espesor total de suelo. Combinando las dos ecuaciones anteriores se puede estimar la saturación relativa del perfil de suelo como:

\(\frac{h}{z_t}=\frac{Qa}{bTsinθ} \)

En cuanto al componente geotécnico, el modelo SHALSTAB está basado en análisis de equilibrio límite con talud infinito y el criterio de falla de Mohr-Coulomb.

\(ρ_sgz_tsinθcosθ=C'+[ρ_s-ρ_w\frac{h}{z_t}]gzcos^2θtanϕ \)

Donde \(ρ_s\) es la densidad del suelo, \(ρw\) es la densidad del agua, \(g\) es la aceleración de la gravedad, \(C’\) es la cohesión efectiva del suelo, y \(φ\) es el ángulo de fricción. Esta ecuación puede ser escrita en términos de la relación h/z como:

\(\frac{h}{z_t}=\frac{ρ_s}{ρ_w}(1-\frac{tanθ}{tanϕ})+\frac{C'}{wgz_tcos^2θtanϕ}\)

Finalmente acoplando el modelo hidrológico con el modelo de estabilidad se obtienen la siguiente ecuación:

\(\frac{a}{b}=\frac{T}{Q}sinθ[\frac{ρ_s}{ρ_w}(1-\frac{tanθ}{tanϕ})+\frac{C'}{wgzcos^2θtanϕ}]\)

A partir de esta ecuación es posible determinar cuatro condiciones de estabilidad para cada celda de análisis. Las celdas donde la relación entre el área de drenaje aferente y la longitud de la celda (a/b) es mayor que la expresión al lado derecho de la ecuación corresponde a celdas inestables, en caso contrario son celdas estables. Las dos condiciones restantes corresponden a condiciones de estabilidad que no dependen de la lluvia. Las celdas estables en condiciones completamente saturadas de todo el perfil de suelo (\(h/z_t=1\)) son denominadas incondicionalmente estables y cumplen la siguiente condición:

\(tanθ < 1-(\frac{ρ_s}{ρ_w})tanϕ+\frac{C'}{ρ_sgz_tcos^2θ} \)

Y las celdas inestables en condiciones secas (\(h/z_t=0\)) se denominan incondicionalmente inestables y cumplen la siguiente condición:

\(tanθ >= tanϕ+\frac{C'}{ρ_sgz_tcos^2θ}\)

Python#

Se deben inicialmente importar las librerías. En este caso se instalará la librería Rasterio

!pip install rasterio
Requirement already satisfied: rasterio in c:\users\edier\miniconda3\lib\site-packages (1.3.9)
Requirement already satisfied: affine in c:\users\edier\miniconda3\lib\site-packages (from rasterio) (2.4.0)
Requirement already satisfied: attrs in c:\users\edier\miniconda3\lib\site-packages (from rasterio) (23.1.0)
Requirement already satisfied: certifi in c:\users\edier\miniconda3\lib\site-packages (from rasterio) (2023.11.17)
Requirement already satisfied: click>=4.0 in c:\users\edier\miniconda3\lib\site-packages (from rasterio) (8.1.6)
Requirement already satisfied: cligj>=0.5 in c:\users\edier\miniconda3\lib\site-packages (from rasterio) (0.7.2)
Requirement already satisfied: numpy in c:\users\edier\miniconda3\lib\site-packages (from rasterio) (1.25.1)
Requirement already satisfied: snuggs>=1.4.1 in c:\users\edier\miniconda3\lib\site-packages (from rasterio) (1.4.7)
Requirement already satisfied: click-plugins in c:\users\edier\miniconda3\lib\site-packages (from rasterio) (1.1.1)
Requirement already satisfied: setuptools in c:\users\edier\miniconda3\lib\site-packages (from rasterio) (67.8.0)
Requirement already satisfied: colorama in c:\users\edier\miniconda3\lib\site-packages (from click>=4.0->rasterio) (0.4.6)
Requirement already satisfied: pyparsing>=2.1.6 in c:\users\edier\miniconda3\lib\site-packages (from snuggs>=1.4.1->rasterio) (3.0.9)
import rasterio as rio
import numpy as np
import matplotlib.pyplot as plt

A continuación se procede a importar los mapas requeridos: cohesion, fricción, permeabilidad, peso unitario del suelo, área acumulada, pendiente y espesor del suelo.

raster = rio.open('https://github.com/edieraristizabal/Libro_cartoGeotecnia/blob/master/data/arenosa/cohesion.asc?raw=true') # (kN/m2)
cohesion=raster.read(1)
cohesion=np.where(cohesion==cohesion[0,0],np.nan,cohesion)
plt.imshow(cohesion)
plt.colorbar();
_images/80bc161d6b48715c2d652319bacd78583f9676b33c49752c4610c98b0770790f.png
raster = rio.open('https://github.com/edieraristizabal/Libro_cartoGeotecnia/blob/master/data/arenosa/friction_rad.asc?raw=true') #(radianes)
friccion=raster.read(1)
friccion=np.where(friccion==friccion[0,0],np.nan,friccion) #ya py lee los nans
plt.imshow(friccion)
plt.colorbar();
_images/6a6cfaf8a4a2e41b4dd27856022c61416f2f93fa039f1debb308fbebe0b4911a.png
raster = rio.open('https://github.com/edieraristizabal/Libro_cartoGeotecnia/blob/master/data/arenosa/ks.asc?raw=true') #(cm/hora)
ks=raster.read(1)
ks[0,0]
ks=np.where(ks==ks[0,0],np.nan,ks)
plt.imshow(ks)
plt.colorbar();
_images/62c80fcf591406928df74b8883fd70f28e2a1a410b53a02d26ac87c100e08520.png
raster = rio.open('https://github.com/edieraristizabal/Libro_cartoGeotecnia/blob/master/data/arenosa/gammas.asc?raw=true') #(kN/m3)
peso=raster.read(1)
peso=np.where(peso==peso[0,0],np.nan,peso)
plt.imshow(peso)
plt.colorbar();
_images/95ee1c0690df09a4320baaf279f4b6e2630c0de68829e928dc14aedaeef37dc3.png
raster = rio.open('https://github.com/edieraristizabal/Libro_cartoGeotecnia/blob/master/data/arenosa/flowacum_m2.asc?raw=true') #(m2)
flujo=raster.read(1)
flujo=np.where(flujo==flujo[0,0],np.nan,flujo)
plt.imshow(flujo)
plt.colorbar();
_images/7de73d5568bc3a481e123104b570e48dad2c3bf20510cf6eb215949938447370.png
raster = rio.open('https://github.com/edieraristizabal/Libro_cartoGeotecnia/blob/master/data/arenosa/slope_rad.asc?raw=true') #(radianes)
pendiente=raster.read(1)
pendiente=np.where(pendiente==pendiente[0,0],np.nan,pendiente)
plt.imshow(pendiente)
plt.colorbar()
<matplotlib.colorbar.Colorbar at 0x24113a32850>
_images/a4ad09f4a09622c2304dabe072c76bd33d5502a140a09a942074c8135245cf61.png
raster = rio.open('https://github.com/edieraristizabal/Libro_cartoGeotecnia/blob/master/data/arenosa/zs.asc?raw=true') #(m)
espesor=raster.read(1)
espesor=np.where(espesor==espesor[0,0],np.nan,espesor)
plt.imshow(espesor)
plt.colorbar()
<matplotlib.colorbar.Colorbar at 0x241111eded0>
_images/2d2d44e8214d94c78fa31f6cbff1fa8034a5809d77bc63b8c9d5ea9c227641ab.png

Para implementar en Python las ecuaciones del modelo Shalstab, se debe primero definir las constantes del peso unitario del agua, y la resolución espacial del DEM utilizado.

En este paso tambien, se incorpora la intensidad de la lluvia que quiere evaluarse. Este valor no modifica las celdas incondicionalmente estables o inestables, sólo modifica el número de celdas potencialmente inestables.

GammaW   = 9.81 #peso unitario del agua [kN/m3]
dx = 12.5 # resolucion espacial del raster [m]
q = 100 #intensidad de la precipitacion [mm/h]

Finalmente se procede a correr el modelo

"STABILITY ANALYSIS"

MatEst = np.zeros(raster.shape)
Matq   = np.zeros(raster.shape)

M4=flujo/dx
M5=((0.01 * ks * (espesor * np.cos(pendiente)) * np.sin(pendiente)) / (0.001*q)) * ((peso / GammaW) * (1 - np.tan(pendiente) / np.tan(friccion)) + (cohesion / (GammaW * espesor * np.cos(pendiente)**2 * np.tan(friccion)))) 
MatEst1=np.where(M4>M5,3,MatEst) #  unstable

MatEst2=np.where(M4<=M5,4,MatEst1) # Stable

M1=np.tan(pendiente)
M2=(1 - (GammaW/peso)) * np.tan(friccion) + (cohesion / (peso * espesor * np.cos(pendiente)**2))
MatEst3  =np.where(M1<M2,1,MatEst2) # Unconditionally stable
	
M3=np.tan(friccion) + (cohesion / (peso * espesor * np.cos(pendiente)**2))
MatEst4  =  np.where(M1>=M3,2,MatEst3) # Unconditionally Unstable

Matq = (1000 * 0.01 * ks * espesor * np.cos(pendiente) * np.sin(pendiente)) * (dx / flujo) * ((peso / GammaW) * (1 - (np.tan(pendiente) / np.tan(friccion))) + cohesion / (GammaW * espesor * np.cos(pendiente)**2 * np.tan(friccion)))
np.nanmin(Matq)

MatEst=np.where(MatEst4 == 0, np.nan,MatEst4)
Matq=np.where(Matq==np.inf,np.nanmax(Matq[Matq!=np.inf]),Matq)

Matq = np.where(M1<M2,-1,np.where(M1>=M3,-2,Matq))
C:\Users\edier\AppData\Local\Temp\ipykernel_31676\3061725108.py:19: RuntimeWarning: divide by zero encountered in divide
  Matq = (1000 * 0.01 * ks * espesor * np.cos(pendiente) * np.sin(pendiente)) * (dx / flujo) * ((peso / GammaW) * (1 - (np.tan(pendiente) / np.tan(friccion))) + cohesion / (GammaW * espesor * np.cos(pendiente)**2 * np.tan(friccion)))

Como resultado se obtienen 2 mapas. El mapa de celdas estables (celdas con valor de 4), inestables (celdas con valor de 3), incondicionalmente inestables (celdas con valor de 2) e incondicionalmente estables (celdas con valor de 1), el cual se denomina MatEst.

plt.imshow(MatEst)
plt.colorbar()
<matplotlib.colorbar.Colorbar at 0x241111ee8d0>
_images/976e78ce15bad64e517524c458aef8fcad52929422942384c439f1e1a79b88fa.png

El segundo mapa que se obtiene es Matq, que representa el valor de lluvia (q) necesario para que fallen las celdas. Este valor no aplica para celdas incondicionalmente inestables o incondicionalmente estables, solo para las celdas potencialmente inestables, es decir que arrojaron valores de 4 y 3 en el mapa MatEst.

plt.imshow(Matq)
plt.colorbar()
<matplotlib.colorbar.Colorbar at 0x24113295e50>
_images/30d67ecf611e13c70fa1d8da2b817d94a68434548696f5d52f9cbd1edaf05922.png

Para exportar los mapas se puede utilizar el siguiente código:

with rio.open('MatEst.tif', 'w', 
              driver='Gtiff',height=raster.shape[0],width=raster.shape[1],count=1,
              dtype='float64',nodata=-999,crs=raster.crs,transform=raster.transform) as dst:
    dst.write(MatEst,1)
with rio.open('Matq.tif', 'w', 
              driver='Gtiff',height=raster.shape[0],width=raster.shape[1],count=1,
              dtype='float64',nodata=-999,crs=raster.crs,transform=raster.transform) as dst:
    dst.write(Matq,1)  

Para conocer el volumen del material desplazado en cada celda se utiliza el siguiente código, donde además se exporta dicho mapa.

Vol = np.where(np.logical_or(MatEst==2, MatEst==3), espesor*12.5*12.5, 0)

with rio.open('Vol.tif', 'w', 
              driver='Gtiff',height=raster.shape[0],width=raster.shape[1],count=1,
              dtype='float64',nodata=-999,crs=raster.crs,transform=raster.transform) as dst:
    dst.write(Vol,1) 

plt.imshow(Vol)
plt.colorbar()
<matplotlib.colorbar.Colorbar at 0x24111063e90>
_images/aa42efbb182b4ba8646e08e30cb21d3f9137daa739270db1c528d0a7c075ac08.png

Modelo TRIGRS#

El modelo TRIGRS (Transient Rainfall Infiltration and Gridbased Slope-Stability) incorpora en su análisis el proceso de infiltración del agua lluvia y el cálculo del factor de seguridad. Este es un modelo con base física que evalúa la distribución temporal y espacial de movimientos en masa superficiales detonados por lluvia a través del cálculo de los cambios transitorios de la presión de poros y su incidencia en la variación del factor de seguridad, debido a la infiltración de la lluvia.

TRIGRS

Fig. 65 Modelo TRIGRS.#

El modelo de infiltración está basado en la solución lineal de las ecuaciones de Richards, y el flujo de agua en el suelo es el resultado de la sumatoria del estado estacionario y el componente transitorio asociado al evento de lluvia modelado. La solución para el caso de frontera basal impermeable a una profundidad finita está dada por:

Donde 𝜓 es la cabeza de presión, 𝑡 el tiempo.

\(𝑍 = 𝑧/𝑐𝑜𝑠 𝛿\)

Donde 𝑍 es la coordenada en dirección vertical (positiva hacia abajo), \(𝑧\) la coordenada en dirección normal al talud y \(δ\) es el ángulo del terreno con la horizontal; \(𝑑\) es la profundidad inicial del nivel en dirección vertical.

\(𝐵 = 𝑐𝑜𝑠2𝛿 − \frac{𝐼𝑍𝐿𝑇}{𝐾_𝑠}\)

Donde \(𝐾_𝑠\) es la conductividad hidráulica saturada en dirección \(𝑍\), \(𝐼𝑍𝐿𝑇\) la tasa de infiltración estacionaria (inicial) en la superficie del suelo. 𝐼𝑛𝑍 es la tasa de infiltración a una intensidad dada para el n-ésimo intervalo de tiempo.

\(𝐷1 = \frac{𝐷_0}{𝑐𝑜𝑠2𝛿}\)

Donde \(𝐷_0\) es la difusividad hidráulica saturada (\(𝐷0 = \frac{𝐾_𝑠}{𝑆_𝑠}\), donde \(𝑆_𝑠\) es el almacenamiento especifico). 𝑁 es el número total de intervalos y \(𝐻(𝑡 − 𝑡_𝑛)\) es la función de paso de Heaviside, donde \(𝑡_𝑛\) es el tiempo en el n-ésimo intervalo en la secuencia de infiltración de lluvia. La función \(𝑖𝑒𝑟𝑓𝑐\) tiene la forma \(𝑖𝑒𝑟𝑐𝑓 (𝜂) = 1√𝜋exp(−𝜂2) − 𝜂 𝑒𝑟𝑓𝑐 (𝜂)\), donde \(𝑒𝑟𝑓𝑐(𝜂)\) es la función de error complementario.

Para estimar los cambios en el factor de seguridad como función de la profundidad (\(Z\)) y el tiempo (\(t\)) se utiliza el modelo de estabilidad de talud infinito de acuerdo con la siguiente ecuación:

\(FS(Z,t)=\frac{tan}{tanδ}+\frac{c'-ψ(Z,t)𝛾_𝑤tanϕ'}{𝛾_sZsenδcosδ}\)

Donde \(𝑐′\) es la cohesión efectiva del suelo, \(𝜙′\) el ángulo de fricción efectivo, \(𝛾𝑤\) el peso unitario del agua, \(𝛾_𝑠\) el peso unitario del suelo y \(𝜓(𝑍,𝑡)\) la cabeza de presión en función de la profundidad y el tiempo \(𝑡\).

modelos

Fig. 66 Modelo TRIGRS.#

Tutorial#

Elaborado por Maria Isabel Arango.

Ingrese a la página web TRIGRS y descargue TRIGRS de acuerdo a su sistema operativo. Descomprima el archivo .zip y ábralo. Adicionalmente, puede descargar el manual del usuario donde encontrará toda la base teórica detrás del modelo, ilustraciones y un tutorial.

Como verá, la carpeta TRIGRS contiene varias subcarpetas y archivos en formato txt. TRIGRS y sus programas complementarios (TopoIndex, GridMatch y UnitConvert) no funcionan con una interfaz, sino que se ejecutan en ventanas de entrada y salida simples. Cada programa funciona a través de un ejecutable (exe), los cuales están ubicados en la subcarpeta bin”. A su vez, cada ejecutable utiliza un archivo de inicialización que contiene los datos básicos necesarios para ejecutar el programa, así como los nombres de otros archivos de entrada.

En este taller, se trabajará sobre una zona de estudio ubicada en el municipio de San Carlos (Antioquia-Colombia) correspondiente a una parte de la cuenca de la quebrada La Arenosa, donde el 21 de septiembre de 1990 se presentó un evento de lluvia extremo que causó cientos de deslizamientos, muchos de los cuales alcanzaron los cauces y causaron una avenida torrencial. En este tutorial se implementará el modelo TRIGRS para identificar las zonas inestables usando los parámetros de la zona de estudio y bajo el escenario de precipitación del 21 de septiembre de 1990. La información requerida para llevar a cabo el análisis de estabilidad con TRIGRS incluye:

  • Modelo de elevación digital (DEM)

  • Mapa de pendientes

  • Mapa de dirección de flujo

  • Mapa de profundidad del suelo

  • Propiedades mecánicas e hidrológicas del suelo. En caso de tener varias zonas con propiedades homogéneas, mapa de zonificación del suelo

  • Escenario de precipitación. En caso de ser una precipitación heterogénea espacialmente, mapa espacialisado de las intensidades.

En todos los casos, los mapas que se utilicen deben ser de tipo raster en formato ASCII y deben cumplir con los siguientes criterios:

  • Deben estar en el mismo sistema de coordenadas

  • Deben tener las mismas dimensiones

  • Deben estar unidades metro-kilogramo-segundo

Preparación de los mapas#

En la carpeta Mapas Base, encontrará los siguientes archivos raster. Para descargar los archivos seleccionelo, cuando Github visualice dicho archivo en la parte superior derecha seleccione Raw y con el boton derecho Save as…. Guardelos como archivos ASCII eliminando la terminación txt que sale por defecto.

  • dem_arenosa_clip: modelo de elevación digital

  • pendiente_arenosa_clip: pendiente en grados

  • flowdir_arenosa_clip: dirección de flujo en formato ESRI

  • profsuelo_arenosa_clip: mapa con la profundidad del suelo en metros

direccion

Fig. 67 Códigos usados por ESRI para la dirección de flujo Profsuelo_Arenosa_Clip: Profundidad del suelo, en metros de la zona de estudio. Suelos_Arenosa_Clip: Zonificación de los tipos de suelo presentes en la zona de estudio.#

Busque el archivo dem_arenosa_clip.asc y ábralo directamente desde el explorador para verlo en el editor de texto. Las primeras 6 líneas forman parte de la información base con la que cuenta todo archivo ASCII: número de columnas, de filas, coordenadas de la primera columna, coordenadas de la primera fila, tamaño de celda y valor para los datos nulos. Todos los archivos tipo ASCII que le ingrese a los archivos ejecutables deben tener valores idénticos de estos parámetros, de lo contrario genera error. Abra los archivos pendiente_arenosa_clip, flowdir_arenosa_clip y profsuelo_arenosa_clip, y compruebe que los parámetros de todos los archivos sean idénticos.

En GIS abra el raster suelos_arenosa_clip.tif. Explore en la tabla de atributos a cuáles unidades de suelo corresponden las categorías. Ahora, es necesario reclasificar el raster para que el modelo lo entienda, asignándole a cada unidad de suelo un número. En caso de utilizar ArcToolbox, vaya a Spatial Analysis > Reclass > Reclassify. En el campo Input raster agregue el raster de suelos, en el campo Reclass Field seleccione Simb_Sue, y reclasifique el raster de la siguiente forma. Guárdelo.

Tabla 1. Codificación de los tipos de suelo.

Valores Antiguos

Valores Nuevos

POc1 (Aluviales)

1

YAe1 (Suelo Residual)

2

YAf2 (Suelo residual)

2

NoData

NoData

De esta forma, las unidades de suelo quedan divididas en dos categorías: depósitos aluviales (1) y suelos residuales del Batolito Antioqueño (2). Cada una de estas unidades tiene unas propiedades mecánicas e hidrológicas homogéneas.

Convierta el raster reclasificado a formato ASCII (Conversion Tools > From Raster > Raster to ASCII). Guárdelo con el nombre suelos_arenosa_reclass_clip.asc. Ábralo con el editor de texto y revise que los parámetros coincidan con los del resto de mapas.

Ejecución de los programas#

Una vez tenga lista la información de entrada, debe ejecutar los programas auxiliares y finalmente TRIGRS, que se encuentran en la carpeta bin. Los programas auxiliares tienen como objetivo preparar la información ingresada y crear capas intermedias que son necesarias para el cálculo final. Estos programas deben ejecutarse siempre en el siguiente orden, ya que normalmente, los archivos de salida de uno son los archivos de entrada del siguiente.

  • Unit Converter (Opcional): En caso de tener archivos en unidades distintas a las aceptadas por el modelo, este programa permite convertirlas (pfs a Pa, pies a metros, etc).

  • Grid Match: Permite corroborar que todos los archivos ASCII contienen raster con el mismo número de celdas, coordenadas y tamaño de celda.

  • Topo Index: Utiliza el DEM y la dirección de flujo para definir el patrón de distribución del flujo, calcular los factores de ponderación para distribuir la escorrentía de la superficie y calcular los tamaños de la matriz que utilizará TRIGRS para sus cálculos.

  • TRIGRS: Una vez preparados los mapas y generado las capas intermedias, se corre el ejecutable final que realiza los cálculos de presión de poros y factor de seguridad para los instantes de tiempo que usted determine.

Cada ejecutable que se encuentra en la carpeta bin contiene un código que se ejecuta cada vez que se dá doble click sobre el archivo .exe. Estos códigos necesitan para su correcto funcionamiento un archivo de entrada, el cual contiene los parámetros con los que debe correr el código, y es el archivo que usted como usuario debe modificar. Los archivos de entrada los encuentra en la carpeta TRIGRS con los siguientes nombres:

  • uc_in.txt. archivo de inicio del Unit Coverter

  • gm_in.txt. Archivo de inicio del Grid Match

  • tpx_in.txt. archivo de inicio de Topo Index

  • tr_in.txt. Archivo deinicio de TRIGRS

Cada vez que se ejecuta los .exe, este buscará su archivo de inicio en la misma carpeta donde está el archivo .exe. Por lo tanto, organice todos los archivos de inicio de la carpeta TRIGRS en la misma carpeta con los archivos .exe. Trate de conservar los archivos de inicio originales en la carpeta TRIGRS siempre, ya que le serán de utilidad para detectar errores en el archivo que está editando. También, evite cambiar los nombres de los archivos de inicio, ya que esto genera fácilmente errores.

Dado que toda la información está en formato metro-kilogramo-segundo, no se ejecutará Unit Converter. Siga los siguientes pasos para ejecutar el resto de programas:

Grid Match#

Grid Match permite asegurarse que todos los archivos tipo ASCII con los que trabajará tienen las mismas dimensiones y no tienen condiciones que generen error. Estas condiciones incluyen también que un mapa tenga una celda nula en la misma ubicación donde otros mapas tienen información.

  • Abra el archivo gm_in.txt. Los archivos de entrada tienen un formato de fácil interacción, donde cada parámetro está explicado en su parte superior. Los archivos ya vienen con información cargada a modo de ejemplo, así que usted debe modificar los parámetros según su caso de estudio. También debe separar con comas los parámetros tal cual están separadas sus explicaciones. (Por ejemplo, fíjese en la línea 3 que dice: rows, columns. Eso significa que en la línea 4 usted debe reemplazar el 10, 10 por las dimensiones del DEM con el que trabajará, con el número de filas separado por una coma del número de columnas. Los espacios son indiferentes, la coma designa la separación entre dos parámetros, y el punto es el separador de decimales).

  • En la línea 2, escriba el número de mapas que desea corroborar. En nuestro caso son 5 (DEM, pendiente, dirección de flujo, profundidad del suelo y tipos de suelo).

  • En la línea 4, escriba las dimensiones de los mapas que va a ingresar. Lea estos valores en cualquiera de los archivos ASCII de entrada como nrows y ncols.

  • En las líneas 6 en adelante, usted debe enumerar uno por uno los archivos ASCII que desea corroborar con su ruta de acceso y cada uno en una línea. Ingrese primero el DEM.

  • Para este ejercicio, seguirá la opción más simple de ingreso de los archivos. Para esto, copie y pegue los archivos ASCII de información base (dem_arenosa_clip.asc, pendiente_arenosa_clip.asc, flowdir_arenosa_clip.asc, profsuelo_arenosa_clip.asc, y suelos_arenosa_reclass_clip.asc) en la carpeta donde tenga los ejecutables .exe(sólo copie y pegue los archivos .asc).

  • su archivo gm_in debe lucir así:

gm

Fig. 68 Archivo gm_in.txt#

  • Una vez termine de editar el archivo gm_in, guarde los cambios.

  • Ejecute gridmatch.exe.

  • Mientras el programa se ejecuta, aparecerá una ventana negra mostrando el avance. Dado que el Grid Match sigue un proceso simple, normalmente su tiempo de ejecución no es mayor a dos segundos, después de los cuales la ventana se cerrará, en caso de que no hayan errores. En caso de que exista algún error con los mapas, la ventana le mostrará de qué se trata el error (le dirá cuál mapa tiene problemas y cuál es el error).

  • Cada vez que un ejecutable termina de realizar su proceso, genera un archivo txt de resumen donde especifica tiempo de ejecución, resultados y ubicación de los archivos de salida. Este archivo resumen se genera en la misma carpeta donde está el ejecutable, y tiene el nombre del ejecutable + Log (GridMatchLog, TopoIndexLog, etc.). Busque y abra el archivo GridMatchLog.txt.

  • Dentro del archivo, usted verá en las líneas 1-6 información sobre la corrida y en las líneas 7 -8 la copia del archivo de inicio. En las líneas 20-57, realiza un resumen de cada uno de los mapas que se revisaron, especificando el número de errores que se encontraron en cada uno. Finalmente, le muestra una confirmación de que Grid Match terminó normalmente, con lo que ya puede estar seguro de que los datos de entrada están libres de errores.

Topo index#

Topo Index permite calcular los tamaños de matriz para TRIGRS y preparar datos de enrutamiento de escorrentía. Utiliza el DEM y las direcciones de flujo para determinar el orden correcto para los cálculos de enrutamiento de escorrentía y para calcular los factores de ponderación que determinan cómo se distribuye el exceso de agua a las celdas de las cuadrículas adyacentes.

  • Abra el archivo txp_in.txt

  • En la línea 2, ingrese el nombre del proyecto (La Arenosa Clip). Este nombre puede incluir letras, números y espacios.

  • Los primeros dos valores de la línea 4 corresponden al número de filas y columnas de los archivos ASCII que se ingresarán (83, 74). El tercer valor corresponde al formato de codificación que tiene el mapa de dirección de flujo. Dado que el mapa que se le entregó está en formato ESRI, el valor correspondiente es 1.

  • El primer valor de la línea 6 corresponde al exponente (ω) que determina el factor de ponderación para la distribución de escorrentía en las celdas adyacentes. Usted puede darle un valor al exponente:

Exponente

Patrón de distribución de la escorrentía

ω < 0

Estrecho, el flujo sólo se distribuye a una o dos celdas adyacentes

ω = 0

En forma de abanico, distribución uniforme a celdas adyacentes

0 < ω <= 20

En forma de abanico, la distribución favorece a celdas inferiores

ω > 20

Estrecho, escorrentía se distribuye por zonas escarpadas

En la página 42 del manual del usuario, usted puede consultar las fórmulas de las ponderaciones según los valores que toma el exponente ω. Para este ejercicio, trabaje con un valor de -1.

  • El segundo valor en la línea 6 especifica el número máximo de iteraciones que debe usar Topo Index para crear el orden para el cálculo de la escorrentía. Cuando su DEM tiene errores hidrológicos, Topo Index arreglará el error y empezará una nueva iteración. Este proceso terminará cuando el enrutamiento converja, o cuando se alcance el número máximo de iteraciones que usted defina. Para este ejercicio, defina el número de iteraciones a 100.

  • Recomendación: Empiece con un número estándar de 100 iteraciones, y si el enrutamiento no converge, llévelas a 200. Si finalmente su archivo no converge, probablemente deba reparar su DEM en ArcGis realizando un Fill.

  • En las líneas 8 y 10 ingrese, respectivamente, la ruta de los archivos del DEM y de la dirección de flujo (siga los lineamientos de la Nota 2)

  • Las líneas 11, 13, 15, 17 y 19 le hace preguntas específicas sobre si desea guardar archivos intermedios del procesamiento. A estas preguntas usted debe contestar con T = Si (True), o F = No (False). Dado que estos archivos son necesarios para ejecutar TRIGRS más adelante, responda con T a todas ellas.

  • En la línea 22, debe ingresar la ruta de la carpeta donde desea guardar los archivos intermedios que se generarán después de correr el Topo Index. Puede crear una subcarpeta dentro de la carpeta bin que se llame intermedios, y especificar la carpeta. Tenga en cuenta asignar al final el backslash ().

  • Finalmente, en la línea 24 usted debe ingresar un código que se asignará a los archivos generados.

  • Recomendación: Asígnele un código numérico que designe el lugar y el número de la corrida. En este caso, asigne el código Arenosa_01

  • Cuando termine de editar el archivo tpx_in, guarde los cambios. Su archivo debe lucir así:

gm

Fig. 69 Archivo gm_in.txt#

  • Ejecute el programa TopoIndex.exe.

  • Cuando el proceso finalice, busque y abra el archivo TopoIndexLog.txt. El archivo resumen incluye, igual que el GridMatchLog, información del tiempo de ejecución y una copia del archivo de inicio. A partir de la línea 34, se puede ver un resumen del proceso de corrección de errores en cada iteración. En este caso, sólo tomaron 3 iteraciones para que el programa calculara las rutas de escorrentía arreglando los errores que fue encontrando. En la línea 55 se encuentra la información del número de celdas en el cálculo de escorrentía, y a partir de la línea 57 se encuentra la confirmación de que no hubo problemas durante el proceso de cálculo.

  • Abra la carpeta “intermedios”. En ella, podrá ver los siguientes archivos:

Archivo

Descripción

TIdsneiList_Arenosa.txt

Lista de celdas receptoras en pendiente descendente en la dirección del flujo D8 (2)

TIdscelGrid_Arenosa.txt

Raster de celdas receptoras en pendiente descendente en la dirección del flujo D8 (1)

TIcelindxGrid_Arenosa.txt

Raster de índice de celda que especifica el orden de cálculo para el enrutamiento de escorrentía en TRIGRS

TIcelindxList_Arenosa.txt

Lista de números de celda y número de índice correspondiente (orden de cálculo)

TIflodirGrid_Arenosa.txt

Raster de dirección de flujo

TIdscelList_Arenosa.txt

Lista de celdas en pendiente descendente para las cuales se han calculado factores de ponderación distintos de cero (3)

TIwfactorList_Arenosa.txt

Lista de factores de ponderación para las celdas receptoras en pendiente descendente (4)

Los números en paréntesis denotan el orden en el que se deberán ingresar estos archivos en el archivo de inicio de TRIGRS, que verá más adelante.

TRIGRS#

TRIGRS se basa en la suposición de la infiltración unidimensional (únicamente hacia abajo) en un campo de flujo constante para determinar los cambios en la presión de poros en un modelo de talud infinito y la estabilidad de los suelos poco profundos durante las tormentas. TRIGRS incluye modelos para condiciones iniciales saturadas y no-saturadas. Para condiciones iniciales saturadas, cada simulación con TRIGRS analiza un único punto en el tiempo especificado por el usuario durante una secuencia de tormenta. Por lo tanto, el usuario ejecuta una serie de simulaciones para diferentes tiempos específicos utilizando TRIGRS para estudiar el historial de tiempo de inestabilidad de la ladera durante una secuencia de lluvia. Para condiciones iniciales no-saturadas, TRIGRS puede analizar una serie de tiempo completa durante una única simulación.

Abra el archivo tr_in.txt. Este archivo contiene varios tipos de información que debemos ingresar:

  • En este taller, se modela con condiciones saturadas y luego con condiciones no saturadas, para compararar los resultados. Primero, realice la modelación usando los parámetros sugeridos en la columna Valor Arenosa Sat. Después, la modelación usando los parámetros sugeridos en la columna Valor Arenosa no Sat

  • En la línea 2, debe ingresar el nombre del proyecto. Defínalo como Arenosa Clip

  • En los primeros cuatro valores de la línea 4 debe ingresar los datos de tamaño de matriz y dimensiones del cálculo de escorrentía. Copie estos valores del archivo TopoIndexLog (Data cells, Rows, Columns, Downslope cells) que corresponden, respectivamente, a los campos imax, row, col, nwf en el archivo de inicio de TRIGRS. Los parámetros tx, nmax, y todos los demas en las lineas siguientes se explican a continuación:

Parámetro

Descripción

Rango de valores

Valor Arenosa Sat

Valor Arenosa No Sat

tx

En combinación con nper, determina cuántos intervalos de tiempo se utilizan en los cálculos. El esfuerzo computacional y los requisitos de almacenamiento aumentan directamente con tx. En condiciones saturadas, se puede establecer en 1 para permitir que TRIGRS calcule una solución lo más rápidamente posible durante un tiempo transcurrido en particular.

tx≥1

1

10

nmax

Número máximo de raíces de la ecuación que calcula la conductividad hidráulica de la zona no saturada (Ecuación 6 del manual del usuario). Si se ingresa un valor negativo, TRIGRS asumirá que trabajará en condiciones saturadas.

nmax≥1

-1

30

nzs

Determina el número de incrementos verticales en que será dividida la capa de suelo para calcular presión de poros y factor de seguridad. Mayor número de incrementos dará como resultado un análisis más detallado en profundidad, pero aumenta de forma considerable el tiempo de cálculo.

nzs≥1

10

10

mmax

Un valor negativo indica el cálculo de presiones de poros siguiendo la Ecuación 1 del manual del usuario, la cual utiliza un modelo de profundidad infinita. Un valor positivo indica el cálculo de presiones de poros siguiendo la Ecuación 2 del manual del usuario, la cual asume una superficie impermeable a una profundidad finita. El valor positivo indica el valor máximo de la Ecuación 2. Ambas ecuaciones se usan únicamente en condiciones saturadas.

\(10<mmax<25\)

20

Indiferente, use cualquier valor

nper

Número de periodos de lluvia

nper>1

7

7

zmin

Mínima profundidad para calcular presiones de poros.

0.001≤zmin≤0.1

0.001

0.001

uww

Peso unitario del agua (9800 kg-m)

uwwr>0

9.8e3

9.8e3

t

Tiempo, en segundos, que desea modelar. En este ejemplo se utilizarán 7 h

nper>0

25200

25200

zones

Número de zonas de suelo con propiedades mecánicas homogéneas

zones>1

2

2

zmax

Profundidad máxima de suelo en que se calcularán presiones de poro y el factor de seguridad. Para no gastar recursos computacionales de forma innecesaria, determine esta profundidad como la de probable falla. Si es un valor homogéneo para toda la zona de estudio, ingrese el valor positivo de la profundidad, en metros. En caso de que sea un valor heterogéneo espacialmente, ingrese un valor negativo y más adelante podrá ingresar un mapa con los valores de la profundidad. En el caso de la Arenosa, esta superficie de falla coincide con la de profundidad del suelo, la cual es heterogénea.

\(czmax> zmin\); si \(zmax<0\), leer valores del raster

-1

-1

depth

Profundidad inicial del nivel freático. Si es un valor homogéneo en toda la zona de estudio, ingrese el valor de la profundidad, en metros. Si es un valor heterogéneo espacialmente, ingrese un valor negativo y más adelante podrá ingresar un mapa con los valores de profundidad. Para la Arenosa, asumirá que en condiciones no saturadas el nivel freático se encontraba a 3.5 metros de profundidad.

\(0<depth≤zmax; si depth<0, leer valores del raster\)

0.5

3.5

rizero

Tasa de infiltración del suelo antes de la tormenta. Si es un valor homogéneo en toda la zona de estudio, ingrese el valor de metros/segundo. Si es un valor heterogéneo espacialmente, ingrese un valor negativo y más adelante podrá ingresar un mapa con los valores.

\(0<rizero, si rizero<0, leer valores del raster\)

6.94e-6

6.94e-6

Min_Slope_Angle

Ángulo de inclinación mínimo; no hay cálculos de presión de poro o factor de seguridad en pendientes inferiores a Min_slope_angle.

\(Min_slope_angle ≥0 \)

0

0

Propiedades físicas del suelo#

Las propiedades de las unidades del suelo se deben ingresar como lo muestran las líneas 9-11 del archivo tr_in. En la primera línea se debe ingresar el número que asignó a la unidad de suelo, y en la segunda y tercera línea están las propiedades del suelo correspondientes. En el archivo de inicio por defecto ya están creados los campos para dos unidades de suelo, sin embargo, si usted está trabajando con una zona de estudio que incluya más unidades, debe crear los campos para cada uno de forma consecutiva en los campos inferiores. En este taller trabajará con las dos unidades de suelo presentes en el raster que se le entregó. Recuerde que a la unidad de depósitos aluviales (POc1) se le asignó el número 1, y a los suelos residuales de Batolito Antioqueño (YAE1-YAF2), se les asignó el número 2. Ingrese los parámetros recomendados a continuación, los cuales no varían en ambas modelaciones (condición saturada y no saturada), excepto por el ultimo parámetro, alpha. Recuerde también conservar siempre las unidades Kg-m-s.

Parámetro

Descripción

POc1

YAE1-YAF2

cohesion

Cohesión del suelo en Pa

1000

5000

phi

Angulo de fricción del suelo en °

34

24

uws

Peso unitario saturado en N/m3

20000

18000

diffus

Difusividad hidráulica (m2/s)

5.33e-4

2.18e-3

K-sat

Conductividad hidráulica vertical del suelo saturado (m/s)

5.33e-6

2.18e-5

Theta-sat

Contenido volumétrico de agua de suelo saturado

0.48

0.46

Thera-res

Contenido volumétrico residual de agua

0.18

0.18

Alpha

Parámetro de distribución del tamaño del suelo (1/m). Un valor positivo indicará a TRIGRS el uso de la solución no saturada, un valor negativo indica el uso de la solución saturada.

2.3 en modelación no saturada, -2.3 en modelación saturada

2.3 en modelación no saturada, -2.3 en modelación saturada

Intensidad del escenario de lluvia e incrementos de tiempo#

La siguiente figura muestra el registro de precipitación obtenido de un pluviómetro ubicado en la cuenca de La Arenosa el día del evento.

lluvia

Fig. 70 Registro de precipitación del evento del 21/09/1990 en la cuenca La Arenosa#

La tormenta comenzó a las 8:00 pm del día del evento y terminó a las 2:00 am del día siguiente. Usted modelará la estabilidad de la cuenca durante las 7 horas de duración de la tormenta (7 horas= 25.200 segundos, que corresponde al valor ingresado en el campo t de los parámetros de modelación). Como verá el registro de precipitación está dado en escala horaria (7 horas de tormenta con registro de intensidad horaria corresponde al valor ingresado en el campo nper de los parámetros de modelación).

  • En la línea 16 del archivo tr_in usted debe ingresar las diferentes intensidades de cada periodo de tiempo registrado, en m/s. (Ojo: Si usted ingresó el número 7 en el campo nper de los parámetros de simulación, debe ingresar exactamente 7 valores de intensidad de lluvia).

  • En la línea 18, debe ingresar los intervalos de tiempo, en segundos, a los que corresponde cada valor de intensidad que ingresó anteriormente. Dado que el primer valor debe ser cero, esta línea debe incluir el número de intervalos + 1 valores. (Si usted ingresó el número 7 en nper, debe incluir 8 valores en los intervalos de tiempo) .

Para el caso de la Arenosa, incluya los siguientes valores de intensidades de lluvia e intervalos de tiempo.

Capt (s)

Cri (m/s)

0

3600

1.58e-05

7200

1.58e-05

10800

2.50e-05

14400

3.60e-06

18000

2.22e-06

21600

1.11e-06

25200

1.11e-06

En el caso de que la intensidad de precipitación no esté como un valor homogéneo en el espacio sino como mapas de intensidades, usted debe ingresar en los campos Cri(n) valores negativos (tantos valores como periodos de tiempo nper) para dar a entender que va a ingresar más adelante mapas de intensidad de precipitación.

Archivos de entrada#

En las siguientes líneas, debe ingresar las rutas de los archivos que entrada que TRIGRS le requiera. Recuerde que en la sección de ingreso de parámetros de computación, cuando una propiedad era heterogénea en el área de estudio, usted podía ingresar un valor negativo, que indica que los valores se deben leer de un raster. En caso de que la propiedad fuera homogénea para toda la zona de estudio, usted podía ingresar un valor positivo, que se tomaba como constante en todo el área.

Ahora, cuando usted indique un valor negativo en la sección de parámetros de computación, debe ingresar la ruta del respectivo raster en la línea donde TRIGRS lo requiera. Cuando usted haya ingresado un valor positivo homogéneo para toda el área de estudio, en el lugar donde TRIGRS le requiera el mapa correspondiente, usted debe ingresar la palabra none, para indicar que es un valor homogéneo.

  • En las líneas 20, 22 y 24 ingrese, respectivamente, las rutas de los archivos pendiente_arenosa_clip.asc , suelos_arenosa_reclass_clip.asc y profsuelo_arenosa_clip.asc.

  • En la línea 26 debe ingresar el mapa de profundidad del nivel freático. Dado que usted asumirá para este ejercicio un nivel freático constante de 0.5 metros para condiciones saturadas y 3.5 metros para condiciones no saturadas, y que ya ingresó los valores positivos en los parámetros de computación, escriba la palabra none. Haga lo mismo en la línea 28, donde debe ingresar el mapa de tasa de infiltración del suelo, la cual usted ya definió con un valor homogéneo.

  • En la línea 29 se le pide que ingrese los mapas de intensidad de precipitación. Dado que usará los valores constantes que ya ingresó, escriba 7 veces la palabra none, cada una en una línea independiente.

  • En las líneas siguientes, debe ingresar las rutas de los archivos generados por TopoIndex. Tenga en cuenta el orden en el que los debe ingresar, el cual está escrito en paréntesis en la Tabla.

Archivos de salida#

En las siguientes líneas, debe especificar las opciones de los archivos de salida de TRIGRS.

  • Debe especificar la carpeta donde desea que TRIGRS guarde los archivos de salida. Deje la línea en blanco si desea colocar los archivos en la misma carpeta donde está el archivo de programa ejecutable TRIGRS. En este caso, cree la carpeta Resultados dentro de la carpeta bin, y especifíquela como la carpeta donde desea que se guarden los archivos de salida. No olvide siempre terminar con backslash ()

  • En la siguiente línea, debe especificar el código que tendrán los archivos de salida, defínalo como Sat para la modelación saturada, y como No_Sat para la modelación no saturada.

  • En las siguientes líneas, le hace preguntas específicas sobre si desea guardar archivos de salida. A estas preguntas usted debe contestar con T = Si (True), o F = No (False). Deje las respuestas por defecto que ya están en el archivo.

  • En una de las preguntas debe especificar si desea que TRIGRS genere una lista detallada de la presión de poros y el factor de seguridad a profundidad para cada celda del raster. Para este taller configure el indicador de salida en -2 para generar una lista completa. Cuando trabaje raster grandes, es mejor desactivar la salida de la lista detallada estableciendo el indicador en 0 para ahorrar espacio en el disco.

  • Después de las preguntas, debe ingresar el número de pasos intermedios en los que desea generar archivos de salida y los intervalos de tiempo en los que deben ser generados. Por ejemplo, si a usted solo le interesa conocer el factor de seguridad y presión de poros al final de la tormenta, debe ingresar el número 1 en el número de veces y en los intervalos de tiempo debe ingresar el tiempo final de la tormenta, en segundos. Si desea conocer al factor de seguridad en la mitad de la tormenta y en el final, ingrese el número 2 y los tiempos de la mitad y final de la tormenta. Para este ejercicio, genere mapas al inicio de la tormenta y después cada hora para conocer las variaciones horarias de la presión de poros y la estabilidad. Ingrese el número 8 en las veces que desea generar mapas, e ingrese los segundos correspondientes a cada hora (0, 3600, 7200, 10800, 14400, 18000, 21600, 25200).

  • En la pregunta de si desea saltarse otros pasos de tiempo, responda F.

Otras opciones de usuario#

Existen otras opciones de usuario explicadas a continuación. Defina los valores para este ejercicio como se indica en paréntesis:

  • Opción de usar una solución analítica para la porosidad rellenable permite al usuario elegir la opción “Solución analítica eficiente” que converge rápidamente para tiempos intermedios y posteriores, o una integración numérica para tiempos muy pequeños (T)

  • Opción de estimar presión de poros positiva en la zona de nivel freático en aumento permite guardar en el disco ya sea la presión de poros positiva debajo del nivel freático ascendente (ingrese T) o los cálculos de la presión de poros negativa sobre el nivel freático inicial (ingresar F). (T)

  • La opción psi0=-1/alpha: Permite a TRIGRS usar los cálculos no-saturados en toda la zona sobre el nivel freático inicial, o solo sobre la zona de ascendencia capilar. El uso de la segunda opción (solo por encima de la franja capilar) tiende a producir un aumento más rápido del nivel freático, debido a que los modelos de infiltración no-saturada operan en una capa más delgada, pero pueden producir una mejor aproximación de la curva característica del suelo-agua. (F)

  • La opción de registrar resultados de balance de masa permite verificar el comportamiento de la solución de zona insaturada para combinaciones particulares de parámetros. Esta opción es principalmente para uso en las primeras etapas del modelado en celdas de prueba. El registro del balance de masa durante las simulaciones en raster grandes ralentizará los cálculos y aumentará en gran medida los requisitos de almacenamiento en disco. (T)

  • La opción de dirección del flujo: Permite al usuario anular la pendiente predeterminada de la línea que determina los valores máximos permitidos de las presiones de poro calculadas. (gener) Cuando termine de ingresar los parámetros del archivo de inicio, guarde los cambios y dé doble clic sobre el ejecutable TRIGRS.exe.

Visualizar los resultados#

Mientras que TRIGRS se ejecuta, aparecerá una ventana negra mostrando el avance de los cálculos. Dependiendo del tamaño del raster que esté analizando, este proceso podría tomar desde minutos hasta días. Cuando el proceso termina, la ventana se cierra, en caso de que no haya errores. En caso de que exista algún error con los mapas o el archivo de inicio, la ventana le mostrará de qué se trata el error y la línea que debe revisar en su archivo de inicio. Cuando terminen los cálculos, abra el archivo TrigrsLog.txt y examínelo.

Abra la carpeta Resultados. Dentro de la carpeta se encuentran los siguientes archivos de salida:

Archivo

Descripción

TRrunoffPer”x”Sat.txt

Raster de escorrentía computada durante un período de tormenta dado; x designa el período de tormenta

TRfs_min_Sat__”x”.txt

Raster del factor de seguridad mínimo calculado en cada celda por tiempo, x denota el periodo de la tormenta.

TRz_at_fs_min_Sat__”x”.txt

Raster de profundidades en la que se calculó el factor de seguridad mínimo para el tiempo del intervalo x.

TRp_at_fs_min_Sat__ ”x.txt

Raster de presión de poros en la profundidad correspondiente al factor mínimo de seguridad en el intervalo x.

TRinfilratPer”x”Sat_.txt

Raster de tasas de infiltración reales durante un período de tormenta dado x.

TRlist_z_p_fs_Sat_.txt

Lista celda por celda de presión y factor de seguridad en cada incremento de profundidad.

TRnon_convrg_SZ_Sat_.txt o TRnon_convrg_UZ_Sat_.txt

Archivos de raster que localizan celdas no convergentes. Se generan automáticamente si los cálculos de series infinitas en cualquier celda no convergen. UZ denota zona insaturada, SZ denota zona saturada.

saturado

Fig. 71 Archivo de inicio de TRIGRS para condiciones saturadas#

Factor detonante sismo#

Existen diferentes métodos para evaluar la estabilidad de laderas ante sismos y estimar los desplazamientos inerciales permanentes. A continuación se presentarán el denominado análisis seudoestático y el análisis de desplazamientos permanentes.

Análisis seudoestático#

El método seudoestático es una generalización del análisis de estabilidad por equilibrio límite, en el cual los efectos de los sismos son representados por una fuerza estática equivalente cuya magnitud es producto del coeficiente sísmico k y el peso de la masa que potencialmente se va a deslizar. Este método evalúa los efectos sísmicos como una fuerza permanente unidireccional, mediante la simulación del incremento de las fuerzas inerciales debido a la aceleración de un sismo, suponiendo que las fuerzas sísmicas son proporcionales al peso de la masa deslizante. Usualmente, solo se adiciona la fuerza sísmica en la componente horizontal en el análisis de equilibrio límite y se consideran los efectos de las fuerzas verticales cercanos a cero. Para calcular el Factor de Seguridad (fs), Matasovic propone la ecuación (1):

\(FS=\frac{\frac{c}{\gamma_sZcos^2\theta}+[1-\frac{\gamma_w(Z-h_w)}{\gamma_s+tan\theta}]-K_htan\theta tan\phi}{K_h+tan\theta}\)

donde \(c\) es la cohesión, \(\phi\) es el ángulo de fricción efectiva, \(\gamma_s\) es el peso unitario de los materiales, \(\theta\) es el ángulo de la pendiente, \(z\) es el espesor del estrato deslizable, \(\gamma_w\) es el peso unitario del agua, \(k_h\) es el coeficiente horizontal seudoestático y \(h_w\) es la distancia al nivel freático medido desde la superficie.

En el título H de la NSR10 se describe el diseño y estabilidad de taludes teniendo en cuenta los efectos sísmicos mediante el análisis seudoestático. Esta metodología plantea determinar el coeficiente sísmico horizontal seudoestático (kh) en función de la aceleración máxima del terreno (Amax), obtenida del espectro de diseño para el periodo de vibración cero establecida por el estudio de microzonificación sísmica del sitio. Para Colombia es posible utilizar el espectro de aceleración superficial para el periodo de retorno de 475 años y amortiguamiento del 5%, tomado de los estudios de microzonificación sísmica donde esten disponibles. La NSR-10 define las aceleraciones máximas de acuerdo con el espectro de diseño, donde presenta la formulación de la aceleración máxima del terreno para un periodo de vibración igual a cero, dada por la ecuación:

\(K_{ST}=k_h=0.8A_max\)

donde \(A_a\) es la aceleración pico efectiva de diseño, \(F_a\) es el coeficiente de amplificación que afecta la aceleración en la zona de períodos cortos e I representa el coeficiente de importancia. Estas variables fueron obtenidas a partir de los espectros de respuesta a nivel de superficie del terreno realizados. El coeficiente de importancia se determina a partir del uso de las edificaciones expuestas; una posible propuesta es utilizar el valor unitario, correspondiente al grupo de importancia I.

Grupo de uso

Coeficiente de importancia I

IV

1.55

III

1.25

II

1.10

I

1.0

A partir de \(A_max\) se obtiene el \(k_h\) con la siguiente ecuación teniendo en cuenta los criterios para el análisis seudoestático dados en la tabla de la NSR-10 [50, H.5.2.5]. POr ejemplo el valor 0,8 representa Suelos, enrocados y macizos rocosos muy fracturados (RQD1<50%):

\(K_{ST}=k_h=0.8A_max\)

Desplazamientos permanentes#

Este método ha dado útiles resultados modelando el comportamiento dinámico de las laderas a nivel regional, aunque no necesariamente predice los desplazamientos medidos, sino que establece un índice útil de cómo una ladera se podría comportar durante un evento sísmico. El método más conocido es el denominado método de Newmark. Múltiples variaciones del método de Newmark han sido propuestas para obtener resultados precisos en términos de desplazamientos al modelar rigurosamente la respuesta dinámica de la ladera. El método requiere evaluar el factor de seguridad, la aceleración crítica del terreno y la aceleración máxima en roca (Peak Ground Acceleration o PGA). Jibson et al. [1] proponen un modelo de equilibrio límite de talud infinito en material friccionante y cohesivo para obtener el factor de seguridad a partir de la ecuación:

\(FS=\frac{c}{\gamma_s Z sen\theta}+\frac{tan\phi}{tan\theta}-\frac{m \gamma_w tan\phi}{\gamma_s tan\theta}\)

donde \(z\) es el espesor del suelo y \(m\) es la proporción del suelo que está saturado. Newmark señala que la aceleración crítica \(a_c\) de un bloque que potencialmente se puede deslizar es función del factor de seguridad estático y de la geometría del bloque expresado en la ecuación:

\(a_c=(FS-1)sen\theta\)

donde \(a_c\) es la aceleración crítica en términos de la aceleración de la gravedad \(g\), y el ángulo que forma la horizontal con el centro de masa del bloque, que generalmente se aproxima a la pendiente. Finalmente, a partir de \(a_c\) y la aceleración máxima en roca (PGA) se calculan los desplazamientos de Newmark a través de una regresión logarítmica estimada por Jibson:

\(LogDn=0.215+log[(1-\frac{a_c}{PGA})^{2.341} (\frac{a_c}{PGA})^{-1.438}]+/-0.510\)

Modelos probabilísticos#

En los enfoques aneriores, determinísticos, los parámetros se suponen constantes obteniendo como resultado único el factor de seguridad como índice de estabilidad. POr le contrario, en los enfoques prababilísticos los parámetros son variables por loq ue el resultado es uan distribución de probabilidad obteniendo como resultado el índice de confiabilidad y la probabilidad de ruptura. Para esto existen métodos como Montecarlo, Estimativas puntuales y el método FOSM. A continuación se describirá el método FOSM.

Modelo FOSM (First Order Second Moment)#

El método estadístico FOSM emplea la expansión de la serie de Taylor de primer orden para derivar el primer y segundo momento de variables de entrada aleatorias. De esta forma, para estimar indirectamente la probabilidad de falla se calcula el Índice de Confiabilidad, dado por la relación entre la media y la desviación estándar de una función de probabilidad que se ajusta al factor de seguridad.

FOSM

Fig. 72 Índice de Confiabilidad#

Las ventajas del modelo FOSM consisten en que los cálculos son simplificados y solo requiere el conocimiento de los valores de los momentos de las distribuciones estadísticas de las variables que forman la función, expresados en la media y la varianza de cada variable, asumiendo una distribución normal tanto para las variables como para el factor de seguridad (FS). De esta manera, para N variables aleatorias no correlacionadas F (\(x_1, x_2,…, x_N\)), conservando solamente los términos del primer orden (lineales) de la serie de Taylor.

Para el uso de N variables independientes con naturaleza probabilística inicialmente se calcula el factor de seguridad media \(E[FS]\) utilizando los valores de los parámetros medios:

\(E[FS]=F(\bar{x}_1, \bar{x}_2, ..., \bar{x}_N)\)

Posteriormente se debe calcular la varianza del factr de seguridad con la siguiente expresión:

\(V[FS]=\sum_{i=1}^{N}(\frac{dFS}{dx_i})^2 V(x_i)\)

Donde \(\bar{x}_𝑖\) y \(V(𝑥_𝑖)\) son la media y varianza de cada variable aleatoria, respectivamente.

Con la varianza, se determina la desviación estandar y finalmente, se obtiene el Índice de Confiabilidad del Factor de Seguridad, calculado por:

\(\beta I=\frac{E[FS]-1}{\sigma[FS]}\)

Donde 𝐸[𝐹𝑆] es el valor esperado del factor de seguridad calculado con los parámetros medios de las variables independientes y 𝜎[𝐹𝑆] es la desviación estándar del Factor de Seguridad (FS) obtenida , teniendo como el FS crítico el valor igual a 1. Este índice expresa la confiabilidad del factor de seguridad en relación con la probabilidad de falla o ruptura. El método FOSM permite evaluar la variabilidad de cualquiera de los parámetros incluidos dentro del análisis.

Ejemplo#

Tomado del Prof. George Fernandes Azevedo de la Universidad Federal do Maranhao (Brasil)

FOSM1

Fig. 73 Ecuación a utilizar y variables aleatorias#

FOSM2

Fig. 74 Factor de seguridad medio#

FOSM3

Fig. 75 .#

FOSM4

Fig. 76 Cálculo de la varianza#

FOSM5

Fig. 77 Cálculo del Indice de Confiabilidad#

Referencias#