SDP - Sonlu Farklar Yöntemi

Yazar

Taygun Bulmus

Yayınlanma Tarihi

26 Ağustos 2024

Sonlu Farklar Yöntemi

Sonlu farklar yöntemini kullanarak sınır değer problemine ait diferansiyel denklem çözeceğiz. Bunun için Taylor açılımını kullanacağız. Taylor açılımında birinci türevi orta noktayı kullanarak yazalım.

\[ y'_{i} = \frac{y_{i+1}-y_{i-1}}{2h} \]

Farz edelim ki bu fonksiyonu \(g_{i}\) için yazmışız.

\[ g'_{i} = \frac{g_{i+1}-g_{i-1}}{2h} \]

Bu \(g_{i}\) fonksiyonu \(y'_{i}\) fonksiyonuna eşit olsun. Yani \(g_{i}=y'_{i}\) olur. Bu durumda

\[ \begin{align*} y''_{i} &= \frac{y'_{i+1}-y'_{i-1}}{2h} \\ &= \frac{1}{2h}\left(\frac{y_{i+2}-y_{i}}{2h}-\frac{y_{i}-y_{i-2}}{2h}\right) \end{align*} \]

Biraz düzenleyelim.

\[ \begin{align*} y''_{i} &= \frac{1}{4h^2}\left(y_{i+2}-y_{i}-y_{i}+y_{i-2}\right) \\ &= \frac{y_{i+2}-2y_{i}+y_{i-2}}{4h^2} \end{align*} \]

Eğer adım aralığını yarıya indirirsek yani \(h \rightarrow h/2\) yaparsak \(y_{i+2}\) noktasından \(y_{i+1}\) noktasına gelebiliriz. Benzer şekilde \(y_{i-2}\) noktasından \(y_{i-1}\) noktasına geliriz. O halde ikinci türeve ait olan denklem şöyle olur.

Not

Bu işlem sadece sabit adım aralığı için geçerlidir.

\[ y''_{i} = \frac{y_{i+1}-2y_{i}+y_{i-1}}{h^2} \]

Bu denklemi kullanarak sınır değer problemlerini çözeceğiz. Yukarıdaki bağıntı, diferansiyel denklemi cebirsel denklem sistemine dönüştürür.

\[ y''_{i} = f\left(x,y_{i},\frac{y_{i+1}-y_{i-1}}{2h}\right)=\frac{y_{i+1}-2y_{i}+y_{i-1}}{h^2} \]

Alıştırma 1

Aşağıdaki denklemi adım adım çözelim.

\[ y''(t)=t, \quad y(0)=1, \quad y(3)=9 \]

Çözüm (Elle)

Adım aralığı \(h=1\) olsun. Çözeceğimiz \(t\)’nin değerleri \(0, 1, 2, 3\) olacaktır. Bu değerlerden \(y\) değerlerini bulalım.

\[ y''(t)=t=\frac{y_{i+1}-2y_{i}+y_{i-1}}{h^2} \]

Çözüme \(y(0)\)’dan başlayalım. Bunun için \(y(0) \equiv y_{i-1}\) almamız gerekecektir. Aksi taktirde köşe noktalar için denklemi çözemeyiz.

Aşağıdaki tabloyu \(h=1\) için yazalım.

\(y_{0}\) \(y_{1}\) \(y_{2}\) \(y_{3}\)
\(y(0)\) \(y(0+h)\) \(y(0+2h)\) \(y(0+3h)\)
\(y(t=0)\) \(y(t=1)\) \(y(t=2)\) \(y(t=3)\)
1 ? ? 9

Örneğin aşağıdaki tabloyu da \(h=0.5\) için yazalım.

\(y_{0}\) \(y_{1}\) \(y_{2}\) \(y_{3}\) \(y_{4}\) \(y_{5}\) \(y_{6}\)
\(y(0)\) \(y(0.5)\) \(y(1)\) \(y(1.5)\) \(y(2)\) \(y(2.5)\) \(y(3)\)
1 ? ? ? ? ? 9

Şimdi yukarıdaki elde ettiğimiz denklemi yazalım.

\[ y''(t)=t=\frac{y_{i+1}-2y_{i}+y_{i-1}}{h^2} \]

\(i=1\) olmalı. Bundan dolayı \(t\) ilk değil bir sonraki adımından yani \(t+h=0+0.5\)’den itibaren başlatmalıyız. \[ 0.5 = \frac{y_{2}-2y_{1}+y_{0}}{h^2} \]

Bu denklemi \(h=0.5\) için adım adım yazalım.

\[ \begin{align*} 0.5\times (0.5)^{2} =& y_{2}-2y_{1}+y_{0}\\ 1\times (0.5)^{2} =& y_{3}-2y_{2}+y_{1}\\ 1.5\times (0.5)^{2} =& y_{4}-2y_{3}+y_{2}\\ 2\times (0.5)^{2} =& y_{5}-2y_{4}+y_{3}\\ 2.5\times (0.5)^{2} =& y_{6}-2y_{5}+y_{4} \end{align*} \]

Yukarıdaki denklemleri matris şeklinde yazalım. Yazarken \(y_{0}=1\) ve \(y_{6}=9\) yazacağız.

\[ \begin{align*} -0.875 =& y_{2}-2y_{1}\\ 0.25 =& y_{3}-2y_{2}+y_{1}\\ 0.375 =& y_{4}-2y_{3}+y_{2}\\ 0.5 =& y_{5}-2y_{4}+y_{3}\\ -8.375 =& -2y_{5}+y_{4} \end{align*} \]

Matris şeklinde yazalım.

\[ \begin{bmatrix} -2 & 1 & 0 & 0 & 0\\ 1 & -2 & 1 & 0 & 0\\ 0 & 1 & -2 & 1 & 0\\ 0 & 0 & 1 & -2 & 1\\ 0 & 0 & 0 & 1 & -2\\ \end{bmatrix} \begin{bmatrix} y_{1}\\ y_{2}\\ y_{3}\\ y_{4}\\ y_{5}\\ \end{bmatrix}= \begin{bmatrix} -0.875\\ 0.25\\ 0.375\\ 0.5\\ -8.375 \end{bmatrix} \]

Bu denklemi çözmeyi biliyoruz. Denklemi çözdüğümüzde \(y_{1}, y_{2}, y_{3}, y_{4}, y_{5}\) değerlerini bulmuş oluruz.

Çözüm (Kod)

import numpy as np
from scipy.linalg import solve
import matplotlib.pyplot as plt
# Katsayı Matrisi
coefMat= np.array([\
      [-2, 1, 0, 0, 0]\
    , [ 1,-2, 1, 0, 0]\
    , [ 0, 1,-2, 1, 0]\
    , [ 0, 0, 1,-2, 1]\
    , [ 0, 0, 0, 1,-2]])
# Sonuç Matrisi
solMat= np.array([-0.875,0.25,0.375,0.5,-8.375])
# Doğrusal Denklemi Çöz
solution=solve(coefMat, solMat)
# Başlangıç Değeri
y0= 1.
# Çözümü Başına Ekle
solution=np.insert(solution, 0, y0)
# Bitiş Değeri
ySon= 9.
# Çözümü Sonuna Ekle
solution=np.append(solution, ySon)
# Çözmek istediğimiz t değerleri
t= np.linspace(0, 3, 7)
# Çiz
plt.plot(t, solution, 'o-', label='Sonlu Farklar')
# --------------------------------------------
# Analitik Çözüm
# y''= t y(0)=1 y(3)=9
# y'= t^2/2  + C1
# y= t^3/6  + C1 t + C2
# y(0)=1 => C2=1
# y(3)=9 => 9=27/6 + 3C1 + 1=9 => C1=27/18
#=> y(t)= t^3/6 + 21t/18 + 1
# --------------------------------------------
plt.plot(t, t**3/6 + 21*t/18 + 1, 'r-', label='Analitik')
plt.legend()
plt.show()

Matris denklemini el ile yazmadan yani başlangıçta biz vermeden çözmeye çalışalım.

import numpy as np
from scipy.linalg import solve
import matplotlib.pyplot as plt
# Fonksiyon
def fonk(t):
    return t
# Sınır Değerleri
y_t0=1
y_t3=9
# Adım aralığı
h=0.5
# Tum t değerleri
t0=0
tSon=3
# adet=(tSon-t0/h)+1
adet= int((tSon-t0)/h)+1
tTum=np.linspace(t0, tSon, adet)
# Yerine koyduğumuzda oluşan denklemler (tabloya bakınız)
# y''(t)=t = (y_{i+1}-2y_{i}+y_{i-1})/h^2
# t * h^2 = y_{i+1}-2y_{i}+y_{i-1}
# ilk denklem için t= h olacak. i=1 olacak
# h * h^2 = y_{2}-2y_{1}+y_{0}
# ilkDenklemSagTaraf= h*(h**2)-y_{0}
# Son denklem için t= tSon-h olacak. i=son-1 olacak
# (tSon-h) * h^2 = y_{son}-2y_{son-1}+y_{son-2}
# sonDenklemSagTaraf=(tSon-h)*(h**2)-y_{son}
ilkDenklemSagTaraf=fonk(h)*(h**2)-y_t0
sonDenklemSagTaraf=fonk(tSon-h)*(h**2)-y_t3
# Ara denklemleri oluştur.
# araDenklemlerSagTaraf[0]=ilkDenklemSagTaraf
araDenklemlerSagTaraf=np.array(ilkDenklemSagTaraf)
for i in range(1, adet-3):
    araDenklemlerSagTaraf= np.append(araDenklemlerSagTaraf, fonk(tTum[i+1])*(h**2))
# Son denklemi ekle
sonucMatrisi= np.append(araDenklemlerSagTaraf,sonDenklemSagTaraf)
# Katsayı Matrisi
katsayiMat= np.zeros((adet-2, adet-2))
# y''(t)=t = (y_{i+1}-2y_{i}+y_{i-1})/h^2
# Denkleminde katsayılar hep aynı olur.
for i in range(adet-2):
    katsayiMat[i, i]= -2
    if not i == adet-3:
        katsayiMat[i, i+1]= 1
    if not i == 0:
        katsayiMat[i, i-1]= 1
# Çöz
print(f"Katsayı Matrisi: \n{katsayiMat}")
print(f"Sonuç Matrisi: \n{sonucMatrisi}")
cozum=solve(katsayiMat, sonucMatrisi)
# İlk Çözümü Ekle
cozum=np.insert(cozum, 0, y_t0)
# Son Çözümü Ekle
cozum=np.append(cozum, y_t3)
# Çiz
plt.plot(tTum, cozum, 'o-', label='Sonlu Farklar')
plt.plot(tTum, tTum**3/6 + 21*tTum/18 + 1, 'r-', label='Analitik')
plt.legend()
plt.show()
Katsayı Matrisi: 
[[-2.  1.  0.  0.  0.]
 [ 1. -2.  1.  0.  0.]
 [ 0.  1. -2.  1.  0.]
 [ 0.  0.  1. -2.  1.]
 [ 0.  0.  0.  1. -2.]]
Sonuç Matrisi: 
[-0.875  0.25   0.375  0.5   -8.375]

Sonlu Farklar Yönteminin Çalışmadığı Örnek

Eğer oluşturduğumuz diferansiyel denklem sistemi doğrusal değilse, yukarıdaki algoritmayı kullanmayız. Aşağıdaki örneği inceleyelim.

\[ y''(t)=tyy'(t) \]

Çözmemiz gereken denklem

\[ y''(t)=tyy'(t)=\frac{y_{i+1}-2y_{i}+y_{i-1}}{h^2} \]

\(i=1\) için yazalım.

\[ ty_{1}y_{1}'=\frac{y_{2}-2y_{1}+y_{0}}{h^2} \]

Türevi de açıkça yazalım.

\[ ty_{1} \frac{y_{2}-y_{0}}{2h}=\frac{y_{2}-2y_{1}+y_{0}}{h^2} \]

Düzenleyelim.

\[ \frac{th}{2} (y_{1}y_{2}-y_{1}y_{0}) -y_{2}+2y_{1}-y_{0} = 0 \]

\(y_{1}y_{2}\) gibi terimlerin varlığından dolayı bu denklem doğrusal olmayan (nonlineer) denklemdir. Klasik yöntemlerle (LU ayrıştırma veya Gauss eleme) çözülemez.

Problemler

Problem 1

Aşağıdaki diferansiyel denklemi sonlu farklar yöntemi kullanarak çözebilir misinz?

\[ \frac{\hbar^{2}}{2 m}\frac{d^{2}}{dx^{2}}\psi(x) = E\psi(x), \qquad \psi(x=0)=0, \qquad \psi(x=L)=0 \text{.} \]

Çözümü \(n=1\) ve \(n=2\) için bulmaya çalışın. Hem bulduğunuz enerji değerlerini hem de dalga fonksiyonunun karesini (olasılık yoğunluğunu, \(|\psi(x)|^{2}\)) çizdirin. \(x\) değerleri \(0\) ile \(L\) arasında \(N=100\) adet olsun.

Denklemdeki katsayıları aşağıdaki gibi alın.

  1. \(\hbar=1\)
  2. \(m=1\) kg
  3. \(L=2\)

Enerji özdeğerinin gerçek değeri aşağıdaki gibi olmalıdır.

\(E=\frac{\pi^{2}\hbar^{2}}{2mL^{2}}n^{2}\) J.

Dalga fonksiyonunun analitik çözümü aşağıdaki gibi olmalıdır.

\[ \psi(x)=\sqrt{\frac{2}{L}}\sin\left(\frac{n\pi x}{L}\right) \]

Tavsiye

Bu problem bir özdeğer problemidir. Diferansiyel operatörü \(\left(\frac{d^{2}}{dx^{2}}\right)\) doğrusal operatörü matris şeklinde yazmalısınız. Bu matrisi sonlu farklar yönteminde yazdık.

Problem 2

Aşağıdaki diferansiyel denklemi sonlu farklar yöntemini kullanarak çözün. Adım aralığını \(h=0.05\) alın. Analitik sonuçla karşılaştırın.

\[ y''(t) = t^{3}+t+5, \quad y(0)=1, \quad y(6)=545.8 \]

Kaynaklar

  1. Numerical Methods in Engineering with Python 3, Jaan Kiusalaas, Cambridge University Press, 2013, Syf:307
  2. Sonlu Farklar Yöntemi, Youtube Videosu
  3. Numerical Solutions of Schrodinger’s Equation, Neill Lambert, 2001