Fuzzy Solution to the Unconﬁned Aquifer Problem

: In this article, the solution to the fuzzy second order unsteady partial differential equation (Boussinesq equation) is examined, for the case of an aquifer recharging from a lake. In the examined problem, there is a sudden rise and subsequent stabilization of the lake’s water level, thus the aquifer is recharging from the lake. The aquifer boundary conditions are fuzzy and create ambiguities to the solution of the problem. Since the aforementioned problem concerns differential equations, the generalized Hukuhara (gH) derivative was used for total derivatives, as well as the extension of this theory concerning partial derivatives. The case studies proved to follow the generalized Hukuhara (gH) derivative conditions and they offer a unique solution. The development of the aquifer water profile was examined, as well as the calculation of the recharging fuzzy water movement profiles, velocity, and volume, and the results were depicted in diagrams. According to presented results, the hydraulic engineer, being specialist in irrigation projects or in water management, could estimate the appropriate water volume quantity with an uncertainty level, given by the α -cuts.


Introduction
The horizontal water flow concerning unconfined aquifers without precipitation is described by the one-dimensional second order unsteady partial differential equation, called Boussinesq equation: . This equation was presented by [1] with the following assumptions: (1) the inertial forces are negligible and (2) the horizontal component of velocity V x does not vary depending on depth, and it is a function of (x, t) In 1904, Boussinesq presented a special solution of this nonlinear equation in the French journal "Journal de Mathématiques Pures et Appliquées". Boussinesq's solution concerned the case of an aquifer overlying an impermeable layer, with boundary conditions that are similar to those of a soil drained by a drain installed in the impermeable substratum. A solution to Boussinesq's equation is published by [2], using the method of small disturbances. In reference [3] a solution is presented to Boussinesq's linear equation concerning inclined and finite-width aquifers. A solution to Boussinesq's linear equation is presented by [4,5], named the Glover-Dumm equation, concerning the issue of drainage. Boussinesq's linear equation in a problem of intermittent instantaneous recharge of the underground aquifer is dealt by [6], with an initial parabolic surface area. A special solution to the aforementioned nonlinear equation is presented by [7], which coincides with the solution given by Boussinesq in 1904. An extensive analysis of the linear equation in the case of an aquifer discharging in a lake was presented by [2,[8][9][10][11][12][13].
Research for similar problem was carried out by [14], who presented an analytical solution of an aquifer with variable boundary conditions at river boundaries. A precise solution to Boussinesq's   (2) In Equation (2), K = the hydraulic conductivity of the aquifer, S = the specific yield of the aquifer or drainable pore space, D = the depth of the lake, h(x,t) = the depth of the aquifer, and x,t = the coordinates (spatial and temporal). The boundary conditions of the problem are: while the initial flow condition is: .
To facilitate calculations non-dimensional variables are introduced: where L = a specific length. The new resulting equation is as follows: with the new boundary conditions: and the initial flow condition:  In Equation (2), K = the hydraulic conductivity of the aquifer, S = the specific yield of the aquifer or drainable pore space, D = the depth of the lake, h(x,t) = the depth of the aquifer, and x,t = the coordinates (spatial and temporal). The boundary conditions of the problem are: while the initial flow condition is: h(x, 0) = H 0 .
To facilitate calculations non-dimensional variables are introduced: where L = a specific length. The new resulting equation is as follows: with the new boundary conditions: H(0, t) = 1, τ > 0, H(x, t) and the initial flow condition: H(ξ, 0) = 0. Fuzzy Case The fuzzy form of Equation (6) is: with the new boundary and initial conditions: where c = the spread and α = the α-cut.
Definition 2. Let X be a Banach space and H be a fuzzy set on X. The α-cuts of H are defined as where cl(supp( H(x))) denotes the closure of the support of H(x).  (9) with the new boundary and initial conditions:

Definition 3. Let
where c = the spread and α = the α-cut.

Mathematical
H is a convex fuzzy set: where c = the spread and α = the α-cut.

Mathematical
H is a convex fuzzy set:  (9) with the new boundary and initial conditions: where c = the spread and α = the α-cut.

Mathematical
H is a convex fuzzy set: Suppose that the functions if and only if one of the following two cases holds:

Fuzzy Case
The fuzzy form of equation (6) is: with the new boundary and initial conditions: where c = the spread and α = the α-cut.

Definition 1. A fuzzy set
H is a convex fuzzy set: Definition 2. Let X be a Banach space and H be a fuzzy set on X. The α-cuts of H are defined as The space of all compact and convex fuzzy sets on X is denoted as Ƒ (X).

