lines 9-151 of file: include/cppad/utility/lu_invert.hpp

{xrst_begin LuInvert}
{xrst_spell
   jp
   permuted
}

Invert an LU Factored Equation
##############################

Syntax
******

   # ``include <cppad/utility/lu_invert.hpp>``

``LuInvert`` ( *ip* , *jp* , *LU* , *X* )

Description
***********
Solves the matrix equation *A* * *X* = *B*
using an LU factorization computed by :ref:`LuFactor-name` .

Include
*******
The file ``cppad/utility/lu_invert.hpp``
is included by ``cppad/cppad.hpp``
but it can also be included separately with out the rest of
the ``CppAD`` routines.

Matrix Storage
**************
All matrices are stored in row major order.
To be specific, if :math:`Y` is a vector
that contains a :math:`p` by :math:`q` matrix,
the size of :math:`Y` must be equal to :math:`p * q` and for
:math:`i = 0 , \ldots , p-1`,
:math:`j = 0 , \ldots , q-1`,

.. math::

   Y_{i,j} = Y[ i * q + j ]

ip
**
The argument *ip* has prototype

   ``const`` *SizeVector* & *ip*

(see description for *SizeVector* in
:ref:`LuFactor<LuFactor@SizeVector>` specifications).
The size of *ip* is referred to as *n* in the
specifications below.
The elements of *ip* determine
the order of the rows in the permuted matrix.

jp
**
The argument *jp* has prototype

   ``const`` *SizeVector* & *jp*

(see description for *SizeVector* in
:ref:`LuFactor<LuFactor@SizeVector>` specifications).
The size of *jp* must be equal to *n* .
The elements of *jp* determine
the order of the columns in the permuted matrix.

LU
**
The argument *LU* has the prototype

   ``const`` *FloatVector* & *LU*

and the size of *LU* must equal :math:`n * n`
(see description for *FloatVector* in
:ref:`LuFactor<LuFactor@FloatVector>` specifications).

L
=
We define the lower triangular matrix *L* in terms 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 :math:`i = 0 , \ldots , n-1` and :math:`j = 0 , \ldots , i`.

U
=
We define the upper triangular matrix *U* in terms 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 :math:`i = 0 , \ldots , n-2` and :math:`j = i+1 , \ldots , n-1`.

P
=
We define the permuted matrix *P* in terms of
the matrix *L* and the matrix *U*
by *P* = *L* * *U* .

A
=
The matrix *A* ,
which defines the linear equations that we are solving, is given by

   *P* ( *i* , *j* ) = *A* [ *ip* [ *i* ] * *n* + *jp* [ *j* ] ]

(Hence
*LU* contains a permuted factorization of the matrix *A* .)

X
*
The argument *X* has prototype

   *FloatVector* & *X*

(see description for *FloatVector* in
:ref:`LuFactor<LuFactor@FloatVector>` specifications).
The matrix *X*
must have the same number of rows as the matrix *A* .
The input value of *X* is the matrix *B* and the
output value solves the matrix equation *A* * *X* = *B* .
{xrst_toc_hidden
   example/utility/lu_invert.cpp
   xrst/lu_invert_hpp.xrst
}
Example
*******
The file :ref:`lu_solve.hpp-name` is a good example usage of
``LuFactor`` with ``LuInvert`` .
The file
:ref:`lu_invert.cpp-name`
contains an example and test of using ``LuInvert`` by itself.

Source
******
The file :ref:`lu_invert.hpp-name` contains the
current source code that implements these specifications.

{xrst_end LuInvert}