Bug: passing around C-style array parameters crashes Clad

[guitargeek]

Using array parameters that are passed to the differentiated function as function arguments crashes Clad.

This is very important to support for the RooFit usecase, where we want to forward arrays of constant values like this, also passing arrays with offsets like obs + 2.

Reproducer:

// Compile with clang++ -std=c++11 -fplugin=/usr/lib/clad.so cladIssue.cxx -o cladIssue

#include "clad/Differentiator/Differentiator.h"

#include <iostream>

inline double foo(double *params, double *obs)
{
   return params[0] * obs[0];
}

double wrapper(double* params, double* obs)
{
   // return params[0] * obs[0]; // works!

   return foo(params, obs); // doesn't work!
}

int main()
{
   auto grad = clad::gradient(wrapper, "params");

   double param = 2.0;
   double obs = 2.0;
   double result = 0.0;
   grad.execute(&param, &obs, &result);

   std::cout << result << std::endl;
}

[vgvassilev]

That is probably related to #735.

[vaithak]

This issue is related to the initialization of derivatives of pointer arguments.
In the above example, the differentiation is w.r.t params only, so the user won't send clad::array_ref<double> _d_obs to us.
How do we initialize _d_obs? We don't have any size information about it.

The error is removed when we also include obs in the independent variables and pass _d_obs as an argument in the execute method.

   auto grad = clad::gradient(wrapper, "params, obs");

   double param = 2.0;
   double obs = 2.0;
   double _d_params = 0.0;
   double _d_obs = 0.0;
   grad.execute(&param, &obs, &_d_params, &_d_obs);

[guitargeek]

Hi @vaithak, thanks for the reply! That would not solve the problem in my use case I think, because doesn't differentiating also with respect to obs introduce a huge overhead? There might be only a few elements in params but thousands in obs, and for all of them I'm not interested in the gradient.

Ans does Clad really have to know the size of obs if it only contains constants? I mean if I don't use obs to forward it to an inner function it also works like expected.

Is there a workaround the problem without generating a different gradient?

Thanks again!

[vaithak]

Hi @guitargeek, I agree with you that having a gradient of obs will cause unnecessary operations. That is why when doing params[0] * obs[0] inside the wrapper function itself, we don't produce gradient operations for obs because we know the user hasn't requested it.
The issue is that when calling a function, Clad currently assumes that all function arguments need to be differentiated, which is wrong. We have a comment regarding this: https://github.com/vgvassilev/clad/blob/master/lib/Differentiator/ReverseModeVisitor.cpp#L1444
@vgvassilev I couldn't find any issue for tracking this, do you know if we already have one?

[guitargeek]

Ah okay, I see. It's not only about function parameters, this constant array can be anywhere. Here I put it in the global scope:

#include "clad/Differentiator/Differentiator.h"

#include <iostream>

double obs[]{1.0};

inline double foo(double *params, double *obs)
{
   return params[0] * obs[0];
}

double wrapper(double* params)
{
   // return params[0] * obs[0]; // works!

   return foo(params, obs); // doesn't work!
}

int main()
{
   auto grad = clad::gradient(wrapper, "params");

   double param = 2.0;
   double result = 0.0;
   grad.execute(&param, &result);

   std::cout << result << std::endl;
}

We need a way to push constant data into nested functions, even if it's an ugly workaround :D I almost found a way:

#include "clad/Differentiator/Differentiator.h"

#include <iostream>
#include <vector>

inline double foo(double *params, double *obs)
{
   return params[0] * obs[0];
}

double wrapper(double * params, double* obs)
{
   double obsCopy[2];
   for(std::size_t i = 0; i < 2; ++i) { obsCopy[i] = obs[i]; }

   // return params[0] * obs[0]; // works!
   //return foo(params, obsCopy); // works

   return foo(params, obsCopy + 1); // doesn't work!
}

int main()
{
   auto grad = clad::gradient(wrapper, "params");

   double param = 2.0;
   std::vector<double> obs{2.0, 3.0};
   double result = 0.0;
   grad.execute(&param, obs.data(), &result);

   std::cout << result << std::endl;
}

So close :) After copying the array to something with know size, I can use it! But only with a pointer to the beginning. I can't do pointer arithmetic :(

The usecase is to inject hundreds of constant arrays like this into the function. The function signature would be way to long if they would all be separate, that's why we works with one concatenated auxiliary array that is passed to the function, and then we use offsetting.

Is supporting the pointer arithmetic maybe more realistic?

edit: what we currently to is to write 100s of for-loops that extract the subarrays. I think that this might be a factor in the long gradient jitting time for RooFit models, which is why I was trying to find a better solution.

[vaithak]

Pointer arithmetic one should be easy to fix, let me try this one out. I will update you once I have a fix for that.

[vaithak]

Hi @guitargeek, the workaround should work after this: #764

[vgvassilev]

Hi @guitargeek, I agree with you that having a gradient of obs will cause unnecessary operations. That is why when doing params[0] * obs[0] inside the wrapper function itself, we don't produce gradient operations for obs because we know the user hasn't requested it. The issue is that when calling a function, Clad currently assumes that all function arguments need to be differentiated, which is wrong.

Don't we have a facility to specify wrt to which parameter we should differentiate against?

We have a comment regarding this: https://github.com/vgvassilev/clad/blob/master/lib/Differentiator/ReverseModeVisitor.cpp#L1444

That went out of date at the time I clicked on it :(

@vgvassilev I couldn't find any issue for tracking this, do you know if we already have one?

Probably worth creating a new issue?

[guitargeek]

Awesome! The workaround with the single array declaration and then offsetting works not, and it reduces jitting time for the gradient of the ATLAS models by a factor 5! Thanks a lot, @vaithak

[vgvassilev]

@guitargeek, glad to hear. Can we close this issue now?

[vaithak]

I think we can close this one, I have already created a new issue for tracking the request for dealing with differentiation with partial arguments in a function call: #765.