Extension of δ_-ziti method in the unit ball: Numerical integration, resolution of Poisson's problem and Heat transfer

1. Introduction

The chemotactic dynamics of a population requires several steps, particularly, aggregation and blow-up. The Keller-Segel model describes this phenomena. It was suggested by Patlak in 1953 [3], Keller-Segel in 1970 [4], which allows for both diffusion and aggregation: depending on the initial data, the solution might exist globally in time or blow up in finite time, depending on the balance of forces between the different parameters involved in the system, the blow-up phenomenon may or may not occur. In fact, the blow-up is a singular behaviour of a Dirac solution. Most of numerical methods ( e.g. Galerkin, Particular method, spectral method..) does not ensure the transition from a regular behaviour to another singular one (i.e. detection of blow-up ). Under certain formulations of the Keller-Segel model, the phenomenon of aggregation has been shown to lead to finite-time blow-up. A large body of works has been devoted to determine when blow-up occurs or whether globally solutions exists: Authors of [5] developed a family of new interior penalty discontinuous Galerkin methods for solving the Keller–Segel chemotaxis model. In [6] they investigated non-negativity of exact and numerical solutions to a generalized Keller–Segel model where this model includes the so-called minimal Keller–Segel model. The main aim of [7] is to study the Keller–Segel model of chemotaxis and to develop a composite particle-grid numerical method with adaptive time stepping which allows to resolve and propagate singular solutions. The purpose of [8] is to formulate a phenomenological model from which the existence and properties of migrating bands can be deduced. Authors of [9] detected the blow-up as $δ$ -function (amoebae aggregation) at the proximity of the origin in dimension one and on a ball in a multidimensional space. Therefore, it was necessary to find a new numerical scheme, which detect this type of singularities easily, without loosing the advantages of classical methods. The $δ_{-} z i t i$ method is on the challenge. It was tested on several type of problems (see [1] and [2]), including the Keller-Segel model, but only in the Cartesian case ( segments, rectangle, cube $\dots$ ).The main goal of $δ_{-} z i t i$ 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 (PDEs and ODEs). This method is based on the classical variation formulation of Galerkin and the most important step, is the construction of our orthonormal family, from the famous function $Φ$ $\in$ $D (Ω)$ with compact support, defined by

