Impacts of different types of measurements on estimating unsaturated flow parameters

Impacts of different types of measurements on estimating unsaturated flow parameters

Journal of Hydrology 524 (2015) 549–561 Contents lists available at ScienceDirect Journal of Hydrology journal homepage: www.elsevier.com/locate/jhy...

1MB Sizes 0 Downloads 19 Views

Journal of Hydrology 524 (2015) 549–561

Contents lists available at ScienceDirect

Journal of Hydrology journal homepage: www.elsevier.com/locate/jhydrol

Impacts of different types of measurements on estimating unsaturated flow parameters Liangsheng Shi a, Xuehang Song a, Juxiu Tong b,⇑, Yan Zhu a, Qiuru Zhang a a b

State Key Laboratory of Water Resources and Hydropower Engineering Sciences (Wuhan University), Wuhan, Hubei 430072, China Collage of Water Resources and Environmental Sciences, China University of Geosciences, Beijing 100083, China

a r t i c l e

i n f o

Article history: Received 29 June 2014 Received in revised form 27 January 2015 Accepted 30 January 2015 Available online 9 February 2015 This manuscript was handled by Corrado Corradini, Editor-in-Chief, with the assistance of Felipe de Barros, Associate Editor Keywords: Data assimilation Ensemble Kalman filter Unsaturated flow Data worth

s u m m a r y This paper assesses the value of different types of measurements for estimating soil hydraulic parameters. A numerical method based on ensemble Kalman filter (EnKF) is presented to solely or jointly assimilate point-scale soil water head data, point-scale soil water content data, surface soil water content data and groundwater level data. This study investigates the performance of EnKF under different types of data, the potential worth contained in these data, and the factors that may affect estimation accuracy. Results show that for all types of data, smaller measurements errors lead to faster convergence to the true values. Higher accuracy measurements are required to improve the parameter estimation if a large number of unknown parameters need to be identified simultaneously. The data worth implied by the surface soil water content data and groundwater level data is prone to corruption by a deviated initial guess. Surface soil moisture data are capable of identifying soil hydraulic parameters for the top layers, but exert less or no influence on deeper layers especially when estimating multiple parameters simultaneously. Groundwater level is one type of valuable information to infer the soil hydraulic parameters. However, based on the approach used in this study, the estimates from groundwater level data may suffer severe degradation if a large number of parameters must be identified. Combined use of two or more types of data is helpful to improve the parameter estimation. Ó 2015 Published by Elsevier B.V.

1. Introduction Unsaturated zone is the link between the land surface and the top of the phreatic zone. It is of great importance in providing water and nutrients to the groundwater environment. In agriculture, the unsaturated zone is essential for optimum water resources management, irrigation and drainage scheduling, fertilizer application, and crop production. Obtaining accurate soil hydraulic parameters is very important but not an easy task. These parameters are usually determined by laboratory experiments with repacked or intact soil cores. However, parameter estimates obtained from laboratory experiments may be considerably different from real parameters due to unavoidable disturbance to soil. There are also many in situ methods for direct estimation of the soil hydraulic parameters, including the crust method, the instantaneous method, and the unit gradient internal drainage (Lazarovitch et al., 2007). The main limitation of in situ methods is that they are time consuming due to the need of adhering to relatively strict initial and boundary conditions (Šimu˚nek and van Genuchten, 1996). Thus, many parameter ⇑ Corresponding author. E-mail address: [email protected] (J. Tong). http://dx.doi.org/10.1016/j.jhydrol.2015.01.078 0022-1694/Ó 2015 Published by Elsevier B.V.

estimation methods have been introduced to identify soil and aquifer hydraulic parameters. Overviews on inverse methods in hydrogeology are presented by Carrera et al. (2005), and recently by Zhou et al. (2014). Different types of observations are adopted for parameter inference, including the outflow flux in soil column experiments (Kool et al., 1985; van Dam et al., 1992, 1994; Toorman et al., 1992), soil water content (Zijlstra and Dane, 1996; Wang et al., 2003; Ritter et al., 2003; Yeh et al., 2005), soil water pressure head (Kool and Parker, 1988; Yeh and Zhang, 1996; Wöhling et al., 2008; Li and Ren, 2011), tracer test data (Mishra and Parker, 1989), cumulative infiltration data (Šimu˚nek and van Genuchten, 1996), and water evapotranspiration flux (Jhorar et al., 2002). In recent years, hydrogeophysical measurement methods, such as electrical resistivity tomography (ERT) and ground penetrating radar (GPR), are also widely used to collect information (Yeh et al., 2002; Kowalsky et al., 2005; Jadoon et al., 2012; Scholer et al., 2013). Furthermore, remote sensing is able to provide data for parameters estimation at larger scale (Vereecken et al., 2008; Montzka et al., 2011). Soil water pressure or soil matric potential is very sensitive to the alternation of soil wetting and drying. Matric potential data, often combined with measurements of flux data (e.g. outflow and

L. Shi et al. / Journal of Hydrology 524 (2015) 549–561

approach was not quite successful and only certain parameters could be identified if only the near-surface soil moisture data were used. With continuous monitored surface soil moisture, Qin et al. (2009) employed particle filter method to retrieve soil moisture profiles and to estimate soil hydraulic parameters simultaneously. Montzka et al. (2011) further investigated the effectiveness of data assimilation in estimating the parameters for different soil types. Although surface soil moisture data are widely used, the measurement accuracy is still a critical issue and its ability to estimate deep soil hydraulic parameters remains challenging. Measuring groundwater level is a common practice in field study. Groundwater level can be continuously measured by an automatic sensing device at a high accuracy and low cost. Groundwater level fluctuation indicates the groundwater recharge or discharge. Thus, it is a valuable data source for groundwater recharge estimation (Crosbie et al., 2005). It was widely used in traditional aquifer testing and parameter estimation in groundwater hydrology (McPherson, 1998). Cirkel et al. (2010) identified the seepage intensities from groundwater level data. However, the groundwater level data is not regarded as an important information source for studies on soil hydraulic parameter estimation. Many studies evaluated the worth of different types of data on the estimation of spatially distributed unsaturated flow parameters. Yeh and Zhang (1996) and Zhang and Yeh (1997) showed usefulness of pressure head and moisture content measurements for estimating unsaturated hydraulic parameters. Li and Yeh (1999) investigated effects of head and tracer concentration measurements on unsaturated flow parameters estimation. More recently, Mao et al. (2013) conducted a rigorous cross-correlation analysis to show the importance of head measurements at different times for estimating unsaturated hydraulic parameters. These studies, however, have not comprehensively compared the data worth of

70

Rainfall 60 50

Rainfall (mm)

evaporation), are the common data to estimate unsaturated parameters at the column scale (Kool et al., 1987; Toorman et al., 1992; Nutzmann et al., 1998). Šimu˚nek et al. (1998) used measured pressure head to estimate unsaturated flow parameters, and Marquardt–Levenberg optimization scheme was used in their study. They found that pressure heads measured near the soil surface are more valuable to parameter estimation than those measured at lower locations. Li and Ren (2011) showed that soil water pressure data can lead to satisfactory parameter estimates, especially for Ks and a in van Genuchten model. In practice, a major obstacle of estimating soil hydraulic parameters with soil water pressure data is the higher measurement error of pressure head. In addition, conventional matric potential measurement only covers the wet part of the moisture retention characteristic (Vereecken et al., 2008). For example, MPS-2 sensor produced by Decagon measures the range of soil water potentials between 9 and 100,000 kPa, with relative accuracy of 25% in the range of 9 kPa to 100 kPa (Decagon Devices Inc., 2014). Furthermore, soil matric potential can only be measured at the point scale, since current techniques do not allow to obtain spatially averaged measurement of matric potential. Unlike the pressure head data, soil moisture content data covers the whole changing range of soil water status. Soil water content can be measured at different scales, depending on the equipment and method used. Point-scale measurement technique is most sophisticated in the real applications. It has been widely used to monitor the vertical soil moisture profile. Some measurement methods (e.g., time domain reflectometry) are able to provide spatially and temporally high resolution measurements at low cost. Vereecken et al. (2008) presented a comprehensive review on the value of soil moisture measurements in vadose zone hydrology at both field and catchment scales. In-situ measurements of point-scale water content at different depths and times during infiltration or drainage experiments are already long used for parameter estimation (Kool et al., 1987). Kool et al. showed that the parameters a and n (in Van Genuchten model) can be solve uniquely using only information on water content profiles during drainage, but simultaneous estimation of three or more parameters requires additional information. Ross (1993) found that besides the measured water content, at least one matric potential was required to obtain good parameter estimates. Ritter et al. (2003) used the measured time series of soil water content to estimate unsaturated flow parameters. An optimization algorithm, global multilevel coordinate search algorithm, was used in their study. The illposedness of their inverse problem was partially due to a large number of parameters to optimize and to insufficient information in the measured data. In other words, although soil moisture data at different depths and times alone can be easily measured, they may be not able to identify the unsaturated parameters. Surface soil moisture is highly variable in space and time. It is defined as the moisture content averaged from the surface down to a given depth. Surface soil moisture can be easily captured over large regions with the rapid development of remote sensing technique. However, remote sensing technique characterizes soil moisture in the shallow soil, usually at a depth between 2 cm and 20 cm (Jackson et al., 1995). For example, the passive microwave L-band (1.4 GHz) sensor has a maximum penetration depth of 5 cm from the soil surface under minimal vegetation cover. Most studies focused on soil moisture profile retrieval by using near-surface soil moisture assimilation scheme (e.g., Walker et al., 2001; Walker and Houser, 2004; Dunne and Entekhabi, 2005). Several studies on estimating soil parameters with near-surface data have been reported in the literature. Ines and Mohanty (2008) employed a genetic algorithm to identify the soil water retention curve and the hydraulic conductivity function by using near-surface soil moisture (0–5 cm) data. For the case of layered soil systems, their

