所以我想为我的跷跷板系统创建一个MPC控制器。所有“繁重的工作”(获取运动方程、状态空间表示等)都已完成,所以我开始在 MATLAB 中进行编码。以下是:
yalmip('clear')
clc; clear all; close all;
%% Parameters
a = 0.116553; % height of center mass of the pendulum, [m]
b = 0; % position of the weight B on the vertical rod (not taken into consideration)
c = 0.180047; % height of the center of mass of the cart, [m]
mA = 4.839; % mass of the pendulum, [kg]
mB = 0; % not taken into consideration
mC = 1; % mass of the cart, [kg]
g = 9.81; % gravity factor, [m/s^2]
kappa = 0.1; % coefficient of the viscous damping in the rotational joint
J = 0.68; % moment of inertia of the pendulum, [kgm^2]
Ke = 0.077; % motor constant of the EM force, [Vs^-1]
Kt = 0.077; % proportional moment motor constant, [NmA^-1]
Ra = 2.6; % electrical resistance, [ohm]
p = 1/3.71; % motor gearbox ratio
r = 7.7*10^(-3); % effective radius of the motor shaft, [m]
A = [0 1 0 0;
(a*mA*g+b*mB*g)/(J+mB*b^2) -kappa/(J+mB*b^2)...
-mC*g/(J+mB*b^2) -Ke*Kt/(mC*(J+mB*b^2)*Ra*p^2*r^2);
0 0 0 1;
(a*mA*g+b*mB*g)/(J+mB*b^2)-g -kappa*c/(J+mB*b^2)...
-c*mC*g/(J+mB*b^2) -(J+mB*b^2+mC*c^2)*Ke*Kt/(mC*(J+mB*b^2)*Ra*p^2*r^2)];
B = [0; Kt/(mC*(J+mB*b^2)*Ra*p^2*r^2); 0; (J+mB*b^2+mC*c^2)*Kt/(mC*(J+mB*b^2)*Ra*p^2*r^2)];
C = [1 0 0 0;
0 0 1 0];
D = 0;
Ts = 15; % sample time
sys = ss(A, B, C, D);
d_sys = c2d(sys,Ts); % discrete time
x0 = [0.1; 0; 0; 0]; % initial values
%% MPC
nx = 2; % number of states
nu = 1; % number of inputs
Q = eye(4);
R = 4;
N = 2;
u = sdpvar(repmat(nu,4,N),repmat(1,4,N));
constraints = [];
objective = 0;
ops = sdpsettings('verbose',2,'solver','quadprog');
xList=[];
uList=[];
SimTime=15;
xList=[xList,x0]
for i=1:SimTime
x = x0;
figure(1)
plot([1,SimTime+N],[-1,-1],'k');
hold on
plot([1,SimTime+N],[1,1],'k');
for k = 1:N
x = A*x + B*u{k};
objective = objective + x'*Q*x + u{k}*R*u{k};
constraints = [constraints , -3<=u{k}<=3, -5<=x<=5];
end
sol=optimize(constraints, objective, ops);
uValues=[];
for j=1:N
uValues=[uValues,value(u{j})];
end
x0=A*x0+B*value(u{1});
xList=[xList,x0]
uList=[uList,value(u{1})]
objective=0;
constraints =[];
pause(0.005)
figure(1)
hold on
stairs(i:i+N-1,uValues,'ro-')
plot(1:i,uList,'c-')
end
figure(2)
plot(xList')
基本上,当我运行它时,模拟开始,当我得到一个“NaN”时它很快就失败了。只有当我把我的初始值都设为 0 时,我才能真正得到一些东西。请记住,我需要在我的系统中设置一些系统中普遍存在的约束(最大电机扭矩,推车可以在跷跷板上移动多远......)。这就是我在连续时间内使用状态空间表示时的全部情况,当我尝试通过“c2d”命令使用离散时间并据此更改其余代码时,模拟甚至无法正常启动。在这种特定情况下,我真的不知道该怎么做。是否有其他方法可以为该系统创建 MPC 调节器?注意:我将 YALMIP 用于 MPC 的优化部分。