IV Jornadas de Ingeniería del Agua

La precipitación y los procesos erosivos

Córdoba, 21 y 22 de Octubre 2015

Modelos Hidráulicos Bidimensionales de

Lluvia Directa aplicados al Cálculo de Eventos
Hidrológicos
Paloma Batanero Akerman
Responsable Soporte y Formación InfoWorks

Eduardo Martínez
Departamento de Ingeniería Civil: Hidráulica y Energética. Universidad P. de Madrid

Ignacio Martínez González

Subdirector de Tecnologías de la Información y Comunicación. INECO

1. Introducción
Una de las principales preocupaciones a la hora de realizar cualquier proyecto de ingeniería
es la adecuada determinación de los hidrogramas de avenida de los cursos de agua que lo
afecten. Hoy en día es posible la simulación de modelos hidrológicos complejos de lluvia
directa sobre malla 2D de elementos finitos, empleando los avances en cálculo vectorial que
ofrece la tecnología CUDA (Arquitectura Unificada de Dispositivos de Cálculo), sin la
necesidad de determinar subcuencas y tránsitos (como en HEC-HMS). Con los datos que
ofrecen los Modelos Digitales del Terreno GRID del PNOA español (Plan Nacional de
Ortofotografía Aérea) con precisión de 5 x 5 m y las Coberturas del Terreno SIG CORINE
(Coordinación de Información del Medio Ambiente), que permite la evaluación de la
rugosidad del suelo, se tiene la información necesaria para construir fácilmente este tipo de
modelos.

2. Metodología seguida
Para demostrar el método, se modelizaron dos cuencas completas mediante lluvia directa
sobre una malla de elementos 2D. Se ensayaron diferentes niveles de precisión (variando el
número de elementos) y se incluyeron dos presas, usando los datos del Modelo Digital del
Terreno del PNOA. También se completó el cálculo hidrológico de ambas mediante HEC-
HMS.
Los dos ejemplos fueron simulados en un PC con procesador Intel Core i7 CPU @2.67 GHz,
equipado con una tarjeta gráfica aceleradora GeForce GTX TITAN Black GPU (2880 núcleos
CUDA) en un SO Windows 7 a 64 bit. Las simulaciones en esta máquina tardaron entre 20 y
2
70 minutos para el río Zapardiel (1,450 km ) y de 6 a 17 horas para la cuenca del Zújar (8,500
2
km - 2 embalses).

B.28.
El programa InfoWorks ICM (Modelo de Cuenca Integrado) fue el software elegido para
realizar los cálculos hidráulicos e hidrológicos en 2D.
ICM es una plataforma de modelización integrada con tecnología de última generación que
puede incorporar cuencas tanto urbanas como fluviales. Este software proporciona una
integración completa de las técnicas de simulación hidrodinámica 1D y 2D, siendo capaz de
simular el desagüe superficial y el enterrado (realizado por el hombre) de las cuencas,
representando con precisión todas las vías de flujo. También permite incorporar la hidrología
y la hidráulica natural y artificial a la vez en un solo modelo. Tiene una total compatibilidad
ArcGIS y utiliza de manera eficiente los avances de hardware disponibles, tales como
procesadores de varios núcleos y tarjetas gráficas de procesamiento dedicadas (GPU).

3. Flujo plano poco profundo

3.1. Ecuaciones
Las Ecuaciones del Movimiento Plano Poco Profundo representan una aproximación del
movimiento de la superficie libre de un fluido sobre un plano horizontal sobre una superficie
inferior dada por la función del lecho zB(x;y). Denominando h al calado del fluido; u y v a las
componentes normalizadas según el calado de la velocidad cartesiana, estas ecuaciones 2D
se pueden escribir de la siguiente manera:

U
  F  Sb  S f [1]
t

Donde los vectores de flujo U, Sb y Sf vienen dados por:

h  0  0 
Z / x ; S   g n u  v
2 2 2
U  hu ; S b   gh  b  f
u 
  [2]
h1 / 3
hv  Z b / y  v 

Y F tiene la expresión:

 hu hv 

F  E , G   h u  g h / 2
2 2
huv  [3]

 huv h u 2  g h2 / 2

