Latest revision as of 08:21, 12 September 2019

Abstract: A meshless method with ridge basis functions for solving the 2D two-flow domain model problem is proposed in this paper, which is using ridge basis functions to construct the approximate functions, and using collocation method to discretize the governing equation. The existence and uniqueness of the solution to the mothod is analyzed. Compared with the traditional finite difference method, the new method reduces the calculation error. Numerical results show that the proposed method is convergent and can improve the error precision.

Keywords: Two-flow domain model, Ridge basis functions, Meshless method, Collocation method.

1. Introduction

The migration of solutes in the soil is an important part of the biosphere cycle. By generalizing the flow velocity characteristics of soil pores, various mathematical models are proposed for description. For undisturbed soil, there are a certain number of lardge pores. The convective dispersion model^[1] cannot explain the early penetration and tailing of the solute breakthrough curve. Therefore, Coats and Smith established a mobile and immobile model^[2] to better describe the characteristics of solute migration. The mobile and immobile model exposed its limitations with the emergence of double peaks in the solute breakthrough curves and preferential flow in the soil. To this end, Skopp proposed a two-flow domain model based on the mobile and immobile model^[3]. The model problem is often a partial differential equation system, and it is generally difficult to obtain the corresponding analytical solution. Thus, the numerical method is currently the most effective method to solve it.

Meshless method, which is based on point approximation and does not need initial division and reconstruction of grid, can completely or partially eliminate the grid. It is an efficient numerical method for solving partial differential equations in recent years^[4]. With the development of meshless method, the current research based on radial basis functions has a series of results^[5-8], but the research based on ridge basis functions is rare in theory and application. Gordon^[9] considered the approximation of the norm space ridge basis functions; Ismailov^[10] gave the best approximation of multivariate functions in L₂ space by linear combination of ridge basis functions on some sets in Rⁿ; After that, Zhang Liwei^[11-12] gave the existence of the interpolation problem of the ridge basis functions in 2D space, estimated its error, and proved ridge basis functions approximation of the linear combination of finite plane waves; Shu Hengmu et al.^[13] constructed an approximation function satisfying the interpolation theory of ridge basis functions, adopted the meshless method of the weighted least squares discrete equation and applied it to the elastic static problem; In recent years, Wang Zhigang et al.^[14] discussed the collocation meshless method of ridge basis functions and gave the existence and uniqueness of the numerical solution. Both ridge basis functions and radial basis functions have the advantages of simple form and isotropy, but the application of the ridge basis functions is more extensive, such as the development equation, the solution of the dynamic system problem. In addition, the selection of the basis function and its parameters, the arrangement of the nodes have direct impact on the interpolation accuracy of radial basis functions and ridge basis functions. However, the direction selection has no effect on radial basis functions, but will affect the interpolation accuracy of ridge basis functions. This is one of the advantages of its application.

This paper mainly considers the 2D two-flow domain model problem, and its equation is expressed

where i is the flow rate area, subscripts A, B are two flow regions, ${\displaystyle \Omega =\left[0,L\right]\times \left[0,L\right]}$ , ${\displaystyle D_{A}}$ , ${\displaystyle D_{B}}$ , ${\displaystyle V_{A}}$ , ${\displaystyle V_{B}}$ , ${\displaystyle {\theta }_{A}}$ , ${\displaystyle {\theta }_{B}}$ , ${\displaystyle \alpha }$ are constants greater than zero.

The remainder of this paper is organized as follows. Section 2 mainly introduces the meshless method with ridge basis functions. Firstly, definitions of the ridge basis function and its interpolation problem are given. Then, the implementation process of the collocation method is described. Finally, the discrete process of 2D two-flow domain model is deduced. Section 3 analyzes the existence and uniqueness of the solution to the proposed method. Section 4 gives numerical examples of 1D and 2D two-flow domain model to verify the effectiveness of the ridge basis function meshless method. Conclusions and future work are given at the end of the paper.

2. Meshless method with ridge basis functions

2.1 Ridge basis functions

Definition 2.1 If there exsit unary function ${\displaystyle \phi :R\rightarrow R}$ and vector ${\displaystyle a\in R^{d}\backslash \left\{0\right\}}$ such that ${\displaystyle f\left(x\right)=\phi \left(a\cdot x\right)}$ holds for any ${\displaystyle x\in R^{d}}$ , then ${\displaystyle f\left(x\right):R^{d}\rightarrow R}$ is a ridge basis function. If ridge basis functions ${\displaystyle f_{1},f_{2},\cdots ,f_{n}}$ are basises of the space, then the space is called the ridge basis function space.

