Ways for faster torque interpolation
5 次查看(过去 30 天)
显示 更早的评论
I am modelling a single stage helical gear transmission, in which I am applying an input torque to the driver gear. However, the torque interpolation in my code takes the longest time and hence a 60 sec simulation takes more than 12 hours to run. Is there a way to interpolate the torque faster or through a more efficient method?
My main code-
function [yp] = reduced_t(t,y,T_a)
beta = 13*(pi/180); %Helix angle (rad)
P = 20*(pi/180); %Pressure angle (rad)
% speed = 1000/60;
Freq = 1000*20/60;
% Common parameters
% e_a = zeros(size(t)); %circumferential displacement of the driver gear (m)
% e_p = zeros(size(t)); %circumferential displacement of the driver gear (m)
R_a = 51.19e-3; %Radius
R_p = 139.763e-3; %Radius
e_a = (pi*2*R_a*tan(beta))/(2*cos(P));
e_p = (pi*2*R_p*tan(beta))/(2*cos(P));
% for i = 2:length(t)
% % theta_a_vec(i) = theta_a_vec(i-1) - speed*2*pi*(t(i) - t(i-1)); % iterative assignment
% e_a(i) = e_a(i-1) + theta_a_vec(i)*R_a;
% e_p(i) = e_p(i-1) - theta_a_vec(i)*R_p;
% end
e = 0;
ks = 0.9e8 + 20000*sin(2*pi*Freq*t); %Shear stiffness
Cs = 0.1 + 0.001*sin(2*pi*Freq*t); %Shear damping
m_a = 0.5;
I_a = 0.5*m_a*(R_a^2); %Moment of inertia of driver gear (Kg.m3)
% Driven gear
m_p = 0.5; %mass of driven gear (kg)
I_p = 0.5*m_p*(R_p^2); %Moment of inertia of driven gear (Kg.m3)
n_a = 20; % no of driver gear teeth
n_p = 60; % no of driven gear teeth
i = n_p/n_a; % gear ratio
% y_ac= e_a + theta_a_vec*R_a; %circumferential displacement at the meshing point on the driver gear (m)
% y_pc = e_p - theta_a_vec*R_p; %circumferential displacement at the meshing point on the driven gear (m)
yp = zeros(4,1); %vector of 4 equations
% Excitation forces
% Fy=ks*cos(beta)*(y_ac-y_pc-e) + Cs*cos(beta)*2*R_a*speed*2*pi; %axial dynamic excitation force at the meshing point (N)
% Fz=ks*sin(beta)*(z_ac-z_pc-z) - Cs*sin(beta)*2*tan(beta)*R_a*yp(3); %circumferential dynamic excitation force at the meshing point (N)
%Time interpolation
time = 0:0.00001:0.06;
Torque = interp1(time,T_a,t);
% Torque = spline(time,T_a,t);
Opp_Torque = 68.853; % Average torque value
% %Time interpolation using parfor loop
% time = 0:0.00001:0.6;
% Torque = zeros(size(t));
% parfor i=1:length(t)
% Torque(i) = interp1(time,T_a,t(i));
% end
% Opp_Torque = 68.853; % Average torque value
%Driver gear equations
yp(1) = y(2);
yp(2) = (Torque- Cs*cos(beta)*R_a*(R_a*y(2) + R_p*y(4)) - ks*cos(beta)*R_a*(R_a*y(1)+R_p*y(3)) -ks*cos(beta)*R_a*(e_a-e_p-e))/I_a;
%Driven gear equations
yp(3) = y(4);
yp(4) = (Opp_Torque*i - Cs*cos(beta)*R_p*(R_a*y(2) + R_p*y(4)) - ks*cos(beta)*R_p*(R_a*y(1)+ R_p*y(3)) -ks*cos(beta)*R_p*(e_a-e_p-e))/I_p;
end
My command window-
tic
A = load("torque_for_Sid_60s.mat");
T_a = A.torque_60s(1:6001); %Torque (Nm)
% speed = 1000/60;
Freq = 1000*20/60;
tspan=0:0.00001:0.06;
y0 = [0;104.719;0;104.719/3];
[t,y] = ode45(@(t,y) reduced_t(t,y,T_a),tspan,y0);
% TE = theta_p*R_p - theta_a_vec*R_a; %transmission error
toc
%Driver gear graphs
nexttile
plot(t,y(:,2))
ylabel('angular velocity (rad/s)')
xlabel('Time')
title('Driver gear angular velocity vs time')
axis padded
grid on
% nexttile
% plot(t,y(:,1))
% ylabel('angular velocity (rad/s)')
% xlabel('Time')
% title('Driver gear angular velocity vs time')
% axis padded
% grid on
% Driven gear graphs
nexttile
plot(t,y(:,4))
ylabel('angular velocity (rad/s )')
xlabel('Time')
title('Driven gear angular velocity vs time')
axis padded
grid on
% nexttile
% plot(t,y(:,3))
% ylabel('angular velocity (rad/s)')
% xlabel('Time')
% title('Driven gear angular velocity vs time')
% axis padded
% grid on
% nexttile
% plot(t,T_a(1:601))
% ylabel('Torque (Nm)')
% xlabel('Time (sec)')
% title('Torque vs time')
% axis padded
% grid on
My torque file - torque_for_Sid_60s.mat
0 个评论
回答(1 个)
Jan
2023-3-9
interp1 is slower than griddedInterpolant. You can created the interpolant outside the integration to save more time.
A sever problem remains: The linear interpolation is not smooth (deifferentable), but Matlab's ODE integrators are defined for smooth functions only. So you code can confuse the stepsize control. On one hand this can reduce the stepsizes dratsically and increase the computation times. On the other hand, the huge number of function evaluations accumulates the rounding errors. In consequence, the result can be dominated by rounding errors, such that the shown code is a poor random number generator.
Use a smooth function (e.g. a cubic interpolation).
4 个评论
Jan
2023-3-9
Does the first do, what you want? You cancompare it with the output of interp1. If so:
time = 0:0.00001:0.06;
torque_interp = griddedInterpolant(time, T_a, 'pchip'); % or 'spline'
[t,y] = ode45(@(t,y) reduced_t(t,y, torque_interp),tspan,y0);
function [yp] = reduced_t(t,y,torque_interp)
...
Torque = torque_interp(t);
end
另请参阅
产品
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!