MyMathLib系列(行列式计算4--向量部分)

1)将向量组进行消元，变换成阶梯矩阵，这是求向量组的极大线性无关组的基本算法。这个方法在前面曾经给出过，但在这里做了改进，目的是为了可以判断是否线性相关：

        /// <summary>
        /// 方程组消元，最后一列为系数,结果就在CoefficientDeterminant里.
        /// 本算法也可以用来求矩阵的秩.
        /// </summary>
        /// <param name="CoefficientDeterminant">方程组系数数组</param>
        public static void EquationsElimination(decimal[,] CoefficientDeterminant)
        {
            var theRowCount = CoefficientDeterminant.GetLength(0);
            var theColCount = CoefficientDeterminant.GetLength(1);
            int theN = theRowCount;
            int theE = theColCount-1;
            //从第1列到第theE-1列,最后一列不用处理.
            for (int i = 0; i < theE; i++)
            {
                //从第theN-1行到第1行，将D[j,i]依次变为0,需要注意的是：
                //如果第j-1行，的左元素全部为0，才能继续交换.
                for (int j = theN - 1; j > 0; j--)
                {
                    //如果为当前值为0，则不处理,继续处理上一行
                    if (CoefficientDeterminant[j, i] == 0)
                    {
                        continue;
                    }
                    //如果左上邻元素[j-1, i-1]以及其左边的元素都为0方可交换
                    //因为当前元素的左边元素已经全部是零，因此如果要交换不能使本行左边产生非零数，
                    //则需要左上邻及其所有元素皆为0.
                    for (int s = i - 1; s >= 0; s--)
                    {
                        if (CoefficientDeterminant[j - 1, s] != 0)
                        {
                            break;
                        }
                    }
                    //如果[j,i]的上一行[j-1, i]的值为0则交换
                    if (CoefficientDeterminant[j - 1, i] == 0)
                    {
                        for (int k = 0; k <= theE; k++)//这里要交换常数项，所以是k <= theE
                        {
                            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 <= theE; k++)//这里要计算常数项，所以是k <= theE
                        //{
                        //    CoefficientDeterminant[j, k] = CoefficientDeterminant[j, k] - CoefficientDeterminant[j - 1, k] * theRate;
                        //}
                        //改进：做乘法，可以避免小数换算带来的误差
                        var theRate2 = CoefficientDeterminant[j, i];
                        var theRate1 = CoefficientDeterminant[j - 1, i];
                        for (int k = 0; k <= theE; k++)//这里要计算常数项，所以是k <= theE
                        {
                            CoefficientDeterminant[j, k] = CoefficientDeterminant[j, k] * theRate1 - CoefficientDeterminant[j - 1, k] * theRate2;
                        }
                    }
                }
            }
        }

后面的矩阵求秩的过程也会用到这种消元法，不过后面的矩阵求秩在这个算法上略有修改，采用了另外一种方式。下面是向量相关的算法：

/*Vector.cs
 * Albert.Tian on 20141225 
 */
using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;

