应《网络安全法》要求，自2017年10月1日起，未进行实名认证将不得使用互联网跟帖服务。为保障您的帐号能够正常使用，请尽快对帐号进行手机号验证，感谢您的理解与支持！

24小时热门版块排行榜    

小木虫论坛-学术科研互动平台 » 计算模拟区 » 仿真模拟 » MATLAB » 拼了。。700金币求高手。最优控制，动态规划 matlab 仿真。。
52/1返回列表
查看: 3066  |  回复: 17
只看楼主 @他人 存档 新回复提醒 (忽略) 收藏    在APP中查看
当前只显示满足指定条件的回帖，点击这里查看本话题的所有回帖

新手菜鸟1818

金虫 (小有名气)

[求助] 拼了。。700金币求高手。最优控制，动态规划 matlab 仿真。。

初次接触matlab。。初次接触最优控制，动态规划。。实在是没有招了。。

炎炎酷暑。。阴面无风无空调无风扇。。蒸桑拿似的。。我和饺子只差一瓶酱油了。。

所有金币700个全扑了。。只要你能解决问题。。金币不是问题。。

我尝试做出 论文"Data-Driven Robust Approximate Optimal Tracking Control for Unknown General Nonlinear Systems Using Adaptive Dynamic Programming Method“中的仿真。。

自己沿着一个师弟的俩个以前的程序做了一下。。如下：

example11fun.m (貌似是被调用的)

function dx=example11fun(t,x)
x1d=sin(t);
%%%%参数设置
gamma5=0.1;eta=1.5;alphac=0.8;alphaa=0.5;Q=1;
R=1; kr=[20 20];zeta=1.2;
%%%%控制%%%%%%%%%%%%%%%%%%
ue=[2*x(26)*x(31)+x(27)*x(32);x(26)*x(32)+2*x(27)*x(33)];
ue1=norm(ue);
cu=[x(11) x(12);x(13) x(14)];
dxd1=x(3)*x(22)+x(5)*x(23)+x(7)*tanh(x(22))+x(9)*tanh(x(23))*ud1+x(15);
dxd2=x(4)*x(22)+x(6)*x(23)+x(8)*tanh(x(22))+x(10)*tanh(x(23))*ud1+x(16)
ud=(inv(cu))'*([dxd1;dxd2]-[x(3) x(5);x(4) x(6)].*[x(22);x(23)]-[x(7) x(9);x(8) x(10)].*[tanh(x(22));tanh(x(23))]-[x(15);x(16)]);
ud1=norm(ud);
ur=kr*[x(26);x(27)]./(x(26)^2+x(27)^2+zeta);
u=ue1+ud1;
uad=u-ur;

%%%%%%%%
sigma1=2*x(26)*(x(3)*x(26)+x(5)*x(27)+x(7)*tanh(x(24))-x(7)*tanh(x(22))-x(9)*tanh(x(25))+x(9)*tanh(x(23))+x(13)*ue1);
sigma2=x(27)*(x(3)*x(26)+x(5)*x(27)+x(7)*tanh(x(24))-x(7)*tanh(x(22))-x(9)*tanh(x(25))+x(9)*tanh(x(23))+x(13)*ue1)+x(26)*(x(4)*x(26)+x(6)*x(27)+x(8)*tanh(x(24))-x(8)*tanh(x(22))-x(10)*tanh(x(25))+x(10)*tanh(x(23))+x(14)*ue1);
sigma3=2*x(27)*(x(4)*x(26)+x(6)*x(27)+x(8)*tanh(x(24))-x(8)*tanh(x(22))-x(10)*tanh(x(25))+x(10)*tanh(x(23))+x(14)*ue1);
sigmac1=sigma1/(sigma1^2+sigma2^2+sigma3^2+1);
sigmac2=sigma2/(sigma1^2+sigma2^2+sigma3^2+1);
sigmac3=sigma3/(sigma1^2+sigma2^2+sigma3^2+1);

