wmlgmatrix.inl

Wild Math Library数值计算库

        }
    }
    else
    {
        Deallocate();
        m_iRows = 0;
        m_iCols = 0;
        m_iQuantity = 0;
        m_afData = 0;
        m_aafEntry = 0;
    }
    return *this;
}
//----------------------------------------------------------------------------
template <class Real>
bool GMatrix<Real>::operator== (const GMatrix& rkM) const
{
    return memcmp(m_afData,rkM.m_afData,m_iQuantity*sizeof(Real)) == 0;
}
//----------------------------------------------------------------------------
template <class Real>
bool GMatrix<Real>::operator!= (const GMatrix& rkM) const
{
    return memcmp(m_afData,rkM.m_afData,m_iQuantity*sizeof(Real)) != 0;
}
//----------------------------------------------------------------------------
template <class Real>
int GMatrix<Real>::CompareArrays (const GMatrix& rkM) const
{
    return memcmp(m_afData,rkM.m_afData,m_iQuantity*sizeof(Real));
}
//----------------------------------------------------------------------------
template <class Real>
bool GMatrix<Real>::operator<  (const GMatrix& rkM) const
{
    return CompareArrays(rkM) < 0;
}
//----------------------------------------------------------------------------
template <class Real>
bool GMatrix<Real>::operator<= (const GMatrix& rkM) const
{
    return CompareArrays(rkM) <= 0;
}
//----------------------------------------------------------------------------
template <class Real>
bool GMatrix<Real>::operator>  (const GMatrix& rkM) const
{
    return CompareArrays(rkM) > 0;
}
//----------------------------------------------------------------------------
template <class Real>
bool GMatrix<Real>::operator>= (const GMatrix& rkM) const
{
    return CompareArrays(rkM) >= 0;
}
//----------------------------------------------------------------------------
template <class Real>
GMatrix<Real> GMatrix<Real>::operator+ (const GMatrix& rkM) const
{
    GMatrix<Real> kSum(rkM.m_iRows,rkM.m_iCols);
    for (int i = 0; i < m_iQuantity; i++)
        kSum.m_afData[i] = m_afData[i] + rkM.m_afData[i];
    return kSum;
}
//----------------------------------------------------------------------------
template <class Real>
GMatrix<Real> GMatrix<Real>::operator- (const GMatrix& rkM) const
{
    GMatrix<Real> kDiff(rkM.m_iRows,rkM.m_iCols);
    for (int i = 0; i < m_iQuantity; i++)
        kDiff.m_afData[i] = m_afData[i] - rkM.m_afData[i];
    return kDiff;
}
//----------------------------------------------------------------------------
template <class Real>
GMatrix<Real> GMatrix<Real>::operator* (const GMatrix& rkM) const
{
    // 'this' is RxN, 'M' is NxC, 'product = this*M' is RxC
    assert( m_iCols == rkM.m_iRows );
    GMatrix<Real> kProd(m_iRows,rkM.m_iCols);
    for (int iRow = 0; iRow < kProd.m_iRows; iRow++)
    {
        for (int iCol = 0; iCol < kProd.m_iCols; iCol++)
        {
            for (int iMid = 0; iMid < m_iCols; iMid++)
            {
                kProd.m_aafEntry[iRow][iCol] += m_aafEntry[iRow][iMid] *
                    rkM.m_aafEntry[iMid][iCol];
            }
        }
    }
    return kProd;
}
//----------------------------------------------------------------------------
template <class Real>
GMatrix<Real> GMatrix<Real>::operator* (Real fScalar) const
{
    GMatrix<Real> kProd(m_iRows,m_iCols);
    for (int i = 0; i < m_iQuantity; i++)
        kProd.m_afData[i] = fScalar*m_afData[i];
    return kProd;
}
//----------------------------------------------------------------------------
template <class Real>
GMatrix<Real> GMatrix<Real>::operator/ (Real fScalar) const
{
    GMatrix<Real> kQuot(m_iRows,m_iCols);
    int i;

    if ( fScalar != (Real)0.0 )
    {
        Real fInvScalar = ((Real)1.0)/fScalar;
        for (i = 0; i < m_iQuantity; i++)
            kQuot.m_afData[i] = fInvScalar*m_afData[i];
    }
    else
    {
        for (i = 0; i < m_iQuantity; i++)
            kQuot.m_afData[i] = Math<Real>::MAX_REAL;
    }

    return kQuot;
}
//----------------------------------------------------------------------------
template <class Real>
GMatrix<Real> GMatrix<Real>::operator- () const
{
    GMatrix<Real> kNeg(m_iRows,m_iCols);
    for (int i = 0; i < m_iQuantity; i++)
        kNeg.m_afData[i] = -m_afData[i];
    return kNeg;
}
//----------------------------------------------------------------------------
template <class Real>
GMatrix<Real> Wml::operator* (Real fScalar, const GMatrix<Real>& rkM)
{
    GMatrix<Real> kProd(rkM.GetRows(),rkM.GetColumns());
    const Real* afMEntry = rkM;
    Real* afPEntry = kProd;
    for (int i = 0; i < rkM.GetQuantity(); i++)
        afPEntry[i] = fScalar*afMEntry[i];
    return kProd;
}
//----------------------------------------------------------------------------
template <class Real>
GMatrix<Real>& GMatrix<Real>::operator+= (const GMatrix& rkM)
{
    for (int i = 0; i < m_iQuantity; i++)
        m_afData[i] += rkM.m_afData[i];
    return *this;
}
//----------------------------------------------------------------------------
template <class Real>
GMatrix<Real>& GMatrix<Real>::operator-= (const GMatrix& rkM)
{
    for (int i = 0; i < m_iQuantity; i++)
        m_afData[i] -= rkM.m_afData[i];
    return *this;
}
//----------------------------------------------------------------------------
template <class Real>
GMatrix<Real>& GMatrix<Real>::operator*= (Real fScalar)
{
    for (int i = 0; i < m_iQuantity; i++)
        m_afData[i] *= fScalar;
    return *this;
}
//----------------------------------------------------------------------------
template <class Real>
GMatrix<Real>& GMatrix<Real>::operator/= (Real fScalar)
{
    int i;

    if ( fScalar != (Real)0.0 )
    {
        Real fInvScalar = ((Real)1.0)/fScalar;
        for (i = 0; i < m_iQuantity; i++)
            m_afData[i] *= fInvScalar;
    }
    else
    {
        for (i = 0; i < m_iQuantity; i++)
            m_afData[i] = Math<Real>::MAX_REAL;
    }

    return *this;
}
//----------------------------------------------------------------------------
template <class Real>
GMatrix<Real> GMatrix<Real>::Transpose () const
{
    GMatrix<Real> kTranspose(m_iCols,m_iRows);
    for (int iRow = 0; iRow < m_iRows; iRow++)
    {
        for (int iCol = 0; iCol < m_iCols; iCol++)
            kTranspose.m_aafEntry[iCol][iRow] = m_aafEntry[iRow][iCol];
    }
    return kTranspose;
}
//----------------------------------------------------------------------------
template <class Real>
GMatrix<Real> GMatrix<Real>::TransposeTimes (const GMatrix& rkM) const
{
    // P = A^T*B, P[r][c] = sum_m A[m][r]*B[m][c]
    assert( m_iRows == rkM.m_iRows );
    GMatrix<Real> kProd(m_iCols,rkM.m_iCols);
    for (int iRow = 0; iRow < kProd.m_iRows; iRow++)
    {
        for (int iCol = 0; iCol < kProd.m_iCols; iCol++)
        {
            for (int iMid = 0; iMid < m_iRows; iMid++)
            {
                kProd.m_aafEntry[iRow][iCol] += m_aafEntry[iMid][iRow] *
                    rkM.m_aafEntry[iMid][iCol];
            }
        }
    }
    return kProd;
}
//----------------------------------------------------------------------------
template <class Real>
GMatrix<Real> GMatrix<Real>::TimesTranspose (const GMatrix& rkM) const
{
    // P = A*B^T, P[r][c] = sum_m A[r][m]*B[c][m]
    assert( m_iCols == rkM.m_iCols );
    GMatrix<Real> kProd(m_iRows,rkM.m_iRows);
    for (int iRow = 0; iRow < kProd.m_iRows; iRow++)
    {
        for (int iCol = 0; iCol < kProd.m_iCols; iCol++)
        {
            for (int iMid = 0; iMid < m_iCols; iMid++)
            {
                kProd.m_aafEntry[iRow][iCol] +=  m_aafEntry[iRow][iMid] *
                    rkM.m_aafEntry[iCol][iRow];
            }
        }
    }
    return kProd;
}
//----------------------------------------------------------------------------
template <class Real>
GVector<Real> GMatrix<Real>::operator* (const GVector<Real>& rkV) const
{
    assert( rkV.GetSize() == m_iCols );
    GVector<Real> kProd(m_iRows);
    for (int iRow = 0; iRow < m_iRows; iRow++)
    {
        for (int iCol = 0; iCol < m_iCols; iCol++)
            kProd[iRow] += m_aafEntry[iRow][iCol]*rkV[iCol];
            
    }
    return kProd;
}
//----------------------------------------------------------------------------
template <class Real>
GVector<Real> Wml::operator* (const GVector<Real>& rkV,
    const GMatrix<Real>& rkM)
{
    assert( rkV.GetSize() == rkM.GetRows() );
    GVector<Real> kProd(rkM.GetColumns());
    Real* afPEntry = kProd;
    for (int iCol = 0; iCol < rkM.GetColumns(); iCol++)
    {
        for (int iRow = 0; iRow < rkM.GetRows(); iRow++)
            afPEntry[iCol] += rkV[iRow]*rkM[iRow][iCol];
    }
    return kProd;
}
//----------------------------------------------------------------------------
template <class Real>
Real GMatrix<Real>::QForm (const GVector<Real>& rkU, const GVector<Real>& rkV)
    const
{
    assert( rkU.GetSize() == m_iRows && rkV.GetSize() == m_iCols );
    return rkU.Dot((*this)*rkV);
}
//----------------------------------------------------------------------------