Donde C () es la capacidad capilar específica de la humedad del suelo (cm³∙cm^-2)  , con C () ≠ 0 porque el medio es no saturado y además representa la curva de succión-humedad:  dado que q es el contenido de agua en e suelo (cm³∙cm^-2), t designa el tiempo (d),      K () es la conductividad hidráulica del suelo no saturado (cm/d) y z la dimensión vertical del suelo (cm), se considera creciente en el sentido hacia arriba.

Para completar el modelo representado por la ecuación de Richards se emplean relaciones constitutivas que facilitan la obtención de , C y K, según van Genuchten (1980).

La relación de retención del agua está dada por donde , describe el contenido del agua en el suelo residual (cm³∙cm^-2), representa la porosidad del suelo (cm³∙cm^-2), y n son parámetros que miden la textura del suelo.

La función de conductividad hidráulica viene dada por , K_s representa la conductividad hidráulica saturada.
Para la implementación matemática del modelo se desarrolla un esquema en Diferencias Finitas por el método Crank-Nicholson [CN]. La discretización en profundidad y en tiempo es lo suficientemente pequeña como para garantizar la eficiencia en las ejecuciones del modelo
 Métodos de asimilación
En el estudio se consideran dos métodos de asimilación de los valores de la humedad superficial en el modelo descrito. El primero de ellos Inserción Directa, y el segundo Filtro de Kalman. A continuación se describen algunos elementos básicos para entender el principio en el que se sustentan. ]]>
Asimilación por Inserción Directa (ID)

En la técnica ID se asume como verdadero el valor de la observación superficial, en detrimento de la estimación del modelo para las profundidades correspondientes. En consecuencia para cada momento de asimilación se sustituye directamente el valor del modelo por el valor de la observación. Esta técnica suele ser ventajosa por su simplicidad, pero requiere que el error de las observaciones sea pequeño. En el caso concreto de esta investigación dicha inserción consiste en actualizar los valores de humedad estimados de los nodos superficiales, con los valores de humedad “observados”.

Asimilación mediante Filtro de Kalman (FK)

En el FK el estado desconocido x_k, con densidad de probabilidad p(x_k), evoluciona en el tiempo k de acuerdo a una probabilidad p(x_k|x_k-1). Las observaciones y_k son independientes dado el estado y pueden ser generadas de acuerdo a la densidad p(y_k|x_k). El sistema dinámico estado-espacio es escrito según:

donde v_k es el ruido del proceso que interviene a través de la función de transición F, y n_k es el ruido de la observación que daña la medición del estado a través de la función de observación H. La entrada exógena del sistema, u_k, es asumida conocida. En principio F y H dependen de los parámetros del modelo w.
Basado en estas consideraciones Kalman (1960) derivó el siguiente esquema recursivo de la actualización óptima de la media y la covarianza del estado:

]]>
 

Los términos en esta recursión están dados por:

donde la predicción de x_k es y corresponde al valor esperado en el tiempo k-1 de una función no lineal de variables aleatorias x_k-1 y (interpretación similar para la predicción de excepto que el valor esperado es tomado en el tiempo k). La ganancia óptima K_k es expresada como una función del valor esperado de la matriz de covarianza del error de predicción del estado y el error de predicción de la observación, con el valor esperado de la matriz de autocorrelación del error de predicción de la observación.

Estudio sintético

