آخرين ارسالهاي تالار

خطا: mod_kunenalatest:كيوننا نسخه 1.7 (يا بالاتر) بر روي سيستم شما نصب نيست!
پیغام
  • Kunena is not installed or the installed Kunena version is not supported. The plug-in has now been disabled. Please install/upgrade Kunena to version 1.7 for the Kunena Discuss Plug-in to function properly.
جاذب لورنز lorenz attractor مشاهده در قالب PDF چاپ فرستادن به ایمیل
نوشته شده توسط saeed   
شنبه, 07 فروردین 1389 ساعت 23:53

جاذب لورنز lorenz attractor

در این بخش برنامه جاذب لورنز را با استفاده از دو روش رانگ کوتای مرتبه4 و ode45 با توجه به توضیحات داده شده مینویسم.

جاذب لورنز یک سیستم آشوبناک است که برای نختسن بار توسط ادوارد لورنز از انستیتو تکنولوژی ماساچوست معرفی شد.

http://mathworld.wolfram.com/LorenzAttractor.html

http://en.wikipedia.org/wiki/Lorenz_attractor

http://ozviz.wasp.uwa.edu.au/~pbourke/fractals/lorenz

سیستم معادلات لورنز عبارتند از

نمی خواهم در مورد نحوه استخراج این معادلات بحث کنم و تنها به این بسنده می­کنم که اساس آنها مباحث فیزیکی در ارتباط با جابجایی اتمسفراست. متغیرهای x، y و z نشان دهنده کمیت ها فیزیکی همچون دما و سرعت سیال هستند و پارامترها معرفی شده مربوط به ویژگیهای جو می باشند که با داده های زیر معین می شوند،

در ابتدا این معادلات را با استفاده از روش رانگ کوتا حل می کنیم. تابع Lorenz_rk اینکار را انجام می دهد. حاصل کار بصورت یک انیمیشن نمایش داده می­شود.

function lorenz_rk

%-----------------

N=5000;

sigma = 10;

rho = 25;

beta = 8/3.;

t=linspace(0,50,N); % time between 0-5s

h=t(2)-t(1);   % time step

f=@(x,y) sigma*(y-x);

g=@(x,y,z) x*(rho-z)-y;

u=@(x,y,z) x*y-beta*z;

x(1)=1;

y(1)=0;

z(1)=0;

for i=1:length(t)-1

kx1=f(x(i),y(i));

ky1=g(x(i),y(i),z(i));

kz1=u(x(i),y(i),z(i));

%------------------------

kx2=f(x(i)+h/2*kx1,y(i)+h/2*ky1);

ky2=g(x(i)+h/2*kx1,y(i)+h/2*ky1,z(i)+h/2*kz1);

kz2=u(x(i)+h/2*kx1,y(i)+h/2*ky1,z(i)+h/2*kz1);

%------------------------

kx3=f(x(i)+h/2*kx2,y(i)+h/2*ky2);

ky3=g(x(i)+h/2*kx2,y(i)+h/2*ky2,z(i)+h/2*kz2);

kz3=u(x(i)+h/2*kx2,y(i)+h/2*ky2,z(i)+h/2*kz2);

%------------------------

kx4=f(x(i)+h*kx3,y(i)+h*ky3);

ky4=g(x(i)+h*kx3,y(i)+h*ky3,z(i)+h*kz3);

kz4=u(x(i)+h*kx3,y(i)+h*ky3,z(i)+h*kz3);

%------------------------

x(i+1)=x(i)+h/6*(kx1+2*kx2+2*kx3+kx4);

y(i+1)=y(i)+h/6*(ky1+2*ky2+2*ky3+ky4);

z(i+1)=z(i)+h/6*(kz1+2*kz2+2*kz3+kz4);

end

body = line( ...

'color','b', ...

'LineStyle','none', ...

'Marker','.', ...

'erase','none', ...

'xdata',[],'ydata',[],'zdata',[]);

axis([ min(y)-eps max(y)+eps...

min(z)-eps max(z)+eps])

for j=2:N

set(body,'XData',y(j),'YData',z(j))

drawnow

end

plot(y,z,'.')

xlabel('y')

ylabel('z')

 

 

 

با استفاده از تابع ode45 این سیستم را میتوانیم بشکل زیر حل کنیم که با نام Lorenz_ode مشخص شده است.

function Lorenz_ode

close all

x0=1;

y0=0;

z0=0;

N=10000;

[t,Y] = ode45(@lorenz_at,linspace(0,50,N),[x0 y0 z0]);

body = line( ...

'color','b', ...

'LineStyle','none', ...

'Marker','.', ...

'erase','none', ...

'xdata',[],'ydata',[],'zdata',[]);

axis([ min(Y(:,2))-eps max(Y(:,2))+eps...

min(Y(:,3))-eps max(Y(:,3))+eps])

for j=2:N

set(body,'XData',Y(j,2),'YData',Y(j,3))

drawnow

end

plot(Y(:,2),Y(:,3),'.')

xlabel('y')

ylabel('z')

function dy=lorenz_at(t,y)

dy=zeros(3,1);

SIGMA = 10;

RHO = 25;

BETA = 8/3.;

dy(1)=SIGMA*(y(2)-y(1));

dy(2)=y(1)*(RHO-y(3))-y(2);

dy(3)=y(1)*y(2)-BETA*y(3);

نظر ها (2)
  • Masoud_Eng
    مثل همیشه عالی !!!
  • ناشناس  - راهنمایی
    سلام. تشکر خسته نباشید. یه سوال. چرا متعییرها در متلب نمایش داده نمی شوند.
تنها کاربران عضو شده می توانند نظر ارسال کنند!
آخرین بروز رسانی در سه شنبه, 25 آبان 1389 ساعت 16:15
 
logo-samandehi