1. Introduction
The difficulties encountered in the supply of drinking water in developing countries (DCs), and particularly in CongoBrazzaville, have led some populations in large urban centers in recent decades to resort to other means of supplying and conserving drinking water (boreholes, wells and springs, rainwater storage systems [1] [2]. Others, on the other hand, are building water storage facilities, such as underground reservoirs. These are structures built in reinforced concrete or masonry, which store and redistribute water to homes by means of a pump, thus increasing water pressure. This type of individual drinking water management appears to be an alternative solution to the problem of untimely water cuts by the national network [1]. According to popular belief, the rate of urban drinking water supply in Congo Brazzaville via the Congolese water distribution company network in the 2010s was just 14%, despite the fact that water resources are quite substantial.
By way of illustration, according to work carried out by Malanda (2014) [1] [3][5], three arrondissements (Moungali, PotoPoto and Ouenzé) in the city of Brazzaville had 256 tanks listed, including 143 reinforced concrete tanks and 113 masonry tanks. These tanks were the subject of a fullscale study on the determining physicochemical parameters of the water stored in these structures. However, it was found that these reinforced concrete tanks, located in wetlands and relatively polluted, are constantly being degraded during operation. These include structural degradation (walls), leading to problems of watertightness and durability of the entire reservoir, to the point of facilitating the migration of active pollutants inside the reservoir, which in turn alters the quality of the stored water, leading to its pollution [6][9]. According to the World Health Organization (WHO) [10], the consumption of polluted water leads to diseases such as cholera and typhoid, which cause the death of at least one and a half million children worldwide every year [5] [6].
This being the case, the description of the diffusion mechanism of polluting substances in the pore space of concrete walls, in the case of water storage in said tanks constitutes under these conditions an interesting research problem. The aim of this research work is to analyse the degree of dissolved ion contamination (${\text{NO}}_{3}^{}$
and Mn^{2+}) considered as potential pollutants, and then to describe the physical phenomenon of pollutant diffusion through this cementitious matrix via its pore network. This phenomenon can also be observed by measuring the durability of the concrete with respect to actions due to physicochemical phenomena, which are expressed in terms of ionic diffusion [3][5]. A concrete reservoir structure of this kind should normally be durable despite the aggressiveness of its immediate environment throughout its operating life. And, according to Demestre and Charron (2012), three modes of transfer can occur in this kind of situation: diffusion (due to a concentration gradient), permeability due to a pressure gradient, and capillary suction (capillary effect). Diffusion transfer generally occurs in structures exposed to harsh environments, where concentrations are very high. Ions pass through the porous matrix of the material. Thus, the rate of transfer depends on the volume fraction of the tortuosity and the connectivity of the pores [11][20].
Figure 1 below illustrates how pollutants present in the environment can enter the reservoir and gradually contaminate the water stored in underground tanks.
These two figures, which are reallife photos taken on site, illustrate the phenomenon of contamination of water stored in concrete tanks, when these water storage structures are located in a polluted, damp environment or in the presence of a water table.
Studies already carried out on this subject reveal that groundwater pollution, considered to be a source of contamination for stored water, is mainly due to domestic discharges or human activities (Figure 2). As a result, water storage structures built in these urban areas are subject to these pollutants when the concrete walls are not perfectly watertight. In other words, the phenomenon occurs as a result of the concrete’s porosity or durability.
Figure 1. Concrete wall and sealing problem, infiltration of groundwater into the reservoir [13].
Figure 2. Landfill site for untreated waste located in an urban area in [13].
For the purposes of this study, pollutant diffusion will be characterized by two physical quantities, representing durability indicators: permeability and diffusivity coefficient.
With regard to durability, Soulie et al. (2006) [21] have indicated that the mechanical behaviour of structures made of cemented granular materials is strongly influenced by the presence of interstitial liquids in the pore space. Under these conditions, when any kind of reagent transfer or phase change occurs, the liquid (solute) can generate structural degradation through dissolution of cemented joints or crystallization of solutes. In this case, the ability of the “concrete” material to be penetrated by aggressive species depends on the porous structure of the cement paste, which is a function of the type of cement and the water/cement ratio. In the mathematical formulation of the problem posed, Malanda (2014) proposed a law modelling pollutant diffusion through concrete, based on Bazant’s law (1971) [22], to reflect the nonlinearity of the phenomenon. Nevertheless, we propose in this work, another approach for the identification this time of the pollutant diffusion coefficient. The aim of this study is to identify the diffusion coefficient of a pollutant through concrete in the case of water storage in underground tanks, with a view to estimating the rheological behaviour of the material subjected to external stresses due to the aggressiveness of the surrounding environment.
This identification of the diffusion coefficient will be achieved through the following specific objectives:
carry out an experimental simulation of the diffusion phenomenon;
develop an approach for modelling the transfer mechanisms of pollutants in the cementitious matrix using partial differential equations (PDEs) in a finite volume numerical scheme;
identify the diffusion coefficient using the optimal control method, which will make it possible to test the chosen mathematical model against the description of the phenomenon of pollutant diffusion through concrete;
numerically simulate the “Ficktype” diffusion phenomenon, highlighting the concentration gradient that favours the migration of pollutants into the tank, using the code developed in Fortran language.
2. Materials and Methods
2.1. Laboratory Experimental SetUp
Two separate samples were analysed for this test phase with a view to determining the concentrations of nitrate ions (${\text{NO}}_{3}^{}$
) and manganese (Mn^{2+}), considered as potential pollutants of the stored water contained in the test tanks. The concrete reserves under consideration are miniaturized to a laboratory scale.
Figure 3 and Figure 4 illustrate the experimental conditions in the laboratory.
The experimental simulation setup consists of:
Figure 3. Plastic tank containing pollutant (manganese ions), Mn^{2+} sand and a concrete tank containing drinking water.
Figure 4. Plastic tank containing pollutant (nitrate ions), ${\text{NO}}_{3}^{}$
sand and concrete tank containing drinking water.
The pollutant solution passes through the concrete walls by the principle of diffusion, thus altering the quality of the stored water.
Figure 5 shows an explicit diagram of this experimental setup.
Figure 5. Schematic diagram of the experimental setup.
Figure 6. UV visible spectrophotometer WPA.
2.2. Experimental Procedure
The equipment used to determine concentrations is the WPA UV visible spectrophotometer (Figure 6) for time and space estimation of concentrations. In practice, nitrate ions ${\text{NO}}_{3}^{}$
(nitrave 5) and manganese ions Mn^{2+} (manganecole n 1) were used as pollutants.
2.3. Numerical Simulations of the Problem Using the Optimal Control Approach
On the basis of measurements obtained on the laboratoryscale phenomenon, the approach consists of approximating the parameters that will enable these measurements to be confirmed. This can be done by numerical simulation, if the model can approximate its physical observations, sometimes deviating from the real parameters. It is therefore necessary to add constraints or a priori factors to narrow down the range of possibilities in order to arrive at a unique solution.
2.3.1. Modelling Diffusive Phenomena
Diffusion takes place according to the law of conservation of matter and thermodynamic equilibrium.
The diffusion equation is the relationship that governs the time and space evolution of particle density. It is the result of a material balance associated with Fick’s law.
According to this hypothesis, particle flow J, in the presence of a concentration gradient in a porous medium, can be modelled by equation 1:
$J=\varnothing D\frac{\partial C}{\partial x}$
(1)
with D as the diffusion coefficient or diffusivity, and ∅ the porosity of the medium.
The concentration C of the solute being deﬁned, as the mass of the solute per unit volume of the solution and ∅, the porosity of the medium; the onedimensional particle balance in the transient regime (flux varying at each point), for an elementary volume dV of containing dN = CSdx particles and a closed surface S, gives:
At x, the total mass of solute entering the element per unit time is:
Me = j (x,t) Sdt (2)
At x + dx abscissa, the total outgoing mass during the elementary displacement is:
$Ms=j\left(x+dx,t\right)Sdt$
(3)
The diﬀerence between incoming and outgoing mass ∆M is:
$\Delta M=MeMs=j\left(x,t\right)Sdtj\left(x+dx,t\right)Sdt=\frac{\partial j\left(x,t\right)}{\partial x}Sdtdx$
(4)
Neglecting adsorption, the diﬀerence between the ﬂux entering and ﬂux leaving the element is equal to the amount of substance dissolved in the element.
The mass exchange rate in the inﬁnitésimal element is given by the following expression:
$dN\left(t+dt\right)dN\left(t\right)=\varnothing \left\{C\left(t+dt,x\right)C\left(t,x\right)\right\}Sdt=\varnothing \frac{\partial C\left(x,t\right)}{\partial t}dt\text{}Sdx$
(5)
Finally, if we assume that there is no process of creation or annihilation of particles in the elementary volume dV = Sdx, for the duration dt, the conservation of matter implies that the diﬀerence of mass per unit time must be equal to the mass exchanged per unit time, i.e.:
$\Delta M=dN\left(t+dt\right)dN\left(t\right)$
(6)
Considering the following equation:
$\frac{\partial C\left(x,t\right)}{\partial t}Sdtdx=\frac{\partial j\left(x,t\right)}{\partial x}Sdtdx$
(7)
which results in:
$\varnothing \frac{\partial C\left(x,t\right)}{\partial t}=\frac{\partial j\left(x,t\right)}{\partial x}$
(8)
According to Fick’s law in porous media:
$j=\varnothing D\frac{\partial C\left(x,t\right)}{\partial x}$
(9)
Hence:
$\varnothing \frac{\partial C\left(x,t\right)}{\partial t}=\frac{\partial}{\partial x}\left(\varnothing D\frac{\partial C\left(x,t\right)}{\partial x}\right)$
(10)
This equation will be completed with boundary and initial conditions corresponding to the experimental conditions we wish to model, for its resolution.
2.3.2. Assumptions and Conditions Adopted
The medium in which the particles are scattered is assumed not to be in motion, in the sense that there is no convection in the medium;
We have limited our field of study essentially to onedimensional (1D) problems, assuming that the concentration of particles varies only in one direction, and their mean displacement also in one direction;
Boundary conditions are of the homogeneous Neumann type. As for the initial condition, the concentration is nonzero both in the concrete and inside the tank;
Porosity and diffusion coefficient are assumed to be homogeneous.
The initial concentration of nitrate ions outside the concrete tank is: 196.8 mg/L and the diffusion coefficient of nitrate ions in water is D0 = 0.17.10 – 8 m^{2}/s ;
The initial concentration of manganese ions outside the concrete tank is 14 mg/L; and the diffusion coefficient of manganese ions in water is D0 = 0.06.10 – 8 m^{2}/s;
The source function is assumed to be zero, i.e. fi,j = 0;
We neglect adsorption.
The problem of contamination of water stored in an underground reservoir located in a saturated and polluted environment can be modelled mathematically by a linear partial differential equation of the parabolic type developed in numerical methods used to solue elliptic of parabolic differential equation:
$\{\begin{array}{c}\varnothing \left(x\right)\frac{\partial C}{\partial t}\left(x,t\right)=\frac{\partial}{\partial x}\left(\varnothing \left(x\right)D\left(x\right)\frac{\partial C}{\partial x}\left(x,t\right)\right),\forall \left(x,t\right)\in \left]0,{L}_{t}\right[\times \left]0,T\right[\\ \frac{\partial C}{\partial x}\left(0,t\right)=\frac{\partial C}{\partial x}\left({L}_{t},t\right)=0,\forall t\in \left]0,T\right[\\ C\left(x,0\right)={C}_{0}\left(x\right),\forall x\in \left]0,{L}_{t}\right[\end{array}$
(11)
in which:
$D\left(x\right)>0$
is the diffusion coefficient; $\varnothing \left(x\right)>0$
is the porosity of the medium; ${C}_{0}$
is the given initial pollutant concentration and C the unknown concentration function defined on $\left]0,{L}_{t}\right[\times \left]0,T\right[$
and with values in $\mathbb{R}$
.
${L}_{t}>0$
(in meters) is the total length of the integration domain of the equation and T (in seconds) is the final time of the experiment.
2.4. Optimal Control Method
2.4.1. Presentation of the Optimal Control Technique
Modelling this control system on the basis of a dynamic parameter (control) involves the use of differential, integral, functional, finitedifference, partialderivative equations, etc. The aim is to stabilize the control system and optimize its performance.
Controls consist of highlighting functions or parameters that are generally subject to constraints. Optimal control, therefore, is a variational method which consists in finding the state of a numerical forecasting model, in order to best fit the observed measurements over a given period of time and make a natural link with the model equations.
Variational analysis relies on the minimization of a “cost” functional J as an additional condition to ensure the uniqueness of the solution. This function is determined by physical or numerical considerations (stability of optimization algorithms) [23][25].
2.4.2. General Formalism
For simplicity’s sake, we’ll use the finitedifference formulation on a rectangular grid with a regular mesh.
Equation (11) discretized on a grid is then written as:
$\frac{dX\left(u,t\right)}{dt}=F\left(X\left(u,t\right),u\right)$
(12)
with $u\u0404{\mathbb{R}}^{l}$
, a parameter to be identified (e.g. diffusion coefficient) and X (u, t), a vector of ${\mathbb{R}}^{N}$
, the discretized concentration. The Cauchy problem with initial condition:
$\{\begin{array}{c}\frac{dX\left(u,t\right)}{dt}=F\left(X\left(u,t\right),u\right)\\ X\left(u,0\right)={X}_{0}\end{array}$
(13)
then admits as a solution $X\left(u,t\right)$
over the interval [0, T].
We assume that the observations distributed in space and time have made it possible to define a timedependent vector Xobs (t) of ${\mathbb{R}}^{N}$
, which we’ll call the observation vector.
The cost function or objective function J(u) is defined by:
$J\left(u\right)=\frac{1}{2}\underset{0}{\overset{T}{{\displaystyle \int}}}<\left(X\left(u,t\right)Xobs\left(t\right)\right),X\left(u,t\right)Xobs\left(t\right)>dt$
(14)
where <· , ·> denotes the inner product in ${\mathbb{R}}^{N}$
.
Since the set of possible states is theoretically all ${\mathbb{R}}^{l}$
, one has an unconstrained optimization problem written as:
$\{\begin{array}{c}u\ast \u0404{\mathbb{R}}^{l}\\ J\left(\text{u}\ast \right)=\mathrm{inf}J\left(u\right)\\ u\text{}\u0404{\mathbb{R}}^{l}\end{array}$
(15)
Variational methods solve the inverse problem by directly minimizing the cost function, usually using a gradient descent method. The gradient is computed by the adjoint method, which comes from the optimal control theory of PDEs [10] [11].
2.4.3. Determining the Optimality Condition, Calculating the Gradient
Let J be a realvalued function defined on a Hilbert space ${\rm H}$
(with inner product <.,.>), of variable u.
To minimize the functional J, we must first perform the derivation operation on it. To do this, we introduce the Gateaux derivative (or directional derivative)$\widehat{J}$
of J in direction h Є ${\rm H}$
defined by:
$\widehat{J}\left(u,h\right)=\underset{\alpha \to 0}{\mathrm{lim}}\frac{J\left(u+\alpha h\right)J\left(u\right)}{\alpha}$
(16)
In the case where $\widehat{J}$
is linearly continuous in h, we can then note the following relationship:
$\widehat{J}\left(u,h\right)=\langle \nabla J\left(u\right),h\rangle $
(17)
$\nabla J\left(u\right)$
is the gradient of J with respect to u.
When J admits continuous partial derivatives for all u in ${\mathbb{R}}^{l}$
, a necessary condition for ${u}^{\ast}$
to be a minimum is that $\nabla j\left({u}^{\ast}\right)=0$
.
In general, problem (15) cannot be solved analytically. However, there are numerical solution algorithms that can be used to determine approximate solutions, including the conjugate gradient method and the quasiNewton method.
Most are iterative descent procedures which, starting from an initial value ${u}_{0}$
, seek the best approximation ${u}_{1}$
of the solution ${u}^{*}$
such that j (${u}_{1}$
) ≤ j (${u}_{0}$
) and iterate the process until a given convergence condition is satisfied.
All these methods use the gradient of the function to be minimized with respect to the control variable (in our case). Hence the importance of determining the gradient of J with respect to u.
A very efficient technique for obtaining the gradient of the cost function with respect to the control variable is the use of the adjoint system that we will present below.
To do this, we consider a vector h belonging to the same space as u and a real number α.
$X\left(u+\alpha h,t\right)$
and $X\left(u,t\right)$
are then respectively solutions of the ordinary differential equations with initial conditions:
$\{\begin{array}{c}\frac{dX\left(u+\alpha h,t\right)}{dt}=F\left(X\left(u+\alpha h\right),t\right)\\ X\left(u+\alpha h,0\right)={X}_{0}\end{array}$
(18)
and:
$\{\begin{array}{c}\frac{dX\left(u,t\right)}{dt}=F\left(X\left(u,t\right),u\right)\\ X\left(u,0\right)={X}_{0}\end{array}$
(19)
Subtracting (18) from (19) then dividing by α and taking the limit when the latter tends to, we obtain that the derivative $\widehat{X}$
of X in the h direction is then solution of the system:
$\{\begin{array}{c}\frac{d\widehat{X}\left(u,h,t\right)}{dt}=\left[\frac{\partial F}{\partial X}\right]\widehat{X}\left(u,h,t\right)+\left[\frac{\partial F}{\partial u}\right]h\\ \widehat{X}\left(u,h,0\right)=0\end{array}$
(20)
where:
$\frac{\partial F}{\partial X}$
is the Jacobian matrix of F with respect to X;
$\frac{\partial F}{\partial U}$
is the Jacobian matrix of F with respect to u.
The Gateaux derivative of J in the h direction is given by
$\widehat{J}\left(u,h\right)=\langle X\left(u,t\right)Xobs\left(t\right),\widehat{X}\left(u,t\right)\rangle $
(21)
To find the linearity of $\widehat{J}$
with respect to h, in equation (22), let’s multiply each by a variable P ${\mathbb{R}}^{\text{N}}$
, and integrate the result over the interval $\left[0,T\right]$
. The result is:
$\underset{0}{\overset{T}{{\displaystyle \int}}}\langle \frac{\partial \widehat{X}}{\partial t},P\rangle dt=\underset{0}{\overset{T}{{\displaystyle \int}}}\langle \left[\frac{\partial F}{\partial X}\right]\widehat{X},P\rangle dt+\underset{0}{\overset{T}{{\displaystyle \int}}}\langle \left[\frac{\partial F}{\partial u}\right]h,P\rangle dt$
(22)
or after integration by parts of the:
$\langle \widehat{X}\left(T\right),P\left(T\right)\rangle \langle \widehat{X}\left(0\right),P\left(0\right)\rangle =\underset{0}{\overset{T}{{\displaystyle \int}}}\widehat{X},\frac{dP}{dt}+{\left[\frac{\partial F}{\partial X}\right]}^{t}P>dt+\underset{0}{\overset{T}{{\displaystyle \int}}}\langle h,\left[\frac{\partial F}{\partial u}\right]P\rangle dt$
(23)
When P is chosen as the solution to the equation:
$\{\begin{array}{c}\frac{dP}{dt}+{\left[\frac{\partial F}{\partial X}\right]}^{t}P=X\left(u,t\right)Xobs\left(t\right)\\ P\left(T\right)=0\end{array}$
(24)
finds:
$\underset{0}{\overset{T}{{\displaystyle \int}}}\langle \widehat{X}\left(u,t\right),X\left(u,t\right)Xobs\left(t\right)\rangle dt=\underset{0}{\overset{T}{{\displaystyle \int}}}\langle h,{\left[\frac{\partial F}{\partial X}\right]}^{t}P\rangle dt$
(25)
$\widehat{J}\left(u,h\right)=\underset{0}{\overset{T}{{\displaystyle \int}}}\langle h,{\left[\frac{\partial F}{\partial u}\right]}^{t}P\rangle dt$
(26)
To obtain the gradient of J with respect to explicitly, starting from, need to h and P in the canonical bases of ${\mathbb{R}}^{l}$
and ${\mathbb{R}}^{N}$
. This gives the components of $\nabla J\left(u\right)$
in the form:
${\left(\nabla J\left(u\right)\right)}_{i}=\underset{0}{\overset{T}{{\displaystyle \int}}}{\left({\left[\frac{\partial F}{\partial u}\right]}^{t}P\right)}_{i}dt1\le i\le l$
(27)
2.4.4. Numerical Resolution of the Diffusion Equation in Dimension 1 (1D)
Different numerical methods can be used for solving diffusion problems.
In this work, the finite volume method was used to discretize equation (11) or (12) in space, and for stability reasons, the implicit Euler method was used for its discretization in time. At each time step, the resulting linear system was solved using the GMRES method.
Indeed, for a regular onedimensional mesh, with a flow conservation approximation, the finite volume method seems suitable. it adapts well to nonlinear conservation equations, as it ensures discrete conservation of extensive quantities [13]. It is a method of decomposition by subdomains (cell or control volumes) in which flows are approximated according to the principle of “finite differences”, but quite different from the partial differential equation (PDE) scheme.
The stability of this method also preserves the divergence problems that can arise with finite elements. Last but not least, this method fits in well with Fick’s law of particle flow diffusion.
Equation (20) was solved using a derivation technique for the discretized equation of (11), while (24) was solved using a transposition technique for the discretized equation of (22).
3. Results and Discussion
3.1. Presentation and Interpretation of Optimal Control Results, Obtained on the Basis of Experimental Results by Malanda (2014)
3.1.1. Verification of Gradient Accuracy
An essential point for the implementation of the method, is the calculation of the gradient as well as its verification. From the direct system, we obtain the tangent system by linearization around a perturbation, then the adjoint system by transposition. To obtain the adjoint, we write the code for the tangent system, then transpose the instructions line by line (Table 1).
To check the accuracy of the gradient obtained, we consider Taylor’s formula:
$J\left(u+\alpha h\right)j\left(u\right)=\alpha \nabla J\left(u\right).h+\frac{{\alpha}^{2}}{2}{\nabla}^{2}j\left(u\right)\left(h,h\right)+0\left({\alpha}^{2}\right)$
(28)
Where ${\nabla}^{2}j\left(u\right)$
denotes the matrix of J’s second derivatives with respect to u. We then have:
$\underset{\propto \to 0}{\mathrm{lim}}f\left(\alpha \right)=\underset{\alpha \to 0}{\mathrm{lim}}\left(\frac{J\left(u+\alpha h\right)J\left(u\right)}{\alpha \nabla J\left(u\right)h}\right)=1$
(29)
Table 1. Illustration of gradient curve verification values.
No. 
Numerical simulation 1 (sample 1) 
No. 
Numerical simulation 2 (sample 2) 

log(α) 
f(α) 

log(α) 
f(α) 
1 
−5.39794 
1.19502 
1 
−5.41854 
1.19443 
2 
−5.69897 
1.10497 
2 
−5.71957 
1.09989 
3 
−6.00 
1.05072 
3 
−6.0206 
1.04828 
4 
−6.30103 
1.02462 
4 
−6.32163 
1.02343 
5 
−6.60206 
1.01187 
5 
−6.62266 
1.01127 
6 
−6.90309 
1.00558 
6 
−6.92369 
1.00525 
7 
−7.20412 
1.00241 
7 
−7.22472 
1.00225 
8 
1.00241 
1.0008 
8 
−7.52575 
1.00076 
9 
−7.80618 
1.0003 
9 
−7.82678 
0.99974 
10 
−8.10721 
0.99957 
10 
−8.12781 
1.00027 
11 
−8.40824 
0.99884 
11 
−8.42884 
0.99949 
12 
−8.70927 
1.00121 
12 
−8.72987 
0.99716 
13 
−9.0103 
0.99756 
13 
−9.0309 
0.99639 
14 
−9.31133 
1.00048 
14 
−9.33193 
0.99328 
15 
−9.61236 
1.01655 
15 
−9.63296 
0.99172 
16 
−9.91339 
1.01655 
16 
−9.93399 
1.00727 
17 
−10.21442 
1.06328 
17 
−10.23502 
0.92022 
18 
−10.51545 
1.0516 
18 
−10.53605 
0.94509 
19 
−10.81648 
0.70107 
19 
−10.83708 
1.04457 
20 
−11.11751 
0.7478 
20 
−11.13811 
1.09431 
21 
−11.41854 
0.93476 
21 
−11.43914 
1.1938 
22 
−11.71957 
1.12171 
22 
−11.74017 
0.79587 
23 
−12.0206 
2.24341 
23 
−12.0412 
1.59173 
24 
−12.32163 
1.49561 
24 
−12.34223 
3.18346 
25 
−12.62266 
0 
25 
−5.41854 
1.19443 
Figure 9 and Figure 10 illustrate the gradient verification curves obtained from numerical simulations 1 and 2 using the optimal control technique. At the same time, they illustrate the instability of the gradient descent directions and show that, for α between ${10}^{1}et{10}^{6}$
, the value of $f\left(\alpha \right)$
is equal to 1. Thus, gradient descent is an optimization algorithm that makes it possible to find the minimum of any function by gradually converging towards it. Similarly, at the current point, a shift is made in the opposite direction to the gradient, so that the cost function decreases.
Figure 9. Variation of f (α) as a function of the logarithm decimal of α for the ammonia pollutant (numerical simulation 1).
Figure 10. Variation of f (α) as a function of the logarithm decimal of α for the ammonia pollutant (numerical simulation 2).
Gradient descent, on the other hand, is an optimization algorithm that finds the minimum of any function by gradually converging towards it. At the current point, a shift is made in the opposite direction to the gradient, so that the cost function decreases.
Figure 11 and Figure 12 summarize the characteristics of numerical simulations 1 and 2, using the optimal control technique for the ammonia pollutant, after six (06) iterations (Table 2).
The gradient optimization algorithm reduces the function at each iteration. The cost function is reduced respectively from 456.7 to 178 for numerical simulation 1 and from 450.6 to 21.3 for numerical simulation 2.
Figure 11. Variation of the cost function as a function of the number of iterations (numerical simulations 1).
Figure 12. Variation of the cost function as a function of the number of iterations. (numerical simulation 2)
Table 2. Values of the cost function as a function of the number of iterations for the NH_{3} pollutant.
Iterations 
Values of the cost function 
Numerical simulations 1 
Numerical simulations 2 
1 
456.6694 
450.681 
2 
27.26329 
28.939 
3 
22.3317 
24.924 
4 
17.95341 
21.464 
5 
17.82243 
21.388 
6 
17.7485 
21.362 
3.1.2. Concentration Values of Numerical Simulations 1 and 2 Obtained by Optimal Control on the Basis of Experimental Simulations 1 and 2 by Malanda (2014)
Table 3 shows the concentration values of numerical simulations 1 and 2 obtained by optimal control on the basis of experimental simulations 1 and 2 by Malanda (2014). Here the results of the numerical simulation tests in this control study are compared with those obtained by Mlalanda in 2014.
Table 3. Values of the cost function as a function of the number of iterations for the NH3 pollutant.
Time (week) 
Experimental measurements 1 Malanda (2014) 
Experimental measurements 2 Malanda (2014) 
Results concentrations numerical simulations 1 optimal control 
Results concentrations numerical simulations 2 optimal control 
1^{st} Week 
1.12 
1.1 
0.9136 
0.8675 
2^{nd} Week 
1.8 
1.86 
1.5947 
0.6151 
3^{rd} Week 
2.25 
2.05 
2.1321 
2.2398 
4^{th} Week 
2.87 
2.55 
2.5819 
2.779 
5^{th} Week 
3.18 
3.5 
2.9713 
3.2551 
Figure 13 and Figure 14 show the concentration variation curves of numerical simulations obtained by the optimal control method from experimental measurements by Malanda (2014).
Interpretation of the figures below reveals the following:
at the end of the first two weeks, the ion concentration results are almost parallel for sample 1, according to the shape of the curves. By the third week, however, the curves are much closer together, and appear spread out at the end of the experiments. This shows that the concentrations of sample 1 show a divergent evolution between the experimental and numerical readings;
from the first week to the second, the curves evolve, but they all cross, and from the middle of the second week onwards, they change shape, showing a fairly rapid evolution of concentrations in sample 2.
With regard to the experimental simulation curves, the experimental results are close at the start of the experiment, and the pattern at the end of the first two weeks may be due to the influence of the granular structure of the cementitious matrix, or even its porosity. The phenomenon of water saturation of the concrete (experimental simulation) accelerates pollutant diffusion as all the pores become interconnected. On the other hand, the delay observed in the first measurements can be attributed to the difficulty for (polluted) water molecules to fill all the voids in the cement matrix during the imbibition phase [12] [14] [15] [17] [18] [26].
However, in the numerical simulation, these parameters have not been taken into account, and the focus has been on diffusion, without taking into account the time required for the (polluted) water to be absorbed by the concrete, so that equilibrium of concentrations in the tank can be established.
In reality, as the phenomenon of tortuosity has not been taken into account, the proposed modelling concerns effective diffusion rather than pollutant diffusion.
Figure 13. In reality, as the phenomenon of tortuosity has not been taken into account, the proposed modelling concerns effective diffusion rather than pollutant diffusion.
Figure 14. Time variations in concentrations for numerical simulation 2 and experimental simulation 2.
3.1.3. Comparison of the Results of Numerical Simulations 1 and 2 by Malanda (2014) with Those Obtained Using the Optimal Control Method
Using experimental measurements 1 and 2 by Malanda (2014), numerical tests were carried out on the basis of the calculation code developed with the Fortan 90 language and the curves were plotted with Origin Pro software.
Table 4 shows the results obtained following the numerical simulation tests.
Table 4. Results of Concentration of numerical simulations 1 and 2 of optimal control and numerical concentrations from Malanda (2014).
Time (week) 
Experimental measurements 1 (Malanda 2014) 
Experimental measurements 2 (Malanda 2014) 
Numerical simulation 1 (Malanda 2014) 
Numerical simulation 2 (Malanda 2014) 
Optimum control concentrations 1 
Concentrationsoptimal control 2 
1^{st} Week 
1.12 
1.1 
0.1741 
0.3574 
0.9136 
0.8675 
2^{nd} Week 
1.8 
1.86 
0.9926 
1.4763 
1.5947 
1.6151 
3^{rd} Week 
2.25 
2.05 
1.8667 
2.4764 
2.1321 
2.2398 
4^{th} Week 
2.87 
2.55 
2.6193 
3.2659 
2.5819 
2.779 
5^{th} Week 
3.18 
3.5 
3.2455 
3.8846 
2.9713 
3.2551 
Similarly, the comparative analysis of the two simulations carried out (experimental and numerical) is shown by assembling the curves (Figures 1518) and from Table 4.
Figure 15. Temporal variations in concentrations from numerical simulation 2 and experimental simulation 2.
Figure 16. Temporal variations in concentrations from numerical simulation 1 and experimental simulation.
Figure 17. Temporal variations in concentrations from numerical simulations 2 and experimental simulation 1.
Figure 18. Temporal variations in concentrations from numerical simulations 1 and experimental simulation 1.
Table 5. Results of identified ammonia diffusion parameters.
Environments 
Identified ammonia diffusion parameters (m^{2}/s) 
Numerical simulation 1 
Numerical simulation 2 
Effective diffusion coefficient D1 of ammonia in sand 
D_{1} = 0.18 × 10^{−8} 
D_{1} = 0.37 × 10^{−8} 
Effective diffusion coefficient D2 of ammonia in concrete 
D_{2} = 0.13 × 10^{−7} 
D_{2} = 0.91 × 10^{−8} 
Diffusion coefficient D0 of ammonia in water 
D_{0} = 6.93 × 10^{−9} 
D_{0} = 6.93 × 10^{−9} 
Interpretation of the figures shown above reveals the following two phases:
in the first phase (at the start of the experiments), the results (concentrations) at the end of the first two weeks for the two types of simulation appear to be spread out;
the second phase (at the end of the experiments) brings these results much closer together and establishes physical convergence.
Nevertheless, the physical convergence of the results shows that the numerical method chosen better simulates the diffusion of ions in solution in the cementitious matrix (concrete) and, consequently, the effective migration of pollutants through it. On the basis of the above, it follows that the proposed model approaches experimental measurements better than that proposed by Malanda (2014). Numerically, this model better predicts the evolution of particle diffusion flow through a cemented granular material (CGM), and ultimately the transfer of dissolved pollutants from one point to another.
The ammonia diffusion parameters identified are given in Table 5.
Generally speaking, this experiment shows that the diffusion coefficient in concrete is higher than in water and sand. This implies ipso facto the influence of the nature of the material, which constitutes an explanatory variable for the differences encountered in the value of the diffusion coefficient. Indeed, the greater the porosity, the greater the proportion of pollutants in the concrete tank.
This observation can also be explained by the fact that the diffusion coefficient in the porous medium increases as the porosity decreases. Consequently, the diffusion coefficient is greater for lower porosity, and therefore inversely proportional.
3.2. Experiments with Nitrate and Manganese Ions
Experimental Laboratory Simulations
In accordance with the experimental setup, these two tanks respectively contained two pollutants (nitrates and manganese), diluted in water and then discharged into fine sand (see Figure 2 and Figure 4). The concentrations of the selected pollutants were determined over time using the abovementioned spectrophotometer.
The results of the eleven (11) concentration measurements obtained over 32 days are given in Table 6.
Table 6. Results of Nitrate and Manganese ion concentration measurements in the concrete tank containing drinking water (experimental simulation).
Experimental simulations 
Time/days 
Nitrate ions (mg/l) sample 1 
Manganese ions (mg/2) sample 2 
0 day 
0.00 
0.00 
5^{th} day 
0.408 
0.316 
7^{th} day 
0.561 
0.427 
9^{th} day 
0.747 
0.547 
12^{th} day 
0.886 
0.608 
14^{th} day 
1.129 
0.710 
16^{th} day 
1.560 
0.721 
19^{th} day 
1.702 
0.761 
21^{st} day 
1.981 
0.768 
23^{rd} day 
2.124 
0.770 
26^{th} day 
2.366 
1.010 
32^{nd} day 
2.401 
1.156 
Overall, there is a temporal variation in the concentration results measured in the two samples. In fact, both samples show gradual changes in the concentrations measured as they diffuse through the concrete walls into drinking water. We also notice :
for sample 1: nitrate ion concentration time values range from 0.00 to 2,401 mg/L; with an initial nitrate ion concentration of 196.8 ml/L outside the concrete tank.
For sample 2: manganese ion concentration time values ranging from 0.00 to 1,156 mg/L; with an initial manganese ion concentration of 14 ml/L outside the concrete tank.
So, although the diffusion of the selected pollutants is clearly established, we note that the measured values are differentiated for each measurement carried out for the two samples.
3.3. Numerical Simulations
3.3.1. Results and Critical Analysis of the Proposed Nitrate Ion Formulation
The results of numerical simulation tests on nitrate ion diffusion are summarized in Table 7.
This produces a curve with some variability in nitrate ion concentrations for the numerical test performed. Figure 19 shows the concentration of nitrate ions in the cement matrix and inside the tank containing stored water after 32 days. Figure 20, on the other hand, illustrates the temporal variations in concentrations obtained from the experimental and numerical simulations.
Table 7. Results of experimental and numerical Nitrate ion concentrations.
Time(days) 
Results of experimental concentration measurements (mg/L) 
Numerical results for concentrations (mg/L) 
Difference(mg/L) exact results 
0 day 
0.00 
0.00 
0 
5^{th} day 
0.408 
0.03213 
0.37587 
7^{th} day 
0.561 
0.09858 
0.46242 
9^{th} day 
0.747 
0.19854 
0.54846 
12^{th} day 
0.886 
0.39353 
0.49247 
14^{th} day 
1.129 
0.5438 
0.5852 
16^{th} day 
1.560 
0.70476 
0.85524 
19^{th} day 
1.702 
0.95949 
0.74251 
21^{st} day 
1.981 
1.1351 
0.8459 
23^{r}^{d} day 
2.124 
1.31366 
0.81034 
26^{th} day 
2.366 
1.58496 
0.78104 
32^{nd} day 
2.401 
2.13273 
0.26827 
Figure 19. Typical graph showing nitrate ion concentration in the cement matrix and inside the tank containing stored water after 32 days.
Figure 20. Temporal variations in nitrate ion concentrations (numerical simulation and experimental simulation).
In addition, Figure 19 above shows the spatial variation in nitrate ion concentration at each point of the chosen onedimensional mesh at a given time, based on the principle of the calculation code developed.
In addition, as the laboratory experiment has not yet reached its stationary state, the values of the effective diffusion coefficients obtained D1 and D2 for nitrate and manganese ions could evolve.
Similarly, the comparative analysis of the two types of simulations (experimental and numerical) is shown by assembling the curves (Figure 20) and from Table 7.
It should be noted that concentrations increase continuously as a function of time only, and the concentration variations observed in these evolutions confirm the random nature of diffusion: it is difficult to predict future concentrations, and the speed with which pollutant particles move is uniformly varied.
Furthermore, Figure 21 shows the temporal variations in nitrate ion concentrations (numerical simulation and experimental simulation), and highlights the following sentences:
In the first phase (at the start of the experiment), concentrations for both types of simulation appear to be spread out;
The second phase (at the end of the experiment) seems to bring these curves closer together.
These results show that the measurements taken in the laboratory experiments and the numerical simulation confirm the temporal evolution of pollutant diffusion. The process is increasing and the equilibrium state has not yet been reached, but the phenomenon seems to be stabilizing at the end of the experiment (11th measurement).
On the basis of these results, we can see that in numerical modelling of this type, the longer the duration of the experimental measurements, the less stringent the stability requirement becomes.
3.3.2. Results and Critical Analysis of the Proposed Manganese Ion Formulation
We thus obtain a curve with a certain variability of manganese ion concentrations for the numerical test carried out. Figure 21 shows the concentration of nitrate ions in the cement matrix and inside the tank containing stored water after 32 days (Table 8).
Similarly, the comparative analysis of the two types of simulation (experimental and numerical) is reflected in the assembly of curves (Figure 21) from Table 8. Figure 22 illustrates the temporal variations in concentrations between the numerical and experimental simulations. Here, the divergence in the shape of 02 curves (experimental and numerical simulations) does not call into question the approved theory, but rather provides evidence that the cementitious matrix (concrete) is not watertight to preserve the quality of the stored water. The results shown in Figure 22 above, for the pollutant used (manganese ions), seem to show that at the start of the experiment, the concentrations for the two types of simulation appear close together, and then the behaviour for the two simulations becomes almost fairly parallel until the end of the experiment (at t = 32 days). This may be explained by the fact that the diffusion time was not optimal, and by the low initial concentration of manganese ions.
The divergence of the curves of the two (02) methods (experimental measurements and numerical simulation) does not call into question the proven theory, but demonstrates that the cement matrix is not watertight enough to preserve the quality of the water in buried receivers in a polluted environment.
Figure 21. Manganese ion concentration in the cement matrix and inside the tank containing stored water after 32 days.
Figure 22. Time variations in manganese ion concentrations (numerical simulation and experimental simulation).
Table 8. Experimental and numerical manganese ion concentration (Mn^{2+}).
Time(days) 
Results of experimental concentration measurements 
Numerical concentration results 
Angle differencesExp and numerical results 
0 day 
0.00 
0.00 
0 
5^{th} day 
0.316 
0.15559 
0.16041 
7^{th} day 
0.427 
0.19635 
0.23065 
9^{th} day 
0.547 
0.23051 
0.31649 
12^{th} day 
0.608 
0.27444 
0.33356 
14^{th} day 
0.710 
0.30043 
0.40957 
16^{th} day 
0.721 
0.32449 
0.39651 
19^{th} day 
0.761 
0.35776 
0.40324 
21^{st} day 
0.768 
0.37842 
0.38958 
23^{rd} day 
0.770 
0.39808 
0.37192 
26^{th} day 
1.010 
0.42596 
0.58404 
32^{nd} day 
1.156 
0.47711 
0.67889 
Figure 23. Time variations of numerical simulations of nitrate ion concentrations with concentrations of numerical simulation 1 and simulation with porosity of 14%, 16% and 18%.
In light of these results, we can appreciate the temporal variation in concentration according to the mathematical model adopted and as a function of concrete porosity. This is because concentrations depend on porosity size. If we try to reduce the porosity to a value between 0.04 and 0.07, it becomes clear that the model has its limits. In other words, when the concrete becomes practically watertight, mathematically, pollutant diffusion stops. This is undoubtedly due to the degeneracy of the diffusion equation, which becomes stationary and changes from parabolic to elliptical (Figure 23).
It may be noted that porosity and diffusion coefficient are inversely proportional, and that the variation in major ion concentrations observed depends fundamentally on the porosity of the concrete (Table 9).
Table 10, which groups together the results of the last three (03) numerical simulations for concrete porosities of 14%, 16% and 18%, shows us that the temporal variation of the diffusion coefficient through the concrete obtained in the tank depends on the porosity of the concrete. These results show that porosity is an explanatory variable in the variation of the diffusion coefficient in the concrete obtained.
Looking at the figure, we can see that the effective porosity and diffusion coefficient are inversely proportional, and that the variation in concentration observed depends fundamentally on the porosity of the concrete.
In addition, the diffusion coefficients identified for the two selected pollutants are given in the table (Table 11).
Table 9. Nitrate Ion concentration values were simulated with concrete porosities of 14%, 16% and 18%.
Time(days) 
Concentration results numerical simulation at 8% 
Concentration results for porosity at 14% 
Concentration results for porosity at 16% 
Concentration results for porosity at 18% 
0 day 
0.00 
0.00 
0.00 
0.00 
5^{th} day 
0.03213 
0.01461 
0.01198 
0.01009 
7^{th} day 
0.09858 
0.05639 
0.04885 
0.04314 
9^{th} day 
0.19854 
0.13128 
0.11785 
0.10734 
12^{th} day 
0.39353 
0.29887 
0.27772 
0.26075 
16^{th} day 
0.5438 
0.43988 
0.41531 
0.39541 
19^{th} day 
0.70476 
0.59863 
0.57225 
0.55082 
21^{st} day 
0.95949 
0.86159 
0.83541 
0.81426 
23^{rd} day 
1.1351 
1.04911 
1.0248 
1.00543 
26^{th} day 
1.31366 
1.24373 
1.2225 
1.20601 
32^{nd} day 
1.58496 
1.5453 
1.5306 
1.52021 
Table 10. Summary of diffusion coefficient values obtained for concrete porosities equal to 14%, 16% and 18%.
Values of concrete porosity 
Effective diffusion coefficient of nitrate ions in sand 
Effective diffusion coefficient of nitrate ions in concrete 
8% 
0.11 × 10^{−}^{8} 
0.45 × 10^{−}^{9} 
14% 
0.12 × 10^{−}^{8} 
0.33 × 10^{−}^{9} 
16% 
0.12 × 10^{−}^{8} 
0.31 × 10^{−}^{9} 
18% 
0.12 × 10^{−}^{8} 
0.29 × 10^{−}^{9} 
Table 11. Results of identified diffusion coefficients for selected pollutants.
Environments 
Diffusion coefficient (m^{2}/s) 
Nitrate ions 
Manganese ions 
D1 effective diffusion coefficient in sand 
D_{1}= 0.11 × 10^{−}^{8} 
D_{1}= 0.43 × 10^{−}^{10} 
D2 diffusion coefficient in concrete 
D_{2}= 0.45 × 10^{−}^{9} 
D_{2}= 0.16 × 10^{−}^{7} 
D0 diffusion coefficient in water 
D_{0} = 0.17 × 10^{−8} 
D_{0}= 0.06 × 10^{−}^{8} 
In addition, several other authors have already worked on the problems of transfer and diffusion through the cementitious matrix that is concrete. The focus is on the porosity parameter, the microstructural analysis of concrete and the corrective aspects of pores. We can point to the work of Hailong Ye et al. (2013), Chen (2011), Tetsuya Tshida et al., (2009), Care er Have (1999), Mainguy (1999), Mchirgui (2013), and Dana (1999) [12] [13] [20]. All these works demonstrate that concrete porosity is its weakest link. Indeed, porosity facilitates ion transfer and migration via the pore network, thanks to its connectivity within the cementitious matrix.
All in all, it has been proven that ${\text{NO}}_{3}^{}$
and Mn^{2+} ions pass through the porous concrete matrix contained in the material zones. Thus, the rate of transfer depends on the volume fraction and connectivity of the soils according to Janz (1997), Ye et al. (2013), and other anthors Hailong ye et al. (2013) [26][30]. By the way, the work carried out by Chen (2011) [15], on the experimental study of the permeability of a sound concrete, under variable thermal and hydric conditions, made it possible to formulate a concrete intended for storage structures for various materials and installed in any type of environment. The results obtained are consistent with those of our research work. Nevertheless, for concrete porosities below 14%, the results can be considered optimum [13].
Also, the cases studied with the nitrate ion pollutant, for porosity values of 14%, 16% and 18%, give us results that appear reasonable. It can therefore be noted that pollutant transport in cementitious materials under the conditions mentioned is an adventitiousdiffusive phenomenon. However, in the modelling, advective transfer due to bulk movement or interstitial solution phase is taken into account, as well as ionic diffusion due to concentration gradients according to [Tetsuya Tshida et al., 2009] [28]. This justifies the fact that the transfer or diffusion coefficients depend on the characteristics of the material, and therefore on the formation of the concrete. However, the porosity of the cementitious matrix subjected to a concentration gradient during the transfer phase also depends on its microstructure. Larger porosity zones, known as “transition halos”, appear around the aggregates and promote transfer [28][31]. It should be noted, however, that the porosity of the material has a particularly strong influence on concrete when it comes to the transfer and diffusion of pollutants.
4. Conclusion
The aim of our work was to predict the diffusion coefficients of ${\text{NO}}_{3}^{}$
and Mn^{2+} ions in concrete, based on the different experimental phases and the chosen protocol. The coefficients obtained from numerical simulations and optimal control were compared with the results of experimental simulations. The results obtained, which are consistent, offer good prospects for predicting the diffusion of pollutants in concrete. Experimental measurements have also highlighted the porous nature of the cementitious matrix in concrete. Molecularscale diffusion of pollutants under a concentration gradient towards the interior of the reservoir resulted in pollution of the stored water. Two simulations (experimental and numerical) were carried out to optimize concentration dependence and characterize diffusion. In particular, these experiments demonstrated that the major ions considered as pollutants migrate through concrete at different concentrations and according to a time and space distribution. The results obtained highlight the physical convergence of experimental and numerical simulations, which clearly explain the problem of pollutant diffusion in a cementitious matrix, and raise the question of the durability of such tanks under the conditions described. This could therefore justify the fact that exchanges can probably take place between the surrounding environment and the interior of the reservoir, highlighting the concentration gradient; by migration of pollutants through the concrete, with the source being dissolved domestic waste in the free water table, if the intrinsic characteristic of porosity is highlighted.