Giovanni Di Grezia
d65fd0e0a9
sebbene funzioni, le trasformazioni sulle matrici sono troppo onerose e servirebbe un algoritmo piu efficiente
100 lines
2.7 KiB
Matlab
100 lines
2.7 KiB
Matlab
function [ area, err, ier ] = integratore_adattivo_globale( f, m, toll , maxite )
|
|
%algoritmo dei trapezi composto, addattivo con controllo di errore globale
|
|
% input
|
|
% f = function handle
|
|
% m = matrice n x 2 ordinata in modo crescente per righe rispetto alla
|
|
% prima colonna. la prima colonna rappresenta un punto x, la seconda
|
|
% colonna la f(x). ogni valore della seconda colonna deve essere uguale a
|
|
% 0 o essere la f(x) del punto x della prima colonna sulla medesima riga.
|
|
% toll = tolleranza assoluta voluta
|
|
% maxite = numero massimo di iterazioni
|
|
% output
|
|
% area = area finale
|
|
% err = errore assoluto
|
|
% ier = indice di errore
|
|
|
|
[x,~] = size(m);
|
|
i=1;
|
|
k=i+1;
|
|
j=x;
|
|
ite = 0;
|
|
ier = 0;
|
|
report=[];
|
|
while (i <= j - 1) && (ite < maxite)
|
|
[~,larea,lerr,m1] = static_quadratic(f,m(i:k,:));
|
|
m=[m(1:i-1,:);m1;m(k+1:j,:)];
|
|
ite = ite + 1;
|
|
j= j + 1;
|
|
report=[report;i,k+1,larea,lerr];
|
|
i = i + 2;
|
|
k = k + 2;
|
|
end
|
|
report = sortrows(report,4);
|
|
report = flipud (report);
|
|
[lenght,~] = size(report);
|
|
area=0;
|
|
err = 0;
|
|
for i=1:lenght
|
|
err= err + report(i,4);
|
|
area = area + report(i,3);
|
|
end
|
|
|
|
while ( err > toll)&& (ite < maxite)
|
|
i=report(1,1);
|
|
k=report(1,2) - 1;
|
|
[lenght,~] = size(report);
|
|
report=report(2:lenght,:);
|
|
[~,larea,lerr,m1] = static_quadratic(f,m(i:k,:));
|
|
m=[m(1:i-1,:);m1;m(k+1:j,:)];
|
|
ite = ite + 1;
|
|
j= j + 1;
|
|
report=[report;i,k+1,larea,lerr];
|
|
report = sortrows(report,1);
|
|
[lenght,~] = size(report);
|
|
for i=1:lenght
|
|
if(report(i,1) >= k+1)
|
|
report(i,1) = report(i,1) + 1;
|
|
report(i,2) = report(i,2) + 1;
|
|
end
|
|
end
|
|
i=k+1;
|
|
k=k+2;
|
|
[~,larea,lerr,m1] = static_quadratic(f,m(i:k,:));
|
|
m=[m(1:i-1,:);m1;m(k+1:j,:)];
|
|
ite = ite + 1;
|
|
j= j + 1;
|
|
report=[report;i,k+1,larea,lerr];
|
|
report = sortrows(report,1);
|
|
[lenght,~] = size(report);
|
|
for i=1:lenght
|
|
if(report(i,1) >= k)
|
|
report(i,1) = report(i,1) + 1;
|
|
report(i,2) = report(i,2) + 1;
|
|
end
|
|
end
|
|
report = sortrows(report,4);
|
|
report = flipud (report);
|
|
[lenght,~] = size(report);
|
|
area=0;
|
|
err = 0;
|
|
for i=1:lenght
|
|
err= err + report(i,4);
|
|
area = area + report(i,3);
|
|
end
|
|
end
|
|
|
|
if ite >= maxite
|
|
disp('Warning: Massimo numero di step raggiunti')
|
|
ier = 1;
|
|
|
|
end
|
|
|
|
|
|
hold on
|
|
fplot (f,[m(1,1) m(j,1)],'k')
|
|
plot (m(:,1),m(:,2),'*r')
|
|
hold off
|
|
|
|
end
|
|
|