【matlab】解线性方程组的迭代方法

（一）矩阵特征多项式、特征值、特征向量，稀疏矩阵

1、测试函数eig,poly,poly2str

>> format short g

>> A=[1 2 3;4 5 6;7 8 9]

A =

     1     2     3

     4     5     6

     7     8     9

>> p=poly(A)

p =

            1          -15          -18  -1.7567e-14

>> poly2str(p,'x')

ans =

    '   x^3 - 15 x^2 - 18 x - 1.7567e-14'

>> d=eig(A)

d =

       16.117

      -1.1168

  -9.7592e-16

1. 测试函数sparse, spdiags

>> A=[0 0 5 0;3 0 3 0;0 0 0 1;0 4 3 0]

A =

     0     0     5     0

     3     0     3     0

     0     0     0     1

     0     4     3     0

>> S=sparse(A)

S =

   (2,1)        3

   (4,2)        4

   (1,3)        5

   (2,3)        3

   (4,3)        3

   (3,4)        1

>> S=sparse([1 2 2],[1 3 2],[6 7 8],5,5)

S =

   (1,1)        6

   (2,2)        8

   (2,3)        7

>> full(S)

ans =

     6     0     0     0     0

     0     8     7     0     0

     0     0     0     0     0

     0     0     0     0     0

     0     0     0     0     0

>> s=sparse(3,4)

s =

   全零稀疏矩阵: 3×4

>> full(s)

ans =

     0     0     0     0

     0     0     0     0

     0     0     0     0

>> B=[1 2 3;4 5 6;7 8 9]

B =

     1     2     3

     4     5     6

     7     8     9

>> d=[-2,1,3] ;

>> S=spdiags(B,d,3,6)

S =

   (3,1)        7

   (1,2)        2

   (2,3)        5

   (1,4)        3

   (3,4)        8

   (2,5)        6

   (3,6)        9

>> full(S)

ans =

     0     2     0     3     0     0

     0     0     5     0     6     0

     7     0     0     8     0     9

（二）专题实验（迭代法求解线性方程组）

1、编写Jacobi迭代格式解线性方程组的函数；

function X_solve=myJacobi(A,b,X0,mytol)

%A为系数矩阵.

%b为右端常向量.

%X0为初始层，默认取0.

%mytol表示允许误差, 要求默认值为1e-6.

%nargin表示输入参数个数

X=diag(A);D=diag(X);L=-tril(A,-1);U=-triu(A,1);X1=ones(size(b));n=1;

if nargin==3   %如果参数=3,

    mytol=1e-6;

elseif nargin==2   %如果参数=2

    X0=zeros(size(b));

    mytol=1e-6;

else    %如果参数<2

    error('输入参数不足')

end

[r,c]=size(A);

if r==c

    if abs(det(A))>=mytol

        if abs(det(D))>=mytol

            B=D^(-1)*(L+U);

            d=max(abs(eig(B))); %特征值绝对值的最大值

            if d<1  %迭代法收敛

                disp('迭代法收敛')

                while max(abs(X0-X1))>=mytol

                    X1=X0;

                    X0=B*X0+D^(-1)*b;

                    n=n+1;

                end

            else

                disp('迭代法不收敛')

            end

        else

            error('主对角线系数为0')

        end

    else

        error('系数矩阵不可逆')

    end

else

    error('系数矩阵不是方阵')  

end

disp(X0)

disp(n)

end

2、编写Gauss-Seidel迭代格式求解线性方程组的函数；

function X_solve=myGS(A,b,X0,mytol)

%A为系数矩阵.

%b为右端常向量.

%X0为初始层，默认取0.

%mytol表示允许误差, 要求默认值为1e-6.

%nargin表示输入参数个数

X=diag(A);D=diag(X);L=-tril(A,-1);U=-triu(A,1);X1=ones(size(b));n=1;

if nargin==3   %如果参数=3,

    mytol=1e-6;

elseif nargin==2   %如果参数=2

    X0=zeros(size(b));

    mytol=1e-6;

else    %如果参数<2

    error('输入参数不足')

end

[r,c]=size(A);

if r==c

    if abs(det(A))>=mytol

        if abs(det(D))>=mytol

            BG=(D-L)^(-1)*U;

            d=max(abs(eig(BG))); %特征值绝对值的最大值

            if d<1  %迭代法收敛

                disp('迭代法收敛')

                fG=(D-L)^(-1)*b;

                while max(abs(X0-X1))>=mytol

                    X1=X0;

                    X0=BG*X0+fG;

                    n=n+1;

                end

            else

                disp('迭代法不收敛')

            end

        else

            error('主对角线系数为0')

        end

    else

        error('系数矩阵不可逆')

    end

else

    error('系数矩阵不是方阵')  

end

disp(X0)

disp(n)

end

3、要求用矩阵形式编写；

4、对下面四个例题进行实验

>> A=[10 -1 -2;-1 10 -2;-1 -1 5];

>> b=[7.2;8.3;4.2];

>> myJacobi(A,b)

迭代法收敛

    1.1000

    1.2000

    1.3000

15

>> A=[10 -1 -2;-1 10 -2;-1 -1 0.5];

>> b=[7.2;8.3;4.2];

>> myJacobi(A,b)

迭代法收敛

   24.5000

   24.6000

   106.6000

   282

>> A=[1 -9 -10;-9 1 5;8 7 1];

>> b=[1;0;4];

>> myGS(A,b)

迭代法不收敛

     0

     0

     0

     1

>> A=[5 -1 -3;-1 2 4;-3 4 15];

>> b=[-1;0;4];

>> myGS(A,b)

迭代法收敛

   -0.0984

   -1.1639

    0.5574

24

程序编写要求：

1. 函数输入形参要求如下，也可以修改得更好、更适合你自己.

function Xsolve=myJacobi(A,b,X0,mytol)

%A为系数矩阵.

%b为右端常向量.

%X0为初始层，默认取0.

%mytol表示允许误差, 要求默认值为1e-6.

%Xsolve表示收敛情况下满足精度的数值解.

%nargin表示输入参数个数，下面这一段是简单的判断，你可以编写的更好.

if nargin<4    %如果参数<4,

mytol=1e-6; 

if nargin<3   %如果参数<3

X0=zeros(size(b));

if nargin<2

                    error(‘输入参数不足.’);

                end

            end

end

2、第一步要判断系数矩阵是否是方阵，如果不是，提示错误；（可以使用size函数）

3、第二步要判断系数矩阵是否可逆，用“abs(det(A))<mytol”或者别的方式判断，不满足要提示；

4、第三步要判断系数矩阵主对角元素是否为零，也是用“绝对值<mytol”的方式判断；

5、第四步用迭代法基本定理（计算谱半径）判断迭代法是否收敛，要在命令窗口输出是否收敛的结论；

6、第五步在收敛的基础上进行迭代格式的程序编写；

7、迭代格式终止循环的判断方式用后一层和前一层计算结果的差的绝对值的最大值小于某个给定的小量mytol，比如1e-6，即，max(abs(A-B))<mytol. 也可以用其它方式终止循环；

8、收敛的情况下，和精确解做简单比较，精确解可以采用Ab的形式；

9、注意：命令窗口中显示的结果不是真实结果，只是在输入格式的条件下进行四舍五入的显示，要想看到小数点后多个位的值，可以采用format long 等形式。

最后

以上就是称心项链为你收集整理的【matlab】解线性方程组的迭代方法的全部内容，希望文章能够帮你解决【matlab】解线性方程组的迭代方法所遇到的程序开发问题。

如果觉得靠谱客网站的内容还不错，欢迎将靠谱客网站推荐给程序员好友。