## AKSIOMA: Jurnal Matematika dan Pendidikan Matematika
Vol. 13, No. 1 April 2022 e-ISSN 2579-7646
Stabilitas metode runge-kutta dengan interpolasi cubic-spline untuk menyelesaikan delay differential equations (dde)
1 Rahma Qudsi, 2 Agus Dahlia
1 Fakultas Keguruan dan Ilmu Pendidikan, Universitas Islam Riau
2 Fakultas Teknik, Universitas Islam Riau Email korespondensi: agus.dahlia@edu.uir.ac.id

## Abstrak
Algoritma  dari  suatu  metode  numerik  terkadang  menghasilkan  hasil  yang tidak konvergen ketika input yang dimasukkan berubah-ubah. Oleh karena itu, dibutuhkan  suatu  sifat  kestabilan dari  metode  numerik  dalam  menyelesaikan Delay Differential Equations (DDE). Menentukan kestabilan dari suatu metode numerik  dalam  menyelesaikan  persamaan  diferensial  merupakan  suatu  hal yang  penting  agar  memperoleh  solusi  yang  diinginkan,  yaitu  nilai  besarnya kesalahan  selalu  stabil  atau  konvergen.  Penelitian  ini  bertujuan  untuk menganalisis  sifat  kestabilan  dari  Runge-Kutta  dengan  Interpolasi  Cubic- Spline  dalam  menentukan  solusi  dari  Delay  Differential  Equation.  Dalam penelitian  ini,  analisa  kestabilan  menggunakan  metode  step-by-step  untuk solusi numerik dari Delay Differential Equation. Penelitian ini dibatasi hanya pada satu persamaan uji Linear Delay Differential Equations yaitu persamaan yang  berbentuk 𝒚 ′ (𝒕) = 𝝀𝒚(𝒕) + 𝝁𝒚(𝒕 − 𝝉) dengan 𝝉 > 𝟎 serta 𝝀 dan 𝝁 adalah  bilangan  kompleks.  Metode  yang  digunakan  oleh  penelitian  ini merupakan kajian literatur dengan menemukan bentuk P-stability dari metode tersebut.  Kestabilan  dari  metode  yang  diajukan  dapat  dilihat  dari  sifat  P- stability  dan  A-  Stability  metode  tersebut.  Dari  hasil  penelitian  diperoleh bahwa  Metode  Runge-Kutta  dengan  interpolasi  Cubic-Spline  memenuhi  sifat- sifat  kestabilan  dari  metode  numerik  dalam  menyelesaikan  Delay  Differential Equations, dimana  P-stability  dan  A-  Stability  dipenuhi  oleh  metode numerik tersebut.  Sehingga  dapat  disimpulkan  bahwa  metode  Runge-Kutta  dengan interpolasi  Cubic-Spline  dapat  digunakan  untuk  menyelesaikan  suatu  Delay Differential Equations (DDE) dalam berbagai kondisi.
Kata kunci: Cubic-Spline; DDE; Runge-Kutta; stabilitas

## Abstract
Algorithms  of  numerical  methods  sometimes  produce  non-converging  results when the inputs are variable. Therefore, a stable numerical method is needed to solve  Delay  Differential  Equations.  Determining  the  stability  of  a  numerical method  in  solving  differential  equations  is  important  to  obtaining  the  desired solution.  The  value  of  the  magnitude  of  the  error  is  always  stable  and convergent.  This  study  aims  to  analyze  the  stability  properties  of  the  Runge- Kutta method with Cubic-Spline Interpolation in determining the solution of the Delay Differential Equation. In this study, the stability analysis uses a step-by- step method for the numerical solution of the Delay Differential Equation. This research  is  limited  to  the  Linear  Delay  Differential  Equations  test  equation  in the form of 𝒚 ′ (𝒕) = 𝝀𝒚(𝒕) + 𝝁𝒚(𝒕 − 𝝉) Where 𝝉 > 𝟎 and 𝝀 and 𝝁 are complex numbers, the method used in this study is a literature review by finding the P- stability form  of  the  method.  The  stability of  the proposed  method  can be  seen from  the  P-stability  and  A-  Stability properties  of  the  method.  From  the  study results,  it  was  found  that  the  Runge-Kutta  method  with  Cubic-Spline interpolation  met  the  stability  properties  of  the  numerical  method  in  solving

## AKSIOMA: Jurnal Matematika dan Pendidikan Matematika
Vol. 13, No. 1 April 2022 e-ISSN 2579-7646
Delay  Differential  Equations,  namely  P-st  A-Stability.  So  it  can  be  concluded that  this  method  can  be  used  to  solve  Delay  Differential  Equations  in  various conditions.
Keywords: Cubic-Spline; DDE; Runge-Kutta; Stability