El trabajo se basa en la generación sintética de valores de humedad (o potencial matricial) mediante un modelo de simulación basado en la ecuación de Richards. En él se adopta como perfil “real” (PR) aquel que resulta de una realización del modelo en el que se toma una condición inicial arbitraria “real”. Una segunda realización, también arbitraria, que considera una humedad inicial “errónea”, bastante alejada de la real, se considera como perfil “erróneo” (PE). Los valores de humedad superficial correspondientes a los perfiles reales más un error aleatorio de los mismos con un ruido blanco fueron considerados como observaciones del modelo.
Mientras, la evolución de los perfiles de humedad que resultan de considerar una condición inicial errónea, pero sujetos al proceso de asimilación, se identifican como perfiles “estimados”. Los parámetros del modelo y las condiciones de frontera en cada una de estas simulaciones son los mismos. El trabajo va dirigido a ilustrar en qué medida el perfil estimado mediante el proceso de asimilación de las observaciones evoluciona en el tiempo hacia el perfil real. Para comparar los resultados de esta investigación se han asumido algunas de las condiciones establecidas por Entekhabi, et al. (1994) y Walker (2001). En las simulaciones, se definen dos condiciones iniciales del perfil de la variable de estado potencial, el real y el erróneo. Sus valores son -50 cm y -300 cm uniformes a lo largo del perfil del suelo, con una profundidad 100 cm respectivamente. Los parámetros del suelo asumidos son:
K_s=25 cm/d, =0,54, =0,2, =0,008. Las condiciones de fronteras que se asumen son las de tipo Neumann dado que se ha considerado un flujo constante de 5,0112∙10-1 cm/d en la parte superior (superficie de suelo) y flujo cero en el inferior del perfil de suelo.

 

]]>  

RESULTADOS Y DISCUSIÓN

La Figura 1 muestra la evolución del perfil estimado mediante el uso de las técnicas de asimilación después de 12 horas, un día y tres días de iniciadas las simulaciones. Para propósitos comparativos se incluye el perfil real y el perfil erróneo.
Mediante el empleo de la asimilación directa el perfil estimado se aproxima de manera progresiva, pero muy lentamente, al perfil real. Aun a los tres días hay un considerable error, cuyo promedio es 7,54∙10-2 cm3∙cm-3 , de modo que no es sino hasta los 20 días que se logra obtener una “correspondencia total” entre el perfil estimado y el real. El método de ID es poco eficiente en aminorar el error de la predicción de la variable de estado en capas inferiores con respecto a la superficie del suelo.
En tanto, el proceso de asimilación usando el FK supera con creces el obtenido mediante el de ID; al cabo de 14 horas el perfil real y el estimado prácticamente coinciden, no siendo así para el caso del uso ID. Para este tiempo los errores absolutos de las profundidades del perfil estimado mediante FK está comprendido entre los valores 2,04∙10-6 cm3∙cm-3 y 3,33∙10-3 cm3∙cm-3 La gran ventaja del FK es que, aun cuando la asimilación se asocia a valores de humedad superficial, este método es capaz de afectar favorablemente el valor de la variable de estado a prácticamente todas las profundidades. En la medida que mayor sea el valor de la covarianza entre los valores de humedad a las menores y mayores profundidades, más efectivo tiende a ser el proceso de asimilación. Cuando se considera ID o el perfil que no considera asimilación (PE) sus errores absolutos correspondientes a las profundidades se hallan en los intervalos y  y,  respectivamente al cabo de 14 horas.

En la Tabla 1 se muestra las cotas inferiores y superiores de los errores relativos de los puntos que conforman los perfiles estimados mediante las técnicas empleadas en los tiempos fijados en la Figura 1. Como se puede apreciar el error relativo mínimo usando el FK es entre cuatro y diez órdenes inferiores al correspondiente usando ID; y este a su vez tiene un orden más pequeño que el obtenido usando el PE. Las diferencias entre los errores relativos máximos considerando ID y PE son prácticamente similares, excepto en el tercer tiempo escogido donde el error relativo máximo considerando ID es casi la mitad con respecto a PE; mientras el KF supera por mucho a ambos con errores próximos a cero, con órdenes que oscilan entre -6 y -2.

La Figura 2 pone de manifiesto la inmediata respuesta del FK en términos de reducción del error. Este método es sin dudas sumamente eficiente en la consecución del objetivo propuesto con el empleo de la asimilación. La curva de variación temporal del error absoluto usando ID se aproxima a una distribución exponencial con parámetro aproximadamente igual a 0,16, es decir, alrededor de seis días y seis horas el error es 2,7183 (número de Euler) veces menor. Aunque la técnica ID es muy inferior al FK, conlleva a una mejoría apreciable en comparación con el PE, cuya distribución de errores no es monótona decreciente.

]]>

