CShap 对称正定矩阵的逆矩阵 源代码

下面是 C# 中计算对称正定矩阵的逆矩阵的源代码：

using System;

class InverseMatrix
{
    static void Main()
    {
        // 定义矩阵大小和元素
        int n = 3;
        double[,] A = { { 3.0, -1.0, 1.0 }, { -1.0, 2.0, 0.0 }, { 1.0, 0.0, 2.0 } };

        // 计算对称正定矩阵的逆矩阵
        double[,] Ainv = Inverse(A, n);

        // 输出结果
        Console.WriteLine("A:");
        PrintMatrix(A, n);
        Console.WriteLine("Ainv:");
        PrintMatrix(Ainv, n);

        // 等待用户输入
        Console.ReadLine();
    }

    // 打印矩阵
    static void PrintMatrix(double[,] A, int n)
    {
        for (int i = 0; i < n; i++)
        {
            for (int j = 0; j < n; j++)
            {
                Console.Write("{0,8:F4}", A[i, j]);
            }
            Console.WriteLine();
        }
        Console.WriteLine();
    }

    // 计算对称正定矩阵的逆矩阵
    static double[,] Inverse(double[,] A, int n)
    {
        double[,] Ainv = new double[n, n];
        double[,] L = new double[n, n];
        double[,] U = new double[n, n];
        double[,] I = new double[n, n];

        // 初始化单位矩阵
        for (int i = 0; i < n; i++)
        {
            for (int j = 0; j < n; j++)
            {
                I[i, j] = (i == j) ? 1.0 : 0.0;
            }
        }

        // LU 分解
        LUDecomposition(A, n, L, U);

        // 解方程
        for (int k = 0; k < n; k++)
        {
            double[] y = new double[n];

            // 前向替换
            for (int i = 0; i < n; i++)
            {
                double sum = 0.0;
                for (int j = 0; j < i; j++)
                {
                    sum += L[i, j] * y[j];
                }
                y[i] = I[i, k] - sum;
            }

            double[] x = new double[n];

            // 后向替换
            for (int i = n - 1; i >= 0; i--)
            {
                double sum = 0.0;
                for (int j = i + 1; j < n; j++)
                {
                    sum += U[i, j] * x[j];
                }
                x[i] = (y[i] - sum) / U[i, i];
            }

            // 将解转换为矩阵形式
            for (int i = 0; i < n; i++)
            {
                Ainv[i, k] = x[i];
            }
        }

        return Ainv;
    }

    // LU 分解
    static void LUDecomposition(double[,] A, int n, double[,] L, double[,] U)
    {
        for (int i = 0; i < n; i++)
        {
            for (int j = 0; j < n; j++)
            {
                if (j < i)
                {
                    double sum = 0.0;
                    for (int k = 0; k < j; k++)
                    {
                        sum += L[i, k] * U[k, j];
                    }
                    L[i, j] = (A[i, j] - sum) / U[j, j];
                    U[i, j] = 0.0;
                }
                else if (j == i)
                {
                    double sum = 0.0;
                    for (int k = 0; k < j; k++)
                    {
                        sum += L[i, k] * U[k, j];
                    }
                    L[i, j] = 1.0;
                    U[i, j] = A[i, j] - sum;
                }
                else
                {
                    double sum = 0.0;
                    for (int k = 0; k < i; k++)
                    {
                        sum += L[i, k] * U[k, j];
                    }
                    L[i, j] = 0.0;
                    U[i, j] = (A[i, j] - sum) / L[i, i];
                }
            }
        }
    }
}

该代码使用 LU 分解来计算对称正定矩阵的逆矩阵。首先，将原矩阵和单位矩阵进行 LU 分解，然后对每一列单位矩阵求解线性方程组，得到对应的列向量，将这些列向量拼接起来形成逆矩阵。其中，前向替换和后向替换分别用于解下三角矩阵和上三角矩阵的线性方程组。