计算机仿真技术实验报告实验三.docx
- 文档编号:17594398
- 上传时间:2023-07-26
- 格式:DOCX
- 页数:22
- 大小:41.79KB
计算机仿真技术实验报告实验三.docx
《计算机仿真技术实验报告实验三.docx》由会员分享,可在线阅读,更多相关《计算机仿真技术实验报告实验三.docx(22页珍藏版)》请在冰点文库上搜索。
计算机仿真技术实验报告实验三
计算机仿真技术实验报告
实验三利用数值积分算法的仿真实验
实验三利用数值积分算法的仿真实验
实验目的
1)熟悉MATLAB勺工作环境;
2)掌握MATLAB勺.M文件编写规则,并在命令窗口调试和运行程序;
3)掌握利用欧拉法、梯形法、二阶显式Adams法及四阶龙格库塔法构建系统仿真模型的方法,并对仿真结果进行分析。
实验内容
系统电路如图2.1所示。
电路元件参数:
直流电压源E=1V,电阻R=10门,电感
L=0.01H,电容C二WF。
电路元件初始值:
电感电流iL(0)=0A,电容电压uc(0)=0V。
系统输出量为电容电压uc(t)。
连续系统输出响应uc(t)的解析解为:
(2-1)
uc(t)二Us(1—e_at(cos■tsinta/J)
其中,
a=
2L
图2.1RLC串联电路
三、要求
1)利用欧拉法、梯形法、二阶显式Adams法及显式四阶Runge-Kutta法构建系统仿真模型,并求出离散系统的输出量响应曲线;
2)对比分析利用欧拉法、梯形法、二阶显式Adams法及显式四阶Runge-Kutta法构建系统仿真模型的仿真精度与模型运行的稳定性问题;
3)分别编写欧拉法、梯形法、二阶显式Adams法及显式四阶Runge-Kutta法的.m函数文
件,并存入磁盘中。
.m函数文件要求输入参数为系统状态方程的系数矩阵、仿真时间及仿真步长。
编写.m命令文件,在该命令文件中调用已经编写完成的上述.m函数文件,完成
仿真实验;
4)subplot和plot函数将输出结果画在同一个窗口中,每个子图加上对应的标题。
4.实验原理
(1)连续系统解析解
连续系统输出响应uc(t)的解析解为:
uc(t)二Us(1-e^t(costsintx/,))
其中,
X页,
(2)原系统的传递函数
根据所示电路图,我们利用电路原理建立系统的传递函数模型,根据系统的传递函数是在零初始条件下输出量的拉普拉斯变换与输入量的拉普拉斯变换之比,可得该系统的传递函数:
G(s)二
Uc(s)
E(s)
1/LC
s2R/Ls1/LC
(3)系统的仿真模型
在连续系统的数字仿真算法中,较常用的有欧拉法、梯形法、二阶显式Adams法及显式四阶Runge-Kutta法等。
欧拉法、梯形法和二阶显式Adams法是利用离散相似原理构造的仿真算法,而显式四阶Runge-Kutta法是利用Taylor级数匹配原理构造的仿真算法。
对于线性系统,其状态方程表达式为:
X(to)=X0
;X(t)=Ax(t)+Bu(t)
W)二Cx(t)Du(t)
其中:
xhX't)X2(t)…xn(t)T是系统的n维状态向量
u(t)=Ui(t)比⑴…Um(t)T是系统的m维输入向量
y(t)二yi(t)y2(t)…yr(t)T是系统的r维输出向量
A为nn阶参数矩阵,又称动态矩阵,B为nm阶输入矩阵,C为rn阶输出矩阵,D为rm阶交联矩阵。
根据图所示电路,系统状态方程模型:
x(t)二Ax(t)BEy(t)二ex(t)
式中,状态变量x=[x1,x2]T二[iL,uC]T,输出变量y(t)=uC,系数矩阵为:
A「―R/L-1/UD刁/Ll厂rA=|,B=|,C=101」。
]1/C0」]0一
(1)欧拉法
利用前向欧拉法构建线性系统的仿真模型为:
xm1二xm*xmh二—Ahxm-hBum
ym1二Cxm-1'DUm-1
式中,h为积分步长,1为单位矩阵。
利用后向欧拉法构建线性系统的仿真模型为:
.-1
Xm・1二xm-xm』二1-AhXghBUm1
.ym・1二Cxm-1'Dum-1
对于前向欧拉法,系数矩阵为:
Az=1hA,Bz二hB,Cz二C,D=0。
对于后向
欧拉法,系数矩阵为:
£=1-hA',Bz二1-hA一1hB,Cz二C。
(2)梯形法
利用梯形法构建线性系统的仿真模型为:
xm单
xm1
y^i
匚h/
rhj
h口1
1—一A
<2丿
1+—A
<2丿
xm+㊁B(Um+U^1)
xm
号Xm-Xm.1
CXm1DUm1
对图所示的系统,利用梯形法构造的系统差分方程具有形式:
Xm1二AtXm2BtE
ym二CtXm-1
二1
其系数矩阵为:
A=H-hA/21hA/2,,
Bt二1-hA/2JhB2,Ct二C,D=0。
(3)二阶显式Adams法
利用二阶显式Adams法构建线性系统的仿真模型为:
h--]
Ix羽=x+一23F—16Fq+5F2
1m12-mm_1m-2
』m+1=Cxm41+DUm+1
Fm=Axm-BUm
式中:
Fm」二AxmdBUmJ
Fm上二Axm/BUm/
二阶显式Adams法为多步计算方法,利用多步计算方法对系统进行仿真时,需要与之具有相同计算精度的单步计算方法辅助计算。
二阶显式Adams法的计算精度为二阶,可以采用梯形法或改进的Euler法等辅助计算。
利用改进的Euler法构建线性系统的仿真模型为:
k1=f(tm,Xm)=Axm+BE
k2二f(tmh,xmhkj=Axmk[hBE
xm
号k1k2)
其中,m=0,1。
由式计算出片和x2后,便可以转入由二阶显式Adams法构造的离散系统模型计算,即系统差分方程。
其计算方程为:
Fm=AxmBE
Fm二Axm_1BE
Fm^=AxmJ2'BE
Xm-1
Xm
—23Fm-16F
12
5F
m_2
ym1二Cxm1
(4)显式四阶Runge-Kutta法
利用显式四阶Runge-Kutta法构建线性系统的仿真模型为:
k[=
f(tm,
xm)
=Axm+BE
k2=
f(tm
h
+—
2
xm
+hk1)
=A(Xm
+k-ih/
2)+
BE
k3=
f(tm
+h,
2
xm
+hk2)
=A(Xm
+k2h/
2)+
BE
k4=
f(tm
+h,
Xm
+hk3)=
A(Xm+
k3h)+
BE
Xm41
=xm
h+-
6
(k1
+2k2+
2k3+k4)
[_ym舟
=Cx
mH!
5.实验过程
1■实验程序
(1)前向欧拉法
function[]=RLC(R丄,C,U,t,h)R=10;
L=0.01;
C=1.0e-6;
U=1;
t=0.01;
h=2.0e-4;
m=fix(t/h);
n=2;
A=[-R/L-1/L;1/C0];
B=[1/L;0];
D=[01];
E=[10;01];
%前向欧拉法%
fori=1:
1:
n
x1(1:
n,1)=0;
end
fork=1:
m
x1(1:
n,k+1)=x1(1:
n,k)+(A*x1(1:
n,k)+B)*h;
end
fork=1:
1:
m
y1(k)=D*x1(1:
n,k);
end
%解析解%
p=R/(2*L);
w=sqrt(1/(L*C)-(R/(2*L)F2);
fork=1:
1:
m
y(k)=U*(1-exp(-p*(k-1)*h)*(cos(w*(k-1)*h)+sin(w*(k-1)*h)*p/w));end
%输出曲线%
fork=1:
1:
m
t(k)=(k-1)*h;
end
subplot(2,3,1),plot(t,y,'g',t,y1,'r')
legend('y解析解,','y1前向欧拉')
title('前向欧拉法')
(2)后向欧拉法
function[]=RLC(R,L,C,U,t,h)
R=10;
L=0.01;
C=1.0e-6;
U=1;
t=0.01;
h=2.0e-4;
m=fix(t/h);
n=2;
A=[-R/L-1/L;1/C0];
B=[1/L;0];
D=[01];
E=[10;01];
%后向欧拉法%
fori=1:
1:
n
x2(1:
n,1)=0;
end
A1=inv(E-A*h);
fork=1:
m
x2(1:
n,k+1)=A1*(x2(1:
n,k)+B*h);
end
fork=1:
1:
m
y2(k)=D*x2(1:
n,k);
end
%解析解%
p=R/(2*L);
w=sqrt(1/(L*C)-(R/(2*L)F2);
fork=1:
1:
m
y(k)=U*(1-exp(-p*(k-1)*h)*(cos(w*(k-1)*h)+sin(w*(k-1)*h)*p/w));end
%输出曲线%
fork=1:
1:
m
t(k)=(k-1)*h;
end
subplot(2,3,2),plot(t,y,'g',t,y2,'r')
legend('y解析解,','y2后向欧拉')title('后向欧拉法')
(3)梯形法
function[]=RLC(R,L,C,U,t,h)
R=10;
L=0.01;
C=1.0e-6;
U=1;
t=0.01;
h=2.0e-4;
m=fix(t/h);
n=2;
A=[-R/L-1/L;1/C0];
B=[1/L;0];
D=[01];
E=[10;01];
%梯形法%
fori=1:
1:
n
x3(1:
n,1)=0;
end
A2=inv(E-A*h/2);
fork=1:
m
x3(1:
n,k+1)=A2*(x3(1:
n,k)+B*h+A*x3(1:
n,k)*h/2);end
fork=1:
1:
m
y3(k)=D*x3(1:
n,k);
end
%解析解%
p=R/(2*L);
w=sqrt(1/(L*C)-(R/(2*L)F2);
fork=1:
1:
m
y(k)=U*(1-exp(-p*(k-1)*h)*(cos(w*(k-1)*h)+sin(w*(k-1)*h)*p/w));end
%输出曲线%
fork=1:
1:
m
t(k)=(k-1)*h;
end
subplot(2,3,3),plot(t,y,'g',t,y3,'r')
legend('y解析解,','y3梯形法')
title('梯形法')
(4)二阶显式Adams法
function[]=RLC(R,L,C,U,t,h)
R=10;
L=0.01;
C=1.0e-6;
U=1;
t=0.01;
h=2.0e-4;
m=fix(t/h);
n=2;
A=[-R/L-1/L;1/C0];
B=[1/L;0];
D=[01];
E=[10;01];
%二阶显示Adams法%
fori=1:
1:
n
x4(1:
n,1)=0;
end
fork=1:
m
x4(1:
n,k+1)=A2*(x4(1:
n,k)+B*h+A*x4(1:
n,k)*h/2);
end
fork=3:
m
fm1=23*(A*x4(1:
n,k)+B);
fm2=-16*(A*x4(1:
n,k-1)+B);
fm3=5*(A*x4(1:
n,k-2)+B);
x4(1:
n,k+1)=x4(1:
n,k)+(fm1+fm2+fm3)*h/12;
end
fork=1:
1:
m
y4(k)=D*x4(1:
n,k);
end
%解析解%
p=R/(2*L);
w=sqrt(1/(L*C)-(R/(2*L)F2);
fork=1:
1:
m
y(k)=U*(1-exp(-p*(k-1)*h)*(cos(w*(k-1)*h)+sin(w*(k-1)*h)*p/w));end
%输出曲线%
fork=1:
1:
m
t(k)=(k-1)*h;
end
subplot(2,3,4),plot(t,y,'g',t,y4,'r')
legend('y解析解,','y4Adams法')title('二阶显式Adams法')
(5)四阶Runge-Kutta法
function[]=RLC(R,L,C,U,t,h)
R=10;
L=0.01;
C=1.0e-6;
U=1;
t=0.01;
h=2.0e-4;
m=fix(t/h);
n=2;
A=[-R/L-1/L;1/C0];
B=[1/L;0];
D=[01];
E=[10;01];
%四阶Runge-Kutta法%
fori=1:
1:
n%状态变量初值
x5(1:
n,1)=0;
end
fork=1:
m
x5(1:
n,k+1)=A2*(x5(1:
n,k)+B*h+A*x5(1:
n,k)*h/2);
end
fork=1:
1:
m
k1=A*x5(1:
n,k+1);
k2=A*(x5(1:
n,k+1)+h*k1/2);
k3=A*(x5(1:
n,k+1)+h*k2/2);
k4=A*(x5(1:
n,k+1)+h*k3);
x5(1:
n,k+1)=x5(1:
n,k+1)+h*(k1+2*k2+2*k3+k4)./6;
end
fork=1:
1:
m
y5(k)=D*x5(1:
n,k);
end
%解析解%
p=R/(2*L);w=sqrt(1/(L*C)-(R/(2*L)F2);
fork=1:
1:
m
y(k)=U*(1-exp(-p*(k-1)*h)*(cos(w*(k-1)*h)+sin(w*(k-1)*h)*p/w));end
%输出曲线%
fork=1:
1:
m
t(k)=(k-1)*h;
endsubplot(2,3,5),plot(t,y,'g',t,y5,'r')
legend('y解析解,','y5Runge-Kutta法')title('显式四阶Runge-Kutta法')
2■仿真图形
取积分步长h=2*10-4s,可以得到以下几个仿真图形:
(1)前向欧拉法
1x10
-1
-2
-3
-4
-5
0
0.001
y解析解,y1前向欧拉
0.002
0.003
i*
0.004
0.005
0.006
r
0.007
0.008
0.009
0.01
(2)后向欧拉法
后向欧拉法
1.6
1.4
1.2
1
0.8
0.6
0.4
0.2
y解析解,
-y2后向欧拉“
0.0010.0020.003
0.0040.005
0.006
0.0070.008
0.0090.01
(3)梯形法
梯形法
2
1.8
1.6
1.4
1.2
y解析解,
y3梯形法
1
0.8
0.6
0.4
0.2
0
0
0.001
0.002
0.003
0.004
0.005
0.006
0.007
0.008
0.009
0.01
(4)二阶显式Adams法
0.5
„„26
X10
二阶显式Adams
y解析解,y4Adams法
-0.5
-1
-1.5
-2
-2.5
-3
x1016
1
前向欧拉法
y解析解,y1前向欧拉
后向欧拉法
-1
-2
-3
-4
-50
0.005
0.01
1.5
0.5
00
0.005
0.001
0.002
0.003
0.004
0.005
y解析解,y2后向欧拉
梯形法
1.5
0.5
0.01
0.005
0.006
0.007
0.008
y解析解,y3梯形法
0.01
0.009
0.01
(5)四阶Runge-Kutta法
四阶Runge-Kutta
y解析解,
y5龙格库塔
显式四阶Runge-Kutta法
y解析解,
y5龙格库塔
-1
0
0.005
1.5
1X1016前向欧y%法析解
y1前向欧拉
后向欧拉法
0
-1
2
2
梯形法
y解析解,
1.5
y2后向欧拉
x10
■■二
阶显
式Adams法
y解析解,
-2
1
y3梯形法
1.5
y解析解,
-3
1
0
y4Adams法
-4
0.5
-1
1
一
-5
0
0.0050.01
-2
0.5
/1_
0
0
0.005-3
0.01
0
一
11I\
0
0.005
0.01
-4
1\F
0
0.002
0.004
0.006
0.0080.01
-0.5
-10
0.001
0.002
0.003
0.004
0.005
0.006
0.007
0.008
0.009
0.01
6.实验结论
1.从仿真的稳定性看,当选取不同的积分步长时,欧拉法稳定性最低,梯形法稳定性其次,而显式四阶Runge-Kutta法、二阶显示Adams法稳定性较好。
2.从仿真的难易性看,欧拉法为单步计算法,用到一个过去的值,计算起来比较简单。
而梯形法则是用两条折线所谓面积来近似,与欧拉法相比较为困难。
二阶显示Adams法需
要知道k个初始值,不能自起步,二次函数很复杂,因此此方法较复杂。
而显式四阶
Runge-Kutta法建模最为复杂,仿真时间也较长。
3.从仿真的精度看,选取不同的积分步长,就可以得出,欧拉法的精度最低,梯形法稍好,显式四阶Runge-Kutta法,二阶显示Adams法的函数曲线与真实曲线较为接近,精度最高。
综上所述,这几种仿真方法各有优劣。
欧拉法较为简单,但精度较低。
Runge-Kutta
法,二阶显示Adams法较为复杂,但精度和稳定性都比较高。
因此,在实际仿真中应根据问题的需求选择不同的方法。
对精度要求较小时,可选用欧拉法,梯形法;当对精度要求较高时,就应选择显式四阶Runge-Kutta法、二阶显示Adams法。
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 计算机仿真 技术 实验 报告