Prev | Next |
Syntax
# include <cppad/utility/lu_factor.hpp>
sign = LuFactor(ip, jp, LU)
A
where
A
is a square matrix.
cppad/utility/lu_factor.hpp
is included by cppad/cppad.hpp
but it can also be included separately with out the rest of
the CppAD
routines.
sign
has prototype
int sign
If
A
is invertible,
sign
is plus or minus one
and is the sign of the permutation corresponding to the row ordering
ip
and column ordering
jp
.
If
A
is not invertible,
sign
is zero.
ip
has prototype
SizeVector &ip
(see description of SizeVector
below).
The size of
ip
is referred to as
n
in the
specifications below.
The input value of the elements of
ip
does not matter.
The output value of the elements of
ip
determine
the order of the rows in the permuted matrix.
jp
has prototype
SizeVector &jp
(see description of SizeVector
below).
The size of
jp
must be equal to
n
.
The input value of the elements of
jp
does not matter.
The output value of the elements of
jp
determine
the order of the columns in the permuted matrix.
LU
has the prototype
FloatVector &LU
and the size of
LU
must equal @(@
n * n
@)@
(see description of FloatVector
below).
A
as the matrix corresponding to the input
value of
LU
.
P
in terms of
A
by
P(i, j) = A[ ip[i] * n + jp[j] ]
L
in terms of the
output value of
LU
.
The matrix
L
is zero above the diagonal
and the rest of the elements are defined by
L(i, j) = LU[ ip[i] * n + jp[j] ]
for @(@
i = 0 , \ldots , n-1
@)@ and @(@
j = 0 , \ldots , i
@)@.
U
in terms of the
output value of
LU
.
The matrix
U
is zero below the diagonal,
one on the diagonal,
and the rest of the elements are defined by
U(i, j) = LU[ ip[i] * n + jp[j] ]
for @(@
i = 0 , \ldots , n-2
@)@ and @(@
j = i+1 , \ldots , n-1
@)@.
sign
is non-zero,
L * U = P
If the return value of
sign
is zero,
the contents of
L
and
U
are not defined.
sign
is zero,
the determinant of
A
is zero.
If
sign
is non-zero,
using the output value of
LU
the determinant of the matrix
A
is equal to
sign * LU[ip[0], jp[0]] * ... * LU[ip[n-1], jp[n-1]]
SizeVector
must be a SimpleVector
class with
elements of type size_t
.
The routine CheckSimpleVector
will generate an error message
if this is not the case.
FloatVector
must be a
simple vector class
.
The routine CheckSimpleVector
will generate an error message
if this is not the case.
FloatVector
.
The type
Float
must satisfy the conditions
for a NumericType
type.
The routine CheckNumericType
will generate an error message
if this is not the case.
In addition, the following operations must be defined for any pair
of
Float
objects
x
and
y
:
Operation | Description |
log(x)
|
returns the logarithm of
x
as a
Float
object
|
lu_factor.hpp
defines the template function
template <class Float>
bool AbsGeq<Float>(const Float &x, const Float &y)
in the CppAD
namespace.
This function returns true if the absolute value of
x
is greater than or equal the absolute value of
y
.
It is used by LuFactor
to choose the pivot elements.
This template function definition uses the operator
<=
to obtain the absolute value for
Float
objects.
If this operator is not defined for your use of
Float
,
you will need to specialize this template so that it works for your
use of LuFactor
.
Complex numbers do not have the operation <=
defined.
The specializations
bool AbsGeq< std::complex<float> >
(const std::complex<float> &x, const std::complex<float> &y)
bool AbsGeq< std::complex<double> >
(const std::complex<double> &x, const std::complex<double> &y)
are define by including lu_factor.hpp
These return true if the sum of the square of the real and imaginary parts
of
x
is greater than or equal the
sum of the square of the real and imaginary parts of
y
.
LuFactor
by itself.
The file lu_solve.hpp
provides a useful example usage of
LuFactor
with LuInvert
.