Algorithms for Solving the Parametric Self-Adjoint 2D Elliptic Boundary-Value Problem Using High-Accuracy Finite Element Method

Full Text

Introduction The adiabatic representation is widely applied for solving multichannel scattering and bound-state problems for systems of several quantum particles in molecular, atomic and nuclear physics [1-4]. Such problems are described by elliptic boundary value problems (BVPs) in a multidimensional domain of the configuration space, solved using the Kantorovich method, i.e., the reduction to a system of self-adjoint ordinary differential equations(SODEs) using the basis of surface functions of an auxiliary BVP depending on the independent variable of the SODEs parametrically. The elements of matrices of variable coefficients of these SODEs including the matrix of the first derivatives are determined by the integrals of products of surface eigenfunctions and/or their first derivatives with respect to the parameter. Thus, the key problem of such method is to develop effective algorithms and programs for calculating with given accuracy the surface eigenfunctions and the corresponding eigenvalues of the auxiliary BVP, together with their derivatives with respect to the parameter, and the corresponding integrals that present the matrix elements of the effective potentials in the SODEs [5-7]. In this paper we propose new calculation schemes and algorithms for the solution of the parametric 2D elliptic boundary-value problem using high-accuracy finite element method (FEM) with rectangular and triangular elements. The algorithms were implemented in a package of programs that calculate with the given accuracy eigenvalues, eigenfunctions and their first derivatives with respect to the parameter of the parametric self-adjoint elliptic differential equations with the boundary conditions of the Dirichlet and/or Neumann type in the 2D finite domain and the integrals of products of the surface eigenfunctions and their first derivatives with respect to the parameter that express the matrix elements of the effective potentials in the SODEs. We also propose a method of constructing the etalon potential in the auxiliary parametric BVP that allows the calculation of the parametric surface basis functions in the analytical form. These functions can be then used for the reduction of the original 2D BVP to the SODEs containing the additional potential matrix, representing the discrepancy between the original potential and the etalon one, averaged with the basis functions. The efficiency of the calculation schemes and algorithms is demonstrated by benchmark calculation of the 2D BVPs describing quadrupole vibrations in the collective nuclear model [4, 8, 9]. The structure of the paper is the following. In Section 2 the Kantorovich method for solving the 2D and 3D BVPs is considered. In Sections 3 and 4 the 2D FEM schemes and algorithms for solving the parametric 2D BVP and calculating derivatives with respect to the parameter together with the corresponding matrix elements are presented. In Section 5 the benchmark calculations of 2D FEM algorithms and programs are analyzed. In Conclusion we discuss the results and perspectives. Kantorovich method with the etalon potential Let us consider the BVP in the domain Ω(x , xs) ⊂ Rn-1 × R1: where D(x ; xs) is a self-adjoint elliptic differential operator in the finite region Ω(x ; xs) Rn-1, E is the spectral parameter, corresponding to the energy of the quantum system, si(xs) > 0, ∂s si(xs) and (x , xs) ∂s (x , xs) are real-valued continuous bounded functions in Ω(x , xs), and Ψ(x , xs) satisfies the Dirichlet and/or Neumann boundary condition (BC) at the boundary ∂Ω ∂Ω(x , xs) of the domain Ω(x , xs) and the orthonormalization conditions The solution Ψ(x , xs)∈W 2(Ω) of the BVP (1) is sought in the form of the Kantorovich expansion [3] using the set of parametric eigenfunctions Φ (x ; xs) ∈ (xs) ∼ W 2(Ω(x ; xs)) of the parametric BVP in the domain Ω(x ; xs) ⊂ Rn-1 For example, D(x ; xs) at n 1 = 2 is determined in a conventional form in Section • To avoid cumbersome notations, below we will use the simplest definition of D(x ; xs) (xmin(xs), xmax(xs))=Ω (xs) and depending on the variable xs ∈ Ω as a parameter. We assume that these functions obey the BCs at the boundary points {xmin(xs), xmax(xs)}=∂Ω (xs), of the interval Ω (xs). The eigenfunctions satisfy the orthonormality condition Here (xs) : 1(xs) < . . . < max (xs) < . . . is the desired set of real eigenvalues. During the simulation the etalon potential 0(x ; xs, g(xs)) in Eq. (5) will be chosen as 0(x ; xs, g(xs)) = (x , xs) or calculated separately for different values of xs Ωs from the conditions If this parametric eigenvalue problem has no analytical solution, then it is solved numerically by the FEM: at n 1 = 1 using the program ODPEVP [6] and at n 1 = 2 using the program, implementing the algorithm presented in Section 3. Substituting the expansion (3) into Eq. (1) with Eqs. (6) and (7) taken into account, we arrive at the set of self-adjoint ODEs for the unknown vector functions (i)(xs, E) ≡ Here I, W(xs) and Q(xs) are the max × max matrices The effective potentials i (xs)=i(xs) are calculated by integrating the difference d (x , xs) = (x , xs)-(x ; xs, g(xs)) with the basis functions: If the etalon potential is equal to the original one, (x ; xs, g(xs))= (x , xs) then i (xs)=0. If the original potential (x , xs) is presented in the tabular form, e.g., as in Ref. [4], then the etalon potential (x ; xs, g(xs)) can be considered as a certain approximation or interpolation of the original one, and the approximation error d (x , xs) = 0. In this case the etalon potential (x ; xs, g(xs)) can be determined minimizing d (x , xs) in accordance with the condition (8) of the least squares method. Then the integrals i (xs) of the approximation error d (x , xs) averaged with the parametric basis functions are included into the final system of ODEs. It means that we take the approximation error of the original potential into account and require the solutions of the 2D BVP (1) to have the given accuracy, see, e.g., Section 5.4. The solutions of the discrete spectrum E : E1 < E2 < . . . < E < . . . obey the BCs at the points xt ={xmin, xmax}=∂Ω , bounding the interval Ω and satisfy the orthonormality conditions 3. FEM algorithm for solving the parametric 2D BVP Let us consider a boundary value problem for the parametric self-adjoint 2D PDE in the domain Ω,y = (xmin, xmax) × (ymin, ymax) with the Dirichlet and/or Neumann boundary conditions where t = min, max. Here z Ω = [zmin, zmax] is a parameter, the functions 1(x) > 0, are continuous and bounded for (x, y) ∈ Ω,y . Also assume that the BVP (15), (16) has only the discrete spectrum, so that (xs) : 1(xs) < . . . < max (xs) < . . . is the desired set of real eigenvalues. The eigenfunctions satisfy the orthonormality conditions The FEM calculation scheme is derived from the variational functional In each • The domain ∆ = [xmin, xmax] × [ymin, ymax] is covered by⋃︀ the⋃︀system of n × subdomain ∆ subdomains ∆i = [xi-1, xi] × [y-1, y ] in such a way that ∆ = and the Lagrange elements {︀ are determined. By means of the Lagrange elements (x) and (y), we define the set of piecewise polynomial functions N (x) and M (y) as follows: for any and The functions {N (x)} and {M (y)} form a basis in the space of polynomials of the -th order. Now, the function Φ(x, y; z) ∈ h ∼ 1(Ωh ,h) is approximated by a finite sum of piecewise polynomial functions N (x) and M (y) The domain Ω(x, y) = ⋃︀Q ∆ , specified as a polygon in the plane (x ≡ z1, y ≡ z2) ∈ 2, is covered with finite elements, the triangles ∆ with the vertices (z11, z21), (z12, z22), (z13, z23) (here zik zik; , i = 1, 2, k = 1, 3, q = 1, Q). On each of the triangles ∆ (the boundary is considered to belong to the triangle) the shape functions (z1, z2) are introduced. For this purpose we divide the sides of the triangle into equal parts and draw three families of parallel straight lines through the partition points. The straight lines of each family are numbered from 0 to , so that the line passing through the side of the triangle has the number 0, and the line passing through the opposite vertex of the triangle has the number . Three straight lines from different families intersect in one point Al ∆ , which will be numbered by the triplet (n1, n2, n3), ni 0, n1 + n2 + n3 = , where n1, n2 and n3 are the numbers of the straight lines passing parallel to the side of the triangle that does not contain the vertex (z11, z21), (z12, z22) and (z13, z23), respectively. The coordinates of this point zl = (z1l , z2l ) are determined by the expression (z1l , z2l ) = (z11, z21)n1/ + (z12, z22)n2/ + (z13, z23)n3/. As shape functions we use the Lagrange triangular polynomials (z1, z2) of the order that satisfy the condition l (z1l′ , z2l′ ) = dll′ , i.e., equal 1 in one of the points Al and zero in the other points. In this method the piecewise polynomial functions N (z1, z2) in the domain Ω are constructed by joining the shape functions (z , z ) in triangle ∆ : and possess the following properties: functions N (z1, z2) are continuous in the domain Ω; the functions N (z , z ) equal 1 in one of the points A and zero in the rest points; Nl (z1l′ , z2l′ ) = dll′ in the entire domain Ω. Here l takes the values l = 1, N . The functions N (x, y) form a basis in the space of polynomials of the -th order. Now, the function Φ(x, y; z) ∈ h ∼ 1(Ωh ,h ) is approximated by a finite sum of piecewise basis functions N (x, y) The vector function has a generalized first-order partial derivative and belongs to the Sobolev space 1(Ωh ,h ) [10]. After substituting the expansion (19) or (20) into the variational functional (18) and minimizing it [10, 11], we obtain the generalized eigenvalue problem (21) Here A is the stiffness matrix; B is the positive definite mass matrix; h is the vector approximating the solution on the finite-element grid; and h is the corresponding eigenvalue. The matrices A and B have the following form: where the local matrices a and b are calculated for rectangular elements as or the matrix elements and b are calculated for triangular elements as where ∆ is a triangular element with the vertices (z11, z21)(z12, z22)(z13, z23), ∆ is the triangle with vertices (0, 0)(1, 0)(0, 1), and is the Jacobian of the transformation from the global frame (x ≡ z1, y ≡ z2) to the local one (z1′ , z2′ ); dz1dz2 = dz1′ dz2′ and (z1, z2) are expressed via (z1′ , z2′ ) by the relations: In this case we have explicit expression for shape functions (z′ , z′ ): Remark. We start from the initial (global) coordinate frame (x z1, y z2) in which the coordinates of vertices of the triangle ∆ are equal to (z11, z21), (z12, z22), (z13, z23), and the local coordinate frame (x ≡ z1′ , y ≡′z2′ ) in which the coordinates of the same vertices of the triangle ∆ are equal to (z11, z2′ 1) = (0, 0), (z1′ 2, z2′ 2) = (1, 0), (z1′ 3, z2′ 3) = (0, 1). We seek the relation between the global and the local coordinates of the triangle vertices in the form Substituting the coordinates of three vertices in the global and the local frame into Eq.(26), we obtain the system of algebraic equations for calculating the coefficients **: Substituting the calculated coefficients into Eq.(26), we arrive at the formula (25) for the transformation of coordinates (z1, z2) → (z1′ , z2′ ), from which we express the inverse transformation of coordinates (z1′ , z2′ ) → (z1, z2) where -1 is the Jacobian of the inverse transformation from the local frame (z1′ , z2′ ) to the global one (x z1, y z2): dz1′ dz2′ = -1dz1dz2. The integrals (22) and (23)-(24) are evaluated using the Gaussian quadrature of the order + 1. Note that in the present approach the maximum half-band width of the matrices A and B is small compared to their dimension and is not greater than In order to solve the generalized eigenvalue problem (21), the subspace iteration method [10, 11] elaborated by Bathe [11] for the solution of large symmetric bandedmatrix eigenvalue problems has been chosen. This method uses the skyline storage mode which stores the components of the matrix column vectors within the banded region of the matrix, and is ideally suited for banded finite-element matrices. The procedure chooses a vector subspace of the full solution space and iterates upon the successive solutions in the subspace (for details, see [11]). The iterations continue until the desired set of solutions in the iteration subspace converges to within the specified tolerance on the Rayleigh quotients for the eigenpairs. If the matrix A in Eq. (21) is not positively defined, the problem (21) is replaced with the following problem: The number a (the shift of the energy spectrum) is chosen such that the matrix A˜ is positive defined. The eigenvector of the problem (27) is the same, and h = ˜h + a. The algorithm for calculating the parametric derivatives of eigenfunctions and the matrices of effective potentials Taking a derivative of the boundary problem (15)-(17) with respect to the parameter z, we find that ∂ Φi(x, y; z) is a solution of the following boundary problem where t = min, max. The parametric BVP (28), (29) has a unique solution, if and only if it satisfies the conditions Below we present an efficient numerical method that allows the calculation of ∂ Φi(x, y; z) with the same accuracy as achieved for the eigenfunctions of the BVP (15)-(17) and the use of it for computing the matrices of the effective potentials defined as The boundary problem (28)-(30) is reduced to the linear system of inhomogeneous algebraic equations with respect to the unknown ∂h/∂z: The normalization condition (17), the condition of orthogonality between the function and its parametric derivative (30), and the additional conditions (29) for the solution of (32) read as Then the potential matrix elements h (z) and Qh (z) (31) corresponding to Eqs. (11)-(12) can be calculated using the formulas Since h is an eigenvalue of (21), the matrix L in Eq. (32) is degenerate. In this case the algorithm for solving Eq. (32) can be written in three steps as follows: Step k1. Calculate the solutions v and w of the auxiliary inhomogeneous systems of algebraic equations (35) with the non-degenerate matrix L¯ and the right-hand sides b¯ and d where is the number of the element of the vector Bh having the greatest absolute value. Step k2. Evaluate the coefficient g (36) Step k3. Evaluate the vector ∂ h (37) From the above consideration it is evident, that the computed derivative has the same accuracy as the calculated eigenfunction. Let D(x, y; z) in Eq. (15) be a continuous and bounded positively defined operator on the space 1 with the energy norm, i(z), Φi(x, y; z) ∈ 2 being the exact solutions of Eqs. (15)-(17), and h(z), Φh(x, y; z) ∈ 1 being the corresponding numerical solutions. Then the following estimates are valid [10, 12] where h is the largest distance between any two points in ∆ (see [12], p. 161), is the order of the finite elements, i is the number of the corresponding solutions, and the constants 1 and 2 are independent of the step h. The following theorem can be formulated. Theorem 1. Let D(x, y; z) in Eq. (15) be a continuous and bounded positively defined operator on the space 1 with the energy norm. Also let ∂ (x, y; z) be continuous and bounded for each value of the parameter z. Then for the exact values of the solutions ∂ i(z), ∂ Φi(x, y; z) ∈ 2, i (z), Qi (z) from (28)-(31) and the corresponding numerical values ∂ h(z), ∂ Φh(x, y; z) ∈ 1, h (z), Qh (z) from (32)-(34) the following estimates are valid: where h is the largest distance between any two points of the finite element ∆, is the order of finite elements, i, are the numbers of the corresponding solutions, and the constants 3, 4, 5 and 6 are independent of the step h. The proof is straightforward following the proof schemes in accordance with [10, 13]. Benchmark calculations The solution of 2D BVP in the triangular domain As a benchmark example we consider the exactly solvable BVP for a membrane in the form of equilateral triangle with side equal to 4/3, in the conventional variables with the Dirichlet or Neumann conditions at the boundary ∂Ω(x, y) of the region Ω(x, y). In both cases the eigenvalues i are integer [14-16]. The BVPs were solved on the uniform finite-element grid composed of n2 equilateral triangles with the side equal to 4/(3n). In Table 1 we present the estimations of the finite-element scheme order depending on the size of elements. Table 1 Discrepancies di (n) = i (n) i of the first 1 = 3 and the fourth 4 = 3 eigenvalues of the BVP(40) with Dirichlet and Neumann boundary conditions, respectively, and the Runge coefficients Rui = log2((di (n) di (2n))/(di (2n) di (4n))) for the schemes of the first = 1 and the second = 2 orders; n is the number of elements n d1(n) d1(2n) d1(4n) Ru1 1 6 0.06914677126132 0.01717343333794 0.00428595766335 2.011 2 3 0.00470310603030 0.00030790240009 0.00001932528605 3.928 n d4(n) d4(2n) d4(4n) Ru4 1 6 0.06914677126135 0.01717343333808 0.00428595766371 2.011 2 3 0.00470310603029 0.00030790240005 0.00001932528593 3.928 They are seen to correspond rigorously to the theoretical estimations (38) of the order (h2) of the eigenvalues approximation. The first eight eigenfunctions of the BVPs (40) with the Neumann and Dirichlet boundary conditions are shown in Figs. 1 and 2. Figure 1. The first eight eigenfunctions of the BVPs (40) in the triangle region with the Neumann boundary conditions. The eigenvalues i = 0, 1, 1, 3, 4, 4, 7, 7 at i = 1, ..., 8. Figure 2. The first eight eigenfunctions of the BVP (40) in the triangle region with the Dirichlet boundary conditions. The eigenvalues i = 3, 7, 7, 12, 13, 13, 19, 19 at i = 1, ..., 8. Solution of the parametric 2D BVP for the oscillator potential We consider the parametric 2D BVP (15)-(17) at 1(x) = 2(x, y) = 3(x) = 4(y) = 5(x, y) = 1, with the potential function (41) which has the known spectrum ={1,2} = (2(1 - 1) + 1) + 3/2(2(2 - 1) + 1) and the eigenfunctions Φ1,2 (z1, z2; z3) = Φ1-1(z1; z3)Φ2-1(z2; 0), where Φ1 (z1; z3) and Φ2-1(z2; 0) are determined in (45) at w1(xs) = 1, z01(xs = z3) = z3 and w2(xs) = 3/2, z02(xs) = 0, respectively. The matrix elements Qi={i1,i2},={1,2} and i={i1,i2},={1,2} are calculated in the analytical form where n= max(i1, 1)-1. Solving the BVP, we use 294 finite elements, i.e. equilateral triangles with the side equal to 1 that form a regular hexagon with the vertices {(7 cos n/3, 7 sin n/3)}5. In each element the interpolation polynomials of the six order were applied. Stiffness and mass matrices with dimensions of 5167 5167 were used. The calculated eigenvalues and the exact ones i are: i = (2.5, 4.5, 5.5, 6.5, 7.5) . The calculated matrices Qh and h , i, = 1, ..., 5 from Eq. (34) of the parametric 2D BVP with the potential function (41) are From the comparison of the calculated and the exact values of i and Qi , i from (42), one can see that the achieved accuracy of results is of the order of 8-9 significant digits. The parametric surface eigenfunctions and their derivatives with respect to the parameter under the dipole splitting are shown in Fig. 3. Figure 3. The first four eigenfunctions of the parametric 2D BVP with the potential function (41) and their derivatives with respect to the parameter at = 0. The solution of 2D BVP for the C3 oscillator potential Let us consider the 2D BVP (15)-(17) for x=22, y=20 with the potential function The quadrupole potential energy is approximated by the quartic potential: We use the set of parameters 1= 120, 2=240, 3=1200, 0=65/16 that provide a crude approximation for the shape of 156Gd92, which has been fitted in the following points1: the minima at (22, 20)=(0, 1/4), (0, 1/4)=0; the maxima at (22, 20)=(0, 0), (0, 0)=65/16; and the saddle point at (22, 20)=(0, 1/5), (0,1/5)=729/400 (see Fig. 4). We choose the mass parameter to be m=2=124. Thus, there are ground and doubly degenerate excited states, localized in three wells. The BVP was solved using three algorithms (and the forth in Section 5.4): The solution of the BVP was calculated using the FEM scheme from Section 3 on the rectangular grid [ 0.4, 0.3, ..., 0.4] [ 0.4, 0.3, ..., 0.4] with Lagrange interpolation polynomials of the order = 12. The first 18 eigenvalues were calculated with 9 significant digits and are presented in Table 2. 1The fitting points in our parametrization are related to those of Ref. [4] as 22 = √2a22, 20 = a20. Figure 4. The potential energy of the C3 oscillator having the quadrupole shape, the eigenfunction Ψ1(22, 20) of the ground state (irr. A1) and the degenerate eigenfunctions Ψ2(22, 20) (irr. E1) and Ψ3(22, 20) (irr. E2) Table 2 The first energy levels of the 2D BVP using triangular △ and rectangular D finite elements and the Kantorovich method h using max=28 parametric basis functions, classified by irrs (Class.) of the C3 point group, Ei (MeV) i △i D i h i Eh=h/2m i i Class. 1 381.754344 381.754355 381.754351 1.53933206 A1 2 3 387.240633 387.240633 387.240644 387.240646 387.240641 387.240641 1.56145419 1.56145419 E1 E2 4 5 617.024951 617.024951 617.024967 617.024989 617.024963 617.024963 2.48800388 2.48800388 E1 E2 6 667.104970 667.105020 667.104992 2.68993948 A2 7 695.166557 695.166590 695.166575 2.80309103 A1 8 9 785.680037 785.680037 785.680100 785.680136 785.680078 785.680078 3.16806483 3.16806483 E1 E2 10 898.045395 898.045497 898.045434 3.62115094 A1 11 12 915.823095 915.823095 915.823200 915.823309 915.823167 915.823167 3.69283535 3.69283535 E1 E2 13 14 993.158636 993.158636 993.158784 993.158872 993.158708 993.158708 4.00467221 4.00467221 E2 E1 15 1063.73690 1063.73709 1063.73692 4.28926178 A1 16 1119.21670 1119.21668 1119.21649 4.51296973 A2 17 18 1174.72840 1174.75531 1174.71177 1174.71183 1174.71166 1174.71166 4.73674057 4.73674057 E1 E2 • In the solution of the BVP we used 54 finite elements in the form of equilateral triangles with the side equal to 1/6, forming a regular hexagon with the vertices {(0.5 cos n/3, 0.5 sin n/3)}5 . In each finite element the interpolation polynomials of the fifth order were applied. The stiffness and the mass 721 721 matrices were used. The calculated eigenvalues are presented in Table 2. • The problem (1)-(43) was also solved using the Kantorovich method implemented in the FEM program KANTBP 2 [5], with max=28 parametric basis functions calculated using the FEM program ODPEVP [6]. The solution was calculated in the domain 2 +2 <1/2 with Dirichlet BCs at the boundary 2 +2 =1/2 using the scheme presented in Section 2. The calculations were performed in the case of (x , xs)= (22, 20) as well as (x , xs)=(20, 22) on the FE grid {-1/2, -1/3, -1/6, 0, 1/6, 1/3, 1/2} with Lagrange interpolation polynomials of the order = 12. The point symmetry group C3 of the problem (1), (43) has four irreducible representations (irrs.) A1, A2, E1 and E2 to classify the solutions, the E-type states being doubly degenerate [17, 18]. The calculated eigenvalues are presented in Table 2. The first eigenfunctions for each of the irrs. A1, E1, E2 are shown in Fig. 4. Parametric surface functions for KM in the analytical form Let us consider the BVP for Eq. (5) with the etalon potential (x ; xs, g(xs)): where g(xs) is the set of parameters, g(xs) = 0(xs), w2(xs), z0(xs) . In the considered case the parametric eigenvalue problem (5)-(7) has an exact solution, i.e., the parametric eigenfunctions Φi (x ; xs) and potential curves i (xs) are expressed in the analytical form The integration in the effective potentials (11) with the basis functions (45) is carried out analytically, which yields the expressions where n= max(i, )-1. The effective potentials are calculated by integrating the difference (x , xs)-(x ; xs) with the basis functions: ∫︁ During the simulation the adiabatic parameters 0(xs), w(xs), z0(xs) of the etalon potential (44) were found from the conditions For the potential (43) from the condition (46) we have We performed the calculations for these parameters in the case when (x , xs)=(22, 20), as well as (x , xs)=(20, 22). The results coincide with those of the calculations of h performed in the previous Section 5.3 and presented in Table 2 with 9 significant digits. 6. Conclusion We elaborated new calculation schemes and algorithms for solving the parametric 2D elliptic BVP using the high-accuracy FEM with rectangular (22) and triangular (23)-(24) elements. The algorithm and the programs calculate with the given accuracy the eigenvalues, the eigenfunctions and their first derivatives with respect to the parameter of the parametric self-adjoint elliptic differential equations with boundary conditions of the Dirichlet and/or Neumann type in the finite 2D domain (15)-(17), (18) and the corresponding inhomogeneous boundary-value problem (28)-(33), obtained by taking a parametric derivative of the original 2D BVP. The program also calculates the potential matrix elements, the integrals of the eigenfunctions multiplied by their first derivatives with respect to the parameter (34). The parametric eigenvalues (potential curves) and the matrix elements computed by the program can be used for solving the bound-state and multi-channel scattering problems for a system of the coupled second-order ODES with using the Kantorovich method. We demonstrated the efficiency of the proposed calculation schemes, algorithms and codes by the example of solving the boundary-value problem of a quadruple vibration collective nuclear model. We proposed the construction of parametric surface functions in the analytical form as eigenfunctions of the etalon equation (44), which provides the solution of the 2D BPV with given accuracy and reduces the expenditures of computer resources compared to the conventional basis, numerically calculated using FEM. One can construct parametric functions using different types of etalon potentials, e.g., that of the two-center problem with harmonic oscillator potentials [19]. This approach can be generalized for the BVP in a multidimensional domain using, e.g., the multistep Kantorovich method [3]. The proposed algorithms and codes can be adapted and applied to the analysis of quantum transparency effect, to the study of resonance three-body scattering problems, the quantum diffusion of molecules, the penetration of micro-clusters through surfaces, the fragmentation mechanism in producing very neutron-rich light nuclei, and heavy-ion collisions, as well as the microscopic study of tetrahedral-symmetric nuclei.

About the authors

A A Gusev

Joint Institute for Nuclear Research

Email: gooseff@jinr.ru
6, Joliot-Curie, Dubna, Moscow region, Russia, 141980

O Chuluunbaatar

Joint Institute for Nuclear Research

Email: chuka@jinr.ru
Institute of Mathematics, National University of Mongolia, Ulaanbaatar, Mongolia 6, Joliot-Curie, Dubna, Moscow region, Russia, 141980

S I Vinitsky

Joint Institute for Nuclear Research

Email: vinitsky@theor.jinr.ru
RUDN University (Peoples’ Friendship University of Russia) 6, Miklukho-Maklaya str., Moscow, Russia, 117198 6, Joliot-Curie, Dubna, Moscow region, Russia, 141980

V L Derbov

Saratov State University

Email: derbov@sgu.ru
Saratov, Russia

A Góźdź

Institute of Physics, University of M. Curie-Skl-odowska

Email: andrzej.gozdz@umcs.pl
Lublin, Poland

References