40 30 20 10 0

1

11

21

31

41

51

61

71

81

91

time (day) 10

Potential evapotraspiration 8

ET (mm)

550

6

4

2

0

1

11

21

31

41

51

61

71

81

time (day) Fig. 1. Synthetic rainfall and potential evapotranspiration.

91

551

L. Shi et al. / Journal of Hydrology 524 (2015) 549–561 Table 1 Soil parameters used in simulation. Ks/md1

Layer number

hs

hr

Layer Layer Layer Layer Layer

0.43 0.41 0.45 0.39 0.43

0.078 0.057 0.067 0.1 0.045

1: 2: 3: 4: 5:

loam loamy sand silt loam sandy clay loam sand

a/m1

n

Mean (l)

SD (r)

Mean (l)

SD (r)

Mean (r)

SD (l)

0.2496 3.5016 0.108 0.3144 7.128

0.4368 2.7264 0.2952 0.6576 3.744

3.6 12.4 2 5.9 14.5

2.1 4.3 1.2 3.8 2.9

1.56 2.28 1.41 1.48 2.68

0.11 0.27 0.12 0.13 0.29

1

7 Layer 1 6

Ks (m/d)

Ks (m/d)

0.8

0.6

0.4

0.2

reference

SWH (5%)

SWH (20%)

SWC(0.02)

SSWC(0.02)

SSWC(0.05)

GWL (0.01m)

0 0

10

20

40

50

Layer 2

4 3

GWL (0.05m)

30

5

2 60

70

80

0

90

10

20

30

40

50

Time (d)

Time (d)

(a)

(b)

60

70

80

90

60

70

80

90

1.2

0.5

Layer 4

Layer 3 1

0.4

Ks (m/d)

Ks (m/d)

0.8

0.3 0.2

0.6 0.4

0.1

0.2

0

0

0

10

20

30

40

50

60

70

80

90

0

10

20

30

40

50

Time (d)

Time (d)

(c)

(d) 12

Ks (m/d)

Layer 5

10

8

6 0

10

20

30

40

50

60

70

80

90

Time (d)

(e) Fig. 2. Estimated Ks values of each layer using different types of data when only estimating Ks.

four most popular types of data (i.e., point-scale soil water head, point-scale soil moisture content, surface soil moisture content and groundwater level) on estimating unsaturated flow parameters.

In this work, we present an ensemble Kalman filter (EnKF) framework that considers these four types of data. The data worth of each data is explored through synthetic examples. In particular, we are

552

L. Shi et al. / Journal of Hydrology 524 (2015) 549–561

8

20

Layer 2

Layer 1 18

α (1/m)

α (1/m)

6

4

16 14

2

reference

SWH (5%)

SWH (20%)

SWC(0.02)

SSWC(0.02)

SSWC(0.05)

GWL (0.01m)

GWL (0.05m)

12

0

10 0

4

10

20

30

40

50

60

70

80

90

0

10

20

30

40

50

Time (d)

Time (d)

(a)

(b) 12

Layer 3

60

70

80

90

60

70

80

90

Layer 4

α (1/m)

α (1/m)

10 3

8

2 6

4

1 0

10

20

30

40

50

60

70

80

0

90

10

20

30

40

50

Time (d)

Time (d)

(c)

(d) 18

Layer 5

α (1/m)

17

16

15

14 0

10

20

30

40

50

60

70

80

90

Time (d)

(e) Fig. 3. Estimated a values of each layer using different types of data when only estimating a.

interested in: (1) How much information regarding the van Genuchten model parameters does different types of data contain? (2) Which factors may affect the accuracy of parameter estimation based on different types of data? 2. Theory 2.1. One-dimensional model of variably saturated flow The vertical one-dimensional soil water flow equation can be written as

   @h @h @ @h KðhÞ 1 þ x Ss þ ¼ @t @t @z @z

where Ss is the specific storage, x is equal to one if the soil is saturated and zero if the soil is unsaturated, K(h) is the unsaturated hydraulic conductivity [L/T], h is the pressure head [L], h is positive or equal to zero when the soil is saturated and negative when the soil is unsaturated, h is the volumetric soil water content, t is the time [T], z denotes the vertical dimension [L], assumed positive upward. Ss is ignored in this study since it is very small. Solution of Richards’ equation requires knowledge of unsaturated conductivity and water content versus hydraulic head. In this work, the van Genuchten model is used to describe these relationships,

hðhÞ ¼ hr þ ð1Þ

hs  hr ½1 þ jahjn 

m

ð2Þ

553

L. Shi et al. / Journal of Hydrology 524 (2015) 549–561

2.7

1.7

Layer 2

Layer 1 2.6

1.65

n

n

2.5

1.6 2.4

1.55

reference

SWH (5%)

SWH (20%)

SWC(0.02)

SSWC(0.02)

SSWC(0.05)

GWL (0.01m)

1.5

0

10

20

30

2.3

