# Numerical Modelling of Clay Seal Maturation in Deep Boreholes with Nuclear Waste

Authors: Ting Yang; Hans Mattsson; Roland Pusch; Jan Laue; Sven Knutsson; Xiaodong Liu
Journal: Advances in Materials Science and Engineering
 (2020)

---

## Abstract

One of the major challenges of high-level nuclear waste (HLW) isolation in deep boreholes is to anticipate the maturation behaviour of swelling clay when the waste, surrounded by dense clay encased in perforated tubes, is submerged into the borehole mud. The ultimate homogeneity of this clay seal acting as a barrier is expected to stabilize the borehole and to prevent possible leakage of radioactivity. In this study, a numerical model for predicting the maturation of the clay barrier has been developed. In the model, the water transport is controlled by the differences in the suction potential and the permeability. The model is able to simulate the maturation process, both the expeditious water transport and the clay migration into the surrounding mud, from beginning to end. Results from laboratory tests of the clay maturation were compared with the predictions made by the model. They are in good agreement, but refinement is proposed by taking the impact of the tube perforation on the maturation rate into more consideration. The proposed numerical model will also be useful in selecting a suitable design for clay barriers in HLW boreholes. Different combinations of clay dimensions and initial densities of mud and dense clay can thus be studied to determine the final homogenization and the end densities.


---

## Body

## 1. Introduction

One concept to store high-level radioactive waste (HLW) is to place it in canisters in the bottom regions of a very deep borehole (VDH); see Figure1. To prevent leakage of radioactivity, the canisters are surrounded by a dense swelling clay. Dense swelling clay is also placed in blocks above the canisters in the boreholes. The canisters and the swelling clay are submerged in the boreholes into a clay mud [1]. The dense swelling clay and clay mud are called clay seals. Maturation of a clay seal is the process in which the dense clay comes in contact with water and swells against the borehole walls, creating a barrier of low permeability to prevent the possible leakage of radionuclides.Figure 1
Schematic illustration of VDH concept. The package shown on the right consists of a supercontainer with a clay-embedded canister containing spent fuel [6, 7].This study focuses on the numerical modelling of the maturation of clay seals in VDHs for HLW disposal. To the authors’ knowledge, this type of model has still not been proposed by anyone until now. However, in a paper by Börgesson et al. [2] a laboratory test on bentonite clay was modelled with focus on radially swelling using a commercial finite element software. Good agreement between experimental and numerical result was obtained with a constitutive model the researchers had developed themselves. The maturation of clay seals was not modelled, but since radial swelling is an important part, this strategy might be successful in modelling the maturation of clay seals too. The numerical model proposed in this study simulates the water transport and the soil migration process during the swelling from the initial heterogeneous state of mud and dense clay to the ultimate homogeneous condition of the entire clay seal. The model could be used to determine the duration of the maturation process for different designs of the storage concept. If the maturation process goes too fast, it will be problematic to install the canisters and the dense clay into the boreholes. If the process goes too slow or stops before a proper barrier is created, leakage of radioactivity is not prevented as intended.To verify the numerical model, laboratory tests of the maturation process were performed and the test results compared with the model predictions. The maturation process in real boreholes is also simulated in the study, as an example of how the numerical model can be utilized in real design.


## 2. Clay-Based Barriers in Very Deep Boreholes with Spent Reactor Fuel

Very deep boreholes repository concepts have gained interest in countries like the US and the UK because of lower costs, faster construction, and simpler and safer placement of high-level radioactive waste [3]. Sweden has 10 nuclear reactors, each with a 1,000 MW capacity, which need to dispose of 40 years’ worth of spent fuel which is currently on intermediate storage. A successive placement of canisters could be made in 20 to 25 years, utilizing 4 to 12 multiple borehole sites with depths between 2.5 and 4 km, and a spacing of 200 m between individual holes to prevent thermal overlap. Further details on the concept can be found in, e.g., the paper by Sandstedt et al. [4].Figure1 illustrates the barriers—the metal canisters and the surrounding dense smectite-rich clay—that prevent the dissemination of radionuclides in the “deployment zone” [5]. Above the canisters, the hole is sealed by placing perforated supercontainers filled solely with dense smectite clay. The containers will be submerged in smectite clay mud. The containers are necessary in order to get the clay in place, and the dense clay will be in contact with the mud when the clay swells through the holes in the supercontainers.In this study, the numerical model is utilized on the part of the borehole where the supercontainers are filled with only dense swelling clay. The heat production caused by radioactive decay in the zone with canisters is not accounted for in the numerical model. The influence of the perforated supercontainers (tubes in laboratory tests) on the maturation process is not included in the numerical model, but the influence was examined by deriving a relationship between the maturation times with and without a tube, based on laboratory tests.


## 3. Laboratory Test of Maturation

Laboratory tests of the maturation process were performed to simulate the in situ situation of a VDH. The test setup is presented in Figure2(c). The components in the setup, from the outside to the inside, are as follows: a 9 mm thick pervious acrylate-stabilized sand tube (sand filter), the 6.5 mm thick soft clay mud, a 2 mm thick steel tube with 50% perforation (representing a supercontainer), and dense swelling clay with the radial thickness 43 mm as sealing material. The 124 mm high steel tubes were closed at the ends by steel plates. The entire setup was submerged in water in order to keep the sand filter continuously saturated during the tests.Figure 2
Instrument in the laboratory test. (a) Setup with clay seals after maturation, with the top lid removed; (b) 96-hour laboratory experiment with clay plugs in perforated steel tubes surrounded by mud; (c) schematic illustration of test setup with clay seals (soft clay mud and dense clay blocks).
(a)(b)(c)The expansive Holmehus clay, produced in Southern Denmark, was the sealing material used in the tests [7]. It is an illite/smectite mixed-layer clay with approximately 60% montmorillonite. The behaviour of swelling clays can be studied in, e.g., the paper by Mitchell and Soga [8]. The clay specimens were manufactured by mixing the air-dry granulate with dry water [9] to become fully water-saturated at uniaxial compaction to a dry density of 1550 kg/m3 [7]. The dense clay blocks in the perforated tubes, surrounded by soft clay mud, were contained in the filter for maturation in 6, 12, 24, 48, and 96 hours and one year. Figure 2(b) shows that, after a certain time of maturation, the dense clay started to swell through the perforation, forming small clay columns that displaced mud. Due to the difference in suction between the mud and the dense clay, and also within the different parts of the dense clay, the mud hardened while the dense clay softened during the maturation process. Figure 2(a) illustrates the laboratory setup at an advanced stage of maturation. After each maturation test, the densified mud and the softened dense clay were examined with respect to the water content and the dry density at different radial distances from the middle, by cutting out and examining hollow cylindrical pieces of the specimen [7].


## 4. Numerical Model of Maturation of Clay Seals

### 4.1. The Proposed Model

