| 39 | template <int SegSizeAtCompileTime> |
| 40 | template <typename BlockScalarVector, typename ScalarVector, typename IndexVector> |
| 41 | EIGEN_DONT_INLINE void LU_kernel_bmod<SegSizeAtCompileTime>::run(const Index segsize, BlockScalarVector& dense, ScalarVector& tempv, ScalarVector& lusup, Index& luptr, const Index lda, |
| 42 | const Index nrow, IndexVector& lsub, const Index lptr, const Index no_zeros) |
| 43 | { |
| 44 | typedef typename ScalarVector::Scalar Scalar; |
| 45 | // First, copy U[*,j] segment from dense(*) to tempv(*) |
| 46 | // The result of triangular solve is in tempv[*]; |
| 47 | // The result of matric-vector update is in dense[*] |
| 48 | Index isub = lptr + no_zeros; |
| 49 | Index i; |
| 50 | Index irow; |
| 51 | for (i = 0; i < ((SegSizeAtCompileTime==Dynamic)?segsize:SegSizeAtCompileTime); i++) |
| 52 | { |
| 53 | irow = lsub(isub); |
| 54 | tempv(i) = dense(irow); |
| 55 | ++isub; |
| 56 | } |
| 57 | // Dense triangular solve -- start effective triangle |
| 58 | luptr += lda * no_zeros + no_zeros; |
| 59 | // Form Eigen matrix and vector |
| 60 | Map<Matrix<Scalar,SegSizeAtCompileTime,SegSizeAtCompileTime, ColMajor>, 0, OuterStride<> > A( &(lusup.data()[luptr]), segsize, segsize, OuterStride<>(lda) ); |
| 61 | Map<Matrix<Scalar,SegSizeAtCompileTime,1> > u(tempv.data(), segsize); |
| 62 | |
| 63 | u = A.template triangularView<UnitLower>().solve(u); |
| 64 | |
| 65 | // Dense matrix-vector product y <-- B*x |
| 66 | luptr += segsize; |
| 67 | const Index PacketSize = internal::packet_traits<Scalar>::size; |
| 68 | Index ldl = internal::first_multiple(nrow, PacketSize); |
| 69 | Map<Matrix<Scalar,Dynamic,SegSizeAtCompileTime, ColMajor>, 0, OuterStride<> > B( &(lusup.data()[luptr]), nrow, segsize, OuterStride<>(lda) ); |
| 70 | Index aligned_offset = internal::first_default_aligned(tempv.data()+segsize, PacketSize); |
| 71 | Index aligned_with_B_offset = (PacketSize-internal::first_default_aligned(B.data(), PacketSize))%PacketSize; |
| 72 | Map<Matrix<Scalar,Dynamic,1>, 0, OuterStride<> > l(tempv.data()+segsize+aligned_offset+aligned_with_B_offset, nrow, OuterStride<>(ldl) ); |
| 73 | |
| 74 | l.noalias() = B * u; |
| 75 | |
| 76 | // Scatter tempv[] into SPA dense[] as a temporary storage |
| 77 | isub = lptr + no_zeros; |
| 78 | for (i = 0; i < ((SegSizeAtCompileTime==Dynamic)?segsize:SegSizeAtCompileTime); i++) |
| 79 | { |
| 80 | irow = lsub(isub++); |
| 81 | dense(irow) = tempv(i); |
| 82 | } |
| 83 | |
| 84 | // Scatter l into SPA dense[] |
| 85 | for (i = 0; i < nrow; i++) |
| 86 | { |
| 87 | irow = lsub(isub++); |
| 88 | dense(irow) -= l(i); |
| 89 | } |
| 90 | } |
| 91 | |
| 92 | template <> struct LU_kernel_bmod<1> |
| 93 | { |
nothing calls this directly
no test coverage detected