ANewNumericalMethod to Solve Some PDEs in the Unit Ball and Comparison with the Finite Element and the Exact Solution

In this paper, we give a new strategy to extend a numerical approximation method for two-dimensional reaction-diffusion problems. We present numerical results for this type of equations with a known analytical solution to qualify errors for the new method. We compare the results obtained using this approach to the standard finite element approach. +e proposed method is adequate even with the singular right-hand side of type Dirac.


Introduction
is work is motivated by the extension of a numerical approximation method constructed in [1,2] and called δ − ziti. Note that when the mathematical model is used to approximate a real concrete problem, most of the numerical schemes have undesirable oscillations, especially near the domain's boundary or near the physical phenomena, such as shock waves, relaxation waves, composite waves, blowup [3][4][5][6], and boundary layer which can exist in mathematical models, depending on the parameters involved, the nonlinearity [7,8], the coupling, the boundary, and the nature of the domain. Note that the link between different parameters can affect the degree of regularity of the solution (local, global, explodes in time, singular, etc.) [9][10][11].
e δ − ziti method has been constructed to exhibit such solutions regardless of their degree of regularity. is method has been tested on several models, which present or do not present singularities, and has been compared with exact solutions and classical methods (finite elements, finite volume, finite differences, particular method, and spectral method), but only in the case of a cubic domain (interval, square, cube, or even Cartesian domain in R d in the general case), especially the differential problems for which the application of other classical methods does not guarantee the transition from a regular behavior to another singular (i.e., giving a solution of the problem, by detecting the singularity when it exists. In this case, the method δ − ziti has shown its efficiency, robustness, and its convergence. It could detect the blowup, the boundary layer, and also it approaches the singular solution avoiding unwanted oscillations [1,2]. In this paper, we construct another strategy cited in [2], but for a domain more complicated than the Cartesian one, for example, a disk in the two-dimensional case, or any ball in the general case. e major goal of our investigation is to apply this numerical method in two types of PDE s , even if the right-hand side could present a singularity of type Dirac. In fact, despite the relatively simple model represented by the transient diffusion equation (especially with bad sign), its numerical solution can still be a serious challenge when the solution diffuses with steep boundary layers, or when we start with very sensitive initial data. e numerical study of the diffusion problems has been the aim of several papers. e third boundary-value problem for a loaded heat equation in a p-dimensional parallelepiped is considered in [12]. Considering that the problem is time-dependent, evaluating the solution process over refined meshes and for many hundreds or thousands of time steps can become a serious numerical burden. is can be a significant numerical challenge in two-dimensional problems and even more so in three dimensions [13][14][15][16]. Our paper mainly completes the investigations of [12][13][14][15][16][17][18], in the case of time-dependent diffusion problems. e difference is the use of a new numerical method, which gives us an admissible approximated solution. Considering that the diffusion problem is time-dependent, evaluating the solution process over refined meshes and for many hundreds or thousands of time steps can become a serious numerical burden. is can be a significant numerical challenge in two-dimensional problems and even more. In this work, we will present a numerical solution. Our objective is to solve such kind of problems with a very simple algorithm, which is a big profit, especially for the CPU time. e main goal of the δ − ziti method is to approach a function with several variables, to integrate it in a given domain, and to resolve numerically Partial and Ordinate Differential Equations (PDE s and ODE s ). is method is based on the classical variation formulation of Galerkin and the most important step is the construction of our orthonormal family from the test function Φ with compact support, defined by where x ∈ R d and d is the space dimension. is function is especially used in numerical analysis, distributions, and functional analysis. In [1,2,19], all the mathematical tools of δ − ziti construction in the multidimensional Cartesian case were detailed. e δ − ziti method is composed of two steps: it starts like the Galerkin (or spectral) method with a variational formulation and hence there is the need for an orthonormal basis (total or not) of an adequate space. en, we use the roots of the base instead of the nodes of the mesh, which will give δ − ziti a resemblance with the particle method. At first glance, the δ − ziti method looks as if it is the finite element method, but after using the roots, the method will follow the technique of the particular one. e main aim of this paper is the construction of the δ − ziti method when the domain is a disk in the twodimensional case (in general, a multidimensional ball). To generalize this method, we opt for two strategies: the first one consists in sweeping all the disk with segments, in the two directions, as shown in Figure 1 and reconstructing our basis functions in every segment. To test this strategy, we apply the resulting tools to calculate numerically the integrals and to solve partial differential equations (two tests will be detailed: an elliptic equation "the Poisson problem" and a parabolic one "the heat equation with bad sign"). e second strategy is a direct use of the polar parametrization of a disk; we will show that this strategy is also efficient and gives a good approximation even when a singularity is forced from the right-hand side term. e outline of this work is as follows: in Section 2, we present our approximation method in the monodimensional case and for several dimensions. In Section 3, we present the mathematical tools of construction, which permits the application of this method in the Cartesian case and the calculation, numerically, of some integrals defined in a disk domain. Section 4 is devoted to the construction of the method's fundamental elements, using polar coordinates. Like the previous section, one of the most important parts is the numerical integration using our method and in the two cases, we will compare the exact value of an integral by the numerical one, obtained by the δ − ziti method. In Section 5, we apply our approach to find the numerical solution of the Poisson problem and the heat equation. Our goal is to compare the solution obtained by the δ − ziti method with a given analytical one defined in a disk domain and to calculate the error in L ∞ (Ω). In the next, we present an approximated solution using the finite element method and we compare it with ours. To see that the method used in the current work is efficient, we treated two types of PDE s with an exact solution which converges to the Dirac mass (exact solution for the Poisson model and initial data for the heat transfer with bad sign). e function used is an approximation of the Dirac mass, and we will see that the δ − method can detect this kind of singularities (always in a disk domain). Finally, in Section 6, we finish with some concluding remarks.

Overview of the Monodimensional Construction
In this section, we will present the necessary mathematical tools for the construction of our method, in the monodimensional case.

Construction of the Linearly Independent Family
From the function Φ defined by (1), we define the set φ ε : where C\coloneq (1/ R Φ(x)dx) is a constant of normalization. is sequence φ ε converges to Dirac measure δ 0 in the sense of distributions. Fundamental results of construction in the monodimensional case are given as follows: First, we take a uniform mesh of [a, b] with the step h � (b − a)/m where m is a given integer such that We construct the family (φ i ) 1≤i≤m+1 as follows: Figure 1: Illustration of the elements of the family (φ i ).
e curves of φ i are shown in Figure 1.

Construction of the Orthonormal Family
Let us consider the Hilbert space L 2 (R) with the usual scalar product (, ). Observe that the family (φ i ) 1≤i≤m+1 is linearly independent, then using the Gram-Schmidt process, we construct a unique orthogonal family, noted as (Ψ i ) satisfying the following relation: Using the orthogonality of the family (Ψ i ) i , we can find the second definition of λ i given by and the recurrence application of definition (4) gives the following formula: Let Ψ i � (Ψ i /‖Ψ i ‖) the normalization of Ψ i . e curves of (Ψ i ) are shown in Figure 2.
We can verify the following assumption. For i � 2, . . . , m + 1, Ψ i admits (i − 1) reel root(s) noted as r k (k � 1, . . . , i − 1) in the interval ]x 1 , x i [, more precisely, e method permits approaching a given function f and its integral by the following relations: e main idea of the δ − ziti method is to use the roots r k of the functions (Ψ i ) which satisfy If we take x � r k in (7a), we obtain Relations (9a)-(9d) will be used frequently to construct the numerical scheme associated with the nonlinear system studied. To reduce the iterations number, the authors of [1] proved the following optimality result: where [.] denotes the integer part. In particular, for ε � 10 − M , we conclude that the parameters λ i are nearly stationary from a certain rank, which reduces considerably the number of iterations. Using r i as a root of Ψ i+1 , we can define the parameter λ i by λ i � − (φ i+1 (r i )/φ i (r i )) (for more details of the construction of the δ − ziti method, see [1,2].)

The First Strategy: Cartesian Coordinates
In this section, we are interested in the extension of the δ − ziti method, when Ω is a disk centered in the origin O � (0, 0) (or the ball in the multidimensional case). Our strategy is inspired by the monodimensional case. We sweep the inside of the domain by a set of intervals horizontally and vertically (see Figures 3 and 4 ). Note that the distance between every two successive chords is constant. We choose a fixed number of nodes to subdivide every interval. As the length of each chord varies according to the position of the segment inside the disk, International Journal of Differential Equations we will adapt our subdivision such that each interval has the same number of nodes already fixed. erefore, for each node, we define the associated basis function (denoted Ψ i ); after that, we compute the roots (denoted by r i ) of every element of this basis (Ψ i (r k ) � 0, for all i ≠ k; see Figure 5).

Remark 1.
Note that, for every fixed vertical level j (resp., the horizontal level i), every internal segment is limited by Figure 3). For simplicity, the step of the horizontal subdivision will be noted h j (resp., the vertical step will be noted h i ).
We present an algorithm to calculate the internal nodes (Algorithm 1).