In the case of maturation of clay seals in a borehole, the key assumptions and criteria in the model are as follows:(1)
Darcy’s law is valid between each element in the numerical algorithm(2)
Both dense clay and mud are initially homogeneous(3)
The total volume inside the hollow cylinder filter, containing fully saturated dense clay and mud, is constant during the entire maturation process(4)
Expansion/contraction and flow take place only in the radial direction(5)
There is no temperature dependence, e.g., from heat production caused by radioactive decayDarcy’s law is expected to give an accurate representation of the flow within a porous medium for small discharge velocities [10]. Water flow is slow in the maturation process, so it is assumed that Darcy’s law is applicable on the elements in the numerical algorithm. Darcy’s law states that q=Aki, where the fluid flux q is related to the cross-sectional area A, the hydraulic conductivity k, and the hydraulic gradient i. The permeated fluid volume Q in time t0 is given by(1)Q=Akit0.With the height of the cylinder specimenH and the radius R, the permeated fluid volume becomes(2)Q=2πRHkdhdRt0,where h is the water head and dh/dR is the hydraulic gradient i.Rearranging and integrating equation (2) gives(3)Q⋅∫R1R21RdR=∫h1h22πkt0Hdh.Thus,(4)Q=2πkHlnR2/R1h2−h1·t0,where R1 and R2 represent the flow path. h2−h1 is the hydraulic head difference causing water transport.The hydraulic gradient and thereby the hydration/dehydration of the dense clay is a function of the position of the clay elements and of the maturation time. Therefore, the dense clay is considered to consist of a sequence of elements (Figure3(a)), subjected to one-dimensional hydration or dehydration, implying only radial flow and clay particle migration. Because of the rotational symmetry, only the radial plane of the clay sample is viewed (Figure 3(b)). If a unit sample with height ΔZ and initial radial length R is divided into N elements (N is a natural number), the initial size of each element and the distance from its centre to the symmetry axis are Δrn and rnn=1,2,3,…,N, respectively. For the mud, the corresponding parameters are Δrmud and rmud, respectively.Figure 3
Schematic view of specimen and maturation. (a) Clay specimen surrounded by water-rich mud, absorbing water in radial direction. (b) Division of dense clay into elements for the numerical modelling.
(a)(b)Theoretically, in the clay seals, water is absorbed from the mud by the nearest clay element, followed by the water movement in the series of elements taking place from the centre of an element to the centre of the adjacent one. Therefore, in the model, the values of permeability and suction should be representative from the centre of the clay element to the centre of the adjacent one. The hydraulic conductivity is proposed to be a weighted average valuek¯nt between the centres of clay elements. If the hydraulic conductivities in the flow path direction have large differences in magnitude, they can lead to numerical problems. To avoid these problems, the weighted average hydraulic conductivity is calculated using equation (18) and the weighted average dry density ρd¯nt, which is the value weighted with respect to the thickness of clay elements.Based on the derived equation (4), the suction difference between elements ΔUent associated with the hydraulic head difference gives the volume of transported water from element n − 1 to element n as(5)Qnt=2πrn−1t⟶rntk¯ntΔZlnrn−1t/rnt·ΔUentγwaterΔt,where γwater is the unit weight of water. Δt is the maturation time in the tth time step, assuming T to be the total time. Thus, t∈1,2,3,…,T/Δt, and t−1=0 means the initial status before maturation. The subscript n−1=0 represents the mud. The parameters hydraulic conductivity knt and suction Uent are a function of the dry density. The hydraulic conductivity depends on the void ratio, which is a function of dry density. Yang et al. [7] have experimentally verified that the suction is a function of dry density. The expressions of knt and Uent as function of dry density are based on laboratory experiments and presented in Section 4.2.The initial dry density of the mud in the laboratory tests was 222 kg/m3, meaning that the clay particles formed a stable clay gel. At the beginning of the maturation process, due to the easily accessible water in the mud, the water transport is from the interface of mud/dense clay to the centre of the first element at distance r1t. The permeability controlling the transport into the first element is determined by the k value of the first element. Therefore, the volume of water flowing into the first element during the initial maturation period is(6)Q1t=2πRt⟶r1tk1tΔZlnRt/r1t·ΔUe1tγwaterΔt.Along with the reduction of water in the mud, the soil particles come in closer contact with each other and from the momentk1t/lnRt/r1t≥k¯1t/lnrmudt/r1t, water is no longer considered to be easily accessible from the mud. From that point, it is assumed that water starts to be transported from the centre of the mud to the centre of the first dense clay element, and the densified mud is considered to be one element among all others. The mathematical expression of the water volume flowing into the first element during this subsequent maturation period is(7)Q1′t=2πrmudt⟶r1tk¯1tΔZlnrmudt/r1t·ΔUe1t/γwaterΔt.For each element, water is both flowing into and out of the element. The volume of water that stays in the element is(8)ΔQnt=Qnt−Qn+1t.Since the clay seals are closed at the top and the bottom, the dense clay swells only in the radial direction. The volume of clay expansion in each element is equal to the volume of water taken up by the element. Thus, the radial expansion of each element during the timeΔt becomes(9)Δrn∗t=−2rnt−1+Δrnt−12+2rnt−1+Δrnt−122+ΔQntΔZ⋅1π.The total radial expansion of all clay elements during the time step is(10)rn∗t=∑i=1NΔri∗t.The new radius of the dense clayRt after swelling becomes(11)Rt=Rt−1+rn∗t.The mud will be compressed with the same volume as the increased volume of the dense clay. The reduced radial thickness of the mudΔrmud∗t is the negative value of the total radial expansion of the dense clay elements.Based on the expansion of clay, the state parameters of clay elements and mud can be updated with time as follows:(12)rnt=rnt−1+12Δrn∗t+∑i=n+1NΔri∗t,(13)Δrnt=Δrnt−1+Δrn∗t,(14)rmudt=rmudt−1−12Δrmud∗t,(15)Δrmudt=Δrmudt−1+Δrmud∗t.Thereby, the dry density of each element becomes(16)ρdnt=ρdnt−1ΔrntrntΔrnt−1rnt−1.And for the mud(17)ρdmudt=ρdmudt−1ΔrmudtrmudtΔrmudt−1rmudt−1.


### 4.2. Soil Properties for the Numerical Modelling

To understand how the physicochemical processes of swelling clay are accounted for in the numerical model, these processes need to be briefly described. Air-dry expansive smectite clay has one thin layer of water molecules between the thin clay lamellae and on the outside of the stacks of lamellae. In humid air and submerged in water, additional water is adsorbed by the strong affinity of the clay to water molecules, which become ordered through a pseudo-crystalline process. The adsorption is driven by this affinity (hydration potential) that causes tension in the water (“suction”) and strives to disjoin adjacent particles by causing a “swelling pressure.” The swelling pressure is logically of the same order of magnitude as the tension and in practice one assumes that tension is a negative pressure (“suction”) while swelling pressure is considered as a positive, real pressure. The affinity to water reduces with increasing numbers of adsorbed water “hulls” associated with the clay particles; when only a few water molecule layers are adsorbed, no further swelling will occur.During the numerical modelling, the hydraulic conductivityk and suction for the selected clay material must be known for each element and time step. The suction is expected to be the negative value of the swelling pressure Ps. The parameters k and Ps can be calculated utilizing the following relationships:(18)k=6⋅10−6exp−0.01ρd,(19)Ps=−2.16114⋅106+0.49613⋅106exp1.53378⋅10−3ρd,ρd≥960kgm3,0,ρd<960kgm3.Based on results obtained from laboratory experiments on Holmehus clay conducted by the authors [7], the relations (18) and (19) are established for quantities in SI-units. Note that the expressions will change with the type of clay, i.e., if the mineral composition and the smectite content are different. Equation (19) is used to calculate the suction difference between elements, a parameter introduced in equation (5).The physicochemical processes of swelling clay are included in the numerical model by the experimentally established relations of hydraulic conductivity (equation (18)) and swelling pressure (equation (19)), both as a function of dry density. These relations should be seen as soil properties needed for the numerical model.


### 4.3. Numerical Simulation Procedure