Φ (x) = {\begin{matrix} exp (\frac{1}{∣ x ∣^{2} - R^{2}}) & si ∣ x ∣< R 0 & otherwise, \end{matrix}

(1)

where R $>$ 0, $Ω \subset R^{n}$ and $x \in R^{n}$ .
This function is used especially in numerical analysis, distributions and functional analysis. It is characterized by giving the best approximation of the Dirac measure. In [1] and [2], the multi-dimensional Cartesian case was detailed. The main aim of this paper is the construction of $δ_{-} z i t i$ method when the domain is a disk in the two-dimensional case, (in general, a multi-dimensional 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 to reconstruct our basis functions in every segment, which means that we inject all the work already done in the mono-dimensional case. To test this strategy, we apply the resulting tools to calculate numerically 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”). The second strategy is a direct use of the polar parametrisation of a disk, we will show that this strategy is also efficient and gives us a good approximation ( Integrals ans two tests of resolving PDEs).

The outline of this paper is as follows. In section 3, we will present the mathematical tools of construction, which permits to apply this method in the Cartesian case, as shown in figure 1, to calculate, numerically, some integrals defined in a disk domain.The 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 $δ_{-} z i t i$ method.In the last section 5, we apply our approach to find the numerical solution of Poisson problem and the Heat equation. Our goal is to compare the solution obtained by $δ_{-} z i t i$ method, with a given analytical one defined in a disk domain and to calculate the error in $L^{\infty} (Ω)$ . By the next, we present an approximated solution using the finite element method and we compare it with our one.

2. Overview of the mono-dimensional construction.

All the results presents in this section, are proved in [1] and [2]. The fundamental results of construction are given as follow:First, we take a uniform mesh of [a,b] with the step $h = \frac{b - a}{N}$ , where N is an integer such that $x_{i} = a + (i - 1) h, \forall i \in [1, N + 1] .$
From the function $Φ$ defined in (1), we define $φ_{ϵ}$ by:

φ_{ϵ} (x) = \frac{c}{ϵ} Φ (\frac{x}{ϵ}), for all ϵ > 0,

where c:= $\frac{1}{\int_{R} Φ (x) d x}$ is the constant of normalization.

This sequence $φ_{ϵ}$ converges to Dirac in the sense of distributions, which is often used to detect singularities.

We construct the family $(φ_{i})_{i = 1 \dots N + 1}$ as follows:

⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ \begin{matrix} φ_{i} (x) = φ_{h} (x - x_{i}) = \frac{}{C} h Φ (\frac{x - x_{i}}{h}), for all x \in [x_{i - 1}, x_{i + 1}], i \in [2, N], φ_{1} (x) = φ_{h} (x - x_{1}) = \frac{C}{h} Φ (\frac{x - x_{1}}{h}), for all x \in [x_{1}, x_{2}], φ_{N + 1} (x) = φ_{h} (x - x_{N + 1}) = \frac{C}{h} Φ (\frac{x - x_{N + 1}}{h}), for all x \in [x_{N}, x_{N + 1}] . \end{matrix}

Let consider the Hilbert space $L^{2} (R)$ , with the usual scalar product $(,)$ . Observe that the family ${(φ_{i})}_{1 \leq i \leq N + 1}$

is linearly independent, then using the Gram-Schmidt process, we construct a unique orthogonal family, noted

({~ Ψ}_{i})

satisfying the following relation

⎧ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎩ \begin{matrix} {~ Ψ}_{i} (x) = φ_{i} (x) + λ_{i - 1} {~ Ψ}_{i - 1} (x), λ_{1} = - \frac{α}{β}, λ_{i + 1} = g (λ_{i}), \end{matrix}

(2)

with

g (X) = \frac{λ_{1}}{2 - λ_{1} X}, α = (φ_{1}, φ_{2}), β = (φ_{1}, φ_{1}) .

The spectral method applied to find the direct formula of the basis functions gives,

λi=−(φi,ψi−1)(ψi,ψi−1)\raisebox2.15pt,

and the recurrence application of the definition given in (2) gives the following formula:

{~ Ψ}_{i} (x) = φ_{i} (x) + λ_{i - 1} φ_{i - 1} (x) + λ_{i - 1} λ_{i - 2} φ_{i - 2} (x) + \dots + 1 \prod k = i - 1 λ_{k} φ_{1} .

Let $Ψ_{i} = \frac{{~ Ψ}_{i}}{| |_{i} | |}$ the normalization of ${~ Ψ}_{i}$ ( for more details see [1] and [2]).
The method permits to approach a given function $f$ and an integral, using the following relations:

\begin{matrix} f (x) ≃ N \sum i = 1 c_{i} Ψ_{i} (x), c_{i} ≃ \int_{a}^{b} f (x) Ψ_{i} (x) d x, \int_{a}^{b} f (x) d x ≃ N \sum i = 1 c_{i} I_{i}, \end{matrix}

(3)

where $I_{i} := \int_{a}^{b} Ψ_{i} (x) d x$ . If we take $x = r_{k}$ in (3), we obtain:

ci≃f(ri)Ψi(ri)\raisebox2.15pt,i=1⋯N−1,cN≃f(b)ΨN(b)\raisebox2.15pt,∫baf(x)Ψi(x)dx≃f(ri)Ψi(ri)\raisebox2.15pt,∫baf(x)dx≃N−1∑i=1f(ri)Ψ2i(ri)+f(b)ΨN(b)2\raisebox2.15pt.

(4)

To reduce the iterations number, in [11] they proved that,

| λ_{i + 1} - λ_{i} | < ϵ as soon as i \geq N_{0} = ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ \frac{ln (\frac{ϵ (2 - λ_{1}^{2})}{λ_{1}^{3} - λ_{1}})}{ln (\frac{λ_{1}}{2 + λ_{1}})^{2}} ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ + 1,

where $[.]$ denotes the floor function. In particular for $ϵ = 10^{- M}$ , we concluded that the parameter $λ_{i}$ is 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} = - \frac{φ_{i + 1} (r_{i})}{φ_{i} (r_{i})} .$

3. The first strategy: Cartesian coordinates.

3.1. Construction of intern nodes

In this section, we are interested in the extension of $δ_{-} z i t i$ method, when $Ω$ is a disk centred in the point $O = (0, 0)$ (or the ball in the multi-dimensional case). As a first step, we start by a general presentation of this new strategy.We present the important steps of the construction, inspired by the mono-dimensional case. For this, suppose that we can 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. (see Figure 1)

