C++ 线性方程组的数值解法.docVIP

线性方程组的数值解法 一、下三角方程组的解法 对于线性方程组： 即： ∴ 写成通式的形式： C++源程序： #include iostream.h const int n=4; void LYB(double L[n][n+1]) { for(int i=0;in;i++) { for(int j=0;j=i-1;j++) L[i][n]-=L[i][j]*L[j][n]; } } void main() { double L[n][n+1]={{1,0,0,0,4},{1,1,0,0,7}, {1,-1,1,0,3},{1,-1,-1,1,0}}; LYB(L); for(int i=0;in;i++) coutY[i]=L[i][n]endl; } 二、上三角方程组的解法 对于线性方程组： 即： ∴ 写成通式的形式： #include iostream.h const int n=4; void UXY(double U[n][n+1]) { for(int i=n-1;i=0;i--) { for(int j=i+1;j=n-1;j++) U[i][n]-=U[i][j]*U[j][n]; U[i][n]/=U[i][i]; } } void main() { double U[n][n+1]={{1,-1,-1,1,0},{0,1,-1,1,3}, {0,0,1,1,3},{0,0,0,1,2}}; UXY(U); for(int i=0;in;i++) coutX[i]=U[i][n]endl; } 三、GAUSS消去法解线性方程组 1、消去目标 2、消去过程（假定在以下消去演算中不会出现分母为0的情况） 为直观起见，我们把方程组写成如下形式： 第0步，实现第0列目标：消去第（1）、（2）、……、（n-1）中的项，（0）式保持不变。 方法： 即： 注意： j=0可以不算，如果要算，其结果肯定是0。 经第一步计算后，方程组变成了如下形式： 第1步：在第0步成果的基础上，实现第1列目标：消去第（2）、（3）、……、（n-1）中的项，（0）、（1）式保持不变。 方法： 即： 经第二步计算后，方程组变成了如下形式： 由第一步、第二步的消去过程，不难归纳出任意一步（第k步）的计算通式： Gauss消去法已经将原方程组变形为上三角方程组，逆序回代容易求出最后结果。 在第k步Gauss消去法中，要用到做分母，这就严格要求的绝对值不能太小。对于事先无法估计的值，解决办法是： 每步消去之前，先考察待消去的几个方程的绝对值大小，并交换方程次序，以保证绝对值最大。 程序设计思路： 设计数组，其中前面的0至n-1列存放方程组的系数矩阵,第n列存放常数项。 设计求列中最大值函数，为选择主元作准备。 设计一个交换数据函数，为交换方程做准备， 设计一个消去函数，其中要调用求列中最大值函数确定主元，要调用交换数据函数来交换方程 调用上三角方程组的求解函数求出最终解。 主函数准备数据、调用有关函数、输出结果 课堂练习：编程实现用Gauss消去法解任意大小的方程组。 #include iostream.h #include math.h const int n=4; /*下面是求列中最大值函数，考察范围是：a[k][k]、a[k+1][k]、……、a[n-1][k],返回结果是主元的行号，为选择主元作准备。*/ int Max(double a[n][n+1],int k) { int t=k; double max=fabs(a[k][k]); for(int i=k+1;i=n-1;i++) { if(fabs(a[i][k])max) { t=i; max=fabs(a[i][k]); } } return t; } /*下面是一个交换第k行和第t行方程的函数，选主元的后续工作*/ void change(double a[n][n+1],int k,int t) { double x; for(int j=k;j=n;j++) { x=a[k][j]; a[k][j]=a[t][j]; a[t][j]=x; } } /*下面是高斯消去法的核心函数，其中调用求主元、交换方程函数*/ void gauss(double a[n][n+1]) { for(int k=0;kn-1;k++) { int t=Max(a,k); if(t!=k) change(a,k,t); for(int i=k+1;in;i++) for(int j=k+1;j=n;j++)