(→lorenz.m) |
|||
| (13 versioni intermedie di uno stesso utente non sono mostrate) | |||
| Riga 1: | Riga 1: | ||
| + | == L'attrattore di Lorenz == | ||
| + | [[Immagine:Lorenz.png|right|thumb|400px|''L'attrattore di Lorenz'']] | ||
| + | Utilizzando i due file per Matlab (testati con la versione 7.0) seguenti è possibile visualizzare una animazione (per t < tf) della posizione nel tempo di N soluzioni all'equazione differenziale di Lorenz. In praticolare le condizioni iniziali a partire dalle quali vengono calcolate le N soluzioni sono contenute in una sfera di raggio 0.1. | ||
| + | |||
| + | Dall'animazione è possibile notare come con l'avanzare del tempo le N soluzioni divergano anche a fronte di una piccolissima variazione delle condizioni iniziali, ma tuttavia mantenendosi confinate sull'attrattore di Lorenz (disegnato in celeste). | ||
| + | |||
| + | ==lorenz.m== | ||
%Equazione differenziale Ordinaria dell'attrattore di lorenz | %Equazione differenziale Ordinaria dell'attrattore di lorenz | ||
%----------------------------------------------------------- | %----------------------------------------------------------- | ||
| − | % dx/dt = -10 x | + | % dx/dt = -10 x + 10 y |
% dy/dt = x - y - (sqrt(27)+x)z | % dy/dt = x - y - (sqrt(27)+x)z | ||
% dz/dt = sqrt(27)*(x+y) - z + xy | % dz/dt = sqrt(27)*(x+y) - z + xy | ||
| Riga 10: | Riga 17: | ||
dy(3) = sqrt(27)*(x(1)+x(2))-x(3)+x(1)*x(2); | dy(3) = sqrt(27)*(x(1)+x(2))-x(3)+x(1)*x(2); | ||
end | end | ||
| + | |||
| + | ==draw_lorenz.m== | ||
| + | function [t,x,y,out] = solve_lorenz(x_0,t_f,n) | ||
| + | y = [ ] ; | ||
| + | x = [ ] ; | ||
| + | plot3(0,0,0); | ||
| + | out = cell(n,1); | ||
| + | hold; | ||
| + | for i=1:n | ||
| + | x_i = rand(1,3)*0.1 + x_0; | ||
| + | x = [x x_i']; | ||
| + | [t,temp]=ode45(@lorenz,[0 t_f], x_i); | ||
| + | out{i} = temp; | ||
| + | y_len = length(temp); | ||
| + | y_f = temp(y_len,:); | ||
| + | y = [ y y_f']; | ||
| + | plot3(x_i(1),x_i(2),x_i(3),'ko'); | ||
| + | plot3(y_f(1),y_f(2),y_f(3),'kx'); | ||
| + | end | ||
| + | min_size = 1e10; | ||
| + | min_x = 1e10 ; | ||
| + | min_y = 1e10 ; | ||
| + | min_z = 1e10 ; | ||
| + | max_x = 0; | ||
| + | max_y = 0; | ||
| + | max_z = 0; | ||
| + | for i=1:n | ||
| + | s = size(out{i}); | ||
| + | x_max = max(out{i}(:,1)); | ||
| + | x_min = min(out{i}(:,1)); | ||
| + | y_max = max(out{i}(:,2)); | ||
| + | y_min = min(out{i}(:,2)); | ||
| + | z_max = max(out{i}(:,3)); | ||
| + | z_min = min(out{i}(:,3)); | ||
| + | if s(1) < min_size | ||
| + | min_size = s(1); | ||
| + | end | ||
| + | if x_max > max_x | ||
| + | max_x = x_max; | ||
| + | end | ||
| + | if y_max > max_y | ||
| + | max_y = y_max; | ||
| + | end | ||
| + | if z_max > max_z | ||
| + | max_z = z_max; | ||
| + | end | ||
| + | if x_min < min_x | ||
| + | min_x = x_min; | ||
| + | end | ||
| + | if y_min < min_y | ||
| + | min_y = y_min; | ||
| + | end | ||
| + | if z_min < min_z | ||
| + | min_z = z_min; | ||
| + | end | ||
| + | end | ||
| + | %Plotta per ogni istante la posizione delle n soluzioni | ||
| + | for j = 1:min_size | ||
| + | clf; | ||
| + | axis([min_x max_x min_y max_y min_z max_z]); | ||
| + | plot3(out{n}(:,1),out{n}(:,2),out{n}(:,3),'c'); | ||
| + | hold on; | ||
| + | for i=1:n | ||
| + | v = out{i}(j,:); | ||
| + | plot3(v(1),v(2),v(3),'ko'); | ||
| + | end | ||
| + | pause(0.1); | ||
| + | drawnow; | ||
| + | hold off; | ||
| + | end | ||
| + | end | ||
| + | |||
| + | == Riferimenti == | ||
| + | |||
| + | [http://mathworld.wolfram.com/LorenzAttractor.html Mathword: Lorenz Attractor] | ||
Utilizzando i due file per Matlab (testati con la versione 7.0) seguenti è possibile visualizzare una animazione (per t < tf) della posizione nel tempo di N soluzioni all'equazione differenziale di Lorenz. In praticolare le condizioni iniziali a partire dalle quali vengono calcolate le N soluzioni sono contenute in una sfera di raggio 0.1.
Dall'animazione è possibile notare come con l'avanzare del tempo le N soluzioni divergano anche a fronte di una piccolissima variazione delle condizioni iniziali, ma tuttavia mantenendosi confinate sull'attrattore di Lorenz (disegnato in celeste).
%Equazione differenziale Ordinaria dell'attrattore di lorenz
%-----------------------------------------------------------
% dx/dt = -10 x + 10 y
% dy/dt = x - y - (sqrt(27)+x)z
% dz/dt = sqrt(27)*(x+y) - z + xy
function dy = lorenz(t,x)
dy = zeros(3,1);
dy(1) = -10*x(1) + 10 * x(2);
dy(2) = x(1) - x(2) -(sqrt(27)+x(1))*x(3);
dy(3) = sqrt(27)*(x(1)+x(2))-x(3)+x(1)*x(2);
end
function [t,x,y,out] = solve_lorenz(x_0,t_f,n)
y = [ ] ;
x = [ ] ;
plot3(0,0,0);
out = cell(n,1);
hold;
for i=1:n
x_i = rand(1,3)*0.1 + x_0;
x = [x x_i'];
[t,temp]=ode45(@lorenz,[0 t_f], x_i);
out{i} = temp;
y_len = length(temp);
y_f = temp(y_len,:);
y = [ y y_f'];
plot3(x_i(1),x_i(2),x_i(3),'ko');
plot3(y_f(1),y_f(2),y_f(3),'kx');
end
min_size = 1e10;
min_x = 1e10 ;
min_y = 1e10 ;
min_z = 1e10 ;
max_x = 0;
max_y = 0;
max_z = 0;
for i=1:n
s = size(out{i});
x_max = max(out{i}(:,1));
x_min = min(out{i}(:,1));
y_max = max(out{i}(:,2));
y_min = min(out{i}(:,2));
z_max = max(out{i}(:,3));
z_min = min(out{i}(:,3));
if s(1) < min_size
min_size = s(1);
end
if x_max > max_x
max_x = x_max;
end
if y_max > max_y
max_y = y_max;
end
if z_max > max_z
max_z = z_max;
end
if x_min < min_x
min_x = x_min;
end
if y_min < min_y
min_y = y_min;
end
if z_min < min_z
min_z = z_min;
end
end
%Plotta per ogni istante la posizione delle n soluzioni
for j = 1:min_size
clf;
axis([min_x max_x min_y max_y min_z max_z]);
plot3(out{n}(:,1),out{n}(:,2),out{n}(:,3),'c');
hold on;
for i=1:n
v = out{i}(j,:);
plot3(v(1),v(2),v(3),'ko');
end
pause(0.1);
drawnow;
hold off;
end
end
%Equazione differenziale Ordinaria dell'attrattore di lorenz
%-----------------------------------------------------------
% dx/dt = -10 x - 10 y
% dy/dt = x - y - (sqrt(27)+x)z
% dz/dt = sqrt(27)*(x+y) - z + xy
function dy = lorenz(t,x)
dy = zeros(3,1);
dy(1) = -10*x(1) + 10 * x(2);
dy(2) = x(1) - x(2) -(sqrt(27)+x(1))*x(3);
dy(3) = sqrt(27)*(x(1)+x(2))-x(3)+x(1)*x(2);
end