GWL (0.05m)

40

50

2.2

60

70

80

90

0

10

20

30

40

50

(a)

70

80

90

60

70

80

90

(b)

1.6

1.7 Layer 3

n

60

Time (d)

Time (d)

Layer 4 1.65

1.5

1.6

n

1.55

1.45

1.55

1.4

1.5

1.35

1.45 1.4

1.3 0

10

20

30

40

50

60

70

80

0

90

10

20

30

40

50

Time (d)

Time (d)

(c)

(d) 3 Layer 5

n

2.9

2.8

2.7

2.6

0

10

20

30

40

50

60

70

80

90

Time (d)

(e) Fig. 4. Estimated n values of each layer using different types of data when only estimating n.

m 2

1=m KðhÞ ¼ K s S0:5 Þ  e ½1  ð1  Se

ð3Þ

where Ks is the saturated hydraulic conductivity [L/T], hs is the saturated water content, hr is the residual water content, a is a measure of the first moment of the pore size density function [1/L], n is an inverse measure of the second moment of the pore size density functions, m = 1–1/n, and Se is the saturation degree, Se ¼ ðh  hr Þ=ðhs  hr Þ. In this study, the soil water flow equation is solved by the popular saturated–unsaturated model HYDRUS-1D (Simunek et al., 2005), which ignores the first term on the left-hand side of Eq. (1). Some new solvers such as Ross method (Zha et al., 2013) can also be employed to solve the variably saturated flow. The model uncertainty is ignored in this study.

2.2. Ensemble Kalman filter to assimilate different types of data The EnKF, first developed by Evensen (1994), has been applied to water resources problems (e.g., Moradkhani et al., 2005; Chen and Zhang, 2006; Liu et al., 2008). The EnKF is essentially a sequential Monte Carlo method based on the Bayes’ theorem, using an ensemble of realizations to approximate the statistics of model independent variables (e.g., saturated hydraulic conductivity) and state variables (e.g., pressure head). The EnKF was traditionally used to update hydraulic conductivity by assimilating hydraulic head data in the groundwater problems. Chen and Zhang (2006) showed that EnKF provides an efficient approach for obtaining satisfactory estimation of the hydraulic conductivity field with dynamic head measurements. Hendricks Franssen (2008)

554

L. Shi et al. / Journal of Hydrology 524 (2015) 549–561

investigated the filter inbreeding problem when hydraulic heads and transmissivities are jointly updated. Li et al. (2012) demonstrated the capability of the EnKF to estimate hydraulic conductivity and porosity by assimilating dynamic head and multiple concentration data. Some recent work focused on improving the computational efficiency of EnKF. Shi et al. (2012) proposed an EnKF based on the mutiscale finite element method to reduce the computational cost in the large-scale groundwater problem. Panzeri et al., 2013 coupled EnKF with stochastic moment equations of transient groundwater flow to circumvent the need for computationally intensive MC simulations. In this study, we consider the soil hydraulic parameters Ks, a, and n as random unknowns, while hr and hs are assumed to be deterministically known. The augmented state vector is defined as, T

ð4Þ

where mk and hk are the vectors of state parameter (i.e., Ks, a, and n) and state variable (i.e., soil water head) at time t k , respectively. Ny is the dimension of yk : Ny = Nm + Nh, where Nh is the node number in the flow system and Nm is the number of the parameters to be estimated. The observation vector at time step k for each ensemble member is,

dj;k ¼ dk þ ej;k

ð5Þ

where dk is the observation at time t k , ej;k are independent white noises for the observations, j is the ensemble number index, dj;k is the head measurement vector at time tk for ensemble index j. The dimension of dj;k is denoted as Nd. For any ensemble member j at a given time t k , the state vector is updated by combining model predictions and observations, f f yaj;k ¼ yj;k þ Kk ðdj;k  Hyj;k Þ

ð6Þ

where superscripts a and f refer to model analysis and model foref is the initially guessed or estimated state cast or initial guess, yj;k

variables and parameter for realization j at time t k based on information at time k-1, yaj;k is their improved estimates for realization j at time t k by conditioning on the observed information at time t k , and H is observation operator which represents the relationship between the state vector and the observation vector. The expression of H depends on the type of data to be discussed below. The Kalman gain K is defined as,

Kk ¼ Ckf HT ½HCkf HT þ CDk 

1

ð7Þ

where Ckf is the covariance matrix of state vector at time tk , and CDk is the error covariance matrix of the observations. Ckf is given by

Ckf 

Ne h ih iT  1 X f f yj;k  hykf i yj;k  hykf i Ne  1 j¼1

hykf i

N

sf 1X li h i L i¼1

hsf ¼

ð9Þ

where hi is the soil water content at node i, li is the control length of PNsf li ; Nsf is the number of node i, L is the observation depth, L ¼ i¼1 nodes contained from the soil surface to the observation depth. If the groundwater level is observed, the forecast groundwater depth zd is obtained through linear interpolation between the nearest node above (node i) and below (node i + 1) the water table,

zd ¼

hiþ1 zi  hi ziþ1 hiþ1  hi

ð10Þ

where hiþ1 P 0, and hi 6 0. In this study, two or more types of data can be used jointly to estimate parameters. For example, when soil water content and groundwater level are introduced together, H maps the simulated heads at the grid nodes on both types of measurement. Although EnKF does not require linearization of nonlinear systems, it may not be suitable to parameter estimation for highly nonlinear systems. Considerable attention has been paid to resolve the problems of nonlinearity in ensemble data assimilation, and various approaches have been developed, such as particle filter method (e.g., Moradkhani, 2005; Montzka et al., 2011), Markov Chain Monte Carlo method (e.g., Zeng et al., 2012), iterative EnKF (e.g., Song et al., 2014), and EnKF based on transformation techniques (Zhou et al., 2011; Schöniger1 et al., 2012). Some nonensemble methods such as successive linear estimator (SLE) are able to handle non-linearity and spatial structure of media (Yeh et al., 1996; Zhang and Yeh, 1997). SLE uses a similar predictioncorrection mathematical approach in a successive manner to include the nonlinear relationships. It has been used to estimate unsaturated flow parameters by assimilating head, soil moisture content, concentration, and electrical resistivity measurements (Zhang and Yeh, 1997; Yeh and Šimu˚nek; 2002; Liu and Yeh, 2004; Berg and Illman, 2012). In this study, EnKF is employed due to its numerical simplicity and good adaptability in assimilating various types of data. 3. Numerical examples Synthetic numerical experiments are conducted to verify the ability of EnKF to estimate soil hydraulic parameters by assimilating different types of data, and to investigate their respective data

ð8Þ 0 Layer 1

ykf ,

where denotes the ensemble mean of Ne is total number of realizations in the ensemble. At the first time step (time zero), the covariance is defined a priori. In this study, the initial heads are assumed to be known and deterministic, and all the parameters to be estimated are independent of each other. The ensembles members of parameters are generated with prior mean and variance. In the case of assimilating point-scale soil water head measurements, H is a Nd  Ny matrix, in which the elements corresponding to the location of head measurement are 1 and zero otherwise. When assimilating point-scale soil water content data, the simulated head at the observation location is first converted to moisture content using Eq. (2), and H is still a matrix with 0 s and 1 s as its elements. If the surface soil moisture data are available, the observed data is the averaged soil water content within

0.2

Groundwater level (m)

yk ¼ ½mk ; hk 

