Out of memory when using ode45(@t,x)

4 次查看(过去 30 天)
function xdot=r2dof(t,x,ths,spec,Kpid)
xdot=zeros(8,1);
%% set-points
th1s=ths(1);
th2s=ths(2);
%% Robot Specifications
M1=spec(3);
M2=spec(4);
L1=spec(1);
L2=spec(2);
g=9.8;
%% Inertia Matrix
b11=(M1+M2)*L1^2+M2*L2^2+2*M2*L1*L2*cos(x(4));
b12=M2*L2^2+M2*L1*L2*cos(x(4));
b21=M2*L2^2+M2*L1*L2*cos(x(4));
b22=M2*L2^2;
Bq=[b11 b12;b21 b22];
%% C Matrix
c1=-M2*L1*L2*sin(x(4))*(2*x(5)*x(6)+x(6)^2);
c2=-M2*L1*L2*sin(x(4))*x(5)*x(6);
Cq=[c1;c2];
%% Gravity Matrix
g1=-(M1+M2)*g*L1*sin(x(3))-M2*g*L2*sin(x(3)+x(4));
g2=-M2*g*L2*sin(x(3)+x(4));
Gq=[g1;g2];
%% PID Control
% PID parameters for theta 1
Kp1=Kpid(1);
Kd1=Kpid(2);
Ki1=Kpid(3);
% PID parameters for theta 2
Kp2=Kpid(4);
Kd2=Kpid(5);
Ki2=Kpid(6);
%decoupled control input
f1=Kp1*(th1s-x(3))-Kd1*x(5)+Ki1*(x(1));
f2=Kp2*(th2s-x(4))-Kd2*x(6)+Ki2*(x(2));
Fhat=[f1;f2];
F=Bq*Fhat; % actual input to the system
%% System states
xdot(1)=(th1s-x(3)); %dummy state of theta1 integration
xdot(2)=(th2s-x(4)); %dummy state of theta2 integration
xdot(3)=x(5); %theta1-dot
xdot(4)=x(6); %theta2-dot
q2dot=inv(Bq)*(-Cq-Gq+F);
xdot(5)=q2dot(1); %theta1-2dot
xdot(6)=q2dot(2); %theta1-2dot
%control input function output to outside computer program
xdot(7)=F(1);
xdot(8)=F(2);
2) System Solution and Simulation (“r2dof_cntrl.m”)
close all
clear all
clc
%% Initilization
th_int=[-pi/2 pi/2]; %initial positions
ths=[pi/2 -pi/2]; %set-points
x0=[0 0 th_int 0 0 0 0]; %states initial values
Ts=[0 20]; %time span
%% Robot Specifications
L1=1; %link 1
L2=1; %link 2
M1=1; %mass 1
M2=1; %mass 2
spec=[L1 L2 M1 M2];
%% PID Parameters
% PID parameters for theta 1
Kp1=15;
Kd1=7;
Ki1=10;
% PID parameters for theta 2
Kp2=15;
Kd2=10;
Ki2=10;
Kpid=[Kp1 Kd1 Ki1 Kp2 Kd2 Ki2];
%% ODE solving
% opt1=odeset('RelTol',1e-10,'AbsTol',1e-20,'NormControl','off');
[T,X] = ode45(@(t,x) r2dof(t,x,ths,spec,Kpid),Ts,x0);
//////////////////////////////////////////////////////
ERROR
Out of memory. The likely cause is an infinite recursion within the program.
Error in r2dof (line 75)
[T,X] = ode45(@(t,x) r2dof(1, [0 0 th_int 0 0 0 0], [pi/2 -pi/2], [1 1 1 1],[1 1 1 1 1 1]),[0 20],1);
  3 个评论
