Newton Step Method for Derivatives of Implicit Functions

@(@\newcommand{\B}[1]{{\bf #1}} \newcommand{\R}[1]{{\rm #1}}@)@Newton Step Method for Derivatives of Implicit Functions
Syntax

implicit_newton newton_ad(

     full_step, num_step, aL_fun, F_fun, solve

)

rk = newton_ad.Forward(k, xk)

dw = newton_ad.Reverse(q, w)

Purpose
Given a function @(@ L : \B{R}^{n \times m} \rightarrow \B{R}^m @)@, we define the implicit function @(@ Y : \B{R}^n \rightarrow \B{R}^m @)@ by @(@ L[ x , Y(x) ] = 0 @)@. The partial @(@ L_y [ x , Y(x) ] @)@ is assumed to be invertible for all @(@ x @)@. We define the reduced function @(@ R(x) = F[ x , Y(x) ] @)@ where @(@ F : \B{R}^{n \times m} \rightarrow \B{R}^p @)@. The object newton_ad can be used to compute derivatives of the reduced function @(@ R(x) @)@. Note that in the special case where @(@ F(x, y) = y @)@, @(@ R(x) = Y(x) @)@.

full_step
This argument has prototype



     bool full_step

If it is true, full Newton steps are used. If it is false, partial Newton steps are used.

num_step
This argument has prototype



     size_t num_step

It is the number of Newton steps in the method.

aL_fun
This argument has prototype



     CppAD::ADFun< CppAD::AD<double> >& aL_fun

and is the CppAD function object corresponding to @(@ L(x, y) @)@. Note that a CppAD::ADFun<double> object could be used for the partial Newton step method. Note that this function is not const because its forward mode is used for calculations. Upon return, there are no forward mode coefficient left in aL_fun .

F_fun
This argument has prototype



     const CppAD::ADFun<double>& F_fun

and is the CppAD function object corresponding to @(@ F(x, y) @)@.

solve
This argument has prototype



     const Solve& solve

The type Solver must support the default constructor and the assignment operator. It must also support the following operations:

solve.function
In the syntax



     y = solve.function(x)

the argument x and the return value y have prototypes



     const VECTOR(double)&  x

     VECTOR(double)         y

where x.size() == n and y.size() == m . The return value satisfies the relation @(@ L(x, y) = 0 @)@.

solve.derivative
In the syntax



     b = solve.derivative(x, y)

the arguments have prototypes



     const VECTOR(Scalar)&  x

     const VECTOR(Scalar)&  y

The return value has prototype



     VECTOR(Scalar) b

This returns the value of @(@ L_y (x, y) @)@ for subsequent calls to solve.linear . Only the elements of @(@ L_y (x, y) @)@ that depend on @(@ (x, y) @)@ need be included in the vector @(@ b @)@. The type Scalar is either double or CppAD::AD<double>.

solve.linear
In the syntax



     u = solve.linear(b, v)

the arguments b , v and the return value u have prototypes



     const VECTOR(CppAD::AD<double>)& b

     const VECTOR(CppAD::AD<double>)& v

     const VECTOR(CppAD::AD<double>)  u

where both vectors have size m . The argument b has prototype



     VECTOR(CppAD::AD<double>)& b

The return value satisfies the equation @[@ u = L_y (x, y)^{-1} v @]@ where @(@ L_y (x, y) @)@ corresponds to b .

k
This is the order for this forward mode calculation and has prototype



     size_t k

Note that computing the k-th order, uses the results of the previous lower order calculations which are internally stored in the function objects.

xk
The argument xk has prototypes



     const VECTOR(double)& xk

its size is m and it is the k-th order Taylor coefficient for x . The results of the calculation are stored internally and used during higher order forward mode calculations. If k > 0 , the previous call to forward must have for order k - 1 or greater.

rk
The return value rk has prototype



     VECTOR(double) rk

its size is m and it is the k-th order Taylor coefficient for R(x) .

q
This argument has prototype



     size_t q

and is the number of Taylor coefficients in the function that we are differentiating.

w
This vector has prototype



     const VECTOR(double)& w

and its size is p*q . For @(@ i = 0 , \ldots , p-1 @)@ and @(@ k = 0 , \ldots , q-1 @)@, let @(@ r_i^k @)@ denote the k-th Taylor coefficients for @(@ R_i (x) @)@ The function we are differentiating is @[@ W(r) = \sum_{i=0}^{p-1} \sum_{k=0}^{q-1} w_{i * q + k} r_i^k @]@

dw
This return value has prototype



     VECTOR(double) dw

and its size is n*q . For @(@ i = 0 , \ldots , n-1 @)@ and @(@ k = 0 , \ldots , q-1 @)@, let @(@ X_i^k @)@ denote the k-th Taylor coefficients for @(@ x @)@. The value dw[i*q+k] is the partial of @(@ W[r(X)]% @)@ with respect to @(@ X_i^k @)@.

Example
The routine test_circle_implicit_newton is a simple example and test using this class. The routine test_control_reduced_objective is a control problem example and test using this class.

Input File: src/implicit_newton.hpp