Prev Next full_newton

@(@\newcommand{\B}[1]{{\bf #1}} \newcommand{\R}[1]{{\rm #1}}@)@
Execute Full Newton Steps For Control Constraint

Syntax
num_itr = control::full_newton(
xy_out, xy_in, L_fun, criteria, max_itr, L_y, work
)

Prototype 

template <class Scalar> size_t
full_newton(
     VECTOR(Scalar)&                                    xy_out        ,
     const VECTOR(Scalar)&                              xy_in         ,
     CppAD::ADFun<Scalar>&                              L_fun         ,
     Scalar                                             criteria      ,
     size_t                                             max_itr       ,
     const CPPAD_SPARSE(Scalar)&                        L_y           ,
     CppAD::sparse_jac_work&                            work          )

Purpose
Let @(@ ( x^0 , y^0 ) @)@ denote the value of @(@ (x, y) @)@ in xy_in . Let @(@ y^k @)@ denote the value of @(@ y @)@ after k Newton steps (not its k-th order Taylor coefficient). The k-th Newton step is @[@ y^k = y^{k-1} - L_y ( x^{k-1} , y^{k-1} )^{-1} L ( x^{k-1}, y^{k-1} ) @]@

xy_in
This vector has size n + m . Its first n components specify x . The other components specify the initial value for y . It is the value of @(@ ( x^0 , y^0 ) @)@ for the first Newton step.

xy_out
This vector has size n + m . Its first n components specify x . The other components specify the initial value for y . It is the value of @(@ ( x^k , y^k ) @)@ for the last Newton step.

L_fun
The operation sequence for the constraint function @(@ L : \B{R}^{2 J} \times \B{R}^{4 J} \rightarrow \B{R}^{4 J} @)@ is stored in L_fun .

criteria
This is the convergence criteria in terms of the Euclidean norm squared. Convergence is accepted when @(@ | L(x, y) |^2 @)@ is less than criteria .

max_itr
This is the maximum number of Newton steps to execute (which must be greater than zero). If criteria = 0 , xy_out will correspond to exactly max_itr Newton steps.

L_y
This argument is the value of @(@ L_y (x, y) @)@ at the value of @(@ (x, y) @)@ corresponding to xy_in .

work
This is a work vector used to reduce the work. It must correspond to L_y and not be empty; see jac_constraint .

num_itr
Is the number of Newton steps executed. If num_itr == max_itr , the convergence criteria may not be satisfied.

Example 

bool test_control_full_newton(void)
{    bool ok = true;
     double eps99 = 99.0 * std::numeric_limits<double>::epsilon();
     //
     // record L_fun
     CppAD::ADFun<double>       L_fun;
     size_t                     J = 10;
     double                     delta_t = 0.1;
     VECTOR(double) p(4);
     for(size_t i = 0; i < 4; i++)
          p[i] = 0.2;
     control::rec_constraint(L_fun, J, delta_t, p);
     size_t m      = L_fun.Range();
     size_t n      = L_fun.Domain() - m;
     //
     // value of x and initial y
     VECTOR(double) xy_in(n + m), xy_out;
     for(size_t i = 0; i < n+m; i++)
          xy_in[i] = 0.0;
     //
     double criteria  = double(n + m) * eps99;
     size_t max_itr   = 10;
     CPPAD_SPARSE(double)    L_y;
     CppAD::sparse_jac_work  work;
     jac_constraint(L_y, L_fun, xy_in, work);
     //
     size_t num_itr = control::full_newton(
          xy_out, xy_in, L_fun, criteria, max_itr, L_y, work
     );
     // max sure it did not require maximum nunber of iterations
     ok &= num_itr < max_itr;
     // make sure xy has not gone to a huge value
     ok &= norm_squared( xy_out ) / (6 * J) < 1.0;
     // check the convergence criteria
     VECTOR(double) L = L_fun.Forward(0, xy_out);
     ok &= norm_squared( L ) < criteria;
     //
     return ok;
}
Input File: src/control.hpp