CAPITULO 4

Modelación Bidimensional 4.1 Objetivos de la Modelación

Dar a conocer la utilización del software River2D, en la modelación numérica del Río Huancabamba.

Lograr un mayor conocimiento del alcance de los modelos numéricos en la

aplicación de la Ingeniería Hidráulica.

Calcular detalladamente la distribución de velocidades, variación de tirantes, número de fraude, así como otros parámetros, en el tramo a estudiar, para evaluación de diferentes caudales.

4.2 Tipos de modelos numéricos El modelo es la representación de algún objeto o sistema3, en un lenguaje o forma de fácil acceso y uso, con el objetivo de poder entenderlo y poder generar sus respuestas para diferentes entradas. Los modelos se clasifican en físicos, analógicos, numéricos e híbridos. Los modelos físicos representan el sistema físicamente a escala o bajo circunstancias específicas. Por otro lado, los modelos analógicos permiten el estudio de un modelo complejo, mediante un sistema más simple con características análogas, explicando sus interrelaciones como una analogía eléctrica y/o mecánica. Los modelos numéricos, que es de importancia en esta tesis, tratan de describir el sistema natural por medio de fórmulas matemáticas y secuencias lógicas. Existen también los modelos híbridos que utilizan una combinación del modelamiento matemático y analógico. Los primeros modelos numéricos se desarrollaron exclusivamente para resolver un problema particular, generalmente el flujo del agua. Actualmente los modelos numéricos son genéricos con capacidad para resolver una gama de condiciones de flujo similares, con módulos adicionales para resolver otros fenómenos como transporte de sedimentos, dispersión de contaminantes, calidad de agua e inclusive modelación de hábitat de peces u otras criaturas acuáticas. Existen muchas maneras de clasificar los modelos numéricos, aquí se darán algunas de ellas según algunos criterios. ---------------------------------------------------------------------------------------------------------- 3 El sistema es un conjunto de procesos físicos, químicos y/o biológicos que actúan sobre unas variables de entrada para convertirlas en variables de salida

82

4.2.1 Clasificación Los modelos para simular los diferentes elementos de una cuenca pueden ser clasificados de acuerdo a diferentes criterios. Rafaelli (2003), los clasifica de la siguiente manera: - En base a la descripción del proceso. - En base a la técnica de resolución. 4.2.1.1 En base a la descripción del proceso Concentrado Son los modelos simplificados que trabajan con valores integrados espacialmente. Distribuido La región de estudio se puede discretizar basándose en la definición de una grilla regular donde se realiza el cálculo distribuido o teniendo como base las líneas de escurrimiento y las curvas de nivel. La tecnología moderna (sensores remotos, sistemas de información geográfica – SIG, procesamiento de imágenes satelitales, modelos digitales de terreno, etc) están permitiendo la evolución de los modelos clásicos del tipo de sistemas concentrados, hacia aquellos más complejos con estructura distribuida. Estocásticos y determinísticos En los modelos estocásticos se considera que las variables son aleatorias con distribución de probabilidades y con dependencia en el tiempo. Se considera los modelos determinísticos que corresponden a una ley definida de certidumbre.

4.2.1.2 Según las técnicas de resolución Modelos de caja negra El análisis se limita simplemente a la correlación de entradas y salidas del sistema. Se conocen también como empíricos o analíticos y se basan en la experimentación. A este grupo pertenecen los modelos basados en la inteligencia artificial como las redes neuronales, los algoritmos genéticos y la lógica difusa. Modelos de caja gris Pueden identificar alguna de las relaciones causa – efectos dentro del sistema. Se conocen también como conceptuales y se representan los procesos físicos que relacionan las entradas con las salidas. Modelos de caja blanca Establecen ecuaciones simplificadas de todos los procesos que intervienen en el fenómeno. En este grupo se encuentran los modelos numéricos que resuelven las ecuaciones

83

diferenciales parciales que explican el flujo superficial y subterráneo. Se conocen también como físicamente basados y consideran las ecuaciones de conservación de masa, energía y movimiento, para describir el movimiento de agua en la superficie y en el subsuelo, zona saturada y no saturada, usando técnicas de discretización numérica como el Método de los Elementos Finitos (MEF), Método de Diferencias Finitas (MDF) o el Método de Volúmenes Finitos (MVF). Estos métodos resuelven las ecuaciones que tienen dominio continuo mediante la solución en un número finitos de puntos discretos en dicho dominio. Cuando los valores buscados (tirante, velocidad, etc.) en dichos puntos discretos son encontrados, la solución en cualquier otro punto puede ser aproximada mediante métodos de interpolación. Decidir cuál de los tres métodos, MDF, MEF o MVF es el mejor es aún un tema en debate. En general, para geometrías simples el MDF es una excelente alternativas pues es muy rápido y simple; sin embargo geometrías complejas exigen a su vez mallas complejas que se modelar mejor con el MEF o MVF.

4.2.2 Dimensionalidad En los modelos unidimensionales (1D) se asume que una de las dimensiones prevalece sobre las otras dos. Esta dimensión es la longitudinal a lo largo del eje del río o canal. La información topográfica e hidráulica se introduce mediante secciones transversales, en las cuales se calculan el tirante y velocidad promedios en toda la sección transversal. Es decir, toda la sección es representada por un único valor medio de velocidad, no considerándose variaciones en la distribución de velocidades tanto horizontal como verticalmente. Asumen por defecto que el flujo es perpendicular a la sección transversal, lo cual es una de sus limitaciones. Estos modelos son aplicables en tramos de ríos y canales muy largos, generalmente mayores a 20 veces el ancho y cuando se busca principalmente determinar el máximo nivel de agua (para diseño de diques por ejemplo).

Figura 4.1: Esquema de modelos hidráulicos unidimensional (X), bidimensional (X,Y) y tridimensional (X,Y,Z).

84

