MyMathLib系列(行列式计算3)

来源:互联网 发布:linux svn 删除仓库 编辑:程序博客网 时间:2024/04/30 04:21

到今天,行列式和线性方程组部分就完成了。

using System;using System.Collections.Generic;using System.Linq;using System.Text;namespace MyMathLib{        /// <summary>    /// 行列式计算,本程序属于MyMathLib的一部分,欢迎使用,参考,提意见。    /// 有时间用函数语言改写,做自己得MathLib,里面的算法经过验证,但没经过    /// 严格测试,如需参考,请慎重.    /// </summary>    public static partial class LinearAlgebra    {#region 线性方程组        /// <summary>        /// 根据拉普拉斯定理计算行列式值。        /// </summary>        /// <param name="Determinants">N阶行列式</param>        /// <returns>计算结果</returns>        public static decimal CalcDeterByLaplaceLaw(decimal[,] Determinants)        {            var n = Determinants.GetLength(0);            //如果阶数小于3,则没必要采用拉普拉斯展开            if (n <= 3)            {                return CalcDeterminantAij(Determinants, false);            }            var theRows = GetLaplaceRowsOdd(n);            return CalcDeterByLaplaceLaw(Determinants, theRows);        }        /// <summary>        /// 求解线性方程组,这里要求N        /// </summary>        /// <param name="CoefficientDeterminant">线性方程组系数行列式</param>        /// <param name="ConstantTerms">常数项</param>        /// <returns></returns>        public static decimal[] LinearEquations(int UnknownElements,decimal[,] CoefficientDeterminant, decimal[] ConstantTerms)        {            var theRowCount = CoefficientDeterminant.GetLength(0);            var theColCount = CoefficientDeterminant.GetLength(1);            if (UnknownElements == theRowCount && theColCount == UnknownElements)            {                var theD = CalcDeterByLaplaceLaw(CoefficientDeterminant);                if(theD==0)                {                    return null;                }                decimal[] theResults = new decimal[UnknownElements];                                for (int i = 1; i <= UnknownElements; i++)                {                    //置换第i列,注意保存原来的值,下次计算时恢复                    var theTemp = new decimal[UnknownElements];                    for (int j = 1; j <= UnknownElements; j++)                    {                        theTemp[j-1] = CoefficientDeterminant[j-1, i-1];                        CoefficientDeterminant[j - 1, i - 1] = ConstantTerms[j - 1];                    }                    var theDi = CalcDeterByLaplaceLaw(CoefficientDeterminant) / theD;                    theResults[i-1] = theDi;                    //复原系数行列式.                    for (int j = 1; j <= UnknownElements; j++)                    {                        CoefficientDeterminant[j - 1, i - 1] = theTemp[j - 1];                    }                }                return theResults;            }            else            {                throw new Exception("参数格式不正确!");            }        }        /// <summary>        /// 求解线性方程组(消元法),这里与化三角求行列式的方法类似,这里要求方程个数和元个数相同。        /// 如果方程个数小于元的个数,消元没问题,但涉及到一般解,会牵扯到符号运算,这里暂不考虑.        /// </summary>        /// <param name="CoefficientDeterminant">线性方程组系数行列式,最右边N+1列是常数项</param>        /// <returns></returns>        public static decimal[] LinearEquationsEM(int UnknownElements,decimal[,] CoefficientDeterminant)        {            var theRowCount = CoefficientDeterminant.GetLength(0);            var theColCount = CoefficientDeterminant.GetLength(1);            if (UnknownElements == theRowCount && theColCount == UnknownElements + 1)            {                decimal[] theResults = new decimal[UnknownElements];                int theN = UnknownElements;                //从第1列到第theN-1列                for (int i = 0; i < theN - 1; i++)                {                    //从第theN-1行到第i+1行,将D[j,i]依次变为0                    for (int j = theN - 1; j > i; j--)                    {                        //如果为当前值为0,则不处理,继续处理上一行                        if (CoefficientDeterminant[j, i] == 0)                        {                            continue;                        }                        //如果[j,i]的上一行[j-1, i]的值为0则交换                        if (CoefficientDeterminant[j - 1, i] == 0)                        {                            for (int k = 0; k <= theN; k++)//这里要交换常数项,所以是k <= theN                            {                                decimal theTmpDec = CoefficientDeterminant[j, k];                                CoefficientDeterminant[j, k] = CoefficientDeterminant[j - 1, k];                                CoefficientDeterminant[j - 1, k] = theTmpDec;                            }                        }                        else                        {                            //将当前行减去上一行与theRate的积。                            var theRate = CoefficientDeterminant[j, i] / CoefficientDeterminant[j - 1, i];                            for (int k = 0; k <= theN; k++)//这里要计算常数项,所以是k <= theN                            {                                CoefficientDeterminant[j, k] = CoefficientDeterminant[j, k] - CoefficientDeterminant[j - 1, k] * theRate;                            }                        }                    }                }                //处理结果                if (CoefficientDeterminant[UnknownElements - 1, UnknownElements - 1] == 0)                {                    if (CoefficientDeterminant[UnknownElements - 1, UnknownElements] == 0)                    {                        throw new Exception("无效方程,有无穷个解!");                    }                    else                    {                        throw new Exception("方程无解!");                    }                }                //结果处理,回代                for (int i = UnknownElements - 1; i >= 0; i--)                {                    //计算已求项                    decimal theTempDec = 0;                    for (int j = i + 1; j < theN; j++)                    {                        theTempDec += CoefficientDeterminant[i, j] * theResults[j];                    }                    //计算结果,如果系数为0,则为无效方程                    if (CoefficientDeterminant[i, i] == 0)                    {                        throw new Exception("无效方程");                    }                    theResults[i] = (CoefficientDeterminant[i, UnknownElements] - theTempDec) / CoefficientDeterminant[i, i];                }                return theResults;            }            else            {                throw new Exception("参数格式不正确!");            }        }                #endregion    }}

0 0
原创粉丝点击