$a_{j} = - \sqrt{1 - y_{j}^{2}}, b_{j} = - a_{j}$ and $c_{i} = - \sqrt{1 - x_{i}^{2}}, d_{i} = - c_{i}$

An algorithm which explains the steps of the construction will be presented in the next part . The main idea is to fix the number of nodes in every interval and to vary the step of subdivision associated to every interval (horizontally and vertically). Every node is noted by $(x_{i}^{j}, y_{j}^{i})$ , (see Figure 1).

Remark 3.1.

Note that, for every fixed vertical level $j$ (respectively the horizontal level $i$ ), every internal segment is limited by $a_{j} = - \sqrt{1 - y_{j}^{2}}$ and $b_{j} = \sqrt{1 - y_{j}^{2}}$ , (respectively, $c_{i} = - \sqrt{1 - x_{i}^{2}}$ and $d_{i} = \sqrt{1 - x_{i}^{2}}$ ).
For simplicity, the step of horizontal subdivision will be noted $h_{j}$ (respectively, the vertical step will be noted $h_{i}$ ).

Here we present an algorithm to calculate the internal nodes.

Data:

N =

The number of nodes, in every interval

Fix an interval [a,b]

h^{1} = \frac{b - a}{N}

;

for $i = 1 \dots N$ do

x_{i}^{1} = a + (i - 1) h^{1}

horizontal nodes for the first intervall.

y_{i}^{1} = a + (i - 1) h^{1}

vertical nodes for the first intervall

end for

for $j = 2 \dots N$ do

for $i = 1 \dots N$ do

a=-

\sqrt{1 - (- 1 + j . h)^{2}}

b=-a

h^{j} = \frac{b - a}{N}

x_{i}^{j} = a + (i - 1) h^{j}

end for

Algorithm 1 Construction of the nodes in the Cartesian case.

3.2. Construction of the orthonormal set

For every node $x_{i}^{j}$ (respectively, $y_{j}^{i}$ ) we associate the function $φ_{i}^{j}$ , ( noted $φ_{i}$ if there is no ambiguity) (respectively ,the family $φ_{j}^{i}$ will be noted $φ_{j}$ ) defined by:

\begin{matrix} φ_{i} (x) := \frac{c}{h_{j}} Φ (\frac{x - x_{i}^{j}}{h^{j}}), \forall i, j = 1 \dots N, φ_{j} (y) := \frac{c}{h_{i}} Φ (\frac{y - y_{j}^{i}}{h^{i}}), \forall j, i = 1 \dots N, \end{matrix}

(5)

where,
$∙$ $h^{j}$ is the step of construction in the horizontal interval of indication $j$ , which describe the distance between the node $x_{i}^{j}$ and $x_{i + 1}^{j}$ .
$∙$ $h_{i}$ is the step of subdivision in the vertical interval of indication $i$ , which describe the distance between $y_{j}^{i}$ and $y_{j + 1}^{i}$ .

It is simple to see that the family $(φ_{i})$ is linearly independent, so we can construct an orthogonal family $(_{i})_{i = 1 \dots 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,

H o r i z o n t a l l y : ⎧ ⎪ ⎪ ⎨ ⎪ ⎪ ⎩ \begin{matrix} {~ Ψ}_{1} (x) & = φ_{1} (x) {~ Ψ}_{i} (x) & = φ_{i} (x) + i - 1 \sum k = 1 λ_{k}^{(i)} {~ Ψ}_{k} (x), for all i = 2, \dots, N, \end{matrix}

V e r t i c a l l y : ⎧ ⎪ ⎪ ⎨ ⎪ ⎪ ⎩ \begin{matrix} {~ Ψ}_{1} (y) & = φ_{1} (y) {~ Ψ}_{i} (y) & = φ_{i} (y) + i - 1 \sum k = 1 λ_{k}^{(i)} {~ Ψ}_{k} (y), for all i = 2, \dots, N, \end{matrix}

which will be reduced in the following theorem, already proved in the mono-dimensional case, (see [1] and [2]).

Theorem 3.1.

The orthogonal family $(_{i})_{i = 1 \dots N}$ (vertically and horizontally), verify the following recurrence relation:

⎧⎪ ⎪ ⎪ ⎪ ⎪⎨⎪ ⎪ ⎪ ⎪ ⎪⎩~Ψ1=φ1~Ψi+1=φi+1+λi~Ψi,    for all  i=1⋯N−1,λi−1=−(φi,φi−1)(~Ψi−1,~Ψi−1)\raisebox2.15pt,

(6)

where $(,)$ is the usual scalar product in the Hilbert space $L^{2} ([a, b])$

Corollary 3.2.

The family $({~ Ψ}_{i})$ and the set $(λ_{i})$ defined in the theorem (3.1), verify the following relations:

\begin{matrix} 1) {~ Ψ}_{i} = φ_{i} + λ_{i - 1} φ_{i - 1} + λ_{i - 1} λ_{i - 2} φ_{i - 2} + \dots + λ_{i - 1} \dots λ_{1} φ_{1} . 2) {~ Ψ}_{i} (x_{i}^{j}) = φ_{i} (x_{i}^{j}) = \frac{c}{h_{j}^{2} e} . 3) {~ Ψ}_{i} (y_{i}^{j}) = φ_{i} (y_{i}^{j}) = \frac{c}{h_{i}^{2} e} . 4) In every fixed level, (φ_{i}, {~ Ψ}_{i - 1}) = (φ_{i}, φ_{i - 1}) . 5) - 1 < λ_{i} = - \frac{(φ_{i}, φ_{i + 1})}{(_{i},_{i})} < 0. \end{matrix}

