【matlab基础知识代码】（十五）线性方程求解（三种情况）

>>  A=[1 2 3 4; 4 3 2 1; 1 3 2 4; 4 1 3 2]; B=[5 1; 4 2; 3 3; 2 4];
 x=inv(A)*B,norm(A*x-B)

x =

   -1.8000    2.4000
    1.8667   -1.2667
    3.8667   -3.2667
   -2.1333    2.7333


ans =

   7.2911e-15

>> y=inv(sym(A))*B,norm(A*y-B)
 
y =
 
[  -9/5,   12/5]
[ 28/15, -19/15]
[ 58/15, -49/15]
[-32/15,  41/15]
 
 
ans =
 
0

>> A=[1,4,0,-1,0,7,-9;2,8,-1,3,9,-13,7;0,0,2,-3,-4,12,-8;-1,-4,2,4,8,-31,37];
 B=[3; 9; 1; 4]; C=[A B]; 
rank(A), rank(C)

ans =

     3


ans =

     3

 通解与检验

>>  Z=null(sym(A)), x0=sym(pinv(A)*B), 
syms a1 a2 a3 a4; 
x=Z*[a1;a2;a3;a4]+x0, E=A*x-B
 
Z =
 
[-4, -2, -1,  3]
[ 1,  0,  0,  0]
[ 0, -1,  3, -5]
[ 0, -2,  6, -6]
[ 0,  1,  0,  0]
[ 0,  0,  1,  0]
[ 0,  0,  0,  1]
 
 
x0 =
 
 92/395
368/395
459/790
 -24/79
347/790
247/790
303/790
 
 
x =
 
3*a4 - 2*a2 - a3 - 4*a1 + 92/395
                    a1 + 368/395
      3*a3 - a2 - 5*a4 + 459/790
      6*a3 - 2*a2 - 6*a4 - 24/79
                    a2 + 347/790
                    a3 + 247/790
                    a4 + 303/790
 
 
E =
 
0
0
0
0

`pinv()` 是计算矩阵的伪逆。`pinv(A)` 会计算矩阵A的伪逆
`syms` 用来声明符号变量。在这里，`a1`、`a2`、`a3` 和 `a4` 被声明为符号变量
`Z=null(sym(A))`：这行代码计算矩阵A的零空间，将其转换为符号表示，并将其赋值给Z。零空间是指所有乘以该矩阵后得到零向量的向量集合。(求基础解系)
`x0=sym(pinv(A)*B)`：这行代码计算矩阵A的伪逆，将其乘以矩阵B，将结果转换为符号表示，并将其赋值给x0。本质上，它解决了方程Ax=B，找到了x的解。（特解）
x=Z*[a1;a2;a3;a4]+x0：这行代码为向量x定义了一个符号表达式。它将Z中的零空间向量与符号变量`[a1;a2;a3;a4]`相乘，然后加上x0的结果。所以，x是矩阵A的零空间向量和方程Ax=B的特解的线性组合。（方程通解）

基本行变换方法

>> C=[A B]; D=rref(C)

D =

     1     4     0     0     2     1    -3     4
     0     0     1     0     1    -3     5     2
     0     0     0     1     2    -6     6     1
     0     0     0     0     0     0     0     0

>> A=[4,7,1,4; 3,7,4,6]; B=[3; 4]; C=[A B]; rank(A), rank(C)

ans =

     2


ans =

     2
>> syms a1 a2 b1 b2; x1=null(sym(A))*[a1; a2]+sym(A\B), A*x1-B,
 
x1 =
 
                 3*a1 + 2*a2
8/21 - (12*a2)/7 - (13*a1)/7
                    a1 + 1/3
                          a2
 
 
ans =
 
0
0
>> a=rref(sym([A B])); x2=[a(:,3:5)*[-b1; -b2; 1]; b1; b2], A*x2-B
 
x2 =
 
          3*b1 + 2*b2 - 1
1 - (12*b2)/7 - (13*b1)/7
                       b1
                       b2
 
 
ans =
 
0
0

>>  A=[1 2 3 4; 2 2 1 1; 2 4 6 8; 4 4 2 2]; B=[1:4]'; C=[A B];rank(A), rank(C)

ans =

     2


ans =

     3

>>  x=pinv(A)*B, norm(A*x-B)


x =

    0.5466
    0.4550
    0.0443
   -0.0473


ans =

    0.4472