Los modelos bidimensionales (2D) consideran las variaciones en las dos dimensiones del plano horizontal. Las variaciones de la velocidad u otra magnitud de interés en la columna vertical de agua se promedian y se asumen como un único valor. Estos modelos son especialmente útiles en flujos muy extendidos (como estuarios, lagos, etc.) donde la variación vertical de velocidad es pequeña, por eso suelen llamarse modelos de aguas someras o poco profundas. Estrictamente no son aplicables a casos en que la variación vertical de la velocidad es apreciable, como por ejemplo el flujo sobre un vertedero o a través de un orificio sumergido; sin embargo mediante uso de expresiones empíricas o similares pueden incorporar estas singularidades dentro de la modelación. Los modelos tridimensionales (3D) representan el estado más avanzado de la modelación. Estos modelos son capaces de calcular las tres componentes espaciales de la velocidad, y por tanto aplicables a cualquier caso práctico. 4.2.3 Método de cálculo El agua es un medio continuo constituido por un número casi infinito de partículas. Determinar la velocidad de cada partícula es virtualmente imposible y poco práctico; por ello este medio continuo se divide en elementos discretos de tamaño finito, que la memoria de una computadora es capaz de manejar. En los casos más simples de flujo 1D, la discretización se realiza a nivel de secciones transversales, siendo el método de cálculo por etapas: se hace el balance energético en una sección y se procede a la siguiente, hasta terminar con todas. Si hay problemas de convergencia en una sección, se asume el tirante crítico y se continúa con la siguiente. Por ello estos programas son robustos, rápidos y numéricamente confiables. Los problemas 2D y 3D requieren resolver las ecuaciones diferenciales y por tanto una discretización más fina mediante el uso de mallas. Los métodos principales de solución, como ya se ha mencionado anteriormente son el Método de Elementos Finitos (MEF), el Método de Diferencias Finitas (MDF) y más recientemente el Método de Volúmenes Finitos (MVF). Los modelos numéricos de 2D y 3D proporcionan mucho detalle y precisión, siendo capaces de manejar condiciones de borde e iniciales complejas, para los cuales no existen en la mayoría de casos soluciones analíticas. Sin embargo, la mayoría de los métodos de solución requieren la solución simultánea de todos los nodos del dominio, lo que incrementa la carga computacional haciéndolos más lentos. También pueden presentarse problemas de convergencia numérica que puede conducir a soluciones irreales, por lo que el buen criterio del ingeniero siempre debe estar presente. 4.2.4 Regímenes de flujo El flujo subcrítico (Números de Froude Fr <1) es simple y todos los programas tienen capacidad para modelarlo. El flujo supercrítico (Número de Froude Fr >1) resulta más complejo para muchos modelos numéricos, por la posibilidad que se forman ondas de choque o saltos hidráulicos en los cambios a régimen subcrítico, lo cual da lugar a un flujo rápidamente variado de difícil modelación.

85

Además, en rigor el flujo supercrítico en tramos largos generalmente se produce en canales artificiales revestidos; ya que en la naturaleza los ríos aluviales ajustan su pendiente para evitar la aparición del flujo supercrítico en tramos muy largos, pues en el fondo se producen escalonamientos que forman una sucesión de pozas y rápidos con una mezcla de flujos subcrítico y supercrítico alternados, que los programas numéricos difícilmente pueden modelan con precisión. 4.2.5 Variación en el tiempo En el flujo permanente se asume que las todas las magnitudes hidráulicas son constantes en el tiempo. Esta es la opción por defecto de todos los modelos numéricos. En el flujo impermanente o transitorio las magnitudes como tirantes o velocidades pueden cambiar con el tiempo, como en el caso del tránsito de una onda de avenidas a través de un río. Algunos programas tienen capacidad para modelación en régimen impermanente, siendo especialmente útiles para estudios de transporte se sedimento, pues la erosión y deposición modifican gradualmente la sección hidráulica. Sin embargo, los flujos impermanentes muy bruscos, como los originados por la rotura de presas, requieren de tratamientos especiales. 4.2.6 Contornos del cauce Las expresiones clásicas de la hidráulica consideran los contornos o fronteras del canal como lecho rígido indeformable, lo cual es cierto para muchas estructuras artificiales; pero no para cauces aluviales, en los cuales se presenta lecho móvil. La modelación de ríos aluviales, especialmente los de lecho fino, requieren casi siempre capacidades de modelación en lecho móvil y flujo impermanente para lograr reproducir los cambios reales del fondo. La deformación del cauce durante una avenida es un fenómeno extremadamente complejo, y a pesar que no existen métodos universalmente aceptados para del transporte de sedimento, si existen modelos numéricos con opciones para calcular erosión y sedimentación, ya sea general o local, aunque casi siempre requieren ser calibrados con datos medidos en campo, lo cual limita mucho su aplicación práctica. 4.3 Modelación numérica de profundidad promedio 4.3.1 Introducción Los adelantos en la capacidad de la computadora y la tecnología del software están haciendo el análisis detallado más rutinario en casi todas ramas de la ingeniería. En canales y ríos bidimensionales (2D), los modelos de profundidad promedio son el inicio para juntar a los modelos unidimensionales en la práctica común. Estos modelos son útiles en los estudios dónde los detalles locales de la distribución de velocidad y profundidad son importantes. Los ejemplos incluyen diseño de puentes, encauzamientos, obras de derivación e incluso la evaluación del hábitat del pez. Esta introducción tiene la intención de dar una breve apreciación general de la modelación de ríos en 2D, resaltando las consideraciones para las aplicaciones prácticas. Con las posibles altas velocidades y grandes pendientes, y relativamente de poca profundidad, los ríos y modelos de flujos presentan una difícil particularidad al desafío

86

computacional. Este hecho es probablemente un factor significante en el retraso de aplicación de modelos de agua poco profundos en ríos comparados a los ríos costeros y a problemas de estuarios. Recientemente, no obstante, se desarrolló originalmente técnicas numéricas para problemas de la mecánica fluidos transitorios que se han sido adaptadas con éxito. Estas técnicas, llamado shock-capturing o esquemas de high-resolution , permite libremente intermezclar el flujo subcrítico y supercrítico. Las aplicaciones en modelos de ríos en 2D normalmente enfocan una limitada magnitud del cauce, típicamente menor de diez veces el ancho del cauce en la longitud. En la mayoría de los casos, la variación de la descarga es relativamente lenta comparado al tiempo de viaje de agua a través del tramo, y las condiciones de estado permanentes siempre son perseguidas. La mayoría de los modelos son diseñados para dar soluciones transitorias, sin embargo, y a menudo la solución permanente deseada es obtenida como una asíntota a la solución transitoria (impermanente) después un largo lapso de tiempo. Dado un conjunto de ecuaciones gobernantes, hay dos pasos esenciales en el desarrollo de un modelo computacional: la discretización y la solución. En la discretización el número infinito de ecuaciones para un infinito el número de incógnitas es reducido a un número finito de ecuaciones en un número finito de malla o de puntos en el espacio y tiempo. A esta fase, los cálculos operacionales son reducidos a operaciones algebraicas. En la solución un esquema o proceso es ideado donde las ecuaciones algebraicas desarrolladas en el primer paso pueden resolverse para los valores nodales desconocidos. El álgebra es reducida a la aritmética que puede ser resuelta dentro del código de la computadora. Hay varias alternativas para cada paso. Métodos comunes de discretización incluyen diferencias finitas, volúmenes finitos y el método de elementos finitos. Los métodos de la solución incluyen soluciones explícitas e implícitas, este último depende de una variedad de métodos de solución de ecuaciones lineales iterativas o directas no lineales. 4.3.2 Modelos de profundidad promedio Hay diversos dominios comerciales y públicos de modelos 2D disponibles. Ellos son basado en una variedad de esquemas numéricos y ofrecen un rango gráfico pre y post procesador de módulos. Los fundamentos físicos son más comunes, sin embargo, todos los modelos 2D resuelven la ecuación básica de conservación de masa y los dos componentes (horizontales) de conservación de velocidad adquirida (impulso). Las salidas del modelo son dos componentes (horizontales) de velocidad y profundidad a cada punto o nodo. Las distribuciones de velocidad en la vertical son asumidas para ser uniformes y las distribuciones de presiones son asumidas para ser hidrostáticas. Como regla, una precaución debe ejercerse en cualquier esfuerzo por resolverse las características del flujo menor que aproximadamente diez veces la profundidad del flujo en la magnitud horizontal. Los diseños del modelo 2D se basaron en diferencias finitas, volúmenes finitos, y también del método de elementos finitos. Cada método tiene sus ventajas y desventajas. Al riesgo de groseramente afirmar que el método de volúmenes finitos ofrece la mejor estabilidad y eficacia mientras que el método de elementos finitos ofrecen la mejor flexibilidad geométrica.

