lines 910-1341 of file: introduction/exp_eps.xrst

{xrst_begin exp_eps_rev2}
{xrst_spell
   preform
   rcll
}

exp_eps: Second Order Reverse Sweep
###################################

Purpose
*******
In general, a second order reverse sweep is given the
:ref:`exp_eps_for1@First Order Expansion`
for all of the variables in an operation sequence.
Given a choice of a particular variable,
it computes the derivative,
of that variables first order expansion coefficient,
with respect to all of the independent variables.

Mathematical Form
*****************
Suppose that we use the algorithm :ref:`exp_eps.hpp-name`
to compute ``exp_eps`` ( *x* , *epsilon* )
with *x* is equal to .5
and *epsilon* is equal to .2.
For this case, the mathematical function for the operation sequence
corresponding to ``exp_eps`` is

.. math::

   f ( x , \varepsilon ) =   1 + x + x^2 / 2

The corresponding derivative of the
partial derivative with respect to :math:`x` is

.. math::
   :nowrap:

   \begin{eqnarray}
   \Dpow{2}{x} f ( x , \varepsilon ) & = &  1
   \\
   \partial_\varepsilon \partial_x f ( x , \varepsilon ) & = &  0
   \end{eqnarray}

epsilon
*******
Since :math:`\varepsilon` is an independent variable,
it could included as an argument to all of the
:math:`f_j` functions below.
The result would be that all the partials with respect to
:math:`\varepsilon` would be zero and hence we drop it to simplify
the presentation.

f_7
***
In reverse mode we choose one dependent variable and
compute its derivative with respect to all the independent variables.
For our example, we chose the value returned by :ref:`exp_eps.hpp-name`
which is :math:`v_7`.
We begin with the function :math:`f_7` where :math:`v_7`
is both an argument and the value of the function; i.e.,

.. math::
   :nowrap:

   \begin{eqnarray}
   f_7 \left(
      v_1^{(0)} , v_1^{(1)} , \ldots , v_7^{(0)} , v_7^{(1)}
   \right)
   & = & v_7^{(1)}
   \\
   \D{f_7}{v_7^{(1)}} & = & 1
   \end{eqnarray}

All the other partial derivatives of :math:`f_7` are zero.

Index 7: f_6
************
The last operation has index 7,

.. math::
   :nowrap:

   \begin{eqnarray}
      v_7^{(0)} & = &   v_4^{(0)} + v_6^{(0)}
      \\
      v_7^{(1)} & = &   v_4^{(1)} + v_6^{(1)}
   \end{eqnarray}

We define the function
:math:`f_6 \left( v_1^{(0)} , \ldots , v_6^{(1)} \right)`
as equal to :math:`f_7`
except that :math:`v_7^{(0)}` and :math:`v_7^{(1)}` are eliminated using
this operation; i.e.

.. math::

   f_6  =
   f_7 \left[ v_1^{(0)} , \ldots , v_6^{(1)} ,
      v_7^{(0)} \left( v_4^{(0)} , v_6^{(0)} \right)  ,
      v_7^{(1)} \left( v_4^{(1)} , v_6^{(1)} \right)
   \right]

It follows that

.. math::

   \begin{array}{rcll}
   \D{f_6}{v_4^{(1)}}
   & = & \D{f_7}{v_4^{(1)}} +
      \D{f_7}{v_7^{(1)}} * \D{v_7^{(1)}}{v_4^{(1)}}
   & = 1
   \\
   \D{f_6}{v_6^{(1)}}
   & = & \D{f_7}{v_6^{(1)}} +
      \D{f_7}{v_7^{(1)}} * \D{v_7^{(1)}}{v_6^{(1)}}
   & = 1
   \end{array}

All the other partial derivatives of :math:`f_6` are zero.

Index 6: f_5
************
The previous operation has index 6,

.. math::
   :nowrap:

   \begin{eqnarray}
      v_6^{(0)} & = & v_5^{(0)} / 2
      \\
      v_6^{(1)} & = & v_5^{(1)} / 2
   \end{eqnarray}

We define the function
:math:`f_5 \left( v_1^{(0)} , \ldots , v_5^{(1)} \right)`
as equal to :math:`f_6`
except that :math:`v_6^{(0)}` and :math:`v_6^{(1)}` are eliminated using
this operation; i.e.

.. math::

   f_5 =
   f_6 \left[ v_1^{(0)} , \ldots , v_5^{(1)} ,
      v_6^{(0)} \left( v_5^{(0)} \right) ,
      v_6^{(1)} \left( v_5^{(1)} \right)
   \right]

It follows that