Construction of the Orthonormal Set.
For every node x j i (resp., y i j ), we associate the function φ j i (noted as φ i if there is no ambiguity) (resp., the family φ i j will be noted as φ j ) defined by where h j is the step of construction in the horizontal interval of indication j, which describes the distance between the nodes x It is simple to see that the family (φ i ) is linearly independent, so we can construct an orthogonal family (Ψ i ) i�1,...,N by using the Gram-Schmidt process, in the space L 2 ([a, b]) (construction in every internal interval of the domain Ω � B(0, 1), horizontally and vertically), verifying the following relation:   International Journal of Differential Equations horizontally: which will be reduced in the following theorem, already proved in the monodimensional case (see [1,2]).

Theorem 1.
e orthogonal family (Ψ i ) i�1,...,N (vertically and horizontally) verifies the following recurrence relation: where (, ) is the usual scalar product in the Hilbert space

Corollary 1.
e family (Ψ i ) and the set (λ i ) defined in eorem 1 verify the following relations:  International Journal of Differential Equations 5 the horizontal and vertical test functions, as well as the roots, verifying the following relations: where Ψ j i (x) are the basis functions in the horizontal dimension (resp., Ψ i j (y) are the basis functions in the vertical dimension) r j i are the roots of Ψ j i (x) (resp., s i j are the roots of Ψ i j (y)) In this section, we are interested in the approximation of a double integral, defined in a disk domain; using (7a), we obtain the following results: Let Ω � B(0, 1), let g be a given function in L 2 (Ω), and let N denote the roots number; therefore, we have the following approximations: where we take r Proof. In [1,2], the authors approximated integral formulas in the one-dimensional case using the δ − ziti method as follows: We recall the principal points of the proof which is inspired by the spectral method: because of the following orthogonality result: Multiplying (18) by Ψ k and integrating it over Using the orthonormal property of (Ψ i ), and, therefore, In the bidimensional case, using the Cartesian strategy, where N is the nodes number in every segment, horizontally and vertically, we will have the following approximation formulas: Here, we take, r j N � b j and s N i � b i .
In Table 1, we present some numerical tests of integration. We compare the exact value noted E x , with the numerical approximation obtained by the δ − ziti method in the Cartesian case. Each interval is subdivided into 50 nodes.