namespace MyMathLib
{
    /// <summary>
    /// 向量相关运算，这里没有给出向量的基本变换。因为
    /// 向量组可以用矩阵表达，因此很多向量的运算最终都归于矩阵运算。
    /// 本算法中的求极大线性无关组，是根据课本算法得出，不是最简单的一种方法，
    /// 在矩阵表达中，可以根据矩阵的初等变换，很容易求到极大线性无关组.
    /// </summary>
    public class Vectors
    {
        /// <summary>
        /// 判断从0到指定结束位的元素是否全是0.
        /// </summary>
        /// <param name="CurrVector">当前向量</param>
        /// <param name="EndIndex">结束索引</param>
        /// <returns></returns>
        public static bool IsZeroVector(decimal[] CurrVector, int EndIndex = 0)
        {
            int theEndIndex = EndIndex;
            if (theEndIndex == 0)
            {
                theEndIndex = CurrVector.Length;
            }

            var theIsZeroVector = true;
            for (int i = 0; i < theEndIndex; i++)
            {
                if (CurrVector[i] != 0)
                {
                    theIsZeroVector = false;
                    break;
                }
            }
            return theIsZeroVector;
        }
        /// <summary>
        /// 判断当前向量是否与前面的向量线性相关，采用的是消元化
        /// </summary>
        /// <param name="CurrVector">当前向量</param>
        /// <param name="PrevMaxLIVectors">已有的线性无关的向量组</param>
        /// <returns></returns>
        private static bool IsLinearCorrelation(decimal[] CurrVector, List<decimal[]> PrevMaxLIVectors)
        {
            //零向量必是线性相关
            var theIsZeroVector = IsZeroVector(CurrVector);
            if (theIsZeroVector)
            {
                return true;
            }
            //如果前面没有向量，则当前非零向量必是线性无关的.
            if (PrevMaxLIVectors.Count <= 0)
            {
                return false;
            }
            //构造方程式，判断是否线性相关
            var theECount = CurrVector.Length;
            var theVCount = PrevMaxLIVectors.Count;
            
            //加入当前向量为常量向量,因此方程组的列为theVCount+1.
            var theEquealGroup = new decimal[theECount, theVCount + 1];
            //导入PrevMaxLIVectors
            for (int i = 0; i < theVCount; i++)
            {
                for (int j = 0; j < theECount; j++)
                {
                    theEquealGroup[j, i] = PrevMaxLIVectors[i][j];
                }
            }

            //加入当前向量作为常量向量
            for (int j = 0; j < theECount; j++)
            {
                theEquealGroup[j, theVCount] = CurrVector[j];
            }
            //消元,这里的消元法变成了求矩阵秩的基本算法.
            LinearAlgebra.EquationsElimination(theEquealGroup);
            //如果theRA<theRA1，则有矛盾方程，必然线性无关，如果theRA>=theRA1，则至少有一个非零解。
            //这里如此判断，主要是因为消元的方法。当然，这在矩阵有证。
            var theRA = 0;
            var theRA1 = 0;
            for (int i = 0; i < theECount; i++)
            {
                if (!IsZeroVector(theEquealGroup.GetVector(i), theVCount))
                {
                    theRA++;
                }
                if (!IsZeroVector(theEquealGroup.GetVector(i)))
                {
                    theRA1++;
                }
            }
            return (theRA >= theRA1);
        }
        /// <summary>
        /// 这个函数暂时无用，判断两个向量的等价性
        /// </summary>
        /// <param name="V1"></param>
        /// <param name="V2"></param>
        /// <returns></returns>
        private static bool IsEquivalent(decimal[] V1, decimal[] V2)
        {
            var theV1IsZeroVector = IsZeroVector(V1);
            var theV2IsZeroVector = IsZeroVector(V2);
            //两个都是0向量
            if (theV1IsZeroVector && theV2IsZeroVector)
            {
                return true;
            }
            //其中一个零向量，另一个是非零向量
            if (theV1IsZeroVector || theV2IsZeroVector)
            {
                return false;
            }
            //一个分量肯定成比例
            if (V1.Length <= 1)
            {
                return true;
            }

            decimal thePreRate = 0;
            var theCount = V1.Length;
            var theIndex = 0;
            for (int i = 0; i < theCount; i++)
            {
                //如果都等于0，则不比较
                if (V1[i] == V2[i] && V2[i] == 0)
                {
                    continue;
                }
                //如果其中一个是0，则必然不成比例,不等价
                if (V1[i] == 0 || V2[i] == 0)
                {
                    return false;
                }
                if (theIndex == 0)
                {
                    thePreRate = V1[i] / V2[i];
                    theIndex++;
                }
                else
                {

                    if (thePreRate != V1[i] / V2[i])
                    {
                        return false;
                    }

                }
            }
            return true;
        }
        /// <summary>
        /// 获取向量组的极大线性无关组，方法是依照教科书写的。
        /// </summary>
        /// <param name="Vectors">向量组</param>
        /// <returns></returns>
        public static List<decimal[]> GetMaxLinearIndependentGroup(decimal[,] Vectors)
        {
            //向量个数
            var theVCount = Vectors.GetLength(0);
            //元数个数
            var theECount = Vectors.GetLength(1);
            List<decimal[]> theTempVectors = new List<decimal[]>();
            for (int i = 0; i < theVCount; i++)
            {
                var theTempVector = Vectors.GetVector(i);
                //如果当前向量不能被前面的向量线性表出，则加入到极大组
                //能不能线性表出，实际上就是求当前向量加入前面的向量组后，
                //还是不是线性相关
                if (!IsLinearCorrelation(theTempVector, theTempVectors))
                {
                    theTempVectors.Add(theTempVector);
                }
            }
            return theTempVectors;
        }
        
    }
}