.. math::

   \begin{array}{rcll}
   \D{f_5}{v_4^{(1)}}
   & = & \D{f_6}{v_4^{(1)}}
   & = 1
   \\
   \D{f_5}{v_5^{(1)}}
   & = & \D{f_6}{v_5} +
      \D{f_6}{v_6^{(1)}} * \D{v_6^{(1)}}{v_5^{(1)}}
   & = 0.5
   \end{array}

All the other partial derivatives of :math:`f_5` are zero.

Index 5: f_4
************
The previous operation has index 5,

.. math::
   :nowrap:

   \begin{eqnarray}
      v_5^{(0)} & = & v_3^{(0)} * v_1^{(0)}
      \\
      v_5^{(1)} & = & v_3^{(1)} * v_1^{(0)} + v_3^{(0)} * v_1^{(1)}
   \end{eqnarray}

We define the function
:math:`f_4 \left( v_1^{(0)} , \ldots , v_4^{(1)} \right)`
as equal to :math:`f_5`
except that :math:`v_5^{(0)}` and :math:`v_5^{(1)}` are eliminated using
this operation; i.e.

.. math::

   f_4 =
   f_5 \left[  v_1^{(0)} , \ldots , v_4^{(1)} ,
      v_5^{(0)} \left( v_1^{(0)}, v_3^{(0)} \right) ,
      v_5^{(1)} \left( v_1^{(0)}, v_1^{(1)}, v_3^{(0)} , v_3^{(1)} \right) ,
   \right]

Given the information from the forward sweep, we have
:math:`v_1^{(0)} =  0.5`,
:math:`v_3^{(0)} =  0.5`,
:math:`v_1^{(1)} =  1`,
:math:`v_3^{(1)} =  1`,
and the fact that the partial of :math:`f_5` with respect to
:math:`v_5^{(0)}` is zero, we have

.. math::

   \begin{array}{rcll}
   \D{f_4}{v_1^{(0)}}
   & = & \D{f_5}{v_1^{(0)}}
     +   \D{f_5}{v_5^{(1)}} * \D{v_5^{(1)}}{v_1^{(0)}}
   & = 0.5
   \\
   \D{f_4}{v_1^{(1)}}
   & = & \D{f_5}{v_1^{(1)}}
     +   \D{f_5}{v_5^{(1)}} * \D{v_5^{(1)}}{v_1^{(1)}}
   & = 0.25
   \\
   \D{f_4}{v_3^{(0)}}
   & = & \D{f_5}{v_3^{(0)}}
     +   \D{f_5}{v_5^{(1)}} * \D{v_5^{(1)}}{v_3^{(0)}}
   & = 0.5
   \\
   \D{f_4}{v_3^{(1)}}
   & = & \D{f_3}{v_1^{(1)}}
     +   \D{f_5}{v_5^{(1)}} * \D{v_5^{(1)}}{v_3^{(1)}}
   & = 0.25
   \\
   \D{f_4}{v_4^{(1)}}
   & = & \D{f_5}{v_4^{(1)}}
   & = 1
   \end{array}

All the other partial derivatives of :math:`f_5` are zero.

Index 4: f_3
************
The previous operation has index 4,

.. math::
   :nowrap:

   \begin{eqnarray}
      v_4^{(0)} = 1 + v_3^{(0)}
      \\
      v_4^{(1)} = v_3^{(1)}
   \end{eqnarray}

We define the function
:math:`f_3 \left( v_1^{(0)} , \ldots , v_3^{(1)} \right)`
as equal to :math:`f_4`
except that :math:`v_4^{(0)}` and :math:`v_4^{(1)}` are eliminated using
this operation; i.e.

.. math::

   f_3 =
   f_4 \left[ v_1^{(0)} , \ldots , v_3^{(1)} ,
      v_4^{(0)} \left( v_3^{(0)} \right) ,
      v_4^{(1)} \left( v_3^{(1)} \right)
   \right]

It follows that

.. math::

   \begin{array}{rcll}
   \D{f_3}{v_1^{(0)}}
   & = & \D{f_4}{v_1^{(0)}}
   & =  0.5
   \\
   \D{f_3}{v_1^{(1)}}
   & = & \D{f_4}{v_1^{(1)}}
   & =  0.25
   \\
   \D{f_3}{v_2^{(0)}}
   & = & \D{f_4}{v_2^{(0)}}
   & = 0
   \\
   \D{f_3}{v_2^{(1)}}
   & = & \D{f_4}{v_2^{(1)}}
   & = 0
   \\
   \D{f_3}{v_3^{(0)}}
   & = & \D{f_4}{v_3^{(0)}}
     +   \D{f_4}{v_4^{(0)}} * \D{v_4^{(0)}}{v_3^{(0)}}
   & = 0.5
   \\
   \D{f_3}{v_3^{(1)}}
   & = & \D{f_4}{v_3^{(1)}}
     +   \D{f_4}{v_4^{(1)}} * \D{v_4^{(1)}}{v_3^{(1)}}
   & = 1.25
   \end{array}