Second Strategy: Polar Coordinates
In this section, we built all the necessary elements for the approximation δ − ziti, using polar coordinates. e domain Ω � B(0, 1) is represented using the polar coordinates, with the following parametrization:

Fundamental Results: Numerical Integration.
To test the previous strategy, we present in Tables 2 and 3 some numerical tests. We compare the exact value with the numerical one, using the δ − ziti method in the polar case.

Elliptic PDE Case: Nonevolutive Diffusion Problem in a
Unit Disk 5.1.1. e Cartesian Case. In this section, let us consider a partial differential equation, which admits an exact solution and we will compare it with the numerical one, using our method in the Cartesian case. Let Ω � B(0, 1). e problem studied is given by and, for simplicity of numerical simulations, we take K � I 2 , with a given analytical solution u ex � 1 − x 2 − y 2 , where the source term takes the value f � 4.

e Strong
Discretization. e first step to approach the previous problem is to multiply equation (7a) by a test function Ψ ij and to integrate the result over the domain Ω � B(0, 1), which gives Using eorem 2, we obtain the following equality: e next step consists of approaching the second derivative, which gives us the following scheme: where u ij denotes the approximation at the root (r j i , s i j ) and N is the nodes number in every internal segment (horizontally and vertically). In the end, we will have a global matrix, with (N − 2) × (N − 2) lines and (N − 2) × (N − 2) columns, defined as follows: where D i is a (N − 2) × (N − 2) tridiagonal matrix, defined by and A i is a (N − 2) × (N − 2) diagonal matrix defined by International Journal of Differential Equations 7 where and, therefore, we should resolve a simple system in the form MX � F, when M is the global matrix defined previously, X is the unknown vector of size (N − 2) × (N − 2), and F is the source term vector of size (N − 2) × (N − 2).

Remark 2.
To complete the resolution of the previous system, we must add boundary conditions (homogeneous Dirichlet in this case).

Numerical Results.
Let Ω � B(0, 1). In this case, we fix the points number in every single segment and we vary the subdivision step. It is clear that the minimum of all the steps is obtained at the first segment (horizontally or vertically) and the maximum is on the segment confused with the diameter of the disk (i.e., for two different intervals, horizontally or vertically, the associated step is not the same). We are interested in the shape of the approximated solution with the δ − ziti scheme, using the Cartesian coordinates and the segments approach. For a fixed node's number in every segment (horizontal or vertical), N � 100, the numerical implementation of the scheme gives us an approximated solution, which is near the exact one, given by u ex (x, y) � 1 − x 2 − y 2 , when f(x, y) � 4. e two solutions are illustrated in Figures 6 and 7.
In Table 4, we present the error between the exact and approximated solution, with different values of nodes number N.In Table 4, We remark that the committed error between the exact solution and the computed one using the δ − ziti method decreases when the number of points on each interval increases, which shows the convergence of the solution.
We present, in the following subsection, a comparison between approximated solutions using the finite element and δ − ziti methods.

