Next: Function
da 21 a 30 Up:Introduzione
Previous: Function
da 1 a 10
Intro Gen: Introduzione
Generale Home: Home
page
% METODO DI EULERO ESPLICITO PER LA RISOLUZIONE DELLE
EQUAZIONI DIFFERENZIALI
function y=eulespl(F,t,y,n,h);
% y
vettore soluzione dell'equazione differenziale
% F
é la funzione a due variabili F(t,y)
% n
é il numero di suddivisioni dell'intervallo [a,b] fissato
% h
é l'ampiezza delle n suddivisioni
for
i=2:n+1
y(i)=y(i-1)+h*feval(F,t(i-1),y(i-1));
end
% METODO DI EULERO IMPLICITO PER LA RISOLUZIONE DELLE
EQUAZIONI DIFFERENZIALI
function [y,risp]=eulimpl(F,G,t,y,n,h,e);
% y
vettore soluzione dell'equazione differenziale
% risp
controlla se il metodo é applicabile
% F
é la funzione a due variabili F(t,y)
% G
é la derivata di F
% t
é la variabile indipendente della funzione y
% n
é il numero di suddivisioni dell'intervallo [a,b] fissato
% h
é l'ampiezza delle n suddivisioni
% e
é l'errore massimo che si deve commettere
risp=0;
i=2;
while i<=n+1 & risp==0
sol1=y(i-1);
wd=feval(G,t(i),sol1);
if
wd==1/h
risp=1;
else
wf=sol1-y(i-1)-h*(feval(F,t(i),sol1));
sol2=sol1-wf/(1-h*wd);
end
while risp==0
& (abs(sol2-sol1))>=e | (abs(wf)>=e)
sol1=sol2;
wd=feval(G,t(i),sol1);
if wd==1/h
risp=1;
else
wf=sol1-y(i-1)-h*feval(F,t(i),sol1);
sol2=sol1-wf/(1-h*wd);
end
end
y(i)=sol2;
i=i+1;
end
% FATTORIZZAZIONE QR
function [Q,R,risp]=fattqr(A,m,n);
% Q ,R
sono le matrici ottenute dalla fattorizzazione
% A
matrice da fattorizzare
% n ,m
dimensioni della matrice
risp=0;
R=A;
Q=eye(m);
if risp==0
for
k=1:n
c=sign(R(k,k));
if c==0
c=1;
end
a=0;
% calcolo della norma del vettore
R(k:m,k)
for i=k:m
a=a+(R(i,k))^2;
end
a=sqrt(a);
%
calcolo delle matrici Q(k) e R(k)
b=a*(abs(R(k,k))+a);
if b==0
risp=1;
else
b=1/b;
u=1;
u(1)=R(k,k)+c*a;
u(2:m-k+1)=R(k+1:m,k);
P=eye(m-k+1,m-k+1)-b*u'*u;
S=0;
if k>1
S(1:k-1,1:k-1)=eye(k-1,k-1);
end
S(k:m,k:m)=P;
R(k:m,k:n)=P*R(k:m,k:n);
Q=S*Q;
end
end
else
Q=0;
R=0;
end
Q=Q';
% METODO DI GAUSS
function [x,risp]=gauss(A,b,n);
% x
vettore soluzione del sistema
% risp
controlla se il metodo é applicabile
% A
vettore dei coefficienti del sistema
% b
vettore dei termini noti
% n
dimensione di A
risp=0;
if b==zeros(n,1);
x=zeros(n,1);
else
k=1;
while k<n
if A(k,k)==0
risp=1;
k=n;
end
if risp==0
for i=k+1:n
for j=k+1:n
L=A(i,k)/A(k,k);
A(i,j)=A(i,j)-L*A(k,j);
end
b(i)=b(i)-L*b(k);
for j=1:k
A(i,j)=0;
end
end
end
k=k+1;
end
% APPLICAZIONE DEL METODO DI GAUSS PER LA RICERCA DEL VETTORE
X:
if risp==0
x(n)=b(n)/A(n,n);
for i=n-1:-1:1
a=0;
for j=i+1:n
a=a+A(i,j)*x(j);
end
x(i)=(b(i)-a)/A(i,i);
end
end
end
if risp==1
x=0;
end
% INTERPOLAZIONE DI HERMITE
function P=intherm(x,y,z,xi,n,m);
% P vettore
soluzione dei punti da interpolare
% x vettore
dei nodi interpolanti
% y valori
della funzione nei nodi x
% z valori
della derivata prima della funzione nei nodi x
% xi vettore
dei punti in cui calcolare la funzione interpolante
% n dimensione
di x , y e z
% m dimensione
di xi
for ki=1:m
p=0;
for k=1:n
a=1;
c=1;
for i=1:k-1
c=c*(x(k)-x(i));
a=a*(xi(ki)-x(i));
end
for i=k+1:n
c=c*(x(k)-x(i));
a=a*(xi(ki)-x(i));
end
a=a/c;
% Calcolo
della derivata prima dei polinomi fondamentali
l=0;
for j=1:n
if j~=k
d=1;
for i=1:n
if
i~=k
if i~=j
d=d*(x(k)-x(i));
end
end
end
l=l+d;
end
end
l=l/c;
p=p+a*a*(y(k)+(z(k)-2*l*y(k))*(xi(ki)-x(k)));
end
P(ki)=p;
end
P=P';
% INTERPOLAZIONE DI LAGRANGE
function L=intlagr(x,y,xi,n,m);
% L vettore
soluzione dei punti da interpolare
% x vettore
dei nodi interpolanti
% y valori
della funzione nei nodi x
% xi vettore
dei punti in cui calcolare la funzione interpolante
% n dimensione
di x , y
% m dimensione
di xi
q=1;
% ordinamento dei nodi
while q==1
q=0;
for i=1:n-1
if x(i)>x(i+1)
u=x(i);
x(i)=x(i+1);
x(i+1)=u;
u=y(i);
y(i)=y(i+1);
y(i+1)=u;
q=1;
end
end
end
for i=1:m
% interpolazione di Lagrange
L(i)=0;
for k=1:n
a=1;
for j=1:k-1
a=a*((xi(i)-x(j))/(x(k)-x(j)));
end
for j=k+1:n
a=a*((xi(i)-x(j))/(x(k)-x(j)));
end
L(i)=L(i)+a*y(k);
end
end
L=L';
asc=x(1):.1:x(n);
% Algoritmo seguente utilizzato
[h,m2]=size(asc);
% unicamente per tracciare il grafico
for i=1:m2
ord(i)=0;
for k=1:n
a=1;
for j=1:k-1
a=a*((asc(i)-x(j))/(x(k)-x(j)));
end
for j=k+1:n
a=a*((asc(i)-x(j))/(x(k)-x(j)));
end
ord(i)=ord(i)+a*y(k);
end
end
plot(asc,ord,'w',xi,L,'+r',x,y,'og');
title('Interpolazione di Lagrange');
% APPROSSIMAZIONE POLINOMIALE AI MINIMI QUADRATI
function [p,risp]=intminqu(x,y,xi,n,m,w);
% p
vettore soluzione dei punti da interpolare
% risp controlla
se il metodo é applicabile
% x
vettore dei nodi interpolanti
% y
valori della funzione nei nodi x
% xi
vettore dei punti in cui calcolare la funzione interpolante
% n
dimensione di x , y
% w
dimensione di xi
% m
grado del polinomio interpolante
B=diag(eye(n));
% determinazione della matrice
for i=1:n
% dei coefficienti del sistema
for j=2:m
% da risolvere
B(i,j)=(x(i))^(j-1);
end
end
[a,f,risp]=minquad(B,y,n,m);
if risp==0
a=a';
for i=1:w
% regola di Horner per valutare
p(i)=0;
% il polinomio interpolante
p(i)=valpol(xi(i),a);
% nei nodi
end
asc=x(1):.1:x(n);
% Algoritmo seguente utilizzato
[h,m2]=size(asc);
% unicamente per tracciare il grafico
for i=1:m2
ord(i)=0;
for j=m:-1:1
ord(i)=ord(i)*asc(i)+a(j);
end
end
plot(asc,ord,'w',xi,p,'+r',x,y,'og');
title('Approssimazione
polinomiale');
else
p=0;
end
% INTERPOLAZIONE CON SCHEMA DI NEVILLE
function z=intnevil(x,y,xi,n,m);
% z
vettore soluzione dei punti da interpolare
% x
vettore dei nodi interpolanti
% y
valori della funzione nei nodi x
% xi vettore
dei punti in cui calcolare la funzione interpolante
% n
dimensione di x , y
% m dimensione
di xi
for k=1:m
L=0;
for i=1:n
L(i,1)=y(i);
end
for j=2:n
for i=j:n
L(i,j)=(L(i,j-1)*(xi(k)-x(i-j+1))-L(i-1,j-1)*(xi(k)-x(i)))/(x(i)-x(i-j+1));
end
end
z(k)=L(n,n);
end
asc=x(1):.1:x(n);
% Algoritmo seguente utilizzato
[h,m2]=size(asc);
% unicamente per tracciare il grafico
for k=1:m2
L=0;
for i=1:n
L(i,1)=y(i);
end
for j=2:n
for i=j:n
L(i,j)=(L(i,j-1)*(asc(k)-x(i-j+1))-L(i-1,j-1)*(asc(k)-x(i)))/(x(i)-x(i-j+1));
end
end
ord(k)=L(n,n);
end
plot(asc,ord,'w',xi,z,'+r',x,y,'og');
title('Interpolazione
con schema di Neville');
% INTERPOLAZIONE DI NEWTON
function L=intnewt(x,y,xi,n,m);
% L
vettore soluzione dei punti da interpolare
% x
vettore dei nodi interpolanti
% y
valori della funzione nei nodi x
% xi vettore
dei punti in cui calcolare la funzione interpolante
% n
dimensione di x , y
% m dimensione
di xi
for i=1:n
F(i,1)=y(i);
end
for k=1:m
for j=2:n
for i=j:n
F(i,j)=(F(i,j-1)-F(i-1,j-1))/(x(i)-x(i-j+1));
end
end
L(k)=F(1,1);
w=1;
for i=2:n
w=w*(xi(k)-x(i-1));
L(k)=L(k)+F(i,i)*w;
end
end
asc=x(1):.1:x(n);
% Algoritmo seguente utilizzato
[h,m2]=size(asc);
% unicamente per tracciare il grafico
for k=1:m2
for j=2:n
for i=j:n
F(i,j)=(F(i,j-1)-F(i-1,j-1))/(x(i)-x(i-j+1));
end
end
ord(k)=F(1,1);
w=1;
for i=2:n
w=w*(asc(k)-x(i-1));
ord(k)=ord(k)+F(i,i)*w;
end
end
plot(asc,ord,'w',xi,L,'+r',x,y,'og');
title('Interpolazione di Newton');
% METODO DI JACOBI PER LA RISOLUZIONE DI SISTEMI LINEARI
function [x,risp]=jacobi(A,b,n);
% x
vettore soluzione del sistema
% risp
controlla se il metodo é applicabile
% A
vettore dei coefficienti del sistema
% b
vettore dei termini noti
% n
dimensione di A
risp=0;
for i=1:n
D(i,i)=1/A(i,i);
end
disp('Per
una verifica molto precisa sulla convergenza del metodo iterativo ');
disp('adottato
introdurre un valore abbastanza elevato del seguente numero ');
D=eye(n)-D*A;
l=autmax(D,n);
if abs(l)>=1
risp=1;
x=0;
else
Kmax=input('Introduci
il numero massimo di iterate per applicare il metodo di Jacobi : ');
err=0;
while err==0
if
Kmax<=0 | round(Kmax)~=Kmax
disp(' ');
disp(' Il valore introdotto non é
corretto');
Kmax=input('Introducilo nuovamente : ');
else
err=1;
end
end
x=diag(eye(n));
for k=1:Kmax
for i=1:n
a=0;
for j=1:i-1
a=a-A(i,j)*x(j);
end
for j=i+1:n
a=a-A(i,j)*x(j);
end
x(i)=(b(i)+a)/A(i,i);
end
end
end
Next: Function
da 21 a 30 Up:Introduzione
Previous: Function
da 1 a 10
Intro Gen: Introduzione
Generale Home: Home
page
Vito Marinelli
8-5-2000