dx=[-x(1)+x(2);  %1 系统方程
    -0.5*x(1)-0.5*x(2)*(1-(cos(2*x(1))+2)^2)+(cos(2*x(1))+2)*u;  %2 系统方程
    (x(18)+0.1*x(19))*(x(1)-x(18)); %%3 自适应调节C11
    (x(18)+0.1*x(19))*(x(2)-x(19)); %%4 自适应调节C12
    (0.1*x(18)+x(19))*(x(1)-x(18)); %%5 自适应调节C21
    (0.1*x(18)+x(19))*(x(2)-x(19)); %%6 自适应调节C22
    (tanh(x(18))+0.2*tanh(x(19)))*(x(1)-x(18)); %%7  自适应调节A11
    (tanh(x(18))+0.2*tanh(x(19)))*(x(2)-x(19)); %%8  自适应调节A12
    (0.2*tanh(x(18))+tanh(x(19)))*(x(1)-x(18)); %%9  自适应调节A21
    (0.2*tanh(x(18))+tanh(x(19)))*(x(2)-x(19)); %%10 自适应调节A22
    (0.1*u)*(x(1)-x(18)); %%11  自适应调节Cu11
    (0.1*u)*(x(2)-x(19)); %%12  自适应调节Cu12
    (u)*(x(1)-x(18));     %%13  自适应调节Cu21
    (u)*(x(2)-x(19));     %%14  自适应调节Cu22
    0.2*(x(1)-x(18));   %%15 自适应调节Au1   
    0.2*(x(2)-x(19));   %%16 自适应调节Au2
    -gamma5*[x(1)-x(18) x(2)-x(19)].*[x(1)-x(18);x(2)-x(19)]/([x(1)-x(18) x(2)-x(19)]*[x(1)-x(18);x(2)-x(19)]+eta); %%17 自适应调节theta
    x(3)*x(18)+x(5)*x(19)+x(7)*tanh(x(18))+x(9)*tanh(x(19))+x(13)*u+x(15)-x(20);  %%18   x1的估计/5式
    x(4)*x(18)+x(6)*x(19)+x(8)*tanh(x(18))+x(10)*tanh(x(19))+x(14)*u+x(16)-x(21); %%19   x2的估计/5式
    -60*(x(1)-x(18))+x(17)*(x(1)-x(18))/([x(1)-x(18) x(2)-x(19)]*[x(1)-x(18);x(2)-x(19)]+eta);  %%20   v1/6式
    -60*(x(2)-x(19))+x(17)*(x(2)-x(19))/([x(1)-x(18) x(2)-x(19)]*[x(1)-x(18);x(2)-x(19)]+eta);  %%21   v2/6式
    x(3)*x(22)+x(5)*x(23)+x(7)*tanh(x(22))+x(9)*tanh(x(23))+x(13)*ud+x(15);    %%22   期望输出xd1/18式
    x(4)*x(22)+x(6)*x(23)+x(8)*tanh(x(22))+x(10)*tanh(x(23))+x(14)*ud+x(16);   %%23   期望输出xd2/18式
    x(3)*x(18)+x(5)*x(19)+x(7)*tanh(x(18))+x(9)*tanh(x(19))+x(13)*u+x(15);  %%24   x1的重写/17式
    x(4)*x(18)+x(6)*x(19)+x(8)*tanh(x(18))+x(10)*tanh(x(19))+x(14)*u+x(16); %%25   x2的重写/17式
    x(24)-x(22);   %%%26 系统误差e1
    x(25)-x(23);   %%%27 系统误差e2
    -alphac*sigmac1*([x(26)^2 x(26)*x(27) x(27)^2]*[x(28);x(29);x(30)]+[x(26) x(27)]*Q*[x(26);x(27)]+R*ue^2); %%28 调节wc1 31式
    -alphac*sigmac2*([x(26)^2 x(26)*x(27) x(27)^2]*[x(28);x(29);x(30)]+[x(26) x(27)]*Q*[x(26);x(27)]+R*ue^2); %%29 调节wc2 31式
    -alphac*sigmac3*([x(26)^2 x(26)*x(27) x(27)^2]*[x(28);x(29);x(30)]+[x(26) x(27)]*Q*[x(26);x(27)]+R*ue^2); %%30 调节wc3 31式
    -alphaa*2*x(26)*(2*x(26)*x(31)+x(27)*x(32)+x(26)*x(11)*x(28)+0.5*(x(27)*x(11)+x(26)*x(12))*x(29)+x(27)*x(12)*x(30));    %%31 调节wa1 38式
    -alphaa*(x(27)*(2*x(26)*x(31)+x(27)*x(32)+x(26)*x(11)*x(28)+0.5*(x(27)*x(11)+x(26)*x(12))*x(29)+x(27)*x(12)*x(30))+x(26)*(x(26)*x(32)+2*x(27)*x(33)+x(26)*x(13)*x(28)+0.5*(x(27)*x(13)+x(26)*x(14))*x(29)+x(27)*x(14)*x(30)));    %%32 调节wa2 38式
    -alphaa*2*x(27)*(x(26)*x(32)+2*x(27)*x(33)+x(26)*x(13)*x(28)+0.5*(x(27)*x(13)+x(26)*x(14))*x(29)+x(27)*x(14)*x(30));    %%33 调节wa3 38式
    ];