La ecuación de la conservación de la energía es:

H h
 v  grad H  g [4]
t t

Donde la energía, H (x,y), se expresa como:

B.28.
 1
H V V  g( Z B  h) [5]
2

La integración de la ecuación (1) en un volumen finito Vk se puede escribir (según F. Alcrudo

y J. Mulet), en la forma:


t V U  dVV  S F  n dS  V (Sb  S f ) dV
k k k
[6]

Los autores usan el método clásico de TVD (Total Variation Dimishing), empleado para
resolver ecuaciones diferenciales parciales a partir de la aproximación Roe de la solución de
Riemann para extrapolar volúmenes finitos de doble paso, con precisión de segundo orden e
integración explicita en el tiempo. Se ha comprobado que el método es esTabla y preciso en
muchas aplicaciones, pudiéndose aplicar en mallas estructuradas y no estructuradas de
elementos tanto triangulares como rectangulares.
La ecuación anterior puede sincretizarse en la forma:

 nwk Δswk  S *bkn  S nfk1/2 

wk nk
Δt Δt
Ukn1/2  Ukn 
2Vk
 Fwk
*n
2
[7]
wk 1

Ukn1  Ukn 
Δt wk nk *n1/2
 F
Vk wk 1 wk

 nwk Δswk  Δt S *bkn1/2  S nfk1  [8]

Donde Uk es el valor medio de las variables conservativas sobre el volumen finito Vk (ya que
en el plano bidimensional, Vk es el área de una superficie), que está limitada por la superficie
Sk (en el plano bidimensional, una línea).
Sk se descompone en muros planos (segmentos de líneas en 2D) con áreas de superficie ∆swk.
El subíndice wk se refiere al muro del plano correspondiente en el que Sk se subdivide. El
vector normal en sentido exterior de wk es nwk.
Los superíndices n, n + 1=2 y n + 1, se refieren a sucesivos pasos temporales, con ∆t,
incremento de tiempo, entre los niveles n y n + 1. F* representa el tensor numérico y Sb* el
vector numérico fuente, que corresponde sólo a la pendiente del lecho.
Las fuerzas de fricción Sf se calculan directamente a partir de la expresión (4), que primero
tiene que ser linealizada para evitar resolver un sistema no lineal, según:
n

U

 S f  n1
S nfk1  S nfk    Uk  Uk
n
 [9]
 k

3.2. Aplicación a Cuencas Naturales

En cuencas naturales las pendientes pueden ser importantes, por lo que las ecuaciones del
Flujo Plano Poco Profundo se ven influenciadas tanto por la componente vertical del calado

B.28.
como por la rugosidad. El calado en el modelo se calcula como diferencia entre la superficie
libre del agua y la cota del terreno, por lo que se mide en la vertical y en consecuencia tanto
la velocidad como el calado y la fricción se ven afectadas por la pendiente, siendo las
relaciones geométricas:

1 L'
Tag a  ; Cos a  ; L  L cos a [10]
p L

Por tanto:

V
d  d cos a; V  a Calado (d)[11]
cos a
Alzado 1 L Velocidad (v)

Calado (d) a
a
p
Alzado 1 L Velocidad (v)
d'

a Alzado proyectado V'

p
L'
Figura 1. Calado real y velocidad en una pendiente y sus proyecciones que serían los resultados
d'
obtenidos por el modelo matemático.
V'
La Alzado proyectado
rugosidad presenta un problema similar, de forma que utilizando la fórmula de Manning,
se necesita emplear un valorL' corregido de n para que compense la diferencia entre la
longitud real de la superficie y la proyectada sobre el plano horizontal, así como la diferencia
de velocidades, de manera que considerando la fórmula de Manning e igualando pérdidas
de energía y asumiendo que Rh  d , se tiene:

n  n cosa 13/ 6 [12]