In the numerical simulation, modelling begins with the selection of material properties, i.e., the initial dry density of mud and dense clay, the geometrical parametersR, Δr, r, and the initial hydraulic conductivity and suction. The modelling was implemented using MATLAB. The procedure is schematically shown in Figure 4 and can be described as follows:(i)
Assume the number of elements to beN and the time for each time step in the computation to be Δt sec. The total computation time is T, which is divided into two parts: t1∈1,2,…,m and t2∈m+1,…,T/Δt. t1 represents the time step numbers for maturation in the initial stage where water moves from the mud/dense clay interface to the first dense clay element, while t2 represents the time step numbers in the second stage implying water movement from the centre of the mud to the adjacent dense clay element (Figure 4, right.).(ii)
Computation starts from the first staget1, followed by the second stage t2. In both of these two stages, the volume of water flowing into the first element and the subsequent elements Qnt is computed. Utilizing the volume of water staying in the elements ΔQnt, the swelling Δrn∗t of each element is computed.(iii)
In each time step, the state parameters of each element are assumed to be constant. As soon asΔrn∗t is computed, all the state parameters R, Δr, r, ρd, k, Ue for both mud and clay elements have to be updated for the next time step (equations (11)–(19)). The iterative simulation process continues until the total time T has been reached, the computation finished, and a solution obtained.(iv)
The procedure above is repeated for the computation withN+a elements with time Δt′ (equation (20), Section 5.1). Here, a∈1,2,…. When the computation is completed, the solution is compared with the solution for N elements. If the solutions converge according to a convergence criterion, the solution is accepted. Otherwise, the number of elements is gradually increased with the associated time interval until convergence is reached (Figure 4, left).Figure 4
Flowchart of numerical solution procedure. Left: overview of the procedure. Right: computation with a certain number of elementsN starting with the first time stage t1, followed by the second time stage t2.


## 4.1. The Proposed Model

In the case of maturation of clay seals in a borehole, the key assumptions and criteria in the model are as follows:(1)
Darcy’s law is valid between each element in the numerical algorithm(2)
Both dense clay and mud are initially homogeneous(3)
The total volume inside the hollow cylinder filter, containing fully saturated dense clay and mud, is constant during the entire maturation process(4)
Expansion/contraction and flow take place only in the radial direction(5)
There is no temperature dependence, e.g., from heat production caused by radioactive decayDarcy’s law is expected to give an accurate representation of the flow within a porous medium for small discharge velocities [10]. Water flow is slow in the maturation process, so it is assumed that Darcy’s law is applicable on the elements in the numerical algorithm. Darcy’s law states that q=Aki, where the fluid flux q is related to the cross-sectional area A, the hydraulic conductivity k, and the hydraulic gradient i. The permeated fluid volume Q in time t0 is given by(1)Q=Akit0.With the height of the cylinder specimenH and the radius R, the permeated fluid volume becomes(2)Q=2πRHkdhdRt0,where h is the water head and dh/dR is the hydraulic gradient i.Rearranging and integrating equation (2) gives(3)Q⋅∫R1R21RdR=∫h1h22πkt0Hdh.Thus,(4)Q=2πkHlnR2/R1h2−h1·t0,where R1 and R2 represent the flow path. h2−h1 is the hydraulic head difference causing water transport.The hydraulic gradient and thereby the hydration/dehydration of the dense clay is a function of the position of the clay elements and of the maturation time. Therefore, the dense clay is considered to consist of a sequence of elements (Figure3(a)), subjected to one-dimensional hydration or dehydration, implying only radial flow and clay particle migration. Because of the rotational symmetry, only the radial plane of the clay sample is viewed (Figure 3(b)). If a unit sample with height ΔZ and initial radial length R is divided into N elements (N is a natural number), the initial size of each element and the distance from its centre to the symmetry axis are Δrn and rnn=1,2,3,…,N, respectively. For the mud, the corresponding parameters are Δrmud and rmud, respectively.Figure 3
Schematic view of specimen and maturation. (a) Clay specimen surrounded by water-rich mud, absorbing water in radial direction. (b) Division of dense clay into elements for the numerical modelling.
(a)(b)Theoretically, in the clay seals, water is absorbed from the mud by the nearest clay element, followed by the water movement in the series of elements taking place from the centre of an element to the centre of the adjacent one. Therefore, in the model, the values of permeability and suction should be representative from the centre of the clay element to the centre of the adjacent one. The hydraulic conductivity is proposed to be a weighted average valuek¯nt between the centres of clay elements. If the hydraulic conductivities in the flow path direction have large differences in magnitude, they can lead to numerical problems. To avoid these problems, the weighted average hydraulic conductivity is calculated using equation (18) and the weighted average dry density ρd¯nt, which is the value weighted with respect to the thickness of clay elements.Based on the derived equation (4), the suction difference between elements ΔUent associated with the hydraulic head difference gives the volume of transported water from element n − 1 to element n as(5)Qnt=2πrn−1t⟶rntk¯ntΔZlnrn−1t/rnt·ΔUentγwaterΔt,where γwater is the unit weight of water. Δt is the maturation time in the tth time step, assuming T to be the total time. Thus, t∈1,2,3,…,T/Δt, and t−1=0 means the initial status before maturation. The subscript n−1=0 represents the mud. The parameters hydraulic conductivity knt and suction Uent are a function of the dry density. The hydraulic conductivity depends on the void ratio, which is a function of dry density. Yang et al. [7] have experimentally verified that the suction is a function of dry density. The expressions of knt and Uent as function of dry density are based on laboratory experiments and presented in Section 4.2.The initial dry density of the mud in the laboratory tests was 222 kg/m3, meaning that the clay particles formed a stable clay gel. At the beginning of the maturation process, due to the easily accessible water in the mud, the water transport is from the interface of mud/dense clay to the centre of the first element at distance r1t. The permeability controlling the transport into the first element is determined by the k value of the first element. Therefore, the volume of water flowing into the first element during the initial maturation period is(6)Q1t=2πRt⟶r1tk1tΔZlnRt/r1t·ΔUe1tγwaterΔt.Along with the reduction of water in the mud, the soil particles come in closer contact with each other and from the momentk1t/lnRt/r1t≥k¯1t/lnrmudt/r1t, water is no longer considered to be easily accessible from the mud. From that point, it is assumed that water starts to be transported from the centre of the mud to the centre of the first dense clay element, and the densified mud is considered to be one element among all others. The mathematical expression of the water volume flowing into the first element during this subsequent maturation period is(7)Q1′t=2πrmudt⟶r1tk¯1tΔZlnrmudt/r1t·ΔUe1t/γwaterΔt.For each element, water is both flowing into and out of the element. The volume of water that stays in the element is(8)ΔQnt=Qnt−Qn+1t.Since the clay seals are closed at the top and the bottom, the dense clay swells only in the radial direction. The volume of clay expansion in each element is equal to the volume of water taken up by the element. Thus, the radial expansion of each element during the timeΔt becomes(9)Δrn∗t=−2rnt−1+Δrnt−12+2rnt−1+Δrnt−122+ΔQntΔZ⋅1π.The total radial expansion of all clay elements during the time step is(10)rn∗t=∑i=1NΔri∗t.The new radius of the dense clayRt after swelling becomes(11)Rt=Rt−1+rn∗t.The mud will be compressed with the same volume as the increased volume of the dense clay. The reduced radial thickness of the mudΔrmud∗t is the negative value of the total radial expansion of the dense clay elements.Based on the expansion of clay, the state parameters of clay elements and mud can be updated with time as follows:(12)rnt=rnt−1+12Δrn∗t+∑i=n+1NΔri∗t,(13)Δrnt=Δrnt−1+Δrn∗t,(14)rmudt=rmudt−1−12Δrmud∗t,(15)Δrmudt=Δrmudt−1+Δrmud∗t.Thereby, the dry density of each element becomes(16)ρdnt=ρdnt−1ΔrntrntΔrnt−1rnt−1.And for the mud(17)ρdmudt=ρdmudt−1ΔrmudtrmudtΔrmudt−1rmudt−1.