the specified observation depth. Forecast surface moisture content hsf is obtained by the arithmetic average of all the moisture contents from the surface to the observation depth,

Case GWL1

0.4

Layer 2

Case GWL2

0.6

Case GWL3

0.8

Layer 3

1 1.2

Layer 4

1.4 1.6

Layer 5

1.8 0

10

20

30

40

50

60

70

80

90

Time (d) Fig. 5. Temporal changes of groundwater level under different rainfall patterns.

555

L. Shi et al. / Journal of Hydrology 524 (2015) 549–561 1.7

8

Ks of Layer 1

α of Layer 1

7 6

0.6

5

0.4

1.65

n of Layer 1

n

0.8

α (1/m)

Ksv (m/d)

1

4 3

1.55

2

0.2

reference Case GWL1

Case GWL2 Case GWL3

1

0 0

10

20

30

40

50

60

70

80

1.5

0

90

0

10

20

30

40

50

60

70

80

0

90

20

6

18

2.6

5

16

2.5

4 K of Layer 2 s 3 2 10

20

30

40

50

60

70

80

n

14

2.4

12

2.3

10

90

4

10

20

30

40

50

60

70

80

0

90

70

80

90

10

20

30

40

50

60

70

80

90

60

70

80

60

70

80

90

60

70

80

90

1.6

n of Layer 3

1.55

3

1.5

n

α (1/m)

60

Time (d)

α of Layer 3

Ks of Layer 3

0.2

50

n of Layer 2

Time (d)

0.3

40

2.2 0

0.4

Ks (m/d)

30

2.7

α of Layer 2

Time (d) 0.5

20

Time (d)

7

0

10

Time (d)

α (1/m)

Ks (m/d)

Time (d)

1.45 1.4

2 0.1

1.35

0 0

10

20

30

40

50

60

70

80

1.3

1

90

0

10

20

30

Time (d)

40

50

60

70

80

0

90

10

20

30

Time (d) 12

1.2

1.6 8

n

α (1/m)

0.6

90

n of Layer 4

1.65

10

0.8

50

1.7

α of Layer 4

Ks of Layer 4

40

Time (d)

1

Ks (m/d)

1.6

1.55 1.5

0.4 6

1.45

0.2 0 0

10

20

30

40

50

60

70

80

1.4

4

90

0

10

20

30

Time (d)

40

50

60

70

80

0

90

30

50

2.9

n of Layer 5

α of Layer 5 16

n

α (1/m)

10

40

3

Ks of Layer 5 17

Ks (m/d)

20

Time (d)

18

12

10

Time (d)

2.8

8 2.7

15

6

0

10

20

30

40

50

Time (d)

60

70

80

90

2.6

14 0

10

20

30

40

50

60

Time (d)

70

80

90

0

10

20

30

40

50

Time (d)

Fig. 6. Estimated Ks, a and n values (from left to right) for each layer (layer 1, 2, 3, 4 and 5, from top to bottom) under different rainfall patterns.

556

L. Shi et al. / Journal of Hydrology 524 (2015) 549–561

8 1

Layer 2 Layer 1

0.8

Ks (m/d)

Ks (m/d)

6 0.6 0.4

4

0.2 0

0

10

20

reference

Case GWL1

Case GWL3

Case GWL4

30

40

50

60

70

80

2

90

0

10

20

30

40

50

60

Time (d)

Time (d)

(a)

(b)

70

80

90

Fig. 7. Estimated Ks values for Layer 1 and 2 by starting data assimilation at different groundwater levels.

worth. All the experiments involve a one-dimensional vertical flow in a soil column of 180 cm. The soil column is devided into five layers: 0–20 cm (top layer, Layer 1), 20–60 cm (Layer 2), 60–100 cm (Layer 3), 100–140 cm (Layer 4), 140–180 cm (Layer 5). The flow domain is discretized into 41 grids. The grids nearest to the interface of two layers have the sizes of 2.5 cm, and the remaining grids have the sizes of 5.0 cm. The total simulation time is 3 months (92 days). We impose the atmospheric boundary as the top boundary condtion and zero flux boundary condition as the bottom boundary condtion. The boundary conditions are deterministic and known in this study. Initial groundwater table is at the depth of 150 cm, located at the Layer 5. Synthetic rainfall and potential evapotranspiration are shown in Fig. 1. The upper boundary conditions for all the cases are based on Fig. 1 unless otherwise specified. Initial total heads for all the nodes are 150 cm, thus there is no water movement at the begining. Synthetic observations including point-scale soil water content, point-scale soil water head, surface soil water content, and groundwater level, are generated by running Hydrus-1D with the ‘‘true’’ parameters, which are given according to Carsel and Parrish (1988). Table 1 lists the representative statistics of the unsaturated flow parameters for different types of soil, where the means are set as the ‘‘true’’ values of the parameters. The generated observations are pertubated by the

2

2 Estimated Ks of Layer 1 from the soil water head data

1 0.5

0.5

0 0

10

20

30

40

50

60

70

80

90

0

10

20

30

50

60

70

80

Estimated Ks of Layer 2 from the soil water head data

Estimated Ks of Layer2 from the surface soil water content data

4 2

6 4

20

30

40

50

Time (d)

10

20

30

60

70

80

90

0

40 50 60 Time (d)

70

80

90

Estimated Ks of Layer 2 from the groundwater level data

6 4 2

2 10

0

8 Ks (m/d)

8

6

0

0

90

10

10

Ks (m/d)

Ks (m/d)

40

Time (d)

10

0

1 0.5

Time (d)

8

Estimated Ks of Layer 1 from the groundwater

1.5 level data Ks (m/d)

reference 0.8μ μ+σ μ+2σ

1

2

Estimated Ks of Layer1 from the surface soil water content data

1.5 Ks (m/d)

Ks (m/d)

1.5

0

observation errors, which are assumed to be Gaussian. The pointscale observations are located at the center of each layer, and the observation depth of surface soil moisture data is 20 cm. In this study, hs and hr are assumed to be known. The ensemble size of EnKF is taken as 500, which has been tested large enough for our numerical investigation. Given the fact that measurement error of soil water matric potential may be very large, two relative measurement errors (5% and 20%) are considered for soil water head data. The associated error of point-scale soil water content data is set at a small value (0.02 m3/m3), while two errors (0.02 m3/m3 and 0.05 m3/m3) are set for surface soil water content to take into account the difficulty of inferring accurate soil water content from remote sensing technique or ground penetrating radar. Absolute errors of 0.01 m and 0.05 m are given to the groundwater level measurements. The data assimilation for all the cases is implemented on a daily basis from time t = 1 d unless otherwise specified. When implementing EnKF, the initial realizations of Ks, a and n are generally generated with the prior knowledge of the unknown parameters. In our synthetic experiments, the initial mean of the parameter to be identified is set as at three levels: l + r, l + 2r, 0.8l, where l and r are presented in Table 1. The goal is to investigate the effect of prior mean on the estimates using different

0

10

20

30

40 50 60 Time (d)

70

80

90

0

0

10

Fig. 8. Estimated Ks values for Layer 1 and 2 under different initial guess values.

20

30

40 50 60 Time (d)

70

80

90

557

L. Shi et al. / Journal of Hydrology 524 (2015) 549–561

types of data. In the following text, all the results are obtained with the initial mean l + r unless otherwise specified.