Los resultados aquí obtenidos coinciden con los obtenidos para el caso de la estimación del perfil de potencial en Hernández y Medina (2011) y con otras literaturas donde consideran otras condiciones (Medina et al., 2008; Walker, 1999; Walter, et al., 2001).
La diferencia principal entre ambas técnicas de asimilación está dada por la naturaleza del FK. En esencia el algoritmo de FK pronostica el nuevo estado a partir de su estimación previa añadiendo un término de corrección proporcional al error de predicción, de forma tal que este último es minimizado estadísticamente, o sea las estimaciones son estadísticamente óptimas en el sentido que ellos minimizan el error que proviene de la estimación. En el caso del algoritmo ID se desprecia dicho error, o sea solamente se tiene en cuenta el valor de la estimación sin realizar un análisis del error que está asociado. Sin el uso de técnicas de asimilación no se tiene en cuenta actualización de las variables de estado, variables que por su naturaleza varían constantemente.
Una ventaja adicional del FK con relación al ID es la estabilidad del algoritmo de asimilación. En el estudio se pudo observar que cuando al humedad inicial estimada es muy diferente de la real, el proceso de asimilación usando ID es susceptible a colapsar. Bajo esa condición inicial la variabilidad de la variable de estado en profundidad es sumamente grande, y a su vez poco probable desde el punto de vista físico. Ello pude conllevar a una situación de intractabilidad para el proceso de asimilación usando ID; en cambio se ha apreciado que aun bajo esta condición el FK es capaz de corregir eficazmente el perfil estimado.
Por último es necesario hacer énfasis en aspectos vinculados al propio modelo de simulación. En el estudio se asume que el modelo de simulación es “perfecto”, tanto desde el punto de vista de su estructura, como de los parámetros que se utilizan. Sin embargo se sabe que, por bueno que sea un modelo de simulación, está siempre sujeto a errores que son regularmente muy difíciles de cuantificar. Un aspecto favorable del FK en el modo empleado es que, aun cuando no permite hacer frente directamente a esta problemática, posibilita optimizar la estimación incorporando el error del modelo, en el error global estimado. Ante esta problemática el ID se ha evidenciado que es un método mucho más frágil.

 

 

CONCLUSIONES

]]> En este trabajo se ha ilustrado las potencialidades de las técnicas de asimilación de datos en la estimación de perfiles de humedad del suelo a partir de la asimilación de valores de humedad superficial en un modelo de simulación. Las técnicas empleadas son ID y FK. El modelo de simulación que se ha asumido es la ecuación de Richards discretizada mediante el esquema en Diferencia Finita Crank-Nicholson y datos sintéticos.
Del estudio realizado resulta que si se emplean técnicas de asimilación se puede apreciar una disminución del error de perfiles estimados con respecto a los reales en el transcurso del tiempo, y con ello evita la posible divergencia de la solución como se muestra en el perfil estimado resultante del modelo. En particular haciendo uso de la técnica de asimilación FK obtiene resultados superiores en menor tiempo que si se considera la técnica de asimilación ID. Luego se recomienda el uso de técnicas de asimilación de datos y específicamente dentro de las dos técnicas aquí vistas, el FK.

 

 

REFERENCIAS BIBLIOGRÁFICAS

1. BACA, R. G., J. N CHUNG & D. J. MULLA: “Mixed transform finite element method for solving the non-linear equation for flow in variably saturated porous media”, International Journal for Numerical Methods en fluids, 24: 441-455, 1997.

2. CELIA, M. A. & E. T. BOULOUTAS: “A general Mass-Conservative Numerical Solution for the Unsaturated Flow Equation”, Water Resources Research, 26(7): 1483−1496, 1990.

]]>

3. ENTEKHABI, D., H. NAKAMURA & E.G. NJOKU: “Solving the inverse problem for soil moisture and temperature profiles by sequential assimilation of multifrecuency remote sensed observations”, IEEE Trans Geosci. Remote Sensing, 32(2): 438-448, 1994.

