Universidad Nacional Agraria La Molina Facultad de Ingeniería Agrícola Programa Maestría en Ingeniería de Recursos Hídricos

Una Introducción al Método de Diferencias Finitas en la Modelación de Aguas Subterráneas

Eusebio Ingol Blanco, Ph.D

From Daene McKinney, Ph.D, UT Texas

Procesos de Modelación

• • • • Identificación del Problema – Sistema y elementos a ser modelados – – Relación e interacción entre ellos Grado de precisión Conceptualización y Desarollo – Representación matemática – – – – – Modelo tipo Método numérico Código de computo Condiciones iniciales Condiciones de frontera Datos de entrada Calibración – Estimación de Parámetros – – – Comparación de modelados con simulados Ajuste de parámetros Análisis de incertidumbre Validación – Usar un grupo de datos independiente al de calibración – Comparar resultados del modelo con observados Identificación del problema y descripción conceptualización desarrollo Calibración Validación y análisis de Sensibilidad Aplicación Presentación resultados

Modelo Conceptual y Matemático

•

Modelo Conceptual

– Idealización y implicación de las condiciones hidrológicas – Representación descriptiva del sistema de agua subterránea que incorpora la interpretación de las condiciones geológicas y hidrológicas – Que procesos son importantes para el modelo – – Cuales son las fronteras Que datos necesitan ser colectados y cuales están disponibles.

Modelo Conceptual y Matematico

•

Modelo Matemático

– Representación matemática de la hidrología subterránea y del transporte de contaminantes – Utiliza las ecuaciones fundamentales del flujo y conservación de masa – Se basa en observaciones reales o extrapolaciones

Métodos

• •

Analítico

– Uso de formulas simples – Solución exacta en el el punto de Calculo – – Modelo homogéneo Se basa en observaciones reales o extrapolaciones

Numérico

– Discretizacion espacial y temporal – Solución aproximada – – Admite heterogeneidades Transforman las Ecuaciones en Derivadas Parciales (PDFs) que gobiernan el flujo en Ecuaciones Diferenciales Ordinarias (ODEs) o ecuaciones algebraicas para su solución

Diferencias Finitas y Elementos Finitos

Que realmente queremos resolver?

• Flujo horizontal en un acuífero confinado ¶ ¶

x

ç

T x

¶

h

¶

x

ö ø ¶

y

æ è

T y

¶

h

¶

y

ö ø

K

' (

h

0

b

' -

h

) ±

W

å

w

= 1

Q w

d (

x

-

x w

) =

S

¶ ¶

t h

flujo filtración fuente/pozo almacenamiento • • • superficie Ecuaciones que gobiernan Carga en acuifero confinado Condiciones iniciales Capa confinada Condiciones de frontera

Q x z y x

confinado

K

basamento

b h

Método de Diferencias Finitas

• Finite-difference method – Replace derivatives in governing equations with Taylor series approximations – Generates set of algebraic equations to solve

Taylor Series

• Expresion en series de Taylor de h(x) en un punto t x+ D x cerca a x • Si se trunca la serie despues del n th termino, el error sera

Primera Derivada- Hacia Adelante

• Considerar una expansión en series de Taylor hacia adelante de una función h(x) cerca al punto x • Resolver para 1 st derivada

h

(

x

)

x

 D

x

D

x x h

(

x

) D

x x

 D

x x

Primera Derivada- Hacia Atras

• Considerar una expansión en series de Taylor hacia atras de una función h(x) cerca al punto x • Resolver para 1 st derivada

h

(

x

)

x

 D

x

D

x x h

(

x

) D

x x

 D

x x

Segunda Derivada- Central

Adicionar y resolver para

h

¢ ¢ (

x

)

Aproximación en Diferencias Finitas

h

(

x

)

h i

1

i

1 D

x i h i

D

x h i

+ 1

i

+ 1

x

Hacia atras 1 st derivada Hacia adelante 1 st derivada Central 2 nd derivada

Grids y Discretización

• • • • Proceso de discretizacion Grid definido para cubrir el dominio

y, j

Objetivo es predecir los valores de carga en los puntos de nodo de la malla – Determinar efectos de bombeo – Flujo de un rio, etc Método D.F

– Popular y fácil de implementar – Atractivo para simple geometría D

