从伺服电机的模型预测控制器生成奖励函数

这个例子使用了:

打开实时脚本

这个例子展示了如何从模型预测控制器对象中定义的成本和约束规范自动生成奖励函数。然后使用生成的奖励函数来训练一个强化学习代理。

简介

您可以使用generateRewardFunction从模型预测控制器中指定的成本和约束开始，生成用于强化学习的奖励函数。由此产生的奖励信号是成本(由目标函数定义)和约束违反惩罚的总和，取决于当前的环境状态。

本例基于具有不可测输出约束的直流伺服电机(模型预测控制工具箱)例如，在电压和轴转矩约束下，为直流伺服机构设计模型预测控制器。方法中定义的成本和约束规范将在此进行转换货币政策委员会对象转换为奖励函数，并使用它来训练代理控制伺服电机。

打开本例的Simu金宝applink模型，该模型基于上述MPC示例，但已针对强化学习进行了修改。

open_system (“rl_motor”）

创建模型预测控制器

创建电机的开环动态模型，定义在植物和最大允许扭矩τ使用辅助函数。

[plant,tau] = mpcmotormodel;

指定MPC控制器的输入和输出信号类型。轴的角位置，被测量为第一次输出。第二个输出，扭矩，是不可测量的。

植物= setmpcsignals(植物，“MV”, 1“莫”, 1“UO”2);

在被操作的变量上指定约束，并定义缩放因子。

MV = struct(“最小值”, -220,“马克斯”, 220,“ScaleFactor”, 440);

在前三个预测步骤中施加扭矩约束，并指定轴位置和扭矩的比例因子。

OV = struct(“最小值”,{负无穷变化,τ,τ;τ;负)},.．.“马克斯”,{正无穷,τ;τ,τ;正]},.．.“ScaleFactor”,{2π* 2 *τ});

为二次代价函数指定权值，实现角位置跟踪。将扭矩的权重设置为零，从而允许它在其约束范围内浮动。

权重= struct(“MV”，0,“MVRate”, 0.1,机汇的, 0.1 [0]);

创建MPC控制器植物模型的采样时间为0．1S，预测视界10步骤，和一个控制地平线2步骤，使用先前为权重、操作变量和输出变量定义的结构。

mpcobj = mpc(plant,0.1,10,2,Weights,MV,OV);

显示控制器规格。

mpcobj

MPC对象(创建于31- 8月2022 05:46:32):---------------------------------------------采样时间:0.1(秒)预测水平:10控制水平:2植物模型:-------------- 1操纵变量(s)—>| 4状态| | |—> 1测量输出(s) 0测量扰动(s)—>| 1输入| | |—> 1未测量输出(s) 0未测量扰动(s)—>| 2输出| --------------指数:(输入矢量)操纵变量:[1](输出矢量)测量输出:[1]未测量输出:[2]扰动和噪声模型:输出扰动模型:默认(类型“getoutdist(mpcobj)”为详细信息)测量噪声模型:默认(单位增益缩放后)权重:操纵变量:0操纵变量率:0.1000输出变量:[0.1000 0]ECR: 10000状态估计:默认卡尔曼滤波器(类型“getEstimator(mpcobj)”为详细信息)约束条件:-220 <= MV1 (V) <= 220, MV1/rate (V)无约束，MO1 (rad)无约束-78.54 <= UO1 (Nm)(t+1) <= 78.54 -78.54 <= UO1 (Nm)(t+2) <= 78.54 -78.54 <= UO1 (Nm)(t+3) <= 78.54 UO1 (Nm)(t+4)无约束

该控制器工作在一个有4个状态的装置上，1个输入(电压)和2个输出信号(角度和扭矩)，并具有以下规格:

被操纵变量、被操纵变量率和输出变量的代价函数权重分别为0、0.1和[0.1 0]。
被操纵变量被限制在-220V和220V之间。
被操纵的可变利率是不受约束的。
在前三个预测时间步中，第一个输出变量(角度)不受约束，但第二个输出变量(扭矩)被约束在-78.54 Nm和78.54 Nm之间，在第四步中不受约束。

注意，对于强化学习，只使用第一个预测时间步骤的约束规范，因为奖励是针对单个时间步骤计算的。

生成奖励函数

的规范中生成奖励函数代码货币政策委员会对象使用generateRewardFunction．代码显示在MATLAB编辑器中。

generateRewardFunction (mpcobj)

生成的奖励函数是奖励设计的起点。您可以使用不同的惩罚函数选择来修改函数并调整权重。对于本例，对生成的代码进行以下更改:

计算原始成本权重Qy，Qmv而且Qmvrate是100倍。
默认的外部惩罚函数方法为一步．将方法更改为二次．

