Een model voor de simulatie van het fysische rijpingsproces van gronden in de IJsselmeerpolders

<strong>SUMMARY</strong><strong>Introduction</strong>Subaqueous sediments, rich in clay and organic matter have a high watercontent and are unaerated and almost impermeable just after emergence. Such a sediment is inaccessible due to the low bearing capacity. Plant growth is restricted to a pioneer vegetation of plants, which are able to provide the oxygen demand of their own roots. Due to evapotranspiration in periods with a rainfall deficit the water content of the sediment decreases irreversibly. By capillary forces the soil particles are drawn to each other, so the bulk density increases. As a result, the surface will subside and crack formation starts. Water transport to open-field drains and tile-drains is then possible via the cracks. This dewatering-process is known as 'ripening'. Ripening can be subdivided into a physical, a chemical and a microbiological part. The physical part of the process, which can be considered as the motor of the entire process, is the topic of this study.Knowledge about this process is necessary to be able to predict the qualities of not yet reclaimed sediments and their rate of subsidence. Since the start of the IJsselmeerpolder project, research has been carried out to increase this knowledge. In the first period emphasis was given to research concerning dewatering-measures that had to be taken to obtain a ripened soil as fast as possible. Also much research was focused on the subsidence-rates. A first attempt for a quantitative integrated approach of the physical ripening process was given by SEGEREN (1966).To improve this quantitative insight, including the process effecting factors a numerical simulation model has been developed, based on a conceptual approach.<strong>Description of the process of physical soil ripening</strong>The soil profile in the IJsselmeerpolders can be roughly characterised as a layer of Holocene, mostly rich in clay and organic matter, with a thickness from nearly nil to over 7 m, underlain by Pleistocene, mostly sandy deposits (figure 3 and 4). The clay fraction is composed of app. <em>60%</em> illite, 20% smectite and 20% kaolinite. The sediments are rich in calciumcarbonate, so they remain basic after reclamation.Just after emergence the water contents of the sediments are very high. At that moment app. 2.2 grams of water are bound by 1 gram of clay. In a completely ripened topsoil just 0.3-0.4 grammes of water are found per gram clay. The water content has been found to be linearly related to the clay and organic matter content during all the ripening stages (Eq. 2.1. ZUUR; 1958).The process of physical soil ripening can be approached by using waterbalances. The progress in ripening is positively influenced by periods with a high rainfall deficit and by discharge of water in periods with a rainfall surplus. During the ripening process the pF-curve, the permeability and the thickness of a layer change considerably. The formation of a waterbalance for a ripening soil is complicated by these changes.Under the climatic conditions in the Netherlands in an average year an evaporation surplus is found in summer and a rainfall excess in winter (figure 10). For the evaporation surplus in summer, which determines the progress in ripening, an average value of 163 mm is found, if the difference between the evaporation of a free water surface and the rainfall is taken into account. In very dry summers a value of 500 mm can be reached.In dry periods water is drawn from the soil by evapotranspiration. The soil layers close to the surface lose their water the most and the first, even if the soil is kept bare.Due to the difference in elevation between the water levels in the polders and in the surrounding lakes, seepage is found locally. Seepage has a negative effect on ripening. The amount of seepage is difficult to determine. Especially the non-discharged part of the seepage, which is responsible for a retardation of the ripening process, is difficult to measure.The permeability of an unripened soil is found to be app. 1.7.10 ^-4m.d ^-1. This value is found for different types of determination. No relation could be found between the permeability and the pore space or the clay and organic matter content. The permeability of a ripening soil increases by crack-formation up to values above 100 m.d ^-1. This 'crack-permeability' is completely determined by the existence of the cracks. The permeability of the compact soil between the cracks remains at the mentioned low value.The depth of the cracks increase during the process of ripening. For that reason the groundwater levels found for a certain amount of discharge will be lowered during ripening. So the progress in ripening can be illustrated by using so called Q-y figures of successive years. In these figures the relation between the discharge Q and the groundwater level y is given (figure 15 and 16).In the IJsselmeerpolders ripening is stimulated by the installation of dewatering-methods, such as ditches, open-field drains and subsurface tile drains and by the exploitation <strong/> of the soils with arable crops. The usual scheme of crop rotation in the IJsselmeerpoiders is: oil seed rape, winter wheat and spring barley in the first three years, followed by oil seed rape or oats in the fourth year. By using crops, which can be harvested in summer, damage to the soil structure is prevented.<strong>Problems in the determination of relevant parameters</strong>Research on the process of physical soil ripening is hampered by several causes. Even on the often rather uniform soil layers in the IJsselmeerpolders an important variation in the results of laboratory measurements is found due to the non-uniformity of the samples. Due to unwanted but inevitable ripening during experiments, some parameters in an unripened soil can not be measured in the right way. The very low permeability of the soil makes it almost impossible to obtain correct values for groundwater levels. An area with unripened soils is difficult to enter, due to the low bearing capacity and to the lack of fixed points, like roads, church towers etc. For the above mentioned reasons it is very difficult to obtain sufficient data over a long period.<strong>Numerical simulation of physical soil ripening</strong>The development of the simulation model has proceeded by the formulation of a number of points of departure. The most important ones are:- insight in this process can only be achieved by a conceptual model- the soil moisture suction ψmust be used as the 'master variable', not only defining the storage and the permeabilities but also the rate of compaction- physical soil ripening can be simulated with a one-dimensional model- a cracked soil must be divided into cracks and compact soil; the parameters of the soil must be related to the compact soil.In the model, the soil profile is divided into layers, thin (0.5-2.0 cm) layers in the topsoil and thicker (up to 5 cm) layers in the subsoil. Every layer is defined by its clay and organic matter content and its bulk density. The content of solid parts in each layer remains the same during the simulation, so the thickness is variable.To allow a simulation of app. 10 years with acceptable costs, a new calculation method has been developed to find the ψ-profile on every time step. The calculation of this ψ- profile and all the ψ-dependant parameters is performed in the following way.A first estimate is chosen for ψin the uppermost layer. Based on this estimated value it is calculated which amount is delivered in this first layer. The delivered amount depends on the suction difference and the possible compaction. The difference between the needed and the delivered amount is transported from or to the second layer. This flux defines the value of ψin this layer, using Darcy's law. This procedure is repeated for all layers. At the last layer the integrated amounts of needed and delivered water must be equal, except for a small inaccuracy, which can be chosen. If the amounts are not equal, the total calculation is repeated, using a new estimate for ψin the first layer. The process of iteration is repeated until the given accuracy is reached. It is proved that the calculated solution for the ψ-profile is the only possible solution.The potential evaporation is calculated by multiplying the open water evaporation by a reduction factor, depending on the soil cover (Eq 4.19). When the soil suction in the topsoil reaches high values, the actual evaporation can have a smaller value than the potential one, due to a decreasing permeability.The potential transpiration is calculated by multiplying the difference between the open water evaporation and the actual evaporation by a reduction factor. Values for this reduction factor, related to the stage of development and the kind of crop for the months of the year, are given (table g). The withdrawal by the plant roots (the sink- term) is calculated, using a method given by HOOGLAND et al. (1981). This method is based on a preference for extraction in the upper layers and on the assumption that a certain maximum amount can be withdrawn from a volume unit of soil in a unit of time. The withdrawal is reduced if the soil suction exceeds values of 500 cm. A method has been developed to calculate pF-curves, using a data-base of pF-curves, determined in the laboratory. The pF-curves of clay soils (>8% clay) are defined by the bulk density and the clay and organic matter content. The pF-curves of sandy soils are only defined by the coarseness of the sand. The advantage of these calculation methods is that no pF-curves have to be delivered as input, an important advantage for ripening soils with changing pF-curves It is proved that a good agreement exists between the calculated and determinated pF-curves (tables 10 and 11).Approximation methods have also been developed for the relation between the soil suction ψand the permeability k. For the sandy soils the relations are obtained using the BROOKS and COREY-method. The saturated permeability is calculated for these soils with the formula given by KOZENY (Eq. 4.39). For the clay soils it is assumed that the saturated permeability remains on a value of 1.7 <strong>^.</strong> 10 ^-4m.d ^-1, and the k-ψrelationship can be given as k = aψ ^-n. A linear relation was found between log a and the exponent n (Eq. 4.41). The value for the exponent n appeared to be related to the clay-content (Eq. 4.42). The approximation methods for both the pF-curves and the permeability are only valid for the compact soil.A new formula has been developed to calculate the rate of compaction (Eq. 4.52). This formula is based on the assumption of an existing relation between the logarithm of the grain stress in the soil and the pore space. Due to the compaction the bulk density will increase. As a result the surface will subside and cracks will be formed. A relation is given to calculate the distribution of the compaction over subsidence and crack-formation (Eq. 4.63). It is indicated that only subsidence can appear as a result of compaction if a soft layer is loaded by the upper layers. In this case developing cracks will be closed, due to this load and the weakness of the soil. A relation between the waterfactor n of the first non-cracked layer and the overburden pressure is given (figure 39).To save costs a simplified method has been developed to calculate the rewetting of the profile. An exact calculation of ψ-profiles during this rewetting-process is not essential for a ripening model.<strong>Validation and use of the model</strong>The developed simulation model is calibrated and tested, using measured data from research areas for subsidence in Zuidelijk Flevoland. Calibration was necessary while the values of two parameters in formula 4.55 could not be found in another way. Although more parameters had to be estimated, the calibration is restricted to these two parameters.Just a few research areas in Flevoland were suitable for testing. Even for these areas, on which many parameters had been measured, several input parameters had to be estimated. The amounts of rainfall and evaporation had not been recorded on the spot and had to be calculated using data of surrounding meteorological stations (Eq. 5.1).The results of the simulations are compared with the data, measured in the field. Most attention is given to a comparison between measured and calculated values for the depth of the cracks and the subsidence. The agreement between the measured and the calculated values is rather good (a.o. figure 52 and 53). A sensitivity analysis has been carried out to investigate the influence on the results of:- the length of the time step- the saturated and unsaturated permeabilities- the maximum rate of water withdrawal by plant roots- the constant K ₂ in the compaction formula 4.55In particular a change in permeability proved to have an impressive influence on the results.Seven variants have been simulated to illustrate the application-possibilities of the developed model. A good ripened soil can be obtained if the soil is kept bare and the extraction of water occurs only by evaporation of the bare soil. A permanent grass vegetation instead of the normal crop rotation shows a slower progress in ripening of the soil. If, next to it, the dewatering depth is also restricted to 0.3 m, the progress in ripening is even slower.Nevertheless the groundwater depths in this variant are found in summer on a low level (figure 57). A very good ripened soil can be obtained if the soil is covered by a forest. Rather high groundwater levels are found in summer if any seepage is present (table 25). Seepage proved to have a negative influence on the progress of ripening. Progress is also diminished if a 4 cm thick sand layer, unpenetratable for plant roots, is assumed on a depth of 40 to 44 cm. The layers below this sand layer remain wet due to a decreased capillary rise, so they subside less (figure 59). A very fast ripening progress is found under dry climatic conditions. A variant has been simulated, using climatic data as found in the Rumanian Danube-delta.<strong>Limitations of the model</strong>A simulation model is defined as a simplified representation of a process in reality. A number of factors, influencing the process in reality, which are not taken into account in the model, are elucidated. The effects of frost, tillage activities and clay transport are not included in the model. It is assumed that a good simulation of the process of physical soil ripening is possible with a model in which these factors are not considered.All relations between parameters are based on research results as found in the IJsselmeerpolders As a result the application of the model in its present form is restricted to that area. Especially the mineralogical composition of the clay fraction will effect the process. Nevertheless the approach of the problem can be used as an example when developing a simulation model for physical soil ripening in other parts of the world.