example11.m （文件）

clear;
clc;
t=[0:0.01:50];
x0=[0;0;0;0.5;0.5;0.5;0.5;0;0;0;0.5;0.5;0.5;0.5;0;0;0;0.5;0.5;0.5;0.5;0;0;0;0.5;0.5;0.5;0.5;0;0;0;0.5;0.5];i=1;
[t,x]=ode45('example11fun',t,x0);
x1d=sin(t);
plot(t,x(:,1),'-',t,x1d,'-.');grid

。
。
。

实在是不行了。。运行部了。。进行下去太困难了。。

也许我最初的思路就是错误的。。。

总之。。求助高手。。

你要是能找此文的作者张老师（控制领域绝对大牛），崔老师（猜测为张老师学生，可能有部分程序），张xin老师或者罗老师 要到仿真程序那是最好不过了。。。

700金币坐等你来拿。。
回复此楼

» 猜你喜欢
1楼 2013-08-10 18:43:05
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

月只蓝

主管区长 (职业作家)

【答案】应助回帖

★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ...
xiegangmai: 金币+2, 谢谢参与！ 2013-08-13 23:40:14
新手菜鸟1818: 金币+700, ★★★★★最佳答案, 就是一个字。。。服！！ 2013-08-14 15:05:03
引用回帖:
13楼: Originally posted by 新手菜鸟1818 at 2013-08-13 16:14:17
example11.m 文件如下：

clear;
clc;
t=;
x0=;i=1;
=ode45('example11fun',t,x0);
x1d=sin(t);
plot(t,x(:,1),'-',t,x1d,'-.');grid

example11fun.m  如下：
function dx=example11fun(t,x)
%%%%期望 ...

考察了一下积分步长和积分区间对计算耗时的影响，发现积分区间是主导因素，以t右端点为3.5为例，积分步长分别为0.001、0.01、0.1时，计算总耗时如下(计算机性能不同，耗时也会不同)：

积分步长       0.001     0.01       0.1      0.5
计算耗时/s         9.50      9.48      9.20     9.08
可见积分步长，对计算耗时影响很小。

但是右端点的数值影响很大(以步长0.01为例)：
右端点数值      0.5    1     1.5    2    2.5    3         3.5      
计算耗时/s         0.45   0.46    0.47    0.48    0.66    1.52      9.48   

可见，计算耗时，从右端点大于3以后，开始剧烈增大，为什么会发生在3左右(很接近于3.14，即π，系统的周期性？)，由于对你的体系具体的情况不了解，猜测可能与你研究的体系的特点有关。个人感觉，在程序不变的情况下，也许积分区间增大到5，就需要几天的计算时间，积分区间增大到10，这是无法想象的天文数字。

