Solving Duffing Equation with the new framework for ODEs
20 views (last 30 days)
Show older comments
Athanasios Paraskevopoulos
on 17 Aug 2024 at 14:34
Commented: Athanasios Paraskevopoulos
on 17 Aug 2024 at 16:09
I solved the duffing equation using the new framework for ODEs. Below is the code.
It works fine and plots as expected. However, I was wondering if there is a way to create a phase portrait in the interval t=[0 3000] without looking like a smudge.
% Define the parameters
delta = 0.1;
alpha = -1;
beta = 1;
gamma = 0.35;
omega = 1.4;
% Define the ODE function for the Duffing equation
duffingODE = @(t, y) [y(2);
-delta*y(2) - alpha*y(1) - beta*y(1)^3 + gamma*cos(omega*t)];
% Initial conditions: [x(0), dx/dt(0)]
initialConditions = [0; 0];
% Create an ode object
F = ode(ODEFcn=duffingODE, InitialTime=0, InitialValue=initialConditions);
% Solve the equation over the interval [0, 3000]
sol = solve(F, 0, 3000);
% Interpolate the solution to get more points for plotting
timeFine = linspace(0, 300, 10000); % Create a fine time vector with 10,000 points
solutionFine = interp1(sol.Time, sol.Solution', timeFine)'; % Interpolate solution
% Plot the interpolated time series solution
figure;
subplot(2, 1, 1);
plot(timeFine, solutionFine(1, :), 'LineWidth', 1.5);
xlabel('Time');
ylabel('Displacement');
title('Interpolated Solution of the Duffing Equation');
grid on;
% Plot the interpolated phase portrait
subplot(2, 1, 2);
plot(solutionFine(1, :), solutionFine(2, :), 'LineWidth', 1.5);
xlabel('Displacement x(t)');
ylabel('Velocity dx/dt');
title('Interpolated Phase Portrait of the Duffing Equation');
grid on;
0 Comments
Answers (1)
Steven Lord
on 17 Aug 2024 at 15:41
FYI rather than interpolating the solution yourself with this code:
% Interpolate the solution to get more points for plotting
timeFine = linspace(0, 300, 10000); % Create a fine time vector with 10,000 points
solutionFine = interp1(sol.Time, sol.Solution', timeFine)'; % Interpolate solution
you can call solve on the ODE object with the timeFine vector directly, rather than just the start and end times. From the solve documentation page: "S = solve(F,t) computes the solution for the ODE represented by F at the specified time values in the vector t."
But to your main question: However, I was wondering if there is a way to create a phase portrait in the interval t=[0 3000] without looking like a smudge.
that's just a plain old line. Do you have a picture that shows the type of plot that you want to see, one that doesn't look "like a smudge"? Do one of the thumbnails on this documentation page look more like what you want?
1 Comment
See Also
Categories
Find more on Ordinary Differential Equations in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!