In this section, only one parameter Ks, a or n is estimated at a time, and the other two parameters are assumed to be known. The true values of the unknown and known parameters are listed as the means (l) in Table 1. The estimated results by using point-scale soil water head, point-scale soil water content, surface soil water content, and groundwater level are denoted as ‘‘SWH’’, ‘‘SWC’’, ‘‘SSWC’’, and ‘‘GWL’’, respectively.

of the estimates deteriorates considerably, especially for parameter n. High observation accuracy seems to be required if groundwater level is used as the conditioning data. Figs. 2 and 3 show that with point-scale soil water content data, parameters a and n of all the layers approach to their expected values very quickly (within 20 days), even when water does not move in the deep layer. It should be pointed out that this is just a numerical coincidence since the initial condition in our simulation is specified by soil water head, some information about a and n has been indicated through Eq. (2) if soil water content at t = 0 is given. This phenomenon implies the importance of jointly using different types of data.

3.1.1. Comparison of results Based on prior mean l + r, Figs. 2–4 respectively show the estimated values of Ks, a and n by assimilating the four types of measurements. The potential impacts of prior mean on the estimates are discussed below. From Fig. 2(a), except the case of using surface soil water content data (with an error of 0.05 m3/m3) and the case of using groundwater level data, most cases reach the reference saturated hydraulic conductivity of Layer 1 within 92 days. Smaller measurement errors result in faster convergence to the true values regardless of the observations types. In this study, by assimilating point-scale soil water head at five observation points, a mild range of observation errors (5–20%) does not significantly influence the final estimates. However, the observation error of 0.05 m3/m3 in surface water content data leads to obviously deteriorated results. Similarly, by increasing the measurement error of groundwater level data to 0.05 m, the quality

3.1.2. Estimates from surface soil moisture data Figs. 2–4 also demonstrate that, while assimilating near-surface soil moisture data with small measurement error (0.02 m3/m3) successfully estimates the soil hydraulic parameters of the top layer, the data assimilation does not work well in Layers 2–5. The estimated a and n values of Layer 2–4 show the trend of converging to the true values, but did not approach stable estimates during the simulation period. The better results for top layers attributed to the stronger relationship between surface soil moisture data and soil hydraulic properties of the top layers. In the work of Ines and Mohanty (2007, 2008), they found that surface soil moisture data alone are not adequate to identify the parameters of the deep layer. Additional soil moisture data from deeper depths or evapotranspiration data were suggested for improving the estimation of soil hydraulic properties of heterogeneous systems. Under a large measurement error (0.05 m3/m3), it is noticed

3.1. Estimating individual parameters, Ks, a or n

8

1

Layer 2

Layer 1 0.8

Ks (m/d)

Ks (m/d)

6

0.6

0.4

4 reference SWH (20%) SSWC(0.02) SWC(0.02)+GWL(0.01m)

0.2

0

0

10

20

30

SWH (5%) SWC(0.02) GWL (0.01m) SSWC(0.02)+GWL(0.01m)

40

50

60

70

80

2

90

0

20

30

40

50

Time (d)

(a)

(b)

60

70

80

90

60

70

80

90

1.2

0.5

Layer 3

Layer 4

0.4

1

0.3

0.8

Ks (m/d)

Ks (m/d)

10

Time (d)

0.2

0.1

0.6

0.4

0

0.2

0

10

20

30

40

50

Time (d)

(c)

60

70

80

90

0

10

20

30

40

50

Time (d)

(d)

Fig. 9. Estimated Ks values of each layer using different types of data when simultaneously estimating Ks, a and n.

558

L. Shi et al. / Journal of Hydrology 524 (2015) 549–561

that, while the estimated Ks of the top layers gradually approaches to the true value, it fails to reach the true value during the simulation period. We expect that it will eventually produce the right estimation if more data are introduced. It is emphasized that the ability of assimilating surface soil moisture data may also depend on the soil type and the observation depth of near-surface measurements (Montzka et al., 2011; Walker et al., 2001). 3.1.3. Estimation using groundwater level data It is seen from Fig. 2 that the groundwater level data (with the measurement error of 0.01 m) produce satisfactory estimates of the saturated hydraulic conductivities of Layers 2–4. Similarly, Figs. 3 and 4 show that the groundwater level data of 92 days contain sufficient information for estimating parameters a and n. From the perspective of groundwater recharge, since the rising rate of groundwater level is determined in part by the soil hydraulic parameters, groundwater level data actually is one type of information that builds a connection with the parameters. EnKF is a useful tool to explore the data worth implied by the groundwater level data. To test the performance of EnKF under different rainfall patterns and to further evaluate the data worth of groundwater level data on parameter estimation, two more numerical experiments (Cases GWL2 and GWL3) are designed. To facilitate analysis, the case of assimilating groundwater level data (with error of 0.01 m) based on the upper boundary condition in Fig. 1 is called Case GWL1. In Case GWL2, the daily rainfall is specified as a constant value (4.635 mm/d) during the time interval from t = 0 to t = 40 d. In Case GWL3, a uniform rainfall rate of 5.54 mm/d is used

7

Layer 2 16

α (1/m)

α (1/m)

6

18

SWH (5%) SWC(0.02) GWL (0.01m) SSWC(0.02)+GWL(0.01m)

reference SWH (20%) SSWC(0.02) SWC(0.02)+GWL(0.01m)

Layer 1

during the whole simulation time. Fig. 5 shows the temporal changes of groundwater level in Case GWL1-3. It is seen that the groundwater level of Case GWL3 starts increasing earlier. However, the increasing rate is slower than that of Cases GWL1 and Case GWL2. The estimated values of Ks, a and n from Case GWL1-3 are compared in Fig. 6. According to these figures, most estimates converge to the same values, regardless of rainfall patterns. Generally, the estimates reach the true values in a slower rate when the rising rate of groundwater table is slower, especially for parameters a and n. As discussed above, the change of groundwater level is determined by all the parameters of the soil layers above the water table. Since different parameter combinations may result in the same groundwater level, parameter estimation based only on the groundwater level data will be difficult if multiple parameters need to be identified simultaneously. Thus for heterogeneous soil and when the water table is deep, the groundwater level data may fail to provide accurate estimates. Two additional scenarios (Cases GWL4 and GWL5) are conducted to investigate the impact of initiating the data assimilation at different groundwater levels. The two cases start data assimilation at t = 45 d (water table located at Layer 4, z = 1.1 m) and t = 48 d (water table located at Layer 3, z = 0.7 m), respectively. Except the starting time of data assimilation, all other configurations for the two cases are the same as those of Case GWL1. For all the cases, the assimilation frequency is one day a time after the data assimilation starts. At the beginning time, the groundwater level data in Cases GWL1, GWL4 and GWL5 indicate the information about five, four and three parameters (Ks,

5

4

14

12

3

10 0

10

20

30

40

50

60

70

80

90

0

10

20

30

40

50

Time (d)

Time (d)

(a)

(b)

60

70

80

90

60

70

80

90

4 Layer 3

Layer 4

11

α (1/m)

α (1/m)

3

2

9

7

5

1 0

10

20

30

40

50

60

70

80

90

0

10

20

30

40

50

Time (d)

Time (d)

(c)

(d)

Fig. 10. Estimated a of each layer using different types of data when simultaneously estimating Ks, a and n.

559

L. Shi et al. / Journal of Hydrology 524 (2015) 549–561

a or n) respectively, depending on the number of layers above the water table. The estimated Ks values of Layers 1–2 are presented in Fig. 7. It is noted that Cases GWL4 and GWL5 produce better results than Case GWL1 does even when only 48 and 44 groundwater level data are used for Cases GWL4 and GWL5, respectively (92 data are used for Case GWL1). This may be because the groundwater level data in Case GWL1 connect with more parameters than in Cases GWL4 and GWL5, and this more complex connection hinders the correct estimation. Although the groundwater level data are collected easily with the lowest cost, their data worth may depend on the complexity between the data and parameters, or conversely, depend on the soil heterogeneity. It is expected that the estimates will deteriorate if over-complex connections between the data and the parameters exist or in the highly heterogeneous soil.