总之，如果积分区间右端点取到3左右就可以的话，这个程序算是成功的。
如果实际要求一定要10的话，建议你好好检查一下程序，是不是有什么地方写错了，参数不对等情况。你可以以0.01的积分步长，区间取到4.5，让MATLAB一直算，把电脑通宵开着，看看需要多久。

另外，我对你的程序整理一下，没有改动任何参数的数值，只是把原本2个m文件，合并为一个，并添加了计算消耗时间。代码参考如下：
CODE:
function solve_example11fun
clear all;clc;
tic,
t=0:0.01:3.5;
x0=[0;0;0;0.5;0.5;0.5;0.5;0;0;0;0.5;0.5;0.5;0.5;0;0;0;0.5;0.5;0.5;0.5;0;0;0;0.5;0.5;0.5;0.5;0;0;0;0.5;0.5];
[t,x]=ode45(@example11fun,t,x0);
toc;
x1d=sin(t);
figure,plot(t,x(:,1),'b-',t,x1d,'r-');grid on,xlabel('t'),legend('x(1)与t关系曲线','x1d与t关系曲线','Location','best');


%--------------------------------------------------------------------------%
function dx=example11fun(t,x)
%-------------------------参数设置----------------------------------------%
gamma5=0.1;
eta=1.5;
alphac=0.8;
alphaa=0.5;
Q=1;
R=1;
kr=[20 20];
zeta=1.2;
%--------------------------控制设计---------------------------------------%

ue=[2*x(26)*x(31)+x(27)*x(32);x(26)*x(32)+2*x(27)*x(33)];
ue1=norm(ue);
cu=[x(11) x(12);x(13) x(14)];
ud=cu'*([x(22);x(23)]-[x(3) x(5);x(4) x(6)]*[x(22);x(23)]-[x(7) x(9);x(8) x(10)]*[tanh(x(22));tanh(x(23))]-[x(15);x(16)]);
ud1=norm(ud);
ur=kr*[x(26)/(x(26)^2+x(27)^2+zeta);x(27)/(x(26)^2+x(27)^2+zeta)];
u=ue1+ud1;
uad=u-ur;

%------------------------------------------------------------------------&

sigma1=2*x(26)*(x(3)*x(26)+x(5)*x(27)+x(7)*tanh(x(24))-x(7)*tanh(x(22))-x(9)*tanh(x(25))+x(9)*tanh(x(23))+x(13)*ue1);
sigma2=x(27)*(x(3)*x(26)+x(5)*x(27)+x(7)*tanh(x(24))-x(7)*tanh(x(22))-x(9)*tanh(x(25))+x(9)*tanh(x(23))+x(13)*ue1)+x(26)*(x(4)*x(26)+x(6)*x(27)+x(8)*tanh(x(24))-x(8)*tanh(x(22))-x(10)*tanh(x(25))+x(10)*tanh(x(23))+x(14)*ue1);
sigma3=2*x(27)*(x(4)*x(26)+x(6)*x(27)+x(8)*tanh(x(24))-x(8)*tanh(x(22))-x(10)*tanh(x(25))+x(10)*tanh(x(23))+x(14)*ue1);
sigmac1=sigma1/(sigma1^2+sigma2^2+sigma3^2+1);
sigmac2=sigma2/(sigma1^2+sigma2^2+sigma3^2+1);
sigmac3=sigma3/(sigma1^2+sigma2^2+sigma3^2+1);

