Рассказать подробно про все методы конечно же очень трудно, но мне эта тема кажется интересной и чрезвычайно важной, поскольку с задачей нахождения решения все сталкиваются достаточно часто. В первой статье Почему Гаусс? был описан метод Гаусса (в том числе с модификацими) и некоторые итерационные методы. Однако, учитывая критику Sinn3r, я решил описать и другие методы.

Начнем сначала

Итак, пусть нам снова нужно решить систему линейных алгебраических уравнений порядка $$. Одним из наиболее ярких численных методов является метод, который использует $$-разложение матрицы.

Это интуитивно понятный метод заключается в следующем. Пусть у нас есть система линейных алгебраических уравнений (записанная в вектором виде):

$$

где $$ — страшная и запутанная невырожденная матрица. Представим ее в виде произведения двух других матриц: $$ — нижняя треугольная матрица, и $$ — верхняя треугольная матрица. Тогда очевидно система принимает вид:

$$

Если обозначить $$, то решение исходной системы находится из решения двух других систем:

$$

$$

Преимущество этого метода в том, что решения двух вспомогательных систем находится очень просто (в силу того, что матрицы $$ и $$ — треугольные). Но легко ли найти эти матрицы? Итак, задача решения системы сводится к задаче построения $$-разложения для матрицы.

Описание алгоритма, применяющего $$-разложение достаточно подробно описано здесь. А мы пока посмотрим на другой метод.

Метод квадратого корня

Наложим дополнитльные условия на матрицу $$. Потребуем, чтобы она была симметрическая, то есть $$ или $$. Такую матрицу можно представить в виде

$$

где $$ — верхняя треугольная матрица, $$ — диагональная вещественная матрица, причем ее диагональные элементы по модулю равны единице. Аналогичным образом исходная система сводится к двум другим:

$$

$$

которые тоже решаются элементарно в силу свойств матриц $$ и $$.

Вывод алгоритмических формул довольно громоздкий и остается за рамками публикации. При желании, их можно найти здесь.

Пример кода на Java, реализующий метод прогонки:

import java.io.FileNotFoundException;
import java.io.FileReader;
import java.io.PrintWriter;
import java.util.Locale;
import java.util.Scanner;

public class SquareRootMethod {

    public static void main(String[] args) throws FileNotFoundException {
        Scanner scanner = new Scanner(new FileReader("input.txt"));
        scanner.useLocale(Locale.US);

        int n = scanner.nextInt();

        double[][] a = new double[n + 1][n + 1];
        double[] b = new double[n + 1];

        for (int i = 1; i <= n; i++) {
            for (int j = 1; j <= n; j++) {
                a[i][j] = scanner.nextDouble();
            }
            b[i] = scanner.nextDouble();
        }

        double[] x = new double[n + 1];

        double[] d = new double[n + 1];
        double[][] s = new double[n + 1][n + 1];

        // for i = 1
        d[1] = Math.signum(a[1][1]);
        s[1][1] = Math.sqrt(Math.abs(a[1][1]));
        for (int j = 2; j <= n; j++) {
            s[1][j] = a[1][j] / (s[1][1] * d[1]);
        }

        // for i > 1
        //searching S and D matrix
        for (int i = 2; i <= n; i++) {
            double sum = 0;
            for (int k = 1; k <= i - 1; k++) {
                sum += s[k][i] * s[k][i] * d[k];
            }
            d[i] = Math.signum(a[i][i] - sum);
            s[i][i] = Math.sqrt(Math.abs(a[i][i] - sum));

            double l = 1 / (s[i][i] * d[i]);
            for (int j = i + 1; j <= n; j++) {
                double SDSsum = 0;
                for (int k = 1; k <= i - 1; k++) {
                    SDSsum += s[k][i] * d[k] * s[k][j];
                }
                s[i][j] = l * (a[i][j] - SDSsum);
            }
        }

        //solve of the system (s^t * d)y = b
        double[] y = new double[n + 1];
        y[1] = b[1] / (s[1][1] * d[1]);

        for (int i = 2; i <= n; i++) {
            double sum = 0;

            for (int j = 1; j <= i - 1; j++) {
                sum += s[j][i] * d[j] * y[j];
            }

            y[i] = (b[i] - sum) / (s[i][i] * d[i]);
        }

        //solve of the system sx = y
        x[n] = y[n] / s[n][n];

        for (int i = n - 1; i >= 1; i--) {
            double sum = 0;

            for (int k = i + 1; k <= n; k++) {
                sum += s[i][k] * x[k];
            }

            x[i] = (y[i] - sum) / s[i][i];
        }

        //output
        PrintWriter pw = new PrintWriter("output.txt");
        for (int i = 1; i <= n; i++) {
            pw.printf(Locale.US, "x%d = %f\n", i, x[i]);
        }
        pw.flush();
        pw.close();
        
    }

}

Метод прогонки

С помощью этого метода можно решать только специфические системы, имеющие не более трех неизвестных в каждой строке. То есть при системе

$$

матрица $$ является трехдиагональной:

$$

Сразу заметим, что имеется связь соседних решений:

$$

$$ — какие-то неизвестные числа. Если мы найдем их и какую-то одну переменную, то сможем найти и все остальные.

Вывод формул

$$

$$

присутствует здесь. Ну и в итоге

$$

Отметим, что в формулах поиска $$ присутствует деление на число $$, которое может оказаться нулем, что нужно отслеживать. Но на самом деле имеет место следующее утверждение, доказательство которого есть здесь: алгоритм прогонки является корректным и устойчивым, если выполняются условия:

$$

$$

Рассмотрим код программы на Java, который решает систему

$$

при $$.

package SystemOfLinearAlgebraicEquations;

import java.io.FileNotFoundException;
import java.io.FileReader;
import java.io.PrintWriter;
import java.util.Locale;
import java.util.Scanner;

public class TridiagonalMatrixAlgorithm {

	public static void main(String[] args) throws FileNotFoundException {

		final int n = 21;
		double[] A = new double[n + 1];
		double[] B = new double[n + 1];
		double[] C = new double[n + 1];
		double[] F = new double[n + 1];
		double xi = 2;
		
		C[1] = 1;
		F[1] = 0;
		for(int i = 2; i <= n-1; i++) {
			A[i] = 1;
			C[i] = xi;
			B[i] = 1;
			F[i] =  - (double) i / n;
		}
		A[n] = 0;
		C[n] = 1;
		F[n] = 1;

		double[] alpha = new double[n + 1];
		double[] beta = new double[n + 1];

		// right stroke
		alpha[1] = B[1] / C[1];
		beta[1] = F[1] / C[1];
		for (int i = 2; i <= n - 1; i++) {
			double denominator = C[i] - A[i] * alpha[i-1];
			alpha[i] = B[i] / denominator;
			beta[i] = (F[i] + A[i] * beta[i - 1]) / denominator;
		}

		double[] x = new double[n + 1];

		// back stroke
		x[n] = (F[n] + A[n] * beta[n - 1]) / (C[n] - A[n] * alpha[n - 1]);
		for (int i = n - 1; i >= 1; i--) {
			x[i] = alpha[i] * x[i + 1] + beta[i];
		}

		// output
		PrintWriter pw = new PrintWriter("output.txt");
		for (int i = 1; i <= n; i = i + 5) {
			pw.printf(Locale.US, "x%d = %f\n", i, x[i]);
		}
		pw.flush();
		pw.close();

	}

}