(7)

3.3. Fundamental results: Numerical integrations

In this paragraph, we are interested by the approximation of integrals, where the domain is the unit disk $Ω = B (0, 1)$ , using the horizontal and vertical test functions, as well as the roots, verifying the following relations:

\begin{matrix} Ψ_{i j} (x, y) = Ψ_{i}^{j} (x) . Ψ_{j}^{i} (y), \forall i, j = 1 \dots N, r_{i j} = (r_{i}^{j}, s_{j}^{i}), \end{matrix}

(8)

where,
$∙$ $Ψ_{i}^{j} (x)$ are the basis functions in the horizontal dimension, (respectively, $Ψ_{j}^{i} (y)$ are the basis functions in the vertical dimension).
$∙$ $r_{i}^{j}$ are the roots of $Ψ_{i}^{j} (x)$ (respectively, $s_{j}^{i}$ are the roots of $Ψ_{j}^{i} (y)$ ).
In this section, we are interested by the approximation of a double integral, defined in a disk domain, using (3), we obtain the following results:

Theorem 3.3.

Let $Ω = B (0, 1)$ , $g$ a given function in $L^{2} (Ω)$ , and N denotes the roots number, therefore, we have the following approximations:

∫Ωg(x,y) dxdy≃N∑i,j=1g(rji,sij)Ψji(rji).Ψii(sij)∫bjajΨji(x)dx∫diciΨij(y)dy,∫Ωg(x,y)Ψij(x,y) dxdy≃g(rji,sij)(Ψji(rji).Ψii(sij))2\raisebox2.15pt.

(9)

here we take, $r_{N}^{j} = b_{j}$ and $s_{i}^{N} = b_{i}$ .

In the table (1), we present some numerical tests of integration. We compare the exact value, with the numerical approximation, obtained by $δ_{-} z i t i$ method, in the Cartesian case.

$\int_{Ω} f (x, y) d x d y$	$2 π .$ Ex	$2 π .$ Ap	$2 π .$ Error
$\int_{Ω} (\frac{1}{x^{2} + y^{2} + 1})^{\frac{1}{4}}$	$0.4545285537$	$0.454496459918650$	$0.0000320937813496069$
$\int_{Ω} \sqrt{(\frac{1}{x^{2} + y^{2} + 1})}$	$0.4142135624$	$0.414440692467526$	$0.0002271300675258380$
$\int_{Ω} exp (\frac{1}{x^{2} + y^{2} + 2})$	$0.7508533738$	$0.750772037320043$	$0.0000813364799567839$
$\int_{Ω} ln (\frac{1}{x^{2} + y^{2} + 2})$	- $0.4547712524$	- $0.453623297839054$	$0.0011479545609460200$

Table 1. Comparison between numerical integration using

δ_{-} z i t i

and the exact value, in the Cartesian case.

4. Second strategy: Polar coordinates.

In this section, we built all the necessary elements for the approximation $δ_{-} z i t i$ , using polar coordinates. The domain $Ω = B (0, 1)$ is represented using the polar coordinates , with the following parametrization:

\forall (x, y) \in B (0, 1), x = r cos (θ), y = r sin (θ), (r, θ) \in [0, 1] \times [0, 2 π] .

