ode45 with discrete data

14 次查看(过去 30 天)
My problem consists of a differential equation in which, the damping coefficient of my system is not linear, it depends on the displacement
where, on the x axis is the displacement of the body and on the y axis is the damping coefficient for these positions.
I did the code in matlab and he did the calculation, however, I'm not sure if he did it correctly, since, in the damping coefficient curve, there is no coefficient for all positions.
An important point is, that I used the slipne function to get the equation, so, I wanted your help to understand a little more about the equation
And my question is more for, I know that the displacement of the system will be within the range that is in the graph of the damping coefficient, but does matlab understand that? why, although I haven't provided all the damping coefficient values, will the matlab integrate only in the range that I provided?
The other point is, I plotted the spline function and, from very high displacements, for my system, the damping coefficient starts to have negative values, so, I want to understand how the matlab can read this function and if it solves my problem in the range provided.
function testando
k1 = 116.0665; % N/m
m1 = 0.06; % kg
w = 25.1327; % rad/s
F1 = 0.01;% N
f1 = @(time) F1.*sin(w*time);
x = xlsread('emf.xlsx','B2:B62');
FT = xlsread('emf.xlsx','J2:J62');
x0 = 0; v0 = 0;
IC = [x0,v0];
% Time
t0 = 0; tf = 10;
tspan = [t0,tf];
%sdot = g(t,s)
sdot1 = @(t,s) [s(2);
(f1(t) - spline(x,FT,s(1)).*s(2) - k1.*s(1))./m1];
% Integração numerica
[time,state_values] = ode45(sdot1,tspan,IC);
displacement = state_values(:,1);
% plots
plot(time,displacement),xlabel('time(s)'),ylabel('displacement(m)')
title('Displacement(m) vs. t')
end

采纳的回答

Ameer Hamza
Ameer Hamza 2020-4-15
For your first question, yes, the spline function only uses the displacement value from vector 'x' the coefficient value from vector 'FT'. As long as the displacement value s(1) lies in the range of 'x', the MATLAB will use values in the range you provided.
For the second question, you are correct. The spline function can overshoot from the actual values. I recommend using interp1, which uses linear interpolation by default. The interpolated value will exactly look like the figure you attached.
sdot1 = @(t,s) [s(2);
(f1(t) - interp1(x,FT,s(1)).*s(2) - k1.*s(1))./m1];
  5 个评论
José Anchieta de Jesus Filho
Thank you very much, your codes helped to understand a little more about ode45.
I tested the first function I did, without any conditions, and then I tested the function using the "if and else" conditions. The result was the same.
Really, matlab understood that only the interval I provided was necessary.
Thank you so much for this.
function testando
k1 = 116.0665; % N/m
m1 = 0.06; % kg
w = 25.1327; % rad/s
F1 = 0.01;% N
f1 = @(time) F1.*sin(w*time);
x = xlsread('emf.xlsx','B2:B62');
FT = xlsread('emf.xlsx','K2:K62');
x0 = 0; v0 = 0;
IC = [x0,v0];
% Time
t0 = 0; tf = 10;
tspan = [t0,tf];
%sdot = g(t,s)
sdot1 = @(t,s) [s(2);
(f1(t) - interp1(x,FT,s(1),'spline').*s(2) - k1.*s(1))./m1];
% Integração numerica
[time,state_values] = ode45(sdot1,tspan,IC);
displacement = state_values(:,1);
% plots
figure(2)
plot(time,displacement),xlabel('time(s)'),ylabel('displacement(m)')
title('Displacement(m) vs. t')
end
The code using the conditions
function sdot = odeFun(t,s,x,Ft)
k1 = 116.0665; % N/m
m1 = 0.06; % kg
w = 25.1327; % rad/s
F1 = 0.01;% N
f1 = @(time) F1.*sin(w*time);
if -0.0302 < s(1) && s(1) < 0.0298
c = spline(x,Ft,s(1));
else
c = 0;
end
sdot(1) = s(2);
sdot(2) = (f1(t) - c.*s(2) - k1.*s(1))./m1;
sdot = sdot';
end
Ameer Hamza
Ameer Hamza 2020-4-15
I am glad to be of help.

请先登录,再进行评论。

更多回答(0 个)

类别

Help CenterFile Exchange 中查找有关 Splines 的更多信息

Community Treasure Hunt

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

Start Hunting!

Translated by