3.1.4. Estimates using different initial guess Fig. 8 plots the estimated Ks values of Layers 1 and 2 based on three different prior means as a function of assimilation steps. According to these figures, all estimates using the initial guess of 0.8l converge to the true values. The initial guess does not appear to have significant impact to the final estimates when assimilating soil water head data at five observation points. Similar phenomenon is observed by assimilating soil water content data at five observation points (results not shown). However, the estimates conditioning on the surface soil water content data and groundwater level data deteriorate with a less accurate initial guess. Apparently, the data worth implied by the surface soil water

content data and groundwater level data is corrupted by inaccurate initial guess. 3.2. Simultaneously estimating Ks, a and n Identification of soil hydraulic parameters of heterogeneous soil is challenging, when there are a large number of parameters to estimate or when information contained in the measured data is not sufficient [Ritter et al., 2003]. In this section, parameters Ks, a and n of five layers (a total of 15 parameters) are determined simultaneously. The initial guesses of all the parameters are l + r. To explore the effect of assimilating two types of data together, two more scenarios are considered to jointly use point-scale soil water content and groundwater level data and to jointly use surface soil water content and groundwater level data. The simulating results are denoted by ‘‘SWC + GWL’’ and ‘‘SSWC + GWL’’, respectively, in Figs. 9–11. 3.2.1. Results comparison Figs. 9–11 show that with the fusion of dynamic soil water head data, most estimates are able to converge to their true values, but the estimated Ks and n of Layer 2 fail to approach the true values when the measurement error increases from 5% to 20%. The estimates based on surface soil water content data with 20% relative error and groundwater level data with 0.05 m absolute error are not shown here. Similarly, the estimates based on these types of data deteriorate with increasing measurement error. Comparing with the cases of estimating single parameter in the previous section, the quality of the estimate degrades more with the increasing

1.7

2.7 Layer 1 Layer 2

2.6 reference SWH (20%) SSWC(0.02) SWC(0.02)+GWL(0.01m)

1.6

SWH (5%) SWC(0.02) GWL (0.01m) SSWC(0.02)+GWL(0.01m)

2.5

n

n

1.65

2.4 2.3

1.55 2.2 1.5

2.1 0

10

20

30

40

50

60

70

80

90

0

10

20

30

40

50

Time (d)

Time (d)

(a)

(b)

1.6

60

70

80

90

60

70

80

90

1.7 Layer 3

Layer 4

1.55 1.6

1.45

n

n

1.5 1.5

1.4 1.4 1.35 1.3

1.3 0

10

20

30

40

50

60

70

80

90

0

10

20

30

40

50

Time (d)

Time (d)

(c)

(d)

Fig. 11. Estimated n of each layer using different types of data when simultaneously estimating Ks, a and n.

560

L. Shi et al. / Journal of Hydrology 524 (2015) 549–561

measurement error when dealing with multiple parameters. In a recent work of Mao et al. (2013), it was shown that the quality of the estimate varies with the accuracy of observations, and the final estimate may deviate significantly from the true values due to the increasing observation error. They also found that the quality of the estimates improves as the number of observations in time and space increases, despite of the presence of observation error. Our numerical experiments indicate the necessity of high accuracy measurement to guarantee estimation accuracy if a large number of parameters must be identified. 3.2.2. Estimates from surface soil moisture data Unlike in Section 3.1, under the situation of estimating multiple parameters simultaneously, surface soil water content data only lead to obvious improvement to the estimates of Layer 1 (Figs. 9–11), but no longer exerts strong influence on the parameter estimation of Layer 2 or deeper layers. Furthermore, larger differences between the final estimated Ks, a and n values of Layer 1 and the true values are observed, comparing to the perfect match in Section 3.1. 3.2.3. Estimates from groundwater level data Figs. 9–11 show that assimilating groundwater level data cannot provide globally optimal estimates for all the 15 parameters. For example, the estimated n of Layer 3 and 4 are close to their corresponding true values, while the estimated a of Layer 2-4 deviate significantly from the true values. Comparing them with the results in Section 3.1, the estimates based on groundwater level data suffer severe degradation if multiple parameters are identified. This may be due to the reason that different parameter combinations can lead to the same model response (i.e., groundwater level), especially when a large number of parameters are involved. As suggested by Hopmans and Šimunek, 1999, an increase in the number of estimated parameters entails the need for further measurements of different types. 3.2.4. Estimates from multiple types of data By jointly using soil water content data and groundwater level data, the produced estimates are generally improved, comparing to the cases of only using one type of measurements. For example, estimated Ks of Layer 4 using soil water content data alone deviate significantly from the true value (Fig. 9a). With additional groundwater level data, the estimate is improved. Similarly, estimated a of Layers 2 and 4 and estimated n of Layer 1-2 based on groundwater level data can be largely improved by adding soil water content data. Although some estimates by conditioning on the two types of data degrade slightly (Figs. 9a and 11d) when comparing them with the results obtained by using either groundwater level data or soil water content data, jointly using more types of data seems to be an effective way to improve the estimation. One more case of assimilating both surface soil moisture data and groundwater level data is conducted to evaluate the potential benefit of using different types of data together. Overall, combination uses of different types of data result in improved estimates of the hydraulic parameters in multi-layered soil, although not all parameter estimates improve (see Figs. 9d, 11c and d). 4. Discussion and conclusion In this paper we present an EnKF framework for considering four most popular types of measurement, i.e., point-scale soil water head, point-scale soil water content, surface soil moisture, and groundwater level. Synthetic experiments are designed to compare the data worth of different types of measurements for estimating soil hydraulic properties. The parameter estimates

