1.掌握方程求根数值方法的基本方法,算法设计,实验验证和结果分析。 2.选择使用二分法、迭代法、牛顿法、割线法等方法对给定的方程进行根的求解; 3. 进行结果判断和误差分析。 二、实验内容和原理(必填) 熟悉使用二分法、迭代法、牛顿法、割线法等方法对给定的方程进行根的求解。选择上述方法中的两种方法求方程: 在[1,1.5]内的一个实根,且要求满足精度 二分法:设区间端点为 x1 和 x2,计算区间的中点 x0,可将区间[x1, x2]分为[x1, x0]和[x0, x2]两部分。这时,为了确定下一次求根搜索的区间,必须判断方程的根在哪 一个区间内,根据函数图像可知,方程的根所在区间的两个端点处的函数值的符号一定是 相反的,也就是说,如果 f(x0)与 f(x1)是异号的,则根一定在左区间[x1, x0]内,否则 根一定在右区间[x0, x2]内。 迭代法:由函数方程 f(x)=0,构造一个等价方程 x=g(x),从某个近似根 x0 出发, 令 xn+1=g(xn),n=0,1,2,…可得到序列{xn},若{Xn}收敛,即 lim xn = x*,只要 g (x)连续,x*是方程(1)的根,也就是 f(x)=0 的根。此时{Xn}就是方程(1)的一个 近似序列,n 越大,xn 的近似程度就越好。若{Xn}发散,则迭代法失败。 三、主要仪器设备(必填) 笔记本 四、操作方法与实验步骤(可选) 二分法 1. #include<bits/stdc++.h> 2. using namespace std; 3. double eps=1e-6; 4. double f(double x){ 5. return (pow(x,4)+2*pow(x,2)-x-4); 6. } 7. int main(){ 8. double x1=1,x2=1.5; 9. double root=(x1+x2)/2; 10. double y=f(root); 11. int round=1; 12. while(fabs(y)>eps){ 13. if(y>0){ 14. x2=root; 15. } 16. else{ 17. x1=root; 18. } 19. root=(x1+x2)/2; 20. y=f(root); 21. round++; 22. } 23. cout<<"root:"<<root<<"round:"<<round<<endl; 24. system("pause"); 25. return 0; 26. } 27. 迭代法 1. #include<bits/stdc++.h> 2. using namespace std; 3. const double eps=1e-6; 4. double g(double x){ 5. return (pow((x+4-2*x*x),0.25)); 6. } 7. int main(){ 8. double x1,x0=1.25; 9. do{ 10. x1=x0;x0=g(x1); 11. }while (fabs(x0-x1)>eps); 12. cout<<"root:"<<x0; 13. system("pause"); 14. return 0; 15. } 16. 五、实验数据记录和处理(可选) 六、实验结果与分析(必填) 七、讨论、心得(可选) 迭代法要建立适当的迭代关系,如果迭代函数选择不恰当结果可能会出现nan(非数) 与inf。迭代法中如果方程无解,近似根序列不会收敛,迭代过程会出现死循环。在迭代前应先判断方程是否有解,同时对循环次数作出限制。 | |
1.掌握线性方程组的直接求解的基本方法,算法设计,实验验证和结果分析; 2.合理利用Gauss消元法、LU分解法、追赶法选择并求解给定的线性方程组,选做两题; 3. 进行结果判断和误差分析。 二、实验内容和原理(必填) 合理利用Gauss消元法、LU分解法、追赶法求解下列方程组: ① ② ③
高斯消元法: 通过始终消去主对角线下方的元素,使整个方程组变成一个三角形似的格式。之后再将已得出的答案一个个地代入已被简化的等式中的未知数中,就可求出其余的答案了。所以我们可以用初等行变换把增广矩阵转换为行阶梯阵,然后回代求出方程的解。这种通过消元再回代的求解方法称为高斯消元法。 列主元素消元: 在第 k 步消元前,先找出 k 行下所有第 k 列元素最大的非零元素,将第 r 行与第 k 行进行整行交换,这样既不影响原方程的解,也可以将绝对值最大的作为主元,放在除数的位置上,尽可能减小引入误差。 全主元素消元法: 全主元消去法与列主元消去法类似,只不过列主元消去法是从第 k 列中选取一个最大的元素,与第 k 行进行交换。而全主元消去法是从第 k 行第 k 列开始的右下角矩阵中所有元素中选取一个最大的元素作为主元,同时交换 r 行与 c 列,从而保证稳定性。 LU分解法: 当系数矩阵 A 满足顺序主子式不为 0 时,可将 A 分解为为一个单位下三角矩阵 L 和一个上三角矩阵 U 的乘积,且分解唯一,然后方程式变:Ly=b,Ux=y,接着先求 y,再求出 x。对于题目给定的第二个问题要格外注意将系数矩阵的第一行与第二行交换位置,才能得到正确的答案。 三、主要仪器设备(必填) 笔记本 四、操作方法与实验步骤(可选) 1. #include <stdio.h> 2. 3. #define N 3 4. 5. void gaussianElimination(float matrix[N][N+1]) { 6. int i, j, k; 7. float factor, temp; 8. 9. for (i = 0; i < N; i++) { 10. if (matrix[i][i] == 0.0) { 11. printf("Error: Division by zero\n"); 12. return; 13. } 14. 15. for (j = i + 1; j < N; j++) { 16. factor = matrix[j][i] / matrix[i][i]; 17. 18. for (k = i; k < N + 1; k++) { 19. matrix[j][k] -= factor * matrix[i][k]; 20. } 21. } 22. } 23. 24. for (i = N - 1; i >= 0; i--) { 25. temp = matrix[i][N]; 26. 27. for (j = i + 1; j < N; j++) { 28. temp -= matrix[i][j] * matrix[j][N]; 29. } 30. 31. matrix[i][N] = temp / matrix[i][i]; 32. } 33. } 34. 35. int main() { 36. float matrix[N][N+1]; 37. int i, j; 38. 39. printf("Enter the augmented matrix:\n"); 40. for (i = 0; i < N; i++) { 41. for (j = 0; j < N + 1; j++) { 42. scanf("%f", &matrix[i][j]); 43. } 44. } 45. 46. gaussianElimination(matrix); 47. 48. printf("Solution:\n"); 49. for (i = 0; i < N; i++) { 50. printf("x%d = %.2f\n", i+1, matrix[i][N]); 51. } 52. 53. return 0; 54. } 55. 1. #include <iostream> 2. #include <vector> 3. 4. using namespace std; 5. 6. // LU分解函数 7. void LU_decomposition(vector<vector<double>>& A, vector<vector<double>>& L, vector<vector<double>>& U) { 8. int n = A.size(); 9. for (int i = 0; i < n; i++) { 10. // 计算U的第一行 11. for (int j = i; j < n; j++) { 12. double sum = 0; 13. for (int k = 0; k < i; k++) { 14. sum += L[i][k] * U[k][j]; 15. } 16. U[i][j] = A[i][j] - sum; 17. } 18. 19. // 计算L的第一列 20. for (int j = i; j < n; j++) { 21. if (i == j) { 22. L[i][i] = 1; 23. } 24. else { 25. double sum = 0; 26. for (int k = 0; k < i; k++) { 27. sum += L[j][k] * U[k][i]; 28. } 29. L[j][i] = (A[j][i] - sum) / U[i][i]; 30. } 31. } 32. } 33. } 34. 35. // 前向替换函数 36. void forward_substitution(vector<vector<double>>& L, vector<double>& b, vector<double>& y) { 37. int n = L.size(); 38. for (int i = 0; i < n; i++) { 39. double sum = 0; 40. for (int j = 0; j < i; j++) { 41. sum += L[i][j] * y[j]; 42. } 43. y[i] = b[i] - sum; 44. } 45. } 46. 47. // 后向替换函数 48. void backward_substitution(vector<vector<double>>& U, vector<double>& y, vector<double>& x) { 49. int n = U.size(); 50. for (int i = n - 1; i >= 0; i--) { 51. double sum = 0; 52. for (int j = i + 1; j < n; j++) { 53. sum += U[i][j] * x[j]; 54. } 55. x[i] = (y[i] - sum) / U[i][i]; 56. } 57. } 58. 59. int main() { 60. vector<vector<double>> A = { {5.291, -6.130, -1, 2}, 61. {3e-16, 59.14, 3, 1}, 62. {11.2, 9, 5, 2}, 63. {1, 2, 1, 1} }; 64. vector<double> b = { 46.78, 59.17, 1, 2 }; 65. int n = A.size(); 66. 67. // 初始化L和U矩阵 68. vector<vector<double>> L(n, vector<double>(n, 0)); 69. vector<vector<double>> U(n, vector<double>(n, 0)); 70. 71. // 进行LU分解 72. LU_decomposition(A, L, U); 73. 74. // 前向替换 75. vector<double> y(n, 0); 76. forward_substitution(L, b, y); 77. 78. // 后向替换 79. vector<double> x(n, 0); 80. backward_substitution(U, y, x); 81. 82. // 输出结果 83. cout << "解向量x:" << endl; 84. for (int i = 0; i < n; i++) { 85. cout << "x[" << i << "] = " << x[i] << endl; 86. } 87. 88. return 0; 89. } 90. 五、实验数据记录和处理(可选) 六、实验结果与分析(必填) 七、讨论、心得(可选) 高斯消元法中如果绝对值很小,由于第k次运算中在分母位置,因此作除数会引起很大的误差,从而影响算法的稳定性。为了减少计算过程中舍入误差对方程组求解的影响,应选择绝对值尽可能大的主元作除数。 LU消元法解答第二题时,得到的结果与实际值存在较大的误差,在排除程序的逻辑错误后,发现系数矩阵的第一行第一列存在一个较小的数3e-16,为了避免误差,应该将系数矩阵的第一行与第二行交换后再输入主程序中,才会与实际值比较接近。另外,不要忘记修改常数向量的顺序。 |