mohammed tifr
mohammed tifr 2022-12-10
first thank you for your answer then
the modified code
[T,X] = ode45(@(t,x) r2dof(t,x, [pi/2 -pi/2], [1 1 1 1],[1 1 1 1 1 1]),[0 20],1);
Errors(see the capture)
Torsten
Torsten 2022-12-10
编辑:Torsten 2022-12-10
Works for me.
But your call to ode45 is wrong.
[T,X] = ode45(@(t,x)r2dof(t,x,[pi/2 -pi/2],[1 1 1 1],[1 1 1 1 1 1 ]),tspan,x0)
would be correct.
r2dof(1,[0 0 [-pi/2 pi/2] 0 0 0 0 ],[pi/2 -pi/2],[1 1 1 1],[1 1 1 1 1 1 ])
ans = 8×1
3.1416 -3.1416 0 0 -6.6584 6.6584 6.2832 0
function xdot=r2dof(t,x,ths,spec,Kpid)
xdot=zeros(8,1);
%% set-points
th1s=ths(1);
th2s=ths(2);
%% Robot Specifications
M1=spec(3);
M2=spec(4);
L1=spec(1);
L2=spec(2);
g=9.8;
%% Inertia Matrix
b11=(M1+M2)*L1^2+M2*L2^2+2*M2*L1*L2*cos(x(4));
b12=M2*L2^2+M2*L1*L2*cos(x(4));
b21=M2*L2^2+M2*L1*L2*cos(x(4));
b22=M2*L2^2;
Bq=[b11 b12;b21 b22];
%% C Matrix
c1=-M2*L1*L2*sin(x(4))*(2*x(5)*x(6)+x(6)^2);
c2=-M2*L1*L2*sin(x(4))*x(5)*x(6);
Cq=[c1;c2];
%% Gravity Matrix
g1=-(M1+M2)*g*L1*sin(x(3))-M2*g*L2*sin(x(3)+x(4));
g2=-M2*g*L2*sin(x(3)+x(4));
Gq=[g1;g2];
%% PID Control
% PID parameters for theta 1
Kp1=Kpid(1);
Kd1=Kpid(2);
Ki1=Kpid(3);
% PID parameters for theta 2
Kp2=Kpid(4);
Kd2=Kpid(5);
Ki2=Kpid(6);
%decoupled control input
f1=Kp1*(th1s-x(3))-Kd1*x(5)+Ki1*(x(1));
f2=Kp2*(th2s-x(4))-Kd2*x(6)+Ki2*(x(2));
Fhat=[f1;f2];
F=Bq*Fhat; % actual input to the system
%% System states
xdot(1)=(th1s-x(3)); %dummy state of theta1 integration
xdot(2)=(th2s-x(4)); %dummy state of theta2 integration
xdot(3)=x(5); %theta1-dot
xdot(4)=x(6); %theta2-dot
q2dot=inv(Bq)*(-Cq-Gq+F);
xdot(5)=q2dot(1); %theta1-2dot
xdot(6)=q2dot(2); %theta1-2dot
%control input function output to outside computer program
xdot(7)=F(1);
xdot(8)=F(2);
end

请先登录,再进行评论。

回答(2 个)

Jan
Jan 2022-12-10
编辑:Jan 2022-12-10
I've copied your code without any changes (but adding en end at the bottom) and it is working:
%% Initilization
th_int=[-pi/2 pi/2]; %initial positions
ths=[pi/2 -pi/2]; %set-points
x0=[0 0 th_int 0 0 0 0]; %states initial values
Ts=[0 20]; %time span
%% Robot Specifications
L1=1; %link 1
L2=1; %link 2
M1=1; %mass 1
M2=1; %mass 2
spec=[L1 L2 M1 M2];
%% PID Parameters
% PID parameters for theta 1
Kp1=15;
Kd1=7;
Ki1=10;
% PID parameters for theta 2
Kp2=15;
Kd2=10;
Ki2=10;
Kpid=[Kp1 Kd1 Ki1 Kp2 Kd2 Ki2];
%% ODE solving
% opt1=odeset('RelTol',1e-10,'AbsTol',1e-20,'NormControl','off');
[T,X] = ode45(@(t,x) r2dof(t,x,ths,spec,Kpid),Ts,x0);
plot(T, X)
function xdot=r2dof(t,x,ths,spec,Kpid)
xdot=zeros(8,1);
%% set-points
th1s=ths(1);
th2s=ths(2);
%% Robot Specifications
M1=spec(3);
M2=spec(4);
L1=spec(1);
L2=spec(2);
g=9.8;
%% Inertia Matrix
b11=(M1+M2)*L1^2+M2*L2^2+2*M2*L1*L2*cos(x(4));
b12=M2*L2^2+M2*L1*L2*cos(x(4));
b21=M2*L2^2+M2*L1*L2*cos(x(4));
b22=M2*L2^2;
Bq=[b11 b12;b21 b22];
%% C Matrix
c1=-M2*L1*L2*sin(x(4))*(2*x(5)*x(6)+x(6)^2);
c2=-M2*L1*L2*sin(x(4))*x(5)*x(6);
Cq=[c1;c2];
%% Gravity Matrix
g1=-(M1+M2)*g*L1*sin(x(3))-M2*g*L2*sin(x(3)+x(4));
g2=-M2*g*L2*sin(x(3)+x(4));
Gq=[g1;g2];
%% PID Control
% PID parameters for theta 1
Kp1=Kpid(1);
Kd1=Kpid(2);
Ki1=Kpid(3);
% PID parameters for theta 2
Kp2=Kpid(4);
Kd2=Kpid(5);
Ki2=Kpid(6);
%decoupled control input
f1=Kp1*(th1s-x(3))-Kd1*x(5)+Ki1*(x(1));
f2=Kp2*(th2s-x(4))-Kd2*x(6)+Ki2*(x(2));
Fhat=[f1;f2];
F=Bq*Fhat; % actual input to the system
%% System states
xdot(1)=(th1s-x(3)); %dummy state of theta1 integration
xdot(2)=(th2s-x(4)); %dummy state of theta2 integration
xdot(3)=x(5); %theta1-dot
xdot(4)=x(6); %theta2-dot
q2dot=inv(Bq)*(-Cq-Gq+F);
xdot(5)=q2dot(1); %theta1-2dot
xdot(6)=q2dot(2); %theta1-2dot
%control input function output to outside computer program
xdot(7)=F(1);
xdot(8)=F(2);
end
So if you get errors, you do not run the posted code.
If you still have problem, post the code you run.

