#SOR-序列超松弛算法

\(\textbf{SOR（Successive Over-Relaxation）}\)算法是一种用于解线性方程组的迭代方法，它通过引入松弛因子来加快收敛速度。SOR算法的基本步骤如下：

1. 将系数矩阵\(A\)分解为\(A=D-L-U\)，其中D是A的对角线元素构成的对角矩阵，\(L\)是\(A\)的下三角部分（不含对角线）构成的矩阵，\(U\)是\(A\)的上三角部分（不含对角线）构成的矩阵。

2. 选择一个适当的松弛因子\(ω\)（通常在\(0\)和\(2\)之间），并且假设初始解向量\(x^{(0)}\)。

3. 对于\(k=0,1,2,...\)进行下面的迭代步骤：

   • 对于\(i=1,2,...,n\)，计算新的解向量中的第i个分量：
     \[x_i^{(k+1)} = (1-\omega)x_i^{(k)} + \left(\frac{\omega}{a_{ii}}\right) \left(b_i - \sum_{j=1}^{i-1}a_{ij}x_j^{(k+1)} - \sum_{j=i+1}^{n}a_{ij}x_j^{(k)}\right) \]

4. 重复步骤3直到满足收敛条件。

SOR算法的优点是相对于\(\textbf{Jacobi或Gauss-Seidel}\)方法，可以更快地收敛到解。然而，选择合适的松弛因子是很关键的，不同的问题可能需要不同的松弛因子才能获得最佳的收敛效果。

c++实现:

#include <bits/stdc++.h>
using namespace std;

vector<double> SOR(vector<vector<double>> A, vector<double> b, double omega, double tol, int max_iter) {
    int n = A.size();
    vector<double> x(n, 0.0); // 初始化解向量
    vector<double> x_new(n, 0.0); // 用于存储新的解向量

    for (int k = 0; k < max_iter; ++k) {
        for (int i = 0; i < n; ++i) {
            double sigma1 = 0.0, sigma2 = 0.0;
            for (int j = 0; j < i; ++j) {
                sigma1 += A[i][j] * x_new[j];
            }
            for (int j = i + 1; j < n; ++j) {
                sigma2 += A[i][j] * x[j];
            }
            x_new[i] = (1 - omega) * x[i] + (omega / A[i][i]) * (b[i] - sigma1 - sigma2);
        }
        // 检查收敛性
        double max_diff = 0.0;
        for (int i = 0; i < n; ++i) {
            max_diff = max(max_diff, abs(x_new[i] - x[i]));
        }
        if (max_diff < tol) {
            return x_new;
        }
        x = x_new;
    }
    return x; // 返回最后的解向量
}

int main() {
    vector<vector<double>> A = {{4, 1, 1}, {1, 3, 2}, {1, 2, 5}};
    vector<double> b = {12, 13, 14};
    double omega = 1.25;
    double tol = 1e-5;
    int max_iter = 100;

    vector<double> x = SOR(A, b, omega, tol, max_iter);

    for (int i = 0; i < x.size(); ++i) {
        cout << "x[" << i << "] = " << x[i] << endl;
    }

    return 0;
}

python实现:

import numpy as np

def SOR(A, b, omega, tol, max_iter):
    n = len(A)
    x = np.zeros(n) # 初始化解向量

    for k in range(max_iter):
        x_new = np.copy(x)
        for i in range(n):
            sigma1 = np.dot(A[i, :i], x_new[:i])
            sigma2 = np.dot(A[i, i+1:], x[i+1:])
            x_new[i] = (1 - omega) * x[i] + (omega / A[i, i]) * (b[i] - sigma1 - sigma2)
        # 检查收敛性
        if np.max(np.abs(x_new - x)) < tol:
            return x_new
        x = x_new
    return x # 返回最后的解向量

A = np.array([[4, 1, 1], [1, 3, 2], [1, 2, 5]])
b = np.array([12, 13, 14])
omega = 1.25
tol = 1e-5
max_iter = 100

x = SOR(A, b, omega, tol, max_iter)

for i in range(len(x)):
    print(f"x[{i}] = {x[i]}")