本文详细介绍了高斯-赛德尔、最速下降法和共轭梯度法在解决线性方程组中的应用，包括原理分析和Java代码实现，展示了三种迭代方法的收敛性和迭代结束条件。
一、题目要求
二、原理分析
2.1 高斯-赛德尔方法原理
设方程组Ax = f 有唯一解，将Ax = b 变形为等价的方程组 x = Bx + b，将公式化简为以下形式：
j=1naijxj= bi(i=1,2…n)
（式1）
由此可以建立迭代公式，如式二所示：
xi(k+1)=1aii[bi-j=1i-1aijxj(k+1)-j=i+1naijxj(k)]
(i=1,2…n;k=1,2,……) （式2）
取初试向量x(0)=[x10,x20…xn(0)]T
，进行迭代计算，高斯-赛德尔迭代法比简单迭代法收敛快，相比简单迭代法，高斯-赛德尔跌打法在计算xik+1
时利用了刚刚迭代出的x1k+1，x2k+1，…，xi-1k+1
。
收敛性：已知定理，如果系数矩阵A对称正定，则高斯赛德尔迭代公式收敛。已知系数矩阵A是Hilbert矩阵，是对称正定的，所以高斯赛德尔迭代公式收敛。
2.2 最速下降法原理
已知A =（aii）n*n
为实对称正定矩阵，b，x∈Rn
，则线性方程组Ax = b的解x满足：解x使得二次函数φx= 12Ax,x
-(b,x)取得极小值。所以求解线性方程组Ax = b的解，即为求解二次函数φx= 12Ax,x
-(b,x)的极小值点。求二次函数φx
极小值点的一般方法是：构造一个向量序列{x(k)
},使得φxk→
minφx
。
方法步骤如下：
1）任取一个初始向量x0
；
2）构造迭代格式xk+1=xk+αkpk
，其中pk
是搜索方向，αk
是搜索步长；
3）选择pk
和αk
使得φxk+1
= φ（xk+αkpk
）<φxk
,则当k→∞
时，有φxk→φx*
= minφx
；
4）计算误差ε
，直到||xk-xk-1
||<ε
或者||rk
||<ε
迭代结束。
在以上步骤中，对于搜索步长αk
，当αk
= rk，pkApk，pk
时，αk
是k到k+1的最优步长。对于搜索方向pk取为函数φx
减少最快的方向，即φx
的负梯度方向，此时pk
=rk
。
最速下降法就是使得目标函数下降的最快的方向，反映了目标函数的一种局部性质，最速下降法是线性收敛的算法。
收敛性：最优步长规则下的最速下降法是线性收敛的算法。
2.3 共轭梯度法原理
共轭梯度法同时考量上一步和当前方向，目前的下降方向是上步下降方向和当前残量的线性组合。pk
=rk+βk-1pk-1
。寻求步长α
，使得φx+αp
= min，其中αk
= rk，pkApk，pk
，对于pk
=rk+βk-1pk-1
，选择βk-1使得φx+αp
= min，可得βk-1= pk-1TArkpk-1TApk-1
，综合以上概念可得到共轭梯度发的基本流程。
1）任取一个初始向量x（0） ；
2）构造迭代格式xk+1=xk+αkpk
，其中pk
是搜索方向，αk
是搜索步长；
3）计算步长αk
= rk，rkApk，pk
；
4）rk+1
= rk
- αkApk
；
5）βk= rk+1，rk+1rk，rk
；
6）pk= rk+1
+βkpk
；
7）计算误差ε
，直到||rk
||<ε
迭代结束。
收敛性：最优步长规则下的共轭梯度法是线性收敛的算法。
三、代码实现
代码用java语言实现
public class Main {
public static void main(String[] args) {
/**
* 数值分析作业
*/
double[][] A = hilbert(12);//设A为一个12阶次的hilbert矩阵
double[] b = {1353.0/436,1295.0/594,1121.0/640,934.0/629,2142.0/1651,1162.0/1005
,1645.0/1574,1545.0/1618,901.0/1024,787.0/964,611.0/802,793.0/1110};
double[] x0 = new double[12];
GS(A,b,x0);//高斯赛德尔方法
GD(A,b,x0);//梯度下降法（x-x0）
GDR(A,b,x0);//梯度下降法（r）
CG(A,b,x0);//共轭梯度法
}
static double[][] hilbert(int n){//定义hilbert矩阵的方法
double[][] A = new double[n][n];
for (int i = 0; i < A.length; i++) {
for (int j = 0; j < A[i].length; j++) {
A[i][j] = (1.0/(i+j+1));
}
}
return A;
}
static void GS(double[][] A,double[] b,double[] x0){//高斯-赛德尔迭代法
double[] x =Arrays.copyOf(x0,x0.length);
int number = 0;
while(true){
for (int i = 0; i < A.length; i++) {//计算迭代后的下一次x
double temp = 0;
for (int j = 0; j < A[i].length; j++) {
if (i != j) {
temp += x[j] * A[i][j];
}
}
x[i] = (b[i] - temp) / A[i][i];
}
number++;
double[] err = new double[x.length];
for(int i = 0; i< x.length; i++){//以x-x0的无穷范数为判断标准
err[i] = Math.abs(x[i] - x0[i]);
}
Arrays.sort(err);
if(err[err.length-1] < 0.0001){//当无穷范数小于误差值输出结果
System.out.println(String.format("高斯-赛德尔迭代法：迭代次数k =%d,其X值为：",number)
+Arrays.toString(x));
break;
}
x0 = Arrays.copyOf(x,x.length);
}
}
static void GD(double[][] A,double[] b,double[] x0){//最速下降法方法
double[] x =Arrays.copyOf(x0,x0.length);
int number = 0;
double[] r = new double[x0.length];
while(true){
for(int i = 0; i<A.length;i++){//计算迭代后的r
double temp = 0;
for(int j = 0 ;j < A[i].length;j++){
temp += x[j] * A[i][j];
}
r[i] = b[i] - temp;
}
double element = 0 ;
for(int i = 0;i<r.length;i++){//计算步长式子的分子项
element += r[i]*r[i];
}
double[] r1 = new double[r.length];
for(int i = 0; i<A.length;i++){
double temp = 0;
for(int j = 0 ;j < A[i].length;j++){
temp += r[j] * A[i][j];
}
r1[i] = temp;
}
double denominator = 0;
for(int i = 0;i<r.length;i++) {//计算步长式子的分母项
denominator += r1[i] * r[i];
}
double a = element/denominator;//计算步长
for(int i = 0;i<r.length;i++){
x[i] = x[i]+a*r[i];
}
number++;
double[] err = new double[x.length];
for(int i = 0; i< x.length; i++){//计算x-x0的绝对值为err值，即无穷范数
err[i] = Math.abs(x[i] - x0[i]);
}
Arrays.sort(err);
if(err[err.length-1] < 0.0001){//当x-x0的绝对范数小于误差
System.out.println(String.format("最速下降迭代法（x-x0）：迭代次数k =%d,其X值为：",number)
+Arrays.toString(x));
break;
}
x0 = Arrays.copyOf(x,x.length);
}
}
static void GDR(double[][] A,double[] b,double[] x0){//最速下降法方法
double[] x =Arrays.copyOf(x0,x0.length);
int number = 0;
double[] r = new double[x0.length];
while(true){
for(int i = 0; i<A.length;i++){//计算迭代后的r
double temp = 0;
for(int j = 0 ;j < A[i].length;j++){
temp += x[j] * A[i][j];
}
r[i] = b[i] - temp;
}
double element = 0 ;
for(int i = 0;i<r.length;i++){//计算步长式子的分子项
element += r[i]*r[i];
}
double[] r1 = new double[r.length];
for(int i = 0; i<A.length;i++){
double temp = 0;
for(int j = 0 ;j < A[i].length;j++){
temp += r[j] * A[i][j];
}
r1[i] = temp;
}
double denominator = 0;
for(int i = 0;i<r.length;i++) {//计算步长式子的分母项
denominator += r1[i] * r[i];
}
double a = element/denominator;//计算步长
for(int i = 0;i<r.length;i++){
x[i] = x[i]+a*r[i];
}
number++;
double[] err = new double[x.length];
for(int i = 0; i< x.length; i++){//计算r的绝对值为err值，即无穷范数
err[i] = Math.abs(r[i]);
}
Arrays.sort(err);
if(err[err.length-1] < 0.0001){//当r的绝对范数小于误差
System.out.println(String.format("最速下降迭代法（r）：迭代次数k =%d,其X值为：",number)
+Arrays.toString(x));
break;
}
x0 = Arrays.copyOf(x,x.length);
}
}
static void CG(double[][] A,double[] b,double[] x0){//共轭梯度法方法
double[] x =Arrays.copyOf(x0,x0.length);//复制内容
int number = 0;
double[] r = new double[x0.length];
for(int i = 0; i<A.length;i++){//求r
double temp = 0;
for(int j = 0 ;j < A[i].length;j++){
temp += x[j] * A[i][j];
}
r[i] = b[i] - temp;
}
double[] p = Arrays.copyOf(r,r.length);//p(0) = r(0) = b
double[] ap = new double[r.length];
while(true){
double[] r0 = Arrays.copyOf(r,r.length);//保留本次r值
for(int i = 0; i<A.length;i++){//求AP
double temp = 0;
for(int j = 0 ;j < A[i].length;j++){
temp += p[j] * A[i][j];
}
ap[i] = temp;
}
double element = 0 ;
for(int i = 0;i<r.length;i++){//求r与r的内积
element += r[i]*r[i];
}
double denominator = 0;
for(int i = 0;i<r.length;i++) {//求ap与p的内积
denominator += p[i] * ap[i];
}
double a = element/denominator;//求a
for(int i = 0;i<r.length;i++){//求x迭代一次的值
x[i] = x[i]+a*p[i];
}
number++;//记录迭代次数
for(int i = 0; i<A.length;i++){//求下一次的r
r[i] = r[i] - a*ap[i];
}
double[] err = new double[x.length];
for(int i = 0; i< x.length; i++){//计算梯度r的无穷范数即无穷范数
err[i] = Math.abs(r[i]);
}
Arrays.sort(err);
if(err[err.length-1]<0.0001){//终止迭代条件，梯度趋于0
System.out.println(String.format("共轭梯度迭代法：迭代次数k =%d,其X值为：",number)
+Arrays.toString(x));
break;
}
x0 = Arrays.copyOf(x,x.length);
double element1 = 0 ;
for(int i = 0;i<r.length;i++){//求r(k+1)与r(k+1)的内积
element1 += r[i]*r[i];
}
double denominator1 = 0;
for(int i = 0;i<r.length;i++) {//求r(k)与r(k)的内积
denominator1 += r0[i] * r0[i];
}
double beta = element1/denominator1;//求beta
for(int i = 0; i<A.length;i++){//求下一次的p
p[i] = r[i] + beta*p[i];
}
}
}
}四、结果展示
最速下降法用||xk-xk-1
||<ε
和||rk
||<ε
两种方式判断迭代结束。
扫一扫
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
