25 #ifndef EIGEN_TRIANGULAR_SOLVER_VECTOR_H
26 #define EIGEN_TRIANGULAR_SOLVER_VECTOR_H
32 template<
typename LhsScalar,
typename RhsScalar,
typename Index,
int Mode,
bool Conjugate,
int StorageOrder>
33 struct triangular_solve_vector<LhsScalar, RhsScalar, Index,
OnTheRight, Mode, Conjugate, StorageOrder>
35 static void run(Index size,
const LhsScalar* _lhs, Index lhsStride, RhsScalar* rhs)
37 triangular_solve_vector<LhsScalar,RhsScalar,Index,
OnTheLeft,
40 >::run(size, _lhs, lhsStride, rhs);
45 template<
typename LhsScalar,
typename RhsScalar,
typename Index,
int Mode,
bool Conjugate>
46 struct triangular_solve_vector<LhsScalar, RhsScalar, Index,
OnTheLeft, Mode, Conjugate,
RowMajor>
51 static void run(Index size,
const LhsScalar* _lhs, Index lhsStride, RhsScalar* rhs)
53 typedef Map<const Matrix<LhsScalar,Dynamic,Dynamic,RowMajor>, 0, OuterStride<> > LhsMap;
54 const LhsMap lhs(_lhs,size,size,OuterStride<>(lhsStride));
55 typename internal::conditional<
57 const CwiseUnaryOp<typename internal::scalar_conjugate_op<LhsScalar>,LhsMap>,
61 for(Index pi=IsLower ? 0 : size;
62 IsLower ? pi<size : pi>0;
63 IsLower ? pi+=PanelWidth : pi-=PanelWidth)
65 Index actualPanelWidth = (std::min)(IsLower ? size - pi : pi, PanelWidth);
67 Index r = IsLower ? pi : size - pi;
73 Index startRow = IsLower ? pi : pi-actualPanelWidth;
74 Index startCol = IsLower ? 0 : pi;
78 &lhs.coeffRef(startRow,startCol), lhsStride,
84 for(Index k=0; k<actualPanelWidth; ++k)
86 Index i = IsLower ? pi+k : pi-k-1;
87 Index s = IsLower ? pi : i+1;
89 rhs[i] -= (cjLhs.row(i).segment(s,k).transpose().cwiseProduct(Map<
const Matrix<RhsScalar,Dynamic,1> >(rhs+s,k))).sum();
99 template<
typename LhsScalar,
typename RhsScalar,
typename Index,
int Mode,
bool Conjugate>
100 struct triangular_solve_vector<LhsScalar, RhsScalar, Index,
OnTheLeft, Mode, Conjugate,
ColMajor>
105 static void run(Index size,
const LhsScalar* _lhs, Index lhsStride, RhsScalar* rhs)
107 typedef Map<const Matrix<LhsScalar,Dynamic,Dynamic,ColMajor>, 0, OuterStride<> > LhsMap;
108 const LhsMap lhs(_lhs,size,size,OuterStride<>(lhsStride));
109 typename internal::conditional<Conjugate,
110 const CwiseUnaryOp<typename internal::scalar_conjugate_op<LhsScalar>,LhsMap>,
115 for(Index pi=IsLower ? 0 : size;
116 IsLower ? pi<size : pi>0;
117 IsLower ? pi+=PanelWidth : pi-=PanelWidth)
119 Index actualPanelWidth = (std::min)(IsLower ? size - pi : pi, PanelWidth);
120 Index startBlock = IsLower ? pi : pi-actualPanelWidth;
121 Index endBlock = IsLower ? pi + actualPanelWidth : 0;
123 for(Index k=0; k<actualPanelWidth; ++k)
125 Index i = IsLower ? pi+k : pi-k-1;
127 rhs[i] /= cjLhs.coeff(i,i);
129 Index r = actualPanelWidth - k - 1;
130 Index s = IsLower ? i+1 : i-r;
132 Map<Matrix<RhsScalar,Dynamic,1> >(rhs+s,r) -= rhs[i] * cjLhs.col(i).segment(s,r);
134 Index r = IsLower ? size - endBlock : startBlock;
142 &lhs.coeffRef(endBlock,startBlock), lhsStride,
144 rhs+endBlock, 1, RhsScalar(-1));
154 #endif // EIGEN_TRIANGULAR_SOLVER_VECTOR_H