Prev Next multi_atomic_three_user

@(@\newcommand{\W}[1]{ \; #1 \; } \newcommand{\R}[1]{ {\rm #1} } \newcommand{\B}[1]{ {\bf #1} } \newcommand{\D}[2]{ \frac{\partial #1}{\partial #2} } \newcommand{\DD}[3]{ \frac{\partial^2 #1}{\partial #2 \partial #3} } \newcommand{\Dpow}[2]{ \frac{\partial^{#1}}{\partial {#2}^{#1}} } \newcommand{\dpow}[2]{ \frac{ {\rm d}^{#1}}{{\rm d}\, {#2}^{#1}} }@)@This is cppad-20221105 documentation. Here is a link to its current documentation .
Defines a atomic_three Operation that Computes Square Root

Syntax
atomic_user a_square_root
 a_square_root(au, ay)

Purpose
This atomic function operation computes a square root using Newton's method. It is meant to be very inefficient in order to demonstrate timing results.

au
This argument has prototype
    const ADvector& au
 where ADvector is a simple vector class with elements of type AD<double>. The size of au is three.

num_itr
We use the notation
    num_itr = size_t( Integer( au[0] ) )
 for the number of Newton iterations in the computation of the square root function. The component au[0] must be a parameter .

y_initial
We use the notation
    y_initial = au[1]
 for the initial value of the Newton iterate.

y_squared
We use the notation
    y_squared = au[2]
 for the value we are taking the square root of.

ay
This argument has prototype
    ADvector& ay
 The size of ay is one and ay[0] is the square root of y_squared .

Limitations
Only zero order forward mode is implements for the atomic_user class.

Source 

// includes used by all source code in multi_atomic_three.cpp file
# include <cppad/cppad.hpp>
# include "multi_atomic_three.hpp"
# include "team_thread.hpp"
//
namespace {
using CppAD::thread_alloc;                // multi-threading memory allocator
using CppAD::vector;                      // uses thread_alloc
typedef CppAD::ad_type_enum ad_type_enum; // constant, dynamic or variable

class atomic_user : public CppAD::atomic_three<double> {
public:
    // ctor
    atomic_user(void)
    : CppAD::atomic_three<double>("atomic_square_root")
    { }
private:
    // for_type
    bool for_type(
        const vector<double>&        parameter_u ,
        const vector<ad_type_enum>&  type_u      ,
        vector<ad_type_enum>&        type_y      ) override
    {   bool ok = parameter_u.size() == 3;
        ok     &= type_u.size() == 3;
        ok     &= type_y.size() == 1;
        if( ! ok )
            return false;
        ok     &= type_u[0] < CppAD::variable_enum;
        if( ! ok )
            return false;
        type_y[0] = std::max( type_u[0], type_u[1] );
        type_y[0] = std::max( type_y[0], type_u[2] );
        //
        return true;
    }
    // forward
    bool forward(
        const vector<double>&        parameter_u ,
        const vector<ad_type_enum>&  type_u      ,
        size_t                       need_y      ,
        size_t                       order_low   ,
        size_t                       order_up    ,
        const vector<double>&        taylor_u    ,
        vector<double>&              taylor_y    ) override
    {
# ifndef NDEBUG
        size_t n = taylor_u.size() / (order_up + 1);
        size_t m = taylor_y.size() / (order_up + 1);
        assert( n == 3 );
        assert( m == 1 );
# endif
        // only implementing zero order forward for this example
        if( order_up != 0 )
            return false;

        // extract components of argument vector
        size_t num_itr    = size_t( taylor_u[0] );
        double y_initial  = taylor_u[1];
        double y_squared  = taylor_u[2];

        // Use Newton's method to solve f(y) = y^2 = y_squared
        double y_itr = y_initial;
        for(size_t itr = 0; itr < num_itr; itr++)
        {   // solve (y - y_itr) * f'(y_itr) = y_squared - y_itr^2
            double fp_itr = 2.0 * y_itr;
            y_itr         = y_itr + (y_squared - y_itr * y_itr) / fp_itr;
        }

        // return the Newton approximation for f(y) = y_squared
        taylor_y[0] = y_itr;
        return true;
    }
};
}
Input File: example/multi_thread/multi_atomic_three.cpp