文章目录
- Jacobi迭代方法介绍
- Gauss-Seidel迭代方法介绍
- 具体代码实现
- 示例题目
- 实现效果
Jacobi迭代方法介绍
Jacobi迭代法是一种简单的迭代求解方法,适用于严格对角占优矩阵。其基本思想是利用当前迭代步的已知解来更新下一个迭代步的解。在C语言实现中,我们首先需要定义系数矩阵A、常数向量b以及迭代向量x。然后,通过循环迭代,不断更新x的值,直到满足收敛条件或达到最大迭代次数。
Jacobi迭代法的优点是简单直观,但收敛速度相对较慢。特别是对于非严格对角占优的矩阵,Jacobi迭代法可能不收敛。因此,在实际应用中,我们需要根据问题的具体情况选择合适的迭代方法。
Gauss-Seidel迭代方法介绍
Gauss-Seidel迭代法是对Jacobi迭代法的一种改进。在Gauss-Seidel迭代法中,我们使用当前迭代步已更新的解来更新下一个解。这种策略可以加速收敛速度,特别是在对角线元素较大的情况下。在C语言实现中,我们需要对系数矩阵A进行适当的分解,以便在迭代过程中方便地获取已更新的解。
与Jacobi迭代法相比,Gauss-Seidel迭代法具有更快的收敛速度。但是,Gauss-Seidel迭代法并不总是收敛的,其收敛性取决于系数矩阵A的性质。因此,在使用Gauss-Seidel迭代法时,我们需要确保系数矩阵A满足一定的条件(如严格对角占优)。
具体代码实现
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <string.h>
#include "time.h"#define N 3 //要求解的数组大小N x N
#define ERROR 1e-6 //误差要求
#define IterMaxNum 50000 //最大迭代次数float current_error(float* last, float* now)
{float max_error = 0;float temp;for (int i = 0; i < N; i++){temp = *(now + i) - *(last + i);if (temp < 0)temp = -temp;if (temp > max_error)max_error = temp;}return max_error;
}//雅可比迭代函数
float* Jacobi(float* A, float* b)
{clock_t start, end;double cpu_time_used;// 记录开始时间start = clock();int i, j;int iter_num = 0; //目前的迭代次数float error_now = 100; //目前迭代下的误差//生成数组用于保存每次迭代得到的x结果并保存上一次迭代结果至x0float* x = (float*)malloc(N * sizeof(float));memset(x, 0, N * sizeof(float));float* x0 = (float*)malloc(N * sizeof(float));memset(x0, 0, N * sizeof(float));while(error_now > ERROR){//将上一次迭代的结果保存到x0中去memcpy(x0, x, N * sizeof(float));for (i = 0; i < N; i++){float sum = 0;for (j = 0; j < N; j++){//printf("%f\t", A[i * N + j]);if (i != j)sum += A[i * N + j] * (x0[j]);}x[i] = (b[i] - sum) / A[i * N + i];}error_now = current_error(x, x0);iter_num++;printf("当前迭代次数为:【%d】,迭代误差为:【%f】\n", iter_num, error_now);printf("迭代结果为:%f,%f,%f\n", x[0], x[1], x[2]);if (iter_num > IterMaxNum){printf("迭代次数超过最大值要求,返回最后一次迭代结果");break;}}// 记录结束时间end = clock();// 计算经过的CPU时间cpu_time_used = ((double) (end - start)) / CLOCKS_PER_SEC;// 打印Jacobi迭代的运行时间printf("\n【Jacobi迭代运行了 %f 秒。】\n", cpu_time_used);free(x0);return x;
}float* GS(float* A, float* b)
{int i, j;int iter_num = 0; //目前的迭代次数float error_now = 100; //目前迭代下的误差//生成数组用于保存每次迭代得到的x结果并保存上一次迭代结果至x0float* x = (float*)malloc(N * sizeof(float));memset(x, 0, N * sizeof(float));float* x0 = (float*)malloc(N * sizeof(float));memset(x0, 0, N * sizeof(float));clock_t start, end;double cpu_time_used;// 记录开始时间start = clock();while (error_now > ERROR){memcpy(x0, x, N * sizeof(float));for (i = 0; i < N; i++){float sum = 0;for (j = 0; j < N; j++){//printf("%f\t", A[i * N + j]);if (i != j)sum += A[i * N + j] * (x[j]);}x[i] = (b[i] - sum) / A[i * N + i];}error_now = current_error(x, x0); //获得当前迭代的误差结果iter_num++; //迭代次数加1printf("当前迭代次数为:【%d】,迭代误差为:【%f】\n", iter_num, error_now);printf("迭代结果为:%f,%f,%f\n", x[0], x[1], x[2]);if (iter_num > IterMaxNum){printf("迭代次数超过最大值要求,返回最后一次迭代结果");break;}}// 记录结束时间end = clock();// 计算经过的CPU时间cpu_time_used = ((double) (end - start)) / CLOCKS_PER_SEC;// 打印Gauss-Seidel迭代的运行时间printf("\n【Gauss-Seidel迭代运行了 %f 秒。】\n", cpu_time_used);free(x0);return x;
}int main()
{// 例子float A[N][N] = {{10,-1,-2},{-1,10,-2},{-1, -1, 5},}; //要求解的矩阵float b[N] = {72, 83, 42};float* x;// x = Jacobi((float*)A, b);x = GS((float*)A, b);printf("\n【最终结果为】:%f, %f, %f\n", x[0], x[1], x[2]);free(x); //释放x
}
示例题目
以如下方程组为例,进行迭代方法的测试
实现效果
Jacobi迭代实现效果,迭代了【17】
达到了收敛
Gauss-Seidel迭代实现效果,而Gauss-Seidel只迭代了【10】
次即可达到收敛。
注意,有些示例中,Jacobi迭代是比Gauss-Seidel收敛的快
的,但更多的时候Gauss-Seidel是收敛更快
的,并且存在有时候Jacobi收敛
,而Gauss-Seidel却发散
,要根据情况而定。