Last Updated on Agustus 18, 2022 by prooffic
Postingan kali ini akan membahas tentang bagaimana melakukan Simulasi Numerik Persamaan Diferensial Biasa Orde 2 dengan menggunakan bahasa pemrograman Python. Pembahasan berikut akan menggunakan contoh-contoh tertentu, meskipun begitu, cara yang digunakan dapat diterapkan ke contoh-contoh lainnya. Pembahasan berikut juga akan dibagi menjadi dua, yaitu untuk persamaan diferensial orde 2 yang dilengkapi dengan syarat awal, sering disebut sebagai Masalah Nilai Awal (MNA), dan persamaan diferensial orde 2 yang dilengkapi dengan syarat batas, sering disebut sebagai Masalah Syarat Batas (MSB).
**Selamat menikmati**
Pembahasan untuk tiap bagian terdiri dari
- Penyelesaian solusi analitik
- Simulasi solusi numerik dan solusi analitik
Persamaan Diferensial Biasa Orde 2 yang dilengkapi dengan syarat awal (Masalah Nilai Awal)
Diberikan persamaan diferensial orde 2 berikut. $$\begin{aligned} y^{\prime\prime} + 2 y^{\prime}- 3 y = \sin t\end{aligned} $$dengan syarat awal $y(0) = 0$ dan $y^{\prime} (0) = 1$.
Solusi Analitik
Terlebih dahulu akan ditentukan solusi analitik dari persamaan diferensial tersebut. Ingat kembali bahwa solusi umum dari suatu persamaan diferensial biasa berbentuk $$y(t) = y_h (t) + y_p (t) $$dengan $y_h (t) $ dan $y_p (t)$ masing- masing adalah solusi homogen dan solusi partikularnya. Untuk menentukan solusi homogennya, tinjau persamaan homogen berikut. $$ y^{\prime\prime} + 2 y^{\prime}- 3 y = 0 $$
Karena koefisien dari $y^{\prime\prime}, y^{\prime}, y$ tidak memuat $t,$ maka kita dapat menggunakan metode karakteristik untuk menyelesaikannya. Persamaan karakteristiknya adalah $$r^2 + 2r – 3 = 0 $$Solusinya adalah $r = -3$ serta $r = 1.$ Dari sini, solusi umumnya adalah $$y_h (t) = A e^{-3t} + B e^{t} $$dengan $A$ dan $B$ adalah sebarang konstanta.
Untuk menentukan solusi partikularnya, terlebih dahulu perhatikan bahwa kita dapat melihat dengan jelas fungsi $\sin t$ bukan merupakan solusi homogen. Dari sini, solusi partikularnya berbentuk $$y_p (t) = C \sin t + D \cos t $$dengan $C$ dan $D$ adalah konstanta yang akan ditentukan kemudian. Dengan bentuk $y_p$ tersebut, maka $$y_p^{\prime} = C \cos t – D \sin t $$dan $$y_p^{\prime\prime} = -C \sin t – D \cos t $$Subtitusikan $y_p^{\prime\prime}, y_p^{\prime}, y_p$ ke persamaan, maka diperoleh $$\left(-C \sin t – D \cos t\right) + 2 \left( C \cos t – D \sin t \right) – 3 \left( C \sin t + D \cos t \right) = \sin t$$
Oleh karena itu, dengan meninjau koefisien dari masing-masing $\cos t$ dan $\sin t,$ maka diperoleh $C = 2 D$ dan $2D + 4C = -1.$ Sehingga, $D= -1/10$ dan $C = -2/10.$ Dari sini, $$y_p (t) = -\frac{2}{10} \sin t – \frac{1}{10} \cos t $$Kemudian, $$y_(t) = A e^{-3t} + Be^{t} -\frac{2}{10} \sin t – \frac{1}{10} \cos t $$
Langkah berikutnya adalah menentukan konstanta $A$ dan $B$ menggunakan syarat awal yang diberikan sebelumnya. Untuk $y(0) = 0$ diperoleh bahwa $$0 = A + B – \frac{1}{10} $$
Sebelumnya, perhatikan bahwa $$y^{\prime} (t) = -3 A e^{-3t} + B e^{t} – \frac{2}{10} \cos t + \frac{1}{10} \sin t$$Maka, untuk $y^{\prime} (0) = 1$ diperoleh $$1 = -3 A + B – \frac{2}{10} $$Kemudian, dengan menyelesaikan sistem persamaan linear dalam $A$ dan $B$ tersebut, diperoleh $A = -\frac{11}{40}$ dan $B = \frac{15}{40}.$ Oleh karena itu, $$y_(t) = -\frac{11}{40} e^{-3t} +\frac{15}{40} e^{t} -\frac{2}{10} \sin t – \frac{1}{10} \cos t$$
Simulasi solusi numerik dan solusi analitik dengan python
Simulasi akan dilakukan dengan menggunakan bahasa pemrograman Python. Lybrary yang digunakan adalah numpy, matplotlib dan scipy. Numpy digunakan untuk mendefinisikan array, dalam hal ini adalah interval fungsi kita dan juga untuk memunculkan fungsi trogonometri dan eksponensial. Maptplotlib digunakan untuk melakukan plot (simulasi) solusi analitik dan solusi numerik. Library scipy.integrate terutama odeint digunakan untuk menyelesaikan masalah nilai awal secara numerik. Kode python berikut telah dilengkapi dengan catatan untuk mempermudah dalam memahami makna tiap baris kode.
import library from matplotlib
import pyplot as plt from scipy.integrate
import odeint import numpy as np
#mendefinisikan persamaan
def f(u,t):
return (u[1], -2*u[1]+3*u[0]+np.sin(t))
#mendefinisikan syarat batas
y0 = [0,1]
#mendefinisikan interval
xs = np.linspace(0,10,200)
#mendefinisikan solusi analitik
y_analitik = (-11/40)*np.exp(-3*xs)+(15/40)*np.exp(xs)-(2/10)*np.sin(xs)-(1/10)*np.cos(xs)
#penyelesaian dengan metode numerik
us = odeint(f,y0,xs)
ys = us[:,0]
# melakukan plot (simulasi) dari solusi analitik dan solusi numerik
plt.xlabel('t')
plt.ylabel('y')
plt.title('$y^{\prime\prime} + 2 y^{\prime}- 3 y = \sin t; y(0) = 0, y^{\prime} (0) =1 $')
plt.plot(xs,ys, 'r*',label='numerik')
plt.plot(xs,y_analitik,'k',label=r'analitik')
plt.legend(loc=0)
plt.show()
Jika dijalankan, maka akan diperoleh plot sebagaimana pada gambar berikut.
Dari gambar tersebut, kita dapat mengatakan bahwa solusi numerik dan solusi analitiknya berimpit. Metode numerik yang digunakan pada library scipy.integrate, utamanya odeint, adalah metode runge kutta.
Persamaan Diferensial Biasa Orde 2 yang dilengkapi dengan Kondisi Batas (Masalah Syarat Batas)
Diberikan persamaan diferensial orde 2 berikut. $$\begin{aligned} y^{\prime\prime} + y = \sin (x)\end{aligned} $$dengan syarat batas $y(0) = 0$ dan $y (\pi/2) = 1$.
Sama seperti cara sebelumnya, tulis $$y(x) = y_h (x) + y_p (t) $$Tinjau persamaan homogen $$y^{\prime\prime} + y =0 $$Solusi umumnya adalah $$y_h (x) = A \sin (x) + B \cos (x) $$dengan $A$ dan $B$ adalah sebarang konstanta. Kemudian, kita dapat melibat bahwa $\sin (x)$ adalah termasuk ke dalam solusi dari persamaan homogennya. Maka, persamaan partikularnya berbentuk $$y_p (x) = x (C \sin (x) + D \cos (x)) $$dengan $C$ dan $D$ adalah konstanta yang nilainya akan ditentukan. Dari bentuk $y_p$ tersebut, maka $$y_p^{\prime} = (C \sin (x) + D \cos (x)) + x (C \cos (x) – D \sin (x)) $$dan $$y_p^{\prime\prime} = C \cos (x) – D \sin (x) + C \cos (x) – D \sin (x) + x (-C \sin (x) – D \cos (x))$$
Dengan melakukan subtitusi ke persamaan diferensial awal, maka diperoleh bahwa $$2 C \cos (x) – 2 D \sin (x) = \sin (x) $$Dari sini, $C = 0$ dan $D = -1/2. $Oleh karena itu, $$y_p (x) = – \frac{1}{2} x \cos (x) $$Sehingga, $$y(x) = A \sin (x) + B \cos (x) – \frac{1}{2} x \cos (x) $$
Langkah berikutnya adalah menentukan $A$ dan $B$ berdasarkan syarat batas yang diberikan. Untuk $y(0) = 0, $diperoleh $ B = 0 .$ Selain itu, untuk $y(\pi/2) = 1,$ diperoleh $A = 1 . $Dari sini, $$y(x) = \sin (x) – \frac{1}{2} x \cos (x) $$
Simulasi solusi numerik dan solusi analitik dengan python
Lybrary yang digunakan adalah numpy, matplotlib dan scipy.integrate, terutama solve_ivp dan scipy.optimize, terutama fsolve dan math. Numpy digunakan untuk mendefinisikan array, dalam hal ini adalah interval fungsi kita dan juga untuk memunculkan fungsi trogonometri. Maptplotlib digunakan untuk melakukan plot (simulasi) solusi analitik dan solusi numerik. Library scipy.integrate terutama solve_ivp digunakan untuk menyelesaikan masalah syarat batas secara numerik. Library scipy.optimize, terutama fsolve, digunakan untuk memberikan solusi numerik bagi nilai $y^{\prime} (0).$ Library math digunakan menuliskan $\pi /2.$ Kode python berikut telah dilengkapi dengan catatan untuk mempermudah dalam memahami makna tiap baris kode.
import numpy as np
import matplotlib.pyplot as plt from scipy.integrate
import solve_ivp plt.style.use('seaborn-poster')
import math from scipy.optimize
import fsolve a = math.pi/2
#menuliskan/mendefinisikan persamaan
F = lambda t, s: \ np.dot(np.array([[0,1],[-1,np.sin(t)/s[1]]]),s)
#mendefinisikan interval
t_span = np.linspace(0, a, 100)
#mendefinisikan solusi analitik
y_analitik = np.sin(t_span) - (1/2)*t_span*np.cos(t_span)
#kondisi batas
y0 = 0
#Tebakan awal untuk y'(0)
v0 = 0.55
t_eval = np.linspace(0, a, 10)
#melakukan optimasi terhadap nilai y'(0) menggunakan metode numerik
def objective(v0):
sol = solve_ivp(F, [0, a], \ [y0, v0], t_eval = t_eval)
y = sol.y[0]
return y[-1] - 1
v0 = fsolve(objective, 10)
#menyelesaikan masalah syarat batas dengan metode numerik
sol = solve_ivp(F, [0, a], \ [y0, v0], t_eval = t_eval)
#plot solusi analitik dan solusi numerik
plt.plot(t_span,y_analitik,'*k',label=r'solusi analitik')
plt.plot(sol.t, sol.y[0],'r',label='numerik')
plt.xlabel('x')
plt.ylabel('y')
plt.title(r'$y^{\prime\prime} + y = \sin (x); y(0) = 0 , y(\pi/2) = 1$')
plt.legend(loc=0)
plt.grid(False)
plt.show()
Metode numerik yang digunakan pada proses tersebut (penentuan solusi numerik dari $y$) adalah shooting method. Jika kode python tersebut dijalankan, maka diperoleh plot solusi analitik dan numerik sebagaimana diperlihatkan pada gambar berikut.
Dari gambar tersebut, kita dapat mengatakan bahwa solusi numerik dan solusi analitiknya berimpit. Salah satu poin penting pada kode tersebut adalah bahwa kita menggunakan metode numerik untuk menentukan nilai dari $y^{\prime} (0)$ dengan terlebih dahulu memberikan nilai awal. Metode numerik yang digunakan adalah metode secan.
Demikian pembahasan kali ini tentang Simulasi Numerik Persamaan Diferensial Biasa Orde 2. Postingan ini termasuk dalam topik persamaan diferensial. Jikan Anda tertarik dengan topik di persamaan diferensial lainnya, silahkan ke sini. Jika Anda tertarik dengan topik/materi lainnya, silahkan ke sini. Semoga membantu. Sekian dan terima kasih.
