CppAD::Independent on smart_pointers or reference_wrappers

[pnik073]

Question

Is there a recommended way to calculate derivatives, while the AD Independent variables are owned by multiple numerical objects and not from just one AD vector?

Problem

Currently, the CppAD::Independent accepts only AD<> variables (herein ax), given that the internal use of ADTape<Base>::Independent accesses directly ax[i].member_data without encapsulated methods. That restricts us from recording operations on vectors of reference_wrapper, smart pointers, etc.

Example

I will try to demonstrate an example of 2 "numerical" objects Object1 and Object2 instantiated in function ad_from_references.

Ideally, I would expect

• to own the AD<> instances on the calculation objects
• consolidate all AD<> Independent variables in a vector<reference_wrapper<AD<double>>> ax_r or similar memory management container
• start the recording Independent(ax_r)
• calculated from the Object1 and Object2 using directly the AD<> data
• get the derivatives from ADFun

Instead, at the moment I achieve the that by

• using a vector<AD<double>> ax
• start the recording Independent(ax)
• passing reference_wrapper<AD<double>> to the calculation objects via a sort of index tracking mechanism
• the calculation objects operate on reference_wrapper.get()

I may be missing something here, but the last 2 steps would imply that the consumer libraries need to create quite some infrastructure to fully adopt CppAD.

#include <gtest/gtest.h>
#include <cppad/cppad.hpp>
#include <vector>
#include <functional>
#include <algorithm>
#include <iterator>
#include <iostream>

using AdType = CppAD::AD<double>;
using AdVector = std::vector<AdType>;
using ReferenceAdVector = std::vector<std::reference_wrapper<AdType>>;

template<class ValueType>
ValueType func_on_ad(std::vector<ValueType> const& x)
{
    size_t n = x.size();
    ValueType res = 1.;
    for(size_t i = 0; i < n; i++)
    {
        res += CppAD::pow( x.at(i), i+1) ;
    }
    return res;
}

template<class ValueType>
ValueType func_on_wrapped_ad(std::vector<std::reference_wrapper<ValueType>> const& x)
{
    size_t n = x.size();
    ValueType res = 1.;
    for(size_t i = 0; i < n; i++)
    {
        res += CppAD::pow( x.at(i).get(), i+1) ;
    }
    return res;
}

struct Object
{
    Object(AdVector ax):ax_(std::move(ax)){}
    AdVector ax_;
    ReferenceAdVector ax_r_;
    /// A method to avoid
    void set_ax_ref(AdVector & ax, size_t start, size_t end)
    {
        std::transform(
            ax.begin()+start,
            ax.begin()+end+1,
            std::back_inserter(ax_r_),
            [](auto& x) {return std::ref(x);}
            );
    }
};
struct Object1 : public Object
{
    Object1(AdVector ax):Object(std::move(ax)){}
    AdType calculate() {return func_on_ad(ax_);}
    AdType calculate_from_refs() {return func_on_wrapped_ad(ax_r_);}
};
struct Object2: public Object
{
    Object2(AdVector ax):Object(std::move(ax)){}
    AdType calculate() {return CppAD::pow(ax_.back(),3);}
    AdType calculate_from_refs() {return CppAD::pow(ax_r_.back().get(),3);}
};

bool ad_from_references()
{
    AdVector ay(1);

    // Assume that the ax Independent values are initially owned by other numerical objects
    AdVector ax_0_1 {1,2};
    AdVector ax_2 {3};
    Object1 obj1(std::move(ax_0_1));
    Object2 obj2(std::move(ax_2));

    // Create a recordable vector
    AdVector ax;
    // copy or move equivalent for concat on AD<>
    // herein we get new instances of AD<>
    std::move(obj1.ax_.begin(),obj1.ax_.end(), std::back_inserter(ax));
    std::move(obj2.ax_.begin(),obj2.ax_.end(), std::back_inserter(ax));

    // However, it would be more flexible to (also) collect
    // The AD<> from their owners and record do CppAD::Independent(ax_r)
    // However, a vector of references/pointers etc. is not recordable via Independent
    ReferenceAdVector ax_r {
        std::ref(ax[0]),
        std::ref(ax[1]),
        std::ref(ax[2]),
    }; // So, this vector is actually of no use at the moment

    // Instead we need to add the references on the numerical objects
    // and calculate with these references.
    // That forces all consumer libs to create unnecessary mapping mechanisms
    // to access the ax elements.
    obj1.set_ax_ref(ax,0,1); // Range mapping required
    obj2.set_ax_ref(ax,2,2); // Range mapping required

    // Initiate recording
    // Herein the ADTape<Base>::Independent
    // accesses ax[i].members without encapsulated methods
    // such that it is not trivial to use smart pointers/reference wrappers of AD<>
    CppAD::Independent(ax);// Start recording

    // Due to the above line, the recording on the AD<>  of the Objects is not recordable,
    // the `AD<>` of `ax` are different instances 

    // >>> THE GOAL: ideally the consumer libs should operate
    // as follows, if we can record CppAD::Independent(ax_r)
    // AdType ay_ = obj1.calculate();
    // ay_ += obj2.calculate();

    // Instead of using func_on_ad we use calculate_from_refs
    // auto ay_ = func_on_ad(ax);
    AdType ay_ = obj1.calculate_from_refs();
    ay_ += obj2.calculate_from_refs();
    ay[0] = ay_;

    CppAD::ADFun<double> f(ax, ay);

    std::vector<double> x_derivative{1.,2.,3.};
    auto jac_final = f.Jacobian(x_derivative);

    std::cout << "AD f(x): " << Value(Var2Par(ay[0])) << std::endl;
    for ( auto i : jac_final ) {
        std::cout << "AD f'(x) " << i << " " << std::endl;
    }

    return true;
}

TEST(AdReferencePointer_Tests, RefWrapperTest) {
    EXPECT_TRUE(ad_from_references());
}

[bradbell]

If your wrapped vector type support all the operations required by SimpleVector ?
https://cppad.readthedocs.io/latest/SimpleVector.html
you shoud be able to use it for your AD vector type.