finding the value at a specific time from second-order ODE

41 次查看（过去 30 天）
Tien Yin
Tien Yin 2024-12-2，17:08
评论： Tien Yin 2024-12-4，2:20
Ran in:
I'm new to MATLAB
I've been working on solving ODE and I now have learnt two ways to plot the 2-order ODE
  1. using ode45
  2. using conv() fucntion
Here's the problem : I want to know the value at a specific time for example f(2.5) when t=2.5
it can be displayed in any way like output in the command window or on the figure
please teach a way to output the value
I'll attach the files I've been working on
if you think the codes are tltr you could just give suggestions though
[time,y] = ode45('ode45_function' , [0 10] , [0,0]) ;
figure ;
subplot(212);plot(time,y(:,1));xlabel('t');ylabel('y(t)');
m = 5;
k = 18;
c = 1.2;
f0 = 100;
wn = sqrt(k/m);
cc = 2*m*wn;
Dr = c/cc;
wd = wn*sqrt(1-Dr^2);
dt = 0.001; % assumed dt
time = 0:dt:10;
f = zeros(1, length(time));
f(0001:3000) = (f0/3) * (0:dt:3-dt); % Linear increase from 0 to 3 seconds
f(3001:5000) = f0; % Constant from 3 to 5 seconds
% f remains 0 after 5 seconds
g = (1/(m*wd))*exp(-Dr*wn*time).*sin(wd*time);
x = conv(f, g)*dt;
x = x(1:length(time)); % Ensure x is the same length as time
% x is the solution function for the 2-order ODE　
figure;
subplot(221); plot(time, f, 'r'); xlabel('t(s)'); ylabel('f(t)');
subplot(222); plot(time, g); xlabel('t(s)'); ylabel('g(t)');
subplot(223); plot(time, x); xlabel('t(s)'); ylabel('位移(m)');
% this is the function which I want to extract the specific value from
% it for ex: f(3), f(2.5) ...
function ydot=ode45_function(t,y)
m=5;
k=18;
c=1.2;
if t<3
f= 100*t/3;
elseif t<5
f=100;
else
f=0;
end
ydot(1)=y(2);
ydot(2)=(1/m)*(f-c*y(2)-k*y(1));
ydot=ydot';
end

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

采纳的回答

Steven Lord
Steven Lord 2024-12-2，18:44
Ran in:
You could call ode45 with just one output argument and then call deval, or you could call ode45 with the desired time as one element of the timespan input argument, or (if you're using a sufficiently recent release of MATLAB) you could use the ode object instead of just ode45.
I also recommend you pass a function handle into ode45 rather than passing the name of the function.
One output with deval
sol = ode45(@ode45_function , [0 10] , [0,0]) ;
solutionAtTime = deval(sol, 2.5)
solutionAtTime = 2×1
5.2273 1.9079
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
Timespan input contains the desired time
[t, y] = ode45(@ode45_function , [0 2.5 10] , [0,0]);
solutionAtTime2 = y(2, :).'
solutionAtTime2 = 2×1
5.2273 1.9079
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
ODE object
odeobj = ode;
odeobj.ODEFcn = @ode45_function;
odeobj.InitialValue = [0 0];
results = solve(odeobj, 2.5); % Just solve at t = 2.5
solutionAtTime3 = results.Solution
solutionAtTime3 = 2×1
5.2267 1.9092
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
  6 个评论
显示 4更早的评论
Steven Lord
Steven Lord 2024-12-3，15:53
Yes, if the first component of y represents position then the second component represents velocity.

请先登录，再进行评论。

更多回答（1 个）

Torsten
Torsten 2024-12-2，18:04
移动：Torsten 2024-12-2，18:04
The easiest (maybe slightly inaccurate) way is to use interpolation:
f_conv = interp1(time,x,2.5)
f_ode45 = interp1(time,y(:,1),2.5)

类别

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

标签

Community Treasure Hunt

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

Start Hunting!

Translated by