The ridge basis function uses the unary function ${\displaystyle \phi }$ to act on the inner product ${\displaystyle a\cdot x}$ . It can be seen that the function space is a linear combination of some plane waves. Define the ridge basis function interpolation as follows:

Definition 2.2 Give ${\displaystyle n}$ different nodes ${\displaystyle {\left\{x_{i}\right\}}_{i=1}^{n}}$ and ${\displaystyle n}$ data points ${\displaystyle {\left\{d_{i}\right\}}_{i=1}^{n}}$ in ${\displaystyle R^{d}}$ , the function ${\displaystyle F(x)}$ statisfying

Definition 2.2 Give ${\displaystyle n}$ different nodes ${\displaystyle {\left\{x_{i}\right\}}_{i=1}^{n}}$ and ${\displaystyle n}$ data points ${\displaystyle {\left\{d_{i}\right\}}_{i=1}^{n}}$ in ${\displaystyle R^{d}}$ , the function ${\displaystyle F(x)}$ statisfying

where ${\displaystyle {\lambda }_{k}}$ is unknown coefficient, ${\displaystyle a_{v}}$ is fixed direction, ${\displaystyle m}$ is the total number of the disk directions, and ${\displaystyle \phi }$ is a ridge basis function, then the ridge basis function interpolation problem becomes the following system of linear equations:

This system has a unique solution if and only if its coefficient matrix is non-singular, which is related with the distribution of nodes and the selection of direction ${\displaystyle a_{v}}$ ^[14].

In the numerical calculation, the ridge basis functions commonly used are

Gaussians: ${\displaystyle \phi \left(r\right)=exp\left(-c{\left(a\cdot r\right)}^{2}\right)}$

Multi-Quadrics(MQ): ${\displaystyle \phi \left(r\right)={\left(c^{2}+{\left(a\cdot r\right)}^{2}\right)}^{\beta }}$

Markoff distribution functions: ${\displaystyle \phi \left(r\right)=exp\left(-c\vert a\cdot r\vert \right)}$

Thin-Plate Spline(TPS): ${\displaystyle \phi \left(r\right)={\left(a\cdot r\right)}^{2\beta }log\left(a\cdot r\right)}$

Inverse Multi-Quadrics(IMQ): ${\displaystyle \phi \left(r\right){\mbox{=}}{\left(c^{2}+{\left(a\cdot r\right)}^{2}\right)}^{-\beta }}$

where ${\displaystyle c}$ and ${\displaystyle \beta }$ are constants ( ${\displaystyle c>0}$ ), ${\displaystyle r=\Vert x-x_{k}\Vert }$ ..

2.2 Collocation method

The collocation method is a true meshless method, which does not need any grid to calculate the integral. Its basic idea is that each node in the domain satisfies the governing equation, and each node on the boundary satisfies the uniform boundary condition. Generally, the governing equation is set as

Dirichlet boundary condition is

Neumann boundary condition is

where ${\displaystyle L}$ , ${\displaystyle B}$ are differential operators, ${\displaystyle \partial \Omega ={\Gamma }_{1}\cup {\Gamma }_{2}=\varnothing }$ is boundary, ${\displaystyle {\overline {u}}}$ is the specified value of the unknown function on the boundary ${\displaystyle {\Gamma }_{1}}$ , ${\displaystyle u\left(x\right)}$ is an unknown function that needs to be solved.

Let the approximate expression of the unknown function ${\displaystyle u\left(x\right)}$ be ${\displaystyle u^{h}\left(x\right)}$ , using the weighted residual value method, Eq.(5)~Eq.(7) can be replaced by the following formula

where the weighting functions ${\displaystyle w_{1}}$ , ${\displaystyle w_{2}}$ , ${\displaystyle w_{3}}$ can have different definitions. Define ${\displaystyle {\delta }_{i}=w_{1}=w_{2}=w_{3}}$ , ${\displaystyle {\delta }_{i}}$ is Dirac ${\displaystyle \delta }$ function, thereby the following collocation type equation can be obtained

where ${\displaystyle n_{1}}$ is the number of nodes on the boundary ${\displaystyle {\Gamma }_{1}}$ , ${\displaystyle n_{2}}$ is the number of nodes on the boundary ${\displaystyle {\Gamma }_{2}}$ , ${\displaystyle n_{0}=N-n_{1}-n_{2}}$ is the number of nodes in the domain ${\displaystyle \Omega }$ , ${\displaystyle N}$ is total number of nodes. Thus, Eq.(8)~Eq.(10) can be expressed in matrix form as

2.3 The construction of meshless method with ridge basis functions