## 4.2. Soil Properties for the Numerical Modelling

To understand how the physicochemical processes of swelling clay are accounted for in the numerical model, these processes need to be briefly described. Air-dry expansive smectite clay has one thin layer of water molecules between the thin clay lamellae and on the outside of the stacks of lamellae. In humid air and submerged in water, additional water is adsorbed by the strong affinity of the clay to water molecules, which become ordered through a pseudo-crystalline process. The adsorption is driven by this affinity (hydration potential) that causes tension in the water (“suction”) and strives to disjoin adjacent particles by causing a “swelling pressure.” The swelling pressure is logically of the same order of magnitude as the tension and in practice one assumes that tension is a negative pressure (“suction”) while swelling pressure is considered as a positive, real pressure. The affinity to water reduces with increasing numbers of adsorbed water “hulls” associated with the clay particles; when only a few water molecule layers are adsorbed, no further swelling will occur.During the numerical modelling, the hydraulic conductivityk and suction for the selected clay material must be known for each element and time step. The suction is expected to be the negative value of the swelling pressure Ps. The parameters k and Ps can be calculated utilizing the following relationships:(18)k=6⋅10−6exp−0.01ρd,(19)Ps=−2.16114⋅106+0.49613⋅106exp1.53378⋅10−3ρd,ρd≥960kgm3,0,ρd<960kgm3.Based on results obtained from laboratory experiments on Holmehus clay conducted by the authors [7], the relations (18) and (19) are established for quantities in SI-units. Note that the expressions will change with the type of clay, i.e., if the mineral composition and the smectite content are different. Equation (19) is used to calculate the suction difference between elements, a parameter introduced in equation (5).The physicochemical processes of swelling clay are included in the numerical model by the experimentally established relations of hydraulic conductivity (equation (18)) and swelling pressure (equation (19)), both as a function of dry density. These relations should be seen as soil properties needed for the numerical model.


## 4.3. Numerical Simulation Procedure

In the numerical simulation, modelling begins with the selection of material properties, i.e., the initial dry density of mud and dense clay, the geometrical parametersR, Δr, r, and the initial hydraulic conductivity and suction. The modelling was implemented using MATLAB. The procedure is schematically shown in Figure 4 and can be described as follows:(i)
Assume the number of elements to beN and the time for each time step in the computation to be Δt sec. The total computation time is T, which is divided into two parts: t1∈1,2,…,m and t2∈m+1,…,T/Δt. t1 represents the time step numbers for maturation in the initial stage where water moves from the mud/dense clay interface to the first dense clay element, while t2 represents the time step numbers in the second stage implying water movement from the centre of the mud to the adjacent dense clay element (Figure 4, right.).(ii)
Computation starts from the first staget1, followed by the second stage t2. In both of these two stages, the volume of water flowing into the first element and the subsequent elements Qnt is computed. Utilizing the volume of water staying in the elements ΔQnt, the swelling Δrn∗t of each element is computed.(iii)
In each time step, the state parameters of each element are assumed to be constant. As soon asΔrn∗t is computed, all the state parameters R, Δr, r, ρd, k, Ue for both mud and clay elements have to be updated for the next time step (equations (11)–(19)). The iterative simulation process continues until the total time T has been reached, the computation finished, and a solution obtained.(iv)
The procedure above is repeated for the computation withN+a elements with time Δt′ (equation (20), Section 5.1). Here, a∈1,2,…. When the computation is completed, the solution is compared with the solution for N elements. If the solutions converge according to a convergence criterion, the solution is accepted. Otherwise, the number of elements is gradually increased with the associated time interval until convergence is reached (Figure 4, left).Figure 4
Flowchart of numerical solution procedure. Left: overview of the procedure. Right: computation with a certain number of elementsN starting with the first time stage t1, followed by the second time stage t2.


## 5. Numerical Modelling of Laboratory Tests

Numerical modelling has been performed using parameters and boundary conditions relevant to the laboratory tests described in Section3. The initial dry density was 222 kg/m3 for the mud and 1550 kg/m3 for the dense clay. The initial radial thickness was 6.5 mm and 43 mm for mud and dense clay, respectively. Computations for 20, 40, 60, 80, 100, 150, 200, 250, 300, 400, 500, and 800 elements were made, and convergence was reached for 200 elements (see Section 5.1).


### 5.1. Simulation Results of Clay Maturation

#### 5.1.1. Dry Density of Clay Seals in Maturation

Figure5 shows that the dry density of the clay components changes with the hydration time in the numerical modelling with 200 elements. Along with the fluid flow into the fully saturated dense clay, microstructural reorganisation was successively generated in each of the elements. The study revealed that the dry density of the dense clay elements decreased from 1550 kg/m3, while the dry density of the mud increased from 222 kg/m3. In the numerical simulation, the lab-scaled clay seals became isotropic and uniform with a dry density of 1225 kg/m3 in around 540 minutes. The corresponding hydraulic conductivity and swelling pressure would be 2.90E − 11 m/s and 1.08 MPa, respectively.Figure 5
Maturation of clay seals in simulation with 200 elements. Plotting with mud and every 20th element.


#### 5.1.2. Solid Mass Conservation

The process of maturation, implying water movement in the clay seal, involves simultaneous swelling of the dense clay and densification of the soft mud in a closed system. The solid mass of both mud and dense clay is conserved throughout the process.


#### 5.1.3. TimeΔt for Each Time Step in the Numerical Computation

The accuracy of the solution is, in general, improved by using a small spacing of elements and a short spacing of time. However, to obtain an accurate and economical solution, the number of elements and the time intervals should not be chosen randomly with respect to each other. In order to simplify these choices, an empirical formula was developed and utilized in the modelling. In a numerical test procedure, the largestΔt for an acceptable solution was determined for each different number of elements. By curve-fitting the results, Δt was obtained as(20)1Δt=−50+60expN⋅6.31522⋅10−4,N≤4500.Assume the total number of elementsN was given.In a computation, a time step shorter than the designatedΔt in equation (20) can also work but is not economical because of the longer computation time. Equation (20) is applicable for numerical simulations of the maturation of the Holmehus clay. The choice of Δt has to be further investigated and verified for other clay types.


#### 5.1.4. Numerical Solution Convergence

The convergence of the numerical solution is evaluated by using the correlation coefficientρXY, implying a certain statistical relationship between two maturation curves X and Y. ρXY is calculated as covariance divided by standard deviation. For ρXY approaching unity, the similarity of the curves is high. For ρXY=1, there is a maximum positive similarity between the two curves, while for ρXY=0, there is no similarity between them.The maturation graphs for 20, 40, 60, 80, 100, 150, 200, 250, 300, 400, 500, and 800 elements (Figure6) show that convergence is reached between 150 and 300 elements. Compared with the result with 800 elements, the values of the correlation coefficient for the maturation curves are calculated and given in Table 1. Since the correlation coefficient is higher than 0.9999, the simulation with 200 elements is considered to have converged.Figure 6
Maturation of the mud in computations with different number of clay elements.Table 1
Correlation coefficient values for multiple number of elements and 800 elements.
Element no., (N)Correlation coefficient,ρXYElement no., (N)Correlation coefficient,ρXY200.9660512000.999911400.9941542500.999953600.9979263000.999974800.9990214000.9999911000.9994415000.9999971500.9998058001.000000


