MODELIZACIÓN DE LA TURBULENCIA EN AGUAS POCO PROFUNDAS Luis Cea Gómez Grupo de Ingeniería del Agua y del Medioambiente (GIAMA) E.T.S. Ingenieros de Caminos, Canales y Puertos. Universidad de A Coruña. Campus de Elviña s/n. 15071 A Coruña
[email protected]
1. Introducción Un gran número de estudios en ingeniería hidráulica implica el análisis de flujos en lámina libre, muchos de los cuales pueden considerarse flujos poco profundos, refiriéndonos con el término poco profundo a una relación entre dimensiones vertical y horizontal pequeña. Prácticamente la totalidad de flujos en lámina libre son turbulentos. En cualquier río pueden observarse pequeños remolinos que aparecen y desaparecen con un movimiento aparentemente caótico, mostrando la complejidad del movimiento turbulento. Estos remolinos turbulentos son los principales responsables de los procesos de mezcla, por lo que juegan un importe papel en la difusión de sustancias solubles, de sólidos en suspensión, etc. La interacción entre corrientes naturales y estructuras puede crear estructuras turbulentas de dimensiones considerables, como muestra la Figura 1.
Figura 1. Estructuras turbulentas en corrientes naturales.
La capacidad de calcular de forma precisa los campos de velocidad, calado y turbulencia es importante para diseñar el emplazamiento de emisarios, evaluar el transporte de sedimentos, o diseñar la geometría de canales y depósitos, simplemente por citar algunos ejemplos. Los modelos numéricos para el cálculo de flujo en lámina
libre son una herramienta cada vez más utilizada en ingeniería hidráulica. Las principales ventajas de los modelos numéricos respecto a los modelos físicos de laboratorio y a las medidas de campo son: (1) proporcionan una cantidad de información muy completa y extensa; (2) se obtienen resultados de forma más rápida y barata con un modelo numérico que con un modelo experimental; (3) se puede modificar fácilmente la geometría y condiciones de contorno para la evaluación de escenarios futuros. En su contra, los resultados numéricos no son tan precisos como pueden ser las medidas experimentales, y en general necesitan de una validación y calibración previa. La elección de un modelo numérico concreto depende del problema considerado. En la actualidad existen modelos 3D para el cálculo del flujo en lámina libre en estructuras hidráulicas con geometría compleja. El coste computacional de estos modelos es muy elevado, especialmente en ingeniería fluvial, en donde el dominio de estudio es muy extenso y la geometría totalmente irregular, por lo que en la actualidad son mucho más utilizados los modelos 2D e incluso 1D. La potencia de los ordenadores actuales permite asumir la utilización de modelos 2D (2D-SWE) en proyectos de ingeniería hidráulica a un coste computacional razonable. Los modelos 2D-SWE han sido utilizados con éxito para simular flujos poco profundos con un nivel de turbulencia elevado como pueden ser zonas de recirculación, canales de toma o escalas de peces, si bien su uso más generalizado se centra en la simulación del flujo en ríos y en regiones costeras. A pesar de que prácticamente todos los flujos en ingeniería hidráulica son turbulentos, en determinados casos la turbulencia no es lo suficientemente alta como para tener una influencia notoria en el campo de velocidad media. Este suele ser el caso de flujo en ríos, estuarios y en general en zonas costeras con una geometría lo suficientemente suave como para que no se produzcan zonas de recirculación en planta. Sin embargo, incluso en este tipo de situaciones es importante realizar una correcta modelización de la turbulencia, ya que esta juega un papel fundamental en los procesos de transporte y mezcla de contaminantes y sedimentos. La difusión de calor, de un soluto, o de un sedimento en suspensión se produce básicamente por turbulencia, excepto en flujo laminar, el cual no suele darse en general en ingeniería hidráulica, y mucho menos en ríos o estuarios. El coeficiente de difusión turbulenta es varios órdenes de magnitud superior al coeficiente de difusión molecular. Por lo tanto es necesario evaluar previamente la energía cinética turbulenta para poder calcular el flujo difusivo. La modelización de la turbulencia en flujos poco profundos no ha sido tratada tan ampliamente como en otras disciplinas en dinámica de fluidos. Se has propuesto algunos modelos promediados en profundidad para las ecuaciones de aguas someras bidimensionales. Son básicamente modelos que derivan directamente de modelos 3D estándar ya existentes, con ligeras modificaciones que tienen en cuenta la turbulencia producida por la tensión tangencial de fondo. Debe mencionarse el modelo k-ε promediado en profundidad de Rastogi y Rodi (1978), el cual fue uno de los primeros propuesto para flujos poco profundos y todavía uno de los más utilizados en la actualidad. En este trabajo se presentan y describen diferentes modelos de turbulencia específicos para las ecuaciones de aguas someras bidimensionales. Para mostrar las capacidades y limitaciones de los principales modelos de turbulencia para aguas someras se presentan
diversas aplicaciones prácticas con validaciones experimentales. Se ha optado por incluir tanto ejemplos en las que la turbulencia tiene una influencia notable en la velocidad media, como ejemplos en los que, aún siendo el flujo turbulento, los efectos de las tensiones turbulentas en el flujo medio son inapreciables. De esta forma se pretende dar una idea general de cuándo es necearlo utilizar un modelo de turbulencia y cuándo puede obviarse. A pesar de ello, es necesario remarcar que es difícil establecer con carácter general cuando sí es necesario el modelo de turbulencia y cuándo no. En la sección 2 se presentan brevemente las ecuaciones de aguas someras bidimensionales, haciendo hincapié en los términos donde aparecen las tensiones turbulentas. Los diferentes modelos de turbulencia para aguas someras se describen en la sección 3, incluyendo el modelo parabólico, el de longitud de mezcla, 3 versiones del modelo k-ε, y un modelo de tensiones algebraicas. La sección 4 muestra algunos resultados numéricos y validaciones experimentales obtenidos con un modelo de aguas someras bidimensional. Finalmente, en la sección 6 se presentan unas breves conclusiones a modo de resumen general.
2. Las ecuaciones de aguas someras bidimensionales Debido a su carácter marcadamente bidimensional, los flujos poco profundos se pueden representar adecuadamente mediante un campo bidimensional de calado y velocidad horizontal promediada en profundidad. Las ecuaciones que rigen la evolución de dichas variables son las ecuaciones de aguas someras promediadas en profundidad, también conocidas como 2D Shallow Water Equations (2D-SWE) o ecuaciones de St.Venant bidimensionales. Las ecuaciones asumen una escala espacial vertical mucho más pequeña que la escala horizontal (Figura 2), lo cual permite asumir una distribución de presión hidrostática. Al mismo tiempo asumen un perfil vertical de velocidades homogéneo en profundidad. La hipótesis de presión hidrostática se cumple razonablemente en el flujo en ríos, así como en las corrientes generadas por la marea en estuarios y zonas costeras. Asimismo, la hipótesis de distribución uniforme de velocidad en profundidad se cumple habitualmente en ríos y estuarios, aunque pueden existir zonas en las que dicha hipótesis no se cumpla debido a flujos locales tridimensionales o a cuñas salinas. En estos casos es necesario estudiar la extensión de dichas zonas y su posible repercusión en los resultados. En la actualidad, los modelos numéricos basados en las 2D-SWE son los más utilizados en estudios de dinámica fluvial y litoral, evaluación de zonas inundables, transporte de sedimentos y contaminantes, etc.
Figura 2. Separación de escalas horizontal (Ls) y vertical (Hs) en flujos poco profundos.
Las ecuaciones 2D-SWE son un conjunto de 3 ecuaciones diferenciales de transporte (conservación de masa y conservación de momento en las 2 direcciones horizontales) con 3 incógnitas (calado y 2 componentes de la velocidad promediada en profundidad), que se pueden expresar como:
∂h ∂hU x ∂hU y + + =0 ∂t ∂x ∂y e ∂z s τ s, x τ b, x ∂hU x ∂hU 2x ∂hU x U y ∂hτ exx ∂hτ xy + + = −gh + − + 2 Ω h sinλ U y + + ∂t ∂x ∂y ∂x ρ ρ ∂x ∂y
∂hU y ∂t Ux =
+
∂hU x U y ∂x
+
∂hU 2y
1 Zs u x dz h ∫Z b
∂y
Uy =
= −gh
∂ hτ exy ∂ hτ eyy ∂z s τ s, y τ b, y + − − 2 Ω h sinλ U x + + ∂y ρ ρ ∂x ∂y
1 Zs u y dz h ∫Z b
h = zs − zb
en donde h es el calado, u x , u y son las velocidades horizontales, Ux, Uy son las velocidades horizontales promediadas en profundidad, g es la aceleración de la gravedad, zb es la elevación del fondo, zs es la elevación de la superficie libre, τs es la fricción en la superficie libre debida al rozamiento producido por el viento, τb es la fricción debido al rozamiento del fondo, ρ es la densidad del agua, Ω es la velocidad angular de rotación de la tierra, λ es la latitud del punto considerado, τexx, τexy, τeyy son las tensiones tangenciales efectivas horizontales. Las ecuaciones pueden además tener en cuenta el efecto de variaciones en la presión atmosférica, así como aportaciones puntuales o distribuidas de caudal. Las tensiones efectivas horizontales que aparecen en las ecuaciones hidrodinámicas incluyen los efectos de las tensiones viscosas, de las tensiones turbulentas y los términos de dispersión debido a la no homogeneidad en profundidad del perfil de velocidad.
τ ije = τ ijv − u'i u' j + D ij en donde τvij son las tensiones viscosas, u’ es la oscilación de la velocidad instantánea respecto a la velocidad media, u'i u' j son las tensiones turbulentas (también llamadas tensiones de Reynolds), y Dij son los términos de dispersión: D ij =
1 Zs (U i − u i ) (U j − u j )dz h ∫Z b
Los términos de dispersión se desprecian en las ecuaciones 2D-SWE (hipótesis de perfil de velocidad uniforme en profundidad), debido a la imposibilidad de calcularlos de forma general. Su importancia será mayor cuanto menos uniforme sea el perfil de velocidad en profundidad. Una situación típica en la que estos términos pueden cobrar importancia es en canales con codos o radios de curvatura pequeños, así como en la confluencia de canales (Figura 3).
Q1
Q
Q
3
2
Figura 3. Flujos secundarios (izquierda) y perfil vertical de velocidad (derecha). Principales causas de los términos de dispersión.
Las tensiones viscosas se calculan a partir de la viscosidad cinemática del fluido (ν) como: ∂U ∂U j τ ijv = ν i + ∂x ∂ x i j
En general, excepto cerca de las paredes, y excepto en flujo laminar, el orden de magnitud de las tensiones viscosas es mucho menor que el del resto de los términos que aparecen en las ecuaciones hidrodinámicas, y pueden por lo tanto despreciarse. Las tensiones turbulentas son mucho más importantes, especialmente en zonas de recirculación, en donde la producción de turbulencia es elevada (Figura 4). En el caso de las ecuaciones de aguas someras bidimensionales las tensiones turbulentas constituyen 3 nuevas incógnitas a calcular, que sumadas al calado y a las velocidades Ux, Uy produce un total de 6 incógnitas. Esto es lo que se conoce como problema de cierre de la turbulencia, porque es necesario resolver un conjunto de 3 ecuaciones con 6 incógnitas. Debido a ello, es necesario utilizar un modelo de turbulencia que permita calcular las tensiones turbulentas. Lo deseable sería tener un modelo que, al igual que la Ley de Fick permite calcular la difusión molecular de una sustancia soluble, permita calcular la difusión turbulenta de forma general. El problema es que no existe un modelo de turbulencia universal, que permita calcular de forma precisa las tensiones turbulentas, por lo que a lo largo del tiempo se han ido desarrollando diferentes modelos de mayor o menor complejidad. A continuación se describen los modelos de turbulencia más utilizados en aguas someras.
Figura 4. Distintos flujos con zonas de recirculación en planta importantes.
3. Modelos de turbulencia promediados en profundidad 3.1.
Escalas de turbulencia en aguas someras
Una de las principales características de los flujos poco profundos es la separación entre escalas horizontales y verticales, debido a que la extensión vertical del fluido (limitada por la profundidad) es mucho menor que su extensión horizontal. Esta separación de escalas es aplicable tanto a la dimensión espacial como a las velocidades, y por lo tanto a la turbulencia. En el caso de la turbulencia, su principal efecto supone una separación entre estructuras turbulentas (remolinos) tridimensionales y estructuras turbulentas bidimensionales. La escala espacial de la turbulencia 3D está limitada por la profundidad, y por lo tanto son estructuras mucho más pequeñas que las asociadas a la turbulencia 2D, las cuales están únicamente limitadas por la escala horizontal. La turbulencia 3D está generada principalmente por el rozamiento del fondo, mientras que la turbulencia 2D está generada por gradientes de velocidad en el plano horizontal. Es importante que el modelo de turbulencia incluya los efectos tanto de la turbulencia 3D, producida por fricción de fondo, como de la turbulencia 2D, producida por gradientes de velocidad horizontales. En los modelos de aguas someras, el carácter bidimensional del flujo está considerado de forma implícita en las ecuaciones de transporte al considerar un perfil de velocidad homogéneo en profundidad, mientras que la producción tridimensional se incluye habitualmente por medio de un término fuente que depende de la tensión tangencial de fondo.
A continuación se presentan varios modelos de turbulencia promediados en profundidad, incluyendo el modelo parabólico, un modelo de longitud de mezcla, tres versiones del modelo k-ε , y un modelo de tensiones algebraicas. 3.2.
La aproximación de Boussinesq en aguas someras
La aproximación de Boussinesq se utiliza en los modelos lineales de viscosidad turbulenta, y relaciona las tensiones de Reynolds con los gradientes de velocidad media mediante la siguiente relación lineal: ∂U ∂U j 2 − k δ ij − uiu j = νt i + ∂x 3 x ∂ i j
La hipótesis de Boussinesq no deja de ser una mera aproximación que funciona mejor o peor en función de las condiciones de flujo. En los modelos de aguas poco profundas, lo que realmente aparece en las ecuaciones es el valor de las tensiones de Reynolds promediado en profundidad. Se puede asumir directamente el uso de una formulación similar para el cálculo de tensiones de Reynolds promediadas en profundidad, utilizando para ello las siguientes expresiones: − u' 2 = 2ν t
∂U y ∂U − u' v' = ν t x + ∂x ∂y
∂U x 2 − kδ ij ∂x 3
− v' 2 = 2ν t
∂U y ∂y
−
2 kδ ij 3
en donde todas las variables hacen referencia a valores promediados en profundidad. Introduciendo la aproximación de Boussinesq en las ecuaciones de aguas someras y simplificando algunos términos de segundo orden, se obtiene la siguiente formulación: ∂h ∂hU x ∂hU y + + =0 ∂t ∂x ∂y τ b, x ∂z ∂hU x ∂hU 2x ∂hU x U y ∂U x ∂ + + = −gh s − + h νe ∂t ∂x ∂y ∂x ρ ∂x ∂x ∂hU y ∂t
+
∂hU x U y ∂x
+
∂hU 2y ∂y
= −gh
∂U x ∂ + h ν e ∂y ∂y
∂U y ∂ ∂U y ∂z s τ b, y ∂ + h νe − + h ν e ∂y ρ ∂x ∂x ∂y ∂y
νe = ν + νt
en donde por simplicidad en la notación se han omitido la fuerza de Coriolis y el rozamiento por viento. La viscosidad efectiva νe incluye la viscosidad laminar ν y la viscosidad turbulenta νt, aunque como se ha dicho anteriormente, en flujos turbulentos la viscosidad turbulenta es varios ordenes de magnitud superior a la viscosidad laminar (en hidráulica fluvial suele ser entre 3 y 4 órdenes de magnitud superior). 3.3.
Viscosidad turbulenta constante
El orden de magnitud de la viscosidad turbulenta se puede fijar de forma aproximada en función del flujo considerado. Existen diferentes publicaciones en las que se proponen valores aproximados de la viscosidad turbulenta en función del flujo considerado. Este enfoque es muy sencillo, y no se puede considerar como un modelo de turbulencia adecuado ni realista en ningún caso, ya que no tiene en cuenta que la viscosidad
turbulenta varía fuertemente de un punto a otro. Es importante remarcar que no es sólo el valor de la viscosidad turbulenta, sino también su variación espacial la que determina el campo de velocidad media. Además las tablas existentes proporcionan únicamente valores aproximados. Por todo ello no se recomienda utilizar este método, ya que puede llevar a resultados con errores considerables, generalmente por utilizar valores excesivamente elevados de viscosidad turbulenta, así como por no considerar su variabilidad espacial. Otro inconveniente importante de este enfoque se produce en la modelización de flujos en régimen no estacionarios, ya que en estos casos la turbulencia varía no sólo en espacio sino también en el tiempo. 3.4.
Perfil parabólico de viscosidad turbulenta
Este modelo asume una distribución parabólica en profundidad de la viscosidad turbulenta, calculándose a partir de dicha distribución una viscosidad promediada en profundidad, la cual viene dada por la siguiente expresión: ν t = 0.068 u f h
en donde h es el calado y uf es la velocidad de fricción del fondo, calculada a partir de la tensión tangencial del fondo como: uf =
τf ρ
Si se utiliza la fórmula de Manning para calcular el rozamiento del fondo se obtiene la siguiente expresión para la viscosidad turbulenta: ν t = 0.068 g n U h 5/6
Es decir, que la viscosidad turbulenta depende localmente del calado, módulo de la velocidad promediado en profundidad y del coeficiente de Manning. Debido a la sencillez de este modelo, a veces se utiliza un coeficiente multiplicador para permitir ajustar mejor el valor de la viscosidad turbulenta. Este coeficiente se fija de forma arbitraria por el usuario. ν t = C m 0.068 g n U h 5/6
3.5.
Modelo de longitud de mezcla
En el modelo de longitud de mezcla para aguas someras, la viscosidad turbulenta se calcula a partir de las características locales del flujo mediante la siguiente expresión: ν t = [min (0.267 κ h, κ d wall )]
2
u 2S ij S ij + 2.34 f κh
2
en donde κ=0.41 es la constante de von Karman. Es un modelo algebraico relativamente sencillo, que permite obtener resultados aceptables en flujos en los que la turbulencia está generada localmente y principalmente por el rozamiento del fondo. Tiene en cuenta la producción de turbulencia debido a gradientes horizontales de velocidad, pero no considera el transporte convectivo ni la disipación de turbulencia. En flujos con zonas de recirculación fuertes los resultados obtenidos con el modelo de longitud de mezcla empeoran. 3.6.
Modelo k-ε de Rastogi y Rodi (1978)
Es un modelo que resuelve una ecuación de transporte para la energía cinética turbulenta k y para la tasa de disipación de energía turbulenta ε. El modelo tiene en cuenta la producción debido al rozamiento del fondo, la producción por gradientes de velocidad, la disipación y el transporte convectivo. Las ecuaciones del modelo k-ε para aguas someras son las siguientes: ∂k ∂U x k ∂U y k ∂ + + = ∂t ∂x ∂y ∂x j ∂ε ∂U x ε ∂U y ε ∂ + + = ∂t ∂x ∂y ∂x j
ν t = cµ
k2 ε
3 ν + ν t ∂k + 2ν t S ij S ij + c k u f − ε σ k ∂x j h
ν ν + t σε
c k = c f−1/2
2 ∂ε u4 + c ε1 ε 2ν t S ij S ij + c ε f − c ε2 ε k k h2 ∂x j
1/2 c ε = 3.6c 3/2 k c ε2 c µ
cf =
τb 1 ρ U2
con las constantes: c µ = 0.09
c ε1 = 1.44
c ε2 = 1.92
σ k = 1.0
σ ε = 1.31
en donde k es la energía cinética turbulenta, ε es la tasa de disipación de turbulencia Sij es el tensor de deformación. Los términos que incluyen la velocidad de fricción de fondo u* son los responsables de modelar la generación de turbulencia por rozamiento de fondo. El modelo k-ε es un modelo relativamente sofisticado. En flujos turbulentos pocoprofundos proporciona resultados relativamente buenos, siendo uno de los modelos más utilizados en dicho ámbito cuando el nivel de turbulencia es importante. No obstante, su grado de complejidad no garantiza resultados correctos en cualquier tipo de flujo. Al igual que cualquier modelo de turbulencia, los resultados obtenidos con el modelo k-ε deben de analizarse y valorarse de forma crítica, para lo cual es fundamental la experiencia del usuario en la modelización de flujos turbulentos.
3.7.
Las constantes modificadas de Booij
Los términos de producción de turbulencia por fricción de fondo en el modelo k-ε de Rastogi y Rodi dependen de 2 coeficientes, ck y cε. El valor de dichos coeficientes determina el nivel de turbulencia en condiciones de flujo uniforme 1D. Booij (1989) propuso la siguiente modificación en el valor de dichas constantes
c Bk =
1 ck 10
c εB =
1 cε 44
El efecto de disminuir los coeficientes ck y cε es que la producción de turbulencia por fricción de fondo se reduce de manera significativa. Si se utilizan las constantes de Booij para evaluar la energía cinética turbulenta y la disipación en flujo uniforme en canal, se obtiene:
k Bu = 0.44k u
ε Bu = 0.10ε u
ν t,Bu = 1.94ν t, u
En donde los valores ku, εu y νt,u son los proporcionados por el modelo de Rastogi y Rodi bajo las mismas condiciones de flujo. Por lo tanto, en condiciones de flujo uniforme en canal, la modificación de constantes propuesta por Booij reduce la energía cinética turbulenta y la disipación de turbulencia, pero incrementa la viscosidad turbulenta en un factor de 1.94. Por otro lado, en flujos dominados por tensiones horizontales, las diferencias con el modelo de Rastogi y Rodi disminuyen, porque en dicho caso la producción por rozamiento de fondo es mucho menor que la producción por gradientes de velocidad horizontal.
3.8.
El modelo k-ε de Babarutsi y Chu
Babarutsi y Chu (1991) propusieron un modelo k-ε promediado en profundidad con dos tipos de turbulencia: tridimensional y bidimensional. La principal diferencia con el modelo de Rastogi es el modo en que se introducen los efectos de la turbulencia tridimensional producida por fricción de fondo. Babarutsi y Chu proponen separar la viscosidad turbulenta total en dos partes, una parte tridimensional y una parte bidimensional.
ν t = ν 2D + ν 3D t t La viscosidad 3D tiene en cuenta los efectos de la turbulencia de pequeña escala generada por fricción de fondo, y se calcula con una ecuación similar a la utilizada en el modelo parabólico. La ecuación utilizada por Babarutsi y Chu es la siguiente:
ν 3D = 0.08 u f h t que por otro lado es el mismo valor dado por el modelo de Rastogi para flujo uniforme en canal. La parte bidimensional se calcula a partir de las variables k’ y ε’, las cuales representan la turbulencia bidimensional de gran escala generada por gradientes horizontales, a partir de la siguiente expresión: ν
2D t
k' 2 = cµ ε'
En este caso la energía cinética turbulenta k’ no incluye la turbulencia tridimensional generada por fricción de fondo. Las ecuaciones que se utilizan para calcular k’ y ε’ son similares a las utilizadas en el modelo de Rastogi, pero sin producción por fricción de fondo, ya que no se considera la turbulencia 3D. Se incluye además un nuevo término
F’ que modela la transferencia de energía entre la turbulencia bidimensional y la turbulencia tridimensional. Las ecuaciones para k’ y ε’ son las siguientes: ∂k' ∂U x k' ∂U y k' ∂ + + = ∂t ∂x ∂y ∂x j
ν ν + t σk
ν ∂ε' ∂U x ε' ∂U y ε' ∂ ν + t + + = ∂t ∂x ∂y ∂x j σε
∂k' ∂x j
+ 2ν 2D t S ij S ij − F'−ε'
∂ε' ε' 2 + c ε1 ε' 2ν 2D S ij S ij − (1 − c ε3 )F' − c ε2 t k' k' ∂x j
(
)
con la constante cε3 = 0.8. Nótese que el término de difusión se modela con la viscosidad turbulenta total, mientras que el término de producción horizontal se modela con la viscosidad turbulenta bidimensional. Todas las constantes del modelo son las mismas que en el modelo de Rastogi. El término fuente F’ está asociado con el trabajo ejercido por las grandes escalas turbulentas horizontales contra la fuerza de fricción, y se calcula como:
F' =
[ (
)
(
c f u' 2 2U 2x + U 2y + 2u' v'U x U y + v' 2 U 2x + 2U 2y
)]
h U 2x + U 2y
en donde las tensiones de Reynolds se calculan mediante la aproximación de Boussinesq utilizando para ello la viscosidad turbulenta bidimensional νt2D. Para flujo uniforme en canal, la producción de turbulencia bidimensional por gradientes de velocidad horizontal es nula, y la ecuación de energía cinética turbulenta 2D se reduce a F’=-ε’, y por consiguiente todas las variables turbulentas de gran escala se anulan:
k 'u = 0
ε 'u = 0
ν 2D t, u = 0
ν t, u = ν 3D t, u = 0.08 u f h
En este caso toda la turbulencia es tridimensional y generada por fricción de fondo. La viscosidad turbulenta total es similar a la proporcionada por el modelo de Rastogi en las mismas condiciones de flujo. Por otro lado, en un flujo dominado por tensiones tangenciales horizontales la velocidad de fricción de fondo tiende a cero, al igual que el término F’ y la viscosidad turbulenta tridimensional. En este caso la turbulencia es predominantemente bidimensional y generada por gradientes de velocidad horizontales. Las ecuaciones para k’ y ε’ se reducen al modelo k-ε bidimensional estándar. 3.9.
Diferencias entre las 3 versiones del modelo k-ε
Las principales diferencias entre las 3 versiones presentadas del modelo k-ε aparecen en flujos dominados por fricción. La modificación de Booij reduce la turbulencia generada por tensión de fondo, y como consecuencia la turbulencia bidimensional de gran escala se incrementa, incrementando al mismo tiempo la difusión turbulenta del modelo. De las tres versiones, la de Babarutsi y Chu es conceptualmente la más interesante, ya que distingue entre 2 tipos de turbulencia diferentes. Las principales diferencias con el modelo de Rastogi aparacen en flujos con recirculaciones dominados por fricción de fondo. Ambos modelos convergen al modelo k-ε bidimensional estándar cuando la
fricción del fondo es despreciable. En el caso contrario, en flujo uniforme en canal donde toda la turbulencia está generada por la tensión de fondo, ambos modelos proporcionan la misma viscosidad turbulenta. Las diferencias entre ambos modelos aparecen cuando existen zonas de recirculación horizontales dominadas por la tensión de fondo. Teniendo en cuenta que en muchos casos los modelos proporcionan resultados similares, y que no existen apenas comparaciones rigurosas ni validaciones experimentales de los resultados proporcionados por los 3 modelos, es difícil afirmar con carácter general cuál de ellos es el más adecuado. Probablemente por su antigüedad y difusión el modelo más utilizado en la práctica es el de Rastogi y Rodi. 3.10. Modelo de tensiones algebraicas promediado en profundidad Uno de los principales problemas de los modelos de viscosidad turbulenta es la aproximación de Boussinesq, la cual asume una viscosidad turbulenta isotrópica para calcular las tensiones de Reynolds. Esta aproximación produce una reducción en la precisión del modelo cuando la turbulencia es fuertemente anisótropa, como por ejemplo en zonas de recirculación o en regiones con gradientes de presión negativos. Un intento de mejorar la modelización de la turbulencia en estas condiciones es la utilización de modelos de tensiones algebraicas. Un posible modelo de tensiones algebraicas promediado en profundidad que tienen en cuenta la producción de turbulencia por fricción de fondo calcula las tensiones de Reynolds promediadas en profundidad como (Cea, 2007): 2 1 − c 2 k 4 ∂U x 2 2 ∂U y 2 2 ∂U y 4 ∂U x 2 1 u' v'+ Puu, V − Pvv,V + k u' + v' + − − c11 ε 3 ∂x 3 ∂y 3 ∂y 3 3 3 3 ∂x 2 1 − c 2 k 2 ∂U x 2 4 ∂U y 2 2 ∂U x 4 ∂U y 2 1 u' v'+ Pvv, V − Puu, V + k v' 2 = u' − v' + − c11 ε 3 ∂x 3 ∂y 3 ∂x 3 3 3 3 ∂y 1 − c 2 k ∂U y 2 4 ∂U x 2 ∂U x ∂U y u' v'+ Puv, V u' v' = u' − v' − + − c11 ε ∂x 3 ∂y ∂y ∂x u' 2 =
con los constantes: c 2 = 0.6
c11 = 0.8 +
Pk ε
La producción de energía cinética turbulenta Pk se calcula utilizando la viscosidad turbulenta como:
∂U 2 1 ∂U ∂U y x x Pk = 2ν t + + 2 ∂y ∂x ∂x
∂U y + ∂y 2
2
Puu, V + Pvv, V + 2
con la viscosidad turbulenta calculada a partir de un modelo k-ε. Los términos Puu,V Pvv,V Puv,V modelan la producción de tensiones turbulentas por fricción de fondo, y se calculan a partir de las siguientes expresiones:
Puu, V = 2
u f2 U 2x h U
Pvv, V = 2
u f2 U 2y h U
Puv, V = −2
u f2 U x U y h U
Es necesario utilizar algún limitador que impida que las tensiones normales de Reynolds tomen valores negativos, ya que no existe nada en el modelo que asegure que su valor permanezca mayor que cero. Para ello es necesario establecer las siguientes restricciones en el modelo:
u' 2 ≥ 0
v'2 ≥ 0
u' 2 + v'2 ≤ 2k
La tercera condición impide que la tensión normal turbulenta en dirección vertical tome valores negativos, al mismo tiempo que evita valores excesivamente elevados de las tensiones normales horizontales. Si el modelo predice alguna tensión normal negativa, su valor se fuerza a cero, y al mismo tiempo se le resta la misma cantidad de energía introducida en la tensión normal vertical a las tensiones normales horizontales, de forma que la energía turbulenta total se mantenga constante. 3.11. Análisis dimensional de las ecuaciones de aguas someras Si se adimensionalizan las ecuaciones de aguas someras se obtienen los siguientes números adimensionales:
F=
U gh
T=
1 H Cf L
Rl =
UL ν
Rt =
UL νt
los cuales hacen referencia respectivamente a la fuerza de presión (F), a la fuerza de rozamiento del fondo (T), a las tensiones tangenciales laminares (Rl) y a las tensiones tangenciales turbulentas (Rt). La importancia relativa de las fuerzas asociadas a cada número adimensional es inversamente proporcional a la magnitud de dicho número, i.e. cuanto mayor sea un número adimensional, menor será la fuerza asociada a dicho número. Así, para un número de Reynolds laminar elevado, el flujo es turbulento y las fuerzas laminares pierden importancia en el desarrollo del flujo. De igual manera, la importancia de las tensiones turbulentas en la velocidad media dependerá de la magnitud del número de Reynolds turbulento, el cual depende de la viscosidad turbulenta. Se puede realizar una estimación del orden de magnitud de la viscosidad turbulenta a partir del modelo parabólico como: νt ≈
1 1 κ u f h ≈ κ g n h 5/6 U ≈ 0.21 n h 5/6 U 6 6
en donde se ha utilizado la fórmula de Manning para estimar la velocidad de fricción del fondo uf. Esta estimación será más precisa en casos en los que la turbulencia esté generada fundamentalmente por fricción de fondo, como puede ser el caso de ríos, y se alejará más del valor real en casos en los que la turbulencia esté generada principalmente por tensiones de corte horizontales, como por ejemplo en zonas de recirculación. En cualquier caso, utilizando dicha aproximación el número de Reynolds turbulento queda como:
Rt =
UL 4.8 L ≈ νt n h 5/6
Esta expresión se puede utilizar en primera instancia para evaluar la importancia de los esfuerzos turbulentos en el campo de velocidad y calado. Por ejemplo, si estamos modelando un tramo de río con calados del orden de 10m, una sección de 400m de anchura, un coeficiente de Manning estimado de 0.025, y una velocidad media de 0.5m/s, se obtiene una viscosidad turbulenta aproximada de 0.02m2/s y un número de Reynolds turbulento igual a Rt ~ 11000, el cual es un valor bastante elevado, por lo que es de esperar que las tensiones turbulentas tengan un efecto menor en el desarrollo del flujo medio.
4. Programas comerciales Como se ha comentado en la introducción, la modelización de la turbulencia en flujos poco profundos no ha sido tratada tan ampliamente como en otras disciplinas en dinámica de fluidos. Existen numerosos modelos numéricos diseñados para la simulación de flujos poco profundos mediante la resolución de las ecuaciones de aguas someras bidimensionales, algunos de los cuales incluyen modelo de turbulencia, otros que únicamente permiten fijar un valor constante de la viscosidad turbulenta, y otros que no incluyen el efecto de la difusión turbulenta en las ecuaciones de flujo. Prácticamente todos ellos han surgido del desarrollo llevado a cabo en grupos de investigación de distintas universidades o centros de investigación. Probablemente los más extendidos en la actualidad sean los modelos SMS (Brigham Young University) y MIKE21 (Danish Hydraulic Institute). Ambos se comercializan a nivel mundial, estando presentes en un gran número de consultorías. Tienen la principal ventaja de ser muy amigables desde el punto de vista del usuario, y de permitir una fácil definición del problema y un post-proceso de datos relativamente sencillo y eficiente. El resto de modelos no se comercializan de manera tan extendida, aunque la mayoría de ellos se pueden conseguir directamente a través de las universidades que los han desarrollado. Tienen el principal inconveniente de no tener un pre-proceso y post-proceso de datos tan amigables para el usuario no habituado a su utilización. Sin embargo, en lo que se refiere al cálculo numérico y resolución de las ecuaciones del modelo, muchos de ellos tienen las mismas capacidades o incluso superiores a los modelos comerciales. Simplemente por citar alguno de ellos, aunque la lista de programas disponibles es más larga, se pueden mencionar los siguientes: TELEMAC (LNHE-EDF), CARPA (Grupo Flumen, Universidad Politécnica de Catalunya), CCHE2D (University of Mississippi), DAMFLOW2D (Universidad de Málaga y Universidad de Sevilla), TURBILLON (Grupo Giama, Universidad de A Coruña). El modelo Surface-Water Modelling System (Brigham Young University) es un conjunto de módulos diseñados para la modelización de flujos en lámina libre. En la versión SMS 8.0 del año 2002, incluye los términos de difusión turbulenta. El módulo RMA2, para calcular flujo en régimen lento, permite al usuario definir una viscosidad turbulenta anisótropa (aunque recomienda definirla de forma isótropa), que además puede tomar valores diferentes en las distintas zonas del dominio modelado, en función de las características esperadas del flujo. Incluye 3 formas de calcular la dispersión
turbulenta: viscosidad turbulenta constante, viscosidad turbulenta basada en el número de Peclet, y el modelo de Smagorinski. Tanto el modelo de Smagorinski como el número de Peclet se utilizan en algunos modelos de aguas poco profundas más bien como formas de estabilizar numéricamente el modelo que como formas de modelar realmente la turbulencia en el flujo. El número de Peclet se define como: P=
U ∆x νt
siendo ∆x el tamaño de malla. Se trata de imponer una viscosidad turbulenta en cada punto de la malla tal que el valor del Peclet esté comprendido entre aproximadamente 15 y 40. Por su parte el modelo de Smagorinski es similar al modelo de longitud de mezcla pero asumiendo que la longitud característica de los remolinos turbulentos son proporcionales al tamaño de malla. ∂U ∂V ∂U ∂V + 2 ν t = C m A + + ∂x ∂y ∂x ∂y
siendo A el área del elemento de la malla, y Cs un coeficiente de Smagorinsky. Como puede observarse no tiene en cuenta el efecto del fondo ni del calado. El manual de usuario de RMA2 proporciona una tabla con recomendaciones para fijar el valor de la viscosidad turbulenta en función del tipo de flujo, los cuales son demasiado elevados (Tabla 1). Esto conlleva todos los inconvenientes comentados anteriormente en la sección 3.3. , lo cual conduce muy habitualmente a resultados erróneos en situaciones donde la turbulencia es importante. De hecho, el manual de usuario ya advierte que los valores de viscosidad turbulenta son utilizados por el programa “para ayudar a la convergencia numérica así como para calibrar el modelo”.
Tipo de problema Flujo homogéneo horizontal alrededor de una isla Flujo homogéneo horizontal en una confluencia Descarga térmica a un río en régimen lento y permanente Flujo de marea en un estuario con zonas pantanosas Flujo lento en un lago/embalse poco profundo
νt (m2/s) 0.5 - 5.0 1.1 - 5.0 1.0 – 50.0 2.5 - 10.0 0.01 – 0.05
Tabla 1. Valores de viscosidad turbulenta recomendados en el modelo RMA2.
También dentro del modelo SMS, el módulo HIVEL para el cálculo de flujo en régimen rápido, incorpora el modelo parabólico, permitiendo al usuario introducir un coeficiente multiplicador de ajuste de la viscosidad turbulenta. Como se ha visto en la sección 3.4. , en este modelo la viscosidad turbulenta se calcula localmente a partir del módulo de la velocidad, del coeficiente de fricción de fondo (Manning) y del calado. El módulo hidrodinámico del modelo MIKE21 resuelve las ecuaciones de St.Venant bidimensionales incluyendo como particularidad los efectos del oleaje en las corrientes. En lo que se refiere a la turbulencia, incluye 2 formas de tratarla: o usar una viscosidad turbulenta constante, y el modelo de Smagorinski. Tanto en SMS como en MIKE21, los términos de dispersión turbulenta se tratan más como una forma de estabilizar
numéricamente el modelo que como una forma de calcular la turbulencia real en el flujo.
5. Aplicaciones A modo de ejemplo se presentan a continuación los resultados obtenidos mediante la modelización matemática con las ecuaciones de aguas someras bidimensionales de diferentes flujos, incluyendo en todos los casos los esfuerzos turbulentos. Todos los resultados numéricos presentados se han obtenido con el código TURBILLON, desarrollado y utilizado actualmente por el Grupo de Ingeniería del Agua y del Medioambiente de la Universidad de A Coruña, el cual incluye todos los modelos de turbulencia presentados en la sección 3. Más que una descripción detallada de cada ejemplo, se pretende poner de relevancia los efectos del modelo de turbulencia utilizado en los resultados numéricos, así como la validación experimental de los mismos. La comparación de resultados se centrará en velocidades medias, ya que es ahí donde la turbulencia presenta una mayor influencia. Los campos de calado son menos sensibles a la modelización de la turbulencia. 5.1.
Flujo de marea en el estuario Crouch (Essex, Reino Unido)
El estuario Crouch tiene una longitud de aproximadamente 25Km y una anchura del cauce en la desembocadura de unos 1Km. Es un estuario que está totalmente dominado por la marea, siendo las aportaciones de agua dulce despreciables frente al caudal de marea. El flujo es por lo tanto no estacionario, y los valores de viscosidad turbulenta variables en función de la fase de marea. Se estudió la influencia del modelo de turbulencia en los resultados de calado y velocidad, concluyendo que en este caso la influencia de la turbulencia es despreciable en el flujo medio. Los resultados obtenidos son independientes del modelo de turbulencia utilizado, e incluso puede no utilizarse ningún modelo en absoluto (νt = 0). Asumiendo un calado de 10m, un coeficiente de Manning de 0.02, y una velocidad media de 1.5m/s, se obtiene una viscosidad turbulenta aproximada de 0.04m2/s. Esta estimación es del mismo orden de magnitud que los valores máximos de viscosidad turbulenta obtenidos con el modelo k-ε, lo cual confirma que la turbulencia está fundamentalmente generada por la fricción del fondo. Asumiendo una escala espacial de 1Km (que es aproximadamente la anchura del estuario en la desembocadura), se obtiene un número de Reynolds turbulento de Rt ~ 37000, confirmando que el efecto de la difusión turbulenta en el flujo medio es pequeño.
Figura 5. Campo de velocidades en el estuario Crouch (marea entrante).
Figura 6. Campo de viscosidad turbulenta con marea entrante (izquierda) y serie temporal de velocidades (numérica con distintos modelos de turbulencia y experimental) en un punto del cauce principal (derecha).
5.2.
Escalas de peces de hendidura vertical
Este es un ejemplo paradójico en donde un modelo de aguas poco profundas proporciona resultados muy aceptables del campo de velocidades y calados a pesar de que el calado y la escala horizontal son del mismo orden de magnitud. El hecho de que un modelo 2D-SWE funcione bien se debe a que el flujo es marcadamente bidimensional, las velocidades verticales son muy pequeñas excepto en la hendidura, y el campo de velocidad horizontal es muy homogéneo en profundidad. Es un flujo extremadamente turbulento, por lo que se hace indispensable utilizar un modelo de turbulencia adecuado, tipo k-ε o similar. La Figura 7 muestra el campo de velocidad y líneas de corriente calculadas en 2 diseños diferentes de piscinas. La correcta evaluación de las zonas de recirculación que aparecen en la piscina depende del modelo de turbulencia utilizado, como se puede apreciar en la Figura 8.
Figura 7. Campos de velocidad en dos tipos de escala de peces diferentes.
Figura 8. Secciones transversales de velocidad longitudinal, experimental y numérica, utilizando varios modelos de turbulencia
5.3.
Flujo en canal con codo de 90º
En este ejemplo se evalúa la zona de recirculación creada aguas abajo de un codo de 90º. El flujo es mucho más lento y menos turbulento que en la escala de peces. A pesar de ello, la presencia de una zona de recirculación hace necesaria la utilización de un modelo de turbulencia, si bien en este caso las diferencias entre los diferentes modelos utilizados son menores. Otra particularidad de este ejemplo es que aparecen flujos de recirculación secundarios aguas abajo del codo (ver Figura 3) cuyo efecto no se incluye en los modelos 2D-SWE, que desprecian los términos de dispersión (sección 2). A pesar de ello la comparación numérico experimental de velocidades es muy buena (Figura 9 y Figura 10).
Figura 9. Campos de velocidad experimental (izquierda) y numérico (derecha). La zona de recirculación generada aguas abajo del codo se muestra en color blanco.
Figura 10. Secciones transversales de velocidad longitudinal experimental y numérica calculados con diferentes modelos de turbulencia.
5.4.
Flujo en un tanque de floculación
Este es otro ejemplo de un flujo lento en el que aparecen grandes zonas de recirculación. Se modela el flujo en un tanque de floculación en el cual aparecen dos remolinos de recirculación de diferente tamaño. La relación entre calado y escala horizontal es del orden de H/L ~ 0.1. Para el cálculo de las tensiones de Reynolds se utilizan los modelos de turbulencia de longitud de mezcla, k-ε, y de tensiones algebraicas para poner de relieve la importancia de una correcta modelización de la turbulencia en los resultados (Figura 11). Aunque no se presentan resultados experimentales de velocidad, las líneas
de corriente más parecidas a las observadas experimentalmente son las obtenidas con los modelos k-ε y DASM. Si no se utiliza ningún modelo de turbulencia los resultados son totalmente incorrectos. Si se utiliza un modelo sencillo, como el de longitud de mezcla, los valores de viscosidad turbulenta son demasiado pequeños, y por lo tanto el tamaño de la zona de recirculación izquierda demasiado grande.
Figura 11. Campo de velocidad numérico. Sin modelo de turbulencia (izquierda), modelo de longitud de mezcla (medio), modelo de tensiones algebraicas DASM (derecha).
Figura 12. Campo experimental de velocidad en un tanque de floculación.
6. Conclusiones Se han presentado los siguientes modelos de turbulencia para aguas someras, por orden creciente de complejidad:
Viscosidad turbulenta constante Modelo parabólico Modelo de longitud de mezcla Modelo k-ε de Rastogi y Rodi (Rastogi y Rodi 1978) Modificación del modelo k-ε propuesta por Booij Modelo k-ε de Babarutsi y Chu Modelo de tensiones algebraicas para aguas someras
Así como otras formas de tratar la turbulencia que aparecen en algunos programas comerciales, como el modelo de Smagorinsky y el número de Peclet. Estos dos últimos,
más que como modelos de turbulencia se utilizan en algunos modelos de aguas poco profundas para estabilizar numéricamente el modelo. La utilización de modelos de turbulencia de diferente complejidad permite seleccionar el más adecuado en cada caso de estudio. En general el modelo de longitud de mezcla proporciona resultados satisfactorios en ríos y estuarios, pudiendo incluso llegar a no ser necesario utilizar ningún modelo de turbulencia en dichos casos. En estructuras hidráulicas como canales en lámina libre con codos pronunciados y zonas de recirculación, suele ser necesario utilizar por lo menos un modelo de longitud de mezcla, pudiendo ser necesario utilizar un modelo k-ε o de tensiones algebraicas. La elección del modelo de turbulencia que mejor se adecua a cada caso se realiza en base a la experiencia del usuario, teniendo siempre en cuenta que cuanto más complejo es el modelo mayor es el tiempo de cálculo y más compleja la resolución de las ecuaciones A la vista de los diferentes modelos de turbulencia presentados se plantea la duda de cuál es el más preciso, el más correcto o el más adecuado para un problema dado. No existe una respuesta absoluta a dichas cuestiones. Cabe recordar que todos ellos son modelos, con ciertas hipótesis y constantes experimentales, y el hecho de que no exista un modelo universal, se debe a que ninguno de ellos es el mejor en términos absolutos. Este es precisamente el problema de la turbulencia. En general los más complejos son más precisos, pero conllevan un coste computacional mayor, por lo cual puede llegar a ser recomendable en algunas situaciones utilizar uno más sencillo. Para modelos de similar nivel de complejidad, como es el caso de las 3 versiones del k-ε y el de tensiones algebraicas, lo más que se puede hacer es resolver el mismo problema con varios modelos y ver la sensibilidad de los resultados al modelo utilizado.
7. Referencias Babarutsi, S., Chu, V.H. (1991). A two-length-scale model for quasi-two-dimensional turbulent shear flows. Proc. 24th Congr. of IAHR, Vol. C, Madrid. IAHR Booij, R. (1989). Depth-averaged k-ε modelling. Proc. 23rd IAHR Congr., Delft, The Netherlands. vol. A, IAHR Cea, L., Pena, L., Puertas, J., Vázquez-Cendón, M.E, Peña, E. (2007). Application of Several Depth Averaged Turbulence Models to Simulate Flow in Vertical Slot Fishways. Journal of Hydraulic Engineering, Vol.133 (2) pp.160-172 Cea, L., French, J.R., Vázquez-Cendón, M.E. (2006). Numerical modelling of tidal flows in complex estuaries including turbulence: An unstructured finite volume solver and experimental validation. International Journal for Numerical Methods in Engineering, Vol.67 (13) pp.1909-1932 Rastogi, A.K., Rodi, W. (1978). Predictions of heat and mass transfer in open channels. Journal of the Hydraulics Division HY3 Rodi, W. (1993). Turbulence models and their application in hydraulics. A state-of-theart review. A.A. Balkema. IAHR