Definition 4. Let
Suppose that the functions if and only if one of the following two cases holds: According to representation theorem of [36] and the theorem of [37], the membership function and the α-cut form of a fuzzy number H, are equivalent and in particular the α-cuts where c = the spread and α = the α-cut.

Definitions
H is a convex fuzzy set: Let X be a Banach space and H be a fuzzy set on X. The α-cuts

Definition 3.
Let Ҡ(X) be the family of all nonempty compact convex su The space of all sets on X is denoted as Ƒ (X).
if and only if one of the following two cases holds: dx . In both of the above cases, H α (x) derivative is a fuzzy number.
Notation2: The first case concerns the Hukuhara differentiability.
Definition 8. The gH-differentiability implies g-differentiability, but the inverse is not true.

Definition 9.
[gH-p] differentiability. A fuzzy-valued function H of two variables is a rule that assigns to each ordered pair of real numbers (x, t), in a set D xt a unique fuzzy number denoted by H(x, t). Let H(x, t): with the new boundary and initial conditions: where c = the spread and α = the α-cut.

Mathematical
H is a convex fuzzy set: : are real valued functions and partial differentiable w.r.t. x. [35,39,40]: Notation. The same is valid for where c = the spread and α = the α-cut.  x (

Systems of Crisp Problems
Our problem is transformed now to the following fuzzy partial differential equation: with the new boundary and initial conditions: Solutions to the fuzzy problem (20) with the boundary and initial conditions (21), utilizing the theory of [35,38,39,41] can be found, translating the above fuzzy problem to a system of second order of crisp boundary value problems, hereafter called corresponding system for the fuzzy problem. Therefore, four crisp BVPs systems are possible for the fuzzy problem, as follows: (1.1)-system: (1.2)-system: (2.1)-system: (2.2)-system: The method of solving the above problem is based on the selection of appropriate type of derivative, according to the theory of [38,39]. With such a selection, the above problem is transformed to a corresponding system of boundary value problems. A domain with the valid solution results, in which the derivatives have level sets according to the above theory of differentiability. ). τ , 0 ( H The method of solving the above problem is based on the selection of appropriate type of derivative, according to the theory of [38,39]. With such a selection, the above problem is transformed to a corresponding system of boundary value problems. A domain with the valid solution results, in which the derivatives have level sets according to the above theory of differentiability.

Solution of the (1.1)-System
Let now F = H − α , and G = H + α . Subsequently, the above system becomes: The Boltzmann transform is introduced now in Equation (22): Equation (22) takes the following form: with boundary conditions: The solution is: Based on Equation (26) the following arise: Second Equation of (1.1) According to the same procedure as above: Based now on Equation (30) the following arises: The final solution of the system is: This solution satisfies the system (1,1), as well as the boundary and initial conditions of the system. The first derivatives with respect to τ are the following, and form a fuzzy number as: Figure 3 illustrates the first derivative with respect to τ.  Figure 3 illustrates the first derivative with respect to τ. The second derivatives with respect to ξ are the following: (36) and form a fuzzy number as: The second derivatives with respect to ξ are the following: and form a fuzzy number as:    From Equations (34) and (36), the following results: Boundaries condition Initial condition

Solution of the (1.2)-System
F = H − α , G = H + α are introduced and the above system takes the form:

Solution
Adding the two equations: Let now: Ω = F + G, so the above equation becomes: and the solution is: The difference of the two equations is taken and it results to: Writing now: Φ = G − F, the above equation becomes: The solution of the new equation is obtained by introducing Boltzmann transformation z = ξ/ √ 4τ. The new equation is: with the new boundary conditions The solution of the above Equation (44) is: Introducing now the boundary conditions formulated by (45): where the function f(z) is equal to: (48) Figure 5 illustrates the function f(z) with respect to z.  Figure 5 illustrates the function f(z) with respect to z.

 
The form of the function f(z) is indicative, since for the calculation of the above series, 43 terms were taken and the convergence criterion is: The solution of the system is: The form of the function f(z) is indicative, since for the calculation of the above series, 43 terms were taken and the convergence criterion is: The solution of the system is: In the above solution an indicative value for z = 4 was selected and in Figure 6 it is shown that after the value z = 1.017, negative values of the function H − α (ξ, τ) for α = 0, 0.5 appear. This means that for the physical problem this solution is not acceptable for z > 1.017. Water 2019, 11, x FOR PEER REVIEW 13 of 22 The first derivative w.r.t. τ The following results to: The second_derivative w.r.t. ξ The following results to: The first derivative with respect to τ is a fuzzy number and the second derivative with respect to ξ is also a fuzzy number according the theory of [38]. The above system (1.2) satisfies the Equation (38) as well as the initial and boundary conditions and has fuzzy derivatives. The solution for the physical problem is valid only for z ≤ 1.017.
It is proved that the solution of the system (2.1) is similar to the solution of system (1.2) 2.2.6. Solution of the (2.2)-System It is proved that the solution of the system (2.2) is similar to the solution of system (1.1).

Darcy Velocity, Water Flow Recharging Volume
The Darcy's velocity is: ] and the following expression arises: The real velocity V α x=0 is equal to: Introducing dimensionless velocity: The volume recharging from the lake will be: The dimensionless volume is introduced:

Parameters of the Problem
For the aquifer of Figure 1, the following are assumed: K = 0.03 cm/s (gravelly coarse sand [13] Note: The cases (1.1) and (2.2) are considered here, since they give a physical solution to the problem of aquifer recharging from the lake. In the following, the Darcy velocity and the recharging volume of the aquifer are presented.
In Figure 7, the change in dimensionless depth is presented as a function of z, while Figure 8 illustrates dimensionless depth profiles as a function of ξ, for times t = 0.25 d, and 4d. Figure 9 illustrates membership functions µ H at ξ = 0.1, for times t = 0.25 d, and 4 d. Figure 7, shows that the waterfront towards the aquifer approaches the position z = 2.125. That front corresponds in real times t~0.25 d, and 4 d, at distances 125 m, 500 m apart from the lake. Figure 10, illustrates the recharging water volume at x = 0, while Figure 11 illustrates the membership function µ √

Remark
Many times we try to find the difference using trapezoidal or Gaussian membership functions. Here, some explanations concerning the above cases are given.

Trapezoidal membership function
A symmetrical trapezoidal membership function with the same fuzziness c and a core equal to ε 2 is considered. Subsequently, the boundary and initial conditions of the problem will be:

Remark
Many times we try to find the difference using trapezoidal or Gaussian membership functions. Here, some explanations concerning the above cases are given.

Trapezoidal Membership Function
A symmetrical trapezoidal membership function with the same fuzziness c and a core equal to 2ε

Remark
Many times we try to find the difference using trapezoidal or Gaussian membership functions. Here, some explanations concerning the above cases are given.

Trapezoidal Membership Function
A symmetrical trapezoidal membership function with the same fuzziness c and a core equal to 2ε is considered. Subsequently, the boundary and initial conditions of the problem will be: (a) This solution with ε → 0 converges to the solution with a triangular membership function. (b) For cases where the core width (2ε) is large enough, the solution with a trapezoidal membership function will give a cut larger than the present solution with triangular membership function.

Gaussian Membership Function
First a fuzzy Gaussian number is considered as follows: µ G = exp[−(x − x) 2 /2σ 2 ]. Afterwards, a triangular fuzzy number is selected with membership function µ T and having the following properties: The base of triangular fuzzy number is equal to: 2c = 2σ √ 2π → σ = c/ √ 2π . In reference [43] it is proved that in this way the two membership functions almost coincide. Accordingly, using this method that transforms the Gaussian number to triangular one, it is proved that, when using the Gaussian membership function, the results remain the same.

Remark Concerning the Uncertainties of the Parameters
In this analysis, the focus was on the ambiguity of the initial conditions of the problem. It should be stressed, of course, that ambiguity is introduced into a linearized Boussinesq equation, which can provide information on solving the non-linear equation of Boussinesq by numerical methods. Also, there are more possible ambiguities on the parameters of hydraulic conductivity and effective porosity, since they are measurable quantities. The effect of these ambiguities should be considered and included for future research.

Conclusions
The [32] theory with the generalized Hukuhara (g-H) derivative, as well as its extension by [38] to partial differential equations, allows researchers to solve practical problems, which is useful in engineering. It is now possible for engineers to take the fuzziness of various sizes involved into consideration when calculating and making a decision on their work.
Boussinesq's linear equation regarding the aquifer case study has a fuzzy solution that is unique, and the function H(ξ, τ) is [(ii)-p] differentiable with respect to ξ and [(i)-pp] with respect to τ.

The function
∂ H(ξ,τ) ∂ξ is [(ii)-p] differentiable with respect to ξ. The fuzzy water volume variations (spreads) regarding their average value, amount to 15% of the average values of all times, and therefore, equal to the initial relative fuzziness in boundary values.
The fuzzy water flow velocity tends asymptotically to zero, while the recharging fuzzy water volume tends asymptotically to infinity.
It is important here to note, that for practical cases (artificial recharge problem, irrigation, and drainage projects), the engineers should take the right decision, knowing the deviations of the crisp value of water volume from the fuzzy ones, which here attains 15%.
Using trapezoidal membership functions, the difference is negligible if the core width is small. For larger values of core width, α-cuts intervals that are larger than the present solution arise. Using the Gaussian membership function, the results remain the same.