### 5.2. Verification of the Numerical Model

In order to validate the numerical model, it is expected to show that the simulation of the laboratory-scaled clay seals agrees with the experimental results. Although the retardation of the maturation process caused by the perforation of the tubes cannot be readily modelled, it is shown in the following that it is still possible to assess its impact.


#### 5.2.1. Maturation Influence by Perforated Tubes

Two series of laboratory tests were performed with the duration of 12 and 48 hours. In each of the series, two tests were performed, one test with and one test without the perforated tube for the same maturation time. The results are presented in Figure7. Considering the area between the maturation curve and the initial horizontal line, one can calculate the volume obtained if the area is rotated around the dry density axis as Vwith−12h,Vwith−48h and Vwithout−12h,Vwithout−48h for maturation in 12 and 48 hours with and without perforated tube, respectively. Building relationships of the volume of rotation as a function of time for Vwith−12handVwith−48h, as well as Vwithout−12handVwithout−48h, one obtains relationships between Vwith with time Twith, and Vwithout with time Twithout. Assuming that maturation is synchronous for the experiment and computation when Vwith is equal to Vwithout, one arrives at the following relationship of the maturation time for the cases with and without a tube:(21)Twithout=0.013112∗Twith1.309506,where Twith and Twithout represent the maturation time in the experiment and numerical simulation in hours, respectively. With the value of maturation time in the experiment, the corresponding time in computation at the same maturation state is calculated utilizing equation (21). The results are collected in Table 2.Figure 7
Maturation of dense clay with and without perforated tubes for 12 and 48 hours, respectively.Table 2
Maturation time in experiment compared with time in numerical simulation to reach the same maturation state.
Time in experiment,Twith,hTime in simulation,Twithout,min681220245048125963101 (year)105


#### 5.2.2. Comparison with Experimental Results

Based on the relationship of maturation timeTwith and Twithout, as well as the values of the time in Table 2, the dry densities of maturated clay in numerical simulations are compared with the corresponding dry densities from experimental results (Figure 8). In Figure 8, the blue curve, which presents the computed result with 200 elements and the mud (regarded as one element), consists of 201 data points. Each point represents one element, providing information on the dry density and the position in the seal for the maturation time studied. To determine the dry density in the experiments, the dense clay sample to be analysed after hydration was divided in the radial direction into hollow cylindrical pieces with the same thickness. Therefore, each point in the experimental red curve represents the mean dry density of one piece of clay in its radial position.Figure 8
Comparison of clay maturation in numerical modelling and laboratory experiments. In Figure 8(f), the simulation results of 1000 minutes were used for comparison since end of maturation was reached. (a) Maturation with 6 h (experiment), 8 min (simulation). (b) Maturation with 12 h (experiment), 20 min (simulation). (c) Maturation with 24 h (experiment), 50 min (simulation). (d) Maturation with 48 h (experiment), 125 min (simulation). (e) Maturation with 96 h (experiment), 310 min (simulation). (f) Maturation with 1 year (experiment), 1000 min (simulation).
(a)(b)(c)(d)(e)(f)The black curve is plotted based on the simulated result. The curve consists of 8 points, each representing the average of 25 computed elements, in order to achieve the same number of data points as in the experimental result. The orange triangle point, on the right on the blue curve, represents mud. In the simulation, the thickness of the mud changes with the maturation of the dense clay elements, so the position of the orange point changes as well. In contrast, the tested mud in the experiments contains some expanded dense clay since the same thickness of mud as the initial thickness was tested without any regard to the maturation process. Therefore, the experimentally determined dry density of the mud becomes increasingly incorrect as the maturation process progresses. The comparison in Figure8 should mostly be focused on the dense clay part to the left.Figure8 indicates that the numerical model works in a good way since the simulated results agree quite well with the experimental results for 6, 12, 24, and 48 hours and one year of maturation. However, for the maturation time of 96 hours in the experimental test in Figure 8(e), the difference between the simulated and the experimental results is considerably large. It seems that the numerically simulated results are ahead of the experimental results in the maturation process. A possible reason for this discrepancy is the relationship equation (21), which is based on the limited number of experimental tests resulting in this equation. The agreement between the simulated and the experimental results in Figure 8, in particular Figure 8(e), is expected to improve if the relationship given in equation (21) is developed further based on more laboratory tests.


### 5.3. Maturation Progress of Clay Seals

The maturation process of clay seals includes both water transport and clay migration. Water transport is generated by suction differences. Water is absorbed in the interlamellar space of the crystals, as a result of the hydration of the exchangeable cations and the osmotic phenomena [11]. The clay migration is generated by the expansion of the crystal layers and the microstructural reorganisation.From the result of the numerical simulation (with 200 elements), the maturation progress of the clay seals can be studied in Figure9. A water content ratio, corresponding to the colour bar in Figures 9(a)–9(c), 9(g)–9(i), 9(m)–9(o), is introduced for visual representation of the water transport in the mud and dense clay by colour shadings. The water content ratio is normalised between 0 and 1 and is calculated by(22)water content ratio=w−wdense clayinitialwmudinitial−wdense clayinitial,where w is the water content in the element, wdense clayinitial is the initial water content of the dense clay, and wmudinitial is the initial water content of the mud. The value 0 in the colour bar thus represents the same water content as the dense clay has initially, while 1 denotes the same water content as the mud has initially. A higher value of the water content ratio means higher water content.Figure 9
Simulated maturation process of the clay seals. (a, b, c, g, h, i, m, n, o) The water transport in clay as a function of maturation time. (d, e, f, j, k l, p, q, r) The clay migration, presenting density changes of different parts of the clay seal as a function of time.
(a)(b)(c)(d)(e)(f)(g)(h)(i)(j)(k)(l)(m)(n)(o)(p)(q)(r)At the beginning of the maturation process, the water content in the area of the interface of mud and dense clay is very high. Therefore, water moves rapidly into the part of the dense clay that is very close to the mud (Figures9(a)–9(c)). With water gradually spreading into more elements, the dense clay becomes wetter but will not be as water-rich as the elements next to the mud in the earliest stage (Figures 9(g)–9(i)). As the maturation progresses, the difference in water content between the elements will be less significant as in the initial stage. The clay seal becomes increasingly more uniform, and the maturation process slows down gradually since the water-driving suction becomes lower with time (Figures 9(m) and 9(n)). After a long time, the water content is evenly distributed in the entire seal, and homogenization is then achieved (Figure 9(o)).In the modelling, the solid mass in mud and clay elements is conserved during the process. Because of the expandability of smectite clay, the variation of the dry density is associated with a change in volume of the elements, which remain fully saturated throughout the maturation process. Initially, the mud loses a substantial volume of water and densifies rapidly until it reaches a density that approaches that of the adjacent elements. During this period, the red mud-dot in Figures9 changes its position gradually and moves closer to the filter boundary (Figures 9(d)–9(f)). Simultaneously, as shown by the blue dots, the dense clay elements that are located close to the mud absorb water, expand, and move towards the mud.Along with more elements closer to the centre of the sample gradually undergoing expansion, the swelling forces further compress the mud and some elements close to the mud. They will temporarily become even denser than some elements in the more central region (Figures9(j)–9(l)). With time, elements in the central part become less dense than the mud and the elements in its vicinity (Figures 9(p) and 9(q)). Water is sucked from the softer parts of the clay by the denser parts. The volumes of mud and elements slightly adjust until they ultimately obtain the same density (Figures 9(r)).