4.1. Construction of the method’s tools using polar coordinates

This first part of the algorithm, compute $λ_{i}$ and $r_{i}$ in $[0, 1] .$
The construction’s algorithm using the two variables $r$ and $θ$ is given as follows::

Data:

N_{r} =

the root’s number in the interval

[0, 1]

h_{r} = \frac{1}{N_{r} - 1}

;

for $i = 1 \dots N_{r}$ do

x_{i} = a + (i - 1) h_{r}

φ_{i} (r) = Φ (\frac{r - x_{i}}{h_{r}})

end for

α = \int_{0}^{1} φ_{1} (r) φ_{2} (r) d r

β = \int_{0}^{1} φ_{1} (r) φ_{1} (r) d r

λ_{1} = - \frac{α}{β}

for $i = 1 \dots N_{r} - 1$ do

λ_{i + 1} = - \frac{λ_{1}}{2 - λ_{1} λ_{i}}

end for

_{1} (r) = φ_{1} (r)

for $i = 2 \dots N_{r}$ do

_{i} (r) = φ_{i} (r) + λ_{i - 1} {~ Ψ}_{i - 1} (r)

end for

for $i = 1 \dots N_{r}$ do

Ψ_{i} (r) = \frac{{~ Ψ}_{i} (r)}{∣ ∣ {~ Ψ}_{i} (r) ∣ ∣}

end for

for $i = 1 \dots N_{r}$ do

Λ_{i} = ln (- λ_{i})

P (y) = Λ_{i} y^{4} - 2 Λ_{i} y^{3} - Λ_{i} y^{2} + 2 (Λ_{i} - 1) y + 1 = 0

P (y_{i}^{*}) = 0

r_{i} = x_{i} + h_{r} . y_{i}^{*}

end for

Algorithm 2 Construction of

δ_{-} z i t i

’s tools, to compute the roots

r_{i}

The second part of the algorithm permits to compute $λ_{i}^{θ}$ and $s_{i}$ in $[0, 2 π] .$

Data:

N_{θ} =

the root’s number in the interval

[0, 2 π]

h_{θ} = \frac{2 π}{N_{θ} - 1}

;

for $i = 1 \dots N_{θ}$ do

θ_{i} = a + (i - 1) h_{θ}

φ_{i} (θ) = Φ (\frac{θ - θ_{i}}{h_{θ}})

end for

α^{θ} = \int_{0}^{2 π} φ_{1} (θ) φ_{2} (θ) d θ

β^{=} \int_{0}^{2 π} φ_{1} (θ) φ_{1} (θ) d θ

λ_{1}^{=} - \frac{α^{θ}}{β^{θ}}

for $i = 1 \dots N_{θ} - 1$ do

λ_{i + 1}^{=} - \frac{λ_{1}^{θ}}{2 - λ_{1}^{θ} λ_{i}^{θ}}

end for

_{1} (θ) = φ_{1} (θ)

for $i = 2 \dots N_{θ}$ do

_{i} (θ) = φ_{i} (θ) + λ_{i - 1}^{θ} {~ Ψ}_{i - 1} (θ)

end for

for $i = 1 \dots N_{θ}$ do

Ψ_{i} (θ) = \frac{{~ Ψ}_{i} (θ)}{∣ ∣ {~ Ψ}_{i} (θ) ∣ ∣}

end for

for $i = 1 \dots N_{θ}$ do

Λ_{i} = ln (- λ_{i}^{θ})

P (y^{θ}) = Λ_{i} y^{4} - 2 Λ_{i} y^{3} - Λ_{i} y^{2} + 2 (Λ_{i} - 1) y + 1 = 0

P (y_{i}^{*^{θ}}) = 0

s_{i} = θ_{i} + h_{θ} . y_{i}^{*^{θ}}

end for

Algorithm 3 Construction of

δ_{-} z i t i

’s tools, to compute the roots

s_{i}

The polar set $(Ψ_{i j} (r, θ))$ is defined by,

Ψ_{i j} (r, θ) := Ψ_{i} (r) . Ψ_{j} (θ) .

(10)

4.2. Fundamental results: Numerical Integration

To test the previous strategy, we present in the following table, some numerical tests. We compare the exact value with the numerical one, using $δ_{-} z i t i$ method in the polar case.

