Julia实现克莱姆法则求解线性方程组

      在实际应用中，有时候我们需要求解一组方程。一般来说，基于线性方程组的解空间理论，线性方程组有唯一解当且仅当有效方程数等于未知数的个数。这时，可以运用多种方法来求出唯一的解。而克莱姆法则（Cramer's Rule）就是一种求解线性方程组的方法。

1 克莱姆法则

      针对一个方程组，首先把它的系数改写成一个行列式，并判断它的值是否为0，如果为0，则不适用于此规则，否则则说明这个线性方程线有解，可以使用克莱姆法则法则。下面给出一个三元方程组的示例，用于描述克莱姆法则的解题过程:

此图来自网站： http://www.saddlebackmath.com/home/cramer-s-rule

以上述方程组为例，首先需要将方程组按照变量进行对齐，即x的放于前，y放中间，z放后。然后将变量的系数改写成一个矩阵A，同时将方程组的值写成一个矩阵B。首先，求出系数的行列式值，用D= det(A)进行求解。其次求x的值时，D1行列式由系数矩阵A的第一列替换为矩阵B的第一列，然后求行列式值，即Dx= det(M1)，即x = Dx / D即可。同理，求y的值，首先需要将系数矩阵A的第二列替换为矩阵B的第一列，然后求行列式值Dy。

2 克莱姆法则Julia实现

      根据百度百科，Julia语言是一个面向科学计算的高性能动态高级程序设计语言，其语法与其他科学计算语言相似。在许多情况下拥有能与编译型语言相媲美的性能。Julia 是个灵活的动态语言，适合科学和数值计算。关于线程代数的知识，可以参考官网: https://docs.julialang.org/en/v1/stdlib/LinearAlgebra

     Julia语言在矩阵的操作上，非常方便，下面给出上述示例的具体代码实现：

julia> using LinearAlgebra

julia> A = [1 -1 2 ; 1 2 3 ; 2 1 1]
3×3 Array{Int64,2}:
 1  -1  2
 1   2  3
 2   1  1

julia> D = det(A)
-12.0

julia> B = [ -3 ; 4 ; -3 ]
3-element Array{Int64,1}:
 -3
  4
 -3

julia> Mx = deepcopy(A)
3×3 Array{Int64,2}:
 1  -1  2
 1   2  3
 2   1  1

julia> Mx[:,1] = B
3-element Array{Int64,1}:
 -3
  4
 -3

julia> Mx
3×3 Array{Int64,2}:
 -3  -1  2
  4   2  3
 -3   1  1

julia> Dx = det(Mx)
36.0

julia> x = Dx / D
-3.0

julia> My = deepcopy(A)
3×3 Array{Int64,2}:
 1  -1  2
 1   2  3
 2   1  1

julia> My[:,2] = B
3-element Array{Int64,1}:
 -3
  4
 -3

julia> My
3×3 Array{Int64,2}:
 1  -3  2
 1   4  3
 2  -3  1

julia> y = det(My) / D
2.0

julia> Mz = deepcopy(A)
3×3 Array{Int64,2}:
 1  -1  2
 1   2  3
 2   1  1

julia> Mz[:,3] = B
3-element Array{Int64,1}:
 -3
  4
 -3

julia> Mz
3×3 Array{Int64,2}:
 1  -1  -3
 1   2   4
 2   1  -3

julia> z = det(Mz) / D
1.0

julia>  

 由此可见，求出的方程解与示例一致。其实，利用Julia语言，我们可以直接进行方程组求解，示例如下:

julia> using LinearAlgebra

julia> A = [1 -1 2 ; 1 2 3 ; 2 1 1]
3×3 Array{Int64,2}:
 1  -1  2
 1   2  3
 2   1  1

julia> B = [ -3 ; 4 ; -3 ]
3-element Array{Int64,1}:
 -3
  4
 -3

julia> A \ B
3-element Array{Float64,1}:
 -3.0
  2.0
  1.0

julia>   

最后，需要说的就是，还有很多其他的矩阵操作，比如转置，求逆矩阵，矩阵乘法等等。示例如下:

julia> using LinearAlgebra

julia> A = [1 -1 2 ; 1 2 3 ; 2 1 1]
3×3 Array{Int64,2}:
 1  -1  2
 1   2  3
 2   1  1

julia> A'
3×3 Adjoint{Int64,Array{Int64,2}}:
  1  1  2
 -1  2  1
  2  3  1

julia> inv(A)
3×3 Array{Float64,2}:
  0.0833333  -0.25   0.583333
 -0.416667    0.25   0.0833333
  0.25        0.25  -0.25

julia> inv(A) * A
3×3 Array{Float64,2}:
 1.0  1.11022e-16  0.0
 0.0  1.0          5.55112e-17
 0.0  0.0          1.0

julia> A' * A
3×3 Array{Int64,2}:
 6  3   7
 3  6   5
 7  5  14

#不是 A + 1
julia> A .+ 1
3×3 Array{Int64,2}:
 2  0  3
 2  3  4
 3  2  2

julia> dot(A,A)
26

julia> ones(3,2)
3×2 Array{Float64,2}:
 1.0  1.0
 1.0  1.0
 1.0  1.0

julia> zeros(3,2)
3×2 Array{Float64,2}:
 0.0  0.0
 0.0  0.0
 0.0  0.0

julia> x = rand(5)
5-element Array{Float64,1}:
 0.09025788794245759
 0.6031193596136584
 0.6252782877091618
 0.43105107636432316
 0.3206800449549978

julia> [2*x[i-1] + 0.5*x[i] + 0.25*x[i+1] for i=2:length(x)-1]
3-element Array{Float64,1}:
 0.6383950276190349
 1.6266406321729785
 1.5462521248392345

julia>