## 5.1. Simulation Results of Clay Maturation

### 5.1.1. Dry Density of Clay Seals in Maturation

Figure5 shows that the dry density of the clay components changes with the hydration time in the numerical modelling with 200 elements. Along with the fluid flow into the fully saturated dense clay, microstructural reorganisation was successively generated in each of the elements. The study revealed that the dry density of the dense clay elements decreased from 1550 kg/m3, while the dry density of the mud increased from 222 kg/m3. In the numerical simulation, the lab-scaled clay seals became isotropic and uniform with a dry density of 1225 kg/m3 in around 540 minutes. The corresponding hydraulic conductivity and swelling pressure would be 2.90E − 11 m/s and 1.08 MPa, respectively.Figure 5
Maturation of clay seals in simulation with 200 elements. Plotting with mud and every 20th element.


### 5.1.2. Solid Mass Conservation

The process of maturation, implying water movement in the clay seal, involves simultaneous swelling of the dense clay and densification of the soft mud in a closed system. The solid mass of both mud and dense clay is conserved throughout the process.


### 5.1.3. TimeΔt for Each Time Step in the Numerical Computation

The accuracy of the solution is, in general, improved by using a small spacing of elements and a short spacing of time. However, to obtain an accurate and economical solution, the number of elements and the time intervals should not be chosen randomly with respect to each other. In order to simplify these choices, an empirical formula was developed and utilized in the modelling. In a numerical test procedure, the largestΔt for an acceptable solution was determined for each different number of elements. By curve-fitting the results, Δt was obtained as(20)1Δt=−50+60expN⋅6.31522⋅10−4,N≤4500.Assume the total number of elementsN was given.In a computation, a time step shorter than the designatedΔt in equation (20) can also work but is not economical because of the longer computation time. Equation (20) is applicable for numerical simulations of the maturation of the Holmehus clay. The choice of Δt has to be further investigated and verified for other clay types.


### 5.1.4. Numerical Solution Convergence

The convergence of the numerical solution is evaluated by using the correlation coefficientρXY, implying a certain statistical relationship between two maturation curves X and Y. ρXY is calculated as covariance divided by standard deviation. For ρXY approaching unity, the similarity of the curves is high. For ρXY=1, there is a maximum positive similarity between the two curves, while for ρXY=0, there is no similarity between them.The maturation graphs for 20, 40, 60, 80, 100, 150, 200, 250, 300, 400, 500, and 800 elements (Figure6) show that convergence is reached between 150 and 300 elements. Compared with the result with 800 elements, the values of the correlation coefficient for the maturation curves are calculated and given in Table 1. Since the correlation coefficient is higher than 0.9999, the simulation with 200 elements is considered to have converged.Figure 6
Maturation of the mud in computations with different number of clay elements.Table 1
Correlation coefficient values for multiple number of elements and 800 elements.
Element no., (N)Correlation coefficient,ρXYElement no., (N)Correlation coefficient,ρXY200.9660512000.999911400.9941542500.999953600.9979263000.999974800.9990214000.9999911000.9994415000.9999971500.9998058001.000000


## 5.1.1. Dry Density of Clay Seals in Maturation

Figure5 shows that the dry density of the clay components changes with the hydration time in the numerical modelling with 200 elements. Along with the fluid flow into the fully saturated dense clay, microstructural reorganisation was successively generated in each of the elements. The study revealed that the dry density of the dense clay elements decreased from 1550 kg/m3, while the dry density of the mud increased from 222 kg/m3. In the numerical simulation, the lab-scaled clay seals became isotropic and uniform with a dry density of 1225 kg/m3 in around 540 minutes. The corresponding hydraulic conductivity and swelling pressure would be 2.90E − 11 m/s and 1.08 MPa, respectively.Figure 5
Maturation of clay seals in simulation with 200 elements. Plotting with mud and every 20th element.


## 5.1.2. Solid Mass Conservation

The process of maturation, implying water movement in the clay seal, involves simultaneous swelling of the dense clay and densification of the soft mud in a closed system. The solid mass of both mud and dense clay is conserved throughout the process.


## 5.1.3. TimeΔt for Each Time Step in the Numerical Computation

The accuracy of the solution is, in general, improved by using a small spacing of elements and a short spacing of time. However, to obtain an accurate and economical solution, the number of elements and the time intervals should not be chosen randomly with respect to each other. In order to simplify these choices, an empirical formula was developed and utilized in the modelling. In a numerical test procedure, the largestΔt for an acceptable solution was determined for each different number of elements. By curve-fitting the results, Δt was obtained as(20)1Δt=−50+60expN⋅6.31522⋅10−4,N≤4500.Assume the total number of elementsN was given.In a computation, a time step shorter than the designatedΔt in equation (20) can also work but is not economical because of the longer computation time. Equation (20) is applicable for numerical simulations of the maturation of the Holmehus clay. The choice of Δt has to be further investigated and verified for other clay types.


## 5.1.4. Numerical Solution Convergence

The convergence of the numerical solution is evaluated by using the correlation coefficientρXY, implying a certain statistical relationship between two maturation curves X and Y. ρXY is calculated as covariance divided by standard deviation. For ρXY approaching unity, the similarity of the curves is high. For ρXY=1, there is a maximum positive similarity between the two curves, while for ρXY=0, there is no similarity between them.The maturation graphs for 20, 40, 60, 80, 100, 150, 200, 250, 300, 400, 500, and 800 elements (Figure6) show that convergence is reached between 150 and 300 elements. Compared with the result with 800 elements, the values of the correlation coefficient for the maturation curves are calculated and given in Table 1. Since the correlation coefficient is higher than 0.9999, the simulation with 200 elements is considered to have converged.Figure 6
Maturation of the mud in computations with different number of clay elements.Table 1
Correlation coefficient values for multiple number of elements and 800 elements.
Element no., (N)Correlation coefficient,ρXYElement no., (N)Correlation coefficient,ρXY200.9660512000.999911400.9941542500.999953600.9979263000.999974800.9990214000.9999911000.9994415000.9999971500.9998058001.000000


## 5.2. Verification of the Numerical Model

In order to validate the numerical model, it is expected to show that the simulation of the laboratory-scaled clay seals agrees with the experimental results. Although the retardation of the maturation process caused by the perforation of the tubes cannot be readily modelled, it is shown in the following that it is still possible to assess its impact.


### 5.2.1. Maturation Influence by Perforated Tubes

Two series of laboratory tests were performed with the duration of 12 and 48 hours. In each of the series, two tests were performed, one test with and one test without the perforated tube for the same maturation time. The results are presented in Figure7. Considering the area between the maturation curve and the initial horizontal line, one can calculate the volume obtained if the area is rotated around the dry density axis as Vwith−12h,Vwith−48h and Vwithout−12h,Vwithout−48h for maturation in 12 and 48 hours with and without perforated tube, respectively. Building relationships of the volume of rotation as a function of time for Vwith−12handVwith−48h, as well as Vwithout−12handVwithout−48h, one obtains relationships between Vwith with time Twith, and Vwithout with time Twithout. Assuming that maturation is synchronous for the experiment and computation when Vwith is equal to Vwithout, one arrives at the following relationship of the maturation time for the cases with and without a tube:(21)Twithout=0.013112∗Twith1.309506,where Twith and Twithout represent the maturation time in the experiment and numerical simulation in hours, respectively. With the value of maturation time in the experiment, the corresponding time in computation at the same maturation state is calculated utilizing equation (21). The results are collected in Table 2.Figure 7
Maturation of dense clay with and without perforated tubes for 12 and 48 hours, respectively.Table 2
Maturation time in experiment compared with time in numerical simulation to reach the same maturation state.
Time in experiment,Twith,hTime in simulation,Twithout,min681220245048125963101 (year)105