%---------------------------------定义常微分方程组------------------------%
dx=[-x(1)+x(2);  %1 系统方程
     -0.5*x(1)-0.5*x(2)*(1-(cos(2*x(1))+2)^2)+(cos(2*x(1))+2)*u;  %2 系统方程
     (x(18)+0.1*x(19))*(x(1)-x(18)); %%3 自适应调节C11
     (x(18)+0.1*x(19))*(x(2)-x(19)); %%4 自适应调节C12
     (0.1*x(18)+x(19))*(x(1)-x(18)); %%5 自适应调节C21
     (0.1*x(18)+x(19))*(x(2)-x(19)); %%6 自适应调节C22
     (tanh(x(18))+0.2*tanh(x(19)))*(x(1)-x(18)); %%7  自适应调节A11
     (tanh(x(18))+0.2*tanh(x(19)))*(x(2)-x(19)); %%8  自适应调节A12
     (0.2*tanh(x(18))+tanh(x(19)))*(x(1)-x(18)); %%9  自适应调节A21
     (0.2*tanh(x(18))+tanh(x(19)))*(x(2)-x(19)); %%10 自适应调节A22
     (0.1*u)*(x(1)-x(18)); %%11  自适应调节Cu11
     (0.1*u)*(x(2)-x(19)); %%12  自适应调节Cu12
     u*(x(1)-x(18));     %%13  自适应调节Cu21
     u*(x(2)-x(19));     %%14  自适应调节Cu22
     0.2*(x(1)-x(18));   %%15  自适应调节Au1   
    0.2*(x(2)-x(19));   %%16  自适应调节Au2
     -gamma5*((x(1)-x(18))^2+(x(2)-x(19))^2)/((x(1)-x(18))^2+(x(2)-x(19))^2+eta); %%17 自适应调节theta
     x(3)*x(18)+x(5)*x(19)+x(7)*tanh(x(18))+x(9)*tanh(x(19))+x(13)*u+x(15)-x(20);  %%18   x1的估计/5式
     x(4)*x(18)+x(6)*x(19)+x(8)*tanh(x(18))+x(10)*tanh(x(19))+x(14)*u+x(16)-x(21); %%19   x2的估计/5式
     -60*(x(1)-x(18))+x(17)*(x(1)-x(18))/((x(1)-x(18))^2+(x(2)-x(19))^2+eta);  %%20   v1/6式
     -60*(x(2)-x(19))+x(17)*(x(2)-x(19))/((x(1)-x(18))^2+(x(2)-x(19))^2+eta);  %%21   v2/6式
     x(3)*x(22)+x(5)*x(23)+x(7)*tanh(x(22))+x(9)*tanh(x(23))+x(13)*ud1+x(15);    %%22   期望输出xd1/18式
     x(4)*x(22)+x(6)*x(23)+x(8)*tanh(x(22))+x(10)*tanh(x(23))+x(14)*ud1+x(16);   %%23   期望输出xd2/18式
     x(3)*x(18)+x(5)*x(19)+x(7)*tanh(x(18))+x(9)*tanh(x(19))+x(13)*u+x(15);  %%24   x1的重写/17式
     x(4)*x(18)+x(6)*x(19)+x(8)*tanh(x(18))+x(10)*tanh(x(19))+x(14)*u+x(16); %%25   x2的重写/17式
     x(24)-x(22);   %%%26 系统误差e1
     x(25)-x(23);   %%%27 系统误差e2
     -alphac*sigmac1*((x(26)^2*x(28)+x(26)*x(27)*x(29)+x(27)^2*x(30))+Q*x(26)^2+Q*x(27)^2+R*ue1^2); %%28 调节wc1 31式
     -alphac*sigmac2*((x(26)^2*x(28)+x(26)*x(27)*x(29)+x(27)^2*x(30))+Q*x(26)^2+Q*x(27)^2+R*ue1^2); %%29 调节wc2 31式
     -alphac*sigmac3*((x(26)^2*x(28)+x(26)*x(27)*x(29)+x(27)^2*x(30))+Q*x(26)^2+Q*x(27)^2+R*ue1^2); %%30 调节wc3 31式
     -alphaa*2*x(26)*(2*x(26)*x(31)+x(27)*x(32)+x(26)*x(11)*x(28)+0.5*(x(27)*x(11)+x(26)*x(12))*x(29)+x(27)*x(12)*x(30));    %%31 调节wa1 38式
     -alphaa*(x(27)*(2*x(26)*x(31)+x(27)*x(32)+x(26)*x(11)*x(28)+0.5*(x(27)*x(11)+x(26)*x(12))*x(29)+x(27)*x(12)*x(30))+x(26)*(x(26)*x(32)+2*x(27)*x(33)+x(26)*x(13)*x(28)+0.5*(x(27)*x(13)+x(26)*x(14))*x(29)+x(27)*x(14)*x(30)));    %%32 调节wa2 38式
     -alphaa*2*x(27)*(x(26)*x(32)+2*x(27)*x(33)+x(26)*x(13)*x(28)+0.5*(x(27)*x(13)+x(26)*x(14))*x(29)+x(27)*x(14)*x(30));
     ];