变更后的成本和罚款规格如下:

Qy = [10 0];Qmv = 0;Qmvrate = 10;Py = Wy * exterorpenalty (y,ymin,ymax，“二次”）;Pmv = Wmv * exteriorPenalty(mv,mvmin,mvmax，“二次”）;Pmvrate = Wmvrate * exteriorPenalty(mv-lastmv,mvratemin,mvratemax, mvratemax)“二次”）;

对于本例，修改后的代码已保存在MATLAB函数文件中rewardFunctionMpc.m．显示生成的奖励函数。

类型rewardFunctionMpc.m

function reward = rewardFunctionMpc(y,refy,mv,refmv,lastmv) % rewardFunctionMpc从MPC规范中生成奖励。% %输入参数描述:% % y:第k步植物的输出变量+ 1% refy:第k步的参考输出变量+ 1% mv:第k步的操作变量% refmv:第k步的参考操作变量% lastmv:第k步的操作变量% % %限制(MPC和NLMPC): % -基于预测区间第一步计算的奖励。因此，信号预览和控制水平设置被忽略。% -不支持在线成本和约束更新。金宝app% -不考虑自定义成本和约束规格。% -不支持随时间变化的成本权重和约束。金宝app% -不考虑混合约束规范(对于MPC情况)。%强化学习工具箱% 02- 6月-2021 16:05:41 %#codegen %%来自MPC对象的规范%在'States'， 'OutputVariables'和% 'ManipulatedVariables'属性中指定的标准线性边界ymin = [-Inf -78.5398163397448];ymax = [Inf 78.5398163397448];Mvmin = -220; mvmax = 220; mvratemin = -Inf; mvratemax = Inf; % Scale factors as specified in 'States', 'OutputVariables', and % 'ManipulatedVariables' properties Sy = [6.28318530717959 157.07963267949]; Smv = 440; % Standard cost weights as specified in 'Weights' property Qy = [0.1 0]; Qmv = 0; Qmvrate = 0.1; %% Compute cost dy = (refy(:)-y(:)) ./ Sy'; dmv = (refmv(:)-mv(:)) ./ Smv'; dmvrate = (mv(:)-lastmv(:)) ./ Smv'; Jy = dy' * diag(Qy.^2) * dy; Jmv = dmv' * diag(Qmv.^2) * dmv; Jmvrate = dmvrate' * diag(Qmvrate.^2) * dmvrate; Cost = Jy + Jmv + Jmvrate; %% Penalty function weight (specify nonnegative) Wy = [1 1]; Wmv = 10; Wmvrate = 10; %% Compute penalty % Penalty is computed for violation of linear bound constraints. % % To compute exterior bound penalty, use the exteriorPenalty function and % specify the penalty method as 'step' or 'quadratic'. % % Alternatively, use the hyperbolicPenalty or barrierPenalty function for % computing hyperbolic and barrier penalties. % % For more information, see help for these functions. % % Set Pmv value to 0 if the RL agent action specification has % appropriate 'LowerLimit' and 'UpperLimit' values. Py = Wy * exteriorPenalty(y,ymin,ymax,'step'); Pmv = Wmv * exteriorPenalty(mv,mvmin,mvmax,'step'); Pmvrate = Wmvrate * exteriorPenalty(mv-lastmv,mvratemin,mvratemax,'step'); Penalty = Py + Pmv + Pmvrate; %% Compute reward reward = -(Cost + Penalty); end

要集成此奖励函数，请在Simulink模型中打开MATLAB函数块。金宝app

open_system (“rl_motor /奖励函数”）

用下面的代码行附加函数并保存模型。

r = rewardFunctionMpc(y,refy,mv,refmv,lastmv);

MATLAB函数块现在将执行rewardFunctionMpc.m在模拟。

对于这个例子，MATLAB函数块已经被修改和保存。

创建一个强化学习环境

环境动力学模型在伺服机构子系统。在这种环境下，

观察结果是最后8个时间步长的参考和实际输出变量(角度和扭矩)。
作用是电压 $V$ 适用于伺服电机。
采样时间为 $T_{年代} ＝ 0 ． 1 年代$ ．
仿真总时间为 $T_{f} ＝ 20. 年代$ ．

指定仿真总时间和采样时间。

Tf = 20;Ts = 0.1;

为环境创建观察和行动规范。

numObs = 32;numAct = 1;oinfo = rlNumericSpec([numObs 1]);ainfo = rlNumericSpec([numAct 1]，“LowerLimit”, -220,“UpperLimit”, 220);

创建强化学习环境rl金宝appSimulinkEnv函数。

环境= rlSi金宝appmulinkEnv(“rl_motor”，“rl_motor / RL代理”、oinfo ainfo);

创建一个强化学习代理

修正随机种子的重现性。

rng (0)

本例中的代理是双延迟深度确定性策略梯度(TD3)代理。

创建两个评论家表示。

%的批评家cnet = [featureInputLayer(numObs，“归一化”，“没有”，“名字”，“状态”) fullyConnectedLayer (128,“名字”，“fc1”) concatenationLayer(1、2、“名字”，“concat”) reluLayer (“名字”，“relu1”) fullyConnectedLayer (64,“名字”，“一个fc3”文件) reluLayer (“名字”，“relu2”) fullyConnectedLayer (1,“名字”，“CriticOutput”));actionPath = [featureInputLayer(numAct，“归一化”，“没有”，“名字”，“行动”) fullyConnectedLayer (8,“名字”，“取得”));criticNetwork = layerGraph(cnet);criticNetwork = addLayers(criticNetwork, actionPath);临界网络= connectLayers(临界网络，“取得”，“concat / in2”）;criticOptions = rlRepresentationOptions(“LearnRate”1 e - 3,“GradientThreshold”1);critic1 = rlQValueRepresentation(criticNetwork,oinfo,ainfo，.．.“观察”, {“状态”}，“行动”, {“行动”}, criticOptions);critic2 = rlQValueRepresentation(criticNetwork,oinfo,ainfo，.．.“观察”, {“状态”}，“行动”, {“行动”}, criticOptions);

创建一个参与者表示。

actorNetwork = [featureInputLayer(numObs，“归一化”，“没有”，“名字”，“状态”) fullyConnectedLayer (128,“名字”，“actorFC1”) reluLayer (“名字”，“relu1”) fullyConnectedLayer (64,“名字”，“actorFC2”) reluLayer (“名字”，“relu2”) fullyConnectedLayer (numAct“名字”，“行动”));actorOptions = rlRepresentationOptions(“LearnRate”1 e - 3,“GradientThreshold”1);actor = rlDeterministicActorRepresentation(actorNetwork,oinfo,ainfo，.．.“观察”, {“状态”}，“行动”, {“行动”}, actorOptions);

使用指定代理选项rlTD3AgentOptions．代理从最大容量的经验缓冲区进行训练1 e6通过随机选择小批量的大小256．的折现因子0.995喜欢长期奖励。

agentOpts = rlTD3AgentOptions(“SampleTime”Ts,.．.“DiscountFactor”, 0.995,.．.“ExperienceBufferLength”1 e6,.．.“MiniBatchSize”, 256);

该TD3代理中的探索模型为高斯模型。噪声模型在训练过程中为动作增加一个均匀的随机值。设置噪声的标准偏差为One hundred.．标准差的衰减速率为1 e-5的最小值0.005．

agentoptions . explorationmodel . standarddeviationmin = 0.005;agentoptions . explorationmodel . standard deviation = 100;agentoptions . explorationmodel . standarddeviationdecayrate = 1e-5;

使用参与者和评论家表示创建TD3代理。有关TD3制剂的更多信息，请参见rlTD3Agent．

agent = rlTD3Agent(actor，[critic1,critic2]，agentOpts);

培训代理人

要训练代理，首先指定使用的训练选项rlTrainingOptions．对于本例，使用以下选项:

每次训练最多运行2000集，每集持续时间最长装天花板(Tf / Ts)时间的步骤。
当智能体收到的平均累积奖励大于－2在20.连续集。此时，代理可以跟踪参考信号。

trainOpts = rlTrainingOptions(.．.“MaxEpisodes”, 2000,.．.“MaxStepsPerEpisode”装天花板(Tf / Ts),.．.“StopTrainingCriteria”，“AverageReward”，.．.“StopTrainingValue”2,.．.“ScoreAveragingWindowLength”, 20);

培训代理使用火车函数。训练这个代理是一个计算密集型的过程，可能需要几分钟才能完成。为了在运行此示例时节省时间，请通过设置加载预训练的代理doTraining来假．要亲自训练特工，请设置doTraining来真正的．

doTraining = false;如果doTraining trainingStats = train(agent,env,trainOpts);其他的负载(“rlDCServomotorTD3Agent.mat”）结束

训练进度快照如下图所示。由于训练过程中固有的随机性，您可以预期不同的结果。

验证控制器响应

为了验证训练代理的性能，模拟模型并在Scope块中查看响应。该强化学习agent能够在满足扭矩和电压约束的情况下跟踪参考角度。

sim卡(“rl_motor”）;

关闭模型。

close_system (“rl_motor”）