Discretize model dx dt = σx − σy, dy dt = xρ − xz − y, dz dt = xy − βz, using Runge-Kutta methods (of order 2 or order 4).
조회 수: 6 (최근 30일)
이전 댓글 표시
Using the code :
clear all
close all
clc
Beta = [10;28;8/3];
x0 = [10;10;10];% I.C
dt = 0.001; % time step size
tspan = 0:dt:20;
figure(1)
options = odeset('RelTol',1e-12,'AbsTol',1e-12*ones(1,3));
[t,x] = ode45(@(t,x)lorenz(t,x,Beta),tspan,x0,options);
plot3(x(:,1),x(:,2),x(:,3),'r','lineWidth',1.5);
set(gca,'color','w','xcolor','r','ycolor','r','zcolor','r')
set(gcf,'color','w')
xlabel('porportional to the rate of convention')
ylabel('horizontal temperature variation')
zlabel('vertical temperature variation')
hold on
options1 = odeset('RelTol',1e-6,'AbsTol',1e-6*ones(1,3));
[t1,x1] = ode45(@(t1,x1)lorenz(t1,x1,Beta),tspan,x0,options1);
plot3(x1(:,1),x1(:,2),x1(:,3),'b','lineWidth',1.5);
set(gca,'color','w','xcolor','b','ycolor','b','zcolor','b')
set(gcf,'color','w')
grid on
legend('1e-12','1e-6')
figure(2)
subplot(3,1,1)
plot(t,x(:,1),'r',t1,x1(:,1),'b')
title('x(t)') %porportional to the rate of convention in respect to time
legend('1e-12','1e-6')
xlabel('time')
ylabel('y(t)') %porportional to the rate of convention
subplot(3,1,2)
plot(t,x(:,2),'r',t1,x1(:,2),'b')
title('y(t)') %horizontal temperature variation in respect to time
ylabel('HA temp var')
subplot(3,1,3)
plot(t,x(:,3),'r',t1,x1(:,3),'b')
title('z(t)') %vertical temperature variation in respect to time
xlabel('time')
ylabel('VA temp var') %vertical temperature variation
find the Runge-Kutta methods (of order 2 or order 4) Use the same parameter values 10;28;8/3 and initial condition x(0) = 10, y(0) = 10, z(0) = 10, run the simulation for 0 ≤ t ≤ 20. Plot the approximated solution using Runge-Kutta methods and the solution using ode45 in the same figure.
댓글 수: 1
답변 (1개)
Karim
2023년 3월 15일
편집: Karim
2023년 3월 15일
After having a look at the wikipedia article (Lorenz system - Wikipedia), it appears that you have some mistakes in your equations that you showed in the title. The equations should be:
- dx/dt = σ *(y-x)
- dy/dt = x*(ρ−z) − y
- dz/dt = x*y − β*z
Based on your script Beta = [10; 28; 8/3], i assume that 'Beta' holds the values for sigma rho and beta? So
- σ = 10
- ρ = 28
- β = 8/3
the part missing from your code is the lorenz function itself. So we need to create a function. An example can be:
% create the anonymous Function
Beta = [10;28;8/3];
lorenz = @(t,x) [ Beta(1)*( x(2)-x(1) ); % this represents dx/dt
x(1)*( Beta(2)-x(3) ) - x(2); % this represents dy/dt
x(1)*x(2) - Beta(3)*x(3) ]; % this represents dz/dt
Now we can solve using the ode45 solver
% initial values
x0 = [1;1;1];
% time span
dt = 1e-3;
tspan = 0:dt:20;
% set the options for the solver
options = odeset('RelTol',1e-12,'AbsTol',1e-12*ones(1,3));
% call the solver using the 'lorenz' function
[t,x] = ode45(lorenz,tspan,x0,options);
% create a figure
figure
plot3(x(:,1),x(:,2),x(:,3),'r','lineWidth',1.5);
xlabel('porportional to the rate of convention')
ylabel('horizontal temperature variation')
zlabel('vertical temperature variation')
grid on
댓글 수: 0
참고 항목
카테고리
Help Center 및 File Exchange에서 Ordinary Differential Equations에 대해 자세히 알아보기
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!