## A. Pendahuluan
Salah satu persamaan matematika yang paling sering digunakan dalam berbagai  pemodelan  masalah  nyata  adalah Delay  Differential  Equations (DDE).  DDE  memiliki  suatu  bentuk  yang  unik,  dimana   persamaan  ini tidak  hanya  bergantung  kepada  waktu  sekarang  tetapi  juga  bergantung kepada  waktu  sebelumnya.  Beberapa  penelitian  mengenai  metode penyelesaian DDE telah banyak diteliti, diantaranya Oberle & Pesch (1981) mengenalkan  metode  modifikasi  Runge-Kutta  dengan  interpolasi  Hermite untuk  menyelesaikan  DDE,  selanjutnya  Ponnammal  (2017)  menggunakan Singly  Diagonally  Implicit Runge-Kutta  untuk  menyelesaikan  DDE. Ditahun  yang  sama,  Yaghoobi et.  al. (2017)    menggunakan Cubic-Spline untuk mengaproksimasi waktu tunda di DDE. Dimana Cubic-Spline adalah suatu  metode  interpolasi  yang  mengaproksimasi  fungsi  berdasarkan kumpulan  data  atau  titik-titik  koordinat,  dimana  setiap  pasangan  titik diinterpolasi dengan persamaan polinomial berorde 3.
Metode Runge-Kutta dengan interpolasi Cubic-Spline merupakan salah satu  metode  numerik  yang  dapat  digunakan  untuk  menyelesaikan  DDE. Persamaan Runge-Kutta yang digunakan tertuang dalam persamaan (1).
{ 𝑌 𝑛+1 𝑖 = 𝑦 𝑛 + ℎ 𝑛 ∑ 𝑎 𝑖𝑗 𝑘 𝑖 𝑠
𝑗=1 ,             𝑖 = 1, … , 𝑠 𝑘 𝑖 = 𝑓 (𝑡 𝑛+1 𝑗 , 𝑌 𝑛+1 𝑗 , 𝜂 (𝑡 𝑛+1 𝑗 − 𝛼(𝑡 𝑛+1 𝑗 , 𝑌 𝑛+1 𝑗 ))) (1)
dengan interpolasi Cubic-Spline ditampilkan dalam persamaan (2).
𝜂 ′′ (𝑡) = (𝑡 𝑖+1 − 𝑡)𝑀 𝑖 + (𝑡 − 𝑡 𝑖 )𝑀 𝑖+1 ℎ 𝑖 ,        𝑖 = 0, 1, … , 𝑛 − 1             (2)
Dalam  penyelesaian  numerik  terdapat  sifat  yang  harus  dipenuhi  oleh persamaan yaitu well posed , yang artinya jika persamaan tersebut dikenai sedikit  gangguan  maka  akan  mempengaruhi  solusi  yang  diperoleh.  Oleh karena  itu,  metode  yang  baik  untuk  menyelesaikan  suatu  persamaan diferensial harus memenuhi tiga sifat utama berikut: stabil, konvergen, dan konsisten.  Namun  dalam  penelitian  ini,  hanya  berfokus  pada  sifat kestabilan dari metode Runge-Kutta dengan interpolasi Cubic-Spline untuk memperoleh solusi dari Delay Differential Equations (DDE).
Adapun bentuk umum dari masalah nilai awal untuk Delay Differential Equations (DDE) diberikan pada persamaan (3).

## AKSIOMA: Jurnal Matematika dan Pendidikan Matematika
Vol. 13, No. 1 April 2022 e-ISSN 2579-7646
𝑦 ′ (𝑡) = 𝑓(𝑡, 𝑦, 𝑦(𝑡 − 𝛼)),        𝑡 ≥ 𝑡 0 𝑦(𝑡 0 ) = 𝑦 0 𝑦(𝑡) = 𝜑(𝑡),         𝑡 0 − 𝛼 ≤ 𝑡 ≤ 𝑡 0 (3)
dengan 𝑦(𝑡) fungsi bernilai vektor- 𝑛 , 𝛼 > 0 adalah konstanta delay , dan 𝜑(𝑡) adalah  fungsi  nilai  awal,  yang  diasumsikan  kontinu  terus  menerus ( piecewise  continuous )  pada  interval 𝑡 0 − 𝛼 ≤ 𝑡 ≤ 𝑡 0 ,  untuk  kasus 𝜑(𝑡) ≠ 𝑡 0 juga termasuk.
Solusi dari persamaan (1) adalah suatu fungsi diferensial piecewise 𝑦(𝑡) yang  terdefinisi  pada  interval 𝑡 0 − 𝛼 ≤ 𝑡 ≤ 𝑡 𝑓 (𝑡 𝑓 > 𝑡 0 ) yang  kontinu  untuk semua 𝑡 > 𝑡 0 dan memenuhi persamaan diferensial. DDE hanya memenuhi limit  kiri-kanan  titik 𝑡 0 + 𝛼 dan  titik 𝜉 + 𝛼 ,  dimana 𝜉 adalah  suatu  titik lompatan  pada 𝜑 .  Namun,  kemulusan  dari  grafik 𝑓 dan 𝜑 tidak mengakibatkan solusinya juga mulus, walaupun 𝜑(𝑡 0 ) = 𝑦 0 dipenuhi.
Definisi 1. DDE (3) dikatakan stabil jika untuk setiap solusi fungsi awal yang  kecil 𝑦(𝑡) mendekati  nol  ketika 𝑡 → ∞ ,  sehingga  untuk  semua konstanta yang sangat kecil 𝛿 > 0 , 𝑡 → ∞ untuk ‖𝑢(𝑡)‖ < 𝛿, −𝑑 ≤ 𝑡 ≤ 𝑡 0 .(Al- Mutib, 1984)
Definisi  2. Suatu  metode  numerik  dikatakan stabil jika  untuk  setiap solusi 𝑦(𝑡) memenuhi 𝑦(𝑡) = 0 jika 𝑡 → ∞ . (Al-Mutib, 1984)
Suatu  sistem  DDE 𝑥 ′ (𝑡) = 𝑓(𝑡, 𝑥 𝑡 ) yang  memenuhi 𝑓(𝑡, 0) = 0, 𝑡 ∈ ℝ sehingga 𝑥(𝑡) = 0 adalah  solusinya.  Solusi 𝑥 = 0 dikatakan stabil jika untuk  setiap 𝜎 ∈ ℝ dan 𝜀 > 0 ,  terdapat 𝛿 = (𝜎, 𝜀) > 0 sehingga 𝜙 ∈ 𝐶 dan ‖𝜙‖ < 𝛿 mengakibatkan ‖𝑥 𝑡 (𝜎, 𝜙)‖ < 𝜀, 𝑡 ≥ 𝜎 dengan 𝑥(𝑡, 𝜎, 𝜙) → 0, 𝑡 → ∞ .
Beberapa  penelitian  telah  dilakukan  untuk  melihat  kestabilan  dari metode  penyelesaian  DDE  diantaranya  Metode  Kutta-Merson,  Metode Trapesium,  dan  Metode  Implisit  Runge-Kutta  untuk  menyelesaikan  DDE telah  diperiksa  kestabilannya  oleh  Al-Mutib  (1984).  Bellen  (2000)  telah menyelesaikan  masalah  stabilitas  metode  Runge-Kutta  untuk  DDE. Kemudian  Liu et.  al. (2006)  telah  memeriksa  stabilitas  dari  metode modifikasi  Runge-Kutta  untuk  DDE,  dan  menghasilkan  bahwa  metode modifikasi Runge-Kutta lebih akurat dan stabil untuk menyelesaikan DDE. Penelitian  ini  akan  melihat  sifat-sifat  dari  metode  numerik  Runge-Kutta yang  dimodifikasi  dengan  interpolasi Cubic-Spline yang  digunakan  untuk menentukan  solusi  DDE.  Sifat  yang  diteliti  adalah  sifat  kestabilan  dari metode yang diperkenalkan.

## B. Metode Penelitian
Penelitian  ini  merupakan  penelitian  kajian  literatur,  untuk menganalisa kestabilan dari metode yang diperkenalkan digunakan metode step-by-step untuk metode numerik dalam menentukan solusi numerik dari Delay  Differential  Equation (DDE).  Metode  ini  akan  memperlihatkan

## AKSIOMA: Jurnal Matematika dan Pendidikan Matematika
Vol. 13, No. 1 April 2022 e-ISSN 2579-7646
apakah  metode  Runge-Kutta  yang  diperkenalkan  memenuhi  atau menunjukkan sifat kestabilan seperti dalam Definisi 1 dan Definisi 2 untuk menyelesaikan DDE, dimana waktu tunda dalam DDE akan diaproksimasi dengan menggunakan interpolasi Cubic Spline .
Untuk  memperoleh  hasil  penelitian,  dilakukan  penurunan-penurunan dari  teorema-teorema  yang  mendukung  sifat  kestabilan  dari  suatu persamaan diferensial, dimana teorema kestabilan yang digunakan adalah teorema P-Stability dan A-Stability yang  juga  berlaku  untuk  Persamaan Diferensial Biasa.

## C. Hasil dan Pembahasan
Menurut Barwell (1975), kestabilan dari metode numerik yang diajukan dapat dilihat dari apakah metode tersebut memenuhi sifat P-Stability .
Definisi  3. Daerah P-stability dari  metode  numerik  untuk  DDE  adalah himpunan 𝑆 𝑃 dari  pasangan  bilangan  kompleks (𝛼, 𝛽) , 𝛼 = ℎ𝜆, 𝛽 = ℎ𝜇 , sehingga solusi numerik {𝑦 𝑛 } 𝑛≥0 dari (3), dapat ditulis dalam bentuk seperti yang ditampilan dalam persamaan (4).
ℎ = 𝜏 𝑚 ⁄ ,        𝑚 ≥ 1,    𝑚 bilangan bulat                   (4)
yang  memenuhi lim 𝑛→∞ 𝑦 𝑛 = 0 untuk  semua  konstanta delay dan  nilai  awal 𝜑(𝑡) .
Salah satu cara untuk menganalisa stabilitas dari solusi numerik yang paling banyak digunakan adalah dengan menggunakan persamaan (5) yang dikenal sebagai persamaan uji
𝑦 ′ (𝑡) = 𝜆𝑦(𝑡) + 𝜇𝑦(𝑡 − 𝜏),      𝑡 ≥ 𝑡 0 ,     𝑦(𝑡) = 𝑦 0 (𝑡),     𝑡 ≤ 0            (5) dengan 𝜏 > 0 adalah delay ,  dan 𝜆, 𝜇 melambangkan  konstanta  yang merupakan bilangan kompleks, dengan 𝑦 0 adalah fungsi bernilai kompleks dan  kontinu.  Dari  Barwell  (1975),  solusi 𝑦(𝑡) dari  persamaan  (5)  akan menuju ke nol dengan 𝑡 → ∞ , untuk semua 𝑦 0 , jika 𝑅𝑒(𝜆) < −|𝜇|
Aplikasikan persamaan (1) ke persamaan (5), diperoleh persamaan (6).
𝑌 𝑛+1 𝑖 = 𝑦 𝑛 + ℎ ∑ 𝑎 𝑖𝑗 (𝜆𝑌 𝑛+1 𝑗 + 𝜇𝜂(𝑡 𝑛−𝑚+1 𝑗 )) 𝑠 𝑗=1 ,      𝑖 = 1, … , 𝑠             (6)
dengan
𝜂(𝑡) = (𝑡 𝑖+1 − 𝑡) 3 𝑀 𝑖 + (𝑡 − 𝑡 𝑖 ) 3 𝑀 𝑖+1 6ℎ 𝑖 + (𝑡 𝑖+1 − 𝑡)𝑓 𝑖 + (𝑡 − 𝑡 𝑖 )𝑓 𝑖+1 ℎ 𝑖 − ℎ 𝑖 6 [(𝑡 𝑖+1 − 𝑡)𝑀 𝑖 + (𝑡 − 𝑡 𝑖 )𝑀 𝑖+𝑡 ]
jika 𝜂(𝑡) disubtitusikan  pada 𝑌 𝑛+1 𝑖 di  persamaan  (6)  maka  diperoleh persamaan (7).

## AKSIOMA: Jurnal Matematika dan Pendidikan Matematika
Vol. 13, No. 1 April 2022 e-ISSN 2579-7646
𝑌 𝑛+1 𝑖 = 𝑦 𝑛 + ℎ ∑ 𝑎 𝑖𝑗 (𝜆𝑌 𝑛+1 𝑗 𝑠 𝑗=1
+ 𝜇 ( (𝑡 𝑖+1 − 𝑡) 3 𝑀 𝑖 + (𝑡 − 𝑡 𝑖 ) 3 𝑀 𝑖+1 6ℎ 𝑖 + (𝑡 𝑖+1 − 𝑡)𝑓 𝑖 + (𝑡 − 𝑡 𝑖 )𝑓 𝑖+1 ℎ 𝑖 − ℎ 𝑖 6 [(𝑡 𝑖+1 − 𝑡)𝑀 𝑖 + (𝑡 − 𝑡 𝑖 )𝑀 𝑖+𝑡 ]))                                                         (7)
𝜂(𝑡) pada persamaan (6) merupakan bagian dari 𝑘 𝑖 pada persamaan Runge- Kutta atau dapat ditulis sebagai berikut:
𝑘 𝑖 = 𝑓 (𝑡 𝑛+1 𝑗 , 𝑌 𝑛+1 𝑗 , 𝜂 (𝑡 𝑛+1 𝑗 − 𝛼(𝑡 𝑛+1 𝑗 , 𝑌 𝑛+1 𝑗 ))) 𝑘 𝑖 = 𝑓 (𝑡 𝑛+1 𝑗 , 𝑌 𝑛+1 𝑗 , (𝑡 𝑖+1 − 𝑡) 3 𝑀 𝑖 + (𝑡 − 𝑡 𝑖 ) 3 𝑀 𝑖+1 6ℎ 𝑖 + (𝑡 𝑖+1 − 𝑡)𝑓 𝑖 + (𝑡 − 𝑡 𝑖 )𝑓 𝑖+1 ℎ 𝑖 − ℎ 𝑖 6 [(𝑡 𝑖+1 − 𝑡)𝑀 𝑖 + (𝑡 − 𝑡 𝑖 )𝑀 𝑖+𝑡 ])                                                           (8)
Selanjutnya,  misal  didefinisikan 𝑢 = (1, … , 1) 𝑇 dan  untuk 𝑛 ≥ 1 , diperoleh
𝑘 𝑛 = (𝑘 𝑛 (1) , … , 𝑘 𝑛 (𝑞) ) 𝑇 𝑎 = (𝑎 1 , 𝑎 2 , … , 𝑎 𝑞 ) 𝑇
𝑡 − 𝑡 𝑛 = (𝑡 − 𝑡 1 , 𝑡 − 𝑡 2 , … , 𝑡 − 𝑡 𝑞 ) 𝑇
𝑡 𝑛+1 − 𝑡 = (𝑡 2 − 𝑡, 𝑡 3 − 𝑡, 𝑡 𝑞+1 − 𝑡) 𝑇
ℎ = (ℎ 1 , ℎ 2 , … , ℎ 𝑞 ) 𝑇
𝑀 = (𝑀 1 , 𝑀 2 , … , 𝑀 𝑞 ) 𝑇
𝐴 = [𝑎 𝑖𝑗 ] 𝑖,𝑗=1 𝑠
untuk 𝑛 ≥ 𝑚 ,  persamaan  (7)  menjadi  persamaan  (9)  dan  persamaan  (8) menjadi persamaan (10).
𝑌 𝑛+1 = 𝑦 𝑛 + ℎ𝐴𝑘 𝑛+1 (9)
𝑘 𝑛+1 = 𝜆(𝑦 𝑛 𝑢 + ℎ𝐴𝑘 𝑛+1 ) + 𝜇 ( (𝑡 𝑖+1 − 𝑡) 3 𝑀 𝑖 + (𝑡 − 𝑡 𝑖 ) 3 𝑀 𝑖+1 6ℎ 𝑖 + (𝑡 𝑖+1 − 𝑡)𝑓 𝑖 + (𝑡 − 𝑡 𝑖 )𝑓 𝑖+1 ℎ 𝑖 − ℎ 𝑖 6 [(𝑡 𝑖+1 − 𝑡)𝑀 𝑖 + (𝑡 − 𝑡 𝑖 )𝑀 𝑖+𝑡 ])                                                                                        (10)

## AKSIOMA: Jurnal Matematika dan Pendidikan Matematika
Vol. 13, No. 1 April 2022 e-ISSN 2579-7646
Dalam Bellen & Zennaro (2003), suatu fungsi rasional dikatakan fungsi A-stability jika  dapat  ditulis  dalam  bentuk  ke  dalam  bentuk  persamaan (11).
𝑅(𝛼) = 1 + 𝛼𝑎 𝑇 (𝐼 − 𝛼𝐴) −1 𝑢                                         (11) atau bentuk ekuivalennya ditampilkan dalam persamaan (12).
𝑅(𝛼) = det[𝐼 − 𝛼𝐴 + 𝛼𝑢𝑎 𝑇 ] det[𝐼 − 𝛼𝐴] (12)
Agar diperoleh bentuk seperti persamaan (11), maka manipulasi bentuk persamaan (10) menjadi persamaan (13) berikut.
𝑘 𝑛+1 = 𝜆𝑦 𝑛 𝑢 + 𝜆ℎ𝐴𝑘 𝑛+1 + 𝜇 ( (𝑡 𝑖+1 − 𝑡) 3 𝑀 𝑖 + (𝑡 − 𝑡 𝑖 ) 3 𝑀 𝑖+1 6ℎ 𝑖 + (𝑡 𝑖+1 − 𝑡)𝑓 𝑖 + (𝑡 − 𝑡 𝑖 )𝑓 𝑖+1 ℎ 𝑖 − ℎ 𝑖 6 [(𝑡 𝑖+1 − 𝑡)𝑀 𝑖 + (𝑡 − 𝑡 𝑖 )𝑀 𝑖+𝑡 ])
𝑘 𝑛+1 − 𝜆ℎ𝐴𝑘 𝑛+1 = 𝜆𝑦 𝑛 𝑢 + 𝜇 ( (𝑡 𝑖+1 − 𝑡) 3 𝑀 𝑖 + (𝑡 − 𝑡 𝑖 ) 3 𝑀 𝑖+1 6ℎ 𝑖 + (𝑡 𝑖+1 − 𝑡)𝑓 𝑖 + (𝑡 − 𝑡 𝑖 )𝑓 𝑖+1 ℎ 𝑖 − ℎ 𝑖 6 [(𝑡 𝑖+1 − 𝑡)𝑀 𝑖 + (𝑡 − 𝑡 𝑖 )𝑀 𝑖+𝑡 ])
(𝐼 − ℎ𝜆𝐴)𝑘 𝑛+1 = 𝜆𝑦 𝑛 𝑢 + 𝜇 ( (𝑡 𝑖+1 − 𝑡) 3 𝑀 𝑖 + (𝑡 − 𝑡 𝑖 ) 3 𝑀 𝑖+1 6ℎ 𝑖 + (𝑡 𝑖+1 − 𝑡)𝑓 𝑖 + (𝑡 − 𝑡 𝑖 )𝑓 𝑖+1 ℎ 𝑖 − ℎ 𝑖 6 [(𝑡 𝑖+1 − 𝑡)𝑀 𝑖 + (𝑡 − 𝑡 𝑖 )𝑀 𝑖+𝑡 ])                                                         (13)
Persamaan (11) dikali dengan ℎ(𝐼 − ℎ𝜆𝐴) −1 maka diperoleh persamaan (14) berikut.
ℎ𝑘 𝑛+1 = ℎ𝜆(𝐼 − ℎ𝜆𝐴) −1 𝑦 𝑛 𝑢
+ ℎ𝜇(𝐼 − ℎ𝜆𝐴) −1 ( (𝑡 𝑖+1 − 𝑡) 3 𝑀 𝑖 + (𝑡 − 𝑡 𝑖 ) 3 𝑀 𝑖+1 6ℎ 𝑖 + (𝑡 𝑖+1 − 𝑡)𝑓 𝑖 + (𝑡 − 𝑡 𝑖 )𝑓 𝑖+1 ℎ 𝑖 − ℎ 𝑖 6 [(𝑡 𝑖+1 − 𝑡)𝑀 𝑖 + (𝑡 − 𝑡 𝑖 )𝑀 𝑖+𝑡 ])   (14)
Ganti  persamaan  (14)  dengan  memisalkan 𝛼 = ℎ𝜆 , 𝛽 = ℎ𝜇 ,  dan 𝜔 = (𝐼 − ℎ𝜆𝐴) −1 , sehingga diperoleh persamaan (15).

## AKSIOMA: Jurnal Matematika dan Pendidikan Matematika
Vol. 13, No. 1 April 2022 e-ISSN 2579-7646
ℎ𝑘 𝑛+1 = 𝛼𝜔𝑦 𝑛 𝑢 + 𝛽𝜔 ( (𝑡 𝑖+1 − 𝑡) 3 𝑀 𝑖 + (𝑡 − 𝑡 𝑖 ) 3 𝑀 𝑖+1 6ℎ 𝑖 + (𝑡 𝑖+1 − 𝑡)𝑓 𝑖 + (𝑡 − 𝑡 𝑖 )𝑓 𝑖+1 ℎ 𝑖 − ℎ 𝑖 6 [(𝑡 𝑖+1 − 𝑡)𝑀 𝑖 + (𝑡 − 𝑡 𝑖 )𝑀 𝑖+𝑡 ])                                                         (15)
selanjutnya  subtitusi  persamaan  (15)  ke  persamaan  (9)  sehingga  didapat persamaan (16).
𝑦 𝑛+1 = 𝑦 𝑛 + ℎ𝑎 𝑇 𝑘 𝑛+1 𝑦 𝑛+1 = 𝑦 𝑛 + 𝑎 𝑇 ℎ𝑘 𝑛+1 = 𝑦 𝑛 + 𝑎 𝑇 (𝛼𝜔𝑦 𝑛 𝑢 + 𝛽𝜔 ( (𝑡 𝑖+1 − 𝑡) 3 𝑀 𝑖 + (𝑡 − 𝑡 𝑖 ) 3 𝑀 𝑖+1 6ℎ 𝑖 + (𝑡 𝑖+1 − 𝑡)𝑓 𝑖 + (𝑡 − 𝑡 𝑖 )𝑓 𝑖+1 ℎ 𝑖 − ℎ 𝑖 6 [(𝑡 𝑖+1 − 𝑡)𝑀 𝑖 + (𝑡 − 𝑡 𝑖 )𝑀 𝑖+𝑡 ])) 𝑦 𝑛+1 = 𝑦 𝑛 + 𝛼𝑎 𝑇 𝜔𝑦 𝑛 𝑢 + 𝛽𝑎 𝑇 𝜔 ( (𝑡 𝑖+1 − 𝑡) 3 𝑀 𝑖 + (𝑡 − 𝑡 𝑖 ) 3 𝑀 𝑖+1 6ℎ 𝑖 + (𝑡 𝑖+1 − 𝑡)𝑓 𝑖 + (𝑡 − 𝑡 𝑖 )𝑓 𝑖+1 ℎ 𝑖 − ℎ 𝑖 6 [(𝑡 𝑖+1 − 𝑡)𝑀 𝑖 + (𝑡 − 𝑡 𝑖 )𝑀 𝑖+𝑡 ])
𝑦 𝑛+1 = 𝑦 𝑛 (1 + 𝛼𝑎 𝑇 𝜔𝑢)
+ 𝛽𝑎 𝑇 𝜔 ( (𝑡 𝑖+1 − 𝑡) 3 𝑀 𝑖 + (𝑡 − 𝑡 𝑖 )
3 𝑀 𝑖+1 6ℎ 𝑖 + (𝑡 𝑖+1 − 𝑡)𝑓 𝑖 + (𝑡 − 𝑡 𝑖 )𝑓 𝑖+1 ℎ 𝑖 − ℎ 𝑖 6 [(𝑡 𝑖+1 − 𝑡)𝑀 𝑖 + (𝑡 − 𝑡 𝑖 )𝑀 𝑖+𝑡 ])                                                         (16)
Berdasarkan persamaan (14), maka fungsi A-stability diperoleh
𝑅(𝛼) = 1 + 𝛼𝑎 𝑇 𝜔𝑢
Untuk melihat P-Stability , tulis ulang kembali persamaan (16) kedalam bentuk matriks, sehingga diperoleh
𝑌 𝑛+1 = ( 𝑦 𝑛+1 ℎ𝑘 𝑛+1 )
Selanjutnya, persamaan (15) dapat ditulis

## AKSIOMA: Jurnal Matematika dan Pendidikan Matematika
Vol. 13, No. 1 April 2022 e-ISSN 2579-7646
( 1 0 0 𝐼 ) 𝑌 𝑛+1 = (1 + 𝛼𝑎 𝑇 𝜔𝑢 0 𝛼𝜔𝑢 𝑢 ) 𝑌 𝑛 + ( 𝛽𝑎 𝑇 𝜔((𝑡 𝑛+1 − 𝑡) 𝑇 ) 3 𝑀 𝑇 + (𝑡 − 𝑡 𝑛 ) 3 𝑀 𝑇 6ℎ 𝑇 0 𝛽𝜔((𝑡 𝑛+1 − 𝑡) 𝑇 ) 3 𝑀 𝑇 + (𝑡 − 𝑡 𝑛 ) 3 𝑀 𝑇 6ℎ 𝑇 0) + ( 𝛽𝑎 𝑇 𝜔 (𝑡 𝑛+1 − 𝑡) 𝑇 ℎ 𝑇 0 𝛽𝜔 (𝑡 𝑛+1 − 𝑡) 𝑇 ℎ 𝑇 0) 𝑌 𝑛 + ( 𝛽𝑎 𝑇 𝜔 (𝑡 𝑛+1 − 𝑡) 𝑇 ℎ 𝑇 0 𝛽𝜔 (𝑡 𝑛+1 − 𝑡) 𝑇 ℎ 𝑇 0) 𝑌 𝑛−1 − ( 𝛽𝑎 𝑇 𝜔 ℎ 𝑇 6 [(𝑡 𝑛+1 − 𝑡) 𝑇 𝑀 𝑇 + (𝑡 − 𝑡 𝑛 ) 𝑇 𝑀 𝑇 ] 0 𝛽𝜔 ℎ 𝑇 6 [(𝑡 𝑛+1 − 𝑡) 𝑇 𝑀 𝑇 + (𝑡 − 𝑡 𝑛 ) 𝑇 𝑀 𝑇 ] 0 )
Kemudian, persamaankarakteristik dari persamaan (5)
𝑆(𝛼, 𝛽, 𝜁) = det
( ( 1 0 0 𝐼 ) 𝜁 𝑛+1 − (1 + 𝛼𝑎 𝑇 𝜔𝑢 0 𝛼𝜔𝑢 𝑢 ) 𝜁 𝑛 − ( 𝛽𝑎 𝑇 𝜔((𝑡 𝑛+1 − 𝑡) 𝑇 ) 3 𝑀 𝑇 + (𝑡 − 𝑡 𝑛 ) 3 𝑀 𝑇 6ℎ 𝑇 0 𝛽𝜔((𝑡 𝑛+1 − 𝑡) 𝑇 ) 3 𝑀 𝑇 + (𝑡 − 𝑡 𝑛 ) 3 𝑀 𝑇 6ℎ 𝑇 𝑢) − ( 𝛽𝑎 𝑇 𝜔 (𝑡 𝑛+1 − 𝑡) 𝑇 ℎ 𝑇 0 𝛽𝜔 (𝑡 𝑛+1 − 𝑡) 𝑇 ℎ 𝑇 0) 𝜁 𝑛 − ( 𝛽𝑎 𝑇 𝜔 (𝑡 𝑛+1 − 𝑡) 𝑇 ℎ 𝑇 0 𝛽𝜔 (𝑡 𝑛+1 − 𝑡) 𝑇 ℎ 𝑇 0) 𝜁 𝑛−1 + ( 𝛽𝑎 𝑇 𝜔 ℎ 𝑇 6 [(𝑡 𝑛+1 − 𝑡) 𝑇 𝑀 𝑇 + (𝑡 − 𝑡 𝑛 ) 𝑇 𝑀 𝑇 ] 0 𝛽𝜔 ℎ 𝑇 6 [(𝑡 𝑛+1 − 𝑡) 𝑇 𝑀 𝑇 + (𝑡 − 𝑡 𝑛 ) 𝑇 𝑀 𝑇 ] 0) ) dengan 𝐼 adalah  matriks  identitas,  maka  relasi  rekurensinya  ditampilkan dalam persamaan (17).
𝑆(𝛼, 𝛽, 𝜁) = 𝜍 𝑛+1 − (1 + 𝛼𝑎 𝑇 𝜔𝑢)𝜍 𝑛 − (𝛽𝑎 𝑇 𝜔((𝑡 𝑛+1 − 𝑡) 𝑇 ) 3 𝑀 𝑇 + 𝛽𝑎 𝑇 𝜔((𝑡 − 𝑡 𝑛 ) 𝑇 ) 3 𝑀 𝑇 ) − 𝛽𝑎 𝑇 𝜔 ( (𝑡 𝑛+1 − 𝑡) 𝑇 ℎ 𝑇 ) 𝜍 𝑛 − 𝛽𝑎 𝑇 𝜔 ( (𝑡 − 𝑡 𝑛 ) 𝑇 ℎ 𝑇 ) 𝜁 𝑛−1 + 𝛽𝑎 𝑇 𝜔 ℎ 𝑇 6 [(𝑡 𝑛+1 − 𝑡) 𝑇 𝑀 𝑇 + (𝑡 − 𝑡 𝑛 ) 𝑇 𝑀 𝑇 ]                                          (17)
Sehingga  dapat  dilihat  bahwa  persamaan  (17)  memenuhi  bentuk  dari fungsi rasional (11), oleh karena itu, persamaan polinomial (17) dikatakan P-Stability .

## AKSIOMA: Jurnal Matematika dan Pendidikan Matematika
Vol. 13, No. 1 April 2022 e-ISSN 2579-7646

## D. Simpulan
Dari hasil analisis kestabilan metode Runge-Kutta dengan Cubic-Spline untuk menyelesaikan DDE diperoleh bahwa metode ini memenuhi sifat 𝑃 − 𝑠𝑡𝑎𝑏𝑖𝑙𝑖𝑡𝑦 ,  sehingga  dapat  digunakan  untuk  menyelesaikan  DDE  dengan nilai  kesalahan  ( error) stabil  atau  solusi  yang  diperoleh  akan  konvergen. Sehingga solusi yang diperoleh dapat digunakan.

## E. Daftar Pustaka
Al-Mutib, A. N. (1984). Stability properties of numerical methods for solving delay  differential  equations. Journal  of  Computational  and  Applied Mathematics , 10 (1), 71 – 79. doi:10.1016/0377-0427(84)90071-2. Barwell,  V.  K.  (1975).  Special  functional  stability  problems  for  equations.
BIT , 15 , 2 – 7.
Bellen,  A.  (2000).  Numerical  methods  for  delay  differential  equations:
Accuracy  and  stability  problems. IFAC  Proceedings  Volumes , 33 (23),
127 – 128. doi:10.1016/s1474-6670(17)36928-8. Bellen,  A.  &  Zennaro,  M.  (2003). Numerical Methods  for  Delay  Differential Equations . Oxford: Clarendon Press.
Liu,  M.  Z.,  Yang,  Z.  W.,  &  Xu,  Y.  (2006).  The  stability  of  modified  Runge- Kutta  methods  for  the  pantograph  equation. Mathematics  of Computation , 75 (255), 1201 – 1216. doi: 10.1090/s0025-5718-06-01844-8. Oberle,  H.  J.  &  Pesch,  H.  J.  (1981).  Numerical  treatment  of  delay differential  equations  by  Hermite  Interpolation. Numerische Mathematik , 37 (2), 235 – 255. doi:10.1007/BF01398255. Ponnammal, K. (2017). Singly Diagonally Implicit Runge-Kutta fourth and fifth order methods for Solving Delay Differential Equations . 6 (7), 182 – 197.
Yaghoobi, S., Moghaddam, B. P., & Ivaz, K. (2017). An efficient cubic spline approximation  for  variable-order  fractional  differential  equations  with time  delay. Nonlinear  Dynamics , 87 (2),  815 – 826.  doi:10.1007/s11071- 016-3079-4.