y i,j

+1

i

-1

,j i,j i+

1

,j i,j

-1 D

x

Notation

h

(

x

,

y

,

z

,

t

) =

h l i

,

j

,

k

Grid cell Mesh Domain Node point

x, i

Grids en Tres Dimensiones

• • • Un sistema acuífero esta dividido en bloque rectangulares por un grid. El grid es organizado por filas (i), columnas (j), y capas (k), y cada bloque es llamado celda.

Tipos de capas – Confinada – – No confinada Convertible

Capas pueden ser de Diferente material

Flujo en Acuífero Confinado 1-D

• • • Medio Homogéneo, isotrópico, 1-D, flujo confinado Ecuación principal ¶ ¶

x

æ è

T

¶

h

¶

x

ö ø

S

¶ ¶

t h

Condición inicial

h

(

x

,0) = 6.1

m

• Condiciones de frontera

h

(0,

t

) = 6.1

m h

(

L

,

t

) = 1.5

m z

Nodo

y x h A i =

0 Superficie del suelo 1 Acuifero 2 Confinada capa 3 4 5 Celda grid 6 7 D

x

8 D

x

= 1

m

,

L

= 10

m

,

b

= 1.5

m h A

= 6.1

m

,

h B K

= 0.5

m

/

d

,

S

= 1.5

m

, = 0.02

9 10

b h B

Aproximación de Derivadas

• Ecuación que rige ¶ 2

h

¶

x

2 =

S

¶

h T

¶

t i

,

l

 1

t

,

l

D

t i

 1 ,

l

• 2nd derivada x

i

,

l i

 1 ,

l

D

x i

,

l

 1 • 1st derivada t

x

,

i

Hacia adelante Hacia atras

Método Explicito

• • • Use toda la información del paso de tiempo anterior para calcular el valor en este paso de tiempo Procede punto por punto a través del dominio Podría ser inestable para largos periodos de tiempo

t

,

l i

 1 ,

l i

,

l

 1

i

,

l

D

t i

 1 ,

l

D

x i

,

l

 1

x

,

i

¶ 2

h

¶

x

2 =

S

¶

h T

¶

t

Aprox. DF

Método Explicito

h l i

+ 1 =

h l i

+ (

l i-

1 2

h l i

+

h l i

+ 1 )

l

+1 nivel tiempo desconocido

l

nivel tiempo conocido

Método Explicito

Superficie del terreno

h l i

+ 1 =

h l i

+

r

(

h l i-

1 2

h l i

+

h l i

+ 1 ) Nodo

h A

Acuifero Confinada layer D

x i =

0 1 2 3 4 5 Celda grid 6 7 8 9 10

b h B

Considerar: r = 0.48

r

= 0.52

D

x

= 1

m

,

L

= 10

m

,

b

= 1.5

m h A

= 6.1

m

,

h B K

= 0.5

m

/

d

,

S

= 1.5

m

, = 0.02

Resultados (

D

t = 18.5 min; r = 0.48 < 0.5)

Resultados (

D

t = 20 min; r = 0.52 > 0.5)

Que pasa aqui?

• • • En el time t = 0  no flujo En el time t > 0  flujo El agua proveniente del almacenamiento en una celda grid sobre el tiempo Dt • El agua fluyendo fuera de la celda en el intervalo D

t h A

Superficie del suelo Acuifero Confinada capa D

x

D

x i =

0 1 2 Celda Grid

i

…

i

-1

i i

+1 … 8 9 10

b h B

r > 0.5

El intervalo de tiempo es demasiado largo Las celdas no contienen la suficiente agua Causa inestabilidad

Método Implícito

• • • Usa información de un punto en el paso de tiempo anterior para calcular el valor en todos los puntos de este paso de tiempo Resuelve para todos los punto en el dominio simultáneamente Es mas estable ¶ 2

h

¶

x

2 =

S

¶

h T

¶

t

Aprox. DF

i

 1 ,

l

 1

i

,

l

 1

t

,

l

D

t i

 1 ,

l i

,

l

D

x i

 1 ,

l

 1

i

,

l

 1

i

 1 ,

l

 1

i

 1 ,

l i

 1 ,

l

 1

x

,

i

Método Implícito

-

rh l

+ 1

i-