Let ${\displaystyle C_{A}^{n}=C_{A}^{n}\left(X\right)}$ be the numerical solution of ${\displaystyle C_{A}\left(X,t_{n}\right)}$ , ${\displaystyle C_{B}^{n}=C_{B}^{n}\left(X\right)}$ be the numerical solution of ${\displaystyle C_{B}\left(X,t_{n}\right)}$ , where ${\displaystyle t_{n}=n\Delta t}$ , ${\displaystyle n=1,2,\cdots ,s}$ , ${\displaystyle \Delta t=T/s}$ , then Eq.(1) can be discretized in the time leyer as

If ${\displaystyle \phi \left(\cdot \right)}$ is not a polynomial, then ridge basis functions can approximate to almost all the functions. Here, we respectively denote functions ${\displaystyle C_{A}^{h}\left(X\right)}$ , ${\displaystyle C_{B}^{h}\left(X\right)}$ as approximations of ${\displaystyle C_{A}\left(X\right)}$ , ${\displaystyle C_{B}\left(X\right)}$ :

where ${\displaystyle X=\left(x,y\right)}$ , ${\displaystyle I_{i}\in \Omega }$ , ${\displaystyle B_{i}\in \partial \Omega }$ , ${\displaystyle N_{i}}$ is the number of nodes in the region, ${\displaystyle N_{b}}$ is the number of nodes on the boundary, ${\displaystyle {\alpha }_{1i}}$ , ${\displaystyle {\beta }_{1i}}$ , ${\displaystyle {\alpha }_{2i}}$ , ${\displaystyle {\beta }_{1i}}$ are unknown coefficients, ${\displaystyle a_{v}}$ is a fixed direction, ${\displaystyle m}$ is the total number of the disk directions, ${\displaystyle \phi }$ is a ridge basis function.

Let ${\displaystyle \varphi \left(\Vert X-I\Vert \right)=\sum _{v=1}^{m}\phi \left(a_{v}\cdot \left(X-\right.\right.}$${\displaystyle \left.\left.I\right)\right)}$ , then

By adopting ridge basis functions collocation method, it is required that Eq.(11) and Eq.(12) satisfy the partial differential equation in the solution domain, which can be obtained

Thus there is

and let

Then Eq.(13) can be expressed in the following matrix equcation

where ${\displaystyle H_{1}}$ , ${\displaystyle H_{2}}$ are ${\displaystyle \left(N_{b}+N_{i}\right)\times \left(N_{b}+N_{i}\right)}$ matrices, ${\displaystyle a}$ , ${\displaystyle b}$ are the undetermined coefficient vectors, ${\displaystyle F_{1}}$ , ${\displaystyle F_{2}}$ are the given right end term vectors, namely

${\displaystyle H_{2}=\left[{\begin{array}{llllll}\varphi \left(\Vert B_{1}-B_{1}\Vert \right)&\cdots &\varphi \left(\Vert B_{1}-B_{N_{b}}\Vert \right)&\varphi \left(\Vert B_{1}-I_{1}\Vert \right)&\cdots &\varphi \left(\Vert B_{1}-I_{N_{i}}\Vert \right)\\{\mbox{          }}{\mbox{ }}\vdots &&{\mbox{          }}{\mbox{ }}\vdots &{\mbox{          }}{\mbox{ }}\vdots &&{\mbox{          }}{\mbox{ }}\vdots \\\varphi \left(\Vert B_{N_{b}}-B_{1}\Vert \right)&\cdots &\varphi \left(\Vert B_{N_{b}}-B_{N_{b}}\Vert \right)&\varphi \left(\Vert B_{N_{b}}-I_{1}\Vert \right)&\cdots &\varphi \left(\Vert B_{N_{b}}-I_{N_{i}}\Vert \right)\\{\psi }_{B}\left(\Vert I_{1}-B_{1}\Vert \right)&\cdots &{\psi }_{B}\left(\Vert I_{1}-B_{N_{b}}\Vert \right)&{\psi }_{B}\left(\Vert I_{1}-I_{1}\Vert \right)&\cdots &{\psi }_{B}\left(\Vert I_{1}-I_{N_{i}}\Vert \right)\\{\mbox{           }}{\mbox{ }}\vdots &&{\mbox{           }}{\mbox{ }}\vdots &{\mbox{           }}{\mbox{ }}\vdots &&{\mbox{           }}{\mbox{ }}\vdots \\{\psi }_{B}\left(\Vert I_{N_{i}}-B_{1}\Vert \right)&\cdots &{\psi }_{B}\left(\Vert I_{N_{i}}-B_{N_{b}}\Vert \right)&{\psi }_{B}\left(\Vert I_{N_{i}}-I_{1}\Vert \right)&\cdots &{\psi }_{B}\left(\Vert I_{N_{i}}-I_{N_{i}}\Vert \right)\end{array}}\right],}$