使用C#计算矩阵零空间：MathNet代码问题及正确实现请求

问题分析 #

你编写的扩展方法能得到正确的行列数，但数值不符合预期，核心问题是未考虑数值计算的容差设置，且对Kernel()返回基向量的处理存在细节偏差，另外Kernel()默认的基向量形式和MatlabNull()的有理/正交基逻辑不一致。

修正后的Kernel()实现版本 #

这个版本基于MathNet的Kernel()方法，修正边界处理并显式设置容差，匹配数值计算需求：

using MathNet.Numerics.LinearAlgebra;

public static class MatrixExtensions
{
    public static Matrix<double> NullMatrix(this Matrix<double> matrix)
    {
        // 处理空矩阵或零空间为空的满秩矩阵
        if (matrix == null || matrix.RowCount == 0 || matrix.ColumnCount == 0)
        {
            return Matrix<double>.Build.Dense(matrix?.ColumnCount ?? 0, 0);
        }

        // 计算容差：基于矩阵范数、机器精度和矩阵维度，避免数值误差干扰
        double tolerance = matrix.L2Norm() * double.Epsilon * Math.Max(matrix.RowCount, matrix.ColumnCount);
        // 获取零空间的基向量集合
        var nullSpaceBasis = matrix.Kernel(tolerance);

        int vectorDimension = matrix.ColumnCount;
        int basisCount = nullSpaceBasis.Count();

        // 构建结果矩阵：行数为原矩阵列数，列数为基向量个数
        var result = Matrix<double>.Build.Dense(vectorDimension, basisCount);

        int columnIndex = 0;
        foreach (var basisVector in nullSpaceBasis)
        {
            // 将每个基向量直接赋值为结果矩阵的对应列
            result.SetColumn(columnIndex, basisVector);
            columnIndex++;
        }

        return result;
    }
}

匹配MatlabNull('svd')的正交基版本 #

如果需要和MatlabNull()默认的正交基行为一致（返回标准正交基），可通过奇异值分解(SVD)实现：

using MathNet.Numerics.LinearAlgebra;

public static class MatrixExtensions
{
    public static Matrix<double> NullMatrixOrthogonal(this Matrix<double> matrix)
    {
        if (matrix == null || matrix.RowCount == 0 || matrix.ColumnCount == 0)
        {
            return Matrix<double>.Build.Dense(matrix?.ColumnCount ?? 0, 0);
        }

        // 执行奇异值分解
        var svd = matrix.Svd(true);
        double tolerance = svd.SingularValues.Max() * double.Epsilon * Math.Max(matrix.RowCount, matrix.ColumnCount);

        // 确定矩阵秩，计算零空间维度
        int rank = svd.SingularValues.Count(s => s > tolerance);
        int nullity = matrix.ColumnCount - rank;

        if (nullity == 0)
        {
            return Matrix<double>.Build.Dense(matrix.ColumnCount, 0);
        }

        // 提取对应右奇异向量，转置后作为零空间的正交基矩阵
        return svd.VT.SubMatrix(rank, nullity, 0, matrix.ColumnCount).Transpose();
    }
}

关键说明 #

• Kernel()方法返回零空间的有理基（非正交），SVD方法返回标准正交基，对应MatlabNull('svd')的行为
• 容差设置直接影响零空间的判定，必须结合矩阵特性和机器精度计算，避免数值误差导致的错误
• 必须处理零空间为空的情况（满秩矩阵），此时返回行数为原矩阵列数、列数为0的空矩阵

内容的提问来源于stack exchange，提问作者Ajinkya Jagtap