$\frac{1}{2 π} \int_{Ω} f (x, y) d x d y$	Ex	Ap	Error
$\int_{Ω} (\frac{1}{x^{2} + y^{2} + 1})^{\frac{1}{4}} d x d y$	$0.4545285537$	$0.454496459918650$	$0.0000320937813496069$
$\int_{Ω} \sqrt{(\frac{1}{x^{2} + y^{2} + 1})} d x d y$	$0.4142135624$	$0.414440692467526$	$0.0002271300675258380$
$\int_{Ω} exp (\frac{1}{x^{2} + y^{2} + 2}) d x d y$	$0.7508533738$	$0.750772037320043$	$0.0000813364799567839$
$\int_{Ω} ln (\frac{1}{x^{2} + y^{2} + 2}) d x d y$	$0.4547712524$	- $0.453623297839054$	$0.0011479545609460200$
$\int_{Ω} x y d x d y$	0	$0.000104893284924629$	$0.0001$
$\int_{Ω} \frac{ln (\sqrt{x^{2} + y^{2}})}{\sqrt{x^{2} + y^{2}}} d x d y$	-1	$- 1.00011299531961$	$0.000112995319609510$

Table 2. Comparison between numerical integration using

δ_{-} z i t i

and the exact value, in the polar case.

Using the polar coordinates $r \in [0, 1]$ and $θ \in [0, 2 π]$ for some types of integrals, the following table shows us the error between exact and approximated solution founded using $δ_{-} z i t i$ method.

$\int_{Ω} f (r, θ) d r d θ$	Ex	Ap	Error
$\int_{Ω} \frac{r . sin (θ)}{(r^{2} + 1) . t^{\frac{1}{3}}}$	$0.2204366348$	$.223358184762906$	$0.00292154976290612$
$\int_{Ω} \frac{sin (θ)}{r^{\frac{1}{2}}}$	$0$	$0.000796466856449210$	$0.000796466856449210$
$\int_{Ω} \frac{1}{\sqrt{2 π r}}$	$\frac{1}{π}$	$\frac{0.994955526251694}{π}$	$0.00160570586425771$
$\int_{Ω} \frac{ln (r)}{\sqrt{2 π}}$	$\frac{- 0.500}{2 π}$	$\frac{- 0.511242350773147}{2 π}$	$0.00357855139409662$

Table 3. Generalised integrals expressed with polar variables

In the previous table, we remark that even we choose a generalised integral, like the example $\int_{Ω} \frac{sin (\frac{θ}{2})}{x^{\frac{1}{2}}}$ , which is in fact an operation of Riemann integral, we found a good approximation using $δ_{-} z i t i$ roots. For the two last examples, $\int_{Ω} \frac{1}{\sqrt{2 π r}}$ and $\int_{Ω} \frac{ln (r)}{\sqrt{2 π}}$ , other approximation methods (e.g. Simpson, Trapeze..) didn’t give any result, which is an important point for our construction.

5. Numerical applications.

5.1. Elliptic PDE case : Poisson problem

5.1.1. The Cartesian case

In this section, let 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)$ . The problem studied is given by:

		$- Δ u = f i n Ω,$			(11a)
		$u (x, y) = 0 i n \partial Ω,$			(11b)

with a given analytical solution $u_{e x} = 1 - x^{2} - y^{2}$ , and the source term function is defined by $f = 4$ .

The strong discretization

The first step to approach the previous problem, is to multiply the equation (11a) by a test function $Ψ_{i j}$ and to integrate the result over the domain $Ω = B (0, 1)$ , which gives,

- \int_{Ω} Δ u (x, y) . Ψ_{i j} (x, y) d x d y = \int_{Ω} f (x, y) . Ψ_{i j} (x, y) d x d y .

(12)

Using the theorem (3.3), we obtain the following scheme:

- \frac{Δ u (r_{i}^{j}, s_{j}^{i})}{Ψ_{i j} (r_{i}^{j}, s_{j}^{i})} = \frac{f (r_{i}^{j}, s_{j}^{i})}{Ψ_{i j} (r_{i}^{j}, s_{j}^{i})}, \forall (r_{i}^{j}, s_{j}^{i}) \in Ω .

(13)

The next step, consists to approach the second derivative, which gives us the following scheme:

⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ \begin{matrix} \frac{u_{i - 1, j} + 2 u_{i j} + u_{i + 1, j}}{(r_{i + 1}^{j} - r_{i}^{j}) (r_{i}^{j} - r_{i - 1}^{j})} + \frac{u_{i, j - 1} + 2 u_{i j} + u_{i, j + 1}}{(s_{j + 1}^{i} - s_{j}^{i}) (s_{j}^{i} - s_{j - 1}^{i})} = f (r_{i}^{j}, s_{j}^{i}), i, j = 2 \dots N - 1, u_{1, j} = u_{M, j} = 0, j = 1 \dots N, u_{i, 1} = u_{i, M} = 0, i = 1 \dots N, \end{matrix}

(14)

where, $N$ is the nods number in every internal segment (horizontally and vertically). At the end, we will have a global matrix, with $(N - 2) \times (N - 2)$ lines and $(N - 2) \times (N - 2)$ columns, defined as follows::

M = ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ \begin{matrix} D^{2} & A^{3} & 0 & \dots & 0 A^{2} & D^{3} & A^{4} & \dots & 0 ⋮ 0 & \dots & A^{n - 3} & D^{n - 2} & A^{n - 1} 0 & 0 & \dots & A^{n - 2} & D^{n - 1} \end{matrix} ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦

where $D^{i}$ is a $(N - 2) \times (N - 2)$ tri-diagonal matrix, defined by:

	$D_{k, k}^{i + 1} = - \frac{2 Ψ_{i, k + 1}}{d x_{k} . d x_{k - 1}} - \frac{2 Ψ_{i, k + 1}}{d y_{k} . d y_{k - 1}},$	$i, k = 1 \dots N - 2,$
	$D_{k, k + 1}^{i + 1} = - \frac{2 Ψ_{i, k + 2}}{d y_{k} . d y_{k - 1}},$	$i, k = 1 \dots N - 2,$
	$D_{k - 1, k}^{i + 1} = - \frac{2 Ψ_{i, k}}{d y_{k} . d y_{k - 1}},$	$i, k = 1 \dots N - 2,$

and $A^{i}$ is a $(N - 2) \times (N - 2)$ diagonal matrix defined by:

A_{k, k}^{i} = - \frac{Ψ_{i, k}}{d x_{k} . d x_{k - 1}}, i = 3 \dots N - 1, a n d k = 2 \dots N - 1,

where,

	$dxk=rjk−rjk−1\raisebox2.15pt,$
	$dyk=rik−rik−1\raisebox2.15pt,$

therefore, we should resolve a simple system in the form $M X = F$ , when $M$ is the global matrix defined previously, $X$

in the unknown vector of size

(N - 2) \times (N - 2)

and

F

is the source term vector of size

(N - 2) \times (N - 2)

Remark 5.1.

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

5.1.2. 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 by the shape of the approximated solution with $δ_{-} z i t i$ 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 to the exact one, given by $u_{e x} (x, y) = 1 - x^{2} - y^{2},$ and $f (x, y) = 4$ .

Figure 3. Approximated solution using $δ_{-} z i t i$

In the table 4, we present the error between exact and approximated solution, with different values of nodes number $N$ .

$N$	$h_{m i n}$	$h_{m a x}$	$E r_{m a x}$	$E r_{m e a n}$
$60$	$0.008534$	$0.03333$	$0.01174$	$0.00758$
$100$	$0.003979$	$0.02000$	$0.0046$	$0.00288$
$150$	$0.00217$	$0.01333$	$0.00204$	$0.0013169$
$200$	$0.001410$	$0.01000$	$0.00167$	$0.0007514823$

Table 4. The committed error using several values of nodes number

N

where,

	$E r_{m a x}$	$:= max i, j (\| u_{e x} (i, j) - u (i, j)) \|,$
	$E r_{m e a n}$	$:= \frac{1}{N} \sum i, j \| u_{e x} (i, j) - u (i, j)) \| .$

We present in the following subsection, a comparison between approximated solutions using Finite Element and $δ_{-} z i t i$ methods.

5.1.3. Comparison with Finite elements method

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 [10] and [12]). The starting point for the FEM is a PDE expressed in variational form. The basic 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 Galerkin analogy and $δ_{-} z i t i$ . In this part, we present the approximated solution of the Poisson’s problem defined in (11a), using Finite Element Method.

Figure 4. Approximated solution using finite element method

5.1.4. The 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 $δ_{-} z i t i$ method in the polar case.

Let $Ω = B (0, 1)$ , the strategy presented consists in using the results of approximation in the mono-dimensional case and taking into consideration the following function basis:

Ψ_{i j} (r, θ) = Ψ_{i} (r) Ψ_{j} (θ), \forall i, j = 1 \dots N .

Figure 5. The domain $Ω$ using polar parametrization

The problem presented in the previous subsection 5.1.1 is equivalent of the polar one, expressed as follows:

		$- Δ u = f i n Ω,$			(15a)
		$u (r = 1, θ) = 0 \forall θ \in [0, 2 π],$			(15b)

with

	$Δ u := \frac{\partial^{2} u}{\partial r^{2}} + \frac{1}{r} \frac{\partial u}{\partial r} + \frac{1}{r^{2}} \frac{\partial^{2} u}{\partial θ^{2}},$
	$f (r, θ) = f (r cos θ, r sin θ) = 4$

where $f$ id the source term defined in Cartesian problem, with an exact polar solution $u_{e x} = 1 - r^{2}$ . Note that, the roots of the basic functions $(Ψ_{i} (r))$ will be noted $r_{i}$ and $θ_{j}$ are those associated with $(Ψ_{j} (θ))$ . To obtain a numerical scheme using $δ_{-} z i t i$ method, we should multiply the equation 15a by a test function $Ψ_{i j} (r, θ)$ and after, we use the strong result of approximation 3.3, which gives:

⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ \begin{matrix} Δ_{i j} u = f (r_{i}, θ_{j}), & i, j = 2 \dots N - 1, u_{1, j} = u_{2, j}, & j = 1 \dots N, u_{M, j} = 0, & j = 1 \dots N, u_{i, 1} = u_{i, 2}, & i = 1 \dots N, u_{i, M} = u_{i,} \end{matrix}

Extension of δ_-ziti method in the unit ball: Numerical integration, resolution of Poisson's problem and Heat transfer

Authors

Numerical Integration on Graphs: where to sample and how to weigh

Efficient Monte Carlo Method for Integral Fractional Laplacian in Multiple Dimensions

Poisson Integrators based on splitting method for Poisson systems

A C^1 Petrov-Galerkin method and Gauss collocation method for 1D general elliptic problems and superconvergence

Energetic Stable Discretization for Non-Isothermal Electrokinetics Model

Numerical computation of Neumann controllers for the heat equation on a finite interval

Geodesics in Heat

1. Introduction

2. Overview of the mono-dimensional construction.

3. The first strategy: Cartesian coordinates.

3.1. Construction of intern nodes

Remark 3.1.

3.2. Construction of the orthonormal set

Theorem 3.1.

Corollary 3.2.

3.3. Fundamental results: Numerical integrations

Theorem 3.3.

4. Second strategy: Polar coordinates.

4.1. Construction of the method’s tools using polar coordinates

4.2. Fundamental results: Numerical Integration

5. Numerical applications.

5.1. Elliptic PDE case : Poisson problem

5.1.1. The Cartesian case

The strong discretization

Remark 5.1.

5.1.2. Numerical results

5.1.3. Comparison with Finite elements method

5.1.4. The polar case

Extension of δ_-ziti method in the unit ball: Numerical integration, resolution of Poisson's problem and Heat transfer

Authors

Related Research

Numerical Integration on Graphs: where to sample and how to weigh

Efficient Monte Carlo Method for Integral Fractional Laplacian in Multiple Dimensions

Poisson Integrators based on splitting method for Poisson systems

A C^1 Petrov-Galerkin method and Gauss collocation method for 1D general elliptic problems and superconvergence

Energetic Stable Discretization for Non-Isothermal Electrokinetics Model

Numerical computation of Neumann controllers for the heat equation on a finite interval

Geodesics in Heat

This week in AI

1. Introduction

2. Overview of the mono-dimensional construction.

3. The first strategy: Cartesian coordinates.

3.1. Construction of intern nodes

Remark 3.1.

3.2. Construction of the orthonormal set

Theorem 3.1.

Corollary 3.2.

3.3. Fundamental results: Numerical integrations

Theorem 3.3.

4. Second strategy: Polar coordinates.

4.1. Construction of the method’s tools using polar coordinates

4.2. Fundamental results: Numerical Integration

5. Numerical applications.

5.1. Elliptic PDE case : Poisson problem

5.1.1. The Cartesian case

The strong discretization

Remark 5.1.

5.1.2. Numerical results

5.1.3. Comparison with Finite elements method

5.1.4. The polar case