based on different types of data show different patterns of temporal evolutions when the data are sequentially assimilated. The estimation accuracy depends on the type of measurements, the boundary condition (such as rainfall pattern), the prior information, and the magnitude of the measurement errors. For all the data types, smaller measurements errors lead to faster convergence to the true values. Successful estimation based on surface soil moisture data and groundwater level data requires high accuracy measurement. All estimates from point-scale data (five observation points in this study), regardless of the initial guesses, converge to the true values. However, the data worth implied by the surface soil water content data and groundwater level data is prone to corruption by an inaccurate initial guess. Surface soil moisture data are capable of identifying soil hydraulic parameters for the top layers, but exert less or no influence to the deeper layers especially when estimating multiple parameters simultaneously. This may be due to the direct relationship between the soil hydraulic properties of the top layers and the surface soil moisture data. The serial data of groundwater level contain valuable information to infer parameters Ks, a and n, and the estimates generally converge to the same values regardless of the rainfall patterns. This is very attractive since groundwater level data at the monitoring wells can be easily collected at low cost. However, in highly heterogeneous soils (or equivalently a large number of unknown parameters involved) the over-complex connection between groundwater level data and soil parameters may reduce its worth in parameter estimation. Our results demonstrate that the estimates based on groundwater level data suffer from severe degradation if a large number of parameters must be identified. An increase in the number of estimated parameters entails the need of assimilating different types of data together. Jointly using point-scale soil water content data and groundwater level data and jointly using surface soil moisture content and groundwater level data are investigated in this study. Overall, assimilating more types of data is helpful to improve the estimates of soil hydraulic parameters. As illustrated in this study, EnKF is a useful approach to estimate unsaturated flow parameters. Different types of data can be easily assimilated in the framework of EnKF. We note that the data worth contained in different types of data may further depend on soil types, locations of observation points, layer structures, and initial conditions. Other types of data (such as flux, tracer, and heat data) can also be employed for the parameter estimation, which however is not investigated in this present study. Our future interest is to develop a numerical framework that can assimilate different types of data at various spatial scales. Acknowledgements This work is supported by Natural Science Foundation of China through Grants 51179132, 51279141, 51328902 and 51209187. We thank the two anonymous reviewers for their valuable comments which helped to improve the quality of this manuscript. References Berg, S.J., Illman, W.A., 2012. Improved predictions of saturated and unsaturated zone drawdowns in a heterogeneous unconfined aquifer via transient hydraulic tomography: laboratory sandbox experiments. J. Hydrol. 470–471, 172–183. Carrera, J., Alcolea, A., Medina, A., Hidalgo, J., Slooten, L.J., 2005. Inverse problem in hydrogeology. Hydrogeol. J. 13 (1), 206–222. Carsel, R.F., Parrish, R.S., 1988. Developing joint probability distributions of soil water retention characteristics. Water Resour. Res. 24, 755–769. Chen, Y., Zhang, D., 2006. Data assimilation for transient flow in geologic formations via ensemble Kalman filter. Adv. Water Resour. 29 (8), 1107–1122. Cirkel, D.G., Witte, J.P.M., van der Zee, S.E.A.T.M., 2010. Estimating seepage intensities from groundwater level time series by inverse modelling: a sensitivity analysis on wet meadow scenarios. J. Hydrol. 385, 132–142.

L. Shi et al. / Journal of Hydrology 524 (2015) 549–561 Crosbie, R.S., Binning, P., Kalma, J.D., 2005. A time series approach to inferring groundwater recharge using the water table fluctuation method. Water Resour. Res. 41, W01008. Decagon Devices Inc., 2014. MPS-2 Dielectric water Potential Sensor Operator’s Manual. Dunne, S., Entekhabi, D., 2005. An ensemble-based reanalysis approach to land data assimilation. Water Resour. Res. 41, W02013. Evensen, G., 1994. Sequential data assimilation with a nonlinear quasi-geostrophic model using Monte Carlo methods to forecast error statistics. J. Geophys. Res. 99 (C5), 10143–10192. Hendricks Franssen, H.J., Kinzelbach, W., 2008. Real-time groundwater flow modeling with the ensemble Kalman filter: joint estimation for states and parameters and the filter inbreeding problem. Water Resour. Res. 44, W09408. Hopmans, J.W., Šimunek, J., 1999. Review of inverse estimation of soil hydraulic properties. In: van Genuchten, M.Th., Leij, F.J., Wu, L. (Eds.), Characterization and measurement of the hydraulic properties of unsaturated porous media. University of California, Riverside, CA, pp. 643–659. Ines, A.V.M., Mohanty, B.P., 2007. Near-Surface soil moisture assimilation for quantifying effective soil hydraulic properties under different hydroclimatic conditions. Vadose Zone J. 7, 39–52. Ines, A.V.M., Mohanty, B.P., 2008. Near-surface soil moisture assimilation for quantifying effective soil hydraulic properties using genetic algorithm: 1. Conceptual modeling. Water Resour. Res. 44, W06422. Jackson, T.J., Levine, D.M., Swift, C.T., Schmugge, T.J., Schiebe, F.R., 1995. Large-area mapping of soil moisture using the ESTAR passive microwave radiometer in Washita92. Remote Sens. Environ. 54, 27–37. Jadoon, K.Z., Weihermüller, L., Scharnagl, B., Kowalsky, M.B., Bechtold, M., Hubbard, S.S., Vereecken, H., Lambot, S., 2012. Estimation of soil hydraulic parameters in the field by integrated hydrogeophysical inversion of time-lapse groundpenetrating radar data. Vadose Zone J, http://dx.doi.org/10.2136/vzj2011.0177. Jhorar, R.K., Bastiaanssen, W.G.M., Feddes, R.A., Van Dam, J.C., 2002. Inversely estimating soil hydraulic functions using evapotranspiration fluxes. J. Hydrol. 258 (1), 198–213. Kool, J.B., Parker, J.C., 1988. Analysis of the inverse problem for transient unsaturated flow. Water Resour. Res. 24, 817–830. Kool, J., Parker, J., Van Genuchten, M.T., 1985. Determining soil hydraulic properties from one-step outflow experiments by parameter estimation: I. Theory and numerical studies. Soil Sci. Soc. Am. J. 49 (6), 1348–1354. Kool, J.B., Parker, J.C., van Genuchten, M.T., 1987. Parameter estimation for unsaturated flow and transport models—a review. J. Hydrol. 91, 255–293. Kowalsky, M.B., Finsterle, S., Peterson, J., Hubbard, S., Rubin, Y., Majer, E., Ward, A., Gee, G., 2005. Estimation of field-scale soil hydraulic and dielectric parameters through joint inversion of GPR and hydrological data. Water Resour. Res. 41, W11425. Lazarovitch, N., Ben-Galb, A., Šimu˚nekc, J., Shani, U., 2007. Uniqueness of soil hydraulic parameters determined by a combined wooding inverse approach. Soil Sci. Soc. Am. J. 71 (3), 860–865. Li, C., Ren, L., 2011. Estimation of unsaturated soil hydraulic parameters using the ensemble Kalman filter. Vadose Zone J. http://dx.doi.org/10.2136/vzj2010.0159. Li, B., Yeh, T.C.J., 1999. Cokriging estimation of the conductivity field under variably saturated flow conditions. Water Resour. Res. 35 (12), 3663–3674. Li, L., Zhou, H., Gómez-Hernández, J.J., Hendricks Franssen, H.J., 2012. Jointly mapping hydraulic conductivity and porosity by assimilating concentration data via ensemble Kalman filter. J. Hydrol. 428, 152–169. Liu, S., Yeh, T.C.J., 2004. An integrative approach for monitoring water movement in the vadose zone. Vadose Zone J. 3, 681–692. Liu, G., Chen, Y., Zhang, D., 2008. Investigation of flow and transport processes at the MADE site using ensemble Kalman filter. Adv. Water Resour. 31 (7), 975–986. Mao, D., Yeh, T.C.J., Wan, L., Hsu, K.C., Lee, C.H., Wen, J.C., 2013. Necessary conditions for inverse modeling of flow through variably saturated porous media. Adv. Water Resour. http://dx.doi.org/10.1016/j.advwatres.2012.08.001. McPherson, B.J.O.L., 1998. Aquifer hydraulics: a comprehensive guide to hydrogeologic data analysis. Eos, Trans. Am. Geophys. Union 79 (pp. 635–635). Mishra, S., Parker, J.C., 1989. Parameter estimation for coupled unsaturated flow and transport. Water Resour. Res. 25, 385–396. Montzka, C., Moradkhani, H., Weihermüller, L., Hendricks Franssen, H.J., Canty, M., Vereecken, H., 2011. Hydraulic parameter estimation by remotely-sensed top soil moisture observations with the particle filter systems. J. Hydrol. 399, 410– 421. Moradkhani, H., 2005. Uncertainty assessment of hydrologic model states and parameters: Sequential data assimilation using the particle filter. Water Resour. Res. 41 (5), 1–17. Moradkhani, H., Sorooshian, S., Gupta, H.V., Houser, P., 2005. Dual state-parameter estimation of hydrological models using ensemble Kalman filter. Adv. Water Resour. 28 (2), 135–147. Nutzmann, G., Thiele, M., Maciejewski, S., Joswig, K., 1998. Inverse modelling techniques for determining hydraulic properties of coarse-textured porous media by transient outflow methods. Adv. Water Resour. 22, 273–284. Panzeri, M., Riva, M., Guadagnini, A., Neuman, S.P., 2013. Data assimilation and parameter estimation via ensemble Kalman filter coupled with stochastic moment equations of transient groundwater flow. Water Resour. Res. 49 (3), 1334–1344.

