1 # ifndef CPPAD_UTILITY_LU_INVERT_HPP
2 # define CPPAD_UTILITY_LU_INVERT_HPP
168 template <
typename SizeVector,
typename FloatVector>
170 const SizeVector &ip,
171 const SizeVector &jp,
172 const FloatVector &LU,
181 CheckNumericType<Float>();
184 CheckSimpleVector<Float, FloatVector>();
185 CheckSimpleVector<size_t, SizeVector>();
189 size_t n = ip.size();
191 size_t(jp.size()) == n,
192 "Error in LuInvert: jp must have size equal to n * n"
195 size_t(LU.size()) == n * n,
196 "Error in LuInvert: Lu must have size equal to n * m"
198 size_t m = size_t(B.size()) / n;
200 size_t(B.size()) == n * m,
201 "Error in LuSolve: B must have size equal to a multiple of n"
208 for(k = 0; k < m; k++)
210 for(p = 0; p < n; p++)
212 etmp = B[ ip[p] * m + k ] / LU[ ip[p] * n + jp[p] ];
213 B[ ip[p] * m + k ] = etmp;
215 for(i = p+1; i < n; i++)
216 B[ ip[i] * m + k ] -=
217 etmp * LU[ ip[i] * n + jp[p] ];
224 etmp = B[ ip[p] * m + k ];
226 for(i = 0; i < p; i++ )
227 B[ ip[i] * m + k ] -=
228 etmp * LU[ ip[i] * n + jp[p] ];
232 for(i = 0; i < n; i++)
#define CPPAD_ASSERT_KNOWN(exp, msg)
Check that exp is true, if not print msg and terminate execution.
Define the CppAD error checking macros (all of which begin with CPPAD_ASSERT_)
void LuInvert(const SizeVector &ip, const SizeVector &jp, const FloatVector &LU, FloatVector &B)