Plotting 2nd order ode with time dependent parameters

4 ビュー (過去 30 日間)
Jake-Anthony Vickers
Jake-Anthony Vickers 2021 年 10 月 12 日
回答済み: Ashutosh Singh Baghel 2021 年 10 月 20 日
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

回答 (1 件)

Ashutosh Singh Baghel
Ashutosh Singh Baghel 2021 年 10 月 20 日
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'.

カテゴリ

Help Center および File ExchangeProgramming についてさらに検索

タグ

製品


リリース

R2020b

Community Treasure Hunt

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

Start Hunting!

Translated by