Duffing equation:Transition to Chaos

4 次查看（过去 30 天）
Athanasios Paraskevopoulos
编辑： Swastik Sarkar 2024-8-29
Ran in:
The Original Equation is the following:
Let . This implies that
Then we rewrite it as a System of First-Order Equations
Using the substitution for , the second-order equation can be transformed into the following system of first-order equations:
Exploring the Effect of γ.
% Define parameters
gamma = 0.338;
alpha = -1;
beta = 1;
delta = 0.1;
omega = 1.4;
% Define the system of equations
odeSystem = @(t, y) [y(2);
-delta*y(2) - alpha*y(1) - beta*y(1)^3 + gamma*cos(omega*t)];
% Initial conditions
y0 = [0; 0]; % x(0) = 0, v(0) = 0
% Time span
tspan = [0 2000];
% Solve the system
[t, y] = ode45(odeSystem, tspan, y0);
% Plot the results
figure;
plot(t, y(:, 1));
xlabel('Time');
ylabel('x(t)');
title('Solution of the nonlinear system');
grid on;
% Plot the phase portrait
figure;
plot(y(:, 1), y(:, 2));
xlabel('x(t)');
ylabel('v(t)');
title('Phase-Plane $$\ddot{x}+\delta \dot{x}+\alpha x+\beta x^3=0$$ for $$\gamma=0.318$$, $$\alpha=-1$$,$$\beta=1$$,$$\delta=0.1$$,$$\omega=1.4$$','Interpreter', 'latex');
grid on;
% Define the tail (e.g., last 10% of the time interval)
tail_start = floor(0.9 * length(t)); % Starting index for the tail
tail_end = length(t); % Ending index for the tail
% Plot the tail of the solution
figure;
plot(y(tail_start:tail_end, 1), y(tail_start:tail_end, 2), 'r', 'LineWidth', 1.5);
xlabel('x(t)');
ylabel('v(t)');
title('Phase-Plane $$\ddot{x}+\delta \dot{x}+\alpha x+\beta x^3=0$$ for $$\gamma=0.318$$, $$\alpha=-1$$,$$\beta=1$$,$$\delta=0.1$$,$$\omega=1.4$$','Interpreter', 'latex');
grid on;
ax = gca;
chart = ax.Children(1);
datatip(chart,0.5581,-0.1126);
Then I used but the results were not that I expected for
. The paper I study is Duffing Equation
My code gives me the following. Any suggestion?
% Define parameters
gamma = 0.35;
alpha = -1;
beta = 1;
delta = 0.1;
omega = 1.4;
% Define the system of equations
odeSystem = @(t, y) [y(2);
-delta*y(2) - alpha*y(1) - beta*y(1)^3 + gamma*cos(omega*t)];
% Initial conditions
y0 = [0; 0]; % x(0) = 0, v(0) = 0
% Time span
tspan = [0 3000];
% Solve the system
[t, y] = ode45(odeSystem, tspan, y0);
% Plot the phase portrait
figure;
plot(y(:, 1), y(:, 2));
xlabel('x(t)');
ylabel('v(t)');
title('Phase-Plane $$\ddot{x}+\delta \dot{x}+\alpha x+\beta x^3=0$$ for $$\gamma=0.318$$, $$\alpha=-1$$,$$\beta=1$$,$$\delta=0.1$$,$$\omega=1.4$$','Interpreter', 'latex');
grid on;
% Define the tail (e.g., last 10% of the time interval)
tail_start = floor(0.9 * length(t)); % Starting index for the tail
tail_end = length(t); % Ending index for the tail
% Plot the tail of the solution
figure;
plot(y(tail_start:tail_end, 1), y(tail_start:tail_end, 2), 'r', 'LineWidth', 1.5);
xlabel('x(t)');
ylabel('v(t)');
title('Phase-Plane $$\ddot{x}+\delta \dot{x}+\alpha x+\beta x^3=0$$ for $$\gamma=0.35$$, $$\alpha=-1$$,$$\beta=1$$,$$\delta=0.1$$,$$\omega=1.4$$','Interpreter', 'latex');
grid on;
ax = gca;
chart = ax.Children(1);
datatip(chart,0.5581,-0.1126);
  3 个评论
