This second Note is dedicated to the discussion of the case of 3-rd type boundary conditions posed on the edges of the modeling region and will present an outline of the algorithm for the account of both the Dirichlet and 3-rd type boundary conditions on different edges.
Before passing to the discussion of that scheme, a few words should be said about another ADI scheme – namely 2D Peaceman – Rachford finite difference scheme. The correct account of the intermediate boundary conditions for both Dirichlet and 3-rd order types has been thoroughly discussed in . One essential difference between the 3D Douglas – Rachford and 2D Peaceman – Rachford schemes is that in the latter scheme the spatial coordinates enter symmetrically, such that the second equation (in 2D Peaceman – Rachford) contains finite difference representations of both spatial derivatives. This fact leads to a cumbersome infer of relation between the values of the unknown function on different edges – for details see .
Surprisingly, although the Douglas – Rachford scheme is designed for a 3D spatial region, it is easier to treat the 3-rd order boundary conditions for it correctly due to the “non-symmetrical” entrance of the spatial coordinates in three equations – see Eq. (1) – (3):
where - is the Laplace operator; - is the unknown function, defined on the following mesh:
где , (5)
i,j,s- - are the nodes numbers of the spatial mesh along the X, Y and Z directions respectively. - are the total numbers of nodes of the spatial mesh along the X, Y and Z directions respectively. - are the steps values along the directions X, Y and Z respectively. - is the value of time sub-step : the full time step is equal to.
Let the following boundary conditions being posed on the faces of a rectangular modeling region (3-rd order boundary conditions on all the faces):
In eq. (7) – (12)the nodes with indexes are fictitious nodes which lay outside the modeling region and are introduced with the purpose to retain the second order accuracy of the Douglas – Rachford scheme.
Let us now re-write the Eqs. (1), (2) and (3) in such a form, that the unknowns , and stand in l.h.s.:
and determine which quantities are not known explicitly.
Now, let us consider Eq. (14). In order to find the values of , one has to solve a system of linear algebraic equations:
… и т.д. .
In its turn, the solution of this system requires the knowledge of and , which can be found from the Eq. (9):
While one can reasonably use Eq. (19), a natural question arises with the Eq. (20) – how one should interpret it for the time step ? - is not the true unknown function u at an intermediate time step !
In order to overcome this difficulty, it would be naturally to take the known boundary condition (9) at time steps k and k+3 and, thus, to obtain (9) at k+2:
where the parameters and T are taken at time moment k+3. This choice is based on the fact that our goal is to find .
After that, Eq. (17) takes the following form:
In absolute analogy, consider Eq. (15). In order to find the values of , one has to solve a system of linear algebraic equations:
… and so on.
Finally, let us consider Eq. (13). In order to find the values of , one has to solve a system of linear algebraic equations:
... adn so on . (28)
In Eq. (27) the values of , , and at can be found from the corresponding 3-rd order boundary conditions. In case of the Dirichlet boundary conditions on one of the faces, the range of indexes in (28) must be changed appropriately (e.g.).
As one can see from the discussion above, the correct account of intermediate boundary conditions in case of the 3-rd order boundary conditions posed of the face of the modeling region is not so cumbersome as in case of the Dirichlet boundary conditions.
Now, let us discuss the account of intermediate boundary conditions in case of both Dirichlet and 3-rd order BC’s on different faces. Let the following boundary conditions being posed on the faces of a rectangular modeling region:
i.e., the Dirichlet boundary conditions are posed on the Y faces and the 3-rd order boundary conditions – on all the other.
In accordance with the above presented discussion, the algorithm for the account of intermediate boundary conditions is the following:
1. Find the values of at using the following relation (Eq. (19) of Note 1):
and at from the boundary conditions (29) and (30).
2. Find the values of при using the following relation (Eq. (18) of Note 1):
This is possible because we know the necessary values of and at .
3. Solve the Eq. (13). One now has everything in order to find at .
4. Solve Eq. (14) by finding solution of a system (16) – (18).
5. Solve Eq. (15) by finding solution of a system (23) – (25).
Thus, following the steps 1 – 5, one finds the solution of the heat equation on the next time step.
One should note that the discussed above method for the account of intermediate boundary conditions has the order of approximation with respect to time variable not higher than . In order to improve this approximation one must use a more sophisticated method for the account of boundary conditions.
The author acknowledges T.I. Lakoba (University of Vermont) for helpful discussions.
1. Lakoba, T. (13 04 2012 r.). Lecture notes on Numerical Differential Equations. Vermont, USA: University of Vermont. – internet resource. http://www.cems.uvm.edu/~tlakoba/