Con el fin de conocer el error que se comete al no considerar las pendientes en la ladera se
han analizado los errores en el cálculo de calados y velocidades y en la estimación del valor
de la n de Manning para pendientes entre el 1 y el 20%.
Pendiente (%) Error en calados y velocidades (%) Error en n (%)
1.0% 0.005% 0.01%
2.0% 0.020% 0.04%
5.0% 0.125% 0.27%
10.0% 0.496% 1.08%
15.0% 1.106% 2.44%
20.0% 1.942% 4.34%
Tabla 1. Errores cometidos en la estimación del calado, velocidad y valor de n, según la pendiente

B.28.
4. Modelo del río Zapardiel
El primer modelo presentado es el río Zapardiel (cuenca del Duero), con una superficie total
2
de 1.450 km . Es una cuenca sin regulación, pero la sección central presenta algunas áreas
de almacenamiento naturales. Se han comparado dos simulaciones: la primera, un modelo
tradicional HEC-HMS y otro de lluvia directa sobre malla 2D mediante InfoWorks ICM.

4.1. Modelo HEC-HMS

El modelo HEC-HMS realizado de comparación tiene 6 subcuencas y 4 tránsitos.

3520_
2

3520_
1 3510_
2
3505_
2

3510_
3505_ 1
1

Figura 2. Modelo Hidrológico de cuenca río Zapardiel HEC-HMS

Se empleó un NC medio de 75, junto con el hidrograma unitario SCS. Los tránsitos fueron
modelados con el método Muskingum-Cunge de 8 puntos. La precipitación, de 172 mm
durante una hora, corresponde a una precipitación efectiva 100 mm en cada subcuenca.
2
Subcuenca Área (km ) Longitud (km)
3505_1 906.1 60.10
3505_2 78.3 9.10
3510_1 128.8 28.50
3510_2 51.7 9.40
3520_1 208.2 16.90
3520_2 81.2 17.12
Tabla 2. Características de los elementos del modelo del río Zapardiel

B.28.
4.2. Modelo de Lluvia Directa de InfoWorks ICM
Simultáneamente se construyó un modelo de precipitación directa sobre malla 2D a partir
de un MDT GRID PNOA 5x5 m. La malla 2D general se probó con diferentes tamaños de
2
elementos (5.000, 10.000 y 20.000 m ). Junto a esta, se definió una malla más precisa de
2
1.000 m en los cursos de agua.

Figura 3. Modelo del terreno GRID empleado en el Zapardiel y malla 2D definida en InfoWorks ICM
Se vertieron sobre la cuenca 100 mm de precipitación total durante una hora, igual que en el
modelo HEC-HMS. El coeficiente de rugosidad de Manning empleado fue de 0.040.

4.3. Comparación de Resultados

Los resultados son muy parecidos en los tres modelos realizados, que varían entre 618,4
3 2 3
m /s obtenidos con un tamaño de elemento de 10.000 m , 619,0 m /s (elementos de 5,000
2 3
m ) y 627,5 m /s (con elementos de 20.000 m2). La rama descendente del hidrograma tiene
una pendiente suave, debido al efecto de laminación esperado en la cuenca y el cauce.
Al comparar los resultados de este modelo con el de lluvia directa, se aprecia una gran
similitud en la forma del hidrograma, pero con la diferencia fundamental del caudal máximo
3 3
(787 m /s frente a los 619 m /s del modelo de precipitación directa con ICM). La razón de
esta diferencia se debe a que el modelo HEC-HMS no reproduce el efecto de laminación en

B.28.
las zonas más llanas de las subcuencas, y por tanto la parte baja del hidrograma tiene una
pendiente más pronunciada.

Caudal
Min Max
...Y 5.000 m2 malla central>Lluvia 100 mm 1 h 0.000 618.924
... 10.000 m2 malla central>Lluvia 100 mm 1 h 0.000 618.381
... 20.000 m2 malla central>Lluvia 100 mm 1 h 0.000 627.475
Figura 4. Resultados de caudales obtenidos en el modelo 2D del río Zapardiel al final de la cuenca

El hidrograma de salida obtenido con el modelo HEC-HMS fue:

Figura 5. Resultados de caudales obtenidos al final del modelo HEC-HMS del río Zapardiel