87

La base para el método de elementos finitos es una técnica más general conocido como el método residual ponderado. La idea es que las ecuaciones que gobiernan pueden ser resueltas aproximadamente por el uso de una “función prueba” que se especifica pero tiene varios ajustes de los grados de libertad. En cierto modo, el proceso es análogo al ajuste de una curva, por una línea recta, a los datos observados. La línea recta es una función específica que tiene dos grados de libertad. Los valores para esos dos parámetros son buscados a que otorguen el menor error. Las ecuaciones de continuidad y de velocidad adquirida de profundidad promedio del flujo pueden escribirse en la forma abreviada como: ( ) 0,, =VUHC 4.1 0),,( =VUHM x 4.2 0),,( =VUHM y 4.3 Si nosotros introducimos las funciones de prueba para las variables H, U, y V; denotado por Hˆ, Uˆ, y Vˆ, y las sustituimos en las ecuaciones, ellas no satisfarán las ecuaciones exactamente. Un “residuo” resultará como:

RcVUHC =∧∧∧

),,( 4.4

RVUHM xx =∧∧∧

),,( 4.5

RVUHM yy =∧∧∧

),,( 4.6 El objetivo es ahora hacer los residuos tan pequeño como sea posible. El aproximar al residual ponderado es multiplicar el residuo por una función ponderada, integrando encima del área entera y colocar el resultado a cero. Para cada grado de libertad en la función de prueba, una separada función ponderada es usada para generar el mismo número de ecuaciones como de incógnitas. La ecuación de continuidad, para el ejemplo, se vuelve:

0),,( =∫∧∧∧

dAVUHCN i 4.7 donde Ni es la i-ésima función de prueba. La integración encima del área sirve para reducir la distribución espacial a sólo una cantidad. El truco con el método residual ponderado es escoger una buena función de prueba suficiente para que trabaje en diferentes situaciones. La idea de interpolación se pone muy útil aquí. Las funciones de interpolación pueden ser escritas en la forma de influencia o de funciones “bases” para cada punto, para que la función de prueba sea una suma de condiciones.

88

NHH j

n

jj∑=

=

∧

1

4.8 El método de elementos finitos toma el método residual ponderado un par de pasos más allá. Primero, el dominio es dividido en discretas áreas, llamados elementos finitos que se definen por las conexiones, entre los puntos nodales. La triangulación que se muestra debajo es un el ejemplo.

Figura 4.2 : Ejemplo de triangulación

Los nodos están en los vértices de los triángulos. La interpolación simple o las funciones bases son fáciles de construir, permitiendo variaciones lineales entre los nodos, por ejemplo. Estas funciones bases son consideradas a ser “locales” en ellos sólo se extienden de un nodo particular cerca a los nodos vecinos. Fuera de esa área, ellos son cero. Las funciones básicas para el nodo j pueden ser visualizadas como una función “tienda” con el polo (de altura la unidad) al nodo j, y fijo todos los nodos que son adjuntados al nodo j por un borde del triángulo. La operación de integración en el método residual ponderado es llevada a cabo como una suma de integrales del elemento en el método de elementos finitos. Por la naturaleza local de las funciones bases, sólo esas funciones con los nodos en el elemento particular necesitan ser integradas. Para un triángulo lineal, con tres nodos, un total de nueve integrales son requeridas (la combinación de tres funciones básicas con tres funciones de prueba), aún cuando hay miles de otros nodos. La ventaja más grande del método del elemento finito es la flexibilidad geométrica. Los elementos pueden cambiar el tamaño y forma prontamente, permitiendo bordes complejos a ser trazados, así como permitiendo el refinamiento de la malla rápidamente en áreas importantes. Cualquier método de discretización, siendo aproximado, engendra algunos errores. Teniendo una idea de la fuente y magnitud de ese error es importante una evaluación computacional de los resultados y tomar pasos para mejorarlos. Para asuntos de esta tesis, el error será definido como la diferencia entre el valor numérico a cualquier punto (no sólo a los nodos) y la solución analítica en ese punto. Para el método de elementos finitos, hay dos fuentes importantes de origen de errores relacionado a los funcionamientos llevados a cabo. El primero es el error en la interpolación. En la figura 4.3 se muestran dos curvas diferentes con interpolaciones lineales; es obvio que uno con más pedazos es mejor.

89

Fig. 4.3 Interpolación del error

El segundo origen del error es más sutil. El método residual ponderado realmente es un proceso promedio ponderado. El valor nodal a cualquier punto particular refleja las condiciones circundantes al punto. Si hay una variación grande en las condiciones, entonces el valor del punto puede ser influenciado de una manera u otra. Hay dos ejemplos comunes de este problema en el cálculo de flujos de profundidad promedio. El primero está cercano a la inclinación del fondo. El nodo al fondo del lecho refleja ambas condiciones del canal (velocidad uniforme relativamente) y el fondo (la velocidad más lenta). El promedio tiende a ser más bajo que la velocidad al punto específico. La solución a ambos tipos de errores de discretización es tener elementos más pequeños. Esto significa más elementos y más nodos para la misma área física. El arte de la modelación computacional es proporcionar suficiente refinamiento en las áreas importantes mientras mantenemos razonables tiempos de solución. El resultado del método de elementos finitos, o cualquier otro método de discretización, es establecer las ecuaciones algebraicas no lineales para todas las profundidades y velocidades desconocidas. El proceso de resolver éstas ecuaciones es lo que está exigiendo habilidad de la computadora. Más modelos de computadora resuelven el flujo de profundidad promedio en condiciones transitorias, aún si resultados de estados permanentes son deseados. Esto es una manera conveniente de proporcionar un control e iteración estable de una arbitraria primera suposición o condición inicial. Dos acercamientos son generalmente usados, referidos a los métodos explícito e implícito. Los métodos explícitos resuelven para las nuevas variables, valores en cualquier nodo basado en los valores de los nodos circundantes en el momento anterior (en el tiempo previo). Este método tiene la tremenda ventaja de calcular cada valor independientemente. Ninguna solución de la matriz es necesaria para que el almacenamiento sea mínimo y la ejecución sea rápida. La desventaja es que el time step es restringido para guardar una solución estable. La condición Courant - Lewy - Freidrichs (CFL) es una guía al maximum time step (máximo paso de tiempo) permitido.

+< ∆∆ gHV

xt min

4.9

90

Esta condición se relaciona al tiempo tomado por una onda de agua poco profunda para viajar de un nodo al próximo. Cuando el espaciamiento de la malla es hecho más pequeño en cualquier área, el time step para el proceso de la solución debe disminuirse. El segundo método es implícito. Aquí los valores de todas las variables en el nuevo tiempo son considerados a depender de otras así como los valores al paso de tiempo anterior. Desde que las incógnitas se interrelacionan, considerables esfuerzos usando métodos algebraicos lineales y no lineales se requieren para resolver las ecuaciones. Típicamente, las matrices grandes son formadas. Si hay N nodos en una malla, hay entonces 3N variables desconocidas para resolver durante cada paso de tiempo. La matriz llena que describe sus interconexiones habría de ser 3N X 3N. Si N es 10,000, entonces la matriz requeriría 7.2 gigabytes de memoria para almacenar. Afortunadamente, la matriz llena no es normalmente requerida. La mayoría de estos valores son ceros. La estimación mas realista de almacenamiento es 3N X 3B donde B es el ancho de banda de la malla. Para problemas de 10,000 nodos, los anchos de banda estarán en el orden de 100. La matriz requerida está entonces menos de 100 Megabytes - grande, pero manejable. El tiempo del cómputo para la solución de la matriz es aproximadamente proporcional a NB2. Si B se suma con la raíz cuadrada de N, entonces el tamaño de la matriz es proporcional a N3/2 y el tiempo de la solución es proporcional a N2. 4.3.3 Requerimiento de datos

Puede ser un cliché decir que un modelo sólo es tan bueno como la entrada de los datos, pero es verdad. Como en los datos de entrada, los modelos hidrodinámicos 2D requieren la topografía del fondo del cauce, rugosidad, así como condiciones de borde, y las condiciones iniciales del flujo. A demás, alguna clase de malla discreta debe diseñarse para capturar las variaciones del flujo. Obteniendo una representación exacta de la topografía del fondo es probablemente lo más crítico, difícil, y el tiempo que toma los ejercicios de modelación 2D es considerable. Los simples estudios de secciones transversales generalmente son inadecuados. Combinando GPS y los sistemas sonoros de profundidad para los ríos grandes y los estudios distribuidos de Estación Total para los ríos más pequeños tienen que ser buscados para ser efectivos. En cualquier evento, debe esperarse a gastar mínimo una semana para la colección de datos de campo para el estudio. Los datos del campo deben procesarse y deben verificarse a través de un modelo de calidad digital del terreno antes de usarse como la entrada para los modelos 2D. La rugosidad del fondo, en la forma de una altura de rugosidad o del valor del “n” de Manning, es un parámetro de entrada menos crítico. Comparado con modelos tradicionales unidimensionales dónde muchos efectos bidimensionales son resumidos en el factor de resistencia, el término de la resistencia bidimensional sólo responde al cortante directo del fondo. Las observaciones del material del fondo y del tamaño de las formas del fondo (dunas, bancos de arena, etc.) son normalmente suficientes para establecer estimaciones de rugosidad razonables. Si los valores de rugosidad requeridos son poco realistas, es probable que haya problemas con la topografía del fondo.

91

Las condiciones del borde normalmente toman la forma de una descarga total en las secciones de entrada y en las secciones de salida las elevaciones de superficie de agua son fijas. Colocando los bordes del flujo a alguna distancia de las áreas de interés es importante para minimizar el efecto de incertidumbre de las condiciones de borde. Las condiciones iniciales son importantes, incluso para el flujo permanente, desde que ellos son normalmente usados como la suposición inicial en el procedimiento de la solución iterativa. Una buena suposición reducirá significativamente el tiempo de ejecución total y puede hacerse la diferencia entre una ejecución permanente y una impermanente. El diseño de la malla o cuadricula es el arte de la modelación 2D. El número total de grados de libertad (el número de nodos, tres incógnitas por nodo) está limitado por la capacidad de la computadora y el tiempo disponible. El desafío es distribuir estos nodos de tal manera que la solución más exacta se obtenga para un propósito particular. Estrechamente los nodos espaciados en las regiones de alto interés o en regiones de variación de flujo, los cambios graduales en el espaciamiento del nodo, y regularidad de elemento o forma de la celda son consideraciones importantes. Para eliminar las preocupaciones sobre los efectos del diseño de la malla en la solución final, las soluciones de más de un diseño deben ser comparadas. Afortunadamente, los pre-procesadores gráficos están disponibles para ayudar al diseño (y después para interpretar los resultados). Desafortunadamente, éstos requieren mucho más esfuerzo de la programación, y por consiguiente puede ser más caro, que los mismos modelos 2D hidrodinámicos. 4.4 Selección del modelo numérico a utilizar Existe una diversidad de modelos numéricos que simulan diferentes fenómenos hidráulicos e hidrológicos. Cada software ha sido creado especialmente para simular fenómenos y eventos específicos y viene con un formato definido y entornos definidos por el diseñador (DOS, WINDOWS, etc). Para la ejecución de estos programas se necesita obtener la licencia de uso. Esta licencia exige un costo el cual es más elevado cuando el software es más sofisticado. Sin embargo existe también programas que se obtienen en forma gratutita y tienen libre distribución ofreciendo también resultados muy confiables para obtener parámetros útiles para el diseño. Los modelos que simulan fenómenos o eventos hidrológicos e hidráulicos, se pueden agrupar en modelos de medio ambiente y contaminación, modelos de agua superficial y modelos de agua subterránea. Los modelos de agua superficial se pueden dividir en modelos de flujo hidrodinámicos, modelos de circulación costera, modelos de circulación en lagos, modelos de inundaciones, modelos de análisis de redes de distribución de agua, modelos de ingeniería sanitaria, modelos de hidrología y modelos de calidad del agua superficial. Por ser de interés en este capítulo se tratará brevemente algunos modelos de agua superficial que simulan flujo hidrodinámico para luego justificar el uso del River 2D en esta tesis.

92

4.4.1 AquaDyn AquaDyn es un modelo hidrodinámico de gran alcance y de fácil uso para la simulación de los recursos de agua, factor de riesgo, y estudios de impacto. AquaDyn permite la descripción y el análisis completo de condiciones hidrodinámicas (caudales y niveles de agua) de canales abiertos tales como ríos, lagos, o estuarios. Los ingenieros especialistas y los responsables pueden utilizar los módulos especializados de la simulación de AquaDyn para predecir impactos en condiciones particulares del flujo del agua. Por ejemplo, AquaDyn proporciona una manera confiable de pronosticar las consecuencias de diversas actividades tales como construcción de diques, embarcaderos de los puentes, y terraplenes. AquaDyn puede modelar flujo permanente e impermanente, tanto para régimen crítico o supercrítico y por lo tanto permiten que el usuario considere y estudie los efectos de vertederos, de contracciones, y de ondas de marea. La licencia de uso cuesta US$ 1975 4.4.2 AQUASEA Este software corre en entorno de Windows y es de gran alcance y fácil uso. Utiliza el método de elementos finitos como método de solución para simular problemas de marea en flujos de ríos, estuarios y en áreas costeras. Además puede simular problemas en la circulación de agua en lagos y problemas que implican transporte de sedimentos suspendidos. AQUASEA contiene un modelo hidrodinámico de flujo y un modelo de transporte-dispersión. El programa puede simular variaciones de flujos y de nivel de agua en respuesta a varias funciones que gobiernan en los ríos, lagos, estuarios, bahías y áreas costeras. El modelo de transporte-dispersión simula la separación de una sustancia en el ambiente bajo influencia del flujo fluido y de los procesos existentes de la dispersión. La sustancia puede ser un agente contaminador de cualquier sal buena, conservadora o no-conservadora, inorgánica u orgánica, sedimento suspendido, oxígeno disuelto, fósforo inorgánico, nitrógeno y de otros. Su costo es de US$ 1800. 4.4.3 HEC-RAS Es un software que tiene un sistema integrado que realiza cálculos hidráulicos unidimensionales (1D) para una red completa de canales naturales y construidos. El sistema se abarca de un interfaz utilizador gráfico (GUI), componentes del análisis, capacidades del almacenaje de datos y de gerencia y gráficos hidráulicos. El HEC-RAS simula para flujo permanente e impermanente. HEC-RAS se ha desarrollado para el cuerpo de ejército de ESTADOS UNIDOS de ingenieros. Sin embargo, el software desarrollado en el centro hidrológico de la ingeniería se pone a disposición el público siempre que sea apropiado. El uso no es restricto y los individuos fuera del cuerpo de ingenieros pueden utilizar el programa sin carga. Se desarrolla en Windows y su distribución es forma libre y gratuita. 4.4.4 HEC-2 El Módulo de Interfaz de HEC-2 proporciona una interfaz a modo de AutoCAD para uso del motor de cálculo de perfiles de la superficie del agua del U.S. Army Corps of Engineers. El Módulo de Interfaz de HEC-2 incluye un motor de cálculo de HEC-2 y suministra todas las herramientas para crear el modelo del perfil de la superficie del agua. BOSS CAD permite al usuario cambiar entre HEC-2 y HEC-RAS, o viceversa, rápidamente. El motor de cálculo de HEC-2 es un modelo bidimensional, para estado de

93

flujo gradualmente variado y constante. Pueden calcularse por separado los perfiles subcríticos y supercríticos. El modelo puede tener en cuenta el remanso creado por puentes, obras de fábrica, azudes y otras estructuras. El modelo puede usarse para valorar las limitaciones de los cauces de avenida, identificar las zonas de riesgo en avenidas, el tratamiento de planas de inundación, valoración de las mejoras del cauce, y determinar los flujos divididos. 4.4.5 HSCTM2D Modela el flujo hidrodinámico, transporte de sedimento y del transporte de contaminantes. HSCTM2D trabaja con elemento finito que modela el sistema para el flujo del agua de dos dimensiones. El sistema que modela consiste en dos módulos, uno para modelar hidrodinámico (HYDRO2D) y el otro para modelar del transporte del sedimento y del contaminante (CS2D). El HSCTM2D, que modela el sistema, se puede utilizar para simular a corto plazo (menos de 1 año) y a largo plazo el problema de la sedimentación y transporte de sedimentos y el problema de transporte de contaminantes en aguas de superficie bien mezcladas verticalmente. Este programa se desarrolla en entorno DOS 16 – bit con un lenguaje de desarrollo FORTRAN. Este software es de alcance para todo el público en general de manera gratuita. 4.4.6 River2D Es un modelo hidrodinámico hecho para una profundidad promedio de dos dimensiones del elemento finito. Ofrece parámetros hidrodinámicos confiables (velocidades, niveles de agua, número de fraude, etc) y adicionalmente ofrece una evaluación del habitat de peces. La habitación modelo de River2D consiste en realmente cuatro programas: R2D_Bed, R2D_Ice, R2D_Mesh y River2D. Estos corren en entorno windows. R2D_Bed fue diseñado para corregir datos de la topografía del terreno. El programa de R2D_Mesh se utiliza para el desarrollo de los acoplamientos de cómputo que serán entrados en última instancia para River2D. Estos programas se utilizan típicamente en sucesión. El proceso de modelación normal implica crear un archivo preliminar de la topografía del terreno (texto) de los datos crudos del campo, en R2D_Bed. Después en el archivo R2D_Mesh se crea una malla la cual puede ser refinada las veces que se desee e incluso se pueden corregir datos de la topografía. También se utiliza para solucionar para las profundidades y las velocidades del agua a través de la discretización. Finalmente, River2D se utiliza para visualizar y para interpretar los resultados y para realizar análisis del flujo y el habitat de los pescados del tipo de PHABSIM. Este programa se obtiene en forma gratuita. 4.4.7 RIVERMorph Este programa fue diseñado para proveer de profesionales de la restauración de la corriente, uso que analiza datos geomorfológicos y facilita el gravamen, supervisando y diseñando con un interfaz fácil de utilizar. La restauración de la corriente y el diseño de canal natural tienen problemas hidráulicos y geomorfológicos complejos. RIVERMorph se puede utilizar para realizarse los cálculos necesarios para analizar diseño de corrientes. El software funciona en Windows y se puede obtener en forma gratuita

94

4.4.8 SED3D Es un modelo numérico de tres dimensiones que simula la hidrodinámica de los flujos de agua y el transporte del sedimento en ríos y especialmente en lagos y los estuarios. SED3D simula el flujo y el transporte del sedimento en los lagos, los estuarios, los puertos, y las aguas costeras. Este programa tiene un sistema que se puede utilizar para simular el flujo y el transporte del sedimento en varios cuerpos del agua bajo fuerzas de vientos, mareas, afluencias de agua dulce, y de gradientes de densidad con la influencia de la aceleración de Coriolis y de la geometría batimétrica. Este modelo generalmente es tridimensional pero puede funcionar bidimensionalmente tanto en la dirección x-y o en la dirección x-z. El límite apropiado y las condiciones iniciales correctas permiten que el modelo pueda computar los componentes de velocidad, temperatura, salinidad, y concentración suspendida del sedimento. Este software corre en el sistema operativo VAX VMS con el lenguaje de desarrollos FORTRAN y su distribución es libre y gratuita. 4.4.9 SMS 9 Este programa modela el agua superficial y tiene un ambiente comprensivo para la hidrodinámica, transporte del agente contaminador y del sedimento. SMS 9 incluye herramientas que modelan en dos dimensiones (2D) con método de elementos finitos, en dos dimensiones(2D) con método de diferencias finitas, en tres dimensiones con método de elementos finitos y en una dimensión (1D) para modelar el remanso. SMS también incluye un interfaz modelo genérico e incluye una variedad de información aplicable para modelar el agua superficial, elevaciones de la superficie del agua y velocidades del flujo. Simula para flujo permanente e impermanente. Los usos adicionales incluyen modelar la intrusión de la salinidad, la dispersión de la energía de la onda, características de la onda (direcciones, magnitudes y amplitudes) y de otras. Su precio es de US$ 9950. 4.4.10 TUFLOW Es un software de 1D y 2D de diferencias finitas que simulan inundaciones y mareas. Simula la hidrodinámica de los cuerpos del agua usando ecuaciones de conservación de masa y de la superficie libre. TUFLOW se desarrolla para determinar patrones del flujo en aguas costeras, estuarios y ríos. Sus fuerzas incluyen una excelente estabilidad y convergencia. Representan en 1D y 2D estructuras hidráulicas como diques y terraplenes. Simula en flujo en régimen crítico y supercrítico. Buena modelación en inundaciones importantes en crecidas de y en hidráulica de estuario y áreas costeras. TUFLOW utiliza GIS para manejar, para manipular y para presentar datos. La licencia cuesta US$ 5995 4.4.11 WMS (Watershed Modeling System) Este software tiene un ambiente que modela gráficamente todas las fases de la hidrología y de la hidráulica de la vertiente de aguas. WMS incluye las herramientas de gran alcance para automatizar modelar procesos tales como delineación automatizada del lavabo, cálculos geométricos del parámetro, cómputos de recubrimiento de GIS (profundidad de la precipitación, coeficientes de la aspereza, etc.), extracción seccionada transversalmente de datos del terreno, etc. El diseño modular del programa permite al usuario seleccionar los módulos en combinaciones de encargo, permitiendo que el usuario elija solamente las capacidades de

95

modelación hidrológicas que se requieran. Los módulos adicionales de WMS se pueden comprar y agregar en cualquier momento. El software se ligará dinámicamente a estos módulos subsecuentes en el tiempo de pasada que agrega automáticamente capacidad que modela adicional al software. Este software cuesta US$ 3700 con todos los módulos completos. En la tabla 4.1 se presenta un resumen de las características principales de los modelos explicados anteriormente. En esa tabla se ha tomado ciertos criterios para poder seleccionar el modelo más indicado. Después del análisis de todos los criterios, se llegó a la conclusión de que el River 2D es el más indicado para la modelación en esta tesis puesto modela para flujos permanentes e impermanentes, es de acceso gratuito (Internet), rápido aprendizaje en el manejo del software y brinda resultados muy confiables.

96

Tabla 4.1 : Características principales de softwares de simulación de flujo de ríos

Criterios Dim. Fenómenos que Tipo y Régimen de Flujo Métodos Entorno Costo modelan de en que Softwares Solución trabajan

AquaDyn 2D Flujo Hidrodinámico Permanente e impermanente MEF Windows US$ 1975 Factor de Riesgo Subcrítico y crítico Estudios de Impacto

AQUASEA 2D Flujo Hidrodinámico Permanente e impermanente MEF Windows US$ 1800 Transp..de sedimentos Subcrítico y crítico Transporte - dispersión

HEC - RAS 1D Flujo Hidrodinámico Permanente e impermanente MEF Windows Gratuito Subcrítico

HEC - 2 2D Flujo Hidrodinámico Permanente e impermanente MEF MS - DOS Gratuito Perfil de sup. del agua Subcrítico y crítico Espolón

HSCTM2D 2D Flujo Hidrodinámico Permanente e impermanente MEF DOS 16 Gratuito Transp. de sedimentos Subcrítico y crítico Transp. contaminantes

River 2D 2D Flujo Hidrodinámico Permanente e impermanente MEF Windows Gratuito Habitat de peces Subcrítico Evaluación del hielo

RIVERMorph 2D Flujo Hidrodinámico Permanente Windows US$ 1000 Perfil de sup. del agua Subcrítico Parámetros geomorf.

SED3D 3D Flujo Hidrodinámico Permanente e impermanente VAX VMS Transp. de sedimentos Subcrítico y crítico

SMS 9 3D,2D, Flujo Hidrodinámico Permanente e impermanente MDF Windows US$ 9950 1D Transp. de sedimentos Subcrítico y crítico Transp. contaminantes

TUFLOW 2D, 1D Flujo Hidrodinámico Permanente e impermanente MDF Windows US$ 5995 Inundaciones costeras Subcrítico y crítico Mareas

WMS 2D,1D Flujo Hidrodinámico Permanente e impermanente MDF Windows US$ 3700 Análisis de cuencas Subcrítico y crítico

97

4.5 River 2D 4.5.1 Generalidades El River2D es un modelo bidimensional, de profundidad promedio hidrodinámica y modelo de hábitat de peces. Se desarrolla específicamente para uso en arroyos naturales y ríos. Es un modelo de Elementos Finitos, basado en una formulación conservadora de Petrov-Galerkin upwinding. El modelo de River2D se desarrolló en la Universidad de Alberta con financiamiento proporcionado por: Concejo de Investigación de Canadá en las Ciencias Naturales e Ingeniería, el Departamento de Pesquerías y Los océanos, el Gobierno de Canadá, Protección Medioambiental de Alberta, y por el Estudio Geológico de los Estados Unidos. El código fue escrito por F. Hicks, A. Ghanem, J. Sandelin, P. Steffler, y J. Blackburn. T Waddle ha proporcionado una guía e informe en el desarrollo del modelo para las aplicaciones del mundo real. Continuado el desarrollo del modelo se ha facilitado y se ha apoyado por C. Katopodis. El derecho de propiedad literaria es retenido por estas personas y por la Universidad de Alberta. El modelo River2D, en la plataforma de un programa ejecutable como el Windows (95/98/2000/ME/NT/XP), está disponible en el dominio público. El programa es proporcionado como es, sin la garantía de integridad o aplicabilidad a cualquier problema en particular. El programa y las utilidades asociadas, los ejemplos de archivos de datos y la documentación pueden copiarse libremente y puede distribuirse con tal de que esta noticia sea incluida y usada en modelo reconocido propiamente. 4.5.2 Formulación El modelo River2D es bidimensional, como ya se sabe, es un modelo de elemento finito de profundidad promedio. Propone ser usado en arroyos naturales y ríos. Tiene características especiales para el ajuste del flujo en transiciones supercríticas o subcríticas, cubiertas de hielo, y en área húmeda variable. Es básicamente un modelo transitorio pero proporciona una acelerada convergencia para las condiciones de régimen permanente. El River2D usa las unidades del Sistema Internacional S.I ( Kg, m, s) para los datos de entrada y de salida. 4.5.2.1 Modelo hidrodinámico El componente hidrodinámico del modelo River2D es basado en dos dimensones, las ecuaciones de profundidad promedio de St. Venant se expresan en forma conservadora. Estas tres ecuaciones representan la conservación de masa de agua y de los dos componentes del vector de velocidad adquirida. Las variables dependientes realmente resueltas son la profundidad y las intensidades de la descarga en la dirección de las dos coordenadas respectivas. Conservación de masa:

0=∂

∂+

∂

∂+

∂

∂

y

y

x

x

t

qqH

4.10

98

Conservación de la velocidad adquirida en X

∂

∂+

∂

∂+−=

∂

∂+

∂

∂+

∂

∂+

∂

∂)(1)(1)(

2)()( 0

2 ττ ρρ xyy

xxx

fyxx

xy

xxt

x HHSSgHHgqVqU

q

4.11 Conservación de la velocidad adquirida en Y

∂

∂+

∂

∂+−=

∂

∂+

∂

∂+

∂

∂+

∂

∂)(1)(1)(

2)()( 0

2 ττ ρρ yyy

yxx

fyyy

yy

yxt

y HHSSgHHgqVqU

q

4.12 donde H es la profundidad de flujo, U y V son la velocidad de profundidad promedio en las direcciones X e Y respectivamente. El qx y los qy son las intensidades de la descarga respectivas a que se relacionan a través de los componentes de la velocidad. qx = HU

4.13 qy = HV

4.14 g es la aceleración debido a la gravedad y ρ es la densidad del agua. S0x y S0y son la pendiente del fondo en las direcciones X e Y; Sfx y Sfy son la correspondiente fricción slopes (pendiente de fricción), los τxx , τxy , τyx y τyy son las componentes horizontales de la fuerza turbulenta. 4.5.2.2 Suposiciones básicas • La distribución de presión en la vertical es hidrostática. Generalmente, esto limita la

exactitud en las áreas de pendientes empinadas y los cambios rápidos de pendiente en el fondo. Hablando bruscamente, características del fondo de tamaño horizontal menor que 10 veces la profundidad (típicamente la duna) no se modelara con precisión. Similarmente, pendientes en la dirección de flujo en exceso del 10% no será modelado correctamente.

• Las distribuciones horizontales de la velocidad por encima de la profundidad son

esencialmente constantes. Una distribución de velocidad supuesta puede ser usada en la interpretación de proporcionar la velocidad de profundidad promedio, pero la distribución es tratada como constante por los cálculos internos.

• Las fuerzas de Coriolis y de viento son asumidas despreciables. Para cuerpos muy

grandes de agua, particularmente para grandes lagos y estuarios, estas fuerzas pueden ser significativas.

99

4.5.3 Métodos de modelo hidrodinámico

4.5.3.1 Método de Elementos Finitos

El Método de Elementos Finitos usado en el modelo hidrodinámico del River2D es basado en la formulación residual ponderada del Perfil aerodinámico Upwind Petrov-Galerkin. En esta técnica, río arriba influye las funciones de prueba a ser usadas para asegurar la estabilidad de la solución bajo el pleno alcance de las condiciones del flujo, incluyendo los flujos subcrítico, el supercrítico, y el transcrítico. Una discretización totalmente conservadora es llevada a cabo para asegurar que ninguna masa flujo está perdida o ganada por encima del dominio modelado. Esto también permite aplicación de condiciones de borde como el flujo natural o condiciones forzadas. Una elevación fija se impone como un equivalente de la fuerza hidrostática, por ejemplo y las descargas de entrada conocidas son usadas directamente. Las discretizaciones ejemplares apoyadas incluyen igual orden lineal, interpolaciones cuadráticas, y cúbicas de todas las variables encima de los elementos triangulares y cuadriláteros. El generador de la malla sólo apoya los triángulos lineales, sin embargo, estos elementos son los más simples posible en dos dimensiones y produce el tiempo de la ejecución mínima para un número dado de nodos. Aunque el modelo realmente es inestable, puede usarse para realizar un análisis transitorio y obtener una solución de estado permanente. Para los resultados del estado permanente, una convergencia acelerada en el procedimiento acelera el proceso a la final realización aumentando sistemáticamente el incremento de tiempo. 4.5.3.2 El Método de Newton Raphson Usando el Método de Elementos Finitos para resolver las ecuaciones hidrodinámicas resultadas de un sistema de ecuaciones no simétrico, las ecuaciones no lineales pueden ser resueltas por métodos explícitos o implícitos. En River2D, un método implícito es tomado para requerir una solución simultánea de un sistema de ecuaciones. Porque el sistema es no lineal, el método de Newton-Raphson reiterativo es empleado. En este método, la solución reiterativa es lograda usando la solución de los valores desconocidos previos al paso de tiempo. El Φn, se usa como la primera suposición de la solución al nuevo paso de tiempo, Φn+1. A menos que las suposiciones son correctas, entonces un residuo resultará, tal que:

{ }Rf mnmn ,1,1 ++=

φ

4.15 donde n y m denotan el paso de tiempo y la iteración respectivamente. Las correcciones ∆Φ, para obtener las próximas suposiciones en el proceso reiterativo son determinadas para resolver el sistema siguiente de ecuaciones. [ ]{ } { }RJ mnmnmn ,11,1,1 ++++ =∆φ 4.16

100

donde [Jn+1,m] se llama la matriz Jacobiana y se define como:

[ ] { }{ }

∂∂

= +

++

φ mn

mnmn RJ ,1

,1,1

4.17 Después de cada m-ésima iteración los valores de las variables de la solución son actualizadas usando:

φφφ 1,1,11,1 +++++ ∆+= mnmnmn 4.18 El Jacobiano puede evaluarse analíticamente, usando la ecuación 4.17, o numéricamente, usando la siguiente aproximación:

{ }{ }

−−

≈

∂∂

++

+

+

δδφφ

φ j

j

mn

jimn

ji

mn

mn RRR,1,1

,1

,1

4.19 donde δ es una perturbación pequeña y los índices i y j representan el número de la ecuación y el número variable respectivamente. En River2D el usuario tiene la opción de evaluar la matriz de Jacobiano analíticamente o numéricamente. El Jacobiano analítico es más rápido al evaluar que el numérico, pero debido a la suposición en la derivación matemática, puede ser menos exacto que el Jacobiano numérico. En particular, el Jacobiano analítico no considera los cambios en la matriz del upwinding o el efecto del cambio de ubicación del borde del agua. Al realizar los análisis del estado permanente, es mejor empezar el proceso de la solución usando el Jacobiano analítico. Si la solución no puede lograr el nivel deseado de convergencia, es aconsejable a cambiar al Jacobiano numérico para el resto de la simulación. El Jacobiano analítico debe ser adecuado para las simulaciones transitorias, como el borde del agua y la matriz del upwinding es asumida para ser fijas por encima del intervalo de tiempo. 4.5.4 Método de solución de ecuaciones

Dentro de cada iteración Newton-Raphson, el sistema de ecuaciones lineales descritas por ecuación 4.16 requiere solución. En River2D, el usuario tiene la opción de escoger dos soluciones de ecuaciones lineales: una solución directa y una iterativa. 4.5.4.1 Solución directa

La primera solución es llamada la solución de ecuación de la zona activo (Stasa 1985), y emplea una solución directa al sistema lineal no simétrico. Es decir, una solución basada en

101

la eliminación de Gauss. En la solución de ecuación de zona activa, la matriz del coeficiente [A] en el sistema lineal es implicada por: [ ]{ } { }bxA = 4.20 es descompuesto matrices triangulares inferiores y superiores. [L] y [U] tal que: [ ] [ ][ ]ULA = 4.21 Por consiguiente, usando la ecuación 4.21, la ecuación 4.20 puede volverse a escribir como [ ][ ]{ } { }bxUL = 4.22 qué es, en el efecto, dos sistemas de ecuaciones: [ ]{ } { }bcL = 4.23 [ ]{ } { }cxU = 4.24 El sistema definido por ecuación 4.20 es luego resuelta por substitución delantera en la ecuación 4.23 para calcular {c} y luego por substitución hacia atrás en la ecuación 4.24 para obtener la solución del vector {x}. En River2D, la matriz del coeficiente es una Matriz Jacobiana muy grande con la dimensión 3N por 3N, dónde N es el número de los nodos computacionales. La Matriz Jacobiana es típicamente esparcida (muchas de sus entradas son 0) y bandeada (la mayoría de ceros se encuentran a lo largo de la diagonal principal). La solución de la ecuación de la zona activa toma ventaja de estas propiedades para reducir los requisitos de almacenamiento usando un método de almacenamiento llamado skyline (Stasa 1985). Esto reduce los requisitos de memoria significativamente de 3N por 3B, dónde B es el ancho de banda de la matriz. Sin embargo, a menos que el ancho de banda es pequeño, los requisitos del almacenamiento pueden ser bastante significantes. El ancho de banda de la matriz es función de la geometría y de los números de nodos de la malla. Los cauces estrechos producirán pequeños anchos de banda mientras los cauces anchos producirán grandes anchos de banda. Desafortunadamente, la velocidad de la solución de la ecuación de la zona activa es una función de los requerimientos de memoria para la Matriz Jacobiana.

102

4.5.4.2 Solución iterativa La segunda solución usa un método reiterativo para resolver sistemas lineales no simétricos llamado Método Generalized Minimal Residual o Método de GMRES (Saad y Schultz 986). Esta solución estaba incorporada en River2D porque no tiene grandes condiciones de almacenamiento como la solución directa. Esto permite por consiguiente que para mallas más grandes puedan ser resueltas usando la solución directa. En el método de GMRES, una suposición inicial de la solución al vector {x0} en la ecuación 4.20 es escogido. Esto resulta en: { } [ ]{ } { }bxAr −= 00

4.25 donde {r0} se llama el vector residual. La formación de un residuo ilustra por qué este método, como todo los métodos reiterativos, requiere de un almacenamiento muy pequeño. No hay ninguna necesidad de guardar la matriz del coeficiente [A] explícitamente, como es necesario en los métodos directos; sólo necesita que se guarde implícitamente en el vector del producto [A]{x}. En River2D, este producto vectorial se obtiene sumando la matriz del producto vectorial al nivel del elemento. El vector residual después es usado para formar una base ortogonal:

{ }vvvV kk ,...,, 21= 4.26

que es llamado sub espacio de Krylov, Kk . De aquí en adelante es llamado un conjunto de vectores Krylov Vk. La solución de GMRES del sistema lineal no simétrico, descrito en ecuación 4.20, es formada usando una combinación lineal de estos vectores Krylov. Eso es: { } { } { } { } { }vyvyvyxx kk...22110 ++=

4.27 donde los yi son los coeficientes de combinación lineal. Una vez que el residuo se minimiza a una tolerancia aceptable, entonces se dice que la solución ha convergido y ningún vector más se requiere. El número de vectores, o el número de iteraciones k, que son requeridas depende en cómo la tolerancia es estricta. Sin embargo, el número de iteraciones no puede exceder el número de incógnitas, N. Esto es debido al hecho que el sub espacio Krylov se define por los vectores de N. Desafortunadamente, el tiempo de ejecución de este algoritmo aumenta geométricamente con cada la iteración adicional. Esto es porque el cálculo del vector vi es basado inclusive en los vectores v1 a vi-1. El algoritmo es reiniciado tantas veces como es exigido para lograr un nivel deseado de convergencia. La versión reiniciada de GMRES, denotado por GMRES(m) es el método reiterativo usado en River2D.

103

En River2D, GMRES (m) es llevado a cabo para que continúe la reiniciación hasta que el residuo se minimice dentro de una tolerancia o hasta que un número máximo de reinicios, o iteraciones, haya ocurrido. Actualmente, el goal tolerance (meta de tolerancia) es comparada a la magnitud del residual al final de cada iteración dividida por la magnitud del vector {b}. Esto es:

{ }{ }

∑

∑

=

==N

ii

N

ii

m

b

r

br

1

2

1

2

4.28

El número de pasos hasta el reinicio, m, y el número máximo de iteraciones, k, también son parámetros definidos por el usuario. El valor predeterminado para estos parámetros es de 10 reinicios y 10 iteraciones. El valor predeterminado para el goal tolerance es 0.01. Estos valores fueron escogidos basado en un pequeño número de pruebas.