Comparison with the Finite Elements Method.
e finite element method (FEM) is a widely used analogy to resolve some types of partial differential equations. A large class of works was already done to resolve the Poisson problem using FEM (see [20,21]). e starting point for the FEM is a PDE expressed in a variational form. e basic Table 2: Comparison between numerical integration using the δ − ziti and the exact value in the polar case.   Table 3: Generalized integrals using the polar parametrization.
2π(x 2 + y 2 ) )dxdy (− 0.500/2π) � 0.1590909091 0.1560909091 0.003 8 International Journal of Differential Equations recipe for turning a PDE into a variational problem is to multiply the equation by a test function v and to integrate the resulting expression over all the domain Ω: it is the common step between the Galerkin analogy and δ − ziti. In this part, we present the approximated solution of Poisson's problem defined in (7a), using the finite element method, which is illustrated in Figure 8.

e Polar
Case. Now, we consider the same partial differential equation defined before, which admits a polar analytical solution; we will compare it with the approximated one founded using the δ − ziti method in the polar case.
Let Ω � B(0, 1); the strategy presented consists in using the results of approximation in the monodimensional case and taking into consideration the following function basis: (37) e problem presented in the previous Subsection 5.1.1 is equivalent to the polar one, expressed as follows: with u(ρ, θ) ≔ u(ρ cos(θ), ρ sin(θ)), where f is the source term defined in Cartesian problem, with an exact polar solution u ex � 1 − ρ 2 . Note that, the roots of the basic functions (Ψ i (ρ)) will be noted; r i and θ j are those associated with (Ψ j (θ)).
To obtain a numerical scheme using the δ − ziti method, we should multiply equation (38a) by a test function Ψ ij (ρ, θ) and after that, we use the strong result of approximation (16), which gives   International Journal of Differential Equations erefore, the goal is to find u(ρ, θ) solution to the polar problem given in (38a) and (38b). Like the Cartesian analogy, we resolve in this case a simple system in the form MX � F; when M is the global matrix defined previously, we should just add the polar terms (1/ρ) and (1/ρ 2 ) in the corresponding terms of the matrices A i and D i , X in the polar unknown vector of size (N − 2) × (N − 2), and F is the source term vector of size (N − 2) × (N − 2).

Numerical Tests.
Let N r � 100 be the root's number for the radius and N θ � 100 the other one for the angles. Two solutions, exact and approximated, are illustrated in Figures 9 and 10 . Table 5 shows us the error between exact and approximated solution, using different N r and N θ . In Figure 11, we present the evolution of the error between exact and approximated solution, with several values of the nodes number N r for the diffusion problem. Table 5 shows how the error is reduced when the number of nodes increases from N r � 60 to N r � 200 (resp. N θ ) which is illustrated in Figure 11. e results obtained in this study suggest that having just 60 radial nodes leads to an efficient solution.

Detection of Blowup in the Stationary Diffusion Problem.
To evaluate the efficiency of the δ − ziti approach, we consider the following diffusion problem: where the source term is given by where C is the normalization constant defined by (1) and Φ is the regular function given by It is easy to see that the analytical solution of (42) and (43) is given by is function converges to Dirac mass in the weak sense when ε goes to 0. We present in Figure 12 the exact solution (resp., Figure 13 presents the approximated solution obtained by the scheme presented in (31)), to see that the δ − ziti method detects this kind of singularity.

