EUROPEAN JOURNAL OF PURE AND APPLIED MATHEMATICS Vol. 16, No. 4, 2023, 2384-2396 ISSN 1307-5543 – ejpam.com Published by New York Business Global Cubic Trigonometric B-spline Method for Solving a Linear System of Second Order Boundary Value Problems Ahmed Salem Heilat Department of Mathematics, Jadara University, P.O. Box(733), 21111 Irbid, Jordan Abstract. This paper introduces a novel trigonometric B-spline collocation method for solving a specific class of second-order boundary value problems. The study showcases the method’s practicality and effectiveness through various numerical examples. Furthermore, it evaluates the technique’s performance by calculating maximum errors for different step sizes in the spatial do- main. The paper also conducts a comparative analysis with alternative methods, demonstrating the superior accuracy of the trigonometric B-spline approach. 2020 Mathematics Subject Classifications: 34A30, 65L10 Key Words and Phrases: Linear system of second order boundary value problems, boundary value problems, cubic trigonometric B-spline basis functions. 1. Introduction System of second order ordinary differential equation associated with boundary condi- tions arises in several applications and have been discussed in [1–22]. Issue of existence of solution to such system has been discussed in [1–3]. The standard numerical approaches to solve second order boundary value problems are to use the shooting and finite differ- ence methods. In this paper, we focus on the use of cubic trigonometric B-spline method (CTBSM). Consider the following system of second order boundary value problems:{ u′′(x) + a1(x)u ′(x) + a2(x)u(x) + a3(x)v ′′(x) + a4(x)v ′(x) + a5(x)v(x) = f1(x) v′′(x) + b1(x)v ′(x) + b2(x)v(x) + b3(x)u ′′(x) + b4(x)u ′(x) + b5(x)u(x) = f2(x). (1) The following boundary conditions are associated with the system u(a) = u(b) = 0, v(a) = v(b) = 0, (2) DOI: https://doi.org/10.29020/nybg.ejpam.v16i4.4947 Email address: ahmed heilat@yahoo.com (Ahmed Salem Heilat) https://www.ejpam.com 2384 © 2023 EJPAM All rights reserved. A. S. Heilat / Eur. J. Pure Appl. Math, 16 (4) (2023), 2384-2396 2385 where a ≤ x ≤ b, f1(x) and f2(x) are continuous functions, ai(x) and bi(x), for i = 1, 2, 3, 4, 5, are sufficiently smooth real-valued functions of x. Many approximate techniques have been suggested to solve linear and nonlinear sys- tems of second order BVPs. In [4], Saadatmandi et al. introduced a homotopy pertur- bation approach aimed at addressing nonlinear second-order boundary value problems. This method generates solutions in the form of convergent series with readily computable or minor perturbations. Ogunlaran and Ademola [5] introduced the Laplace homotopy analysis method as an alternative for solving the same system. This technique obtains solutions without resorting to discretization or imposing restrictive assumptions. The re- sulting solution manifests in the form of a rapidly convergent series. Lu [6] introduced a variational iteration approach to solve systems analogous to problem (1)-(2). [7, 8] pre- sented a method to obtain the analytical and approximate solutions in the form of series in the reproducing kernal space for linear and nonlinear system of second order BVPs. In [9, 10], Dehghan, Saadatmandi, and El-Gamel utilized a numerical approach employing the Sinc-collocation method for solving second-order nonlinear systems. This method is acknowledged as effective for problems with potential singularities, infinite domains, or boundary layers. Furthermore, it is applied to streamline the computation of solutions for the system represented by equations (1)-(2) into a set of algebraic equations. Local ra- dial basis function based differential quadrature collocation method [11] is also suggested to solve system of second order BVPs, this method is easy to implementation and the results is more efficient. Saadatmandi and J. Askari [12] solved similar systems using the Chebyshev finite element method. Moreover, the continuous genetic algorithm, as described in reference [13], has proven to be successful in addressing our problem. Addi- tionally, a comprehensive analysis of the convergence and sensitivity of genetic operators and control parameters for the algorithm has been conducted. Finally, [14–16] are also proposed different methods to solve the system (1)-(2). The aim of this research is to use CTBSM for solving linear system of second order BVPs. The system had already been treated using spline collocation approach [17], cubic B-spline scaling functions [18], B-spline method [19], Extended cubic B-spline method [20], Hybrid cubic B-spline method [21]. CTBSM was developed to solve dynamical systems in [22], the results are promising. In our present work, a new trigonometric B-spline technique is described and presented for solving a linear system of second order BVPs. The technique is based on the cubic trigonometric B-spline functions. CTBSM is used as an interpolating function in the space dimension. The regular system is attempted using a trigonometric B-spline collocation approach constructed over uniform mesh. The efficiency and applicability of the technique are demonstrated by applying the scheme to solve several examples. The numerical results demonstrate that this method is superior as it yields more accurate solutions. The outline of this study is as follows: In section 2, the CTBSM is utilized as an interpo- lating function in the space dimension for solving a linear system of second order BVPs. A numerical solution of linear system of second order BVPs is presented in section 3. Numerical results are also considered in section 4 to show the achievable of the proposed method. Finally, the concluding remarks of this study are in section 5. A. S. Heilat / Eur. J. Pure Appl. Math, 16 (4) (2023), 2384-2396 2386 2. Cubic Trigonometric B-Spline Method In this section, the CTBSM is applied for the numerical solution of second order BVPs (1)-(2). Let a = x0 < x1 < ... < xn = b be the mesh point in the interval [a, b] such that xi = a+ ih, i = 0, 1, 2, ..., n with h = b−a n . Let the approximate solution U(x) and V (x) to the exact solution u(x) and v(x) at the point xi respectively, and can be defined as:  U(x) = n−1∑ i=−3 αiT 4 i (x), x ∈ [x0, xn], αi ∈ R, V (x) = n−1∑ i=−3 βiT 4 i (x), x ∈ [x0, xn], βi ∈ R. (3) where αi, βi are real coefficient to be determined for the approximated solutions U(x) and V (x) respectively and Ti(x) are trigonometric B-spline basis functions which have C2 continuity at each knots and defined over the mesh by [21–23] T 4 i (x) = 1 ω  p3(xi), x ∈ [xi, xi+1], p(xi)(p(xi)q(xi+2) + q(xi+3)p(xi+1)) + q(xi+4)p 2(xi+1), x ∈ [xi+1, xi+2], q(xi+4)(p(xi+1)q(xi+3) + q(xi+4)p(xi+2)) + p(xi)q 2(xi+3), x ∈ [xi+2, xi+3], q3(xi+4), x ∈ [xi+3, xi+4]. (4) where, p(xi) = sin(x−xi 2 ), q(xi) = sin(xi−x 2 ), ω = sin(h2 ) sin(h) sin( 3h 2 ). For the numerical solution of proposed problem, the values of Ti(x) and its derivatives T ′ i (x), T ′′ i (x) at nodal points are need and recorded in Table. 1. Table 1: Coefficient of Ti(x), T ′ i (x) and T ′′ i (x) x xi−1 xi xi+1 Ti(x) m1(x) m2(x) m1(x) T ′ i (x) m3(x) 0 m4(x) T ′′ i (x) m5(x) m6(x) m5(x) where m1(x) = sin2(h2 ) sin(h) sin(3h2 ) , m2(x) = 2 1 + 2 cos(h) , m3(x) = 3 4 sin(3h2 ) , m4(x) = −3 4 sin(3h2 ) , m5(x) = 3(1 + 3 cos(h)) 16 sin2(h2 )(2 cos( h 2 ) + cos(3h2 )) , m6(x) = −3 cot2(h2 ) 2 + 4 cos(h) . A. S. Heilat / Eur. J. Pure Appl. Math, 16 (4) (2023), 2384-2396 2387 From equation (3) and (4), the values of U(x), U ′ (x), U ′′ (x), V (x), V ′ (x), and V ′′ (x) at the knots xi are determined in the terms of αi and βi as U(xi) = m1(x)αi−3 +m2(x)αi−2 +m1(x)αi−1, U ′ (xi) = m3(x)αi−3 +m4(x)αi−1, U ′′ (xi) = m5(x)αi−3 +m6(x)αi−2 +m5(x)αi−1. (5) and V (xi) = m1(x)βi−3 +m2(x)βi−2 +m1(x)βi−1, V ′ (xi) = m3(x)βi−3 +m4(x)βi−1, V ′′ (xi) = m5(x)βi−3 +m6(x)βi−2 +m5(x)βi−1. (6) 3. Solution of System of Second Order Boundary Value Problems In this section, a numerical solution of a class of system of second order BVPs (1)-(2) is obtained using collocation approach based on cubic trigonometric basis functions. It is necessary to satisfy the proposed system of second order differential equation (1)-(2) by putting the approximate (3) in (1) at x = xi, it follows as: n−1∑ i=−3 αi[T ′′ i (xi) + a1(xi)T ′ i (xi) + a2(xi)Ti(xi)]+ n−1∑ i=−3 βi[a3(xi)T ′′ i (xi) + a4(xi)T ′ i (xi) + a5(xi)Ti(xi)] = f1(xi), i = 0, 1, ..., n (7) n−1∑ i=−3 βi[T ′′ i (xi) + b1(xi)T ′ i (xi) + b2(xi)Ti(xi)]+ n−1∑ i=−3 αi[b3(xi)T ′′ i (xi) + b4(xi)T ′ i (xi) + b5(xi)Ti(xi)] = f2(xi), i = 0, 1, ..., n (8) The boundary conditions (2) can be written as: n−1∑ i=−3 αiTi(xi) = 0, x = 0, n (9) n−1∑ i=−3 βiTi(xi) = 0, x = 0, n (10) The values and derivatives of spline functions at knots xi are determined using Table 1 and putting in (7)-(10). A system of equations with unknown α−3, α−2, α−1, ..., αn−1, β−3, β−2, β−1, ..., βn−1 A. S. Heilat / Eur. J. Pure Appl. Math, 16 (4) (2023), 2384-2396 2388 of order (2n+3)× (2n+3) is obtained and can be solved by Thompson algorithm [24] for the trigonometric spline solution of proposed problem. This system can be written in the matrix-vector form as follows AB = C (11) whereB = [α−3, α−2, ..., αn−1, β−3, β−2, ..., βn−1] T , C = [0, f1(x0), ..., f1(xn), 0, 0, f2(x0), ..., f2(xn), 0] T , and A is a 2(n+ 3)× 2(n+ 3) matrix given by A = M1 | M2 · · · · · · · · · M4 | M3 . where MK ,K = 1, 2, 3, 4 are four submatrices, can be defined as: M1 =  m1(x) m2(x) m1(x) 0 · · · 0 0 ϵ1(x0) ρ1(x0) τ1(x0) 0 · · · 0 0 0 ϵ1(x1) ρ1(x1) τ1(x1) 0 · · · 0 ... ... ... ... ... ... ... 0 · · · 0 0 ϵ1(xn) ρ1(xn) τ1(xn) · · · · m1(x) m2(x) m1(x)  (n+3)×(n+3) M2 =  0 0 0 0 · · · 0 0 ϵ2(x0) ρ2(x0) τ2(x0) 0 · · · 0 0 0 ϵ2(x1) ρ2(x1) τ2(x1) 0 · · · 0 ... ... ... ... ... ... ... 0 · · · 0 0 ϵ2(xn) ρ2(xn) τ2(xn) · · · · 0 0 0  (n+3)×(n+3) M3 =  m1(x) m2(x) m1(x) 0 · · · 0 0 ϵ3(x0) ρ3(x0) τ3(x0) 0 · · · 0 0 0 ϵ3(x1) ρ3(x1) τ3(x1) 0 · · · 0 ... ... ... ... ... ... ... 0 · · · 0 0 ϵ3(xn) ρ3(xn) τ3(xn) · · · · m1(x) m2(x) m1(x)  (n+3)×(n+3) M4 =  0 0 0 0 · · · 0 0 ϵ4(x0) ρ4(x0) τ4(x0) 0 · · · 0 0 0 ϵ4(x1) ρ4(x1) τ4(x1) 0 · · · 0 ... ... ... ... ... ... ... 0 · · · 0 0 ϵ4(xn) ρ4(xn) τ4(xn) · · · · 0 0 0  (n+3)×(n+3) The input of M1, M2, M3, and M4 are given below for i = 0, 1, ..., n. ϵ1(xi) = 3(1 + 3 cos(h)) 16 sin2(h2 )(2 cos( h 2 ) + cos(3h2 )) + a1(xi) −3 4 sin(3h2 ) + a2(xi) sin2(h2 ) sin(h) sin(3h2 ) A. S. Heilat / Eur. J. Pure Appl. Math, 16 (4) (2023), 2384-2396 2389 ρ1(xi) = −3 cot2(h2 ) 2 + 4 cos(h) + a1(xi)(0) + a2(xi) 2 1 + 2 cos(h) τ1(xi) = 3(1 + 3 cos(h)) 16 sin2(h2 )(2 cos( h 2 ) + cos(3h2 )) + a1(xi) 3 4 sin(3h2 ) + a2(xi) sin2(h2 ) sin(h) sin(3h2 ) ϵ2(xi) = a3(xi) 3(1 + 3 cos(h)) 16 sin2(h2 )(2 cos( h 2 ) + cos(3h2 )) + a4(xi) −3 4 sin(3h2 ) + a5(xi) sin2(h2 ) sin(h) sin(3h2 ) ρ2(xi) = a3(xi) −3 cot2(h2 ) 2 + 4 cos(h) + a4(xi)(0) + a5(xi) 2 1 + 2 cos(h) τ2(xi) = a3(xi) 3(1 + 3 cos(h)) 16 sin2(h2 )(2 cos( h 2 ) + cos(3h2 )) + a4(xi) 3 4 sin(3h2 ) + a5(xi) sin2(h2 ) sin(h) sin(3h2 ) ϵ3(xi) = 3(1 + 3 cos(h)) 16 sin2(h2 )(2 cos( h 2 ) + cos(3h2 )) + b1(xi) −3 4 sin(3h2 ) + b2(xi) sin2(h2 ) sin(h) sin(3h2 ) ρ3(xi) = −3 cot2(h2 ) 2 + 4 cos(h) + b1(xi)(0) + b2(xi) 2 1 + 2 cos(h) τ3(xi) = 3(1 + 3 cos(h)) 16 sin2(h2 )(2 cos( h 2 ) + cos(3h2 )) + b1(xi) 3 4 sin(3h2 ) + b2(xi) sin2(h2 ) sin(h) sin(3h2 ) ϵ4(xi) = b3(xi) 3(1 + 3 cos(h)) 16 sin2(h2 )(2 cos( h 2 ) + cos(3h2 )) + b4(xi) −3 4 sin(3h2 ) + b5(xi) sin2(h2 ) sin(h) sin(3h2 ) ρ4(xi) = b3(xi) −3 cot2(h2 ) 2 + 4 cos(h) + b4(xi)(0) + b5(xi) 2 1 + 2 cos(h) τ4(xi) = b3(xi) 3(1 + 3 cos(h)) 16 sin2(h2 )(2 cos( h 2 ) + cos(3h2 )) + b4(xi) 3 4 sin(3h2 ) + b5(xi) sin2(h2 ) sin(h) sin(3h2 ) Error analysis of the boundary value problem (1) and implementation on method is dis- cussed in the following section through three numerical examples. 4. Numerical Results and Discussions Some numerical examples are considered in this section to demonstrate the competency of the proposed trigonometric spline method. Numerical results found by the method are compared with existing methods in the literature such as [6] and [19] and with the ana- lytical solution at knots xi using different N . It was establish that proposed method in contrast with these methods is more accurate. Absolute errors of the proposed method is calculated numerically by the following formula [25–27] L∞ = max i | u(xi)− U(xi) | or L∞ = max i | v(xi)− V (xi) | A. S. Heilat / Eur. J. Pure Appl. Math, 16 (4) (2023), 2384-2396 2390 L2 = √√√√ N∑ i=1 (u(xi)− U(xi)) 2 or L2 = √√√√ N∑ i=1 (v(xi)− V (xi)) 2 Example 1.Consider the following system [6]  u ′′ (x) + (2x− 1)u ′ (x) + cos(πx)v ′ (x) = f1(x) v ′′ (x) + xu(x) = f2(x) u(0) = u(1) = 0, v(0) = v(1) = 0, (12) where 0