5. RÍO ZÚJAR
El segundo realizado corresponde al río Zújar, de la Cuenca del Río Guadiana, con una
2 3
superficie total de 8.500 km , incluyendo dos grandes presas: La Serena, de 3.200 hm y el
3
Embalse del Zújar, situado aguas abajo de la primera, con un volumen útil de 300 hm . El
tamaño de este río puede clasificarse como medio-alto, y supone un buen ejemplo de la
metodología utilizada cuando hay regulación de embalses en la cuenca.
A partir de los datos de elevación PNOA, se generó un MDT de toda la cuenca, con tamaño
de celda de 5 x 5 m.
Los modelos hidráulicos bidimensionales definen el terreno a partir de mallas de superficie
formadas por elementos de diferentes tamaños y formas, que disponen sus vértices en el
terreno, de los cuales obtienen la cota de cada del elemento. Si se ignoran otros
condicionantes, el tamaño de los elementos especifica la precisión del modelo, así como el
número total de éstos y a su vez, el tiempo total de cálculo. Por tanto se debe buscar un

B.28.
compromiso entre la precisión requerida y el tiempo de cálculo; definiendo diferentes zonas
de mallado con diferentes tamaños de elementos.
En el caso de la hidrología de precipitación directa, es necesario proporcionar una precisión
suficiente en los cauces, especialmente en aquellos tramos con secciones más estrechas.
Esto se debe hacer con el fin de evitar obstrucciones de flujo que se comportarán como
zonas de almacenamiento, reduciendo al mínimo el flujo total aguas abajo.

La Serena

Embalse del Zújar

Figura 6. Modelo Digital del Terreno de la cuenca del río Zújar y resultados iniciales de inundación

Se ha desarrollado un método sencillo para identificar los principales cauces de una cuenca,
donde se debe disponer una mayor precisión en los elementos del modelo 2D. Para ello se
parte de un modelo inicial con elementos del mismo tamaño en toda su extensión. A partir
de una primera simulación, se localizan los cauces de la cuenca en aquellas áreas con
resultados de inundación mayores de un determinado calado. A partir de estos resultados,
se lleva a cabo un procesamiento de datos con la ayuda de un SIG para delimitar las zonas
que definen los cursos de agua y que se mallan con mayor precisión en el siguiente paso.
Repitiendo este proceso unas pocas veces, se alcanza un punto óptimo, en el que las zonas
de flujo principal se hayan delimitadas adecuadamente. Se puede dar por terminado el
proceso cuando haya unos 4 ó 5 elementos de malla en el interior del lecho del río.

B.28.
Los resultados presentan un incremento gradual en el flujo al disminuir el tamaño del
Simulation Plot Realizado por PC (31/10/2014 14:08:03) Page 2 of 2
Simulación: >Simulaciones>Simulacion_Cauces_v05>Lluvia 33.7 mm (80mm) 6h (30/10/2014 19:21:26)
elemento, hasta un valor en el que se estabilizan los resultados. Los gráficos siguientes
Simulación: >Simulaciones>Simulacion_Cauces_v02>Lluvia 33.7 mm
Simulación: >Simulaciones>Simulacion_Cauces_v03>Lluvia 33.7 mm
(80mm)
(80mm)
6h
6h
(30/10/2014 18:34:37)
(30/10/2014 18:48:19)
muestran resultados de caudales finales en afluentes.
Simulación: >Simulaciones>Simulacion_Cauces_v04>Lluvia 33.7 mm
Lista de Selecciones: Custom Selection
(80mm) 6h (30/10/2014 18:58:36)

Caudal
Min Max
Tamaño general del elemento 5.000 m2 -0.257 1382.335
Tamaño general del elemento 10.000 m2 -0.357 1378.366
Tamaño general del elemento 20.000 m2 -0.353 1290.516
Single Simulation
TamañoPlot Realizado
general del por PC (31/10/2014
elemento 50.000 m214:20:06) Page 1 of 1 -0.362 1201.711
Simulación: >Simulaciones>Simulacion_Cauces_v05>Lluvia 33.7 mm (80mm) 6 h (30/10/2014 19:21:26)

Figura 12: Caudales en el río Guadalefra (afluente del Zújar); tamaños de malla 5000-50000 m
Lista de Selecciones: Custom Selection

