Metoda Strzału

Równania różniczkowe i całkowe. Równania różnicowe. Transformata Laplace'a i Fouriera oraz ich zastosowanie w równaniach różniczkowych.
janusz47
Użytkownik
Użytkownik
Posty: 7910
Rejestracja: 18 mar 2009, o 16:24
Płeć: Mężczyzna
Podziękował: 30 razy
Pomógł: 1670 razy

Metoda Strzału

Post autor: janusz47 »

Metoda Strzału

Rozpatrujemy zagadnienie brzegowe:

\(\displaystyle{ \begin{cases} y''+ p(x)y' +q(x)y = r(x) \\ y(a) = \alpha, \ \ y(b)=\beta \ \ (1) \end{cases} }\)

Rozwiązanie zagadnienia \(\displaystyle{ (1) }\) Metodą Strzału polega na zamianie zagadnienia brzegowego na kombinację liniową dwóch zagadnień początkowych drugiego rzędu:

\(\displaystyle{ u''+ p(x)u' +q(x)u = r(x) \\ u(a) = \alpha, \ \ u'(\alpha)= 0 \ \ (2)}\)

\(\displaystyle{ v''+ p(x)v' +q(x)v = r(x) \\ v(a) = 0, \ \ v'(a)= 1 \ \ (3) }\)

\(\displaystyle{ y(x) = u(x)+A v(x) \ \ (4) }\)

Na podstawie \(\displaystyle{ (2), (3), (4), (1) }\)

\(\displaystyle{ u''(x)+ Av''(x)+p(x)[u'(x)+Av"(x)] + q(x)[u(x)+ Av(x)]=r(x) }\)

\(\displaystyle{ u''(x) +p(x)u'(x)+q(x)u(x) +A[v''(x)+p(x)v'(x)+q(x)v(x)]=r(x) }\)

Stałą \(\displaystyle{ A }\) znajdujemy z warunków brzegowych:

\(\displaystyle{ y(a) = u(a)+Av(a) = \alpha +0 = \alpha, }\)

\(\displaystyle{ y(b) = u(b) +Av(b) = \beta }\)

\(\displaystyle{ A = \frac{\beta - u(b)}{v(b)}.}\)

Jeśli \(\displaystyle{ v(b) \neq 0 }\) to rozwiązanie zagadnienia brzegowego \(\displaystyle{ (1) }\) ma postać:

\(\displaystyle{ y(x) = u(x) + \frac{\beta - u(b)}{v(b)}v(x).}\)

Do rozwiązań zagadnień początkowych \(\displaystyle{ (2), (3) }\) stosujemy Metodę Runge-Kutta 4 rzędu.

Z warunkami istnienia rozwiązań tej metody można zapoznać się w podręczniku:

DAVID KINCAID WARD CHENEY analiza numeryczna. WNT Warszawa 2006.

Program w MATLAB

Kod: Zaznacz cały

 
function mstrzal(f1,f2,a,b,alfa,beta,n) 
%Rozwiązanie zagadnienia początkowego drugiego rzędu Metodą Strzałów. 
% f1= -p(x)*u-q(x)*v+r(x), f2= -p(x)*u-q(x)*v. 
h=(b-a)/n; 
y1=alfa; 
y2=0; 
u=0; 
v=alfa; 
for i=1:n 
% Runge-Kutta rzędu 4 
 x=a+(i-1)*h; 
 k1=feval(f1,x,u,v); 
 c1=u; 
 k2=feval(f1,x+h/2,u+h*k1/2,v+h*c1/2); 
 c2=u+h/2*k1; 
 k3=feval(f1,x+h/2,u+h*k2/2,v+h*c2/2); 
 c3=u+h/2*k2; 
 k4=feval(f1,x+h,u+h*k3,v+h*c3); 
 c4=u+h*k3; 
 u=u+h*(k1+2*k2+2*k3+k4)/6; 
 v=v+h*(c1+2*c2+2*c3+c4)/6; 
 y1(i+1)=v; 
end 
u=1; 
v=0; 
for i=1:n 
% Runge-Kutta rzędu 4 
 x=a+(i-1)*h; 
 k1=feval(f2,x,u,v); 
 c1=u; 
 k2=feval(f2,x+h/2,u+h*k1/2,v+h*c1/2); 
 c2=u+h/2*k1; 
 k3=feval(f2,x+h/2,u+h*k2/2,v+h*c2/2); 
 c3=u+h/2*k2; 
 k4=feval(f2,x+h,u+h*k3,v+h*c3); 
 c4=u+h*k3; 
 u=u+h*(k1+2*k2+2*k3+k4)/6; 
 v=v+h*(c1+2*c2+2*c3+c4)/6; 
 y2(i+1)=v; 
 end 
fprintf('\n') 
disp(' Metoda Strzału ') 
fprintf('\n') 
disp([' u(b) = ',num2str(y1(n+1)),' ']) 
disp([' v(b) = ',num2str(y2(n+1)),' ']) 
fprintf('\n') 
disp('_______________________________________________') 
disp('   xi         ui          vi           yi      ') 
disp('_______________________________________________') 
for i=1:n+1 
 x=a+(i-1)*h; 
 w=y1(i)+(beta-y1(n+1))/y2(n+1)*y2(i); 
  s='n'; 
 if (s=='n') 
 fprintf('%6.2f %12.6f %12.6f %12.6f\n',x,y1(i),y2(i),w) 
 else 
 err=abs(w-s); 
 fprintf('%6.2f %12.6f %12.6f %12.6f %12.6f %10.2e\n', x, y1(i),y2(i), w, s, err) 
 end 
end 
Przykład

Kod: Zaznacz cały

function f=f1(x,u,v)
f=3/x*u-3/x.^2*v+2*x.^2*exp(x);
 
function f=f2(x,u,v)
f=3/x*u-3/x.^2*v;

mstrzal('f1','f2',1,2,0,4*exp(2),10)

 Metoda Strzału

 u(b) = 13.2463 
 v(b) = 3 

_______________________________________________
   xi         ui          vi           yi      
_______________________________________________
  1.00     0.000000     0.000000     0.000000
  1.10     0.032983     0.115500     0.660913
  1.20     0.158380     0.264001     1.593649
  1.30     0.423715     0.448502     2.862042
  1.40     0.888400     0.672003     4.541813
  1.50     1.625684     0.937504     6.722522
  1.60     2.724900     1.248006     9.509811
  1.70     4.294042     1.606508    13.027984
  1.80     6.462735     2.016010    17.422977
  1.90     9.385635     2.479512    22.865754
  2.00    13.246340     3.000015    29.556224
Ostatnio zmieniony 13 sie 2022, o 13:31 przez Jan Kraszewski, łącznie zmieniany 1 raz.
Powód: Poprawa wiadomości.
ODPOWIEDZ