Index 3: f_2
************
The previous operation has index 3,

.. math::
   :nowrap:

   \begin{eqnarray}
      v_3^{(0)} & = & v_2^{(0)} / 1
      \\
      v_3^{(1)} & = & v_2^{(1)} / 1
   \end{eqnarray}

We define the function
:math:`f_2 \left( v_1^{(0)} , \ldots , v_2^{(1)} \right)`
as equal to :math:`f_3`
except that :math:`v_3^{(0)}` and :math:`v_3^{(1)}` are eliminated using
this operation; i.e.

.. math::

   f_2 =
   f_3 \left[ v_1^{(0)} , \ldots , v_2^{(1)} ,
      v_3^{(0)} \left( v_2^{(0)} \right) ,
      v_3^{(1)} \left( v_2^{(1)} \right)
   \right]

It follows that

.. math::

   \begin{array}{rcll}
   \D{f_2}{v_1^{(0)}}
   & = & \D{f_3}{v_1^{(0)}}
   & =  0.5
   \\
   \D{f_2}{v_1^{(1)}}
   & = & \D{f_3}{v_1^{(1)}}
   & =  0.25
   \\
   \D{f_2}{v_2^{(0)}}
   & = & \D{f_3}{v_2^{(0)}}
     +   \D{f_3}{v_3^{(0)}} * \D{v_3^{(0)}}{v_2^{(0)}}
   & = 0.5
   \\
   \D{f_2}{v_2^{(1)}}
   & = & \D{f_3}{v_2^{(1)}}
     +   \D{f_3}{v_3^{(1)}} * \D{v_3^{(1)}}{v_2^{(0)}}
   & = 1.25
   \end{array}

Index 2: f_1
************
The previous operation has index 1,

.. math::
   :nowrap:

   \begin{eqnarray}
      v_2^{(0)} & = & 1 * v_1^{(0)}
      \\
      v_2^{(1)} & = & 1 * v_1^{(1)}
   \end{eqnarray}

We define the function
:math:`f_1 \left( v_1^{(0)} , v_1^{(1)} \right)`
as equal to :math:`f_2`
except that :math:`v_2^{(0)}` and :math:`v_2^{(1)}` are eliminated using
this operation; i.e.

.. math::

   f_1 =
   f_2 \left[  v_1^{(0)} , v_1^{(1)} ,
      v_2^{(0)} \left( v_1^{(0)} \right)  ,
      v_2^{(1)} \left( v_1^{(1)} \right)
   \right]

It follows that

.. math::

   \begin{array}{rcll}
   \D{f_1}{v_1^{(0)}}
   & = & \D{f_2}{v_1^{(0)}}
     +   \D{f_2}{v_2^{(0)}} * \D{v_2^{(0)}}{v_1^{(0)}}
   & =  1
   \\
   \D{f_1}{v_1^{(1)}}
   & = & \D{f_2}{v_1^{(1)}}
     +   \D{f_2}{v_2^{(1)}} * \D{v_2^{(1)}}{v_1^{(1)}}
   & = 1.5
   \end{array}

Note that :math:`v_1` is equal to :math:`x`,
so the second partial derivative of
``exp_eps`` ( *x* , *epsilon* )
at *x* equal to .5 and *epsilon* equal .2 is

.. math::

   \Dpow{2}{x} v_7^{(0)}
   = \D{v_7^{(1)}}{x}
   = \D{f_1}{v_1^{(0)}}
   = 1

There is a theorem about algorithmic differentiation that explains why
the other partial of :math:`f_1` is equal to the first partial of
``exp_eps`` ( *x* , *epsilon* )
with respect to :math:`x`.
{xrst_toc_hidden
   introduction/exp_eps_rev2.cpp
}

Verification
************
The file :ref:`exp_eps_rev2.cpp-name` contains a routine
that verifies the values computed above.
It only tests the partial derivatives of
:math:`f_j` that might not be equal to the corresponding
partials of :math:`f_{j+1}`; i.e., the
other partials of :math:`f_j` must be equal to the corresponding
partials of :math:`f_{j+1}`.

Exercises
*********

#. Consider the case where :math:`x = .1`
   and we first preform a zero order forward mode sweep
   for the operation sequence used above (in reverse order).
   What are the results of a
   first order reverse mode sweep; i.e.,
   what are the corresponding values for
   :math:`\D{f_j}{v_k}` for all :math:`j, k` such that
   :math:`\D{f_j}{v_k} \neq 0`.
#. Create a modified version of
   :ref:`exp_eps_rev2.cpp-name`
   that verifies the values you obtained for the previous exercise.
   Also create and run a main program that reports the result
   of calling the modified version of
   :ref:`exp_eps_rev2.cpp-name` .

{xrst_end exp_eps_rev2}