Sam Chak
Sam Chak 2022-12-11
The desired control torques cannot be integrated directly in and , because the actual torques are delivered by the motors. This, some minor modifications are made to assume and include the motor dynamics in and , as well as joint dynamics and .
However, I didn't check if the kinematics and are true for large angles. Swinging ° to ° are relatively large motion. Logically, there should be some constraints (or maybe singularity), like the elbow and shoulder.
%% Initilization
th_int = [-pi/2 pi/2]; % initial positions
ths = [ pi/2 -pi/2]; % setpoints
x0 = [0 0 th_int 0 0 0 0]; % states initial values
Ts = linspace(0, 10, 1001); % time span
%% Robot Specifications
L1 = 1; % link 1
L2 = 1; % link 2
M1 = 1; % mass 1
M2 = 1; % mass 2
spec = [L1 L2 M1 M2];
%% PID Parameters
% PID parameters for theta 1
Kp1 = 15;
Kd1 = 7;
Ki1 = 10;
% PID parameters for theta 2
Kp2 = 15;
Kd2 = 10;
Ki2 = 10;
Kpid = [Kp1 Kd1 Ki1 Kp2 Kd2 Ki2];
%% ODE solving
% opt1=odeset('RelTol',1e-10,'AbsTol',1e-20,'NormControl','off');
[t, x] = ode45(@(t,x) r2dof(t, x, ths, spec, Kpid), Ts, x0);
plot(t, x(:,3:4)*180/pi), grid on, ylabel('Angle, [\circ]')
legend('\theta_1', '\theta_2', 'fontsize', 16)
yticks(-150:60:150)
plot(t, x(:,7:8)), grid on, ylabel('Torque, [Nm]')
legend('\tau_1', '\tau_2', 'fontsize', 16)
function xdot = r2dof(t, x, ths, spec, Kpid)
xdot = zeros(8, 1);
%% set-points
th1s = ths(1);
th2s = ths(2);
%% Robot Specifications
M1 = spec(3);
M2 = spec(4);
L1 = spec(1);
L2 = spec(2);
g = 9.8;
%% Inertia Matrix
b11 = (M1 + M2)*L1^2 + M2*L2^2 + 2*M2*L1*L2*cos(x(4));
b12 = M2*L2^2 + M2*L1*L2*cos(x(4));
b21 = M2*L2^2 + M2*L1*L2*cos(x(4));
b22 = M2*L2^2;
Bq = [b11 b12; b21 b22];
%% C Matrix
c1 = - M2*L1*L2*sin(x(4))*(2*x(5)*x(6) + x(6)^2);
c2 = - M2*L1*L2*sin(x(4))*x(5)*x(6);
Cq = [c1; c2];
%% Gravity Matrix
g1 = - (M1 + M2)*g*L1*sin(x(3)) - M2*g*L2*sin(x(3) + x(4));
g2 = - M2*g*L2*sin(x(3) + x(4));
Gq = [g1; g2];
%% PID Control
% PID parameters for theta 1
Kp1 = Kpid(1);
Kd1 = Kpid(2);
Ki1 = Kpid(3);
% PID parameters for theta 2
Kp2 = Kpid(4);
Kd2 = Kpid(5);
Ki2 = Kpid(6);
% decoupled control input
f1 = - Kp1*(x(3) - th1s) - Kd1*x(5) - Ki1*(x(1));
f2 = - Kp2*(x(4) - th2s) - Kd2*x(6) - Ki2*(x(2));
Fhat = [f1; f2];
F = Bq*Fhat; % desired torques computed by formulas to the system, x7 & x8 are the actual torques delivered by the motors to the joints
%% System states
xdot(1) = x(3) - th1s; % dummy state of theta1 integration
xdot(2) = x(4) - th2s; % dummy state of theta2 integration
xdot(3) = x(5); % theta1-dot
xdot(4) = x(6); % theta2-dot
q2dot = inv(Bq)*(- Cq - Gq + [x(7); x(8)]); % actuated torques are delivered by the motors
xdot(5) = q2dot(1); % theta1-2dot
xdot(6) = q2dot(2); % theta1-2dot
% control input function output to outside computer program
% (assume the motor is fast enough to be approximated by a 1st-order function)
xdot(7) = - 1e2*(x(7) - F(1)); % can modify the motor constant
xdot(8) = - 1e2*(x(8) - F(2));
end

类别

Help CenterFile Exchange 中查找有关 PID Controller Tuning 的更多信息

Community Treasure Hunt

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

Start Hunting!

Translated by