(30 intermediate revisions by 2 users not shown)
Line 12: Line 12:
  
 
In this paper, it is presented a formulation of a generalized finite difference scheme to solve the Motz problem. It is based on a general difference scheme defined by an optimality condition, which has been developed to solve Poisson-like equations whose domains are approximated by a wide variety of grids over general regions. Numerical examples showing second-order accuracy of the calculated solutions are presented.
 
In this paper, it is presented a formulation of a generalized finite difference scheme to solve the Motz problem. It is based on a general difference scheme defined by an optimality condition, which has been developed to solve Poisson-like equations whose domains are approximated by a wide variety of grids over general regions. Numerical examples showing second-order accuracy of the calculated solutions are presented.
 +
 +
'''Keywords''': Generalized finite difference, Motz problem
  
 
==1. The Motz problem==
 
==1. The Motz problem==
  
The Motz problem is a benchmark for the Laplace equation problem that is often used for testing various special numerical methods; it was proposed in the literature for the solution of elliptic boundary value problems with boundary singularities. The boundary value problem is stated as follows (Figure [[#img-1|1]])
+
The Motz problem is a benchmark for the Laplace equation problem that is often used for testing various special numerical methods; it was proposed in the literature for the solution of elliptic boundary value problems with boundary singularities. The boundary value problem is stated as follows ([[#img-1|Figure 1]])
  
 
{| class="formulaSCP" style="width: 100%; text-align: left;"  
 
{| class="formulaSCP" style="width: 100%; text-align: left;"  
Line 57: Line 59:
 
where <math display="inline">r</math> and <math display="inline">\theta </math> are the polar coordinates; its convergence ratio is 2, and therefore valid over the whole domain  (Figure [[#img-1|1]]).
 
where <math display="inline">r</math> and <math display="inline">\theta </math> are the polar coordinates; its convergence ratio is 2, and therefore valid over the whole domain  (Figure [[#img-1|1]]).
  
Several papers discuss the accurate numerical calculation of the coefficients <math display="inline">A_i</math>. Lu <span id='citeF-1'></span>[[#cite-1|[1]]] calculated them with 17 significant digits using a collocation Treftz method. Li ''et al.'' <span id='citeF-2'></span>[[#cite-2|[2]]] used a special  boundary approximation method to derive the leading coefficients of the asymptotic solution expansion. Similar approaches were followed by Li ''et al.''  <span id='citeF-3'></span>[[#cite-3|[3]]] and Arad ''et al.'' <span id='citeF-4'></span>[[#cite-4|[4]]], who  employed boundary methods combined with least-squares and penalty method techniques to handle the boundary information, and by Georgiou  <span id='citeF-5'></span>[[#cite-5|[5]]], who employed Lagrange multipliers.
+
Several papers discuss the accurate numerical calculation of the coefficients <math display="inline">A_i</math>. Lu <span id='citeF-1'></span>[[#cite-1|[1]]] calculated them with 17 significant digits using a collocation Treftz method. Li ''et al.'' <span id='citeF-2'></span>[[#cite-2|[2]]] used a special  boundary approximation method to derive the leading coefficients of the asymptotic solution expansion. Similar approaches were followed by Li <span id='citeF-3'></span>[[#cite-3|[3]]] and Arad ''et al.'' <span id='citeF-4'></span>[[#cite-4|[4]]], who  employed boundary methods combined with least-squares and penalty method techniques to handle the boundary information, and by Georgiou  <span id='citeF-5'></span>[[#cite-5|[5]]], who employed Lagrange multipliers.
  
Some other authors focus on solving the Motz problem using different numerical techniques. Bernal and Kindelan <span id='citeF-6'></span>[[#cite-6|[6]]] thoroughly discussed the performance of the radial basis function method in the solution of the Motz problem, including the global radial basis function method due to Kansa <span id='citeF-7'></span>[[#cite-7|[7]]], as well as  Shu's local  differential quadrature method based on radial basis functions <span id='citeF-8'></span>[[#cite-8|[8]]]. They showed that the exponential convergence obtained by increasing resolution and increasing the shape parameter is lost due to the presence of the singularity. The accuracy of the solution can be increased by using collocation at some boundary nodes to restore the convergence of the method, although it is necessary to use functions similar to the terms in equation ([[#eq-1|(1)]]), which are capable of capturing the behaviour of the solution near the discontinuity.
+
Some other authors focus on solving the Motz problem using different numerical techniques. Bernal and Kindelan <span id='citeF-6'></span>[[#cite-6|[6]]] thoroughly discussed the performance of the radial basis function method in the solution of the Motz problem, including the global radial basis function method due to Kansa <span id='citeF-7'></span>[[#cite-7|[7]]], as well as  Shu's local  differential quadrature method based on radial basis functions <span id='citeF-8'></span>[[#cite-8|[8]]]. They showed that the exponential convergence obtained by increasing resolution and increasing the shape parameter is lost due to the presence of the singularity. The accuracy of the solution can be increased by using collocation at some boundary nodes to restore the convergence of the method, although it is necessary to use functions similar to the terms in equation ([[#eq-1|1]]), which are capable of capturing the behaviour of the solution near the discontinuity.
  
 
To produce an accurate solution based on the flexibility of meshfree methods, Hong <span id='citeF-9'></span>[[#cite-9|[9]]] addressed a singular basis function enrichment technique in a meshless method, using the leading terms of the local series expansion at the boundary singularity,  as enrichment functions for a local approximation space constructed to calculate the solution. The enrichment technique is derived from the Generalized Finite Element Method (GFEM) and produces accurate results, although it requires smooth global support basis functions and basis enrichment to  produce the numerical solution.
 
To produce an accurate solution based on the flexibility of meshfree methods, Hong <span id='citeF-9'></span>[[#cite-9|[9]]] addressed a singular basis function enrichment technique in a meshless method, using the leading terms of the local series expansion at the boundary singularity,  as enrichment functions for a local approximation space constructed to calculate the solution. The enrichment technique is derived from the Generalized Finite Element Method (GFEM) and produces accurate results, although it requires smooth global support basis functions and basis enrichment to  produce the numerical solution.
Line 86: Line 88:
 
|-
 
|-
 
| style="text-align: center;" | <math>L(u,A(x,y),B(x,y),C(x,y),D(x,y),E(x,y))=</math>
 
| style="text-align: center;" | <math>L(u,A(x,y),B(x,y),C(x,y),D(x,y),E(x,y))=</math>
| style="width: 5px;text-align: right;white-space: nowrap;" | (2)
 
 
|-
 
|-
 
| style="text-align: center;" | <math> A(x,y)\frac{\partial ^{2}u}{\partial x^{2}}+B(x,y)\frac{\partial ^{2}u}{\partial y^{2}}+C(x,y)\frac{\partial u}{\partial x}+D(x,y)\frac{\partial u}{\partial y}+E(x,y)u </math>
 
| style="text-align: center;" | <math> A(x,y)\frac{\partial ^{2}u}{\partial x^{2}}+B(x,y)\frac{\partial ^{2}u}{\partial y^{2}}+C(x,y)\frac{\partial u}{\partial x}+D(x,y)\frac{\partial u}{\partial y}+E(x,y)u </math>
 
|}
 
|}
 +
| style="width: 5px;text-align: right;white-space: nowrap;" |(2)
 
|}
 
|}
  
Line 161: Line 163:
 
{| style="text-align: left; margin:auto;width: 100%;"  
 
{| style="text-align: left; margin:auto;width: 100%;"  
 
|-
 
|-
| style="text-align: center;" | <math>\left( \begin{array}{rrrrr}1&  &1 &... & 1\\ 0&  &\Delta x_1 &... &\Delta x_q\\ 0&  &\Delta y_1 &... &\Delta y_q\\ 0&  &(\Delta x_1)^2 &... &(\Delta x_q)^2\\ 0&  &\Delta x_1\Delta y_1 &... &\Delta y_q\Delta z_q\\ 0&  &(\Delta y_1)^2 &... &(\Delta y_q)^2\\ \end{array} \right) \left( \begin{array}{lclcl}\Gamma _0\\ \Gamma _1\\ \Gamma _2\\ .\\ .\\ .\\ \Gamma _q\\ \end{array} \right) = \left( \begin{array}{rclcl}E(x_0,y_0)\\ C(x_0,y_0)\\ D(x_0,y_0)\\ 2A(x_0,y_0)\\ 0\\ 2B(x_0,y_0)\\ \end{array} \right). </math>
+
| style="text-align: center;" | <math>\left( \begin{array}{rrrrr}1&  &1 &... & 1\\ 0&  &\Delta x_1 &... &\Delta x_q\\ 0&  &\Delta y_1 &... &\Delta y_q\\ 0&  &(\Delta x_1)^2 &... &(\Delta x_q)^2\\ 0&  &\Delta x_1\Delta y_1 &... &\Delta y_q\Delta z_q\\ 0&  &(\Delta y_1)^2 &... &(\Delta y_q)^2\\ \end{array} \right) \left( \begin{array}{c}\Gamma _0\\ \Gamma _1\\ \Gamma _2\\ \cdot\\ \cdot\\ \cdot\\ \Gamma _q\\ \end{array} \right) = \left( \begin{array}{c}E(x_0,y_0)\\ C(x_0,y_0)\\ D(x_0,y_0)\\ 2A(x_0,y_0)\\ 0\\ 2B(x_0,y_0)\\ \end{array} \right). </math>
 
|}
 
|}
 
| style="width: 5px;text-align: right;white-space: nowrap;" | (5)
 
| style="width: 5px;text-align: right;white-space: nowrap;" | (5)
 
|}
 
|}
  
One must note that, in general, this system is not well-determined.  This poses the calculation of the coefficients in terms of a well studied algebraic problem which can be solved robustly.  Hence, the main question is how to produce a solution of equation ([[#eq-5|(5)]]) which provides a consistent second-order scheme. Several alternatives are possible, for example, in <span id='citeF-10'></span>[[#cite-10|[10]]], an efficient heuristic calculation to solve numerically the diffusion equation was proposed. It was based on an unconstrained optimization problem for every inner grid point. To apply it, first, we separate the first equation of the matrix system ([[#eq-5|(5)]])
+
One must note that, in general, this system is not well-determined.  This poses the calculation of the coefficients in terms of a well studied algebraic problem which can be solved robustly.  Hence, the main question is how to produce a solution of equation ([[#eq-5|5]]) which provides a consistent second-order scheme. Several alternatives are possible, for example, in <span id='citeF-10'></span>[[#cite-10|[10]]], an efficient heuristic calculation to solve numerically the diffusion equation was proposed. It was based on an unconstrained optimization problem for every inner grid point. To apply it, first, we separate the first equation of the matrix system ([[#eq-5|5]])
  
 
<span id="eq-6"></span>
 
<span id="eq-6"></span>
Line 209: Line 211:
 
{| style="text-align: left; margin:auto;width: 100%;"  
 
{| style="text-align: left; margin:auto;width: 100%;"  
 
|-
 
|-
| style="text-align: center;" | <math>\Gamma =\left( \begin{array}{lclcl}\Gamma _1\\ \Gamma _2\\ .\\ .\\ .\\ \Gamma _q\\ \end{array} \right),\quad  \beta = \left( \begin{array}{rclcl}C(x_0,y_0)\\ D(x_0,y_0)\\ 2A(x_0,y_0)\\ 0\\ 2B(x_0,y_0)\\ \end{array} \right). </math>
+
| style="text-align: center;" | <math>\Gamma =\left( \begin{array}{c}\Gamma _1\\ \Gamma _2\\ .\\ .\\ .\\ \Gamma _q\\ \end{array} \right),\quad  \beta = \left( \begin{array}{c}C(x_0,y_0)\\ D(x_0,y_0)\\ 2A(x_0,y_0)\\ 0\\ 2B(x_0,y_0)\\ \end{array} \right). </math>
 
|}
 
|}
 
| style="width: 5px;text-align: right;white-space: nowrap;" | (9)
 
| style="width: 5px;text-align: right;white-space: nowrap;" | (9)
Line 216: Line 218:
 
Next, <math display="inline">\Gamma _0</math> is obtained from equation ([[#eq-6|6]]).
 
Next, <math display="inline">\Gamma _0</math> is obtained from equation ([[#eq-6|6]]).
  
Equation  ([[#eq-7|(7)]]) defines the generalized finite differences in general grids. To solve it is straightforward when the grid is rectangular or nearly rectangular, and a sufficient number of nodes or symmetries are available.  
+
Equation  ([[#eq-7|7]]) defines the generalized finite differences in general grids. To solve it is straightforward when the grid is rectangular or nearly rectangular, and a sufficient number of nodes or symmetries are available.  
  
 
In irregular grids or clouds of points,  which lack symmetry, the system defined by <math display="inline">M</math> is underdetermined if <math display="inline">q>6</math>, which means that, in that case, at least seven neighbor nodes around <math display="inline">p_0</math> are required to get local consistency.  
 
In irregular grids or clouds of points,  which lack symmetry, the system defined by <math display="inline">M</math> is underdetermined if <math display="inline">q>6</math>, which means that, in that case, at least seven neighbor nodes around <math display="inline">p_0</math> are required to get local consistency.  
  
On the other hand, in a rather regular structured -or mapping- grid, symmetries reduce the number of required neighbors. An important consequence is that the symmetric location of nodes around the central node avoids rank deficiency of <math display="inline">M</math> in the system ([[#eq-7|(7)]]).
+
On the other hand, in a rather regular structured -or mapping- grid, symmetries reduce the number of required neighbors. An important consequence is that the symmetric location of nodes around the central node avoids rank deficiency of <math display="inline">M</math> in the system ([[#eq-7|7]]).
  
 
==3. Grids and stencils==
 
==3. Grids and stencils==
Line 226: Line 228:
 
The main issue to consider when solving the Motz problem is the singularity in the partial derivative <math display="inline">\frac{\partial u}{\partial x}</math> at the origin of coordinates. The usual approach to produce an accurate approximation at that point is to refine the grid around it. However, this causes not only an increase in the computational cost, but often the results are not satisfactory due to the high gradient values of the solution on the positive <math display="inline">x</math>-axis; even with robust techniques as those based on radial basis functions <span id='citeF-6'></span>[[#cite-6|[6]]], a special treatment is required near the origin. Some other approaches, like that of Hong, require a rather sophisticated process of patch mapping to produce accurate results <span id='citeF-9'></span>[[#cite-9|[9]]].
 
The main issue to consider when solving the Motz problem is the singularity in the partial derivative <math display="inline">\frac{\partial u}{\partial x}</math> at the origin of coordinates. The usual approach to produce an accurate approximation at that point is to refine the grid around it. However, this causes not only an increase in the computational cost, but often the results are not satisfactory due to the high gradient values of the solution on the positive <math display="inline">x</math>-axis; even with robust techniques as those based on radial basis functions <span id='citeF-6'></span>[[#cite-6|[6]]], a special treatment is required near the origin. Some other approaches, like that of Hong, require a rather sophisticated process of patch mapping to produce accurate results <span id='citeF-9'></span>[[#cite-9|[9]]].
  
When finite differences are considered, it is important to emphasize that the method has been successfully applied to problems with strong boundary layers, for instance, by Ivanenko  <span id='citeF-24'></span>[[#cite-24|[24]]] and Pereyra <span id='citeF-25'></span>[[#cite-25|[25]]], when the equidistribution principle is addressed. In other words, the key is to produce a suitably adapted grid through a simple procedure. In this paper, to derive the rule that controls the refinement around the origin, the clue is given by the shape of the significative terms in the closed-form solution series of equation ([[#eq-1|(1)]]).  
+
When finite differences are considered, it is important to emphasize that the method has been successfully applied to problems with strong boundary layers, for instance, by Ivanenko  <span id='citeF-24'></span>[[#cite-24|[24]]] and Pereyra <span id='citeF-25'></span>[[#cite-25|[25]]], when the equidistribution principle is addressed. In other words, the key is to produce a suitably adapted grid through a simple procedure. In this paper, to derive the rule that controls the refinement around the origin, the clue is given by the shape of the significative terms in the closed-form solution series of equation ([[#eq-1|1]]).  
  
 
Bearing in mind the equidistribution principle as discussed by Ivanenko,  since the graph of the solution <math display="inline">u(x)</math>  is  monotonic along the positive <math display="inline">x</math>-axis, and resembles that of  a boundary layer, there is a mapping <math display="inline">\xi \mapsto x=\xi ^\alpha </math> such that the expression
 
Bearing in mind the equidistribution principle as discussed by Ivanenko,  since the graph of the solution <math display="inline">u(x)</math>  is  monotonic along the positive <math display="inline">x</math>-axis, and resembles that of  a boundary layer, there is a mapping <math display="inline">\xi \mapsto x=\xi ^\alpha </math> such that the expression
Line 265: Line 267:
 
|}
 
|}
  
According to the estimations presented by Lu <span id='citeF-1'></span>[[#cite-1|[1]]], the first five coefficients in the equation ([[#eq-1|(1)]]) are the  significative ones: Those terms correspond to <math display="inline">i=0,1,...,5</math>. Since the value <math display="inline">i>1</math> in equation ([[#eq-11|(11)]]) represents an increase in the distance between the grid nodes around the origin, only the value <math display="inline">i=0</math>, which corresponds to <math display="inline">\alpha=2</math>, is considered.
+
According to the estimations presented by Lu <span id='citeF-1'></span>[[#cite-1|[1]]], the first five coefficients in the equation ([[#eq-1|1]]) are the  significative ones: Those terms correspond to <math display="inline">i=0,1,...,5</math>. Since the value <math display="inline">i>1</math> in equation ([[#eq-11|11]]) represents an increase in the distance between the grid nodes around the origin, only the value <math display="inline">i=0</math>, which corresponds to <math display="inline">\alpha=2</math>, is considered.
  
 
Thus, the mesh adaptation procedure is as follows:
 
Thus, the mesh adaptation procedure is as follows:
Line 271: Line 273:
 
<li>The grid for the whole domain of the problem is generated in the first step as a regular rectangular grid. </li>
 
<li>The grid for the whole domain of the problem is generated in the first step as a regular rectangular grid. </li>
 
<li>Every cell is subdivided into two triangles in the second step to define six neighbours in the stencil for every inner grid node. </li>
 
<li>Every cell is subdivided into two triangles in the second step to define six neighbours in the stencil for every inner grid node. </li>
<li>Finally, the grid is adapted by modifying the grid coordinates following the mapping given by equation([[#eq-11|(11)]]). </li>
+
<li>Finally, the grid is adapted by modifying the grid coordinates following the mapping given by equation([[#eq-11|11]]). </li>
 
+
 
</ol>
 
</ol>
  
In the numerical tests, the interval <math display="inline">[1.5,2.5]</math> was considered for the exponent <math display="inline">\alpha </math>. For comparison purposes, the mapping was applied to the <math display="inline">x</math>-axis (Type I grids), to the <math display="inline">y</math>-axis (Type II grids), and both axes (Type III). An example of the grids corresponding to the different distributions is shown in Figure [[#img-2|2]].
+
 
 +
In the numerical tests, the interval <math display="inline">[1.5,2.5]</math> was considered for the exponent <math display="inline">\alpha </math>. For comparison purposes, the mapping was applied to the <math display="inline">x</math>-axis (Type I grids), to the <math display="inline">y</math>-axis (Type II grids), and both axes (Type III). An example of the grids corresponding to the different distributions is shown in [[#img-2| Figure 2]].
  
 
<div id='img-2'></div>
 
<div id='img-2'></div>
Line 288: Line 290:
  
  
The three kinds of grids define a regular location of neighbors at the inner nodes: star-shaped subgrids  with 6 neighbours that can be used to approximate the derivatives in the left-hand side of equation ([[#eq-6|6]]); an example of such subgrid can be seen in Figure [[#img-3|3]]. The Laplacian operator is approximated at the point <math display="inline">p_0</math>, and the corresponding neighbours are <math display="inline">p_1</math>, <math display="inline">p_2</math>, <math display="inline">p_3</math>, <math display="inline">p_4</math>, <math display="inline">p_5</math>, and <math display="inline">p_6</math>. The values <math display="inline">h_1</math>, <math display="inline">h_2</math>, <math display="inline">k_1</math>, and <math display="inline">k_2</math> are the step size lengths in the stencil. Bearing in mind equation ([[#eq-5|(5)]]), this poses a problem with 6 equations and <math display="inline">k=7</math> unknowns; since the values <math display="inline">h_1</math>, <math display="inline">h_2</math>, <math display="inline">k_1</math>, and <math display="inline">k_2</math> are non zero, after preconditioning of equation ([[#eq-6|(6)]]) by scaling, the rank of <math display="inline">M</math> is equal to five,  which allows a robust calculation of the coefficients and leaves one degree of freedom available.  
+
The three kinds of grids define a regular location of neighbors at the inner nodes: star-shaped subgrids  with 6 neighbours that can be used to approximate the derivatives in the left-hand side of equation ([[#eq-6|6]]); an example of such subgrid can be seen in [[#img-3| Figure 3]]. The Laplacian operator is approximated at the point <math display="inline">p_0</math>, and the corresponding neighbours are <math display="inline">p_1</math>, <math display="inline">p_2</math>, <math display="inline">p_3</math>, <math display="inline">p_4</math>, <math display="inline">p_5</math>, and <math display="inline">p_6</math>. The values <math display="inline">h_1</math>, <math display="inline">h_2</math>, <math display="inline">k_1</math>, and <math display="inline">k_2</math> are the step size lengths in the stencil. Bearing in mind equation ([[#eq-5|5]]), this poses a problem with 6 equations and <math display="inline">k=7</math> unknowns; since the values <math display="inline">h_1</math>, <math display="inline">h_2</math>, <math display="inline">k_1</math>, and <math display="inline">k_2</math> are non zero, after preconditioning of equation ([[#eq-6|6]]) by scaling, the rank of <math display="inline">M</math> is equal to five,  which allows a robust calculation of the coefficients and leaves one degree of freedom available.  
  
 
<div id='img-3'></div>
 
<div id='img-3'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
+
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
 
|-
 
|-
|[[Image:Draft_Dominguez Mota_636399092-inner.png|498px|Selection of neighbours for an inner grid node]]
+
|style="padding:10px;"| [[Image:Draft_Dominguez Mota_636399092-inner.png|498px|Selection of neighbours for an inner grid node]]
 
|- style="text-align: center; font-size: 75%;"
 
|- style="text-align: center; font-size: 75%;"
| colspan="1" | '''Figure 3:''' Selection of neighbours for an inner grid node
+
| colspan="1" style="padding:10px;"| '''Figure 3'''. Selection of neighbours for an inner grid node
 
|}
 
|}
 +
 +
 
Regarding the boundary condition, a straightforward discretization using  the GFDM method was addressed in <span id='citeF-26'></span>[[#cite-26|[26]]]: at a boundary node of a structured grid, equation ([[#eq-7|7]]) is evaluated at its five neighbors, taking the right-hand side <math display="inline">\beta </math> to define the directional derivative. However, since the neighbors are often not balanced around the boundary node, such selection might produce first approximations which could cause loss of order in the Motz problem with the grids considered here.
 
Regarding the boundary condition, a straightforward discretization using  the GFDM method was addressed in <span id='citeF-26'></span>[[#cite-26|[26]]]: at a boundary node of a structured grid, equation ([[#eq-7|7]]) is evaluated at its five neighbors, taking the right-hand side <math display="inline">\beta </math> to define the directional derivative. However, since the neighbors are often not balanced around the boundary node, such selection might produce first approximations which could cause loss of order in the Motz problem with the grids considered here.
  
Line 310: Line 314:
 
|}
 
|}
  
where <math display="inline">n=(l_x,l_y)</math> is the outer normal vector at <math display="inline">p_0</math>, is approximated using the scheme of equation ([[#eq-5|5]]) with <math display="inline">A=B=E=0</math> and <math display="inline">C=l_x</math>, <math display="inline">D=l_y</math> using 7 points: a ghost point <math display="inline">p_g</math>, and the nodes <math display="inline">p_0</math>, <math display="inline">p_1</math>, <math display="inline">p_2</math>, <math display="inline">p_3</math>, <math display="inline">p_4</math> and <math display="inline">p_5</math> (see fig. [[#img-4|4]]), which yields the equation ([[#eq-12|12]])
+
where <math display="inline">n=(l_x,l_y)</math> is the outer normal vector at <math display="inline">p_0</math>, is approximated using the scheme of equation ([[#eq-5|5]]) with <math display="inline">A=B=E=0</math> and <math display="inline">C=l_x</math>, <math display="inline">D=l_y</math> using 7 points: a ghost point <math display="inline">p_g</math>, and the nodes <math display="inline">p_0</math>, <math display="inline">p_1</math>, <math display="inline">p_2</math>, <math display="inline">p_3</math>, <math display="inline">p_4</math> and <math display="inline">p_5</math> ([[#img-4|Figure 4]]), which yields the equation ([[#eq-12|12]])
  
 
<span id="eq-12"></span>
 
<span id="eq-12"></span>
Line 336: Line 340:
 
|}
 
|}
  
For convenience, the ghost point <math display="inline">p_g</math> is calculated in such a way that <math display="inline">p_0</math> is the midpoint of <math display="inline">p_g</math> and <math display="inline">p_c</math>, the latter being the centroid of the nodes <math display="inline">p_0</math>, <math display="inline">p_1</math>, <math display="inline">p_2</math>, <math display="inline">p_3</math>, <math display="inline">p_4</math>, and <math display="inline">p_5</math> (see figure [[#img-4|4]]).
+
For convenience, the ghost point <math display="inline">p_g</math> is calculated in such a way that <math display="inline">p_0</math> is the midpoint of <math display="inline">p_g</math> and <math display="inline">p_c</math>, the latter being the centroid of the nodes <math display="inline">p_0</math>, <math display="inline">p_1</math>, <math display="inline">p_2</math>, <math display="inline">p_3</math>, <math display="inline">p_4</math>, and <math display="inline">p_5</math> ([[#img-4|Figure 4]]).
  
 
The value of <math display="inline">u(p_g)</math> given by the equation ([[#eq-13|13]]) is substituted in the approximation to the Laplacian (equation ([[#eq-7|7]]) at <math display="inline">p_0</math> using the same 7 points <math display="inline">p_g</math>, <math display="inline">p_0</math>, <math display="inline">p_1</math>, <math display="inline">p_2</math>, <math display="inline">p_3</math>, <math display="inline">p_4</math>, and <math display="inline">p_5</math>, which eliminates the value of <math display="inline">u(p_g)</math>.
 
The value of <math display="inline">u(p_g)</math> given by the equation ([[#eq-13|13]]) is substituted in the approximation to the Laplacian (equation ([[#eq-7|7]]) at <math display="inline">p_0</math> using the same 7 points <math display="inline">p_g</math>, <math display="inline">p_0</math>, <math display="inline">p_1</math>, <math display="inline">p_2</math>, <math display="inline">p_3</math>, <math display="inline">p_4</math>, and <math display="inline">p_5</math>, which eliminates the value of <math display="inline">u(p_g)</math>.
  
Te assembled system for the whole grid becomes sparse as the number of grid points per side increases; it was solved using a sparse lu decomposition to produces the approximation to the solution of the problem.  <div id='img-4'></div>
+
Te assembled system for the whole grid becomes sparse as the number of grid points per side increases; it was solved using a sparse lu decomposition to produces the approximation to the solution of the problem.   
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
+
 
 +
<div id='img-4'></div>
 +
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 45%;"
 
|-
 
|-
|[[Image:Draft_Dominguez Mota_636399092-ghost.png|498px|Selection of ghost point at a boundary node of a general structured grid with Neumann condition]]
+
|style="padding:10px;"| [[File:Chavez_et_al_2020a_2003_498px-Draft_Dominguez_Mota_636399092-ghost.png]]
 +
<!--style="padding:10px;"| [[Image:Draft_Dominguez Mota_636399092-ghost.png|498px|Selection of ghost point at a boundary node of a general structured grid with Neumann condition]]-->
 
|- style="text-align: center; font-size: 75%;"
 
|- style="text-align: center; font-size: 75%;"
| colspan="1" | '''Figure 4:''' Selection of ghost point at a boundary node of a general structured grid with Neumann condition
+
| colspan="1" style="padding:10px;"| '''Figure 4'''. Selection of ghost point at a boundary node of a general structured grid with Neumann condition
 
|}
 
|}
  
 
==4 Numerical tests==
 
==4 Numerical tests==
  
Several tests were performed  with <math display="inline">m</math> and <math display="inline">2m-1</math> nodes per side on the <math display="inline">y</math> and <math display="inline">x</math> directions, respectively, varying <math display="inline">m</math> from 31 to 351. Grids of the three types in the previous section were generated, and the Motz problem was solved using the proposed scheme for every grid size and distribution. An example of the finite difference discretization described in the previous section is sketched in figure [[#img-5|5]].
+
Several tests were performed  with <math display="inline">m</math> and <math display="inline">2m-1</math> nodes per side on the <math display="inline">y</math> and <math display="inline">x</math> directions, respectively, varying <math display="inline">m</math> from 31 to 351. Grids of the three types in the previous section were generated, and the Motz problem was solved using the proposed scheme for every grid size and distribution. An example of the finite difference discretization described in the previous section is sketched in [[#img-5|Figure 5]].
  
 
<div id='img-5'></div>
 
<div id='img-5'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
+
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
 
|-
 
|-
|[[Image:Draft_Dominguez Mota_636399092-malla2.png|396px|Grid and calculated solution]]
+
|style="padding:10px;"| [[Image:Draft_Dominguez Mota_636399092-malla2.png|396px|Grid and calculated solution]]
 
|- style="text-align: center; font-size: 75%;"
 
|- style="text-align: center; font-size: 75%;"
| colspan="1" | '''Figure 5:''' Grid and calculated solution
+
| colspan="1" style="padding:10px;"| '''Figure 5'''. Grid and calculated solution
 
|}
 
|}
 +
 +
 
To assess the accuracy of the numerical solution, since  the approximation was calculated using structured grids, the quadratic error
 
To assess the accuracy of the numerical solution, since  the approximation was calculated using structured grids, the quadratic error
  
Line 371: Line 380:
 
|}
 
|}
  
was used, where <math display="inline">u_{i,j}</math> is the solution calculated using the difference scheme at the point <math display="inline">p_{i,j}</math>, <math display="inline">u</math> is the exact solution (equation [[#eq-1|1]]) at the same point, and <math display="inline">A_{i,j}</math> is the area of the rhombus defined by the nodes <math display="inline">p_{i+1,j}, p_{i-1,j}, p_{i,j+1}, p_{i,j-1}</math>.
+
was used, where <math display="inline">u_{i,j}</math> is the solution calculated using the difference scheme at the point <math display="inline">p_{i,j}</math>, <math display="inline">u</math> is the exact solution (equation [[#eq-1|(1)]]) at the same point, and <math display="inline">A_{i,j}</math> is the area of the rhombus defined by the nodes <math display="inline">p_{i+1,j}, p_{i-1,j}, p_{i,j+1}, p_{i,j-1}</math>.
 +
 
 +
The measure of the roughness of the grid is given by the meh norm <math display="inline">h</math>; the log graph of the  mean square errors versus the grid norm for the type I and II grids are shown in Figures [[#img-6|6]] and [[#img-7|7]].  One can note that there is a clear convergence tendency on all the tests despite the boundary singularity. But, when only a nonuniform distribution either along the <math display="inline">x</math> or the <math display="inline">y</math>-axis is used,  the second-order approximation is lost, in a similar way as reported by Bernal using radial basis functions. The best results are obtained with the grids which are non-uniform on <math display="inline">x</math> and <math display="inline">y</math>, as can be seen in  [[#img-7|Figure 7]].  
  
 
<div id='img-6'></div>
 
<div id='img-6'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
+
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
 
|-
 
|-
|[[Image:Draft_Dominguez Mota_636399092-figErrorX20.png|300px|]]
+
|style="padding:10px;"| [[Image:Draft_Dominguez Mota_636399092-figErrorX20.png|300px|]]
|[[Image:Draft_Dominguez Mota_636399092-figErrorY20.png|300px|Log values of the quadratic errors, nonuniform node distribution on x or y]]
+
|style="padding:10px;"| [[Image:Draft_Dominguez Mota_636399092-figErrorY20.png|300px|Log values of the quadratic errors, nonuniform node distribution on x or y]]
 
|- style="text-align: center; font-size: 75%;"
 
|- style="text-align: center; font-size: 75%;"
| colspan="2" | '''Figure 6:''' Log values of the quadratic errors, nonuniform node distribution on <math>x</math> or <math>y</math>
+
| colspan="2" style="padding:10px;"| '''Figure 6'''. Log values of the quadratic errors, nonuniform node distribution on <math>x</math> or <math>y</math>
 
|}
 
|}
The measure of the roughness of the grid is given by the meh norm <math display="inline">h</math>; the log graph of the  mean square errors versus the grid norm for the type I and II grids are shown in figures [[#img-6|6]] and [[#img-7|7]].  One can note that there is a clear convergence tendency on all the tests despite the boundary singularity. But, when only a nonuniform distribution either along the <math display="inline">x</math> or the <math display="inline">y</math>-axis is used,  the second-order approximation is lost, in a similar way as reported by Bernal using radial basis functions. The best results are obtained with the grids which are non-uniform on <math display="inline">x</math> and <math display="inline">y</math>, as can be seen in figure [[#img-7|7]].
 
  
The values of the errors, grid norms, and estimated convergence order are shown in tables  [[#table-1|1]] to [[#table-3|3]]. The values in the last table show that, indeed, quadratic convergence is obtained when the same nonuniform node distribution along  both axes is used, even in the presence of the boundary singularity.
 
  
 +
<div id='img-7'></div>
 +
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
 +
|-
 +
|style="padding:10px;"|[[Image:Draft_Dominguez Mota_636399092-figErrorXY20.png|500px|Log values of the quadratic errors, nonuniform node distribution on x and y]]
 +
|- style="text-align: center; font-size: 75%;"
 +
| colspan="1" style="padding:10px;"| '''Figure 7'''. Log values of the quadratic errors, nonuniform node distribution on <math>x</math> and <math>y</math>
 +
|}
 +
 +
 +
The values of the errors, grid norms, and estimated convergence order are shown in Tables  [[#table-1|1]] to [[#table-3|3]]. The values in the last table show that, indeed, quadratic convergence is obtained when the same nonuniform node distribution along  both axes is used, even in the presence of the boundary singularity.
 +
 +
<div class="center" style="font-size: 75%;">'''Table 1'''. Quadratic errors and approximation order for type I grids</div>
  
{| class="floating_tableSCP wikitable" style="text-align: right; margin: 1em auto;min-width:50%;"
+
<div id='table-1'></div>
|+ style="font-size: 75%;" |<span id='table-1'></span>Table. 1 Quadratic errors and approximation order for type I grids
+
{| class="wikitable" style="margin: 1em auto 0.1em auto;border-collapse: collapse;font-size:85%;width:auto;"  
|- style="border-top: 2px solid;"
+
|-style="text-align:center"
| <math display="inline">\alpha </math>  
+
! <math display="inline">\alpha </math> !!  <math>m</math> !! <math>n</math> !!  <math>h</math> !! <math>\epsilon </math> !!  Order  
| <math>m</math>
+
|- style="text-align:center"
| <math>n</math>
+
| <math>h</math>
+
| <math>\epsilon </math>
+
| Order  
+
|- style="border-top: 2px solid;"
+
 
|  1.00  
 
|  1.00  
 
| 31  
 
| 31  
Line 402: Line 418:
 
| 6.4941E-03  
 
| 6.4941E-03  
 
| &#8211;&#8211;  
 
| &#8211;&#8211;  
|-
+
|-style="text-align:center"
 
| 1.00  
 
| 1.00  
 
| 71  
 
| 71  
Line 409: Line 425:
 
| 2.7556E-03  
 
| 2.7556E-03  
 
| 1.01  
 
| 1.01  
|-
+
|-style="text-align:center"
 
| 1.00  
 
| 1.00  
 
| 111  
 
| 111  
Line 416: Line 432:
 
| 1.7489E-03  
 
| 1.7489E-03  
 
| 1.01  
 
| 1.01  
|-
+
|-style="text-align:center"
 
| 1.00  
 
| 1.00  
 
| 151  
 
| 151  
Line 423: Line 439:
 
| 1.2810E-03  
 
| 1.2810E-03  
 
| 1.00  
 
| 1.00  
|-
+
|-style="text-align:center"
 
| 1.00  
 
| 1.00  
 
| 191  
 
| 191  
Line 430: Line 446:
 
| 1.0106E-03  
 
| 1.0106E-03  
 
| 1.00  
 
| 1.00  
|-
+
|-style="text-align:center"
 
| 1.00  
 
| 1.00  
 
| 231  
 
| 231  
Line 437: Line 453:
 
| 8.3446E-04  
 
| 8.3446E-04  
 
| 1.00  
 
| 1.00  
|-
+
|-style="text-align:center"
 
| 1.00  
 
| 1.00  
 
| 271  
 
| 271  
Line 444: Line 460:
 
| 7.1060E-04  
 
| 7.1060E-04  
 
| 1.00  
 
| 1.00  
|-
+
|-style="text-align:center"
 
| 1.00  
 
| 1.00  
 
| 311  
 
| 311  
Line 451: Line 467:
 
| 6.1877E-04  
 
| 6.1877E-04  
 
| 1.00  
 
| 1.00  
|-
+
|-style="text-align:center"
 
| 1.25  
 
| 1.25  
 
| 31  
 
| 31  
Line 458: Line 474:
 
| 6.4057E-03  
 
| 6.4057E-03  
 
| &#8211;&#8211;  
 
| &#8211;&#8211;  
|-
+
|-style="text-align:center"
 
| 1.25  
 
| 1.25  
 
| 71  
 
| 71  
Line 465: Line 481:
 
| 2.7333E-03  
 
| 2.7333E-03  
 
| 1.01  
 
| 1.01  
|-
+
|-style="text-align:center"
 
| 1.25  
 
| 1.25  
 
| 111  
 
| 111  
Line 472: Line 488:
 
| 1.7503E-03  
 
| 1.7503E-03  
 
| 0.99  
 
| 0.99  
|-
+
|-style="text-align:center"
 
| 1.25  
 
| 1.25  
 
| 151  
 
| 151  
Line 479: Line 495:
 
| 1.2920E-03  
 
| 1.2920E-03  
 
| 0.98  
 
| 0.98  
|-
+
|-style="text-align:center"
 
| 1.25  
 
| 1.25  
 
| 191  
 
| 191  
Line 486: Line 502:
 
| 1.0260E-03  
 
| 1.0260E-03  
 
| 0.98  
 
| 0.98  
|-
+
|-style="text-align:center"
 
| 1.25  
 
| 1.25  
 
| 231  
 
| 231  
Line 493: Line 509:
 
| 8.5207E-04  
 
| 8.5207E-04  
 
| 0.97  
 
| 0.97  
|-
+
|-style="text-align:center"
 
| 1.25  
 
| 1.25  
 
| 271  
 
| 271  
Line 500: Line 516:
 
| 7.2925E-04  
 
| 7.2925E-04  
 
| 0.97  
 
| 0.97  
|- style="border-bottom: 2px solid;"
+
|- style="text-align:center"
 
| 1.25  
 
| 1.25  
 
| 311  
 
| 311  
Line 507: Line 523:
 
| 6.3784E-04  
 
| 6.3784E-04  
 
| 0.97  
 
| 0.97  
 
 
|}
 
|}
  
  
{|  class="floating_tableSCP wikitable" style="text-align: right; margin: 1em auto;min-width:50%;"
+
<div class="center" style="font-size: 75%;">'''Table 2'''. Quadratic errors and approximation order for type II grids</div>
|+ style="font-size: 75%;" |<span id='table-2'></span>Table. 2 Quadratic errors and approximation order for type II grids  
+
 
|- style="border-top: 2px solid;"
+
<div id='table-2'></div>
| <math display="inline">\alpha </math>  
+
{| class="wikitable" style="margin: 1em auto 0.1em auto;border-collapse: collapse;font-size:85%;width:auto;"  
| <math>m</math>
+
|-style="text-align:center"
| <math>n</math>
+
! <math display="inline">\alpha </math> !!  <math>m</math> !! <math>n</math> !!  <math>h</math> !! <math>\epsilon </math> !!  Order  
| <math>h</math>
+
|- style="text-align:center"
| <math>\epsilon </math>
+
| Order  
+
|- style="border-top: 2px solid;"
+
 
|  1.00  
 
|  1.00  
 
| 31  
 
| 31  
Line 527: Line 539:
 
| 0.0064941  
 
| 0.0064941  
 
| &#8211;&#8211;  
 
| &#8211;&#8211;  
|-
+
|-style="text-align:center"
 
| 1.00  
 
| 1.00  
 
| 71  
 
| 71  
Line 534: Line 546:
 
| 0.0027556  
 
| 0.0027556  
 
| 1.01  
 
| 1.01  
|-
+
|-style="text-align:center"
 
| 1.00  
 
| 1.00  
 
| 111  
 
| 111  
Line 541: Line 553:
 
| 0.0017489  
 
| 0.0017489  
 
| 1.01  
 
| 1.01  
|-
+
|-style="text-align:center"
 
| 1.00  
 
| 1.00  
 
| 151  
 
| 151  
Line 548: Line 560:
 
| 0.001281  
 
| 0.001281  
 
| 1.00  
 
| 1.00  
|-
+
|-style="text-align:center"
 
| 1.00  
 
| 1.00  
 
| 191  
 
| 191  
Line 555: Line 567:
 
| 0.0010106  
 
| 0.0010106  
 
| 1.00  
 
| 1.00  
|-
+
|-style="text-align:center"
 
| 1.00  
 
| 1.00  
 
| 231  
 
| 231  
Line 562: Line 574:
 
| 0.00083446  
 
| 0.00083446  
 
| 1.00  
 
| 1.00  
|-
+
|-style="text-align:center"
 
| 1.00  
 
| 1.00  
 
| 271  
 
| 271  
Line 569: Line 581:
 
| 0.0007106  
 
| 0.0007106  
 
| 1.00  
 
| 1.00  
|-
+
|-style="text-align:center"
 
| 1.00  
 
| 1.00  
 
| 311  
 
| 311  
Line 576: Line 588:
 
| 0.00061877  
 
| 0.00061877  
 
| 1.00  
 
| 1.00  
|-
+
|-style="text-align:center"
 
| 1.25  
 
| 1.25  
 
| 31  
 
| 31  
Line 583: Line 595:
 
| 0.0060306  
 
| 0.0060306  
 
| &#8211;&#8211;  
 
| &#8211;&#8211;  
|-
+
|-style="text-align:center"
 
| 1.25  
 
| 1.25  
 
| 71  
 
| 71  
Line 590: Line 602:
 
| 0.0026189  
 
| 0.0026189  
 
| 0.99  
 
| 0.99  
|-
+
|-style="text-align:center"
 
| 1.25  
 
| 1.25  
 
| 111  
 
| 111  
Line 597: Line 609:
 
| 0.0016827  
 
| 0.0016827  
 
| 0.98  
 
| 0.98  
|-
+
|-style="text-align:center"
 
| 1.25  
 
| 1.25  
 
| 151  
 
| 151  
Line 604: Line 616:
 
| 0.001243  
 
| 0.001243  
 
| 0.98  
 
| 0.98  
|-
+
|-style="text-align:center"
 
| 1.25  
 
| 1.25  
 
| 191  
 
| 191  
Line 611: Line 623:
 
| 0.00098714  
 
| 0.00098714  
 
| 0.98  
 
| 0.98  
|-
+
|-style="text-align:center"
 
| 1.25  
 
| 1.25  
 
| 231  
 
| 231  
Line 618: Line 630:
 
| 0.00081951  
 
| 0.00081951  
 
| 0.97  
 
| 0.97  
|-
+
|-style="text-align:center"
 
| 1.25  
 
| 1.25  
 
| 271  
 
| 271  
Line 625: Line 637:
 
| 0.00070109  
 
| 0.00070109  
 
| 0.97  
 
| 0.97  
|- style="border-bottom: 2px solid;"
+
|- style="text-align:center"
 
| 1.25  
 
| 1.25  
 
| 311  
 
| 311  
Line 632: Line 644:
 
| 0.00061294  
 
| 0.00061294  
 
| 0.97  
 
| 0.97  
 
 
|}
 
|}
  
  
{|  class="floating_tableSCP wikitable" style="text-align: right; margin: 1em auto;min-width:50%;"
+
<div class="center" style="font-size: 75%;">'''Table 3'''. Quadratic errors and approximation order for type III grids</div>
|+ style="font-size: 75%;" |<span id='table-3'></span>Table. 3 Quadratic errors and approximation order for type III grids  
+
 
|- style="border-top: 2px solid;"
+
<div id='table-2'></div>
| <math display="inline">\alpha </math>  
+
{| class="wikitable" style="margin: 1em auto 0.1em auto;border-collapse: collapse;font-size:85%;width:auto;"  
| <math>m</math>
+
|-style="text-align:center"
| <math>n</math>
+
! <math display="inline">\alpha </math> !!  <math>m</math> !! <math>n</math> !!  <math>h</math> !! <math>\epsilon </math> !!  Order  
| <math>h</math>
+
|- style="text-align:center"
| <math>\epsilon </math>
+
| Order  
+
|- style="border-top: 2px solid;"
+
 
|  2.00  
 
|  2.00  
 
| 71  
 
| 71  
Line 652: Line 660:
 
| 9.3949E-04  
 
| 9.3949E-04  
 
| 1.77  
 
| 1.77  
|-
+
|-style="text-align:center"
 
| 2.00  
 
| 2.00  
 
| 111  
 
| 111  
Line 659: Line 667:
 
| 4.1865E-04  
 
| 4.1865E-04  
 
| 1.80  
 
| 1.80  
|-
+
|-style="text-align:center"
 
| 2.00  
 
| 2.00  
 
| 151  
 
| 151  
Line 666: Line 674:
 
| 2.3926E-04  
 
| 2.3926E-04  
 
| 1.81  
 
| 1.81  
|-
+
|-style="text-align:center"
 
| 2.00  
 
| 2.00  
 
| 191  
 
| 191  
Line 673: Line 681:
 
| 1.5584E-04  
 
| 1.5584E-04  
 
| 1.82  
 
| 1.82  
|-
+
|-style="text-align:center"
 
| 2.00  
 
| 2.00  
 
| 231  
 
| 231  
Line 680: Line 688:
 
| 1.1006E-04  
 
| 1.1006E-04  
 
| 1.82  
 
| 1.82  
|-
+
|-style="text-align:center"
 
| 2.00  
 
| 2.00  
 
| 271  
 
| 271  
Line 687: Line 695:
 
| 8.2131E-05  
 
| 8.2131E-05  
 
| 1.83  
 
| 1.83  
|-
+
|-style="text-align:center"
 
| 2.00  
 
| 2.00  
 
| 311  
 
| 311  
Line 694: Line 702:
 
| 6.3786E-05  
 
| 6.3786E-05  
 
| 1.83  
 
| 1.83  
|-
+
|-style="text-align:center"
 
| 2.25  
 
| 2.25  
 
| 31  
 
| 31  
Line 701: Line 709:
 
| 4.5265E-03  
 
| 4.5265E-03  
 
| &#8211;&#8211;  
 
| &#8211;&#8211;  
|-
+
|-style="text-align:center"
 
| 2.25  
 
| 2.25  
 
| 71  
 
| 71  
Line 708: Line 716:
 
| 9.2020E-04  
 
| 9.2020E-04  
 
| 1.91  
 
| 1.91  
|-
+
|-style="text-align:center"
 
| 2.25  
 
| 2.25  
 
| 111  
 
| 111  
Line 715: Line 723:
 
| 3.8741E-04  
 
| 3.8741E-04  
 
| 1.93  
 
| 1.93  
|-
+
|-style="text-align:center"
 
| 2.25  
 
| 2.25  
 
| 151  
 
| 151  
Line 722: Line 730:
 
| 2.1315E-04  
 
| 2.1315E-04  
 
| 1.94  
 
| 1.94  
|-
+
|-style="text-align:center"
 
| 2.25  
 
| 2.25  
 
| 191  
 
| 191  
Line 729: Line 737:
 
| 1.3496E-04  
 
| 1.3496E-04  
 
| 1.94  
 
| 1.94  
|-
+
|-style="text-align:center"
 
| 2.25  
 
| 2.25  
 
| 231  
 
| 231  
Line 736: Line 744:
 
| 9.3190E-05  
 
| 9.3190E-05  
 
| 1.94  
 
| 1.94  
|-
+
|-style="text-align:center"
 
| 2.25  
 
| 2.25  
 
| 271  
 
| 271  
Line 743: Line 751:
 
| 6.8259E-05  
 
| 6.8259E-05  
 
| 1.95  
 
| 1.95  
|-
+
|-style="text-align:center"
 
| 2.25  
 
| 2.25  
 
| 311  
 
| 311  
Line 750: Line 758:
 
| 5.2179E-05  
 
| 5.2179E-05  
 
| 1.95  
 
| 1.95  
|-
+
|-style="text-align:center"
 
| 2.50  
 
| 2.50  
 
| 31  
 
| 31  
Line 757: Line 765:
 
| 4.9361E-03  
 
| 4.9361E-03  
 
| &#8211;&#8211;  
 
| &#8211;&#8211;  
|-
+
|-style="text-align:center"
 
| 2.50  
 
| 2.50  
 
| 71  
 
| 71  
Line 764: Line 772:
 
| 9.2007E-04  
 
| 9.2007E-04  
 
| 2.02  
 
| 2.02  
|-
+
|-style="text-align:center"
 
| 2.50  
 
| 2.50  
 
| 111  
 
| 111  
Line 771: Line 779:
 
| 3.7194E-04  
 
| 3.7194E-04  
 
| 2.02  
 
| 2.02  
|-
+
|-style="text-align:center"
 
| 2.50  
 
| 2.50  
 
| 151  
 
| 151  
Line 778: Line 786:
 
| 1.9951E-04  
 
| 1.9951E-04  
 
| 2.02  
 
| 2.02  
|-
+
|-style="text-align:center"
 
| 2.50  
 
| 2.50  
 
| 191  
 
| 191  
Line 785: Line 793:
 
| 1.2407E-04  
 
| 1.2407E-04  
 
| 2.02  
 
| 2.02  
|-
+
|-style="text-align:center"
 
| 2.50  
 
| 2.50  
 
| 231  
 
| 231  
Line 792: Line 800:
 
| 8.4506E-05  
 
| 8.4506E-05  
 
| 2.02  
 
| 2.02  
|-
+
|-style="text-align:center"
 
| 2.50  
 
| 2.50  
 
| 271  
 
| 271  
Line 799: Line 807:
 
| 6.1224E-05  
 
| 6.1224E-05  
 
| 2.02  
 
| 2.02  
|-
+
|-style="text-align:center"
 
| 2.50  
 
| 2.50  
 
| 311  
 
| 311  
Line 806: Line 814:
 
| 4.6380E-05  
 
| 4.6380E-05  
 
| 2.02  
 
| 2.02  
|-
+
|-style="text-align:center"
 
| 2.75  
 
| 2.75  
 
| 31  
 
| 31  
Line 813: Line 821:
 
| 5.3606E-03  
 
| 5.3606E-03  
 
| &#8211;&#8211;  
 
| &#8211;&#8211;  
|-
+
|-style="text-align:center"
 
| 2.75  
 
| 2.75  
 
| 71  
 
| 71  
Line 820: Line 828:
 
| 9.3491E-04  
 
| 9.3491E-04  
 
| 2.10  
 
| 2.10  
|-
+
|-style="text-align:center"
 
| 2.75  
 
| 2.75  
 
| 111  
 
| 111  
Line 827: Line 835:
 
| 3.6799E-04  
 
| 3.6799E-04  
 
| 2.08  
 
| 2.08  
|-
+
|-style="text-align:center"
 
| 2.75  
 
| 2.75  
 
| 151  
 
| 151  
Line 834: Line 842:
 
| 1.9446E-04  
 
| 1.9446E-04  
 
| 2.07  
 
| 2.07  
|-
+
|-style="text-align:center"
 
| 2.75  
 
| 2.75  
 
| 191  
 
| 191  
Line 841: Line 849:
 
| 1.1975E-04  
 
| 1.1975E-04  
 
| 2.06  
 
| 2.06  
|-
+
|-style="text-align:center"
 
| 2.75  
 
| 2.75  
 
| 231  
 
| 231  
Line 848: Line 856:
 
| 8.0998E-05  
 
| 8.0998E-05  
 
| 2.06  
 
| 2.06  
|-
+
|-style="text-align:center"
 
| 2.75  
 
| 2.75  
 
| 271  
 
| 271  
Line 855: Line 863:
 
| 5.8375E-05  
 
| 5.8375E-05  
 
| 2.05  
 
| 2.05  
|-
+
|-style="text-align:center"
 
| 2.75  
 
| 2.75  
 
| 311  
 
| 311  
Line 862: Line 870:
 
| 4.4042E-05  
 
| 4.4042E-05  
 
| 2.05  
 
| 2.05  
|-
+
|-style="text-align:center"
 
| 3.00  
 
| 3.00  
 
| 31  
 
| 31  
Line 869: Line 877:
 
| 5.7955E-03  
 
| 5.7955E-03  
 
| &#8211;&#8211;  
 
| &#8211;&#8211;  
|-
+
|-style="text-align:center"
 
| 3.00  
 
| 3.00  
 
| 71  
 
| 71  
Line 876: Line 884:
 
| 9.6029E-04  
 
| 9.6029E-04  
 
| 2.17  
 
| 2.17  
|-
+
|-style="text-align:center"
 
| 3.00  
 
| 3.00  
 
| 111  
 
| 111  
Line 883: Line 891:
 
| 3.7142E-04  
 
| 3.7142E-04  
 
| 2.13  
 
| 2.13  
|-
+
|-style="text-align:center"
 
| 3.00  
 
| 3.00  
 
| 151  
 
| 151  
Line 890: Line 898:
 
| 1.9454E-04  
 
| 1.9454E-04  
 
| 2.10  
 
| 2.10  
|-
+
|-style="text-align:center"
 
| 3.00  
 
| 3.00  
 
| 191  
 
| 191  
Line 897: Line 905:
 
| 1.1915E-04  
 
| 1.1915E-04  
 
| 2.09  
 
| 2.09  
|-
+
|-style="text-align:center"
 
| 3.00  
 
| 3.00  
 
| 231  
 
| 231  
Line 904: Line 912:
 
| 8.0299E-05  
 
| 8.0299E-05  
 
| 2.08  
 
| 2.08  
|-
+
|-style="text-align:center"
 
| 3.00  
 
| 3.00  
 
| 271  
 
| 271  
Line 911: Line 919:
 
| 5.7715E-05  
 
| 5.7715E-05  
 
| 2.07  
 
| 2.07  
|- style="border-bottom: 2px solid;"
+
|- style="text-align:center"
 
| 3.00  
 
| 3.00  
 
| 311  
 
| 311  
Line 918: Line 926:
 
| 4.3453E-05  
 
| 4.3453E-05  
 
| 2.06  
 
| 2.06  
 
 
|}
 
|}
  
==5 Results and discussion==
+
==5. Results and discussion==
  
 
Several assertions can be drawn from the numerical results. First, despite the presence of the boundary singularity, the results are accurate enough. This is certainly a strength of the proposed finite difference scheme, which is rather simple. Since the main differences between the scheme discussed in this paper and the standard five-point scheme for the Laplacian are the addition of two extra “diagonal” nodes in the inner grid nodes and the use of a suitable grid spacing, this suggests that the proposed stencils are suitable to produce accurate results.
 
Several assertions can be drawn from the numerical results. First, despite the presence of the boundary singularity, the results are accurate enough. This is certainly a strength of the proposed finite difference scheme, which is rather simple. Since the main differences between the scheme discussed in this paper and the standard five-point scheme for the Laplacian are the addition of two extra “diagonal” nodes in the inner grid nodes and the use of a suitable grid spacing, this suggests that the proposed stencils are suitable to produce accurate results.
  
<div id='img-7'></div>
 
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
 
|-
 
|[[Image:Draft_Dominguez Mota_636399092-figErrorXY20.png|600px|Log values of the quadratic errors, nonuniform node distribution on x and y]]
 
|- style="text-align: center; font-size: 75%;"
 
| colspan="1" | '''Figure 7:''' Log values of the quadratic errors, nonuniform node distribution on <math>x</math> and <math>y</math>
 
|}
 
 
The convergence orders in table [[#table-3|3]] must be highlighted. Having in mind the presence of the singularity and the simplicity of the addressed difference stencils, this is indeed a remarkable result. Since the  boundary node redistribution was defined by a very simple rule,  this suggests that it is unnecessary to consider another kind of grids or more complex refinements near the singularity to solve the Motz problem, although, it must be acknowledged that the mesh adapting strategy was proposed ''ad hoc'' for it. Naturally, the numerical treatment of different singularities in other problems require a different application of the equidistribution principle.
 
The convergence orders in table [[#table-3|3]] must be highlighted. Having in mind the presence of the singularity and the simplicity of the addressed difference stencils, this is indeed a remarkable result. Since the  boundary node redistribution was defined by a very simple rule,  this suggests that it is unnecessary to consider another kind of grids or more complex refinements near the singularity to solve the Motz problem, although, it must be acknowledged that the mesh adapting strategy was proposed ''ad hoc'' for it. Naturally, the numerical treatment of different singularities in other problems require a different application of the equidistribution principle.
  
==6 Conclusions==
+
==6. Conclusions==
  
 
Indeed, generalized finite differences can be used to solve the Motz problem numerically. The numerical results show that generalized finite differences  produce accurate solutions to a complicated benchmark problem, and the most noteworthy fact is that this can be done even with very simple structured grids. This is an important feature of the GFDM method, considering the difficulty of the problem. Thinking about the solution to several linear and nonlinear problems with this method which appear in the literature <span id='citeF-27'></span><span id='citeF-10'></span><span id='citeF-22'></span><span id='citeF-13'></span>[[#cite-27|[27,10,22,13]]], it is certainly a reliable numerical tool.  
 
Indeed, generalized finite differences can be used to solve the Motz problem numerically. The numerical results show that generalized finite differences  produce accurate solutions to a complicated benchmark problem, and the most noteworthy fact is that this can be done even with very simple structured grids. This is an important feature of the GFDM method, considering the difficulty of the problem. Thinking about the solution to several linear and nonlinear problems with this method which appear in the literature <span id='citeF-27'></span><span id='citeF-10'></span><span id='citeF-22'></span><span id='citeF-13'></span>[[#cite-27|[27,10,22,13]]], it is certainly a reliable numerical tool.  
Line 944: Line 944:
 
We want to thank  CIC-UMSNH  9.16,  and Aula-CIMNE Morelia for the financial support for this work.
 
We want to thank  CIC-UMSNH  9.16,  and Aula-CIMNE Morelia for the financial support for this work.
  
===BIBLIOGRAPHY===
+
==References==
  
 
<div id="cite-1"></div>
 
<div id="cite-1"></div>
'''[[#citeF-1|[1]]]''' Lu, T. T. and Hu, H. Y. and Li, Z. C. (2004) "Highly accurate solutions of Motz's and the cracked beam problems", Volume 28. Engineering Analysis with Boundary Elements 1387&#8211;1403
+
[[#citeF-1|[1]]] Lu T.T., Hu H.Y., Li Z.C. Highly accurate solutions of Motz's and the cracked beam problems. Engineering Analysis with Boundary Elements, 28:1387–1403, 2004.
  
 
<div id="cite-2"></div>
 
<div id="cite-2"></div>
'''[[#citeF-2|[2]]]''' Li, Z. C. and Mathon, R. and Sermer, P. (1987) "Boundary methods for solving elliptic problems with singularities and interfaces", Volume 24. SIAM Journal of Numerical Analaysis 487&#8211;498
+
[[#citeF-2|[2]]] Li Z.C., Mathon R., Sermer P. Boundary methods for solving elliptic problems with singularities and interfaces. SIAM Journal of Numerical Analysis, 24:487–498, 1987.
  
 
<div id="cite-3"></div>
 
<div id="cite-3"></div>
'''[[#citeF-3|[3]]]''' Li, Z. C. and Mathon, R. and Sermer, P. (1987) "Boundary methods for solving elliptic problems with singularities and interfaces", Volume 24. SIAM Journal on Numerical Analysis 3 487-498
+
[[#citeF-3|[3]]] Li Z.C. An approach for combining the Ritz-Galerkin and Finite-Element Methods. Journal of Approximation Theory, 39:132-152, 1983.
  
 
<div id="cite-4"></div>
 
<div id="cite-4"></div>
'''[[#citeF-4|[4]]]''' Yosibash, Z. and Ben-Dor, G. and Yakhot, A. (1998) "Computing flux intensity factors by a boundary method for elliptic equations with singularities", Volume 14. Communications in Numerical Methods in Engineering 657&#8211;670
+
[[#citeF-4|[4]]] Arad M., Yosibash Z., Ben-Dor G.Yakhot A. Computing flux intensity factors by a boundary method for elliptic equations with singularities. Communications in Numerical Methods in Engineering, 14:657–670, 1998.
  
 
<div id="cite-5"></div>
 
<div id="cite-5"></div>
'''[[#citeF-5|[5]]]''' Georgiou, G. C. and Olson, L. and Smyrlis, Y. S. (1996) "A singular function boundary integral method for the laplace equation", Volume 12. Communications in Numerical Methods in Engineering 2 127-134
+
[[#citeF-5|[5]]] Georgiou G.C., Olson L., Smyrlis Y.S. A singular function boundary integral method for the Laplace equation. Communications in Numerical Methods in Engineering, 12:127-134, 1996.
  
 
<div id="cite-6"></div>
 
<div id="cite-6"></div>
'''[[#citeF-6|[6]]]''' Bernal, Francisco and Kindelan, Manuel. (2010) "Radial basis function solution of the Motz problem", Volume 27. Engineering Computations 606&#8211;620
+
[[#citeF-6|[6]]] Bernal F., Kindelan, M. Radial basis function solution of the Motz problem. Engineering Computations, 27:606–620, 2010.
  
 
<div id="cite-7"></div>
 
<div id="cite-7"></div>
'''[[#citeF-7|[7]]]''' Kansa, E. J. (1990) "Multiquadrics, a scattered data approximation scheme with applications to computational fluid dynamics II: solutions to parabolic, hyperbolic and elliptic partial differential equations", Volume 19. Computers and Mathematics with Applications 147&#8211;161
+
[[#citeF-7|[7]]] Kansa E.J. Multiquadrics, a scattered data approximation scheme with applications to computational fluid dynamics II: solutions to parabolic, hyperbolic and elliptic partial differential equations. Computers and Mathematics with Applications, 19:147–161, 1990.
  
 
<div id="cite-8"></div>
 
<div id="cite-8"></div>
'''[[#citeF-8|[8]]]''' Shu, C. and Ding, H. and Yeo, K. S. (2003) "Local radial basis function-based differential quadrature method and its application to solve two-dimensional incompressible Navier-Stokes equations", Volume 192. Computer Methods in Applied Mechanics and Engineering 3 941&#8211;954
+
[[#citeF-8|[8]]] Shu C., Ding H., Yeo K.S. Local radial basis function-based differential quadrature method and its application to solve two-dimensional incompressible Navier-Stokes equations. Computer Methods in Applied Mechanics and Engineering, 192:941–954, 2003.
  
 
<div id="cite-9"></div>
 
<div id="cite-9"></div>
'''[[#citeF-9|[9]]]''' Hong, Won-Tak. (2016) "Enriched Meshfree Method for an Accurate Numerical Solution of the Motz Problem"
+
[[#citeF-9|[9]]] Hong W.T. Enriched meshfree method for an accurate numerical solution of the Motz problem. Advances in Mathematical Physics, 2016:1-12, 2016.
  
 
<div id="cite-10"></div>
 
<div id="cite-10"></div>
'''[[#citeF-10|[10]]]''' Domínguez-Mota, F. J. and Armenta, S. M. and Tinoco-Guerrero, G. and Tinoco-Ruiz, J. G. (2014) "Finite difference schemes satisfying an optimality condition for the unsteady heat equation", Volume 106. Mathematics and Computers in Simulation 76&#8211;83
+
[[#citeF-10|[10]]] Domínguez-Mota F.J., Mendoza-Armenta S., Tinoco-Guerrero G., Tinoco-Ruiz J.G. Finite difference schemes satisfying an optimality condition for the unsteady heat equation. Mathematics and Computers in Simulation, 106:76–83, 2014.
  
 
<div id="cite-11"></div>
 
<div id="cite-11"></div>
'''[[#citeF-11|[11]]]''' Chávez-Negrete, C and Domínguez-M, F and Santana-Quinteros, D. (2018) "Numerical solution of Richards' equation of water flow by generalized finite differences", Volume 101. Computers and Geotechnics 1 168&#8211;175
+
[[#citeF-11|[11]]] Chávez-Negrete C., Domínguez-Mota F.J., Santana-Quinteros D. Numerical solution of Richards' equation of water flow by generalized finite differences. Computers and Geotechnics, 101:168–175, 2018.
  
 
<div id="cite-12"></div>
 
<div id="cite-12"></div>
'''[[#citeF-12|[12]]]''' Tinoco-Guerrero, G. and Domínguez-Mota, F. J. and Gaona-Arias, A. and Ruiz-Zavala, M. L. and Tinoco-Ruiz, J. G. (2018) "A stability analysis for a generalized finite-difference scheme applied to the pure advection equation", Volume 147. Mathematics and Computers in Simulation 293&#8211;300
+
[[#citeF-12|[12]]] Tinoco-Guerrero G., Domínguez-Mota F.J., Gaona-Arias A., Ruiz-Zavala M.L., Tinoco-Ruiz J.G. A stability analysis for a generalized finite-difference scheme applied to the pure advection equation. Mathematics and Computers in Simulation, 147:293–300, 2018.
  
 
<div id="cite-13"></div>
 
<div id="cite-13"></div>
'''[[#citeF-13|[13]]]''' Muñoz, J. J. B. and Ureña-Prieto, F. and Corvinos, L. A. G. (2007) "Solving parabolic and hyperbolic equations by the generalized finite difference method", Volume 209. Journal of Computational and Applied Mathematics 2 208&#8211;233
+
[[#citeF-13|[13]]] Benito J.J., Ureña F., Gavete L. Solving parabolic and hyperbolic equations by the generalized finite difference method. Journal of Computational and Applied Mathematics, 209:208–233, 2007.
  
 
<div id="cite-14"></div>
 
<div id="cite-14"></div>
'''[[#citeF-14|[14]]]''' Ureña, F. and Benito, J. J. and Gavete, L. (2011) "Application of the generalized finite difference method to solve the advection-diffusion equation", Volume 235. Journal of Computational and Applied Mathematics 1849&#8211;1855
+
[[#citeF-14|[14]]] Benito J.J., Ureña F., Gavete, L. Application of the generalized finite difference method to solve the advection-diffusion equation. Journal of Computational and Applied Mathematics, 235:1849–1855, 2011.
  
 
<div id="cite-15"></div>
 
<div id="cite-15"></div>
'''[[#citeF-15|[15]]]''' Ureña, F. and Benito, J. J. and Salete, E. and Gavete, L. (2012) "A note on the application of the generalized finite difference method to seismic wave propagation in 2D", Volume 236. Journal of Computational and Applied Mathematics 3016&#8211;3025
+
[[#citeF-15|[15]]] Benito J.J., Ureña F., Salete E., Gavete, L. A note on the application of the generalized finite difference method to seismic wave propagation in 2D. Journal of Computational and Applied Mathematics, 236:3016–3025, 2012.
  
 
<div id="cite-16"></div>
 
<div id="cite-16"></div>
'''[[#citeF-16|[16]]]''' Muñoz, J. J. B. and Ureña-Prieto, F. and Corvinos, L. and Gavete, A. (2017) "Solving second order non-linear elliptic partial differential equations using generalized finite difference method", Volume 318. Journal of Computational and Applied Mathematics 378&#8211;387
+
[[#citeF-16|[16]]] Gavete L., Ureña F., Benito J.J., García A., Ureña M., Salete E. Solving second order non-linear elliptic partial differential equations using generalized finite difference method. Journal of Computational and Applied Mathematics, 318:378–387, 2017.
  
 
<div id="cite-17"></div>
 
<div id="cite-17"></div>
'''[[#citeF-17|[17]]]''' Chan, Hsin-Fang and Wei Li, Po. (2013) "Generalized finite difference method for solving two-dimensional Burgers' equations", Volume 79. Procedia Engineering 55&#8211;60
+
[[#citeF-17|[17]]] Fan C.M., Li P.W. Generalized finite difference method for solving two-dimensional Burgers' equations. Procedia Engineering, 79:55–60, 2014.
  
 
<div id="cite-18"></div>
 
<div id="cite-18"></div>
'''[[#citeF-18|[18]]]''' Chan, Hsin-Fang and Ming Fan, Chia and Kuo, Chia-Wen. (2013) "Generalized finite difference method for solving two-dimensional non-linear obstacle problem", Volume 37. Engineering Analysis with Boundary Elements 1189&#8211;1196
+
[[#citeF-18|[18]]] Chan H.F., Fan C.M., Kuo C.W. Generalized finite difference method for solving two-dimensional non-linear obstacle problem. Engineering Analysis with Boundary Elements, 37:1189–1196, 2013.
  
 
<div id="cite-19"></div>
 
<div id="cite-19"></div>
'''[[#citeF-19|[19]]]''' Wei Li, Po and Ming Fan, Chia. (2017) "Generalized finite difference method for two-dimensional shallow water equations", Volume 80. Engineering Analysis with Boundary Elements 58&#8211;71
+
[[#citeF-19|[19]]] Li P.W., Fan C.M. Generalized finite difference method for two-dimensional shallow water equations. Engineering Analysis with Boundary Elements, 80:58–71, 2017.
  
 
<div id="cite-20"></div>
 
<div id="cite-20"></div>
'''[[#citeF-20|[20]]]''' Chan, Hsin-Fang and Ming Fan, Chia and Kuo, Chia-Wen. (2018) "Generalized finite difference method for solving the double-diffusive natural convection in fluid-saturated porous media", Volume 95. Engineering Analysis with Boundary Elements 175&#8211;186
+
[[#citeF-20|[20]]] Li P.W. , Chen W., Fu Z.J., Fan C.M. Generalized finite difference method for solving the double-diffusive natural convection in fluid-saturated porous media. Engineering Analysis with Boundary Elements 95:175–186, 2018.
  
 
<div id="cite-21"></div>
 
<div id="cite-21"></div>
'''[[#citeF-21|[21]]]''' Liszka, T. and Orkisz, J. (1980) "The finite difference method at arbitrary irregular grids and its application in applied mechanics", Volume 11. Computers and Structures 1 83&#8211;95
+
[[#citeF-21|[21]]] Liszka T., Orkisz J. The finite difference method at arbitrary irregular grids and its application in applied mechanics. Computers and Structures, 11:83–95, 1980.
  
 
<div id="cite-22"></div>
 
<div id="cite-22"></div>
'''[[#citeF-22|[22]]]''' Muñoz, J. J. B. and Prieto, J. B. and Corvinos, L. A. G. and Alonso, B. (2008) "A posteriori error estimator and indicator in generalized finite differences. Application to improve the approximated solution of elliptic pde's.", Volume 85. International Journal of Computer Mathematics 3-4 359&#8211;370
+
[[#citeF-22|[22]]] Benito J.J., Ureña F., Gavete L., Alonso B. A posteriori error estimator and indicator in generalized finite differences. Application to improve the approximated solution of elliptic PDEs. International Journal of Computer Mathematics, 85:359–370, 2008.
  
 
<div id="cite-23"></div>
 
<div id="cite-23"></div>
'''[[#citeF-23|[23]]]''' Celia, M. and Gray, W. (1992) "Numerical Methods for Differential Equations". Prentice-Hall
+
[[#citeF-23|[23]]] Celia, M., Gray, W. Numerical Methods for Differential Equations: fundamental concepts for scientific and engineering applications. Prentice-Hall, 1992.
  
 
<div id="cite-24"></div>
 
<div id="cite-24"></div>
'''[[#citeF-24|[24]]]''' Ivanenko, S. (2003) "Selected chapters on grid generation and applications. Azarenok, B. (Ed.)". Russian academy of sciences
+
[[#citeF-24|[24]]] Ivanenko, S. In Selected chapters on grid generation and applications. B. Azarenok (Ed.), Russian Academy of Sciences, 2009.
  
 
<div id="cite-25"></div>
 
<div id="cite-25"></div>
'''[[#citeF-25|[25]]]''' Pereyra, V. and Sewell, G. (1975) "A singular function boundary integral method for the laplace equation", Volume 23. Numer. Math. 261-268
+
[[#citeF-25|[25]]] Pereyra V., Sewell G. A singular function boundary integral method for the Laplace equation. Numer. Math., 23:261-268, 1975.
  
 
<div id="cite-26"></div>
 
<div id="cite-26"></div>
'''[[#citeF-26|[26]]]''' Francisco-Javier Domínguez-Mota and Sanzon Mendoza-Armenta and José-Gerardo Tinoco-Ruiz and Gerardo Tinoco-Guerrero. (2011) "Numerical Solution of Poisson-like Equations with Robin Boundary Conditions Using a Finite Difference Scheme Dfined by an Optimality Condition", Volume 17. MASCOT11 Proceedings 101-110
+
[[#citeF-26|[26]]] Domínguez-Mota F.J., Mendoza S., Tinoco-Ruiz J.G, Tinoco-Guerrero G. Numerical solution of Poisson-like equations with Robin boundary conditions using a finite difference scheme defined by an optimality condition. Proceedings of MASCOT 2011, 17, 101-110, 2011.
  
 
<div id="cite-27"></div>
 
<div id="cite-27"></div>
'''[[#citeF-27|[27]]]''' Domínguez-M, F. J. and Equihua, M. and Mendoza, S. and Tinoco-R., J. G. (2010) "Solución de ecuaciones diferenciales elípticas en regiones planas irregulares usando mallas convexas generadas por métodos variacionales empleando elementos finitos", Volume 26. Rev. Int. Mét. Num. Cálc. Dis. Ing. 187&#8211;194
+
[[#citeF-27|[27]]] Domínguez-Mota F.J., Equihua M., Mendoza S., Tinoco-Ruiz J.G. Solución de ecuaciones diferenciales elípticas en regiones planas irregulares usando mallas convexas generadas por métodos variacionales empleando elementos finitos. Rev. Int. Mét. Num. Cálc. Dis. Ing., 26:187–194, 2010.

Latest revision as of 14:36, 20 May 2021

Abstract

In this paper, it is presented a formulation of a generalized finite difference scheme to solve the Motz problem. It is based on a general difference scheme defined by an optimality condition, which has been developed to solve Poisson-like equations whose domains are approximated by a wide variety of grids over general regions. Numerical examples showing second-order accuracy of the calculated solutions are presented.

Keywords: Generalized finite difference, Motz problem

1. The Motz problem

The Motz problem is a benchmark for the Laplace equation problem that is often used for testing various special numerical methods; it was proposed in the literature for the solution of elliptic boundary value problems with boundary singularities. The boundary value problem is stated as follows (Figure 1)

A singularity similar to that at the origin appears analogously in several hydrological applications, often related to an abrupt change in the boundary condition, particularly in groundwater infiltration in the alternation of pervious and impervious boundary zones.

Motz problem domain and boundary conditions.
Figure 1. Motz problem domain and boundary conditions


The closed-form solution of Motz problem is given as a series

(1)

where and are the polar coordinates; its convergence ratio is 2, and therefore valid over the whole domain (Figure 1).

Several papers discuss the accurate numerical calculation of the coefficients . Lu [1] calculated them with 17 significant digits using a collocation Treftz method. Li et al. [2] used a special boundary approximation method to derive the leading coefficients of the asymptotic solution expansion. Similar approaches were followed by Li [3] and Arad et al. [4], who employed boundary methods combined with least-squares and penalty method techniques to handle the boundary information, and by Georgiou [5], who employed Lagrange multipliers.

Some other authors focus on solving the Motz problem using different numerical techniques. Bernal and Kindelan [6] thoroughly discussed the performance of the radial basis function method in the solution of the Motz problem, including the global radial basis function method due to Kansa [7], as well as Shu's local differential quadrature method based on radial basis functions [8]. They showed that the exponential convergence obtained by increasing resolution and increasing the shape parameter is lost due to the presence of the singularity. The accuracy of the solution can be increased by using collocation at some boundary nodes to restore the convergence of the method, although it is necessary to use functions similar to the terms in equation (1), which are capable of capturing the behaviour of the solution near the discontinuity.

To produce an accurate solution based on the flexibility of meshfree methods, Hong [9] addressed a singular basis function enrichment technique in a meshless method, using the leading terms of the local series expansion at the boundary singularity, as enrichment functions for a local approximation space constructed to calculate the solution. The enrichment technique is derived from the Generalized Finite Element Method (GFEM) and produces accurate results, although it requires smooth global support basis functions and basis enrichment to produce the numerical solution.

The aforementioned papers give an insight into the different approaches considered in the problem. Finite elements and radial basis functions can be certainly used to produce accurate results near the singularity, although it is also possible to consider the little restrictive and more versatile alternative given by the Generalized Finite Differences Method (GFDM), which has recently been used to produce satisfactory approximations to solve differential equations in a wide number of applications modeled over irregular regions; for instance, on structured grids, a Generalized Crank Nicolson scheme is presented in [10]. Chávez et al showed a stable calculation of the solution of nonlinear Richards equation which avoids oscillations in high gradient regions in [11], and Tinoco et al addressed the conditions for stability in general schemes in [12].

The GFDM is not limited to a grid structure, and can also be applied on unstructured grids and even as a meshless method in clouds of points: Benito et al studied its use on parabolic and hyperbolic equations [13], advection-diffusion equations [14], seismic wave problems [15] and nonlinear elliptic equations [16] . Fan et al have applied general finite differences on Burgers equation [17], obstacle problems [18], shallow water equations [19] and natural double convection problems [20].

In this paper, bearing in mind the engineering applications where the use of GFDM has been successfully applied, this paper shows an accurate numerical approximation to the solution of Motz problem.

It is organized as follows: generalized finite differences are presented in Section 2. The description of the meshes used is given in Section 3. Numerical tests are described in Section 4; discussion and conclusions in Section 5, followed by conclusions in Section 6.

2. General finite differences

A very common and widely used discretization procedure to approximate the solutions of differential equations is the finite difference method. This is because it can be easily applied in block-rectangular regions.

The well-known classic schemes for the first and second-order derivatives can be easily deduced from basic analytical considerations or by straightforward application of the Taylor theorem. Despite the basic idea of the theorem is quite simple, its use leads to more complicated schemes on non-rectangular regions when the geometry associated with the scheme changes. Nevertheless, when symmetry is involved, very interesting approaches can be considered; like those due to Liszka [21] for regular hexagonal geometries and to Benito [22,13] for slightly perturbed clouds of nodes.

In the following paragraphs, we discuss the use of a general difference scheme that can be calculated by considering a set of nodes where the solution of a differential equation is evaluated or approximated. For these nodes it is required to find coefficients in such a way that consistency condition is fulfilled. In other words, the local truncation error to approximate the second-order linear operator

(2)

at the grid point employing the combination

(3)

must satisfy the condition

(4)

as [23].

The coefficients may change from one inner grid point to another in an irregular region, since they depend on the positions of the nodes relative to . Nevertheless, in the case of nearly rectangular grids and the Laplace equation, since the coefficients are invariant under rotations and translations, essentially one set of coefficients is required.

For the sake of brevity, in the following paragraphs the coefficient

will be denoted as .

Let and be the and components of , respectively. Thus, the local truncation error (up to second-order) yields

Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): \left(D(x_0,y_0)-\sum _{i=1}^k\Gamma _i\Delta y_i\right)\frac{\partial u}{\partial y}(p_0)+\left(A(x_0,y_0)-\sum _{i=1}^k\frac{\Gamma _i (\Delta x_i)^2}{2}\right)\frac{\partial ^{2}u}{\partial x^{2}}(p_0)+
Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): \left(-\sum _{i=1}^k\Gamma _i\Delta x_i\Delta y_i\right)\frac{\partial ^{2}u}{\partial x \partial y}(p_0)+\left(B(x_0,y_0)-\sum _{i=1}^k\frac{\Gamma _i (\Delta y_i)^2}{2}\right)\frac{\partial ^{2}u}{\partial y^{2}}(p_0)+
Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): \mathcal{O}\left(\max \{ \Delta x_i, \Delta y_i\} \right)^3.

Therefore, assuming smoothness of Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle u} , to get second-order local consistency the coefficients of the derivatives of Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle u}

at the right-hand side of the expansion must satisfy the linear system

Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): \left( \begin{array}{rrrrr}1& &1 &... & 1\\ 0& &\Delta x_1 &... &\Delta x_q\\ 0& &\Delta y_1 &... &\Delta y_q\\ 0& &(\Delta x_1)^2 &... &(\Delta x_q)^2\\ 0& &\Delta x_1\Delta y_1 &... &\Delta y_q\Delta z_q\\ 0& &(\Delta y_1)^2 &... &(\Delta y_q)^2\\ \end{array} \right) \left( \begin{array}{c}\Gamma _0\\ \Gamma _1\\ \Gamma _2\\ \cdot\\ \cdot\\ \cdot\\ \Gamma _q\\ \end{array} \right) = \left( \begin{array}{c}E(x_0,y_0)\\ C(x_0,y_0)\\ D(x_0,y_0)\\ 2A(x_0,y_0)\\ 0\\ 2B(x_0,y_0)\\ \end{array} \right).
(5)

One must note that, in general, this system is not well-determined. This poses the calculation of the coefficients in terms of a well studied algebraic problem which can be solved robustly. Hence, the main question is how to produce a solution of equation (5) which provides a consistent second-order scheme. Several alternatives are possible, for example, in [10], an efficient heuristic calculation to solve numerically the diffusion equation was proposed. It was based on an unconstrained optimization problem for every inner grid point. To apply it, first, we separate the first equation of the matrix system (5)

Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): \sum _{i=0}^q \Gamma _i=E(x_0,y_0)
(6)

and then we solve the system

Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): M\Gamma =\beta ,
(7)

where

Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): M=\left( \begin{array}{clcl}\Delta x_1 &... &\Delta x_q\\ \Delta y_1 &... &\Delta y_q\\ (\Delta x_1)^2 &... &(\Delta x_q)^2\\ \Delta x_1\Delta y_1 &... &\Delta x_q\Delta y_q\\ (\Delta y_1)^2 &... &(\Delta y_q)^2\\ \end{array} \right),
(8)
Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): \Gamma =\left( \begin{array}{c}\Gamma _1\\ \Gamma _2\\ .\\ .\\ .\\ \Gamma _q\\ \end{array} \right),\quad \beta = \left( \begin{array}{c}C(x_0,y_0)\\ D(x_0,y_0)\\ 2A(x_0,y_0)\\ 0\\ 2B(x_0,y_0)\\ \end{array} \right).
(9)

Next, Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle \Gamma _0}

is obtained from equation (6).

Equation (7) defines the generalized finite differences in general grids. To solve it is straightforward when the grid is rectangular or nearly rectangular, and a sufficient number of nodes or symmetries are available.

In irregular grids or clouds of points, which lack symmetry, the system defined by Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle M}

is underdetermined if Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle q>6}

, which means that, in that case, at least seven neighbor nodes around Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_0}

are required to get local consistency. 

On the other hand, in a rather regular structured -or mapping- grid, symmetries reduce the number of required neighbors. An important consequence is that the symmetric location of nodes around the central node avoids rank deficiency of Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle M}

in the system (7).

3. Grids and stencils

The main issue to consider when solving the Motz problem is the singularity in the partial derivative Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle \frac{\partial u}{\partial x}}

at the origin of coordinates. The usual approach to produce an accurate approximation at that point is to refine the grid around it. However, this causes not only an increase in the computational cost, but often the results are not satisfactory due to the high gradient values of the solution on the positive Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle x}

-axis; even with robust techniques as those based on radial basis functions [6], a special treatment is required near the origin. Some other approaches, like that of Hong, require a rather sophisticated process of patch mapping to produce accurate results [9].

When finite differences are considered, it is important to emphasize that the method has been successfully applied to problems with strong boundary layers, for instance, by Ivanenko [24] and Pereyra [25], when the equidistribution principle is addressed. In other words, the key is to produce a suitably adapted grid through a simple procedure. In this paper, to derive the rule that controls the refinement around the origin, the clue is given by the shape of the significative terms in the closed-form solution series of equation (1).

Bearing in mind the equidistribution principle as discussed by Ivanenko, since the graph of the solution Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle u(x)}

 is  monotonic along the positive Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle x}

-axis, and resembles that of a boundary layer, there is a mapping Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle \xi \mapsto x=\xi ^\alpha }

such that the expression

Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): \frac{\partial u}{\partial \xi }=\frac{\partial u}{\partial x}\frac{\partial x}{\partial \xi }
(10)

is constant. According to eq. (1), each term in the solution along such axis equals Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle x^{(i+\frac{1}{2})}} , Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle i=0,1,...} . Therefore, the constant right-hand side of (10) takes the form

Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): \frac{\partial u}{\partial x}\frac{\partial x}{\partial \xi }=\left[\frac{\partial x}{\partial \xi }\right]\left[\left(i+\frac{i}{2}\right)x^{i-\frac{1}{2}}\right],\qquad i=0,1,...

This implies that the partial derivative Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle \frac{\partial x}{\partial \xi }}

must have the form Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle x^{(\frac{1}{2}-i)}}
and, in consequence, the grid refinement is of the form

Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): \xi \mapsto x=\xi ^{\frac{2}{2i+1}},\qquad i=0,1,...
(11)

According to the estimations presented by Lu [1], the first five coefficients in the equation (1) are the significative ones: Those terms correspond to Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle i=0,1,...,5} . Since the value Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle i>1}

in equation (11) represents an increase in the distance between the grid nodes around the origin, only the value Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle i=0}

, which corresponds to Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle \alpha=2} , is considered.

Thus, the mesh adaptation procedure is as follows:

  1. The grid for the whole domain of the problem is generated in the first step as a regular rectangular grid.
  2. Every cell is subdivided into two triangles in the second step to define six neighbours in the stencil for every inner grid node.
  3. Finally, the grid is adapted by modifying the grid coordinates following the mapping given by equation(11).


In the numerical tests, the interval Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle [1.5,2.5]}

was considered for the exponent Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle \alpha }

. For comparison purposes, the mapping was applied to the Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle x} -axis (Type I grids), to the Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle y} -axis (Type II grids), and both axes (Type III). An example of the grids corresponding to the different distributions is shown in Figure 2.

Draft Dominguez Mota 636399092-quadx.png Draft Dominguez Mota 636399092-quady.png Examples of the kinds of grids used in the tests
Figure 2. Examples of the kinds of grids used in the tests


The three kinds of grids define a regular location of neighbors at the inner nodes: star-shaped subgrids with 6 neighbours that can be used to approximate the derivatives in the left-hand side of equation (6); an example of such subgrid can be seen in Figure 3. The Laplacian operator is approximated at the point Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_0} , and the corresponding neighbours are Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_1} , Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_2} , Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_3} , Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_4} , Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_5} , and Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_6} . The values Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle h_1} , Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle h_2} , Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle k_1} , and Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle k_2}

are the step size lengths in the stencil. Bearing in mind equation (5), this poses a problem with 6 equations and Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle k=7}
unknowns; since the values Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle h_1}

, Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle h_2} , Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle k_1} , and Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle k_2}

are non zero, after preconditioning of equation (6) by scaling, the rank of Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle M}
is equal to five,  which allows a robust calculation of the coefficients and leaves one degree of freedom available. 
Selection of neighbours for an inner grid node
Figure 3. Selection of neighbours for an inner grid node


Regarding the boundary condition, a straightforward discretization using the GFDM method was addressed in [26]: at a boundary node of a structured grid, equation (7) is evaluated at its five neighbors, taking the right-hand side Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle \beta }

to define the directional derivative. However, since the neighbors are often not balanced around the boundary node, such selection might produce first approximations which could cause loss of order in the Motz problem with the grids considered here.

In this paper, we define an ad hoc ghost point to avoid loss of accuracy. The boundary condition

Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): \frac{\partial u}{\partial n}(p_0)=l_x\frac{\partial u}{\partial x}+l_y\frac{\partial u}{\partial y}

where Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle n=(l_x,l_y)}

is the outer normal vector at Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_0}

, is approximated using the scheme of equation (5) with Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle A=B=E=0}

and Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle C=l_x}

, Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle D=l_y}

using 7 points: a ghost point Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_g}

, and the nodes Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_0} , Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_1} , Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_2} , Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_3} , Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_4}

and Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_5}
(Figure 4), which yields the equation (12)

Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): l_x\frac{\partial u}{\partial x}+l_y\frac{\partial u}{\partial y}\approx \hat{\Gamma }_g u(p_g)+\sum _{l=0}^k\hat{\Gamma }_l u(p_l),
(12)

or, in other words,

Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): u(p_g)=\frac{1}{\hat{\Gamma }_g}\frac{\partial u}{\partial n}(p_0)-\sum _{l=0}^k\frac{\hat{\Gamma }_l}{{\hat{\Gamma }_g}} u(p_l).
(13)

For convenience, the ghost point Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_g}

is calculated in such a way that Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_0}
is the midpoint of Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_g}
and Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_c}

, the latter being the centroid of the nodes Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_0} , Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_1} , Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_2} , Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_3} , Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_4} , and Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_5}

(Figure 4).

The value of Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle u(p_g)}

given by the equation (13) is substituted in the approximation to the Laplacian (equation (7) at Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_0}
using the same 7 points Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_g}

, Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_0} , Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://mathoid.scipedia.com/localhost/v1/":): {\textstyle p_1} , , , , and , which eliminates the value of .

Te assembled system for the whole grid becomes sparse as the number of grid points per side increases; it was solved using a sparse lu decomposition to produces the approximation to the solution of the problem.

Chavez et al 2020a 2003 498px-Draft Dominguez Mota 636399092-ghost.png
Figure 4. Selection of ghost point at a boundary node of a general structured grid with Neumann condition

4 Numerical tests

Several tests were performed with and nodes per side on the and directions, respectively, varying from 31 to 351. Grids of the three types in the previous section were generated, and the Motz problem was solved using the proposed scheme for every grid size and distribution. An example of the finite difference discretization described in the previous section is sketched in Figure 5.

Grid and calculated solution
Figure 5. Grid and calculated solution


To assess the accuracy of the numerical solution, since the approximation was calculated using structured grids, the quadratic error

(14)

was used, where is the solution calculated using the difference scheme at the point , is the exact solution (equation (1)) at the same point, and is the area of the rhombus defined by the nodes .

The measure of the roughness of the grid is given by the meh norm ; the log graph of the mean square errors versus the grid norm for the type I and II grids are shown in Figures 6 and 7. One can note that there is a clear convergence tendency on all the tests despite the boundary singularity. But, when only a nonuniform distribution either along the or the -axis is used, the second-order approximation is lost, in a similar way as reported by Bernal using radial basis functions. The best results are obtained with the grids which are non-uniform on and , as can be seen in Figure 7.

Draft Dominguez Mota 636399092-figErrorX20.png Log values of the quadratic errors, nonuniform node distribution on x or y
Figure 6. Log values of the quadratic errors, nonuniform node distribution on or


Log values of the quadratic errors, nonuniform node distribution on x and y
Figure 7. Log values of the quadratic errors, nonuniform node distribution on and


The values of the errors, grid norms, and estimated convergence order are shown in Tables 1 to 3. The values in the last table show that, indeed, quadratic convergence is obtained when the same nonuniform node distribution along both axes is used, even in the presence of the boundary singularity.

Table 1. Quadratic errors and approximation order for type I grids
Order
1.00 31 61 3.3333E-02 6.4941E-03 ––
1.00 71 141 1.4286E-02 2.7556E-03 1.01
1.00 111 221 9.0909E-03 1.7489E-03 1.01
1.00 151 301 6.6667E-03 1.2810E-03 1.00
1.00 191 381 5.2632E-03 1.0106E-03 1.00
1.00 231 461 4.3478E-03 8.3446E-04 1.00
1.00 271 541 3.7037E-03 7.1060E-04 1.00
1.00 311 621 3.2258E-03 6.1877E-04 1.00
1.25 31 61 4.1492E-02 6.4057E-03 ––
1.25 71 141 1.7825E-02 2.7333E-03 1.01
1.25 111 221 1.1351E-02 1.7503E-03 0.99
1.25 151 301 8.3264E-03 1.2920E-03 0.98
1.25 191 381 6.5746E-03 1.0260E-03 0.98
1.25 231 461 5.4318E-03 8.5207E-04 0.97
1.25 271 541 4.6275E-03 7.2925E-04 0.97
1.25 311 621 4.0306E-03 6.3784E-04 0.97


Table 2. Quadratic errors and approximation order for type II grids
Order
1.00 31 61 0.033333 0.0064941 ––
1.00 71 141 0.014286 0.0027556 1.01
1.00 111 221 0.0090909 0.0017489 1.01
1.00 151 301 0.0066667 0.001281 1.00
1.00 191 381 0.0052632 0.0010106 1.00
1.00 231 461 0.0043478 0.00083446 1.00
1.00 271 541 0.0037037 0.0007106 1.00
1.00 311 621 0.0032258 0.00061877 1.00
1.25 31 61 0.041492 0.0060306 ––
1.25 71 141 0.017825 0.0026189 0.99
1.25 111 221 0.011351 0.0016827 0.98
1.25 151 301 0.0083264 0.001243 0.98
1.25 191 381 0.0065746 0.00098714 0.98
1.25 231 461 0.0054318 0.00081951 0.97
1.25 271 541 0.0046275 0.00070109 0.97
1.25 311 621 0.0040306 0.00061294 0.97


Table 3. Quadratic errors and approximation order for type III grids
Order
2.00 71 141 2.8367E-02 9.3949E-04 1.77
2.00 111 221 1.8099E-02 4.1865E-04 1.80
2.00 151 301 1.3289E-02 2.3926E-04 1.81
2.00 191 381 1.0499E-02 1.5584E-04 1.82
2.00 231 461 8.6767E-03 1.1006E-04 1.82
2.00 271 541 7.3937E-03 8.2131E-05 1.83
2.00 311 621 6.4412E-03 6.3786E-05 1.83
2.25 31 61 7.3442E-02 4.5265E-03 ––
2.25 71 141 3.1856E-02 9.2020E-04 1.91
2.25 111 221 2.0338E-02 3.8741E-04 1.93
2.25 151 301 1.4938E-02 2.1315E-04 1.94
2.25 191 381 1.1803E-02 1.3496E-04 1.94
2.25 231 461 9.7560E-03 9.3190E-05 1.94
2.25 271 541 8.3140E-03 6.8259E-05 1.95
2.25 311 621 7.2434E-03 5.2179E-05 1.95
2.50 31 61 8.1262E-02 4.9361E-03 ––
2.50 71 141 3.5333E-02 9.2007E-04 2.02
2.50 111 221 2.2573E-02 3.7194E-04 2.02
2.50 151 301 1.6583E-02 1.9951E-04 2.02
2.50 191 381 1.3106E-02 1.2407E-04 2.02
2.50 231 461 1.0834E-02 8.4506E-05 2.02
2.50 271 541 9.2336E-03 6.1224E-05 2.02
2.50 311 621 8.0450E-03 4.6380E-05 2.02
2.75 31 61 8.9015E-02 5.3606E-03 ––
2.75 71 141 3.8796E-02 9.3491E-04 2.10
2.75 111 221 2.4802E-02 3.6799E-04 2.08
2.75 151 301 1.8227E-02 1.9446E-04 2.07
2.75 191 381 1.4407E-02 1.1975E-04 2.06
2.75 231 461 1.1911E-02 8.0998E-05 2.06
2.75 271 541 1.0152E-02 5.8375E-05 2.05
2.75 311 621 8.8459E-03 4.4042E-05 2.05
3.00 31 61 9.6704E-02 5.7955E-03 ––
3.00 71 141 4.2248E-02 9.6029E-04 2.17
3.00 111 221 2.7026E-02 3.7142E-04 2.13
3.00 151 301 1.9867E-02 1.9454E-04 2.10
3.00 191 381 1.5707E-02 1.1915E-04 2.09
3.00 231 461 1.2987E-02 8.0299E-05 2.08
3.00 271 541 1.1070E-02 5.7715E-05 2.07
3.00 311 621 9.6462E-03 4.3453E-05 2.06

5. Results and discussion

Several assertions can be drawn from the numerical results. First, despite the presence of the boundary singularity, the results are accurate enough. This is certainly a strength of the proposed finite difference scheme, which is rather simple. Since the main differences between the scheme discussed in this paper and the standard five-point scheme for the Laplacian are the addition of two extra “diagonal” nodes in the inner grid nodes and the use of a suitable grid spacing, this suggests that the proposed stencils are suitable to produce accurate results.

The convergence orders in table 3 must be highlighted. Having in mind the presence of the singularity and the simplicity of the addressed difference stencils, this is indeed a remarkable result. Since the boundary node redistribution was defined by a very simple rule, this suggests that it is unnecessary to consider another kind of grids or more complex refinements near the singularity to solve the Motz problem, although, it must be acknowledged that the mesh adapting strategy was proposed ad hoc for it. Naturally, the numerical treatment of different singularities in other problems require a different application of the equidistribution principle.

6. Conclusions

Indeed, generalized finite differences can be used to solve the Motz problem numerically. The numerical results show that generalized finite differences produce accurate solutions to a complicated benchmark problem, and the most noteworthy fact is that this can be done even with very simple structured grids. This is an important feature of the GFDM method, considering the difficulty of the problem. Thinking about the solution to several linear and nonlinear problems with this method which appear in the literature [27,10,22,13], it is certainly a reliable numerical tool.

The results are very encouraging since boundary configurations similar to those of Motz problem appear in several relevant applications, for instance, in interfaces between different media in infiltration problems, or advection-diffusion problems with sources and sinks along the boundary. Both problems are currently under investigation and advances in this sense will be reported in a future paper.

Acknowledgments

We want to thank CIC-UMSNH 9.16, and Aula-CIMNE Morelia for the financial support for this work.

References

[1] Lu T.T., Hu H.Y., Li Z.C. Highly accurate solutions of Motz's and the cracked beam problems. Engineering Analysis with Boundary Elements, 28:1387–1403, 2004.

[2] Li Z.C., Mathon R., Sermer P. Boundary methods for solving elliptic problems with singularities and interfaces. SIAM Journal of Numerical Analysis, 24:487–498, 1987.

[3] Li Z.C. An approach for combining the Ritz-Galerkin and Finite-Element Methods. Journal of Approximation Theory, 39:132-152, 1983.

[4] Arad M., Yosibash Z., Ben-Dor G., Yakhot A. Computing flux intensity factors by a boundary method for elliptic equations with singularities. Communications in Numerical Methods in Engineering, 14:657–670, 1998.

[5] Georgiou G.C., Olson L., Smyrlis Y.S. A singular function boundary integral method for the Laplace equation. Communications in Numerical Methods in Engineering, 12:127-134, 1996.

[6] Bernal F., Kindelan, M. Radial basis function solution of the Motz problem. Engineering Computations, 27:606–620, 2010.

[7] Kansa E.J. Multiquadrics, a scattered data approximation scheme with applications to computational fluid dynamics II: solutions to parabolic, hyperbolic and elliptic partial differential equations. Computers and Mathematics with Applications, 19:147–161, 1990.

[8] Shu C., Ding H., Yeo K.S. Local radial basis function-based differential quadrature method and its application to solve two-dimensional incompressible Navier-Stokes equations. Computer Methods in Applied Mechanics and Engineering, 192:941–954, 2003.

[9] Hong W.T. Enriched meshfree method for an accurate numerical solution of the Motz problem. Advances in Mathematical Physics, 2016:1-12, 2016.

[10] Domínguez-Mota F.J., Mendoza-Armenta S., Tinoco-Guerrero G., Tinoco-Ruiz J.G. Finite difference schemes satisfying an optimality condition for the unsteady heat equation. Mathematics and Computers in Simulation, 106:76–83, 2014.

[11] Chávez-Negrete C., Domínguez-Mota F.J., Santana-Quinteros D. Numerical solution of Richards' equation of water flow by generalized finite differences. Computers and Geotechnics, 101:168–175, 2018.

[12] Tinoco-Guerrero G., Domínguez-Mota F.J., Gaona-Arias A., Ruiz-Zavala M.L., Tinoco-Ruiz J.G. A stability analysis for a generalized finite-difference scheme applied to the pure advection equation. Mathematics and Computers in Simulation, 147:293–300, 2018.

[13] Benito J.J., Ureña F., Gavete L. Solving parabolic and hyperbolic equations by the generalized finite difference method. Journal of Computational and Applied Mathematics, 209:208–233, 2007.

[14] Benito J.J., Ureña F., Gavete, L. Application of the generalized finite difference method to solve the advection-diffusion equation. Journal of Computational and Applied Mathematics, 235:1849–1855, 2011.

[15] Benito J.J., Ureña F., Salete E., Gavete, L. A note on the application of the generalized finite difference method to seismic wave propagation in 2D. Journal of Computational and Applied Mathematics, 236:3016–3025, 2012.

[16] Gavete L., Ureña F., Benito J.J., García A., Ureña M., Salete E. Solving second order non-linear elliptic partial differential equations using generalized finite difference method. Journal of Computational and Applied Mathematics, 318:378–387, 2017.

[17] Fan C.M., Li P.W. Generalized finite difference method for solving two-dimensional Burgers' equations. Procedia Engineering, 79:55–60, 2014.

[18] Chan H.F., Fan C.M., Kuo C.W. Generalized finite difference method for solving two-dimensional non-linear obstacle problem. Engineering Analysis with Boundary Elements, 37:1189–1196, 2013.

[19] Li P.W., Fan C.M. Generalized finite difference method for two-dimensional shallow water equations. Engineering Analysis with Boundary Elements, 80:58–71, 2017.

[20] Li P.W. , Chen W., Fu Z.J., Fan C.M. Generalized finite difference method for solving the double-diffusive natural convection in fluid-saturated porous media. Engineering Analysis with Boundary Elements 95:175–186, 2018.

[21] Liszka T., Orkisz J. The finite difference method at arbitrary irregular grids and its application in applied mechanics. Computers and Structures, 11:83–95, 1980.

[22] Benito J.J., Ureña F., Gavete L., Alonso B. A posteriori error estimator and indicator in generalized finite differences. Application to improve the approximated solution of elliptic PDEs. International Journal of Computer Mathematics, 85:359–370, 2008.

[23] Celia, M., Gray, W. Numerical Methods for Differential Equations: fundamental concepts for scientific and engineering applications. Prentice-Hall, 1992.

[24] Ivanenko, S. In Selected chapters on grid generation and applications. B. Azarenok (Ed.), Russian Academy of Sciences, 2009.

[25] Pereyra V., Sewell G. A singular function boundary integral method for the Laplace equation. Numer. Math., 23:261-268, 1975.

[26] Domínguez-Mota F.J., Mendoza S., Tinoco-Ruiz J.G, Tinoco-Guerrero G. Numerical solution of Poisson-like equations with Robin boundary conditions using a finite difference scheme defined by an optimality condition. Proceedings of MASCOT 2011, 17, 101-110, 2011.

[27] Domínguez-Mota F.J., Equihua M., Mendoza S., Tinoco-Ruiz J.G. Solución de ecuaciones diferenciales elípticas en regiones planas irregulares usando mallas convexas generadas por métodos variacionales empleando elementos finitos. Rev. Int. Mét. Num. Cálc. Dis. Ing., 26:187–194, 2010.

Back to Top

Document information

Published on 14/01/21
Accepted on 02/01/21
Submitted on 12/04/20

Volume 37, Issue 1, 2021
DOI: 10.23967/j.rimni.2021.01.004
Licence: CC BY-NC-SA license

Document Score

0

Views 197
Recommendations 0

Share this document