Gauss-Seidel 迭代法线性方程组求解.doc

、编程实现以下科学计算算法，并举一例应用之。 （参考书籍精通科学计算 ，王正林等著，电子工业出版社，年）“Gauss-Seidel 迭代法线性方程组求解”（1）迭代解法的基本思想：根据给定方程组，设计出一个迭代公式，构造一数组的序列xi0,代入迭代公式，计算出xi1,在代入迭代公式，经过k次迭代运算后得到xik,若xik收敛于某一极限数组xi，则xi就是方程组的近似解。迭代过程本质上就是计算极限的过程，一般不能得到精确解。但迭代的优点是程序简单，适合大型方程组求解，然而，缺点是要判断迭代是否收敛和收敛速度的问题。(2)算法说明：Gauss-Seidel迭代法与简单迭代法类似，只是迭代公式有所改进。简单迭代法：；Gauss-Seidel迭代法：；设方程组Axb，其中A和b中的元素都为常数，且A为非奇异，则A分可写成：A=D-L-U。其中D上网意义同Jacobi迭代法，L为下三角矩阵，U为上三角矩阵，他的迭代公式为：在MATLAB中编程实现的Gauss-Seidel迭代法函数为：gauseidel。功能：用Gauss-Seidel迭代求线性方程组ax=b的解。调用格式：x,n=gauseidel(A,b,x0,eps,M).其中，A为线性方程组的系数矩阵；b为线性方程组中的常数向量；x0为迭代初始向量；eps为解的精度控制（此参数可选）；M为迭代步数控制（此参数可选）；x为线性方程组的解；n为求出所需精度的解实际迭代步数。(3)Gauss-Seidel迭代法的MATLAB程序代码如下:function x,n=gauseidel(A,b,x0,eps,M)%采用Gauss-Seidel迭代法求线性方程组Ax=b的解%线性方程组的系数矩阵：A%线性方程组的常数向量：b%迭代初始向量：x0%解的精度控制：eps%迭代步数控制：M%线性方程组的解：x%求出所需精度的解实际的迭代步数：nif nargin =3 eps=1.0e-6; %eps表示迭代精度 M=200; %M表示迭代步数的限制值elseif nargin =4 M=200;elseif nargin eps x0=x; x=G*x0+f; n=n+1; %n为最终求出解时的迭代步数 if n=M disp(Warning:迭代次数太多，可能不收敛！); endend(4)进行实例分析：A=1.4449 0.7948 0.8801;0.6946 1.9568 0.1730;0.6213 0.5226 1.9797; b=1 0 1; x0=zeros(3,1); x,n=gauseidel(A,b,x0)x = %输出结果 0.5929 -0.2444 0.3836n = 11 %输出迭代次数n(5)运行图以即流程图：运行图：流程图：开始最大步数M为200读取数据否nargin=3？是eps=10-6 ，最大步数M=200n=1x=G*x0+f否x0=x;n=n+1norm(x-x0)eps？是是n=200？否输出结果Warning:迭代次数太多，可能不收敛结束例题流程图：输入系数矩阵A输入常数向量b和初始向量x0x,n=gauseidel(A,b,x0)输出计算结果输出迭代次数一、分析电路：（RC充电电路）（1） 当t0时，ic（0+）=-Uc（0+）/(R2*R3/（R2+R3）)，达到稳态后，电容中将无电流icf=0A，电流源的全部电流将在两个电阻之间分配，保证端电压相同(也就是电容上的终电压)。即：UR3= Ucf=12V。（3） 由一阶响应电路可用三要素法得到电压公式：，时间常数。 ，且Req=12*6/18=4，=4(s)。同时，以相同的方法得到电流公式：，=4(s)。（4） 可以得到电阻R2的表达式：。二、源程序设计： r1=3; r2=12; r3=6; us=18; is=3; c=1; uc0=-12; %电容C的初始电压值 ucf=12; %电容C的最终稳态电压值 icf=0; %电容C的最终稳态电电流值 T=r2*r3/(r2+r3); %电路时间常数 t=0:0.01:25; %时间变量 uct=ucf+(uc0-ucf)*exp(-t./T); %计算电容C的响应电压 ic0=-uc0/(r2*r3/(r2+r3); %电容C的初始电流值 ict=icf+(ic0-icf)*exp(-t./T); %计算电容C的响应电流 ir20=3;%电阻R2的初始稳态电流 ir2t=(is-ict)*r3/(r2+r3); %计算电阻R2的响应电流 subplot(3,1,1); plot(t,uct,-r); title(Uc(t); subplot(3,1,2); plot(t,ict,-g); title(Ic(t); subplot(3,1,3); plot(t,ir2t,-b); title(Ir2(t);%输出Uc（t）和IR2（t）的响应图像三、波形图：四、流程图：开始读取已知数据计算稳态电流电压数据运用响应电路三要素构造函数uct=ucf+(uc0-ucf)*exp(-t./T)、ict=icf+(ic0-icf)*exp(-t./T)和ir2t=(is-ict)*r3/(r2+r3)结束输出Uc（t）、Ic(t)和IR2（t）的响应图像1.分别用五点和九点等距插值：Lagrange插值多项式原理：n次插值就是利用n+1个插值节点构造的n次的插值多项式，根据差值函数条件，可以得到：设有n+1个互异节点，且，i为整数，构造Ln(x),使Ln(xj)=yj，j为整数。定义，若n次多项式lj(x)在n+1个节点上满足条件，为1时，j=k;为0时，j=k;则称这n+1个n次多项式l0(x),l1(x),ln(x)为节点x0,x1,xn上的n次插值基函数。即为：，k=(0,1,2,n);故满足插值条件多项式为：，称其为Lagrange插值多项式。 n-1阶拉格朗日插值matlab源程序：function y=lagrange(x0,y0,x)% 计算n次lagrange多项式插值% x0,y0输入插值节点向量% x以向量形式输入的插值点% y输出的插值点函数数值n=length(x0);m=length(x);for i=1:m z=x(i); s=0.0; for k=1:n p=1.0; for j=1:n if j=k p=p*(z-x0(j)/(x0(k)-x0(j); end end s=p*y0(k)+s; end y(i)=s;end流程图：输入：z输入：x0, y0开始s=0k=1p=1j=1p=p*(z-x0(j)/(x0(k)-x0(j)(j=k)j=j+1j=n?s=p*y0(k)+sk=n?=n?输出：z, s结束k=k+1是是否否（1）五点等距插值的四阶Lagrange多项式源程序： d=2*pi/4; x=0:d:2*pi; y=sin(x); x0=0:0.1:2*pi; y0=lagrange(x,y,x0); y1=sin(x0); plot(x0,y1,-r,x0,y0,-b,x,y,o)legend(sin(x),四阶插值多项式,lagrange插值点)运行图：流程图：输入：x,y=sin(x)得到插值点用Lagrange进行四阶五点差值计算输入：等距插值点步长d=2*pi/4开始绘图比较sin(x)与Lagrange插值图像比较两图误差结束（2）九点等距插值的八阶Lagrange多项式源程序： d=2*pi/8; x=0:d:2*pi; y=sin(x); x0=0:0.1:2*pi; y0=lagrange(x,y,x0); y1=sin(x0); plot(x0,y1,-r,x0,y0,-b,x,y,+)legend(sin(x),八阶插值多项式,lagrange插值点)运行图：流程图：结束绘图比较sin(x)与Lagrange插值图像比较两图误差输入：等距插值点步长d=2*pi/8输入：x,y=sin(x)得到插值点用Lagrange进行八阶九点差值计算开始2.四阶Lagrange插值多项式近似：（1）源程序：d=5/4;x=0:d:5;y=(1+x)./(1+2*x+3*x.2);x0=0:0.05:5;y0=lagrange(x,y,x0);y1=(1+x0)./(1+2*x0+3*x0.2);x1=0:0.2:5;y2=(1+x1)./(1+2*x1+3*x1.2);y3=lagrange(x,y,x1);e=y2-y3;plot(x,y,o,x0,y1,-r,x0,y0,-b,x1,e,+)legend(插值点,原函数图象,四阶Lagrange插值多项式,误差分布图)（2）插值点：x = 0 1.2500 2.5000 3.7500 5.0000y = 1.0000 0.2748 0.1414 0.0937 0.0698（3）Lagrange差值多项式：Ln(x)=l0(x)*y0+ l1(x)*y1+ l2(x)*y2+ l3(x)*y3+ l4(x)*y4（4）间隔0.2计算的误差值：e = Columns 1 through 9 0 -0.0209 -0.0405 -0.0406 -0.0294 -0.0152 -0.0026 0.0063 0.0111 Columns 10 through 18 0.0122 0.0104 0.0067 0.0023 -0.0021 -0.0057 -0.0078 -0.0082 -0.0066 Columns 19 through 26 -0.0034 0.0012 0.0065 0.0115 0.0150 0.0155 0.0112 0（5）函