4. GALANTOWICZ, J. F., D. ENTEKHABI & E. G. NJOKU: “Test of sequential data assimilation for retrieving profile soil moisture and temperature from observed L-band radiobrightness”. IEEE Trans. Geosci. Remote Sensing 37(4): 1860–1870, 1999.

5. GIACOMELLI, A., U. BACCHIEGA, P.A. TROCH & M. MANCINI: “Evaluation of surface soil moisture distribution by means of SAR remote sensing techniques and conceptual hydrological modelling”, J. Hydrol. 166: 445-459, 1995.

6. HANSON, J.D., L.R. AHUJA, M.D. SHAFFER, K.W. ROJAS, D.G. DECOURSEY, H. FARAHANI & K. JOHNSON: “RZWQM: simulating the effects of management on water quality and crop production”. Agric. Systems 57: 161–195, 1998.

7. HANSON, J.D., K.W. ROJAS & M. J. SCHAFFER: “Calibrating the root zone water quality model”. Agronomy Journal, 91: 171–177, 1999.

8. HERNÁNDEZ PEREIRA, Y. y H. MEDINA: “Evaluación de técnicas de asimilación de datos para estimar la humedad del suelo”, Reporte de Investigación del Instituto de Cibernética, Matemática y Física, No 635, ISSN 0138-8916, La Habana, Cuba, 2011.

]]>

9. HILLEL, D.: “Fundamentals of Soil Physics” [s.n.]. Academic, San Diego, Calif., USA, 1980.

10. JACKSON, T.J., T.J. SCHMUGGE & E.T. ENGMAN: “Remote sensing applications to hydrology: soil moisture”, J. Hydrol., 41(4): 517-530, 1996.

11. JAZWINSKI, A. H.: “Stochastic processes and filtering theory, Mathematics in Science and Engineering”, edited by Richard Bellman, University of Southern California, Vol. 64, USA, 1970.

12. KALMAN, R. E.: “A New Approach to Linear Filtering and Prediction Problems”, Trans. ASME, Ser. D: J. Basic Eng., 82: 35-45 A960, 1960.

13. MA, L., H.D. SCOTT, M.D. SHAFFER & L.R. AHUJA: “RZWQM simulations of water and nitrate movement in a tall fescue field”. Soil Sci., 163: 259–270, 1998.

]]>

14. MEDINA, H., G.B. CHIRICO & N. ROMANO: “A new assimilation scheme of near-surface soil moisture to identify effective soil hydraulic parameters”, EGU2008-A-07120, EGU Assembly 2008, Vienna, 13-18 April 2008.

15.MERWE, R. VAN DER: “Sigma-Point Kalman Filters for Probabilistic Inference in Dynamic State-Space Models”, PhD dissertation. University of Washington, USA, 2004. 

16. ROMANO, N., B. BRUNONE & A. SANTINI: “Numerical analysis of one-dimensional unsaturated flow in layered soils”, Adv. Water Resour., 21: 315-324, 1998.

17. SORENSON, H. W.: “Least-squares estimation: from Gauss to Kalman”, Copyright (c) IEEE. Reprinted from IEEE Spectrum, (1970), 7: 63-68, July 2004.

18. TROCH, P.A., C. PANICONI & D. MACLAUGHLIN: “Catchment-scales hydrological modeling and data assimilation”. Adv. Water Res., 26: 131-135, 2003.

19.VAN GENUCHTEN, M.T.: “A closed-form equation for predicting the hydraulic conductivity of unsaturated soils”, Soil Sci. Soc. Am. J. 44: 892-898, 1980.

]]>

20.WALKER, J. P.: Estimating soil moisture profile dynamics from near-surface soil moisture measurements and standard meteorological data, PhD Thesis, the University of Newcastle, New South Wales, Australia, 1999.

21. WALKER, J. P., G. R. WILLGOOSE & J. D. KALMAN: “One-dimensional soil moisture profile retrieval by assimilation of near-surface observations: a comparison of retrieval algorithms”, Advances in Water Resources 24: 631-650, 2001.