### 5.2.2. Comparison with Experimental Results

Based on the relationship of maturation timeTwith and Twithout, as well as the values of the time in Table 2, the dry densities of maturated clay in numerical simulations are compared with the corresponding dry densities from experimental results (Figure 8). In Figure 8, the blue curve, which presents the computed result with 200 elements and the mud (regarded as one element), consists of 201 data points. Each point represents one element, providing information on the dry density and the position in the seal for the maturation time studied. To determine the dry density in the experiments, the dense clay sample to be analysed after hydration was divided in the radial direction into hollow cylindrical pieces with the same thickness. Therefore, each point in the experimental red curve represents the mean dry density of one piece of clay in its radial position.Figure 8
Comparison of clay maturation in numerical modelling and laboratory experiments. In Figure 8(f), the simulation results of 1000 minutes were used for comparison since end of maturation was reached. (a) Maturation with 6 h (experiment), 8 min (simulation). (b) Maturation with 12 h (experiment), 20 min (simulation). (c) Maturation with 24 h (experiment), 50 min (simulation). (d) Maturation with 48 h (experiment), 125 min (simulation). (e) Maturation with 96 h (experiment), 310 min (simulation). (f) Maturation with 1 year (experiment), 1000 min (simulation).
(a)(b)(c)(d)(e)(f)The black curve is plotted based on the simulated result. The curve consists of 8 points, each representing the average of 25 computed elements, in order to achieve the same number of data points as in the experimental result. The orange triangle point, on the right on the blue curve, represents mud. In the simulation, the thickness of the mud changes with the maturation of the dense clay elements, so the position of the orange point changes as well. In contrast, the tested mud in the experiments contains some expanded dense clay since the same thickness of mud as the initial thickness was tested without any regard to the maturation process. Therefore, the experimentally determined dry density of the mud becomes increasingly incorrect as the maturation process progresses. The comparison in Figure8 should mostly be focused on the dense clay part to the left.Figure8 indicates that the numerical model works in a good way since the simulated results agree quite well with the experimental results for 6, 12, 24, and 48 hours and one year of maturation. However, for the maturation time of 96 hours in the experimental test in Figure 8(e), the difference between the simulated and the experimental results is considerably large. It seems that the numerically simulated results are ahead of the experimental results in the maturation process. A possible reason for this discrepancy is the relationship equation (21), which is based on the limited number of experimental tests resulting in this equation. The agreement between the simulated and the experimental results in Figure 8, in particular Figure 8(e), is expected to improve if the relationship given in equation (21) is developed further based on more laboratory tests.


## 5.2.1. Maturation Influence by Perforated Tubes

Two series of laboratory tests were performed with the duration of 12 and 48 hours. In each of the series, two tests were performed, one test with and one test without the perforated tube for the same maturation time. The results are presented in Figure7. Considering the area between the maturation curve and the initial horizontal line, one can calculate the volume obtained if the area is rotated around the dry density axis as Vwith−12h,Vwith−48h and Vwithout−12h,Vwithout−48h for maturation in 12 and 48 hours with and without perforated tube, respectively. Building relationships of the volume of rotation as a function of time for Vwith−12handVwith−48h, as well as Vwithout−12handVwithout−48h, one obtains relationships between Vwith with time Twith, and Vwithout with time Twithout. Assuming that maturation is synchronous for the experiment and computation when Vwith is equal to Vwithout, one arrives at the following relationship of the maturation time for the cases with and without a tube:(21)Twithout=0.013112∗Twith1.309506,where Twith and Twithout represent the maturation time in the experiment and numerical simulation in hours, respectively. With the value of maturation time in the experiment, the corresponding time in computation at the same maturation state is calculated utilizing equation (21). The results are collected in Table 2.Figure 7
Maturation of dense clay with and without perforated tubes for 12 and 48 hours, respectively.Table 2
Maturation time in experiment compared with time in numerical simulation to reach the same maturation state.
Time in experiment,Twith,hTime in simulation,Twithout,min681220245048125963101 (year)105


## 5.2.2. Comparison with Experimental Results

Based on the relationship of maturation timeTwith and Twithout, as well as the values of the time in Table 2, the dry densities of maturated clay in numerical simulations are compared with the corresponding dry densities from experimental results (Figure 8). In Figure 8, the blue curve, which presents the computed result with 200 elements and the mud (regarded as one element), consists of 201 data points. Each point represents one element, providing information on the dry density and the position in the seal for the maturation time studied. To determine the dry density in the experiments, the dense clay sample to be analysed after hydration was divided in the radial direction into hollow cylindrical pieces with the same thickness. Therefore, each point in the experimental red curve represents the mean dry density of one piece of clay in its radial position.Figure 8
Comparison of clay maturation in numerical modelling and laboratory experiments. In Figure 8(f), the simulation results of 1000 minutes were used for comparison since end of maturation was reached. (a) Maturation with 6 h (experiment), 8 min (simulation). (b) Maturation with 12 h (experiment), 20 min (simulation). (c) Maturation with 24 h (experiment), 50 min (simulation). (d) Maturation with 48 h (experiment), 125 min (simulation). (e) Maturation with 96 h (experiment), 310 min (simulation). (f) Maturation with 1 year (experiment), 1000 min (simulation).
(a)(b)(c)(d)(e)(f)The black curve is plotted based on the simulated result. The curve consists of 8 points, each representing the average of 25 computed elements, in order to achieve the same number of data points as in the experimental result. The orange triangle point, on the right on the blue curve, represents mud. In the simulation, the thickness of the mud changes with the maturation of the dense clay elements, so the position of the orange point changes as well. In contrast, the tested mud in the experiments contains some expanded dense clay since the same thickness of mud as the initial thickness was tested without any regard to the maturation process. Therefore, the experimentally determined dry density of the mud becomes increasingly incorrect as the maturation process progresses. The comparison in Figure8 should mostly be focused on the dense clay part to the left.Figure8 indicates that the numerical model works in a good way since the simulated results agree quite well with the experimental results for 6, 12, 24, and 48 hours and one year of maturation. However, for the maturation time of 96 hours in the experimental test in Figure 8(e), the difference between the simulated and the experimental results is considerably large. It seems that the numerically simulated results are ahead of the experimental results in the maturation process. A possible reason for this discrepancy is the relationship equation (21), which is based on the limited number of experimental tests resulting in this equation. The agreement between the simulated and the experimental results in Figure 8, in particular Figure 8(e), is expected to improve if the relationship given in equation (21) is developed further based on more laboratory tests.


## 5.3. Maturation Progress of Clay Seals

