CC BY
16
Математическое моделирование
UDC 519.632.4
Solution of the Boundary-Value Problem for a Systems of ODEs of Large Dimension: Benchmark Calculations in the Framework
of Kantorovich Method
A. A. Gusev*, O. Chuluunbaatar*, S. I. Vinitsky*+, V. L. Derbov*
* Joint Institute for Nuclear Research, Dubna, Russia ^ RUDN University, Moscow, Russia * Saratov State University, Saratov, Russia
We present benchmark calculations of the boundary-value problem (BVP) for a systems of second order ODEs of large dimension with help of KANTBP program using a finite element method. In practice, for solving the BVPs with the long-range potentials and a large number of open channels there is a necessity of solving boundary value problems of the large-scale systems of differential equations that require further investigation of convergence and stability of the algorithms and programs. With this aim we solve here the eigenvalue problem for an elliptic differential equation in a two-dimensional domain with Dirichlet boundary conditions. The solution is sought in the form of Kantorovich expansion over the parametric basis functions of one of the independent variables with the second variable treated as a parameter. The basis functions are calculated in an analytical form as solutions of the auxiliary parametric Sturm-Lioville problem for a second-order ODE. As a result, the two-dimensional problem is reduced to a boundary-value problem for a set of self-adjoint second-order ODEs for functions of the second independent variable. The discrete formulation of the problem is implemented using the finite element method. The efficiency, stability and convergence of the calculation scheme is shown by benchmark calculations for a triangle membrane with a degenerate spectrum.
Key words and phrases: benchmark calculations, boundary-value problem, large-scale systems of ODEs, Kantorovich method, finite element method
1. Introduction
The solving quantum tunneling problem or calculations of spectral and optical properties of electronic states in axially symmetric quantum dots and Helium-like atom (system of two-electron in the Coulomb field) is reduced to the solution of boundary-value problems (BVP) for elliptic differential equations with nonseparable variables in a finite domain [1-3]. One of the ways to solve these problems is implemented as a set of programs ODPEVP-POTHEA-KANTBP [4-6] basing on the Kantorovich method (KM) that provides the reduction of the initial problem to a set of self-adjoint second-order ODEs [7] with further discretization by the finite element method (FEM) [8]. In practice, for solving problems with the long-range potential and a large number of open channels there is a necessity of solving boundary value problems of the large-scale systems of the ODEs that require further investigation of convergence and stability of the algorithms and programs.
Testing such approach, algorithms and programs for the solution of two-dimensional BVPs and large-scale systems of the ODEs is the aim of the present work. We present a computational scheme for solving the eigenvalue problem for an elliptic differential equation in a two-dimensional finite domain with Dirichlet boundary conditions. The solution is sought in the form of Kantorovich expansion over the basis functions of one of the independent variables with the second variable treated as a parameter. The basis functions are calculated as a solution of the parametric eigenvalue problem for an ordinary
Received 12th July, 2016. The work was partially supported by the Russian Foundation for Basic Research (RFBR) (grant 14-01-00420). The reported study was partially funded within the Agreement No. 02.a03.21.0008 dated 24.11.2016 between the Ministry of Education and Science of the Russian Federation and RUDN University.
second-order differential equation. Finally, the initial problem is reduced to a BVP for a set of self-adjoint second-order differential equations for functions of the second independent variable. The discretization of the problems is carried out using the FEM with Lagrange interpolating polynomials. The result is used to formulate a generalized algebraic eigenvalue problem for higher-order matrices. We demonstrate the efficiency of the KANTBP program for solving the boundary-value problem for a systems of the ODEs of large dimension in benchmark calculations for the exactly solvable eigenvalue problem of a triangle membrane with the degenerate spectrum.
2. Kantorovich Method
Let us consider the 2D BVP in the two-dimensional domain Q(xf ,xs) C R2: ( d2 d2 \
(-fcT - W/ + V(x,.x.) - e) »(*/.x.)=0. (1)
where V(xf ,xs) is a real-valued function and ^(xf ,xs) satisfies the Dirichlet condition at the boundary dQ(xf ,xs) of the domain Q(xf ,xs)
^(Xf ,Xs)
= 0. (2)
( Xf ,Xs)edQ(Xf ,xs)
The solution ^(xf,xs)gW2(^) of the BVP (1)-(2) is sought in the form of Kantorovich expansion [7]
Jmax
^i(Xf .xs) = ^2§j(xf; xs)x) (xs). (3)
3 = 1
using the set of eigenfunctions of the parametric BVP
( d2 \
+ V(xf ,xs) - ej (xsU (xf; xs) = 0, (4)
V 9xf
defined in the interval Xf G (^min(^s), ^max(^s)) = QXf (xs) and depending on the variable xs G (^min,^max) = as a parameter. These functions obey the boundary conditions
(xfn(xs); xs) = 0, (xmax(xs); xs ) = 0 (5)
at the boundary points {xmin(x5), xmax(^5)} = dQxf (%s), of the interval Qxf (xs).
The eigenfunctions satisfy the orthonormality condition in the same interval Xf G
^X f (^5):
xm (x s)
= J $i(Xf; xs)^j(Xf; xs)dxf = % (6)
Xmin (x s )
In Eq. (4) €\(xs) < ... < tjmax(xs) < ... is the desired set of real-valued eigenvalues. If this parametric eigenvalue problem has no analytical solution, then it is solved numerically using the ODPEVP program [4] or in the case of two variables POTHEA program [5].
Substituting the expansion (3) into Eq. (1) with Eqs. (5) and (6) taken into account, we arrive at the set of self-adjoint ODEs for the unknown vector functions x(t)(xs,E) =
(
d2 , TT/ > 771 t dQ(xs) d
-I+ U(®a) - EI + Q(Xs)1 x(i) (Xs) = 0. (7)
Here U(œa) and Q(xs) are matrices of the dimension jmax x jmax
Uij (xs ) = ei(xs)ôij + Hij (xs), (8)
it ( \ u < \ f 9^(xf; x°) (xf; x°) J
HtJ (xs) = HJt(xs)= J ---^-dxf, (9)
œmin(œs)
# m )
f -, , (xf;xs) Qij (xs) = -Qji(xs) = - ®i(xf ; zs)-^-ax/.
The solutions of the discrete spectrum E : Ei < E2 < ... < Ev < ... that obey the boundary conditions at the points xls = {^min,^max} = &, bounding the interval QXs,
X(p)(xD = 0, x\ = xfn,xmax (10)
and the orthonormality conditions
£max
J (x(t)(xs))TX(J)(xs)dxs = %, (11)
^min
are calculated by means of the KANTBP program [6].
3. Benchmark Calculation: Triangular Membrane
As a benchmark example we consider the exactly solvable BVP for a triangular membrane in conventional variables (xf ,xs) G fi(®/,xs)
id2 d2 \ - ^ - E) *{x> •I«> = 0 (12)
with the Dirichlet conditions at the boundary dQ(x,y) of the region Q(xf ,xs)
V((xf ,xa) G dn(xf ,xs)) = 0. (13)
In the considered case the parametric eigenvalue problem (4)-(6) has an exact solution, i.e., the parametric eigenfunctions (xf ; xs) and potential curves ei (xs) are expressed in the analytical form
2 ./ Ki(xf-xfln(xB)) \
n2i2 S m\ xrx(xs)-xfn(xs) ei (xs) = --;—--—t-tt, (xf; xs) =-. J . (14)
' v in-max/™ ^ „mi^ \\2 1 1 v f ' 0' 3 " i 3 T v '
(xf (X^) - xf (X°)) J^max(xs) - xYn(xs)
With the basis functions (14) the integration in the effective potentials (9) can be carried out analytically, which yields the expressions
ff 1) i+jdxfax(xB) dx?n(xs) \ \( ) dxs dxs J
Q. .(X) =___dxxi_dx- / i = i
Q (Xs) i2 - j2 Xjax(xs ) -xfn(xs) ' J= '
(, 1)i +j dxfax(xB) dxfin(xB) \ / dxfax(xB) dxfn(xB) \
n ( )_ 4ij(i2 + j2) \( 1) dxB dxB J ^ dxB dxB J
i Xs) _ (i2 - j2)2 (xmax(xs) -Xf[n(Xs))2
(dxfax(xB)\2 / dxfax(xB)\ (dxfin(xB)\ / dxfin(xB)\2 TT , , ■K2f{ dxs I + I dxB I { dxB I + I dxB I
H-ii (xs) _------7-;-r^—^--------+
1s' o i ^max ^ \ „mm / „ W9
(Xfax(Xs) -xfn(xs))2
i dxfax(x3) dxfn(xB) x 2 11 dxs dxs
+ 4 (xmax(x s) - xfn(xs))2
In the symmetric case x—ax(xs) = -x—in(xs) the matrix elements H^ and Qij between even and odd indexes equal zero and one can solve the BVP for even (e) and odd (o) solutions separately.
As a domain we chose the equilateral triangle with side equal to 4i/3, in this case the eigenvalues Ei = + u2 + = 3, 7, 7,12,13,13,19,19,21, 21, 27,..., where = 1,2,..., are integer [9].
Case 1, Xf is paralleled to a triangle side and xs belong to a triangle height:
xfax(xs) = 2ir/3 -xs/V3, x—in(x.s) = -2n/3 + xs/V3, x—in = 0, x—ax = 2i/V3.
Case 2, xs is paralleled to a triangle side and xf belong to a triangle height: x—ax(xs) = 2i/V3 -V3|xs|, x—in(xs) = 0, x—in = -2i/3, x—ax = 2i/3.
In both cases taking into account the symmetry properties of the equilateral triangle, we apply the FEM for discretization of the BVP (7)-(11) using finite element grid QXs = (0(2)3u/4(2)u), v = x—ax — 0.002, where the number of finite elements in each subinterval is presented in parentheses, and Lagrangian interpolation polynomials of p = 12th order, which provides the accuracy 0(hp +1) of the vector-eigenfunctions X(i Kxs,E) = x(i)(xs) = (Xi^ (xs),..., X^L,. (xs))T and 0(h2p) of the eigenvalues Ei, where h = 3 v/8 is the maximal element length [8].
The numerical calculations of eigenvalue problem (7)-(11) were carried out till j—ax = 280 using the new version of the program KANTBP 2.0 implemented in Fortran. In Fig. 1 some typical examples of profiles of the eigenfunctions are presented, corresponding to the exact doubly degenerate eigenvalues E|=Ef=7, E|=E|=13, and E|=Eg=19.
Achieved the discrepancy 5Ef = Ef — Ei of the order of 10-8 for the eigenvalues that is shown in the Table 1. One can see from the table that the convergence rate of the Kantorovich expansion (3) is the order j—^, which corresponds to the theoretical estimations given by the perturbation theory. Similar rate of convergence takes place also in solving of the parametric 2D BVP for a Helium atom [5] and 2D BVP for quadratic membrane [10].
13 ^ 13 19 4 19 M
< »14
Figure 1. Eigenfunctions x,y) of bound states of the 2D boundary-value problem (12), (13) composed by the components Xj )(xs) of the eigenfunctions of the BVP for system of ODEs (7)—(11) and parametric functions (xf; xs) from (14)
Table 1
The discrepancy 6Ef = £f;calc - Ef, a = e, o, vs a number jmax of even (e) and odd (o) basis functions (14) of Kantorovich expansion (3)
Jmax
SEf
SE%
SEI
SEI
SEf0
5E°
SE°
SE°
6 13 28 60 130 280
1.36(-4) 1.37(-5) 1.39(-6) 1.42(-7) 1.41(-8) 1.43(-9)
5.44(-4) 5.41(-5) 5.49(-6) 5.62(-7) 5.56(-8) 5.54(-9)
2.45(-3) 2.35(-4) 2.37(-5) 2.42(-6) 2.39(-7) 2.41(-8)
1.29(-3) 1.22(-4) 1.22(-5) 1.25(-6) 1.23(-7) 1.24(-8)
2.32(-2) 2.02(-3) 1.99(-4) 2.03(-5) 2.01(-6) 2.05(-7)
6.67(-4) 7.44(-5) 7.98(-6) 8.41(-7) 8.42(-8) 8.56(-9)
2.59(-3) 2.85(-4) 3.05(-5) 3.21(-6) 3.22(-7) 3.26(-8)
case 1
2.47(-2)
2.35(-3)
2.44(-4)
2.56(-5)
2.56(-6)
2.59(-7)
6 13 28 60 130 280
8.69(-4) 1.01(-4) 1.13(-5) 1.21(-6) 1.24(-7) 1.29(-8)
8.21(-3) 8.93(-4) 9.82(-5) 1.05(-5) 1.07(-6) 1.09(-7)
1.83(-2) 1.79(-3) 1.93(-4) 2.07(-5) 2.10(-6) 2.13(-7)
2.43(-2) 2.33(-3) 2.52(-4) 2.70(-5) 2.73(-6) 2.78(-7)
0.95
3.57(-2)
3.48(-3)
3.65(-4)
3.68(-5)
3.74(-6)
1.13(-3) 1.34(-4) 1.50(-5) 1.62(-6) 1.64(-7) 1.68(-8)
8.79(-3) 9.97(-4) 1.10(-4) 1.18(-5) 1.20(-6) 1.22(-7)
case 2
4.84(-2)
3.10(-3)
3.27(-4)
3.48(-5)
3.52(-6)
3.58(-7)
exact
=3
EI=7
Ef=12
Ef = 13
£f0=37
E°=7
El=13
E°=37
For the number jmax of the parametric basis functions increased to 280, that requires more RAM and computer time are needed. The dimension of the mass and stiffness matrices and their half-width are following: (12 ■ 4 + 1)jmax x (12 ■ 4 + 1)jmax and (12 ■ 2 + 1)jmax: 294 x 294 and 150 for jmax = 6, 2940 x 2940 and 1500 for jmax = 60, 13720 x 13720 and 7000 for jmax = 280. The calculation time was about 1 seconds for jmax = 6,15 seconds for jmax = 60 and 455 seconds for jmax = 280 in the double precision of Fortran-77 using the PC Intel Core i5 3.33GHz, 4Gb, 64 bit Windows 7.
4. Conclusion
We show and estimate the rate of convergence of Kantorovich expansion (3) in benchmark calculations for the exactly solvable eigenvalue problem of a triangle membrane with the degenerate spectrum, together with the efficiency and stability of the KANTBP program for solving the boundary-value problem for systems of the ODEs of a large dimension.
The proposed benchmark model can be used for testing of algorithms and programs for solving the BVPs for systems of the ODEs or generalized algebraic eigenvalue problems of a large dimension.
References
1. A. A. Gusev, S. I. Vinitsky, O. Chuluunbaatar, V. L. Derbov, A. GoZdZ, P. M. Krassovitskiy, Metastable States of a Composite System Tunneling Through Repulsive Barriers, Theoretical and Mathematical Physics 186 (2016) 21-40.
2. A. A. Gusev, O. Chuluunbaatar, V. P. Gerdt, V. A. Rostovtsev, S. I. Vinitsky, V. L. Derbov, V. V. Serov, Symbolic-Numeric Algorithms for Computer Analysis of Spheroidal Quantum Dot Models, Lecture Notes in Computer Science 6244 (2010) 106-122.
3. S. I. Vinitsky, A. A. Gusev, O. Chuluunbaatar, V. L. Derbov, A. S. Zotkina, On Calculations of Two-Electron Atoms in Spheroidal Coordinates Mapping on Hyper-sphere, Proc. SPIE 9917 (2016) 99172Z.
4. O. Chuluunbaatar, A. A. Gusev, S. I. Vinitsky, A. G. Abrashkevich, ODPEVP: A Program for Computing Eigenvalues and Eigenfunctions and their First Derivatives with Respect to the Parameter of the Parametric Self-Adjoined Sturm-Liouville Problem, Comput. Phys. Commun. 180 (2009) 1358-1375.
5. A. A. Gusev, O. Chuluunbaatar, S. I. Vinitsky, A. G. Abrashkevich, POTHEA: A Program for Computing Eigenvalues and Eigenfunctions and Their First Derivatives with Respect to the Parameter of the Parametric Self-Adjoined 2D Elliptic Partial Differential Equation, Comput. Phys. Commun. 185 (2014) 2636-2654.
6. O. Chuluunbaatar, A. A. Gusev, S. I. Vinitsky, A. G. Abrashkevich, KANTBP 2.0: New Version of a Program for Computing Energy Levels, Reaction Matrix and Radial Wave Functions in the Coupled-Channel Hyperspherical Adiabatic Approach, Comput. Phys. Commun. 179 (2008) 685-693.
7. L. V. Kantorovich, V. I. Krylov, Approximate Methods of Higher Analysis, Wiley, New York, 1964.
8. G. Strang, G. J. Fix, An Analysis of the Finite Element Method, Prentice-Hall, Englewood Cliffs, New York, 1973.
9. F. Pockels, Uber die partielle Differential-Gleichung Au+k2u = 0 und deren auftreten in der mathematischen physik, B. G. Teubner, Leipzig, 1891.
10. A. A. Gusev, L. L. Hai, O. Chuluunbaatar, S. I. Vinitsky, V. L. Derbov, Solution of Boundary-Value Problems using Kantorovich Method, EPJ Web of Conferences 108 (2016) 02026.
УДК 519.632.4
Решение краевых задач для систем ОДУ большой размерности: эталонные расчеты в рамках метода
Канторовича
А. А. Гусев*, О. Чулуунбаатар*, С. И. Виницкий*^,
В. Л. Дербов*
* Объединённый институт ядерных исследований, г. Дубна
^ Российский университет дружбы народов, г. Москва * Саратовский государственный университет, г. Саратов
Представлены эталонные расчеты краевой задачи для систем ОДУ второго порядка большой размерности с помощью программы KANTBP с использованием метода конечных элементов. На практике для решения краевых задач с дальнодействующими потенциалами и
большого числа открытых каналов необходимо решать краевые задачи для систем дифференциальных уравнений большой размерности, которые также требуют изучения сходимости и устойчивости алгоритмов и программ. С этой целью в данной работе решена задача на собственные значения для эллиптического дифференциального уравнения в двумерной области с граничными условиями Дирихле. Решение ищется в виде разложения Канторовича по параметрическим базисным функциям одной из независимых переменных, при этом вторая независимая переменная рассматривается как параметр. Базисные функции вычисляются в аналитическом виде как решения вспомогательной параметрической задачи Штурма-Лиувилля для ОДУ второго порядка. В результате, двумерная задача сводится к краевой задаче для самосопряжённой системы ОДУ второго порядка относительно второй независимой переменной. Дискретизация задачи выполнена в рамках метода конечных элементов. Эффективность, устойчивость и сходимость вычислительной схемы продемонстрирована эталонными расчетами для треугольной мембраны с вырожденным спектром.
Ключевые слова: тестовые расчеты, краевая задача, системы ОДУ большой размерности, метод Канторовича, метод конечных элементов
Литература
1. Metastable States of a Composite System Tunneling Through Repulsive Barriers / A. A. Gusev, S. I. Vinitsky, O. Chuluunbaatar, V. L. Derbov, A. GoZdZ, P. M. Krasso-vitskiy // Theoretical and Mathematical Physics. — 2016. — Vol. 186. — Pp. 21-40.
2. Symbolic-Numeric Algorithms for Computer Analysis of Spheroidal Quantum Dot Models / A. A. Gusev, O. Chuluunbaatar, V. P. Gerdt, V. A. Rostovtsev, S. I. Vinitsky, V. L. Derbov, V. V. Serov // Lecture Notes in Computer Science. — 2010. — Vol. 6244. — Pp. 106-122.
3. On Calculations of Two-Electron Atoms in Spheroidal Coordinates Mapping on Hy-persphere / S. I. Vinitsky, A. A. Gusev, O. Chuluunbaatar, V. L. Derbov, A. S. Zotk-ina // Proc. SPIE. — 2016. — Vol. 9917. — P. 99172Z.
4. ODPEVP: A Program for Computing Eigenvalues and Eigenfunctions and their First Derivatives with Respect to the Parameter of the Parametric Self-Adjoined Sturm-Liouville Problem / O. Chuluunbaatar, A. A. Gusev, S. I. Vinitsky, A. G. Abrashke-vich // Comput. Phys. Commun. — 2009. — Vol. 180. — Pp. 1358-1375.
5. POTHEA: A Program for Computing Eigenvalues and Eigenfunctions and Their First Derivatives with Respect to the Parameter of the Parametric Self-Adjoined 2D Elliptic Partial Differential Equation / A. A. Gusev, O. Chuluunbaatar, S. I. Vinitsky, A. G. Abrashkevich // Comput. Phys. Commun. — 2014. — Vol. 185. — Pp. 26362654.
6. KANTBP 2.0: New Version of a Program for Computing Energy Levels, Reaction Matrix and Radial Wave Functions in the Coupled-Channel Hyperspherical Adiabatic Approach / O. Chuluunbaatar, A. A. Gusev, S. I. Vinitsky, A. G. Abrashkevich // Comput. Phys. Commun. — 2008. — Vol. 179. — Pp. 685-693.
7. Kantorovich L. V., Krylov V. I. Approximate Methods of Higher Analysis. — New York: Wiley, 1964.
8. Strang G., Fix G. J. An Analysis of the Finite Element Method. — New York: Prentice-Hall, Englewood Cliffs, 1973.
9. Pockels F. Uber die partielle Differential-Gleichung Au+k2и = 0 und deren auftreten in der mathematischen physik. — Leipzig: B. G. Teubner, 1891.
10. Solution of Boundary-Value Problems using Kantorovich Method / A. A. Gusev, L. L. Hai, O. Chuluunbaatar, S. I. Vinitsky, V. L. Derbov // EPJ Web of Conferences. — 2016. — Vol. 108. — P. 02026.
© Gusev A.A., ChuluunbaatarO., Vinitsky S.I., Derbov V.L., 2016