显示 1更早的评论
Athanasios Paraskevopoulos
Athanasios Paraskevopoulos 2024-8-16
Ran in:
Hello @Neelanshu
I tried
% Time span with more points for better resolution
tspan = linspace(0, 200,3000); % Increase the number of points
I feel that I get closer but still I can not create the other two plots yet
% Define parameters
gamma = 0.35;
alpha = -1;
beta = 1;
delta = 0.1;
omega = 1.4;
% Define the system of equations
odeSystem = @(t, y) [y(2);
-delta*y(2) - alpha*y(1) - beta*y(1)^3 + gamma*cos(omega*t)];
% Initial conditions
y0 = [0; 0]; % x(0) = 0, v(0) = 0
% Time span with more points for better resolution
tspan = linspace(0, 200,3000); % Increase the number of points
% Solve the system with increased output refinement
options = odeset('Refine', 4); % Additional refinement of the output
[t, y] = ode45(odeSystem, tspan, y0, options);
% Plot the results
figure;
plot(t, y(:, 1));
xlabel('Time');
ylabel('x(t)');
title('Solution of the nonlinear system');
grid on;
% Plot the phase portrait
figure;
plot(y(:, 1), y(:, 2));
xlabel('x(t)');
ylabel('v(t)');
title('Phase-Plane $$\ddot{x}+\delta \dot{x}+\alpha x+\beta x^3=0$$ for $$\gamma=0.318$$, $$\alpha=-1$$,$$\beta=1$$,$$\delta=0.1$$,$$\omega=1.4$$','Interpreter', 'latex');
grid on;
% Define the tail (e.g., last 10% of the time interval)
tail_start = floor(0.9 * length(t)); % Starting index for the tail
tail_end = length(t); % Ending index for the tail
% Plot the tail of the solution
figure;
plot(y(tail_start:tail_end, 1), y(tail_start:tail_end, 2), 'r', 'LineWidth', 1.5);
xlabel('x(t)');
ylabel('v(t)');
title('Phase-Plane $$\ddot{x}+\delta \dot{x}+\alpha x+\beta x^3=0$$ for $$\gamma=0.318$$, $$\alpha=-1$$,$$\beta=1$$,$$\delta=0.1$$,$$\omega=1.4$$','Interpreter', 'latex');
grid on;
% Adding a datatip for visualization (optional)
ax = gca;
chart = ax.Children(1);
datatip(chart,0.5581,-0.1126);

请先登录，再进行评论。

请先登录，再回答此问题。

采纳的回答

Swastik Sarkar
Swastik Sarkar 2024-8-29
编辑：Swastik Sarkar 2024-8-29
Ran in:
Hi @Athanasios Paraskevopoulos,
Upon reviewing the paper linked in the question (https://www.colorado.edu/amath/sites/default/files/attached-files/good_sample_project_0.pdf), it specifies a time range of [0, 3000] and focuses on the initial and final few points for its plots, whereas the current solution plots the range [0, 200] with 3000 points, as indicated usinglinspace(0, 200, 3000) upon the entire duration. So, changing the code to the below one will give you desired plots.
% Define parameters
gamma = 0.35;
alpha = -1;
beta = 1;
delta = 0.1;
omega = 1.4;
% Define the system of equations
odeSystem = @(t, y) [y(2);
-delta*y(2) - alpha*y(1) - beta*y(1)^3 + gamma*cos(omega*t)];
% Initial conditions
y0 = [0; 0]; % x(0) = 0, v(0) = 0
tspan = [0 3000]; % Increase the number of points
% Solve the system with increased output refinement
[t, y] = ode45(odeSystem, tspan, y0);
% Define the tail (e.g., last 10% of the time interval)
extreme_tail_start = floor(0.9966 * length(t)); % Starting index for the tail
% Define the tail (e.g., last 10% of the time interval)
tail_start = floor(0.989 * length(t)); % Starting index for the tail
tail_end = length(t); % Ending index for the tail
% Define the head (e.g., first 10% of the time interval)
head_start = 1; % Starting index for the tail
head_end = floor(0.1 * length(t)); % Ending index for the tail
% Plot the phase portrait
figure
plot(y(head_start:head_end, 1), y(head_start:head_end, 2), 'LineWidth', 1.5);
xlabel('x(t)');
ylabel('v(t)');
title('Phase-Plane $$\ddot{x}+\delta \dot{x}+\alpha x+\beta x^3=0$$ for $$\gamma=0.318$$, $$\alpha=-1$$,$$\beta=1$$,$$\delta=0.1$$,$$\omega=1.4$$','Interpreter', 'latex');
grid on;
% Plot the tail of the solution
figure
plot(y(tail_start:tail_end, 1), y(tail_start:tail_end, 2), 'LineWidth', 1.5);
xlabel('x(t)');
ylabel('v(t)');
title('Phase-Plane $$\ddot{x}+\delta \dot{x}+\alpha x+\beta x^3=0$$ for $$\gamma=0.318$$, $$\alpha=-1$$,$$\beta=1$$,$$\delta=0.1$$,$$\omega=1.4$$','Interpreter', 'latex');
grid on;
% Plot the extreme tail of the solution
figure
plot(y(extreme_tail_start:tail_end, 1), y(extreme_tail_start:tail_end, 2), 'LineWidth', 1.5);
xlabel('x(t)');
ylabel('v(t)');
title('Phase-Plane $$\ddot{x}+\delta \dot{x}+\alpha x+\beta x^3=0$$ for $$\gamma=0.318$$, $$\alpha=-1$$,$$\beta=1$$,$$\delta=0.1$$,$$\omega=1.4$$','Interpreter', 'latex');
grid on;
% Adding a datatip for visualization (optional)
ax = gca;
chart = ax.Children(1);
datatip(chart,0.5581,-0.1126);

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Line Plots 的更多信息

标签

产品


版本

R2024a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by