CC BY
11
Discrete & Continuous Models
& Applied Computational Science__2024, 32 (1) 74-85
ISSN 2658-7149 (Online), 2658-4670 (Print) http://journals.rudn.ru/miph
Research article
UDC 519.6:004.94
PACS 07.05.Tp, 02.60.Pn, 02.70.Bf
DOI: 10.22363/2658-4670-2024-32-1-74-85 EDN: BUEBFE
Application of the Chebyshev collocation method to solve boundary value problems of heat conduction
Konstantin P. Lovetskiy1, Stepan V. Sergeev1, Dmitry S. Kulyabov1,2, Leonid A. Sevastianov1,2
1RUDN University, 6 Miklukho-Maklaya St, Moscow, 117198, Russian Federation 2 Joint Institute for Nuclear Research, 6 Joliot-Curie St, Dubna, 141980, Russian Federation
(received: December 3, 2023;revised: January 12, 2024;accepted: January 25, 2024)
Abstract. For one-dimensional inhomogeneous (with respect to the spatial variable) linear parabolic equations, a combined approach is used, dividing the original problem into two subproblems. The first of them is an inhomogeneous one-dimensional Poisson problem with Dirichlet-Robin boundary conditions, the search for a solution of which is based on the Chebyshev collocation method. The method was developed based on previously published algorithms for solving ordinary differential equations, in which the solution is sought in the form of an expansion in Chebyshev polynomials of the 1st kind on Gauss-Lobatto grids, which allows the use of discrete orthogonality of polynomials. This approach turns out to be very economical and stable compared to traditional methods, which often lead to the solution of poorly defined systems of linear algebraic equations. In the described approach, the successful use of integration matrices allows complete elimination of the need to deal with ill-conditioned matrices.
The second, homogeneous problem of thermal conductivity is solved by the method of separation of variables. In this case, finding the expansion coefficients of the desired solution in the complete set of solutions to the corresponding Sturm-Liouville problem is reduced to calculating integrals of known functions. A simple technique for constructing Chebyshev interpolants of integrands allows to calculate the integrals by summing interpolation coefficients.
Key words and phrases: initial boundary problems, pseudo spectral collocation method, Chebyshev polynomials, Gauss-Lobatto sets, numerical stability, separation of variables
1. Introduction
Many important physics problems that involve two or more independent variables are solved using mathematical models that include partial differential equations, not limited to models based on ordinary differential equations. Models of this kind also include the description of heat propagation in solids using the heat equation with various boundary and initial conditions. An important method for solving partial differential equations known as the separation of variables method will be discussed below. Its essential feature is the reduction of the original partial differential equation to a system of simpler ordinary differential equations, which can be successfully solved based on given initial or boundary conditions.
The desired solution to a partial differential equation is expressed as an infinite series, which is the sum of solutions of individual ordinary differential equations. In many cases, it is convenient to represent the required solutions in the form of a series of sines, cosines, or polynomial functions, for example Chebyshev polynomials. This approach allows an effective use of the collocation method -a projection method for solving both integral and differential equations. Therefore, the first part of the present paper is devoted to a discussion of the Chebyshev interpolation method for solving the one-dimensional heat equation.
The Chebyshev collocation method has proven itself in solving a wide class of problems [1-5]. Particularly, in [6-9] its effectiveness was demonstrated in solving ODEs and problems of restoring functions from known first- and second-order derivatives. In Ref. [10] stable spectral methods for
©Lovetskiy K.P., Sergeev S.V., Kulyabov D.S., Sevastianov L.A., 2024
This work is licensed under a Creative Commons Attribution 4.0 International License https://creativecommons.Org/licenses/by-nc/4.0/legalcode
solving the Poisson equation with Dirichlet-Dirichlet, Dirichlet-Neumann and Neumann-Neumann boundary conditions were analyzed in detail. The present paper describes an algorithm for solving the 1-D Poisson equation with Dirichlet-Robin boundary conditions. The algorithm is based on the developed new method of spectral collocation and illustrates the effectiveness of the variable separation method in solving inhomogeneous heat conduction problems.
Spectral methods have proven themselves to be excellent in solving homogeneous boundary value problems for a wide class of partial differential equations using the method of separation of variables. In cases of inhomogeneous problems, methods for separating variables are not directly applicable. However, in this paper we show how the Chebyshev collocation method can be effectively applied in a two-stage solution scheme for a certain class of inhomogeneous boundary value problems for a 1-D linear parabolic equation.
2. Mathematical model of heat conduction
Thermal conductivity is the property of a material to conduct heat, which is assessed primarily from the point of view of Fourier's law of thermal conductivity. Heat conduction, also called diffusion, is the direct microscopic exchange of kinetic energy of particles (such as molecules) or quasiparticles (such as lattice waves) across a boundary between two systems. On a microscopic scale, thermal conduction occurs when hot, fast-moving, or vibrating atoms and molecules interact with neighboring atoms and molecules, transferring some of their energy (heat) to those neighboring particles. In other words, heat is transferred by conduction when neighboring atoms vibrate relative to each other or when electrons move from one atom to another.
When an object has a different temperature than another body or its surroundings, heat flows so that the body and surroundings reach the same temperature, at which point they are in thermal equilibrium. This spontaneous transfer of heat always occurs from a region of high temperature to another region of lower temperature, as postulated by the second law of thermodynamics. Thermodynamic and mechanical heat transfer are calculated using the heat transfer coefficient - the proportionality between heat flow and the thermodynamic driving force of heat flux. Heat flux is a quantitative vector representation of the movement of heat through a surface [11]. In an engineering context, the term "heat" is perceived as synonymous with thermal energy.
The heat conduction equation models diffusion processes [12], including thermal energy in solids, solutes in liquids, and biological populations. We will consider the heat conduction equation describing the temperature change in a one-dimensional rod of a finite length. Let us also consider several possible types of boundary conditions that can be used when modeling temperature changes.
A commonly used method for solving the heat conduction equation is the complete separation of variables method, which results in the solution of two ordinary differential equations generated by the separation of variables method. To solve one of the emerging subproblems, a uniform approach to solving the heat conduction equation for almost any of the frequently used (Dirichlet-Neumann-Robin) sets of boundary conditions is considered. A technique is proposed for constructing a general solution to the inhomogeneous Poisson equation - the heat conduction equation - regardless of the type of boundary conditions. Concretization of the solution - the determination of a pair of missing coefficients of expansion of the solution according to the selected polynomial basis occurs at the second stage, considering the specified (distinct types and combinations) boundary conditions.
3. Inhomogeneous boundary value problems
Let us consider the solution of an inhomogeneous initial-boundary value problem for a one-dimensional parabolic equation, including a time-independent inhomogeneous part of the equation and time-independent boundary conditions. Two-sided Dirichlet-Dirichlet conditions
32u du
k--(x,t)---(x,t) = -F(x), 0<x<L, t>0, ox2 at
u(0, t) = u0, u(L, t) = u1, t> 0, (1)
u(x, 0) = f(x), 0 < x <L.
The boundary conditions of problem (1) mean that the left and right edges of the rod have different constant temperatures due to ideal contact with heat baths having respectively temperatures u0 and m1.
Recall that a boundary value problem (BVP) is called inhomogeneous if either the partial differential equation or the boundary conditions are inhomogeneous. The well-known method of separating variables is not applicable to such inhomogeneous boundary value problems. However, in some cases it is possible to change the variables in such a way that the inhomogeneous boundary value problem transforms into two problems. One of which is a relatively simple inhomogeneous BVP for an ordinary differential equation (ODE), and the other is a homogeneous BVP for a partial differential equation (PDE).
Assume a function F(x) to describe the intensity of heat generation inside the rod, u0 and u1 being constant. Replacing the unknown u(x, t) with a new variable by means of substitution u(x, t) = v(x, t) + $(x), we reduce the solution of problem (1) to a subsequent solution of two subproblems:
Problem A, an inhomogeneous problem with two-sided Dirichlet-Dirichlet conditions.
kip"(x) + F(x) = 0, jp(0) = u0, ip(L) = u1. (2)
Problem B, a homogeneous boundary value problem for a partial differential equation.
(d2vr . dv. ,
v(0, t) = 0, v(L, t) = 0, (3)
\,v(x,0) = f(x)-ij)(x).
It is important to note that problem A is a simple one-dimensional Poissonproblem — an ordinary differential equation of the second order with given inhomogeneous boundary value Dirichlet conditions. At the same time, this problem is a problem of restoring a function from its known second-order derivative F(x) with two additional conditions, in this case, the Dirichlet boundary conditions.
The second auxiliary subproblem B is a homogeneous BVP that can be solved based on the traditional separation of variables method. The solution to the original problem (1) will be the "sum" of the solutions to problems A and B.
Let us illustrate the considered approach to the analytical solution of the inhomogeneous boundary value problem (1) using a pair of examples.
Example 1. Assume that F(x) = r = const > 0. It is required to solve the problem (1) under the following boundary and initial conditions for L = 1:
u(0,t) = 0, u(1,t) = u0, t > 0, u(x, 0) = f(x), 0 < x <1.
In our example, both the partial differential equation and the boundary condition at point x= 1 are inhomogeneous. Let us perform the change of variables u(x, t) = v(x, t) + tp(x), then
N d2vr . ,". , , duf . dv, .
dxi(x,t)=dxi(x,t) + xp"(x) and Tt(x't)=Tt(x't)-
We substitute these expressions into Eq. (1), so that the equation takes the form
k-—(x, t) + k(x) + r= -—(x, t). ax2 at
Let the function satisfy the equation
kip"(x) + r = 0 or ip "(x) = -j, (4)
with the boundary conditions
tp (0) = 0 and tp (1) = u0. (5)
Then for v(x, t) we obtain a problem of solving the homogeneous parabolic equation
, d2v , . dv .
(-,^)-Tt(x,t) = °, (6)
with the boundary conditions
v(0, t) = 0 and v(1, t) = 0 (7)
and the initial condition
v(x,0) = f(x)-^(x). (8)
Equation (4) is integrated two times, as a result of which we arrive at its solution in the general form
r
Tp(x) = -2^2 + Clx + c2. (9)
Then to determine the constants c1 and c2, we use the boundary conditions
tp(0) = 0 and tp(1) = u0,
distributed over the summands of the sought solution. Substituting these values into the general solution (9), we calculate the values of constants c1 and c2:
C2 = 0, q = 2^ + Uo.
Therefore, the final partial solution of the Dirichlet-Dirichlet problem (4)-(5) for the simplest Poisson equation (4) has the form:
^x) = -kx2 + (-k+uo)x- (10)
From the initial condition u(x,0) = v(x,0) + ip(x) it follows that v(x,0) = u(x,0) — ip(x) = f(x)-ip(x). Therefore, to determine v(x, t), we solve a new boundary value problem
, d2v. . dv. .
k-^(x,t)=-z-(x,t), 0<x<1, t>0, ox2 at
v(0, t) = 0, v(1, t) = 0, t>0, (11)
v(x,0) = f(x)+2^x2 -(^+Uo)x, 0<x<1,
using the method of separation of variables. This method yields a solution to this problem in the form
TO
v(x, t)=Yj Ane-kn2x2t sin nnx, (12)
n=l
where the coefficients An are calculated by the formula
f1 r r
An = 2 I [f(x) + -~^x2 — (-¿¡¿ + u°) x\ sin nnxdx. (13)
0
As a result, the solution to the original problem (1) is obtained by summing the solutions $(x) and v(x, t) of the homogeneous boundary value problem (4)-(5) and the boundary value problem (11)
TO
u(x,t) = -2fcx2 + (2k+u0)x+ ZAne-knVt sin nnx. (14)
n=1
Note that u(x, t) ^ $(x) as t ^ œ in expression (14). That is why the function $(x) as a part of the solution of the heat conduction equation that does not change with time is called a steady solution. Whereas the function v(x, t) ^ 0 as t ^ œ is therefore called a transition solution.
The determination of the coefficients An according to formula (13) can be implemented using the integration technique based on the interpolation of the integrand with Chebyshev polynomials of the first kind on Gauss-Lobatto grids [13].
Example 2. In the next example, we consider an approach to solving an inhomogeneous boundary value problem when the Dirichlet condition is specified at the left end, i.e. the left end is in perfect contact with the heat bath (welded or screwed to a massive holder having a constant temperature u0), and at the right end of the interval the Robin condition is imposed, i.e. this end exchanges heat with the environment at temperature um (hangs freely in the environment). Robin boundary conditions arise, e.g., when the ends are immersed in some liquid or gaseous medium. The initial temperature distribution f(x) along the rod length 0 < x <L is also assumed to be given.
In accordance with the decomposition method, it is proposed to represent the solution to this inhomogeneous boundary value problem under given boundary and initial conditions
S2U du
k--(x,t)—^(x,t) = —F(x), 0<x<L, t>0, ox2 at
u(0, t) = u0, u0 = const,
(15)
= -h(u(L, t) — um), h > 0 and um = const,
x=L
u(x, 0) = f(x), 0 < x < L,
du dx
as a combination of two terms
u(x, t) = v(x, t) + Tp(x), (16)
each of them being a solution to a separate boundary value problem, respectively:
Problem A2. The function of the spatial variable $(x) is a solution to an inhomogeneous ordinary differential equation with the Dirichlet boundary condition at the left end of the interval and the Robin boundary condition at the right end, namely:
d2é(x)
k 7; )(x) + F(x) = 0, 0 <x<L, t>0, ox2
ip(0) = u0, u0 = const, (17)
dip
dx
+ hip(L) = hui, h > 0 and Ui = const.
x=L
Problem B2. The second term, the function of two variables v(x, t), is a solution to a homogeneous boundary value problem with zero Dirichlet-Robin boundary conditions and given initial condition:
, d2v , . dv .
kj~2(x,t)- -¡^(x,t) = 0, 0<x<L, t>0, v(0, t) = 0,
dv\ (18)
+ hv(L,t) = 0, h > 0,
dxlx=L
v(x, 0) = f(x) - Tp(x), 0 <x<L. Let us consider sequentially the methods for solving each of these problems.
4. Solving Problem A2 by the Chebyshev collocation method
By analogy with the method of approximate solution of the Poisson problem with various boundary conditions, such as the Dirichlet conditions on both ends of the interval u(a) = a, u(b) = ¡3, Neumann-Dirichlet condition u'(a) = a, u(b) = ¡3, or Dirichlet-Neumann condition u(a) = a, u'(b) = ¡3 thoroughly studied in Ref. [13], let us consider the solution of problem (17) based on the Chebyshev collocation method and, therefore, polynomial interpolation of the solution.
The polynomial interpolation based on using the basis of Chebyshev polynomials, leads to a necessity to formulate the problem in the interval [-1,1] instead of the initial interval [a, b]. In
the new interval, the polynomial interpolant of the function y(x) is specified in the form of a series expansion y(x) ~ ckTk(x) in Chebyshev polynomials of the first kind Tk(x) with the domain of definition on the segment x e [-1,1].
The transition x^ t to new arguments, from x e [a, b] to t e [-1,1], is implemented using a linear transformation
_ 2x — (b +a) b — a
and, if necessary, the inverse transformation
_ t(b — a) + (b + a) x_ 2 •
In this case, the values of the function, the integrals and derivatives are recalculated using the formulas
f(x)of(^t + ^), xe[a,b], tei—1,1],
f"f(x)dx_ f xe[a,b], te[—1,1],
ja J-i \ '
dx b — a dt _ 2 ,
fi(x) o (t)/ i-j, xe[a,b], te [—1,1], and [14] the upper estimate of the interpolation error has the form
y(x) - 2 ckTk(x) k=0
1 \b-a\n+1
—r max \yn+1(^)\.
Thus, Chebyshev interpolation provides an almost optimal approximation in the sense of the Lœ norm and almost optimal to the L2 norm. In addition, the use of Gauss-Lobatto nodes as interpolation nodes leads to optimal integration formulas.
As in Ref. [6], the method for approximate solution of problem (17) for a second-order ODE consists of sequential solution of several subproblems.
1. calculation of spectral coefficients of polynomial Chebyshev interpolation of the second derivative of the solution - the function of the right-hand side (17) on the Gaussian-Lobatto grid - interpolation of F(x) in the basis of Chebyshev polynomials of the first kind;
2. calculation of those coefficients of the desired solution (except for the first two) that are determined from the differential conditions of the problem (allowing the solution to satisfy the differential conditions) - multiplication of the inverse (with respect to the matrix of spectral Chebyshev differentiation) matrix by the vector of interpolation coefficients;
3. additional determination of the first few coefficients of the solution based on boundary (or other independent additional) conditions.
Let us represent the approximate solution in the form of a finite series of orthogonal Chebyshev polynomials
n
p(x)=2ckTk(x), xe[-1,1]. (19)
k=0
Let us differentiate function (19) twice. The expression for the second derivative has the form
n n
p"(x) = 2 ckT{"(x) = 2 bkTk(x), X G [-1,1]. (20)
k=0 k=0
Using recurrence relations satisfied by Chebyshev polynomials of the first kind and their derivatives [5], [15] and equating the coefficients for identical polynomials in (20), we arrive [5]
at the following dependence of the coefficients q, i = 2,3,... ,n, bk, k = 0,...,n:
D+D+b = c
(21)
where D+ is a generalized inverse matrix with respect to the Chebyshev differentiation matrix in the spectral space [15], [16].
D+D+b
_1 _2
0 1 6 0
— 0
0 -1/2 (n-3) 0
1/2 0 -1/2
(n-2) (n-2)
0 1/2 (n-1) 0
0 0 1/2
n
0 0
0
0
0
0
0 -1/2 (n-1) 0
' b0 ' Co
b1
b2 C2
h C3
b4 = C4
bn-3 cn-3
bn-2 cn-2
bn-1 cn-1
- bn . - Cn -
(22)
Hence, the vector of coefficients {c2, c3,..., cn} is the result of double multiplication of a simple tridiagonal matrix D+ (inverse of the differentiation matrix) by vector {b0, b1,..., bn}.
At the third stage of solving the problem, the first two coefficients of the expansion of the desired solution in Chebyshev polynomials are determined.
5. Dirichlet-Robin boundary conditions
For a one-dimensional problem considered on the interval [-1,1], the Dirichlet-Robin conditions look as follows
ccp(-1) = g(-1l ¡3p(1) + yp'(1) = g(1). ( '
Here a, ¡3, y are given constants. The sign in front of the term with the derivative at the right boundary point is positive, since the outer normal to the domain of definition at the right boundary point is directed to = to, i.e., in the positive direction.
Let us take into account that the derivatives of Chebyshev polynomials of the first kind are simply expressed in terms of polynomials of the second kind:
dT" „TT
and, in addition, the relations
T„(-1) = i-1)n, Tn(1) = 1, Un(-1) = (-1)n(n + 1), Un(1) = (n + 1),
are valid. In this case, the system of equations for calculating the unknown expansion coefficients of the solution has the form
>_q + Zck(_i)k) = s(_i)
k=2
p(c0 + C1 + ^ ck\ + yici + ^ ckk2\ = g(1). \ k=2 ) \ k=2 )
2
X
Let us introduce the notation
n
a(c0 - a) = g(-1) -a^ ck(-1)k = G(-1),
k=2
(25)
(26)
(27)
(28)
¡3(C0 + Ci) + yci = g(1) -P2°k -7 2 ^ = G(1),
k=2 k=2
or
a(co - ci) = G(-1), fro + (£ + r)ci = G(1).
Remove the brackets
ac0 - ac1 = G(-1), ¡3co + fa + yci = G(1) and multiply the first equation by fi and the second equation by a:
fiac0 - ¡3ac1 = @G(-1), Pac0 + + yac1 = aG(1).
To calculate the coefficient c1, subtract the first equation from the second one:
2^ac1 + yac1 = aG(1) - @G(-1), (29)
from which it follows that
„ = aGiD-iSG-)
C1 = <213 + y) . (30)
To calculate the coefficient c0, we substitute the calculated value of c1 into the first equation (27):
= G(-1) = Gj-1) aG(1) - j3G(-1) °0 = a + C1 = a + cc(2j3 + y) '
Q3+y)G(-1) + aG(1) Co = . (31)
6. Solving the problem B2. Method of separation of variables
Recall the formulation of problem (18). It is required to find a solution of the homogeneous boundary value problem with zero Dirichlet-Robin conditions and a given initial condition:
k--(x,t)-lr(x,t) = 0, 0 < x < L, t>0, ox2 at
v(0, t) = 0,
dv \
+hv(L,t) = 0, h > 0,
dx\X=L
v(x, 0) = f(x) - Tp(x), 0 <x<L. where k is the heat conductivity coefficient depending on the material properties.
The method of separation of variables yields a solution to this problem (see, e.g., [17]) in the form
TO
v(x, t)= 2 Ane-kÀnt sin JTnLx, (32)
n=1
where the coefficients An are calculated as
An = 21 U(x) + $(x)] sin ^X^Lxdx, Jo
(33)
where Xn >0 is the set of positive solutions of the equation —h tan(V^), n g N.
Finally, the solution of the initial problem 2 (15) is obtained by adding the solution ^(x) of the boundary value problem (17) and the solution v(x, t) of the homogeneous boundary value problem (18):
œ
u(x, t) = xp{x) + 2 Ane-kÀnt sin JXnLx. (34)
n=1
The determination of coefficients An by formula (33) can be implemented using the technique of integration based on the interpolation of the integrand by the Chebyshev polynomials of the first kind on Gauss-Lobatto grids. In this case, the calculation of the coefficients An reduces to elementary summation of the weighted even coefficients of the interpolant.
7. Conclusion
Among the numerical algorithms for solving initial and boundary value problems for linear ODEs of the first and second order, there are many methods that use the initial approximation (boundary conditions) as the initially active condition that determines all further solution of the problem. These are methods such as Euler, Adams-Bashforth, Runge-Kutta, etc. [18]. Other methods, based on approximation of the solution using global functions [1-5], are based on the construction of systems of equations that simultaneously include both initial (boundary) conditions and conditions that specify the behavior of the derivatives of the desired solution.
The solution to the main inhomogeneous initial-boundary value problem for a one-dimensional parabolic equation is presented in the form of a sequential solution of several subproblems. As a method for solving one of the subproblems - an inhomogeneous ordinary differential equation with Dirichlet-Robin boundary conditions - it is proposed to use the stable and efficient spectral method of Chebyshev collocation.
Polynomial interpolation of the desired solution by Chebyshev polynomials is carried out in several stages. At the first stage, a general solution is identified, i.e., a set of solutions that satisfies the differential equation, but does not necessarily satisfy the initial (boundary) conditions. Considering the initial (boundary) conditions is carried out at the last stage of solving the problem and reduces to solving a linear equation with two unknown coefficients.
The search for a general solution to an inhomogeneous ODE reduces to multiplying the transposed matrix of values of Chebyshev functions on the Gauss-Lobatto grid by the vector of the function values that specifies the right-hand side of the original differential equation to determine the interpolation coefficients for the expansion of the solution derivative. Next, multiplying the bi-diagonal integration matrix [5], [15] by the vector of these coefficients leads to obtaining all the coefficients of the desired solution, except for first ones. At the final stage, the use of the initial (boundary) condition makes it possible to determine the first two coefficients of the polynomial expansion of the solution.
The solution of the second homogeneous subproblem is carried out by the traditional method of separation of variables. In this case, to calculate the coefficients of expansion of the solution according to the basis of the Sturm-Liouville problem, an effective and stable spectral method of Chebyshev collocation is used. Thus, the authors expand the scope of applicability of the developed 2-stage Chebyshev collocation method.
Author Contributions: Conceptualization, software: Konstantin P. Lovetskiy; software: Stepan V. Sergeev; formal analysis: Dmitry S. Kulyabov;theoretical proof of the consistency of Chebyshev collocation method in the case of complex potentials: Leonid A. Sevastianov. All authors have read and agreed to the published version of the manuscript.
Funding: This research was funded by the RUDN University Scientific Projects Grant System, project No 021934-0-000 (Konstantin P. Lovetskiy). This research was supported by the RUDN University Strategic Academic Leadership Program (Dmitry S. Kulyabov). The work of Leonid A. Sevastianov was supported by the Russian Science Foundation (grant No. 20-11-20257). Data Availability Statement: Data sharing is not applicable.
Conflicts of Interest: The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data;in the writing of the manuscript;or in the decision to publish the results.
References
1. Boyd, J. P. Chebyshev and Fourier spectral methods second (Dover Books on Mathematics, 2013).
2. Mason, J. C. & Handscomb, D. C. Chebyshev polynomials doi:10 . 1201 / 9781420036114 (Chapman and Hall/CRC Press, New York, 2002).
3. Greengard, L. Spectral integration and two-point boundary value problems. SIAM Journal on Numerical Analysis 28, 1071-1080. doi:10.1137/0728057 (1991).
4. Shen, J., Tang, T. & Wang, L. Spectral methods doi:10. 1007/978-3-540-71041-7 (Springer Berlin Heidelberg, Berlin, Heidelberg, 2011).
5. Fornberg, B. A practical guide to pseudospectral methods (Cambridge University Press, New York, 1996).
6. Sevastianov, L. A., Lovetskiy, K. P. & Kulyabov, D. S. Multistage collocation pseudo-spectral method for the solution of the first order linear ODE in VIII International Conference on Information Technology and Nanotechnology (ITNT) (2022), 1-6. doi:10.1109/ITNT55410.2022.9848731.
7. Sevastianov, L. A., Lovetskiy, K. P. & Kulyabov, D. S. A new approach to the formation of systems of linear algebraic equations for solving ordinary differential equations by the collocation method. Izvestiya of Saratov University. Mathematics. Mechanics. Informatics 23. in Russian, 3647. doi:10.18500/1816-9791-2023-23-1-36-47 (2023).
8. Lovetskiy, K. P., Kulyabov, D. S. & Hissein, W. Multistage pseudo-spectral method (method of collocations) for the approximate solution of an ordinary differential equation of the first order. Discrete and Continuous Models and Applied Computational Science 30,127-138. doi:10. 22363/2658-4670-2022-30-2-127-138 (2022).
9. Sevastianov, L. A., Lovetskiy, K. P. & Kulyabov, D. S. Numerical integrating of highly oscillating functions: effective stable algorithms in case of linear phase 2021. doi:10.48550/arXiv. 2104. 03653.
10. Lovetskiy, K. P., Kulyabov, D. S., Sevastianov, L. A. & Sergeev, S. V. Chebyshev collocation method for solving second order ODEs using integration matrices. Discrete and Continuous Models and Applied Computational Science 31, 150-163. doi:10.22363/2658-4670-2023-31-2-150-163 (2023).
11. Lienhard, J. H. I. & Lienhard, J. H. V. A heat transfer textbookfifth edition 2020.
12. Tikhonov, A. N. & Samarskii, A. A. Equations of Mathematical Physics (Nauka, M., 2004).
13. Sevastianov, L. A., Lovetskiy, K. P. & Kulyabov, D. S. A new approach to the formation of systems of linear algebraic equations for solving ordinary differential equations by the collocation method. Izvestiya of Saratov University. Mathematics. Mechanics. Informatics 23, 36-47. doi:10. 18500/1816-9791-2023-23-1-36-47 (2023).
14. Stewart, G. W. Afternotes on numerical analysis (Society for Industrial and Applied Mathematics, USA, 1996).
15. Amiraslani, A., Corless, R. M. & Gunasingam, M. Differentiation matrices for univariate polynomials. Numerical Algorithms 83,1-31. doi:10.1007/s11075-019-00668-z (2020).
16. Rezaei, F., Hadizadeh, M., Corless, R. & Amiraslani, A. Structural analysis of matrix integration operators in polynomial bases. Banach Journal of Mathematical Analysis 16, 5. doi:10.1007/ s43037-021-00156-4 (2022).
17. Boyce, W. E. & DiPrima, R. C. Elementary differential equations and boundary value problems 9th Edition (Wiley, New York, 2009).
18. Planitz, M. et al. Numerical recipes: the art of scientific computing 3rd Edition (Cambridge University Press, New York, 2007).
To cite: Lovetskiy K. P., Sergeev S. V., Kulyabov D. S., Sevastianov L. A., Application of the Chebyshev collocation method to solve boundary value problems of heat conduction, Discrete and Continuous Models and Applied Computational Science 32 (1)(2024)74-85.D01:10.22363/2658-4670-2024-32-1-74-85.
Information about the authors
Lovetskiy, Konstantin P.—Candidate of Sciences in Physics and Mathematics, Associate Professor of Department of Computational Mathematics and Artificial Intelligence of Peoples' Friendship University of Russia named after Patrice Lumumba (RUDN University) (e-mail: [email protected], phone: +7(495)952-25-72, ORCID: https://orcid.org/0000-0002-3645-1060) Kulyabov, Dmitry S.—Professor, Doctor of Sciences in Physics and Mathematics, Professor at the Department of Probability Theory and Cyber Security of Peoples' Friendship University of Russia named after Patrice Lumumba (RUDN University);Senior Researcher of Laboratory of Information Technologies, Joint Institute for Nuclear Research (e-mail: [email protected],
phone: +7(495)952-02-50, ORCID: https://orcid.org/0000-0002-0877-7063)
Sevastianov, Leonid A.—Professor, Doctor of Sciences in Physics and Mathematics, Professor at the Department of Computational Mathematics and Artificial Intelligence of Peoples' Friendship University of Russia named after Patrice Lumumba (RUDN University), Leading Researcher of Bogoliubov Laboratory of Theoretical Physics, Joint Institute for Nuclear Research (e-mail: [email protected], phone: +7(495)952-25-72, ORCID: https://orcid.org/0000-0002-1856-4643) Sergeev, Stepan V.—PhD student of Department of Computational Mathematics and Artificial Intelligence of Peoples' Friendship University of Russia named after Patrice Lumumba (RUDN University) (e-mail: [email protected], ORCID: https://orcid.org/0009-0004-1159-4745)
УДК 519.6:004.94
PACS 07.05.Tp, 02.60.Pn, 02.70.Bf
DOI: 10.22363/2658-4670-2024-32-1-74-85 EDN: BUEBFE
Применение метода коллокации Чебышева для решения граничных задач теплопроводности
К. П. Ловецкий1, С. В. Сергеев1, Д. С. Кулябов1,2, Л. А. Севастьянов1,2
1 Российский университет дружбы народов,
ул. Миклухо-Маклая, д. 6, Москва, 117198, Российская Федерация
2 Объединённый институт ядерных исследований,
ул. Жолио-Кюри, д. 6, Дубна, Московская область, 141980, Российская Федерация
Аннотация. Для одномерных неоднородных (по пространственной переменной) линейных параболических уравнений используется комбинированный подход, разбивающий исходную задачу на две подзадачи. Первая из них - неоднородная одномерная задача Пуассона с граничными условиями Дирихле-Робена, поиск решения которой основан на методе чебышевской коллокации. Метод разработан на основе ранее опубликованных алгоритмов решения обыкновенных дифференциальных уравнений, в которых решение ищется в виде разложения по полиномам Чебышева 1-го рода на сетках Гаусса-Лобатто, что позволяет использовать дискретную ортогональность полиномов. Такой подход оказывается весьма экономичным и стабильным по сравнению с традиционными методами, приводящими к решению часто плохо определенных систем линейных алгебраических уравнений. В описываемом подходе удачное применение матриц интегрирования позволяет вообще избавиться от необходимости работы с плохо обусловленными матрицами.
Вторая, однородная задача теплопроводности решается методом разделения переменных. При этом отыскание коэффициентов разложения искомого решения по полному набору решений соответствующей задачи Штурма-Лиувилля сводится к вычислению интегралов от известных функций. Простая методика построения чебышевских интерполянтов подынтегральных функций позволяет вычислять интегралы суммированием интерполяционных коэффициентов.
Ключевые слова: начально-краевые задачи, псевдоспектральный метод коллокации, полиномы Чебышева, множества Гаусса-Лобатто, численная устойчивость, разделение переменных