The maturation process of clay seals includes both water transport and clay migration. Water transport is generated by suction differences. Water is absorbed in the interlamellar space of the crystals, as a result of the hydration of the exchangeable cations and the osmotic phenomena [11]. The clay migration is generated by the expansion of the crystal layers and the microstructural reorganisation.From the result of the numerical simulation (with 200 elements), the maturation progress of the clay seals can be studied in Figure9. A water content ratio, corresponding to the colour bar in Figures 9(a)–9(c), 9(g)–9(i), 9(m)–9(o), is introduced for visual representation of the water transport in the mud and dense clay by colour shadings. The water content ratio is normalised between 0 and 1 and is calculated by(22)water content ratio=w−wdense clayinitialwmudinitial−wdense clayinitial,where w is the water content in the element, wdense clayinitial is the initial water content of the dense clay, and wmudinitial is the initial water content of the mud. The value 0 in the colour bar thus represents the same water content as the dense clay has initially, while 1 denotes the same water content as the mud has initially. A higher value of the water content ratio means higher water content.Figure 9
Simulated maturation process of the clay seals. (a, b, c, g, h, i, m, n, o) The water transport in clay as a function of maturation time. (d, e, f, j, k l, p, q, r) The clay migration, presenting density changes of different parts of the clay seal as a function of time.
(a)(b)(c)(d)(e)(f)(g)(h)(i)(j)(k)(l)(m)(n)(o)(p)(q)(r)At the beginning of the maturation process, the water content in the area of the interface of mud and dense clay is very high. Therefore, water moves rapidly into the part of the dense clay that is very close to the mud (Figures9(a)–9(c)). With water gradually spreading into more elements, the dense clay becomes wetter but will not be as water-rich as the elements next to the mud in the earliest stage (Figures 9(g)–9(i)). As the maturation progresses, the difference in water content between the elements will be less significant as in the initial stage. The clay seal becomes increasingly more uniform, and the maturation process slows down gradually since the water-driving suction becomes lower with time (Figures 9(m) and 9(n)). After a long time, the water content is evenly distributed in the entire seal, and homogenization is then achieved (Figure 9(o)).In the modelling, the solid mass in mud and clay elements is conserved during the process. Because of the expandability of smectite clay, the variation of the dry density is associated with a change in volume of the elements, which remain fully saturated throughout the maturation process. Initially, the mud loses a substantial volume of water and densifies rapidly until it reaches a density that approaches that of the adjacent elements. During this period, the red mud-dot in Figures9 changes its position gradually and moves closer to the filter boundary (Figures 9(d)–9(f)). Simultaneously, as shown by the blue dots, the dense clay elements that are located close to the mud absorb water, expand, and move towards the mud.Along with more elements closer to the centre of the sample gradually undergoing expansion, the swelling forces further compress the mud and some elements close to the mud. They will temporarily become even denser than some elements in the more central region (Figures9(j)–9(l)). With time, elements in the central part become less dense than the mud and the elements in its vicinity (Figures 9(p) and 9(q)). Water is sucked from the softer parts of the clay by the denser parts. The volumes of mud and elements slightly adjust until they ultimately obtain the same density (Figures 9(r)).


## 6. Numerical Modelling of the Evolution of Clay Seals in a Real Borehole

In a VDH for HLW disposal, the upper 500 m to 2000 m part (Figure1) is tightly sealed with smectite-rich clay which is also the case for the lower 2000 m part, where the containers are filled with clay and nuclear waste canisters. In this section, the proposed numerical model is applied for real borehole cases. In the numerical model, no attention is paid to the impact of the heat production caused by radioactive decay. This implies that the numerical model in its present form might only be valid for the upper part of the borehole.Juhlin and Sandstedt [12] and Pusch [6] have proposed a VDH diameter of 80 cm for the borehole and 60 cm for the tubes with clay-embedded HLW canisters, or alternatively 50 cm for the tubes if a borehole casing is required. The maturation of clay seals in a full-scale borehole with the diameter of 80 cm/60 cm for borehole/clay plugs, termed “real Case 1,” and 80 cm/50 cm, termed “real Case 2,” has been numerically simulated with the model. The result from the simulations is presented in Table 3 and Figure 10.Table 3
Numerical simulation of maturation in reference to real VDH with Holmehus clay.
CaseInitial dry density (dense clay/mud) (kg/m3)End of processCaseInitial dry density (dense clay/mud) (kg/m3)End of processDry density (clay seals) (kg/m3)Time (days)km/sUeMPaDry density (dense clay/mud) (kg/m3)Time (days)11550/222969784.06E − 100.00221650/222(910–960)/450431650/2221025742.12E − 100.2292000/222(878–960)/960751800/22211121048.89E − 110.5702250/2221014/10141371900/22211661555.18E − 110.8061900/245(891–960)/933571900/24511781404.59E − 110.8612000/245(879–960)/960781900/49012832771.61E − 111.3891900/4901042/104285Figure 10
Maturation of clay seals in real Case 2. Plotting with mud and every 100th element.For “real Case 1” (Table3), with the 80 cm/60 cm design of the borehole repository, a suitable initial dry density for Holmehus clay plugs is more than 1800 kg/m3, assuming the initial dry density of the mud is 222 kg/m3. The maturation process will take more than 100 days. If the initial dry density of the mud is high, e.g., around 490 kg/m3, the time to reach homogenization will be much longer.For “real Case 2,” with an initial dry density of 222 kg/m3 for the mud and 1650 kg/m3 for the dense Holmehus clay, the migration of solids is not fast enough to densify the mud early because of the high ratio between the thickness of mud and dense clay. As shown in Figure 10, the ultimate dry density of the mud is around 450 kg/m3, while that of the dense clay is 910–960 kg/m3.The maturation progress ends when the dry density of Holmehus clay reaches 960 kg/m3 since the water suction potential of the clay with this dry density is extremely low. The heterogeneous condition lasts until the initial dry density of dense clay is as high as 2000 kg/m3, in which case the mud and clay plugs become relatively uniform with a dry density around 960 kg/m3 after about 75 days. With a higher initial dry density of 490 kg/m3 for the mud and 1900 kg/m3 for the dense clay, the dry density of the clay seals is predicted to be 1042 kg/m3 after 85 days. However, with the high density of the mud, a significant increase in shear resistance for installing clay packages is expected. To speed up the maturation of the entire clay system and to achieve ultimate homogeneous sealing status, it would be preferable to use a more smectite-rich clay than Holmehus clay or to install it in an incompletely saturated form.


## 7. Concluding Remarks

The goal of this paper is to present a numerical model for the simulation of the maturation of clay seals. It can be concluded that the proposed numerical model works well in describing the maturation process that takes place in the smectite clay barrier in the storage of spent reactor fuel in very deep boreholes.The numerical model takes into account neither the influence of the perforated tube on the maturation process nor the supercontainer surrounding the dense clay. In the present study, this influence was considered by formulating a relation between the maturation times with and without a perforated tube. This relation was based on two series of laboratory tests. It is believed that this relation could be improved by increasing the number of laboratory tests.The numerical model is based on the assumption that both the dense clay and the mud are fully saturated during the maturation process. This is probably the case in most practical situations because of the high pore pressures that exist at the great depths in the storage concept. In any case, if non-saturated conditions are to be modelled, the numerical model needs to be modified.By performing numerical simulations for various clay types, as well as using different combinations of dimensions and initial densities of mud and dense clay, the maturation process and the end result could be studied theoretically. In the storage concept, it is important that the boreholes become properly sealed to avoid radioactive waste leakage. It is also important that it is possible to install the supercontainers through the soft mud in the deep boreholes in a practical way. Whether the boreholes become properly sealed depends on reaching full homogenization on the end densities in the seals and on the maturation time. To be able to install the supercontainers, it is important that the maturation process does not go too fast. If that is the case, the clay in the seals could expand against the borehole walls, which would increase the shear resistance, and make the installation of the supercontainers difficult. All practical issues mentioned above serve as examples of what can be further studied by use of the numerical model.



---

*Source: 4014185-2020-10-14.xml*