Compaction parameter estimation using surface movement data in Southern Flevoland

The Southern part of the Flevopolder has shown considerable subsidence since its reclamation in 1967. We have set up an integrated method to use subsidence data, water level data and forward models for compaction, oxidation and the resulting subsidence to estimate the driving parameters. Our procedure, an Ensemble Smoother with Multiple Data Assimilation, is very fast and gives insight into the variability of the estimated parameters and the correlations between them. We used two forward models: the Koppejan model and the Bjerrum model. In first instance, the Bjerrum model seems to perform better than the Koppejan model. This must, however, be corroborated with more elaborate parameter estimation exercises in which in particular the water level development is taken into account.


Introduction
The Southern part of the Flevopolder (the Netherlands) was reclaimed in 1967-1968.After lowering the groundwater level, the sediments compacted, resulting in surface subsidence.The compaction is a combination of elastic and viscoplastic compaction of clay and peat, and of peat oxidation.A good prediction of future compaction and an extrapolation to other areas naturally requires knowledge of the subsurface lithology, the compaction parameters and the phreatic water level development.The present paper describes an exercise of parameter estimation using a limited number of borehole data in the Southern Flevopolder area.An Ensemble Smoother with Multiple Data Assimilation was employed to this end.

Available data
The Flevopolder is a large polder in the Netherlands, which was reclaimed from the IJsselmeer in 1967 (Fig. 1).The surface has an elevation of about 3 m below sea level.The composition of the subsurface had been mapped on a number of places prior to the reclamation.It is a variable layered structure which we categorized into four generic types: clay, humic clay, peat, and sand.For the present study, 10 measure-ment locations were selected, at which the surface level had been monitored yearly between 1967 and 1993 by geodetic levelling.
Data about the phreatic groundwater level are scarce.One of the few measurements covering the complete period is represented in Fig. 2. It shows the fast lowering of the water level at reclamation of the polder, followed by a gradual decrease during the next 7 years.This can be explained by the gradual adjustment of the phreatic level at the measurement location to the drainage level in the ditches.

Forward model
For the subsidence prediction we employed a combination of compaction in the clay, the humic clay and the peat, and oxidation in the humic clay and peat layers above the phreatic level.We used models for primary and secondary settlement for stresses larger than the pre-consolidation stress: the Koppejan model and the Bjerrum model, both with an instantaneous lowering of the groundwater level to a fixed value below the original surface level (Koppejan, 1948;Bjerrum, 1967).In combination with peat oxidation, the two models then predict the following compaction strain and associated Published by Copernicus Publications on behalf of the International Association of Hydrological Sciences.2.0 ± 0.5 1.40 ± 0.06 2.0 ± 0.5 1.39 ± 0.05 2.0 ± 0.5 1.45 ± 0.12 subsidence as a function of time: Here, σ and σ 0 are the actual and original effective vertical stresses; C p and C s are the primary and secondary compression coefficients above pre-consolidation pressure in the Koppejan model; v ox is the oxidation rate; CR and C α are the virgin compressibility (above pre-consolidation pressure) and the creep parameter or coefficient of secondary compression in the Bjerrum model.

Inverse model
The present study aimed at the assessment of the uncertainty of the model parameters in order to improve the his-tory match of the ground level displacement and the reliability of the associated predictions.We considered as uncertain parameters the primary and secondary compression coefficients, the oxidation rates and the phreatic water levels.We assumed the measured lithology at every location as fixed.
Of the four lithologies identified, we assumed the sand to be incompressible.Both models were thus left with 18 uncertain parameters: primary and secondary creep parameters for clay, for humic clay, and for peat (6 parameters); oxidation rates for humic clay and for peat (2 parameters); and groundwater levels for each location (10 parameters).
We define the vector m as the collection of adjustable model parameters; the vector d as the collection of data (surface level measurements vs. time) and the functional G(m), working on the model parameters, as the forward model.The inverse problem is then formulated as the task of estimating the vector m for which G (m) approaches the data vector d best.With additional information present in the form of a prior model (m 0 ) and covariance matrices of the measurements (C d ) and of the prior model (C m ), the conventional least-squares solution is obtained by maximizing the objective function J given by Tarantola (2005) (or by minimizing the exponent in the expression, − log(J )): We used an ensemble approach in which the mean and the covariance of the model vector m 0 are mapped on an ensemble of N e vectors.Different approaches exist to obtain an estimate of the model vector.A global update of the model using all available data can be achieved in a single step; this procedure is called an Ensemble Smoother (Emerick and Reynolds, 2012) -it is, however, suboptimal for non-linear problems as ours.When both data and parameters are assimilated to newly incoming data in subsequent time steps, the procedure is called a filter.EnKF is an example of a filter, but we did not use EnKF because all data were available at all times.Emerick and Reynolds (2013) and Tavakoli et al. (2013) state that the best approach for nonlinear systems such as ours is to use an Ensemble Smoother with Multiple Data Assimilation (ES-MDA).In ES-MDA, the ensemble smoother is applied iteratively multiple times.
In this way, correlations between parameters that result from the smoother are retained in subsequent steps.The advantage of ES-MDA above EnKF is that the data are used as given and that the procedure is computationally less demanding.