Caudal (m3/s)
Min Max Volume (m3)
Línea Resultados Red (2D) 0.000 3957.350 342590682.205

3
Figura 13: Caudales en el río y sección final de la cuenca completa del río Zújar ≈ 4,000 m /s

6. CONCLUSIONES
En este artículo se presenta un nuevo método para la determinación de caudales máximos
de avenida, dando una alternativa válida al análisis tradicional de cuencas y tránsitos
unidimensionales de modelos hidrológicos como HEC-HMS, donde los resultados dependen
de la habilidad del ingeniero para definir los parámetros de la cuenca.
Gracias a la disponibilidad de datos públicos de GRID PNOA y cobertura del terreno CORINE,
actualmente es posible construir modelos de mallas de elementos 2D para simular lluvia
directa sobre cuenca sin necesidad de costosos levantamientos topográficos. El parámetro
más importante en este tipo de modelos es el tamaño medio de los elementos 2D en los
cursos de agua principales. El método propuesto permite una mejor definición de estas
áreas después de unas pocas iteraciones, acelerando por tanto el proceso de construcción
del modelo sin comprometer la calidad de los resultados.
Con el hardware disponible hoy en día en ordenadores de uso general, se puede establecer
un límite superior razonable para mallas 2D de 4 millones de elementos, lo que permite la
2
correcta simulación de la lluvia directa en cuencas de hasta 10.000 km . Para cuencas más

B.28.
grandes, este proceso se podría aplicar de forma repetida a varias zonas más pequeñas o en
sucesivas secciones de la cuenca.
Durante los últimos 5 años, la aparición de los procesadores multinúcleo y las tarjetas
gráficas basadas en el procesamiento paralelo CUDA ha resultado en un aumento de cinco
veces en la velocidad de cálculo matemático. Proyectando esta tendencia, se podría suponer
que en los próximos cinco años será posible la simulación de lluvia directa en 2D en mallas
de más de 20.000.000 de elementos en plazos de tiempo razonables y en PCs de uso
estándar.
Por ello, es razonable pensar que el método propuesto puede ser una herramienta de
referencia para obtener hidrogramas de avenida que podrían ser utilizados para prevenir
eventos de sobrevertido e inundaciones no deseados.

Referencias
Alcrudo, F. and Mulet-Marti, J., 2005. Urban inundation models based on the Shallow Water
equations. Numerical and practical issues. Proceedings of Finite Volumes for Complex
Applications IV. Problems and Perspectives. Pages 1-12. Editado por F. Benkhaldoun, D. Ouazar,
S. Raghay. Hermes Science Pub.
Cleveland, T. G., Garcia, A., He, X., Fang, X., and Thompson, D. B. 2005. Comparison of physical
characteristics for selected small watersheds in Texas as determined by automated and manual
methods. Proc., Texas ASCE Section Fall Meeting, El Paso, Tex.
Fang, Xing, Thompson, David B., Cleveland, Theodore G., Pradhan, P. and Malla Ranjit, 2008.
Time of Concentration Estimated Using Watershed Parameters Determined by Automated and
Manual Methods. Journal of Irrigation and Drainage Engineering, ASCE, March/April, 2008.
Hill, J., and Neary, V. S., Factors affecting estimates of average watershed slope. J. Hydrol. Eng.,
10(2), 130–140. 2005.
Roussel, M. C., Thompson, D. B., Fang, X., Cleveland, T. G., and Garcia, A. C., 2005. Timing
parameter estimation for applicable Texas watersheds. Research Rep. No. 0-4696-2, Texas
Department of Transportation, Austin, Tex.
Su, D. H., and Fang, X., 2004. Estimating traveling time of flat terrain by 2D overland flow
model. Shallow flows, G. H. Jirka and W. S. J. Uijttewaal, eds., Balkema, Rotterdam, Países
Bajos, 629–635.
U. S. Army Corps of Engineers. Hydrologic Engineering Center. HEC-HMS, Applications Guide.
December 2000.
U. S. Army Corps of Engineers. Hydrologic Engineering Center. HEC-HMS; Hydrologic modelling
system user’s manual. 2000.

B.28.