4 解方程组
4.1普通方程组
解普通方程组可以用函数solve(),solve()的基本用法是solve(A,b),其中,A为方程组的系数矩阵,b为方程组的右端。例如:
已知方程组:
2x1+2x3=1
2x1+x2+2x3=2
2x1+x2=3
解法如下:
> A
[,1] [,2] [,3]
[1,] 2 0 2
[2,] 2 1 2
[3,] 2 1 0
> b=1:3
>b
[1] 1 2 3
> solve(A,b)
[1] 1.0 1.0 -0.5
即x1=1,x2=1,x3=-0.5。
4.2 特殊方程组
对于系数矩阵是上三角矩阵和下三角矩阵的方程组。R中提供了backsolve()和fowardsolve()来解决这个问题。
backsolve(r, x, k=ncol(r), upper.tri=TRUE, transpose=FALSE)
forwardsolve(l, x, k=ncol(l), upper.tri=FALSE, transpose=FALSE)
这两个函数都是符合操作的函数,大致可以分为三个步骤:
①通过将系数矩阵的上三角或者下三角变为0的到新的系数矩阵,这通过upper.tri参数来实现,若upper.tri=TRUR,上三角不为0。
②通过将对步骤1中得到的新系数矩阵进行转置得到新的系数矩阵,这通过transpose参数实现,若transpose=FALSE,则步骤1中得到的系数矩阵将被转置。
③根据步骤2得到的系数矩阵来解方程组。
X1+4X2+7X3=1
2X1+5X2+8X3=2
3X1+6X2+9X3=3
方程组的系数矩阵为:
> A
[,1] [,2] [,3]
[1,] 1 4 7
[2,] 2 5 8
[3,] 3 6 9
> b
[1] 1 2 3
> backsolve(A,b,upper.tri=T,transpose=F)
[1] -0.8000000 -0.1333333 0.3333333
过程分解:
①upper.tri=T,说明系数矩阵的上三角不为0。
> B=A
> B[lower.tri(B)]=0
> B
[,1] [,2] [,3]
[1,] 1 4 7
[2,] 0 5 8
[3,] 0 0 9
②transpose=F说明系数矩阵未被转置。
③解方程:
> solve(B,b)
[1] -0.8000000 -0.1333333 0.3333333