4.机电系统动态性能的计算机仿真 4.1概述 机电系统计算机仿真是目前对复杂机电系统进行分析的重要手段与方法。在进 行机电系统分析综合与设计工作过程中,除了需要进行理论分析外,还要对系统的 特性进行实验研究。系统性能指标与参数是否达到预期的要求?它的经济性能如 何?这些都需要在系统设计中给出明确的结论。对于那些在实际调试过程中存在很 大风险或实验费用昂贵的系统,一般不允许对设计好的系统直接进行实验,然而没 有经过实验研究是不能将设计好的系统直接放到生产实际中去的,因此就必须对其 进行模拟实验硏究。当然在有些情况下可以构造一套物理模拟装置来进行实验,但 这种方法十分费时而且费用又高,而在有的情况下物理模拟几乎是不可能的。近年 来随着计算机的迅速发展,采用计算机对机电系统进行数学仿真的方法已被人们采 纳。所谓机电系统计算机仿真就是以机电系统的数学模型为基础,借助计算机对机 电系统的动静态过程进行实验研究。这里讲的机电系统计算机仿真是指借助数字计 算机实现对机电系统的仿真分析。这种实验研究的特点是:将实际系统的运动规律 用数学表达式加以描述,它通常是一组常微分方程或差分方程,然后利用计算机来 求解这一数学模型,以达到对系统进行分析研究的目的 对机电系统进行计算机仿真的基本过程包括:首先建立系统的数学模型,因为 数学模型是系统仿真的基本依据,所以数学模型极为重要。然后根据系统的数学模 型建立相应的仿真模型,一般需要通过一定的算法或数值积分方法对原系统的数学 模型进行离散化处理,从而建立起相应的仿真模型,这是进行机电系统仿真分析的 关键步骤;最后根据系统的仿真模型编制相应的仿真程序,在计算机上进行仿真实 验研究并对仿真结果加以分析。 机电系统计算机仿真的应用与发展已经过了近40年的历程,进入20世纪80 年代以来,随着微型计算机技术以及软件技术的飞速发展与广泛应用,使得机电系 统计算机仿真获得了实质性的发展,并使其走进广大的机电系统生产、设计、研究 的第一线。目前在该领域应用最为广泛的软件包就是美国 Mathworks公司开发的 MATLAB语言软件。虽然它最初并非是为机电系统仿真与设计开发的,但是它的强 大的矩阵运算功能与图形处理及绘制能力,以及在 MATLAB平台上开发出来的各种 工具箱,和面向结构图的 SIMULINK系统分析环境,为机电系统计算机仿真提供了 强有力的软件工具,从而为机电系统计算机仿真与辅助设计开辟新的方法与途径提 供了充分的条件。 考虑到部分读者不具备 MATLAB基础知识,因此,在本书的附录一中对 119
119 4．机电系统动态性能的计算机仿真 4．1 概述 机电系统计算机仿真是目前对复杂机电系统进行分析的重要手段与方法。在进 行机电系统分析综合与设计工作过程中，除了需要进行理论分析外，还要对系统的 特性进行实验研究。系统性能指标与参数是否达到预期的要求？它的经济性能如 何？这些都需要在系统设计中给出明确的结论。对于那些在实际调试过程中存在很 大风险或实验费用昂贵的系统，一般不允许对设计好的系统直接进行实验，然而没 有经过实验研究是不能将设计好的系统直接放到生产实际中去的，因此就必须对其 进行模拟实验研究。当然在有些情况下可以构造一套物理模拟装置来进行实验，但 这种方法十分费时而且费用又高，而在有的情况下物理模拟几乎是不可能的。近年 来随着计算机的迅速发展，采用计算机对机电系统进行数学仿真的方法已被人们采 纳。所谓机电系统计算机仿真就是以机电系统的数学模型为基础，借助计算机对机 电系统的动静态过程进行实验研究。这里讲的机电系统计算机仿真是指借助数字计 算机实现对机电系统的仿真分析。这种实验研究的特点是：将实际系统的运动规律 用数学表达式加以描述，它通常是一组常微分方程或差分方程，然后利用计算机来 求解这一数学模型，以达到对系统进行分析研究的目的。 对机电系统进行计算机仿真的基本过程包括：首先建立系统的数学模型，因为 数学模型是系统仿真的基本依据，所以数学模型极为重要。然后根据系统的数学模 型建立相应的仿真模型，一般需要通过一定的算法或数值积分方法对原系统的数学 模型进行离散化处理，从而建立起相应的仿真模型，这是进行机电系统仿真分析的 关键步骤；最后根据系统的仿真模型编制相应的仿真程序，在计算机上进行仿真实 验研究并对仿真结果加以分析。 机电系统计算机仿真的应用与发展已经过了近 40 年的历程，进入 20 世纪 80 年代以来，随着微型计算机技术以及软件技术的飞速发展与广泛应用，使得机电系 统计算机仿真获得了实质性的发展，并使其走进广大的机电系统生产、设计、研究 的第一线。目前在该领域应用最为广泛的软件包就是美国 Mathworks 公司开发的 MATLAB 语言软件。虽然它最初并非是为机电系统仿真与设计开发的，但是它的强 大的矩阵运算功能与图形处理及绘制能力，以及在 MATLAB 平台上开发出来的各种 工具箱，和面向结构图的 SIMULINK 系统分析环境，为机电系统计算机仿真提供了 强有力的软件工具，从而为机电系统计算机仿真与辅助设计开辟新的方法与途径提 供了充分的条件。 考虑到部分读者不具备 MATLAB 基础知识，因此，在本书的附录一中对
MATLAB语言基础知识作了简单的介绍 4.2面向微分方程的连续系统仿真与 MATLAB实现 对于控制系统的最基本数学描述就是微分方程,因此有必要在介绍机电系统计 算机仿真时对这种方法加以介绍。实际上在前面章节已经介绍了几种典型数值积分 方法以及它们的 MATLAB实现。数值积分的数学基础就是微分方程(状态方程), 如果一个系统以微分方程的形式加以描述,就可以利用前面的知识对该系统进行仿 真研究,有关这部分内容本章就不再加以讨论了。这里将介绍如何利用 MATLAB中 为用户所提供的常微分方程解函数ode和对微分方程进行求解计算,以及在面向微 分方程的系统仿真研究中如何应用M函数 4.2.1基于ode函数的面向微分方程的系统仿真 1.常微分方程解函数ode的基本格式 (1)ode函数的基本分类。在 MATLAB5x中ode函数可以分为两类五种形式。 类是应用于非刚性微分方程求解,如ode45:另一类属于求解刚性微分方程的ode 解函数,如ode23s他们的具体形式如下 lodε45函数用于求解非刚性微分方程.它属于中阶龙格一库塔法 2ode23函数用于求解非刚性微分方程,它属于底阶龙格-库塔法 3 odell3函数用于求解非刚性微分方程的变阶法 4odel5s函数用于求解刚性微分方程的变阶法 5ode23s函数用于求解刚性微分方程的低阶法 (2)ode函数的基本格式。上还这五种ode函数的格式基本上一样.这里主要 以最常用的ode45函数为例加以介绍,其结果多数也可用于其他ode函数。ode45的 基本格式主要有以下几种 O [t, yI=ode45('F,tspan, yO 其中:F是ode文件名字符串 tspan=[10 TFINAL]为一个向量,它指定了积分的启始时刻和结束时刻 y0为积分运算的初始条件。 该函数对微分方程系统y=F(y从m0到TNAL时间段进行积分。函数F(y) 返回一个列向量。在解向量Y中的每一行与列向量T的返回时间相对应。为了在指 定的时刻点T0,T1,…, TFINAL得到解,则使用 tspan=[T0,Tl,…, TFINAL] 2 t, y]=ode45('F, tspan, yO, options) 120
120 MATLAB 语言基础知识作了简单的介绍。 4．2 面向微分方程的连续系统仿真与 MATLAB 实现 对于控制系统的最基本数学描述就是微分方程，因此有必要在介绍机电系统计 算机仿真时对这种方法加以介绍。实际上在前面章节已经介绍了几种典型数值积分 方法以及它们的 MATLAB 实现。数值积分的数学基础就是微分方程（状态方程）， 如果一个系统以微分方程的形式加以描述，就可以利用前面的知识对该系统进行仿 真研究，有关这部分内容本章就不再加以讨论了。这里将介绍如何利用 MATLAB 中 为用户所提供的常微分方程解函数 ode 和对微分方程进行求解计算，以及在面向微 分方程的系统仿真研究中如何应用 M 函数。 4．2．1 基于 ode 函数的面向微分方程的系统仿真 1. 常微分方程解函数 ode 的基本格式 （l）ode 函数的基本分类。在 MATLAB 5.x 中 ode 函数可以分为两类五种形式。 一类是应用于非刚性微分方程求解，如 ode45；另一类属于求解刚性微分方程的 ode 解函数，如 ode23s.他们的具体形式如下： 1 ode45 函数用于求解非刚性微分方程．它属于中阶龙格一库塔法； 2 ode23 函数用于求解非刚性微分方程，它属于底阶龙格-库塔法； 3 ode113 函数用于求解非刚性微分方程的变阶法； 4 ode15s 函数用于求解刚性微分方程的变阶法； 5 ode23s 函数用于求解刚性微分方程的低阶法。 （2）ode 函数的基本格式。上还这五种 ode 函数的格式基本上一样．这里主要 以最常用的 ode45 函数为例加以介绍，其结果多数也可用于其他 ode 函数。ode45 的 基本格式主要有以下几种。 ① [t，y]=ode45('F'，tspan，y0) 其中：'F'是 ode 文件名字符串； tspan= [T0 TFINAL］为一个向量，它指定了积分的启始时刻和结束时刻； y0 为积分运算的初始条件。 该函数对微分方程系统 y = F(t, y) 从 T0 到 TFINAL 时间段进行积分。函数 F(t, y) 返回一个列向量。在解向量 Y 中的每一行与列向量 T 的返回时间相对应。为了在指 定的时刻点 T0，T1，…，TFINAL 得到解，则使用 tspan=[T0，T1，…，TFINAL]。 ② [t , y]＝ode45（'F'， tspan， y0，options）
该函数的解如同前一格式,只是在这里用选项 options中的参数来代替前一格式中的 缺省值, options中的参数变量由函数 odeset予以建立。最常使用的选项参数是标量 相对容差 RelTo"(缺省值为le-3)与绝对容差向量' Abs To'(缺省的全部向量元素均 为1e-6) (3) odeset函数的基本格式。 odeset函数用来建立或改变 options选项结构,其 基本格式为 ptions=odeset (,, valuel, 'name2,, value 通过该 odeset函数建立积分器选项 options的结构,在 options中所命名的属性具有 指定值,对于任何未指定的属性均为缺省值。 options=odeset (oldopts, 'namel,, valuel,.) 通过该 odeset函数改变现存选项结构 adopts 通过该 odeset函数将现有选项结构 adopts与新选项结构 newts结合在一起,新选 项中的属性将覆盖对应原有属性 2.机电系统仿真的ode函数实现 由od函数可知,当已知系统由微分方程y=F(y)加以描述时,并建立以F 作为ODE文件名,则可应用 It, yI=ode45(F, tspan, yo) 实现对系统的求解分析。这里所讲的ODE文件实际就是M函数文件,因此对于以 高阶微分方程加以描述的系统,首先应将其写成一组一阶微分方程的形式,当然如 果系统是以状态方程描述的,即可直接对其编写M函数文件。函数文件的特征就是 文件第一行一定是以 function开始,而且在函数文件中的变量均为局部变量,同时 由函数文件构造的功能函数可以像 MAT-LAB的其他函数一样进行调用。实际上在 ode函数[t,y]=ode45(F, tspan,y0)中,对以F作为文件名的ODE文件就是这 类函数文件。下面通过实例说明对系统仿真过程中的具体应用。 例4-1之已知一三阶系统,其状态方程表达式为 X=AX+Bl 其中 4910 A=12-910B=4 24-24-18 12 输入为u=100,试应用ode函数求解该系统
121 该函数的解如同前一格式，只是在这里用选项 options 中的参数来代替前一格式中的 缺省值，options 中的参数变量由函数 odeset 予以建立。最常使用的选项参数是标量 相对容差'RelTol'（缺省值为 1e-3）与绝对容差向量'AbsTol'（缺省的全部向量元素均 为 1e-6）。 （3）odeset 函数的基本格式。odeset 函数用来建立或改变 options 选项结构，其 基本格式为 options=odeset（'namel'，valuel，'name2'，value2，…） 通过该 odeset 函数建立积分器选项 options 的结构，在 options 中所命名的属性具有 指定值，对于任何未指定的属性均为缺省值。 options＝odeset（oldopts，'namel'， valuel，…） 通过该 odeset 函数改变现存选项结构 oldopts。 options=odeset（oldopt， newopts） 通过该 odeset 函数将现有选项结构 oldopts 与新选项结构 newopts 结合在一起，新选 项中的属性将覆盖对应原有属性。 2．机电系统仿真的 ode 函数实现 由 ode 函数可知，当已知系统由微分方程 y = F(t, y) 加以描述时，并建立以'F' 作为 ODE 文件名，则可应用 [t，y]=ode45('F'，tspan，y0) 实现对系统的求解分析。这里所讲的 ODE 文件实际就是 M 函数文件，因此对于以 高阶微分方程加以描述的系统，首先应将其写成一组一阶微分方程的形式，当然如 果系统是以状态方程描述的，即可直接对其编写 M 函数文件。函数文件的特征就是 文件第一行一定是以 function 开始，而且在函数文件中的变量均为局部变量，同时 由函数文件构造的功能函数可以像 MAT-LAB 的其他函数一样进行调用。实际上在 ode 函数[t，y]＝ode45（'F'，tspan，y0）中，对以'F'作为文件名的 ODE 文件就是这 类函数文件。下面通过实例说明对系统仿真过程中的具体应用。 例 4-1 之已知一三阶系统，其状态方程表达式为 X = AX + BU 其中 − − − − = 24 24 18 12 9 10 4 9 10 A ， = 12 4 2.8 B 输入为 u= 100，试应用 ode 函数求解该系统
解:由于该系统是以状态方程形式给出的,因此即可直接对其编写M函数文 件 +10*x(3)+4*u;24*x(1)-24*x(2)-18*x(3)+12*山 该函数文件对所研究系统的状态方程加以描述并形成rm文件以便在ode解函数中 使用后即可利用ode解函数对系统进行求解分析。考虑状态的零初值条件,即 0.5 1.5 2.5 图4-1系统状态曲线 [t,x]=ode45(r,[03],xO) plot(t, x) 对其同样可存一M文件,如li32m文件,这样即可在命令窗口执行该文件。在运行 过程中它会自动调用相关m函数文件,实现预定的工作,图41给出了相应的运行 结果 根据ode函数用途可知,该函数应用于求解常微分方程,因此对于含有输出方 程的系统状态空间表达式就不能直接利用ode函数求解系统输出响应。考虑到ode 函数返回变量的向量格式,可以根据输出方程的具体要求,由状态变量的线性组合 直接得到所要求的系统输出。利用ode函数可以方便地对以微分方程形式描述的系 统进行仿真分析,然而ode解函数并非是万能的,在应用ode函数时对其进行干预 般是不方便的,因此在很多情况下,对以微分方程形式给出的系统,还必须采用 前面所介绍的各种数值积分方法进行仿真分析。 122
122 解： 由于该系统是以状态方程形式给出的，因此即可直接对其编写 M 函数文 件。 function xd= rr（t，x） u＝100； xd=[-14*x(1)+9*x(2)+10*x(3)+2.8*u；12*x(1)-9*x(2) +10*x(3)+4*u；24*x(1)-24*x(2)-18*x(3)+12*u]； 该函数文件对所研究系统的状态方程加以描述并形成 rr.m 文件以便在 ode 解函数中 使用后即可利用 ode 解函数对系统进行求解分析。考虑状态的零初值条件，即 图 4-1 系统状态曲线 clear x0=[0；0；0]； [t，x]=ode45('rr'，[0 3]，x0)； plot(t，x) 对其同样可存一 M 文件，如 li32.m 文件，这样即可在命令窗口执行该文件。在运行 过程中它会自动调用相关 rr 函数文件，实现预定的工作，图 4-1 给出了相应的运行 结果。 根据 ode 函数用途可知，该函数应用于求解常微分方程，因此对于含有输出方 程的系统状态空间表达式就不能直接利用 ode 函数求解系统输出响应。考虑到 ode 函数返回变量的向量格式，可以根据输出方程的具体要求，由状态变量的线性组合 直接得到所要求的系统输出。利用 ode 函数可以方便地对以微分方程形式描述的系 统进行仿真分析，然而 ode 解函数并非是万能的，在应用 ode 函数时对其进行干预 一般是不方便的，因此在很多情况下，对以微分方程形式给出的系统，还必须采用 前面所介绍的各种数值积分方法进行仿真分析
4.2.2基于M函数的面向微分方程的系统仿真 M函数是 MATLAB中所特有的一类函数,前面介绍ode函数时已经使用了M 函数。这里将介绍采用M函数建立描述系统的数学模型并在此基础上编制相应的程 序。根据需要调用该M文件,实现所要求的运算。下面通过实例加以说明 例42已知系统状态方程为 05572-07814x 0.7814 0 应用数值积分方法(预估校正梯形法)实现对系统状态的仿真研究 解:首先建立描述该系统数学模型的M函数文件,相应的程序如下 function dfrrl (X, u) A=[0.5572-0.7814:0.78140]: df=A*X+B*u 将该函数文件存盘取作rl.m文件。然后根据预估校正梯形公式编制相应的程序。 h=00l;t0=0;tf=15 X=[00];X1=00 T=(tf-tO)h for i=l: T KI=r(X,1); X= X+h*KI K2=rrl (X, 1): X=X+h*(K1+K2)/2 XI=[XI: XI plot(XI) 在该程序段中两次调用了rr1m函数文件,在第一次调用rlm后,对各状态的 计算将产生新值,因此K2就是经过预估计算所得到的相关值,从而即可应用梯形法 实现对系统状态的仿真分析。图4-2给出了系统状态响应仿真曲线。 123
123 4．2．2 基于 M 函数的面向微分方程的系统仿真 M 函数是 MATLAB 中所特有的一类函数，前面介绍 ode 函数时已经使用了 M 函数。这里将介绍采用 M 函数建立描述系统的数学模型并在此基础上编制相应的程 序。根据需要调用该 M 文件，实现所要求的运算。下面通过实例加以说明。 例 4-2 已知系统状态方程为 u x x x x + − − = 0 1 0 7814 0 0 5572 0 7814 2 1 2 1 . . . 应用数值积分方法（预估校正梯形法）实现对系统状态的仿真研究。 解：首先建立描述该系统数学模型的 M 函数文件，相应的程序如下。 function df=rrl（X，u） A= ［-0.5572 -0.7814； 0.7814 0］； B=[1；0] U= 20； df=A*X+B*u 将该函数文件存盘取作 rrl．m 文件。然后根据预估校正梯形公式编制相应的程序。 h=0.0l；t0=0；tf=15； X=[0 0]'； X1=[0 0]； T=(tf-t0)/h for i=l：T K1=rrl（X，1）； X= X+h*K1 K2=rrl（X，1）； X= X+h*（K1+K2）/2 X1=[X1；X']； end plot（X1） 在该程序段中两次调用了 rr1.m 函数文件，在第一次调用 rrl.m 后，对各状态的 计算将产生新值，因此 K2 就是经过预估计算所得到的相关值,从而即可应用梯形法 实现对系统状态的仿真分析。图 4-2 给出了系统状态响应仿真曲线