Plotting 2nd order ode with time dependent parameters

4 visualizzazioni (ultimi 30 giorni)
I need to use the ODE45 function to plot the following ODE (shown below) with the initial conditions𝜃(0)=𝜋/2 rad and 𝜃(𝑡)=0.1 rad/s and 0 ≤ 𝑡 ≤ 10 seconds. One plot for𝜃(𝑡) when 𝛺=𝜋 and another when 𝛺=4𝜋.
theta'' + theta' + 10*cos(2pi*theta) = sin(omega*t)
Im having difficulty understanding how to create a function that can solve this using ode45. I was able to solve earlier ones when it was a simple linear 2nd ode- Homogeneous. Below is how I solved the other 2nd order.
function = second_order
%% y'' + 6y' + 9y = 0 initial function
yinitial = [4, 0];
t = 0:0.001:3;
[T, STATES] = ode45(@ode2_hw8, t, yinitial);
Y = STATES(:,1);
Ydot = STATES(:,2);
figure
subplot(1,2,1)
plot(T,Y)
subplot(1,2,2)
plot(T, Ydot)
end
function [out2] = ode2_hw8 (t, state)
b = 6;
k = 9;
y = state(1);
ydot = state(2);
yddot = -b*ydot - k*y;
out2 = [ydot; yddot];
end

Risposte (1)

Ashutosh Singh Baghel
Ashutosh Singh Baghel il 20 Ott 2021
Hi Jake,
I believe you wish to perform a second order differential ODE with time dependent parameter. Please refer to the following example -
% Second Order ODE with Time-Dependent Terms
% $$ y''(t) +y'(t) +10*cos(2*pi*y) = sin(w*t);
% w = omega = pi, 4*pi, etc
% t = (0,10)
% The initial condition is $y_0 = [pi*0.5 , 0.1]$.
omega = pi;
gt = linspace(0, 10, 50);
g = sin(pi*gt);
tspan = [0 10];
y_0 = [pi*0.5 0.1];
[t,y] = ode45(@(t, y) myode1(t, y, gt, g), tspan, y_0);
plot(t, y(:,1), '-r ', t, y(:,2), '--g');
legend('y(1)', 'y(2)');
xlabel('t');
ylabel('y Solution');
function dydt = myode1(t, y, gt, g)
g = interp1(gt, g, t); % Interpolate the data set (gt, g) at time t
dydt = [y(2); -(y(2)+10*cos(2*pi*y(1))) + g]; % Evaluate ODE at time
end
Please refer to the following MATLAB Documentation page on the 'ode45' function and relevant information on the 'Differential Equation'.

Categorie

Scopri di più su Programming in Help Center e File Exchange

Prodotti


Release

R2020b

Community Treasure Hunt

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

Start Hunting!

Translated by