赞一下 回复此楼
MATLAB、MS小问题、普通问题请发帖求助！时间精力有限，恕不接受无偿私信求助。
17楼2013-08-13 20:35:12
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
查看全部 18 个回答

月只蓝

主管区长 (职业作家)

【答案】应助回帖

★ ★ ★
感谢参与，应助指数 +1
xiegangmai: 金币+3, 鼓励讨论交流！ 2013-08-12 10:53:07
整个程序看似挺长挺复杂，其实本质上是求解一个普通含有33个方程的常微分方程组的初值问题，求解的方法是用四阶精度的龙格库塔法。程序运行不出来，主要原因在于，构建模型时，需要提供的参数、变量方程没有正确给出。在附图1中，编号为1,2,3,4的行，先来看1行和2行，1行和2行为了定义dxd1和dxd2，需要用到ud1，我们自然要去找ud1在哪里被定义了，发现是在4行，4行为了定义ud1，又需要用到ud，我们又去找ud在哪里，发现是在3行，这时候问题就出来了，3行为了定义ud，用到了1行和2行的dx1和dx2，也就是说，1行和2行的dx1和dx2是隐含的，不是显式。这样，在程序运行的时候，必然会首先提示，1行的ud1没有被定义。
建议你仔细推敲一下，你的带求解模型各已知参数、变量的关系，模型构思正确了，程序的正确性才有保证。
拼了。。700金币求高手。最优控制，动态规划 matlab 仿真。。
附图1.jpg
赞一下 回复此楼
MATLAB、MS小问题、普通问题请发帖求助！时间精力有限，恕不接受无偿私信求助。
5楼2013-08-12 10:14:23
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

月只蓝

主管区长 (职业作家)
上一楼的附图太小，重传了张大点的图。
拼了。。700金币求高手。最优控制，动态规划 matlab 仿真。。-1
附图1.jpg
赞一下(1人) 回复此楼

» 本帖已获得的红花（最新10朵）

新手菜鸟1818
MATLAB、MS小问题、普通问题请发帖求助！时间精力有限，恕不接受无偿私信求助。
6楼2013-08-12 10:16:09
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

新手菜鸟1818

金虫 (小有名气)
送红花一朵
引用回帖:
6楼: Originally posted by 月只蓝 at 2013-08-12 10:16:09
上一楼的附图太小，重传了张大点的图。

附图1.jpg

非常感谢您的指点。送花一朵。。

按照您的建议，我把那几行的程序重新调整了一下。。但是还是运行部了。。显示矩阵维数不对。。但是我检查了好几遍，感觉程序里面的矩阵的维数是对的。。请您指点。程序如下：
cu=[x(11) x(12);x(13) x(14)];
ud=cu'*([x(22);x(23)]-[x(3) x(5);x(4) x(6)]*[x(22);x(23)]-[x(7) x(9);x(8) x(10)]*[tanh(x(22));tanh(x(23))]-[x(15);x(16)]);
ud1=norm(ud);
ur=kr*[x(26);x(27)]./(x(26)^2+x(27)^2+zeta);
ur1=norm(ur);
u=ue1+ud1;
uad=u-ur1;

非常感谢您！
赞一下 回复此楼
7楼2013-08-12 16:25:30
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
查看全部 18 个回答
普通表情 龙 兔 虎 猫
信息提示
关闭
请填处理意见
关闭