MOAB  4.9.3pre
Eigen::internal::ldlt_inplace< Lower > Struct Template Reference

#include <LDLT.h>

List of all members.

Static Public Member Functions

template<typename MatrixType , typename TranspositionType , typename Workspace >
static bool unblocked (MatrixType &mat, TranspositionType &transpositions, Workspace &temp, SignMatrix &sign)
template<typename MatrixType , typename WDerived >
static bool updateInPlace (MatrixType &mat, MatrixBase< WDerived > &w, const typename MatrixType::RealScalar &sigma=1)
template<typename MatrixType , typename TranspositionType , typename Workspace , typename WType >
static bool update (MatrixType &mat, const TranspositionType &transpositions, Workspace &tmp, const WType &w, const typename MatrixType::RealScalar &sigma=1)

Detailed Description

template<>
struct Eigen::internal::ldlt_inplace< Lower >

Definition at line 254 of file LDLT.h.


Member Function Documentation

template<typename MatrixType , typename TranspositionType , typename Workspace >
static bool Eigen::internal::ldlt_inplace< Lower >::unblocked ( MatrixType &  mat,
TranspositionType &  transpositions,
Workspace &  temp,
SignMatrix sign 
) [inline, static]

Definition at line 257 of file LDLT.h.

  {
    using std::abs;
    typedef typename MatrixType::Scalar Scalar;
    typedef typename MatrixType::RealScalar RealScalar;
    typedef typename TranspositionType::StorageIndex IndexType;
    eigen_assert(mat.rows()==mat.cols());
    const Index size = mat.rows();

    if (size <= 1)
    {
      transpositions.setIdentity();
      if (numext::real(mat.coeff(0,0)) > 0) sign = PositiveSemiDef;
      else if (numext::real(mat.coeff(0,0)) < 0) sign = NegativeSemiDef;
      else sign = ZeroSign;
      return true;
    }

    for (Index k = 0; k < size; ++k)
    {
      // Find largest diagonal element
      Index index_of_biggest_in_corner;
      mat.diagonal().tail(size-k).cwiseAbs().maxCoeff(&index_of_biggest_in_corner);
      index_of_biggest_in_corner += k;

      transpositions.coeffRef(k) = IndexType(index_of_biggest_in_corner);
      if(k != index_of_biggest_in_corner)
      {
        // apply the transposition while taking care to consider only
        // the lower triangular part
        Index s = size-index_of_biggest_in_corner-1; // trailing size after the biggest element
        mat.row(k).head(k).swap(mat.row(index_of_biggest_in_corner).head(k));
        mat.col(k).tail(s).swap(mat.col(index_of_biggest_in_corner).tail(s));
        std::swap(mat.coeffRef(k,k),mat.coeffRef(index_of_biggest_in_corner,index_of_biggest_in_corner));
        for(Index i=k+1;i<index_of_biggest_in_corner;++i)
        {
          Scalar tmp = mat.coeffRef(i,k);
          mat.coeffRef(i,k) = numext::conj(mat.coeffRef(index_of_biggest_in_corner,i));
          mat.coeffRef(index_of_biggest_in_corner,i) = numext::conj(tmp);
        }
        if(NumTraits<Scalar>::IsComplex)
          mat.coeffRef(index_of_biggest_in_corner,k) = numext::conj(mat.coeff(index_of_biggest_in_corner,k));
      }

      // partition the matrix:
      //       A00 |  -  |  -
      // lu  = A10 | A11 |  -
      //       A20 | A21 | A22
      Index rs = size - k - 1;
      Block<MatrixType,Dynamic,1> A21(mat,k+1,k,rs,1);
      Block<MatrixType,1,Dynamic> A10(mat,k,0,1,k);
      Block<MatrixType,Dynamic,Dynamic> A20(mat,k+1,0,rs,k);

      if(k>0)
      {
        temp.head(k) = mat.diagonal().real().head(k).asDiagonal() * A10.adjoint();
        mat.coeffRef(k,k) -= (A10 * temp.head(k)).value();
        if(rs>0)
          A21.noalias() -= A20 * temp.head(k);
      }
      
      // In some previous versions of Eigen (e.g., 3.2.1), the scaling was omitted if the pivot
      // was smaller than the cutoff value. However, since LDLT is not rank-revealing
      // we should only make sure that we do not introduce INF or NaN values.
      // Remark that LAPACK also uses 0 as the cutoff value.
      RealScalar realAkk = numext::real(mat.coeffRef(k,k));
      if((rs>0) && (abs(realAkk) > RealScalar(0)))
        A21 /= realAkk;

      if (sign == PositiveSemiDef) {
        if (realAkk < 0) sign = Indefinite;
      } else if (sign == NegativeSemiDef) {
        if (realAkk > 0) sign = Indefinite;
      } else if (sign == ZeroSign) {
        if (realAkk > 0) sign = PositiveSemiDef;
        else if (realAkk < 0) sign = NegativeSemiDef;
      }
    }

    return true;
  }
template<typename MatrixType , typename TranspositionType , typename Workspace , typename WType >
static bool Eigen::internal::ldlt_inplace< Lower >::update ( MatrixType &  mat,
const TranspositionType &  transpositions,
Workspace &  tmp,
const WType &  w,
const typename MatrixType::RealScalar &  sigma = 1 
) [inline, static]

Definition at line 385 of file LDLT.h.

  {
    // Apply the permutation to the input w
    tmp = transpositions * w;

    return ldlt_inplace<Lower>::updateInPlace(mat,tmp,sigma);
  }
template<typename MatrixType , typename WDerived >
static bool Eigen::internal::ldlt_inplace< Lower >::updateInPlace ( MatrixType &  mat,
MatrixBase< WDerived > &  w,
const typename MatrixType::RealScalar &  sigma = 1 
) [inline, static]

Definition at line 347 of file LDLT.h.

  {
    using numext::isfinite;
    typedef typename MatrixType::Scalar Scalar;
    typedef typename MatrixType::RealScalar RealScalar;

    const Index size = mat.rows();
    eigen_assert(mat.cols() == size && w.size()==size);

    RealScalar alpha = 1;

    // Apply the update
    for (Index j = 0; j < size; j++)
    {
      // Check for termination due to an original decomposition of low-rank
      if (!(isfinite)(alpha))
        break;

      // Update the diagonal terms
      RealScalar dj = numext::real(mat.coeff(j,j));
      Scalar wj = w.coeff(j);
      RealScalar swj2 = sigma*numext::abs2(wj);
      RealScalar gamma = dj*alpha + swj2;

      mat.coeffRef(j,j) += swj2/alpha;
      alpha += swj2/dj;


      // Update the terms of L
      Index rs = size-j-1;
      w.tail(rs) -= wj * mat.col(j).tail(rs);
      if(gamma != 0)
        mat.col(j).tail(rs) += (sigma*numext::conj(wj)/gamma)*w.tail(rs);
    }
    return true;
  }

The documentation for this struct was generated from the following file:
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines