计算方法-方程求根
目录
一、 引言…………………………………… 3
二、 算法描述……………………………… 3
三、 试验及分析…………………………… 4
四、 自己的特色…………………………… 6
五、 总结…………………………………… 9
六、 程序清单……………………………… 9
PPT：https://wenku.baidu.com/view/31888818a48da0116c175f0e7cd184254a351b38
一、引言
在我们的实际生活和生成当中会出现大量的方程组，在这些方程组需要求解，如果是要用人工的去求解的时候无疑是一个庞大的工作，而且很能避免错误的出现，这样会浪费大量的人力物力和时间。于是我们今天在前辈们关于方程组的求解总结出来的高斯消去法和直接三角分解法（LU分解）的基础之上借助编程过程中的循环来让计算机对方程组经行求解。这样可以节省大量的时间同时避免了错误的出现。下面将是我们本着这个目的去实现高斯消去法和直接三角分解法（LU分解）的过程。虽然有很多的难题还没有解决，但是这是我们小组共同努力的成果。
二、算法描述
1、高斯列主元消去法的原理
Gauss消去法的基本思想是一次用前面的方程消去后面的未知数，从而将方程组化为等价形式：
这个过程就是消元，然后再回代就好了。具体过程如下：
对于 ，若 依次计算
然后将其回代得到：
以上是高斯消去。
但是高斯消去法在消元的过程中有可能会出现 的情况，这时消元就无法进行了，即使主元数 但是很小时，其做除数，也会导致其他元素数量级的严重增长和舍入误差的扩散。因此，为了减少误差，每次消元选取系数矩阵的某列中绝对值最大的元素作为主元素。然后换行使之变到主元位置上，再进行销元计算。即高斯列主元消去法。
2、直接三角分解法（LU分解）的原理
先将矩阵A直接分解为 则求解方程组的问题就等价于求解两个三角形方程组。
直接利用矩阵乘法，得到矩阵的三角分解计算公式为：
由上面的式子得到矩阵A的LU分解后，求解Ux=y的计算公式为
以上为LU分解法。
三、实验及分析
高斯
3 4
1 2 3 6
2 2 3 7
1 4 4 9
答案：1 1 1
3 4
0 0 1 3
3 0 0 2
0 2 0 1
答案：0.667 0.5 3
LU算法：
中文书本P103例7:
实例: 系数矩阵0 1
1 1 结果矩阵 1 1
四、自己的特色
1用高斯算法解一般的方程组
3 5
10 200 300 4 233
22 22 44 1 22
1 20 30 0.4 233
4 6
100 34 345 334 11 552
12 0.2 0.12 12 23 344
122 53 12 43 23 45
123 234 0.3 34 345 34
3 4
1 22 676 345
12 24 124 466
1 22 676 100
4 5
21 345 67 15 456
0.12 35 02.4 234 3455
122 34 13 345 234
12 12 0.1 24 24
5 6
1 2 3 4 6 7
1 2 5 12 4 56
12 2 5 11 2 122
32 11 2 4 6 23
2 4 8 16 10 63
4 3
1 3 45
32 12 24
42 0.1 23
33 15 60
3 3
1 2 4
1 3 5
2 5 9
4 4
1 4 7 3
2 8 14 6
3 12 21 9
3 3 6 2
2
五、总结
通过这次实习我们学到了很多同样也暴露出了很多的问题：
1、这从的方程组求解的两种方法不像前一次迭代的那样好实现，所以在很多的方面大家写的代码又点不完善的地方，有几个同学本来的基础就是有点弱，所以全部的程序没有写完，有的写完了但是会出现些错误，这样暴露出来在我们编程能力不一样，由于上课大家在一起讲解的时候时间有限，这样我们思考要是还是一个人单独作战的时候是不是会使大家的水平更加的层次不起。所以我们决定进行小组内再分组的方案。有两个人组成一组，形成一个变成能力比较好的人带一个比较差的人，这样大家在课下编程的时候要一起行动，相互的帮助可能比我们在课上汇报讨论的时候更能解决实质问题。
2、由于这次大家都对原理比较清楚，在私下也有点交流，所以在两种方法的主要实现过程中没有太大的问题大大加速了汇报的进度。在最后的时候我们进行的问题的集中解决，并且在最后的时候请别的小组的实现的比较好的同学给我们指导，并且相互交流。明白别的小组对我们小组出现的问题是怎么解决的。然后在最后得出了一些大家都有的一些共有的缺陷，进行的商讨。
3、在即上一次的分工之后，我们这次决定相互的交换一下工作，然后让一些编程能力相对比较弱的同学参与到了讨论后的代码修改当中。我们最终的程序或许没有解决我们最初讨论出的问题，但是我想大家参与了就是一种最大的收获。这次的分工是：由彭亚琼、刘宇、方正、陈俊桦、韩学武五位同学负责修改程序，王艳琴和向胜男负责细化许多上台汇报的工作，张钟文和杨耀鹏写word文档，肖大军和梅旭负责制作PPT。
六、程序清单
源码：
#include<stdio.h> #include<math.h> #define zero 1e-5 //允许的误差范围 double x[10000]; //保存每一次的计算值 double fr(double x) { return x-exp(-x); //测试函数f(x) } double fd(double x) { return exp(-x); //迭代法迭代式 } double fn(double x) //牛顿法迭代公式 { return x-(x-exp(-x))/(1+x); } int erfen(double a,double b) { int i=0,j; double x1,x2,x3,y; printf("二分法:\n"); printf("误差 收敛度\n"); x1=a; x2=b; x3 = (x1+x2)/2; while(fabs(x2-x1)>=zero) { x[i++]=x3; //printf("%lf %lf %lf\n",x1,x2,x3); if(fr(x1)*fr(x3)>0 )x1=x3; else if(fr(x1)*fr(x3)<0)x2=x3; else { printf("方程的根为%.7lf\n",x3); break; } x3 = (x1+x2)/2; } for(j=0;j<i;j++) { printf("%lf",x[j]-x3); if(j==0)printf("\n"); (j)&&printf(" %lf\n",x[j]/x[j-1]); //if(i-j>1) 短路特性 } printf("\n经过%d运算,方程的根为%.7lf\n",i,x3); return 0; } int diedai(double a) { int i=0,j; double x1,x2; printf("\n迭代法:\n"); printf("误差 收敛度\n"); //x1=a+1; //x1的初始化 x2=a; do { x[i++]=x2; x1 = x2; //printf("%.7lf\n",x1); x2 = fd(x1); }while(fabs(x1-x2)>=zero); for(j=0;j<i;j++) { printf("%lf ",x[j]-x2); if(j==0)printf("\n"); (j)&&printf("%lf\n",x[j]/x[j-1]); } printf("\n经过%d运算,方程的根为%.7lf\n",i,x2); return 0; } int Newton(double a) { int i=0,j; double x1,x2; printf("\n牛顿法:\n"); printf("误差 收敛度\n"); //x1=a*2; x2=a; do { x[i++]=x2; x1 = x2; //printf("%.7lf\n",x1); x2 = fn(x1); //printf("%.7lf\n",x2); }while(fabs(x1-x2)>=zero); for(j=0;j<i;j++) { printf("%lf ",x[j]-x2); if(j==0)printf("\n"); (j)&&printf("%lf\n",x[j]/x[j-1]); } printf("\n经过%d运算,方程的根为%.7lf\n",i,x2); return 0; } int main() { double a,b; printf("方程求根:\n"); printf("例式：x = e(-x)\n"); printf("请输入区间【0.4 7】："); scanf("%lf %lf",&a,&b); erfen(a,b); diedai(a); Newton(a); return 0; }
#include<stdio.h> #include<math.h> #define zero 1e-5 //允许的误差范围 #define chuzhi 0.5 double x[10000]; //保存每一次的计算值 double fr(double x) { return x-exp(-x); //测试函数f(x) } double fd(double x) { return exp(-x); //迭代法迭代式 } double fn(double x) //牛顿法迭代公式 { return x-(x-exp(-x))/(1+x); } int erfen(double a,double b) { int i=0,j; double x1,x2,x3,y; x1=a; x2=b; x3 = (x1+x2)/2; while(fabs(x2-x1)>=chuzhi) { x[i++]=x3; //printf("%lf %lf %lf\n",x1,x2,x3); if(fr(x1)*fr(x3)>0 )x1=x3; else if(fr(x1)*fr(x3)<0)x2=x3; else { printf("方程的根为%.7lf\n",x3); break; } x3 = (x1+x2)/2; } for(j=0;j<i;j++) { printf("%lf",x[j]-x3); if(j==0)printf("\n"); (j)&&printf(" %lf\n",x[j]/x[j-1]); //if(j!=0) 短路特性 } printf("\n经过%d运算,提供初始值：%.7lf\n",i,x3); a=x3; Newton(a); return 0; } int Newton(double a) { int i=0,j; double x1,x2; printf("\n牛顿法:\n"); printf("误差 收敛度\n"); x2=a; do { x[i++]=x2; x1 = x2; x2 = fn(x1); }while(fabs(x1-x2)>=zero); for(j=0;j<i;j++) { printf("%lf ",x[j]-x2); if(j==0)printf("\n"); (j)&&printf("%lf\n",x[j]/x[j-1]); } printf("\n经过%d运算,方程的根为%.7lf\n",i,x2); return 0; } int main() { double a,b; printf("——————方程求根————————\n"); printf("************例式：x = e(-x)********\n"); printf("请输入区间,包含区间【0.5 0.6】："); //[0.5 7] scanf("%lf %lf",&a,&b); erfen(a,b); return 0; }
