Writing reliable templated code for both `var` and `double` instantiations

While working on fractional BesselK, I am having pains to properly compile templated functions for both double and var template params. A simplified example:

namespace stan {
namespace math {  
  template <typename T> 
    T test(const T& x) {
      return fabs(x);
    }
  

  double xx() {
    return test<double>(5);
  }
}
}

gets me:

./stan/math/rev/arr/fun/imports_debug.hpp: In instantiation of 'T stan::math::test(const T&) [with T = double]':
./stan/math/rev/arr/fun/imports_debug.hpp:15:26:   required from here
./stan/math/rev/arr/fun/imports_debug.hpp:10:20: error: cannot convert 'stan::math::var' to 'double' in return
       return fabs(x);

However, changing fabs to exp works. Additionally, if I add

namespace stan {
namespace math {  
  double fabs(double x) {
    return std::fabs(x);
  }
}
}

It compiles once again.

My conclusion is that this happens only for functions where stan::math does not contain the double version (e.g. fabs, pow), but not for functions where stan::math contains the double version (e.g. log, exp).

So my question is: why some functions have a double version in stan::math and some do not? And is there a “correct” solution that will let me use templating to write both double and var versions of a function using fabs or pow ?

EDIT: I am on Windows using gcc 4.7.0 as supplied with RTools.

One more observation: moving the function outside stan::math makes the problem go away, e.g. this compiles well:

template <typename T> 
  T test(const T& x) {
    return fabs(x);
  }

namespace stan {
namespace math {  
  double xx() {
    return test<double>(5);
  }

  var xx2() {
    return test<var>(var(5));
  }

}
}

EDIT: And just to vent a bit of frustration: My actual code now passes my unit test, but not header tests. Adding includes for rev and prim versions of base functions, i.e.:

#include <stan/math/rev/scal/fun/exp.hpp>
#include <stan/math/prim/scal/fun/exp.hpp>

Makes the header tests pass, but it changes what is computed, i.e. with those includes in place, my unit test compiles but my function gives different results for some of the gradients. This is getting positively nightmarish. My best guess is that I only had the rev version of the functions which caused some auto-conversion from double to var in my code which now does not happen.

EDIT2: I fixed the nightmarish problem, but I have no idea why it works. I kept the headers and I changed exp to std::exp in this piece of code:

namespace stan { namespace math {

template <typename T_v>
var log_modified_bessel_second_kind_frac(const T_v &v, const var &z) {
  double value_vm1;
  T_v value;
  //Details omitted
  double gradient_dz
      = -std::exp(value_vm1 - value_of(value)) - value_of(v) / z.val();
  //Details omitted 
}

}} //namespaces

The call to exp is obviously with a double argument. Why would stan::math::exp do anything different than std::exp here?
I have similar issues with log where being explicit about which one I call changes the results, but the example is not so clear.

1 Like

I think this is partly getting fixed now that the C99 functions have been rolled into C++11. Before that, we didn’t have a reliable way to get the functions that were in C99 but not in C++03.

Now, I think there are some error handling issues where we still want to use boost (lgamma, specifically). I can’t recall what the final resolution on all this is, but @syclik should know. There was a PR floating around that is going to or already has fixed a lot of this (by replacing Boost versions with std:: versions).

Thanks for the reply, but I think it does not really answer my concern - why it is so hard to get proper resolution of basic functions in templated code? I might be missing something, but I don’t see how that is connected to boost vs. std:: as my problem is in resolving std:: vs. stan::math::.

I find it highly plausible that some basic misconception of how C++ actually works is biting me here and doing something properly will resolve it. Or do people just not write templated code in math and always separate the double and var variants?

I’m not exactly sure what you’re asking. The way resolution works is that if you have a call, you look at everything function whose symbol’s been imported, find the one with the best matching base template, and fail if there is no such. Then, after finding best matching base template, it needs to find best specialization, and will also fail if there is no such.

The problem we run into in our math code is that we want to write something like exp(y) rather than std::exp(y) and let it figure out which exp to use based on argument-dependent lookup. But that can get complicated when the compilers bring in their own versions of ::exp(y) in the top-level namespace or in std:: if we’ve been too eager importing std:: symbols. This is why it’s bad to have using at the namespace level—just too much room for ambiguity.

If you know the type of argument, you can just fully specify.

The title of the post was possibly confusing, sorry. I’ve edited it to more closely match my inquiry.

I am asking, what is the best practice to write templated code (as in the first post) for the math library that reliably works for both var and double instantiations. In particular, how to correctly use functions that reside both in stan::math:: and std:: to make sure the correct version is called.

Additionally if you have any guesses why changing exp to std::exp may change the numerical values of results (for double instantiations), as in the second post. But for this case, I don’t have a short reproducible example yet, so I understand it might be hard (but if you want the full ~300 lines examples, you may have it).

Thanks for the clarification. What you’re looking for is something called argument dependent lookup (link to cppreference).

From the original example, it would work if you did this:

namespace stan { 
  namespace math { 
    template <typename T>
    T test(const T& x) { 
       using std::fabs;
       return fabs(x); 
    }
  }
}

Edit: this is not true. Please read the next comment for an explanation. If you were outside the stan::math:: namespace, you’d also want to add a using stan::math::fabs;. ADL will help choose the right function call for you since you put in the using statements.

It’s working for exp because I bet there’s an exp in the global namespace that’s brought in using one of your includes.

1 Like

No, you don’t want to do this.

What you need is a using statement to bring in any primitive definition that’s needed.

ADL will take care of the non-primitive definitions by bringing in the definition of fabs from whatever namespace x is in. All you need to do is include the appropriate deifnitions.

2 Likes

Yikes! Yes, @Bob_Carpenter is absolutely right! ADL will take care of the non-primitive definitions.

(I’ll edit my post striking out the mistake so it’s clear it’s not true.)

Thanks for the clarifications and patience! I’ve tried using using std::fabs before and it seemed to not work, but it turns out it didn’t work because of other dependencies in my code and when I do it with more care, it resolves my problems! Thanks!

2 Likes