Parabolic PDE Case: Heat Equation.
is section is devoted to the application of the δ − ziti method on a diffusion equation, in a domain Σ � [0, T] × Ω, where Ω � B(0, 1). e heat equation describes the distribution of heat (or variation in temperature) in a given region over time. For a function, u(t, x, y) (resp., u(t, r, θ)) of two spatial variables (x, y) in the Cartesian case ((r, θ) in the polar case) and the time variable t, the heat equation is given by where [0, T] is a given time interval and Ω � B(0, 1) is the unit disk with boundary zΩ. Here, x � (x, y) t denotes the space variables, t is the time variable, D is the diffusion coefficient, and f represents the effect of internal source terms. Using the same analogy applied in the previous sections, we multiply equation (46a) by a test function Ψ ij , after we integrate over the domain Ω. It remains just the direct application of our approximations formulas given in (3.3). First, the time domain [0, t max [ is divided into N t uniform subintervals [t n , t n+1 ] with duration dt � t n+1 − t n for n � 0, . . . , N t . To denote the value of u in the root r ij at time t n , we use the notation u n i,j for i, j � 1, . . . , N. At t � t n , we find directly the following explicit scheme: where Our goal is to find an approximated solution, near the exact one, which verifies the boundary conditions. e function is an exact solution of the heat equation, with (50) e above example provides a starting point for numerically solving diffusion problems on the disk. ese ideas can also be used on the sphere. Figures 14 and 15 show the allure of the exact and approximated solution, at a given finite time. Note that, the stability condition CFL is numerically well verified.
In Table 6, we present the error between exact and approximated global solution at a given finite time. Note that the CPU time necessary to obtain an approximated solution at t � 0.01 (when the number of nodes in every interval, horizontally and vertically, is 200) was around 150 seconds, which means that the method is fast and robust and gives us the numerical results in a very short period of time. Table 6 shows how the error is reduced when we refine the mesh from a number of nodes N � 60 to N � 200.
For checking the stability of our numerical scheme, we use the following strategy:   12 International Journal of Differential Equations in which u n (:, :) is the matrix of numerical solution at n th iteration. Errors obtained corresponding to E n dt using the δ − ziti method with N � 200 for the heat transfer problem are reported in Table 7.

Comparison with Finite Element Method.
Since the low order FEM is often used in industrial applications to solve heat transfer problems, comparisons between the linear FEM and the δ − ziti are also included in this section. We present the approximated solution given by the finite element method. In this direction, a large body of works was already done; see [20]. e solution of the heat transfer using the finite element method at t � 1 is illustrated in Figure 16.

Detection of Blowup in the Heat Transfer with Bad Sign.
In the current study, we are interested in the two-dimensional transient diffusion problem with bad sign (negative diffusion coefficient), governed by the following boundary-value equation: where [0, T] is a given time interval and Ω � B(0, 1) is the unit disk with boundary zΩ. Here, x � (x, y) t denotes the space variables, t is the time variable, D is a negative diffusion coefficient, and f represents the effect of internal source terms, given by the following equality:  We can verify that the exact solution is given by which converges to Dirac in the weak sense when ε goes to 0 as shown in Figure 17. e objective is to follow the behavior of the approximated solution, obtained by the same scheme presented in (47)-(48), which is illustrated in Figure 18.

Conclusion
As a conclusion, the current study is an extension of the δ − ziti method when the domain of consideration is a disk, to be able to generalize it in the case of a multidimensional ball in R d . We started by presenting the fundamental tools of construction using two strategies: the Cartesian case, which is based on the hypothesis of being able to sweep the inside of the domain by a set of intervals, horizontally and vertically; therefore, all the work resides in the construction of the nodes in every interval. e second strategy is the polar one, using the polar coordinates in all the steps of construction. e first application in this work is the generalized integration; we obtained a good approximation even if some of the classical methods did not give any result. We applied our method for two types of PDE s . Where the exact solution is regular, we also compared our results with the finite element method, and the result was impressive. On the other hand, the right-hand side is obtained such that the exact solution presents a singularity of type Dirac, where the parameter of regularization ε is near 0 and less than the step of discretization. For the stationary diffusion problem, we compare the exact solution which is a very sensitive function approximating the Dirac mass, with the approximated one which goes also to Dirac. e numerical problem faced in this case is how to find a numerical scheme which detects this type of singularity and gives a solution having the same regular behavior of the exact one: our method can achieve this objective. For the heat transfer, we treated two types of singularities: the negative diffusion parameter (bad sign) and the case when the initial data is in function of φ ε to test the ability to detect such types of singularity. e case of a negative diffusion presents a great difficulty numerically, since most of the numerical schemes tend to explode, even if the domain is Cartesian and therefore the solution is lost. e Improved CFL condition saves time for advanced previsions (CFL around 0.9, which requires a fewer number of iterations to reach the final time). We can conclude that, with the extension and development of δ − ziti, we have reached our goal: to approach differential problem defined in a disk, to find the results of the classical methods in the best conditions (speed, flexibility), even in the presence of highly singular terms (the numerical scheme detect the singularity, and solution continues regularly). Note that the maximum errors are calculated, finding that the errors are small. Also, a comparison of numerical and analytical solutions is made, finding that the proposed scheme has good accuracy.

Data Availability
No data were used to support this study.

Conflicts of Interest
e authors declare that they have no conflicts of interest.