561

Qin, J., Liang, S.L., Yang, K., Kaihotsu, I., Liu, R.G., Koike, T., 2009. Simultaneous estimation of both soil moisture and model parameters using particle filtering method through the assimilation of microwave signal. J. Geophys. Res. 114, D15103. Ritter, A., Hupet, F., Muñoz-Carpena, R., Lambot, S., Vanclooster, M., 2003. Using inverse methods for estimating soil hydraulic properties from field data as an alternative to direct methods. Agr. Water Manage. 59, 77–96. Ross, P.J., 1993. A method of deriving soil hydraulic properties from field water contents for application in water-balance studies. J. Hydrol. 144, 143–153. Scholer, M., Irving, J., Looms, M.C., Nielsen, L., Holliger, K., 2013. Examining the information content of time-lapse crosshole GPR data collected under different infiltration conditions to estimate unsaturated soil hydraulic properties. Adv. Water Resour. 54, 38–56. Schöniger1, A., Nowak1, W., HendricksFranssen, H.J., 2012. Parameter estimation by ensemble Kalman filters with transformed data: approach and application to hydraulic tomography. Water Resour. Res. 48, W04502. Shi, L., Zeng, L., Zhang, D., Yang, J., 2012. Multiscale-finite-element-based ensemble Kalman filter for large-scale groundwater flow. J. Hydrol. 468–469 (10), 22–34. Šimu˚nek, J., van Genuchten, M.T., 1996. Estimating unsaturated soil hydraulic properties from tension disc infiltrometer data by numerical inversion. Water Resour. Res. 32 (9), 2683–2696. Šimu˚nek, J., Wendroth, O., van Genuchten, M.T., 1998. Parameter estimation analysis of the evaporation method for determining soil hydraulic properties. Soil Sci. Soc. Am. J. 62, 894–905. Simunek, J., Sejna, M., van Genuchten, M.T., 2005. The HYDRUS-1D software package for simulating the one-dimensional movement of water, heat, and multiple solutes in variably-saturated media. University of California, Riverside, Research reports: p. 240. Song, X., Shi, L., Ye, M., Yang, J., 2014. Numerical comparison of iterative ensemble Kalman filters for unsaturated flow inverse modeling. Vadose Zone J. http:// dx.doi.org/10.2136/ vzj2013.05.0083v. Toorman, A.F., Wierenga, P.J., Hills, R.G., 1992. Parameter estimation of hydraulic properties from one-step outflow data. Water Resour. Res. 28, 3021–3028. van Dam, J.C., Stricker, J.N.M., Droogers, P., 1992. Inverse method for determining soil hydraulic functions from one-step outflow experiments. Soil Sci. Soc. Am. J. 56, 1042–1050. van Dam, J.C., Stricker, J.N.M., Droogers, P., 1994. Inverse method to determine soil hydraulic functions from multistep outflow experiments. Soil Sci. Soc. Am. J. 58, 646–652. Vereecken, H., Huisman, J., Bogena, H., Vanderborght, J., Vrugt, J., Hopmans, J.W., 2008. On the value of soil moisture measurements in vadose zone hydrology: a review. Water Resour. Res. 44. http://dx.doi.org/10.1029/2008WR006829. Walker, J.P., Houser, P.R., 2004. Requirements of a global near surface soil moisture satellite mission: accuracy, repeat time, and spatial resolution. Adv. Water Resour. 27, 785–801. Walker, J.P., Willgoose, G.R., Kalma, J.D., 2001. One-dimensional soil moisture profile retrieval by assimilation of near-surface observations: a comparison of retrieval algorithms. Adv. Water Resour. 24, 631–650. Wang, W., Neuman, S.P., Yao, T., Wierenga, P.J., 2003. Simulation of large-scale field infiltration experiments using a hierarchy of models based on public, Generic, and Site Data. Vadose Zone J. 2, 297–312. Wöhling, T., Vrugt, J., Barkle, G.F., Woehling, T., 2008. Comparison of three multiobjective optimization algorithms for inverse modeling of vadose zone hydraulic properties. Soil Sci. Soc. Am. J. 72, 305–319. Yeh, T.C.J., Šimu˚nek, J., 2002. Stochastic fusion of information for characterizing and monitoring the vadose zone. Vadose Zone J. 1, 207–221. Yeh, T.C.J., Zhang, J., 1996. A geostatistical inverse method for variably saturated flow in the vadose zone. Water Resour. Res. 32, 2757–2766. Yeh, T.C.J., Jin, M., Hanna, S., 1996. An iterative stochastic inverse method: Conditional effective transmissivity and hydraulic head fields. Water Resour. Res. 32 (1), 85–92. Yeh, T.C.J., Liu, S., Glass, R.J., Baker, K., Brainard, J.R., Alumbaugh, D., LaBrecque, D., 2002. A geostatistically based inverse model for electrical resistivity surveys and its applications to vadose zone hydrology. Water Resour. Res. 38 (12), 14-1– 14-13. Yeh, T.C.J., Ye, M., Khaleel, R., 2005. Estimation of effective unsaturated hydraulic conductivity tensor using spatial moments of observed moisture plume. Water Resour. Res. 41, W03014. http://dx.doi.org/10.1029/2004WR003736. Zeng, L., Shi, L., Zhang, D., Wu, L., 2012. A sparse grid based Bayesian method for contaminant source identification. Adv. Water Resour. 37, 1–9. Zha, Y., Shi, L., Ye, M., Yang, J., 2013. A generalized Ross method for two- and threedimensional variably saturated flow. Adv. Water Resour. 54, 67–77. Zhang, J., Yeh, T.C.J., 1997. An iterative geostatistical inverse method for steady flow in the vadose zone. Water Resour. Res. 33 (1), 63–71. Zhou, H., Gómez-Hernández, J.J., Hendricks Franssen, H.J., Li, L., 2011. An approach to handling non-Gaussianity of parameters and state variables in ensemble Kalman filtering. Adv. Water Resour. 34, 844–864. Zhou, H., Gómez-Hernández, J.J., Li, L., 2014. Inverse methods in hydrogeology: Evolution and recent trend. Adv. Water Resour. 63, 22–37. Zijlstra, J., Dane, J.H., 1996. Identification of hydraulic parameters in layered soils based on a quasi-Newton method. J. Hydrol. 181, 233–250.