0 x L E rh a n C o¸s k u n , K a ra d e n iz Te k n ik M a t e m a tik , e rh a n @ k t u .e d u .tr (2)ba¸slang¬ç de¼ger problemi kapsam¬nda ifade edilebilen difüzyon veya ¬s¬denk- lemi ve daha sonra ise utt = F (x

(1)

B ¨ol ¨um 8

Ba¸slang¬ç-S¬n¬r-De¼ger Problemleri Için Sonlu Fark Yöntemleri·

Bu bölümde amac¬m¬z verilen bir ba¸slang¬ç-s¬n¬r de¼ger problemine ait en uygun sonlu fark yönteminin seçilerek, uygun ad¬m uzunluklar¬ile uygulana- bilmesi için gereki bilgi ve deneyimin kazand¬r¬lmas¬n¬sa¼glamakt¬r.

Bu amaçla tipik baz¬tek boyutlu ba¸slang¬ç-s¬n¬r de¼ger problemlerini in- celeyerek, sözkonusu problemler için

bir aç¬k yöntemin farkl¬ s¬n¬r ¸sartlar¬yla türetilmesi, uygulanmas¬ ve hatas¬,

yönteme ait kararl¬l¬k kriterinin elde edili¸si ve sa¼glanma zorunlulu¼gu, aç¬k yöntemlere k¬yasla daha büyük ad¬m uzunlu¼gu kullan¬m¬na izin veren kapal¬bir yöntemin türetilmesi ve uygulanmas¬ile

benzer alternatif yöntemleri pratik problemler üzerinde MATLAB/Octave ortam¬nda "vektör tabanl¬" kodlar yard¬m¬yla inceliyoruz.

8.1 Giri¸s

Bu bölümde öncelikle tek boyutlu

u_t = F (x; t; u; u_x; u_xx); 0 < x < L; t > 0 (8.1) u(x; 0) = f (x); 0 x L

E rh a n C o¸s k u n , K a ra d e n iz Te k n ik M a t e m a tik , e rh a n @ k t u .e d u .tr

(2)

ba¸slang¬ç de¼ger problemi kapsam¬nda ifade edilebilen difüzyon veya ¬s¬denklemi ve daha sonra ise

u_tt = F (x; t; u; u_t; u_x; u_xx); 0 < x < L; t > 0 (8.2) u(x; 0) = f (x);

u_t(x; 0) = g(x); 0 x L (8.3)

problemi kapsam¬nda ifade edilebilen dalga denklemi gibi tipik baz¬modelleri reel katsay¬l¬

a₁₁u(0; t) + a₁₂u_x(0; t) = (t); t > 0 (8.4) a21u(L; t) + a22ux(L; t) = (t); t > 0

s¬n¬r ¸sartlar¬ile birlikte gözönüne alarak, aç¬k ve kapal¬sonlu fark denklemlerini üretiyor ve ilgili fark yöntemleri yard¬m¬yla tipik baz¬ problemlerin say¬sal çözümlerini elde ediyoruz. Burada a²₁₁ + a²₁₂ > 0 ve a²₂₁ + a²₂₂ > 0 oldu¼gunu, di¼ger bir deyimle a11 ve a12 katsay¬lar¬ndan en az birinin ve benzer biçimde a21ve a22katsay¬lar¬ndan da en az birinin s¬f¬rdan farkl¬oldu¼gunu kabul ediyoruz.

Öncelikle (8.1) biçiminde yaz¬labilen bir parabolik problemi inceleyelim.

8.2 Parabolik problemler için aç¬k bir yön- tem

Difüzyon denklemi veya yanal yüzeyleri yal¬t¬lm¬¸s olan L uzunluklu iletken bir çubuktaki s¬cakl¬k da¼g¬l¬m¬n¬belirleyen

u_t= cu_xx+ g(x; t) (8.5)

denklemini

=f(x; t)j0 < x < L; t > 0g üzerinde

u(x; 0) = f (x); 0 x L (8.6)

ba¸slang¬ç ve

u(0; t) = (t); u(L; t) = (t); t > 0 (8.7)

(3)

Dirichlet s¬n¬r de¼gerleri ile birlikte gözönüne alal¬m. Burada t zaman; x yer de¼gi¸skeni ve u(x; t) ise x noktas¬nda ve t an¬nda çubu¼gun s¬cakl¬¼g¬n¬göster- mektedir. Ayr¬ca g(x; t) ise herhangi bir nedenle(elektrik ak¬m¬na kar¸s¬direnç vb) olu¸sabilen iç ¬s¬kayna¼g¬n¬temsil etmekte ve c ise iletkenlik katsay¬s¬n¬da içeren difüzyon parametresidir.

(8.5)-(8.7) probleminin de¼gi¸skenlere ay¬rma yöntemi ad¬ verilen yöntem yard¬m¬yla sonsuz seri biçiminde ifade edilebilen analitik çözümü genelde elde edilebilir. Ancak elde edilen seri çözümlerinden de pratik sonuçlar elde edebilmek için sadece sonlu say¬da terimin toplam¬dikkate al¬nmal¬d¬r. Do- lay¬s¬yla elde etti¼gimiz sonuçlar, gerçek de¼gerler yerine yine belirli hatalar içeren yakla¸s¬mlar olacakt¬r.

Bu durumda alternatif bir yakla¸s¬m ise problemin say¬sal çözümünü elde etmektir. Öncelikle zaman de¼gi¸skenine göre ileri fark ve yer de¼gi¸skenine göre merkezi fark ayr¬kla¸st¬rma yöntemini uygulayarak olu¸san ve aç¬k yöntem olarak bilinen fark denklemlerini yazal¬m:

[0; L]aral¬¼g¬n¬N adet e¸sit uzunluklu alt aral¬¼ga bölelim ve elde edilen alt aral¬klar¬n uç noktalar¬n¬s¬ras¬yla

x₁ = 0; x₂ = x; ; x_{N +1} = L

ile gösterelim. Bu durumda aral¬klar¬n uzunluklar¬ x = L=N olur. Her- hangi T > 0 için problemi [0; T ] sonlu zaman aral¬¼g¬nda inceleyelim. Bu aral¬¼g¬ M adet t = T =M uzunluklu alt aral¬¼ga bölelim ve elde edilen alt aral¬klar¬n uç noktalar¬n¬s¬ras¬yla

t₁ = 0; t₂ = t; ; t_{M +1}= T ile gösterelim(¸Sekil 8.1).

Bu durumda (8.5),(8.6) ve (8.7) problemin s¬n¬rlar¬içeren tan¬m kümesi

= [0; L] [0; T ] dir.

Bu problem için olu¸sturulacak fark denkleminin tan¬m kümesini ise _ ile gösterelim. _ , ¸Sekil 8.1 deki yatay ve dü¸sey do¼grular¬n dü¼güm ad¬

verilen arakesit noktalar¬n¬n kümesidir:

Küme notasyonu ile

_ = f(xⁱ; t_j)2 : i = 1 : N + 1; j = 1 : M + 1g

E rh a n C o¸s k u n , K a ra d e n iz Te k n ik M a te m a tik , e rh a n @ k t u .e d u .tr

(4)

¸

Sekil 8.1: Aç¬k yöntemde kullan¬lan dü¼güm noktalar¬

dir. _ kümesinin s¬n¬rlar¬n¬

_sol =f(0; t^j)j j = 2 : M + 1g; _^sag =f(L; t^j)j j = 2 : M + 1g olmak üzere _ = _sol[ _^sag ile gösterelim.

¸

Sekil 8.1 de belirtilen (xi; tj)2 _ n _ noktas¬nda t zaman de¼gi¸skenine göre ileri fark ve x yer de¼gi¸skenine göre merkezi fark türev yakla¸s¬mlar¬kullan¬l¬rsa verilen problem için

(u(x_i; t_j+1) u(x_i; t_j))= t + O( t)

= c(u(x_{i 1}; t_j) 2u(x_i; t_j) + u(x_i+1; t_j))=( x)²

+O(( x)²) + g(x_i; t_j) (8.8)

yakla¸s¬m¬elde edilir. (xi; t_j+1) noktas¬ndaki yakla¸s¬m¬hesaplamak için ¸Sekil 8.1 de tj sat¬r¬nda y¬ld¬zlarla gösterilen dü¼güm noktalar¬ndaki de¼gerlerin kullan¬ld¬¼g¬na dikkat edelim. (8.8) de O( t) ve O(( x)²)terimlerinin ihmal edilmesiyle (xi; t_j)2 _ noktas¬nda fark denklemini sa¼glayan say¬sal yakla¸s¬m¬

U_i;j = u(xi; t_j)

ile gösterelim ve ayr¬ca notasyonel kolayl¬k aç¬s¬ndan Gi;j := g(x_i; t_j) olarak tan¬mlayal¬m.

(5)

_ ayr¬k bölgesinin indislerinden olu¸san

__ij =f(i; j) : i = 1 : N + 1; j = 1 : M + 1g kümesine a¼g ve

__j = __1;j[ _^{N +1;j} =f(1; j); (N + 1; j) : j = 2 : M + 1g kümesine ise a¼g s¬n¬r¬ad¬verelim. S¬n¬rlar d¬¸s¬nda kalan

( _n _ )^ij := _ijn _^ij=f(i; j) : i = 2 : N; j = 1 : Mg indisler kümesine ise hesaplama a¼g¬ad¬n¬verelim.

u(x; t)analitik çözümü üzerinde; u(xi; t_j)de¼gerleri _ üzerinde ve bu de¼ger için yakla¸s¬m olan Uij ler ise _ij üzerinde tan¬ml¬d¬r.

Herhangi (i; j) 2 ( _ n _ )^ij için

(U_i;j+1 U_i;j)= t = c(U_{i 1;j} 2U_i;j+ U_i+1;j)=( x)²+ G_i;j (8.9) elde ederiz. r = c t=( x)² olarak tan¬mlayarak, (8:9) sistemini Ui;j+1

e göre çözmek suretiyle her (i; j) 2 ( _ n _ )^ij için

U_i;j+1= U_i;j+ r(U_{i 1;j} 2U_i;j + U_i+1;j) + tG_i;j (8.10) fark denklemini elde ederiz.

u(x; 0) = f (x) ba¸slang¬ç ¸sart¬ise s¬ras¬yla _ ve _i1 üzerinde

u(x_i; 0) = f (x_i)ve U_i;1 = f (x_i); i = 1 : N + 1 (8.11) olarak ifade edilebilir.

u(0; t) = (t) sol s¬n¬r ¸sart¬ise _sol ve _1j üzerinde s¬ras¬yla u(0; t_j) = (t_j) ve U_1;j = (t_j); j = 2 : M + 1

ve u(L; t) = (t) sa¼g s¬n¬r ¸sart¬ise _sag ve _N +1j üzerinde s¬ras¬yla u(L; t_j) = (t_j) ve U_{N +1;j} = (t_j); j = 2 : M + 1 (8.12) olarak ifade edilebilir.

(6)

(8.10),(8.11) ve (8.12) denklemlerine verilen sürekli problemin, yani (8.5),(8.6) ve (8.7) probleminin fark veya daha genel bir ifade ile sonlu say¬da fark içerdi¼gi için sonlu fark denklemleri ad¬verilmektedir.

Yöntemin hatas¬: Baya¼g¬diferensiyel denklemler k¬sm¬ndan hat¬rla- yaca¼g¬m¬z üzere, yöntemin kesme hatas¬ gerçek çözümün fark denklemini sa¼glamad¬¼g¬miktard¬r. (8.8) ve (8.9) kar¸s¬la¸st¬r¬ld¬¼g¬nda

E_kesme(t; t; x) = O( t) + O( x²); t! 0; x ! 0 (8.13) oldu¼gu kolayca görülür. Benzer biçimde yöntemin yerel hatas¬

Eyerel(t; t; x) = O( t²) + O( x³); t! 0; x ! 0 (8.14) ve kümülatif hatas¬n¬n da

E_kum(t; t; x) = O( t) + O( x²); t! 0; x ! 0 (8.15) oldu¼gu baya¼g¬diferensiyel denklemler k¬sm¬ndaki hata analizi yöntem- leri yard¬m¬yla görülebilir(Al¬¸st¬rma 3).

ÖRNEK 8.1. (Aç¬k yöntemin Dirichlet S¬n¬r ¸Sartlar¬ile uygulan¬¸s¬) u_t = 2u_xx; 0 < x < 1; t > 0

u(x; 0) = sin( x); 0 x 1 u(0; t) = t; u(1; t) = 0; t > 0

ile = f(x; t)j0 < x < 1; t > 0g üzerinde tan¬mlanan ba¸slang¬ç-s¬n¬r de¼ger problemi verilmi¸s olsun.

(a) x; t ad¬m uzunluklar¬ile verilen probleme _ üzerinde aç¬k yöntem ile kar¸s¬l¬k gelen fark denklemlerini ve _ij deki kar¸s¬l¬klar¬n¬ifade ediniz.

(b) x = 1=2 ve t = 1=16 için verilen ba¸slang¬ç-s¬n¬r de¼ger probleminin (x; t) = (1=2; 1=16) noktas¬ndaki say¬sal çözümü bilgisayar kullanmadan hesaplay¬n¬z.

(c) Problemin analitik çözümünün

u(x; t) = 1= ³+ (1 + 1= ³)e ² ²^t sin( x) +

X1 n=2

1

3n³ e ² ²ⁿ²^t 1 sin(n x) + t(1 x)

(7)

olarak verildi¼gini kontrol ediniz(Verilen ifadenin ilgili k¬smi diferensiyel ile ba¸slang¬ç ve s¬n¬r de¼gerlerini sa¼glad¬¼g¬n¬gösteriniz(ipucu:

x 1 = 2

X1 n=1

1

n sin(n x); 0 < x < 1 Fourier aç¬l¬m¬n¬gözönüne al¬n¬z, ayr¬ca bknz Al¬¸st¬rma 9).

(d) kümesinde x = 1=10 ad¬m uzunlu¼gu ve

r = c t=( x)² = 2 t=(0:1)² = 200 t = 0:5

de¼gerine kar¸s¬l¬k gelen t = 0:0025 ile ilgili say¬sal yakla¸s¬mlar¬ elde ediniz.

Say¬sal yakla¸s¬mlar¬ gra…ksel olarak, T = 0:2 için [0; T] zaman aral¬¼g¬nda gösteriniz.

(e) Elde etti¼giniz say¬sal yakla¸s¬mlardan

_ = f(x; t)jx = 0 : 0:2 : 1; t = 0:1; 0:15; 0:2g kümesine ait olanlar¬tablo halinde listeleyiniz.

(f ) _ da kümülatif hata de¼gerlerini elde ederek tablo halinde listeleyiniz.

(g) x = 1=10 ad¬m uzunlu¼gu için r = 0:55 de¼gerine kar¸s¬l¬k gelen t = 0:00275 ile olu¸san _ üzerindeki yakla¸s¬mlar¬ gra…ksel olarak gösteriniz.

Sonuçlar¬n¬z …ziksel olarak anlaml¬m¬d¬r? Elde etti¼giniz sonuçlar parabolik problemler için maximum prensibi ile uyumlu mudur?

Çözüm.

(a) (xi; t_j)2 _ n _ için aç¬k yöntem ile probleme ait fark denklemlerinia¸sa-

¼g¬daki gibi elde ederiz:

u(x_i; t_j+1) u(x_i; t_j+1)

t + O( t) = 2u(x_{i 1}; t_j) 2u(x_i; t_j) + u(x_i+1; t_j) x²

+O( x²):

Yukar¬da ifade edilen fark denklemlerinin ( _ n _ )^ij deki kar¸s¬l¬¼g¬ise U_i;j+1 U_i;j

t = 2U_i+1;j 2U_i;j + U_{i 1;j} x²

(8)

veya

U_i;j+1 = U_i;j+ r(U_i+1;j 2U_i;j+ U_{i 1;j}); r = 2 t x² olarak elde edilir.

u(x; 0) = sin( x) ba¸slang¬ç ¸sart¬ise _ ve _i1 üzerinde s¬ras¬yla u(x_i; 0) = sin( x_i) ve U_i;1 = sin( x_i); i = 1 : N + 1 olarak ifade edilebilir.

u(0; t) = t sol s¬n¬r ¸sart¬ise _sol ve _1j üzerinde s¬ras¬yla u(0; t_j) = t_j ve U_1;j = t_j; j = 2 : M + 1

ve u(1; t) = 0 s¬n¬r ¸sart¬ise _sag ve _N +1;j üzerinde s¬ras¬yla u(1; t_j) = 0 ve U_{N +1;j} = 0; j = 2 : M + 1

olarak ifade edilebilir.

(b) x = 1=2olmas¬[0; 1] aral¬¼g¬n¬n iki alt parçaya bölünmesini gerektirir.

O halde alt aral¬klar¬n uç noktalar¬x1 = 0; x2 = 1=2; x3 = 1olarak elde edilir. Benzer biçimde t = 1=16 için t1 = 0 ve t2 = 1=16 elde ederiz.

O halde _ ayr¬k kümesi

(x₁; t₂) = (0; 1=16) (x₂; t₂) = (1=2; 1=16) (x₃; t₂) = (1; 1=16) (x₁; t₁) = (0; 0) (x₂; t₁) = (1=2; 0) (x₃; t₁) = (1; 0)

noktalar¬ndan olu¸sur. (xi; t_j) noktas¬ndaki yakla¸s¬m Uij olmak üzere, verilen problemin ba¸slang¬ç ve s¬n¬r ¸sartlar¬ndan a¸sa¼g¬daki de¼gerleri yazabiliriz.

Sol s¬n¬r ¸sart¬ Sa¼g s¬n¬r ¸sart¬

#(U(0; t) = t) #(U(1; t) = 0) t₂ = 1=16 U₁₂= 1=16 U₂₂=? U₃₂= 0 t₁ = 0 U₁₁= 0 U₂₁= 1 U₃₁= 0 x₁ = 0 x₂ = 1=2 x₃ = 1

Bu örnek için bilinmeyen tek de¼ger olan U22 de¼gerini elde etmek için (8.10) da i = 2; j = 1 yazarak (Gij = 0 oldu¼guna dikkat edelim)

U₂₂= U₂₁+ r(U₁₁ 2U₂₁+ U₃₁)

(9)

elde ederiz.

r = c t=( x)² = 2(1=16)=(1=2)² = 1=2

için yukar¬daki tabloda belirtilen de¼gerleri yerine yazmak suretiyle U₂₂ = 1 + 1=2(0 2 1 + 0) = 0

elde ederiz.

(c) Verilen ifadenin ba¸slang¬ç ve s¬n¬r ¸sartlar¬n¬ sa¼glad¬¼g¬ aç¬kça görülür.

Ayr¬ca

a_n(t) = e ² ²ⁿ²^tsin(n x)

fonksiyonu her n 1için verilen k¬smi diferensiyel denklemi sa¼glad¬¼g¬n- dan hareketle kolayca görülebilir.

(d) [0; 1] aral¬¼g¬nda x = 1=10 olmas¬, aral¬¼g¬n N = 10 alt aral¬¼ga bölün- dü¼günü ifade etmektedir. t = 0:0025 almak suretiyle

M = [[T =dt]] = 80 elde edilir.

x = 1=10 için _ ayr¬k bölgesi için x1 = 0; x₂ = 0:1; :::; x₁₁ = 1 elde edilir. Benzer biçimde t = 0:0025 için

t₁ = 0; t₂ = t; :::; t₈₁= 1=5 elde edilir.

Yönteme ait yakla¸s¬mlar Program 8.1 ile elde edilmi¸stir.

Elde edilen say¬sal çözüm gra…ksel olarak ¸Sekil 8.2 de sunulmaktad¬r.

(e) ¸Sekil 8.2 de sunulan çözüme seçilen baz¬(xi, tj)2 _ de¼gerleri için elde edilen say¬sal de¼gerler Tablo 8.1 de listelenmektedir:

(f) Analitik çözümün _ üzerindeki de¼gerlerini belirlemek amac¬ylaa¸sa¼g¬- daki KDDAnalitik program parças¬n¬kullanabiliriz(Bknz Program 8.2).

Çözümde yer alan serinin ilk 100 teriminin toplam¬n¬ analitik çözüm olarak kabul ediyoruz.

(10)

0 0.2

0.4 0.6

0.8 1

0 0.1 0.2 0.3 0.4 0.5

0 0.2 0.4 0.6 0.8 1

t x

¸

Sekil 8.2: Örnek 8.1 in say¬sal çözümü

tnx 0 0.2 0.4 0.6 0.8 1

0 0 0.5878 0.9511 0.9511 0.5878 0

0.1000 0.1000 0.1375 0.1439 0.1439 0.0855 0 0.1500 0.1500 0.1259 0.0803 0.0803 0.0439 0 0.2000 0.2000 0.1470 0.0350 0.0697 0.0350 0

Tablo 8.1: Örnek 8.1 in say¬sal çözümünden baz¬degerler, t = 0:0025 _ üzerinde t = 0; 0:1; 0:15 ve 0:2 de¼gerlerine kar¸s¬l¬k gelen analitik çözüm tablosu Tablo 8.2 da verilmektedir:

_ üzerinde t = 0; 0:1; 0:15 ve 0:2 de¼gerlerine kar¸s¬l¬k gelen hata e_ij = u(x_i; t_j) U_ij

tablosu Tablo 8.3 de verilmektedir: Ba¸slang¬ç do¼grusu ve s¬n¬rlar üze- rinde hatan¬n s¬f¬ra e¸sit oldu¼guna dikkat edelim.

(g) t = 0:00275 ad¬m uzunlu¼gu için elde edilen yakla¸s¬mlar ¸Sekil 8.3 de gra…ksel olarak sunulmaktad¬r.

t = 0:00275ad¬m uzunlu¼gu ile olu¸san a¼gda seçilen baz¬x ve t de¼gerleri için elde edilen yakla¸s¬mlar ise Tablo 8.4 de verilmektedir.

Gerek Tablo 8.4 ve gerekse ¸Sekil 8.3 den görülece¼gi üzere elde edilen sonuçlar belirli bir tj an¬ndan itibaren problemin …ziksel an- lam¬yla ba¼gda¸smamaktad¬r. Çünkü ba¸slang¬çta uygulanan negatif olmayan s¬cakl¬k de¼gerlerine ra¼gmen yanal yüzeyleri izole edilen

(11)

tnx 0 0:2 0:4 0:6 0:8 1 0 0 0:5878 0:9511 0:9511 0:5878 0:0000 0:1 0:1 0:1667 0:1934 0:1694 0:0987 0:0000 0:15 0:15 0:1555 0:1406 0:1066 0:0506 0:0000 0:2 0:2 0:1764 0:1397 0:0957 0:0484 0:0000

Tablo 8.2: Analitik çözüm de¼gerleri

tnx 0 0:2 0:4 0:6 0:8 1

0 0 0 0 0 0 0:0000

0:1 0 0:029 2 0:049 5 0:025 5 0:013 2 0:0000 0:15 0 0:029 6 0:060 3 0:026 3 0:006 7 0:0000 0:2 0 0:029 4 0:104 7 0:026 0:013 4 0:0000

Tablo 8.3: Hata tablosu

çubuk üzerindeki s¬cakl¬k de¼gerleri negatif olmaktad¬r. Bu sonuç

…ziksel olarak anlaml¬de¼gildir.

Öte yandan elde edilen sonuçlar matematiksel olarak ta anlaml¬

de¼gildir. Bunu görmek için öncelikle parabolik problemler için maximum prensibini ana hatlar¬yla hat¬rlayal¬m:

Maksimum Prensibi:

u_t = cu_xx; 0 < x < L; t > 0 u(x; 0) = f (x); 0 x L

u(0; t) = (t); u(L; t) = (t); t > 0

parabolik ba¸slang¬ç-s¬n¬r de¼ger problemi denklemi maksimum de-

¼gerini t = 0 da ba¸slang¬ç de¼geriyle veya x = 0, x = L, t > 0 s¬n¬rlar¬ üzerinde al¬r. Benzer sonuç minimum için de geçerlidir.

Yani çözüm minimum de¼gerini de yine ba¸slang¬ç do¼grusu veya s¬n¬rlar boyunca al¬r.

Tablo 8.4 de¼gerlerine bakt¬¼g¬m¬zda minimum de¼gerlerin ba¸slang¬ç do¼grusu veya x = 0, x = L yerine T = 0:1925 do¼grusu üzerinde olu¸smaktad¬r. O halde elde edilen sonuçlar matematiksel olarak ta anlaml¬de¼gildir.

(g) de görüldü¼gü üzere matematiksel teoriye veya …ziksel kural- lara uygun olmayan bu tür sonuçlar ilgili say¬sal yöntemin seçilen

(12)

¸

Sekil 8.3: Örnek 8.1 için r = 0:55 > 1=2 kararl¬l¬k kriterine uymayan yakla¸s¬mlar

tnx 0 0:2000 0:4000 0:6000 0:8000 1:000

0 0 0:5878 0:9511 0:9511 0:5878 0

0:0550 0:0550 0:2204 0:3253 0:3183 0:1954 0 0:1100 0:1100 0:1275 0:1368 0:1188 0:0696 0 0:1650 0:1650 0:0871 0:0332 0:0043 0:0038 0 0:1925 0:1925 0:0248 0:1670 0:2015 0:1323 0

Tablo 8.4: Örnek 8.1 in say¬sal çözümünden baz¬degerler, t = 0:00275

( t; x)ad¬m uzunluklar¬ile karars¬z olmas¬ndan kaynaklanmak- tad¬r. Kararl¬l¬k kavram¬ve kararl¬l¬k için ad¬m uzunlu¼gu üzerin- deki k¬s¬tlamay¬Bölüm 13.2.1 de inceliyoruz.

ÖRNEK 8.2. (Aç¬k yöntemin Neumann S¬n¬r ¸Sartl¬problemlere uygulan¬¸s¬) u_t = u_xx; (x; t)2 =f(x; t)j0 < x < 1; t > 0g

u(x; 0) = sin( x); 0 x 1 u_x(0; t) = 1; u_x(1; t) = 0; t 0 ba¸slang¬ç de¼ger problemi verilsin.

(13)

(a) Verilen problemin analitik çözümün u(x; t) = t + 2 1

3+ 2

2e ²^tcos( x) 2

X1 m=2

m²( 1)^m+ ( 1)m²+ 1

2m²(m² 1) e ^m² ²^tcos(m x) +x(1 1

2x) (8.16)

biçiminde ifade edilebildi¼gini kontrol ediniz ve

x =f0; 0:2; 0:4; 0:6; 0:8; 1g ; t = 0:2

noktalar¬ndaki çözüm de¼gerlerini verilen serinin ilk 100 terimi ile hesaplay¬n¬z.

(b) x ekseni yönünde N ve t ekseni yönünde M adet alt aral¬k ile olu¸sturulan _ ayr¬k kümesi üzerinde aç¬k yöntem ile elde edilen fark denklemlerini ve bu denklemlerin _ijdeki kar¸s¬l¬klar¬n¬ yaz¬n¬z.

(c) x = 1=4; t = 1=32 için (xi; 1=32) noktas¬ndaki Ui;2; i = 1; 2; 3; 4; 5 say¬sal çözüm de¼gerlerini (b) de elde etti¼giniz fark denklemleri yard¬m¬yla ve bilgisayar kullanmadan hesaplay¬n¬z.

(d) x = 0:1; t = 0:0025için elde edilen _ kümesi üzerindeki say¬sal çözüm- leri belirleyiniz. Ayr¬ca

x =f0; 0:2; 0:4; 0:6; 0:8; 1g ve

t = f0; 0:06; 0:12; 0:18; 0:2g

de¼gerlerine kar¸s¬l¬k gelen say¬sal çözüm tablolar¬n¬olu¸sturunuz.

(e)

x = 0:1; 0:05; 0:0250; 0:0125; 0:00625 de¼gerlerine kar¸s¬l¬k

r = t=( x)² = 1=4 ba¼g¬nt¬s¬n¬sa¼glayan t de¼gerleri için

x =f0:2; 0:4; 0:6; 0:8g ; t = 0:2

noktalar¬ndaki say¬sal yakla¸s¬mlar¬tablo halinde kar¸s¬la¸st¬r¬n¬z.

(14)

t = 0:2 x = 0 0:131296554560213 x = 0:2 0:306016680909306 x = 0:4 0:432097632322202 x = 0:6 0:514700420692093 x = 0:8 0:560470187983312 x = 1 0:574997991968442

(f ) (e) de listeledi¼giniz sonuçlar ve analitik çözüm de¼gerlerinden x = 0:4; t = 0:2 noktas¬ndaki hata de¼gerlerini farkl¬ x ler için kar¸s¬la¸st¬r¬n¬z. Elde etti¼giniz sonuçlar yöntemin (8.15) ile verilen kümülatif hatas¬ile uyumlu mudur?

(g) x = 0:1; t = 0:0025 ile elde edilen say¬sayal çözüm gra…¼gini [0; 0:2]

aral¬¼g¬nda çizdiriniz.

Çözüm.

(a) Verilen çözüm serisinin denklemi sa¼glad¬¼g¬aç¬kça görülmektedir. S¬n¬r

¸sartlar¬n¬n da sa¼gland¬¼g¬kolayca görülmektedir.

u(x; 0) = sin( x) ba¸slang¬ç de¼gerinin sa¼gland¬¼g¬n¬ kontrol etmek için analitik çözümden t = 0 için elde edilen serinin sin( x) x(1 ¹₂x)fonksiyonunun (0; 1) aral¬¼g¬ndaki Fourier cosin•usaç¬l¬m¬

oldu¼gunu kontrol ediniz.

Belirtilen noktalardaki analitik çözüm de¼gerlerini hesaplayan Prog- ram 8.3 bölüm sonunda verilmektedir.

Analitik çözümden belirtilen noktalarda elde edilen say¬sal çözüm de¼gerleria-

¸sa¼g¬da verilmektedir.

(b) x ekseni yönünde N ve t ekseni yönünde M adet alt aral¬k ile olu¸sturulan ayr¬k küme

_ = f(xⁱ; tj)jx¹ = 0; x₂ = x; : : : ; xN +1 = 1; t₁ = 0; t₂ = t; g olarak elde edilebilir. S¬n¬r ¸sartlar¬n¬fark denklemlerine yans¬tmak için

bölgesinin d¬¸s¬nda

x₀ = x₁ x; x_{N +2}= x_{N +1}+ x sanal noktalar¬n¬n var oldu¼gunu kabul edelim.

(15)

ux(0; tj) = ux(x1; tj) = 1 s¬n¬r ¸sart¬, merkezi fark yakla¸s¬m¬yla (u(x₂; tj) u(x₀; tj))=2 x + O( x²) = 1

olarak ifade edilebilir ve ilgili say¬sal yakla¸s¬mlar U_2;j U_0;j = 2 x

ba¼g¬nt¬s¬n¬sa¼glar ve dolay¬s¬yla U0;j yapay s¬n¬r de¼gerleri U_0;j = U_2;j 2 x

olarak elde edilir.

u_x(1; t_j) = u_x(x_{N +1}; t_j) = 0 s¬n¬r ¸sart¬, merkezi fark yakla¸s¬m¬yla (u(x_{N +2}; t_j) u(x_N; t_j))=2 x + O( x²) = 0

olarak ifade edilebilir ve ilgili yakla¸s¬mlar ise U_{N +2}; j U_N;j = 0

ba¼g¬nt¬s¬n¬sa¼glar ve dolay¬s¬yla, UN +2; j yapay s¬n¬r de¼gerleri U_{N +2}; j = U_N;j

olarak elde edilir.

Fakat (0; tj)sol s¬n¬r¬üzerinde ut= u_xx denklemine aç¬k yöntemle kar¸s¬l¬k gelen fark denkleminden

(U_1;j+1 U_1;j)= t = (U_0;j 2U_1;j + U_2;j)=( x)² elde ederiz. S¬n¬r ¸sart¬ndan elde etti¼gimiz

U_0;j = U_2;j 2 x

de¼gerini bu denklemde yerine yazarak U0;j yapay s¬n¬r de¼gerlerini yok edebiliriz:

(U_1;j+1 U_1;j)= t = ( 2U_1;j+ 2U_2;j 2 x)=( x)² veya r = t=( x)² ile

U_1;j+1 = U_1;j + 2r(U_2;j U_1;j x) (8.17)

elde ederiz. O halde sol s¬n¬rda bulunan (x1; t_j+1) noktas¬ndaki yakla¸s¬m¬hesaplamak için (x1; t_j); (x₂; t_j)noktalar¬ndaki yakla¸s¬m- lar¬bilmek yeterlidir.

(16)

Benzer biçimde sa¼g s¬n¬r üzerindeki

(U_{N +1;j+1} U_{N +1;j})= t = (U_N;j 2U_{N +1;j}+ U_{N +2;j})=( x)² fark denkleminde

U_{N +2;j} = U_N;j oldu¼gunu kullanarak

(UN +1;j+1 UN +1;j)= t = 2(UN;j UN +1;j)=( x)² veya

UN +1;j+1 = UN +1;j + 2r(UN;j UN +1;j) (8.18) elde ederiz. O halde sa¼g s¬n¬rda (xN +1; t_j+1) noktas¬ndaki yakla¸s¬m¬ hesaplamak için (xN +1; tj) ve (xN; tj) noktalar¬ndaki yakla¸s¬mlar kullan¬lmaktad¬r.

S¬n¬r nokta içermeyen i = 2; : : : ; N ; j = 1; 2; : : : ; M noktalar¬nda ise

U_i;j+1= rU_{i 1;j} + (1 2r)U_i;j+ rU_i+1;j (8.19) yakla¸s¬mlar¬geçerlidir.

(c) x = 1=4; t = 1=32 için r = t=( x)² = 1=2 elde ederiz. N = 4 adet alt aral¬¼g¬n uç noktalar¬

x₁ = 0; x₂ = 1=4; x₃ = 1=2; x₄ = 3=4; x₅ = 1 dir.

U_i1 = sin( x_i); i = 1; 2; :::; 5 den U₁₁ = 0; U₂₁= ¹₂p

2; U₃₁ = 1; U₄₁= ¹₂p

2; U₅₁ = 0 elde ederiz.

i = 1; j = 1 ve r = 1=2 için (8.17) dan

U₁₂ = U₁₁+ 2r(U₂₁ U₁₁ x)

= U21 x = 1 2

p2 1

4 = 2p

2 1

4 elde ederiz.

(17)

r = 1=2 için (8.18) dan

U_{N +1;j+1} = U_{N +1;j}+ 2r(U_N;j U_{N +1;j})

= U_N;j ve N = 4; j = 1 için

U₅₂= U₄₁= 1 2

p2

dir. Böylecea¸sa¼g¬daki tabloda belirtilen de¼gerleri elde etmi¸s olu- ruz:

U₁₂ = ²^p^{2 1}₄ U₂₂ U₃₂ U₄₂ U₅₂ = ¹₂p 2 U₁₁ = 0 U₂₁ = ¹₂p

2 U₃₁= 1 U₄₁= ¹₂p

2 U₅₁ = 0 x₁ = 0 x₂ = 1=4 x₃ = 1=2 x₄ = 3=4 x₅ = 1

Geriye kalan de¼gerleri ise (8.19) yard¬m¬yla hesaplayabiliriz:

r = 1=2 için (8.19) den

U_i;j+1 = 1

2(U_{i 1;j} + U_i+1;j) elde ederiz. Bu ba¼g¬nt¬ya göre

i = 2; j = 1 için U₂₂ = 1

2(U₁₁+ U₃₁) = 1

2(0 + 1) = 1 2; i = 3; j = 1 için

U₃₂ = 1

2(U₂₁+ U₄₁) = 1 2

1 2

p2 + 1 2

p2 = 1 2

p2;

i = 4; j = 1 için U₄₂= 1

2(U₃₁+ U₅₁) = 1

2(1 + 0) = 1 2 de¼gerlerini elde ederiz.

(18)

t_jnxⁱ 0 0:2 0:4 0:6 0:8 1 0 0 0:5878 0:9511 0:9511 0:5878 0:0000 0:06 0:3173 0:4994 0:6215 0:6520 0:6171 0:5916 0:12 0:2380 0:4082 0:5206 0:5822 0:6076 0:6136 0:18 0:1529 0:3266 0:4497 0:5285 0:5710 0:5843 0:2 0:1269 0:3016 0:4277 0:5103 0:5561 0:5706 Tablo 8.5: x = 0:1; t = 0:0025; r = 1=4 için yakla¸s¬mlar

x_i 0 0:2 0:4 0:6 0:8 1:0

x = 0:1 0:1269 0:3016 0:4277 0:5103 0:5561 0:5706 x = 0:05 0:1302 0:3049 0:4310 0:5136 0:5594 0:5739 x = 0:0250 0:1310 0:3057 0:4318 0:5144 0:5602 0:5747 x = 0:0125 0:1312 0:3059 0:4320 0:5146 0:5604 0:5749 x = 0:00625 0:1313 0:3060 0:4321 0:5147 0:5605 0:5750 Tablo 8.6: Farkl¬ t ve x (r = t=( x)² = 1=4) ile t = 0:2 noktas¬nda yakla¸s¬mlar

(d) x = 0:1; t = 0:0025 ve dolay¬s¬yla r = t=( x)² = 1=4 için elde edilen yakla¸s¬mlardan seçilen baz¬ (xi; t_j) noktalar¬na kar¸s¬l¬k gelenler Tablo 8.5 de verilmektedir

(e) r = t=( x)² = 1=4sabit oran¬yla fakat farkl¬ t ve xler için t = 0:2 noktas¬nda elde edilen yakla¸s¬mlar ise Tablo 8.6 da sunulmaktad¬r.

Tablo 8.6 dan görülece¼gi üzere t; x ler t=( x)² = 1=4 olmak

¸sart¬yla s¬f¬ra yakla¸st¬kça belirtilen xi noktalar¬ndaki çözümler belirli de¼gerlere yak¬nsamaktad¬rlar. Ancak küçük t de¼geri ile t = 0:2 noktas¬na ula¸smak için çok say¬da iterasyon gerekmektedir. Örne¼gin Tablo 8.6 daki son sat¬r yakla¸s¬mlar¬n¬elde etmek için yakla¸s¬k olarak

T

t = T

( x)²=4 = 0:2

(0:000625)²=4 = 20480 iterasyon gereklidir.

(f) x = 0:4 noktas¬nda t = ( x²)=4 seçene¼gi ile farkl¬( x; t) de¼gerlerine kar¸s¬l¬k gelen hata de¼gerleri Tablo 8.7 de verilmektedir. Tablodan

(19)

x; t x = 0:4; t = 0:2 oran 0:1; 0:5 10 ¹ 0:0044

0:05; 2:5 10 ² 0:0011 ^0:00444_0:00111 = 4

0:025; 1:2 5 10 ² 0:2999 10 ³ _{0:29999 10}^0:00111 3 = 3: 700 1 0:0125; 6: 25 10 ³ 0:0999 10 ³ ^{0:29999 10}_{0:09999 10} ³3 = 3: 322 2 0:00625; 3: 125 10 ³ 0

Tablo 8.7: x = 0:4 noktas¬nda olu¸san hata de¼gerleri ve hata oranlar¬

¸

Sekil 8.4: Örnek 8.2 için say¬sal çözüm gra…¼gi

t = ( x²)=4 seçene¼gi ile ( x=2 , t=2) de¼gerleriyle elde edilen yakla¸s¬mlarda olu¸san hatan¬n ( x, t) ile elde edilen hatan¬n yakla¸s¬k olarak dörtte biri kadar oldu¼gu görülmektedir. Bu sonuç t = ( x²)=4 için

E_kum(t; t; x) = O( t) + O( x²) = O( x²) kümülatif hatas¬ile uyumludur.

(g) ( x = 0:1; t = 0:025)ile elde edilen say¬sayal çözüm gra…¼gi ¸Sekil 8.4 te verilmektedir.

ux(0; t) = 1; ux(1; t) = 0

s¬n¬r ¸sartlar¬n¬n çözüm üzerindeki etkisine dikkat edelim. ux(1; t) = 0

¸sart¬, sa¼g uç noktas¬n¬n ¬s¬yal¬t¬ml¬oldu¼gunu ifade etmektedir. ux(0; t) =

(20)

1s¬n¬r ¸sart¬ise sol uç noktadan ¬s¬kayb¬oldu¼gunu ifade etmektedir. Bu durum ¸Sekil 8.4 den de aç¬kça görülmektedir.

8.2.1 Parabolik denklem için aç¬k yöntemin kararl¬l¬¼g¬

Örnek 8.1(f) deki sonucun aksine, matematiksel teoriye uygun ve …ziksel anlaml¬çözümler için tad¬m uzunlu¼gunun en fazla ne kadar olmas¬gerekti¼gini ara¸st¬ral¬m.

Baya¼g¬ diferensiyel denklemler ile olu¸sturulan ba¸slang¬ç de¼ger problem- lerinin kararl¬l¬¼g¬n¬

y0 = ay; y(t¹) = y₁; a > 0 (8.20) model problemi üzerinde çal¬¸st¬¼g¬m¬z¬hat¬rlayal¬m. Problemin analitik çözü- münün t ! 1 için s¬f¬ra yakla¸st¬¼g¬n¬ve say¬sal çözümün de benzer davran¬¸s¬

sergilemesi gerekti¼ginden hareketle, örne¼gin ileri Euler yönteminde ad¬m uzunlu¼gunun h < 2=a k¬s¬tlamas¬n¬sa¼glamas¬gerekti¼gini gözlemlemi¸stik. Benzer biçimde

Y0 = f(Y ) 2 Rⁿ; Y (0) = Y₀; n 2

biçimindeki sistemler için ise f (Y ) = 0 ¸sart¬n¬sa¼glayan asimtotik kararl¬Y denge noktas¬kom¸sulu¼gunda elde edilen

U⁰ = J U; U (0) = U₀

lineer sistemini model problem olarak inceliyoruz. Burada J sistemin denge noktas¬ndaki Jacobien matrisidir ve özde¼gerlerinin reel k¬s¬mlar¬ negatiftir.

Denge noktas¬n¬n yak¬n kom¸sulu¼gundaki çözümlerin s¬f¬ra yakla¸sma özel- li¼gini, yani U^(j) ! 0; j ! 1; özelli¼gini say¬sal çözümden de bekleriz. Bu durumda uygulanan herhangi bir say¬sal yöntem ile elde edilen

U^(j+1) = AU^(j); U⁽¹⁾ = U (0); j = 1; 2; :::

iterasyonunda A iterasyon matrisi için Gershgorin yöntemiyle tahmin edilen (A) < 1 kriterinin gerekli ve yeterli oldu¼gunu hat¬rlayal¬m.

K¬smi diferensiyel denklemler için uygulanan sonlu fark yöntemlerinin kararl¬l¬k analizi farkl¬ yöntemlerle yap¬labilmektedir. Öncelikle biz (8.5)- (8.7) problemine, baya¼g¬diferensiyel denklem sistemlerindeki kararl¬l¬k teo- risini uygulamak istiyoruz. Bu amaçla (8.5) denklemini (8.7) s¬n¬r ¸sart¬n¬

(21)

da dikkate alarak sadece yer de¼gi¸skenine göre merkezi fark yöntemiyle ayr¬k- la¸st¬rmak suretiyle

dU dt = c

x²AU + G(t) + S(t) (8.21)

ba¸slang¬ç-de¼ger sistemini elde ederiz. Burada

A = 0 BB BB

@

2 1 0 ::: 0

1 2 1 0

::: ::: :::

1 2 1

0 ::: 0 1 2

1 CC CC

A; G(t) = 2 66 66 64

G₂(t) G₃(t)

... G_N ₁(t)

G_N(t) 3 77 77 75

; S(t) = 2 66 66 64

c (t)= x² 0

... 0 c (t)= x²

3 77 77 75

ve Gi(t) = g(x_i; t) dir.

Yukar¬da tan¬mlanan aç¬k yöntem (8.21) sistemine ileri Euler yönteminin uygulanmas¬yla elde edilir:

U^(j+1) = U^(j)+ t c

x²AU^(j)+ G(tj) + S(tj) ; j = 1; 2; :::

= (I + rA)U^(j)+ t(G(t_j) + S(t_j)) (8.22) U⁽¹⁾ = U (0) := [U₂(0); U₃(0); :::; U_N(0)]^T;

burada

U_i(0) = U (x_i; 0); i = 2; 3; :::; N

dir. Di¼ger bir deyimle, (8.10),(8.11) ve (8.12) sonlu fark sistemi ve (8.22) sistemi birbirine denktir. Zaman de¼gi¸skenine göre uygulanacak olan fark denkleminin kararl¬l¬¼g¬, herhangi bir t an¬nda olu¸san hatan¬n kontrollü birikimini sa¼glamay¬amaçlayan bir kriterdir.

Bu amaçla hata içeren bir gU⁽¹⁾ = ]U (0) ba¸slang¬ç ¸sart¬ile j inci ad¬mda elde edilen yakla¸s¬m¬ gU^(j) ile gösterelim. j inci ad¬mdaki kümülatif hatay¬

E^(j)= U^(j) Ug^(j)

ile gösterelim. gU^(j) ler (8.22) sistemini gU⁽¹⁾ ba¸slang¬ç ¸sart¬yla sa¼glarlar:

U^^(j+1) = Ug^(j)+ t c

x²A gU^(j)+ G(t_j) + S(t_j) ; j = 1; 2; :::

= (I + rA) gU^(j)+ t(G(t_j) + S(t_j)) (8.23) Ug⁽¹⁾ = ]U (0)