IAuto软件快速进行二阶(或高阶)微分方程组的仿真及精度控制

在上一篇文章中，介绍了如何通过IAuto求解一阶微分方程组的仿真曲线或数据，但一阶微分方程组相对比较简单，场景也比较单一，在更加复杂的场景中经常会遇到二阶甚至三阶微分方程组的场景，如下案例：

其中最高为二阶，x1、x2、x3分别为函数f1(t)、f2(t)、f3(t)，这里需要输出y的仿真数据曲线，我们用Matlab也可以搞定求出y的数据，但本着”能不用Matlab就不用Matlab，国产软件IAuto也很香”的原则，在IAuto快速求解：

第一步：先新建一个.iax文件（记住一定要优先保存文件路径），并创建一个数学曲线空白页；

第二步：在表达式内容区定义一个dx，并在更多设置中设置x的属性，x范围[0,10]，曲线点数1000；

注意这里变量dx为（x结束值-x的开始值）/曲线点数；

第三步：在表达式内容中添加如下代码，用来定义初始条件和输出值，其中自定义的全局变量IA.X1、IA.X2、IA.X3分别对应微分方程组的x1、x2、x3；IA.Dx1、IA.Dx2、IA.Dx3分别对应x1、x2、x3的一阶导数，IA.DDx1、IA.DDx2、IA.DDx3分别对应x1、x2、x3的二阶导数；