Results
The results of the estimation exercise with both models are presented in Table 1.For the Koppejan model we report re-sults that were achieved after starting from two initial distributions -the first one with closer bounds than the second.The results show that the main adjustments for the Koppejan model are in the oxidation rates and in the phreatic water levels.These parameters have the largest impact on the surface movement rates.The adjustments of the oxidation rates are larger for the less constrained case, as should be expected.For the phreatic water levels, the adjustments are similar.There is quite a large variability in the resulting estimate of the phreatic water level among the different locations -this is related to the types of lithology present around the phreatic water level: peat oxidation has the largest influence on subsidence and boreholes with peat will result in a better constrained water level because the sensitivity to it is larger.
For the Bjerrum model, the largest effect is also related to the peat oxidation rates and the water levels, but the compaction parameters are also significantly constrained.For these parameters, the values of CR and C α are anti-correlated -an increase in one of them can be partially compensated by a decrease in the other.
Figure 3 visualizes the data, the priors and the inversion results.While a large degree of scatter is present for the prior distribution of parameters, the smoother succeeds in constraining this uncertainty and obtaining a reasonable fit to the data.The assimilation is better when using the Bjerrum model.

Discussion
We wish to address some issues here that are striking from the results as presented above.The first is the relative insensitivity to the compaction parameters in the Koppejan model; the second is the appearance of unphysical negative numbers for the immediate compaction parameter CR in the Bjerrum model; the third is the high oxidation rates resulting for the humic clay and for the peat -the former even being the largest.In our opinion, these issues are related to a fourth issue: the fact that the match of the curves using estimated parameters is suboptimal in the first 6 or 7 years.We feel there is a strong relation with the fact that the parameter estimation with the present dataset was accomplished under the assumption of an immediate drop of the water level to a further constant value.This is a simplification that must be addressed in a follow-up study.Although only limited data on the water level are available, Fig. 2 shows that this assumption is not generally applicable.A first-order approximation could be to distribute a 1 m level drop to the final value over the first 7 years after reclamation.We are currently building such an extension into our models.

Conclusions
We estimated compaction parameters for clay, humic clay, and peat, oxidation rates for peat and for humic clay, proc-iahs.net/372/183/2015/Proc.IAHS, 372, 183-187, 2015  and phreatic water levels for 10 locations in the Southern Flevopolder, using surface movement parameters between the moment of reclamation in 1967 and 1993.The procedure, using an Ensemble Smoother with Multiple Data Assimilation, is very fast and gives insight into the variability of the estimated parameters and the correlations between them.
In first instance, the Bjerrum model seems to perform better than the Koppejan model.This must, however, be corroborated with more elaborate parameter estimation exercises in which in particular the water level development is taken into account.

Figure 1 .
Figure 1.Map of the SE Flevopolder with the selected locations.

Figure 2 .
Figure 2. Hydraulic head in well B26E00030001, located on the boundary between the previously reclaimed Eastern Flevoland polder and the subsequently reclaimed Southern Flevoland polder.The reclamation of the Southern Flevoland polder in 1967 is followed by a gradual decrease of the level during approximately 7 years; after which the level remains virtually constant until 1995.

Figure 3 .
Figure 3. Surface movement data and predictions before and after data assimilation.Subsequent series are shifted downward by 0.5 m.Colored circles: data; colored dashed lines: average prior estimate; colored solid lines: average assimilated estimate.Light grey lines: prior ensemble predictions.Dark grey lines: posterior ensemble predictions.Case 2, with less constrained prior gives larger prior spread and better estimates than Case 1.