1 + (1 + 2

r

)

h l i

+ 1 -

rh l i

+ 1 1 =

h l i l

+1 nivel tiempo desconocido

l

nivel tiempo conocido

Flujo en Estado Estable 2-D

• • Ecuación que rige

y

,

j

¶ ¶

x

ç

T x

¶

h

¶

x

÷ + ¶ ¶

y

æ

T y

¶

h

¶

y

ö ø

W

å

w

= 1

Q w

=

S

¶ ¶

t h

Homogéneos, acuífero isotrópico, no pozos D

y

¶ 2

h

¶

x

2 + ¶ 2

y h

2 = 0 Nodo No.

( 1 , 5 ) ( 1 , 5 ) Cargas desconocidas ( 2 , 5 ) ( 5 , 4 ) Cargas conocidas ( 3 , 5 ) ( 4 , 5 ) ( 0 , 4 ) ( 1 , 4 ) ( 2 , 4 ) ( 3 , 4 ) ( 4 , 4 ) ( 5 , 4 ) ( 0 , 3 ) ( 1 , 3 ) ( 2 , 3 ) ( 3 , 3 ) ( 4 , 3 ) ( 5 , 3 ) ( 0 , 2 ) ( 1 , 2 ) ( 2 , 2 ) ( 3 , 2 ) ( 4 , 2 ) ( 5 , 2 ) ( 0 , 1 ) ( 1 , 1 ) ( 2 , 1 ) ( 3 , 1 ) ( 4 , 1 ) ( 5 , 1 ) ( 1 , 0 ) ( 2 , 0 ) ( 3 , 0 ) ( 4 , 0 ) D

x x

,

i

• Igual espaciamiento (promedio de celdas)

h i,j

=

h i-

1

,j

+

h i+

1

,j

+

h i,j-

1 +

h i,j+

1 4

Flujo Anisotropico Heterogéneo 2-D

y

node (

i,j

)

i-1/2 i+1/2

¶ ¶

x

ç

T x

¶

h

¶

x

÷ + ¶ ¶

y

æ

T y

¶

h

ö ¶

y

0

j+1

D

x

cell (

i,j

)

Q y,j+1/2 j+1/2 T x

and T and y

y

son transmisividades en x

j Q x,i-1/2 Q x,i+1/2

D

y j-1/2 j-1 i-1 i Q y,j-1/2 i+1 x

Flujo Anisotropico Heterogéneo 2-D

• Promedio armonico para T o K

T x i+

1

/

2

,j

=

T x i+

1

,j T i,j x

2

T x i+

1

,j +T i,j x A i,j h i

+ 1

,j

+

B i,j h i

1

,j

+

C i,j h i,j

+ 1 +

D i,j h i,j

1 +

E i,j h i,j

= 0

Problema Transitorio

¶ ¶

x

ç

T x

¶

h

¶

x

ö ø ¶ ¶

y

æ

T y

¶

h

ö ¶

y S

¶ ¶

t h A i,j h l

+ 1

i+

1

,j

+

B i,j h l

+ 1

i-

1

,j

+

C i,j h l

+ 1

i,j+

1 +

D i,j h l

+ 1

i,j-

1 + ¢

h l

+ 1

i,j

=

F i,j h l i,j

MODFLOW

• Modelo Matematico de la USGS • • • • http://water.usgs.gov/nrp/gwsoftware/modflow.html

USGS desarrollo el modelo matematico Usa el elemento de Diferencias Finitas Varias versiones – MODFLOW 88, 96, 2000, 2005 Interfaces graficas para MODFLOW – USGS interface, MPI – GWV ( www.groundwater-vistas.com

) – – GMS ( www.ems-i.com

) PMWIN ( www.ifu.ethz.ch/publications/software/pmwin/index_EN )

Que Simula MODFLOW?

1.

2.

3.

4.

5.

6.

7.

8.

9.

10.

11.

12.

Unconfined and confined aquifers Faults and other barriers Fine-grained confining units and interbeds Confining unit - Ground-water flow and storage changes River – aquifer water exchange Discharge of water from drains and springs Ephemeral stream - aquifer water exchange Reservoir - aquifer water exchange Recharge from precipitation and irrigation Evapotranspiration Withdrawal or recharge wells Seawater intrusion