let dx = 10/1000;if(x==0){  IA.X1 = 1;//X1的初始值  IA.Dx1 = 0;//X1的一阶初始值  IA.X2 = 1;//X1的初始值  IA.x2 = 0;//X1的一阶初始值  IA.X3 = 3;//X1的初始值  IA.x3 = 1;//X1的一阶初始值}else{//核心业务代码区};let y = IA.X1+2*IA.X2+IA.X3;y;//输出y

第四步：把微分方程组中每个方程的最高阶放置左边，其他项放在右边，之后将方程中的参数用自定义的全局变量替代：

 IA.DDx1=3-IA.X3-IA.X1-IA.Dx2-2*IA.Dx1;//第一个方程 IA.Dx2 = IA.X2-2-IA.Dx3;//第二个方程 IA.DDx3 = IA.Dx2+IA.X2-1-IA.X3;//第三个方程

将上述三个方程的脚步放入else{}代码段中：

let dx = 10/1000;if(x==0){  IA.X1 = 1;//X1的初始值  IA.Dx1 = 0;//X1的一阶初始值  IA.X2 = 1;//X1的初始值  IA.Dx2 = 0;//X1的一阶初始值  IA.X3 = 3;//X1的初始值  IA.Dx3 = 1;//X1的一阶初始值}else{    IA.DDx1=3-IA.X3-IA.X1-IA.Dx2-2*IA.Dx1;//第一个方程  IA.Dx2 = IA.X2-2-IA.Dx3;//第二个方程  IA.DDx3 = IA.Dx2+IA.X2-1-IA.X3;//第三个方程};let y = IA.X1+2*IA.X2+IA.X3;y;//输出y

第五步：分别对第四步中的IA.DDx1(x1的二阶)、IA.Dx2(x2的一阶)、IA.DDx3(x3的二阶)进行二次积分或者一次积分得到IA.X1、IA.X2、IA.X3;

这里积分的次数取决导数的深度，比如二阶需要两次积分才能得到原函数，三阶需要三次积分才能得到原函数；

 //对xx1(t)二次积分得到X1  IA.Dx1 += IA.DDx1*dx;  IA.X1 += IA.Dx1*dx;   //对x2(t)一次积分得到X2  IA.X2 += IA.Dx2*dx;  //对xx3(t)一次积分得到X3  IA.Dx3 += IA.DDx3*dx;  IA.X3 += IA.Dx3*dx;

将上述积分的代码添加在第四步的脚步代码中如下：

let dx = 10/1000;if(x==0){  IA.X1 = 1;//X1的初始值  IA.Dx1 = 0;//X1的一阶初始值  IA.X2 = 1;//X1的初始值  IA.Dx2 = 0;//X1的一阶初始值  IA.X3 = 3;//X1的初始值  IA.Dx3 = 1;//X1的一阶初始值}else{  IA.DDx1=3-IA.X3-IA.X1-IA.Dx2-2*IA.Dx1;//第一个方程  IA.Dx2 = IA.X2-2-IA.Dx3;//第二个方程  IA.DDx3 = IA.Dx2+IA.X2-1-IA.X3;//第三个方程  //对xx1(t)二次积分得到X1  IA.Dx1 += IA.DDx1*dx;  IA.X1 += IA.Dx1*dx;   //对x2(t)一次积分得到X2  IA.X2 += IA.Dx2*dx;  //对xx3(t)一次积分得到X3  IA.Dx3 += IA.DDx3*dx;  IA.X3 += IA.Dx3*dx;};let y = IA.X1+2*IA.X2+IA.X3;y;//输出y

第六步：点击“绘制静态图表”按钮，即可得到y的仿真数据和仿真曲线；

这里值得注意的是dx的仿真精度，因为每个dx段为10/1000=0.01，所以有一定的仿真误差，如果我们需要将dx的精度变为0.0000001，只需将第五步代码修改如下，其核心业务代码放置在一个for循环中，积分运算时的dx(0.01)由dxn(0.0000001)替代：

let dx = 10/1000;let n = 10000;let dxn = dx/n;if(x==0){  IA.X1 = 1;//X1的初始值  IA.Dx1 = 0;//X1的一阶初始值  IA.X2 = 1;//X1的初始值  IA.Dx2 = 0;//X1的一阶初始值  IA.X3 = 3;//X1的初始值  IA.Dx3 = 1;//X1的一阶初始值}else{for(let i=0;i<n;i++){//控制运算精度  IA.DDx1=3-IA.X3-IA.X1-IA.Dx2-2*IA.Dx1;  IA.Dx2 = IA.X2-2-IA.Dx3;  IA.DDx3 = IA.Dx2+IA.X2-1-IA.X3;  //对xx1(t)二次积分得到X1  IA.Dx1 += IA.DDx1*dxn;  IA.X1 += IA.Dx1*dxn;   //对x2(t)一次积分得到X2  IA.X2 += IA.Dx2*dxn;  //对xx3(t)一次积分得到X3  IA.Dx3 += IA.DDx3*dxn;  IA.X3 += IA.Dx3*dxn;}};let y = IA.X1+2*IA.X2+IA.X3;y;//输出y

除了以上将dx变小为dxn提高精度的方式，我们还可以通过另一种变步长积分的方式提高精度(龙泽-库塔法)，相比dxn这种需要运算n次而言仅需要迭代运算4次，既降低了运算量同时也能达到较高精度效果，其代码如下：

let dx = 10/1000;if(x==0){  IA.X1 = 1;//X1的初始值  IA.Dx1 = 0;//X1的一阶初始值  IA.X2 = 1;//X2的初始值  IA.X3 = 3;//X3的初始值  IA.Dx3 = 1;//X3的一阶初始值}else{let A=[0,1/2,1/2,1];for(let i=0;i<4;i++){  IA.DDx1=3-IA.X3-IA.X1-IA.Dx2-2*IA.Dx1;  IA.Dx2 = IA.X2-2-IA.Dx3;  IA.DDx3 = IA.Dx2+IA.X2-1-IA.X3;  //对xx1(t)二次积分得到X1  IA.Dx1 += IA.DDx1*A[i]*dx/2;  IA.X1 += IA.Dx1*dx/4;   //对x2(t)一次积分得到X2  IA.X2 += IA.Dx2*A[i]*dx/2;  //对xx3(t)一次积分得到X3  IA.Dx3 += IA.DDx3*A[i]*dx/2;  IA.X3 += IA.Dx3*dx/4;}};IA.X1+2*IA.X2+IA.X3;

这两种方式都可以提高精度，第一种运算量更大些，第二种运算量更小些（推荐）；

除了这些，IAuto多功能绘图软件还可以用于各类基础办公绘图场景中；

关注小编，即可免费学习更多关于国产IAuto软件的使用教程！