两物体之间的干摩擦:使用多重实验数据的参数估计
这个例子展示了如何利用多个实验数据估计非线性灰盒模型的参数。一个显示两个固体之间干摩擦的系统将被用作讨论的基础。在该系统中，一个物体是固定的，而另一个物体由于外生力在固定物体上向前和向后移动，如图1所示。
图1:两体系统的示意图。
两个物体之间的干摩擦建模
利用牛顿第三运动定律，运动物体的运动描述为:
F_tot(t) = m*a(t) = m*d/dt v(t) = m*d²/dt²s(t)
其中F_tot(t)等于外生力F(t)减去两个物体之间接触产生的摩擦力。摩擦力假定为滑动摩擦力F_s(t)和干摩擦力F_d(t)之和。前者通常被建模为速度的线性函数，即F_s(t) = k*v(t)，其中k是一个未知的滑动摩擦参数。另一方面，干摩擦是一种相当复杂的现象。在论文中:
克拉维尔、索林、张q。钢板弹簧系统的建模与识别。2001年第三届IFAC汽车控制先进研讨会。
它由一个常微分方程建模:
d/dt F_d(t) = -1/e*|v(t)|*F_d(t) + f/e*v(t)
其中e和f是两个未知参数，分别为距离和力的维度。
表示输入信号u(t) = F(t) [N]，引入状态为:
x1(t) = s(t)运动体位置[m]。x2(t) = v(t)运动物体的速度[m/s]。x3(t) = F_d(t)体间干摩擦力[N]。
模型参数为:
m运动体质量[m]。k滑动摩擦力系数[kg/s]。e距离相关干摩擦参数[m]。f力相关干摩擦参数[N]。
得到以下状态空间模型结构:
d / dt x1 (t) = x2 (t) d / dt x2 (t) = 1 / m * (u (t) - k * x2 (t) - x3 (t)) d / dt x3 (t) = 1 / e * (- | x2 (t) | * x3 (t) + f * x2 (t))
Y (t) = x1(t)
这些方程被输入到一个C-MEX模型文件twobodies_c.c中。它的状态和输出更新方程compute_dx和compute_y如下:
/*状态方程。*/ void compute_dx(double *dx, double t, double *x, double *u, double *p, const mxArray *auxvar){/*检索模型参数。*/ double *m， *k， *e， *f;M = p[0];/*运动物体的质量。*/ k = p[1];/*滑动摩擦力系数。*/ e = p[2];/*距离相关干摩擦参数。*/ f = p[3]; /* Force-related dry friction parameter. */
/* x[0]:位置。*/ /* x[1]:速率。*/ /* x[2]:干摩擦力。*/ dx[0] = x[1];Dx [1] = (u[0]-k[0]*x[1]-x[2])/m[0];Dx [2] = (-fabs(x[1])*x[2]+f[0]*x[1])/e[0];｝
/*输出方程。*/ void compute_y(双*y，双t，双*x，双*u，双*p, const mxArray *auxvar) {/* y[0]:位置。*/ y[0] = x[0];｝
写完描述模型结构的文件后，下一步是创建一个反映建模情况的IDNLGREY对象。我们还添加了关于模型结构的输入、输出、状态和模型参数的名称和单元的信息。注意，Parameters和InitialStates在这里被指定为向量，默认情况下，这意味着在调用NLGREYEST时，将估计所有模型参数，而没有初始状态向量。
文件名=“twobodies_c”；描述模型结构的文件。Order = [1 1 3];%模型订单[ny nu nx]。参数= [380;2200;0.00012;1900);初始参数向量。InitialStates = [0;0;0);初始状态。Ts = 0;%时间连续系统。nlgr = idnlgrey(文件名，顺序，参数，InitialStates, Ts，.．.“名字”，“两体系统”，.．.“InputName”，“外生力量”，.．.“InputUnit”，“N”，.．.“OutputName”，“移动物体的位置”，.．.“OutputUnit”，“米”，.．.“TimeUnit”，“年代”）;NLGR = setinit(NLGR，“名字”, {“移动物体的位置”.．.“运动物体的速度”.．.“物体之间的干摩擦力”});NLGR = setinit(NLGR，“单位”, {“米”“米/秒”“N”});NLGR = setpar(NLGR，“名字”, {运动物体的质量.．.滑动摩擦力系数，.．.“与距离相关的干摩擦参数”.．.与力相关的干摩擦参数});NLGR = setpar(NLGR，“单位”, {“米”公斤/年代”“米”“N”});
输入输出数据
此时，我们加载可用的(模拟的)输入-输出数据。该文件包含来自三个不同(模拟)测试运行的数据，每个测试运行包含1000个噪声损坏的输入输出样本，这些样本使用0.005秒的采样率(Ts)生成。输入u(t)是作用在运动物体上的外生力[N]。在实验中，输入是一个对称的锯齿形信号，其中所有实验的波形重复频率是相同的，但最大信号振幅在测试运行之间有所变化。输出y(t)是移动物体的位置[m](相对于固定物体)。出于建模的目的，我们创建了一个IDDATA对象，其中包含三个不同的实验:
负载(fullfile (matlabroot“工具箱”，“识别”，“iddemos”，“数据”，“twobodiesdata”));Z = merge(iddata(y1, u1, 0.005)， iddata(y2, u2, 0.005)， iddata(y3, u3, 0.005));z.Name =“两体系统”；z.InputName = nlgr.InputName;z.InputUnit = nlgr.InputUnit;z.OutputName = nlgr.OutputName;z.OutputUnit = nlgr.OutputUnit;z.Tstart = 0;z.TimeUnit = nlgr.TimeUnit;
用于后续识别实验的输入-输出数据显示在一个图窗口中。
图(“名字”[z。的名字':输入输出数据']);为I = 1:z。Ne zi = getexp(z, i);次要情节(z。Ne, 2,2 *i-1);%的输入。情节(子。SamplingInstants zi.InputData);标题([z。ExperimentName{我}“:”zi.InputName {1}));如果(i < z.Ne) xlabel(”）;其他的包含([z。域“(”子。TimeUnit“)”]);结束轴(“紧”）;次要情节(z。Ne, 2,2 *i);%输出。情节(子。SamplingInstants zi.OutputData);标题([z。ExperimentName{我}“:”zi.OutputName {1}));如果(i < z.Ne) xlabel(”）;其他的包含([z。域“(”子。TimeUnit“)”]);结束轴(“紧”）;结束
图2:两体系统的输入输出数据。
初始二体模型的性能
在估计四个未知参数之前，我们使用初始参数值来模拟系统。我们使用默认的微分方程求解器(ode45)和默认的求解器选项。当只有两个输入参数调用时，COMPARE将估计完整的初始状态向量，在这种情况下，每个实验一个，即3个实验每个3 × 1的状态向量意味着总共估计9个初始状态。模拟输出和真实输出显示在一个绘图窗口中，可以看到，拟合是不错的，但不如期望的那么好。
CLF比较(z, nlgr);
图3:初始两体模型的真实输出与模拟输出的比较。
参数估计
为了提高拟合，接下来对四个参数进行估计。我们在估计阶段使用来自第一个和最后一个实验的数据，并保留来自第二个实验的数据用于纯粹的验证目的。
opt = nlgreyestOptions(“显示”，“上”）;NLGR = nlgreyest(getexp(z， [1 3])， NLGR, opt);
估计二体模型的性能
为了研究估计模型的性能，最后对其进行了仿真。通过裁剪nlgr的初始状态结构数组，可以完全指定每个实验中估计哪些状态，例如COMPARE。让我们在这里定义和使用一个结构，其中初始状态x1(0)和x2(0)在实验1中估计，x2(0)在实验2中估计，x3(0)在实验3中估计。通过这种修改，测量结果和模型输出之间的比较将显示在一个绘图窗口中。
nlgr。InitialStates = struct(“名字”getinit (nlgr“名字”)，.．.“单位”getinit (nlgr“单位”)，.．.“价值”， 0 (1,3)，“最低”， -Inf(1,3)，.．.“最大”， Inf(1,3)，“固定”，.．.{[假假真];[真假真];[true true false]});比较(z, nlgr, compareOptions(“InitialCondition”，“模型”));
图4:两体模型的真实输出与模拟输出的比较。
特别有趣的是第二次实验的数据结果，这些数据没有用于参数估计。真实系统的动力学在所有实验中都得到了很好的模拟。估计的参数也非常接近用于生成实验数据的参数:
disp (' True估计参数向量'）;
估计参数向量
Ptrue = [400;2 e3;0.0001;1700);流(' %10.5f %10.5f\n'[ptrue ';getpvec (nlgr) '));
400.00000 402.11398 2000.00000 1986.78107 0.00010 0.00010 1700.00000 1705.01863
通过最后使用PRESENT命令，我们可以获得关于估计模型的附加信息:
礼物(nlgr);
nlgr =由'twobodies_c' (MEX-file)定义的连续时间非线性灰盒模型:dx/dt = F(t, u(t)， x(t)， p1，…， p4) y(t) = H(t, u(t)， x(t)， p1，…，p4) + e(t) with 1 input(s), 3 state(s), 1 output(s), and 4 free parameter(s) (out of 4). Inputs: u(1) Exogenous force(t) [N] States: Initial value x(1) Position of moving body(t) [m] xinit@exp1 0 (estimated) in [-Inf, Inf] xinit@exp2 0 (estimated) in [-Inf, Inf] xinit@exp3 0 (fixed) in [-Inf, Inf] x(2) Velocity of moving body(t) [m/s] xinit@exp1 0 (fixed) in [-Inf, Inf] xinit@exp2 0 (estimated) in [-Inf, Inf] xinit@exp3 0 (fixed) in [-Inf, Inf] x(3) Dry friction force between the bodies(t) [N] xinit@exp1 0 (fixed) in [-Inf, Inf] xinit@exp2 0 (fixed) in [-Inf, Inf] xinit@exp3 0 (estimated) in [-Inf, Inf] Outputs: y(1) Position of moving body(t) [m] Parameters: Value p1 Mass of the moving body [m] 402.114 (estimated) in [-Inf, Inf] p2 Sliding friction force coefficient [kg/s] 1986.78 (estimated) in [-Inf, Inf] p3 Distance-related dry friction parameter [m] 0.000104699 (estimated) in [-Inf, Inf] p4 Force-related dry friction parameter [N] 1705.02 (estimated) in [-Inf, Inf] Name: Two body system Status: Model modified after estimation. More information in model's "Report" property. Model Properties
结论
本例描述了在进行IDNLGREY建模时如何使用多个实验数据。可以使用任意数量的实验，并且对于每个这样的实验，都可以在NLGREYEST、COMPARE、PREDICT等中完全指定要估计的初始状态。
额外的信息
有关使用系统识别工具箱™识别动态系统的更多信息，请访问系统识别工具箱产品信息页面。
您也可以从以下列表中选择一个网站: