第八章-RLC电路与常微分方程的ppt课件

1、第八章 RLC电路与常微分方程的解法,8.1 RC电路与常微分方程的欧拉解法 RC电路: K 2 1 R C 先把开关K接通“1” 端,电容C充满电后再把开关K接通“2”端,则这时电容C放电过程满足方程: 即电容C上的电量是时间t的函数,满足以上微分方程,如果设: =RC, t=0时刻电容所带电量为Q0 则有: 考虑数值微分问题: 已知: 求f(x) 在xn 点的导数. 可以: 或,微分方程化为一般形式: 把时间t等间隔离散化: 其中: 做如下近似: 由方程得,即: 记: 则得到解微分方程的欧拉法递推公式: 对于RC电路: 令,得到: 方程的解析解,微分方程化为一般形式: 把时间 t 等间隔离

2、散化: 其中: 欧拉(Euler) 差分公式: 由方程得,即: 记: 则得到解微分方程的欧拉法递推公式: 对于RC电路: 例如,得到,rc(1,6,1,10,欧拉法也可解释为Q(t)在tn处的泰勒展开: 取线性部分: 欧拉方法的截断误差,例: 写出解如下一阶常微分方程的欧拉公式: 得,8.2 RLC电路和改进的欧拉近似法 RLC 电路图: L R C Va K 根据基尔霍夫定律: 由于,得: 由于: 所有: 欧拉法: 把二阶微分方程化成一阶微分方程组,其中t是自变量,Q和I随着t的改变而改变,function Q,I,tt=rlc(Q0,I0,con,T,dt) % RLC电路欧拉解法 Q(1

3、)=Q0;I(1)=I0; R=con(1);L=con(2);C=con(3);V=con(4); tt=0:dt:T; for n=1:length(tt)-1 Q(n+1)=Q(n)+dt*I(n); I(n+1)=I(n)+dt*(V-R*I(n)-Q(n)/C)/L; end plot(tt,Q,r,tt,I,b,rlc(1,0,1,1,1,5,15,0.1,rlc(1,0,1,5,1,5,50,0.1,2. 向后的欧拉方法 方法分为两步: 预估: (一步)校正,或者(k+1步)校正,function Q,I,tt=rlc1(Q0,I0,con,T,dt) % RLC电路向后欧拉解法

4、 Q(1)=Q0;I(1)=I0; R=con(1);L=con(2);C=con(3);V=con(4); tt=0:dt:T; for n=1:length(tt)-1 Q1=Q(n)+dt*I(n); I1=I(n)+dt*(V-R*I(n)-Q(n)/C)/L; Q(n+1)=Q(n)+dt*I1; I(n+1)=I(n)+dt*(V-R*I1-Q1/C)/L; end plot(tt,Q,r-,tt,I,b,rlc1(1,0,1,1,1,5,15,0.1); hold on rlc(1,0,1,1,1,5,15,0.1,3. 改进的欧拉法 方法分两步: 预估: (一步)校正,或(k+

5、1步)校正,function Q,I,tt=rlc2(Q0,I0,con,T,dt) % RLC电路改进欧拉解法 Q(1)=Q0;I(1)=I0; R=con(1);L=con(2);C=con(3);V=con(4); tt=0:dt:T; for n=1:length(tt)-1 Q1=Q(n)+dt*I(n); I1=I(n)+dt*(V-R*I(n)-Q(n)/C)/L; Q(n+1)=Q(n)+dt*(I1+I(n)/2; I(n+1)=I(n)+dt*(V-R*(I1+I(n)/2- (Q1+Q(n)/2/C)/L; end plot(tt,Q,r:,tt,I,b:,RC电路: 向

6、后的欧拉法: 预估: 校正: 改进的欧拉法: 预估: 校正,function Q1,Q2,Q3,tt=rc3(Q0,T,dt,tao) % RC电路欧拉解法 Q1(1)=Q0;Q2(1)=Q0;Q3(1)=Q0; tt=0:dt:T; for n=1:length(tt)-1 Q1(n+1)=Q1(n)-dt*Q1(n)/tao; end for n=1:length(tt)-1 Q=Q2(n)-dt*Q2(n)/tao; Q2(n+1)=Q2(n)-dt*Q/tao; end for n=1:length(tt)-1 Q=Q3(n)-dt*Q3(n)/tao; Q3(n+1)=Q3(n)-d

7、t*(Q+Q3(n)/2/tao; end Qa=Q0*exp(-tt/tao); plot(tt,Qa,b,tt,Q1,r-,tt,Q2,r-,tt,Q3,r:,rc3(1,6,1,10,一般微分方程: 向后的欧拉法: 改进的欧拉法,8.3 龙格-库塔(R-K)方法 对于微分方程: 根据微分中值定理: 即,Q(t) tn tn+1 用tn处Q(t)的导数代替处导数 f(,Q(),则为欧拉法: 用tn+1处Q(t)的导数的估计值代替处导数 f(,Q(),则为向后的欧拉法,即: 用tn和tn+1处Q(t)的导数的估计值的平均代替处导数 f(,Q(),则为改进的欧拉法: 若取多点处斜率(即导数)的

8、加权平均会使误差更小,称为龙格-库塔法,最常用的四阶龙格-库塔法,例1: 求解方程: 梯形法: 四阶龙格-库塔法,有,function y1,y2,xx=rk1(y0,X,dx) % 矩形法和四阶龙格-库塔法 y1(1)=y0;y2(1)=y0; xx=0:dx:X; for n=1:length(xx)-1 y=y1(n)+dx*y1(n); y1(n+1)=y1(n)+dx*(y1(n)+y)/2; end for n=1:length(xx)-1 k1=y2(n); k2=y2(n)+dx*k1/2; k3=y2(n)+dx*k2/2; k4=y2(n)+dx*k3; y2(n+1)=y

9、2(n)+dx*(k1+2*k2+2*k3+k4)/6; end y=exp(xx); plot(xx,y,b,xx,y1,r-,xx,y2,r:,rk1(1,5,1,例2： 四阶的龙格-库塔公式,例3： 四阶的龙格-库塔公式,微分方程组,龙格-库塔公式,例4: 求解阻尼振动方程 首先把它转化为一阶微分方程组,四阶龙格-库塔公式,function x,v,tt=rk2(m,k,c,x0,v0,T,dt) % 四阶龙格-库塔法 x(1)=x0;v(1)=v0; tt=0:dt:T; for n=1:length(tt)-1 k1=v(n); l1=-(c*v(n)+k*x(n)/m; k2=v(

10、n)+dt*l1/2; l2=-(c*(v(n)+dt*l1/2)+k*(x(n)+dt*k1/2)/m; k3=v(n)+dt*k2/2; l3=-(c*(v(n)+dt*l2/2)+k*(x(n)+dt*k2/2)/m,k4=v(n)+dt*k3; l4=-(c*(v(n)+dt*l3)+k*(x(n)+dt*k3)/m; x(n+1)=x(n)+dt*(k1+2*k2+2*k3+k4)/6; v(n+1)=v(n)+dt*(l1+2*l2+2*l3+l4)/6; end plot(tt,x); h=line(Color,1 0 0, Marker,.,MarkerSize,20,EraseMode,xor); for i=1:length(tt) set(h,Xdata,tt(i),Ydata,x(i); pause(dt); end,x,v,tt=rk2(10,10,2