fix apply_scalar_unary and abs so that abs(complex<T>) return T

[SteveBronder]

Summary

Fix for #2737 so that abs(complex<T>) returns the type T and not complex<T>.

@WardBrian can you try this branch out with your stanc3 branch? This should have the stuff in it to allow the full vectorization.

@bob-carpenter this does also touch apply_scalar_unary, do we want abs(int) to return int or double? Right now it returns double but I modified apply_scalar_unary to return the same type as the output of the function we pass it (so sin(int) would still return an double but abs(int) can return an int. Is that alright?

Tests

A static assert is added to the abs() test to make sure that if we have complex<T> as the input we get back T as the output.

Side Effects

The apply_scalar_unary change is one we should discuss to make sure it's still doing what we would like.

Release notes

Fix #2737

Checklist

• Math issue allows complex abs() to be vectorized #2737

• Copyright holder: Steve Bronder

  The copyright holder is typically you or your assignee, such as a university or company. By submitting this pull request, the copyright holder is agreeing to the license the submitted work under the following licenses:
  - Code: BSD 3-clause (https://opensource.org/licenses/BSD-3-Clause)
  - Documentation: CC-BY 4.0 (https://creativecommons.org/licenses/by/4.0/)

• the basic tests are passing

  • unit tests pass (to run, use: ./runTests.py test/unit)
  • header checks pass, (make test-headers)
  • dependencies checks pass, (make test-math-dependencies)
  • docs build, (make doxygen)
  • code passes the built in C++ standards checks (make cpplint)

• the code is written in idiomatic C++ and changes are documented in the doxygen

• the new changes are tested

[rok-cesnovar]

can you try this branch out with your stanc3 branch? This should have the stuff in it to allow the full vectorization.

Just ran it: https://jenkins.flatironinstitute.org/blue/organizations/jenkins/Stan%2FStanc3/detail/PR-1195/11/pipeline

EDIT: It passed.

do we want abs(int) to return int or double? Right now it returns double but I modified apply_scalar_unary to return the same type as the output of the function we pass it (so sin(int) would still return an double but abs(int) can return an int. Is that alright?

unary Stan functions defined in stanc3 that have int as an input and expect int as an output (excluding _rng functions):

abs(int) => int
int_step(int) => int
size(int) => int
logical_negation(int) => int
minus(int) => int
plus(int) => int

The latter 3 are just how a = !b, a = -b and a = +b are named internally. So yeah, what you did is the correct behaviour in my opinion. But we should definitely wait for @bob-carpenter to confirm that this was the intention.

For completeness, these are the unary functions that return real on int input:

Phi(int) => real
Phi_approx(int) => real
acos(int) => real
acosh(int) => real
asin(int) => real
asinh(int) => real
atan(int) => real
atanh(int) => real
cbrt(int) => real
ceil(int) => real
cos(int) => real
cosh(int) => real
digamma(int) => real
erf(int) => real
erfc(int) => real
exp(int) => real
exp2(int) => real
expm1(int) => real
fabs(int) => real
floor(int) => real
inv(int) => real
inv_Phi(int) => real
inv_cloglog(int) => real
inv_erfc(int) => real
inv_logit(int) => real
inv_sqrt(int) => real
inv_square(int) => real
lambert_w0(int) => real
lambert_wm1(int) => real
lgamma(int) => real
log(int) => real
log10(int) => real
log1m(int) => real
log1m_exp(int) => real
log1m_inv_logit(int) => real
log1p(int) => real
log1p_exp(int) => real
log2(int) => real
log_inv_logit(int) => real
logit(int) => real
round(int) => real
sin(int) => real
sinh(int) => real
sqrt(int) => real
square(int) => real
tan(int) => real
tanh(int) => real
tgamma(int) => real
trigamma(int) => real
trunc(int) => real

[rok-cesnovar]

Just one minor cleanup. If Bob approves that the return types are correct, this should be good to merge.

[WardBrian]

I can confirm we want abs(int) to return int, and abs of containers of ints to return containers of ints. We don't necessarily want to change any other functions which use the vectorization framework though, so the rest of the signatures Rok listed should be unchanged

[bob-carpenter]

This is why we want to write clear issues. abs has been problematic because its return type was initially messed up in the language. Here's a hopefully clear spec on exactly what the signatures should be.

int abs(int)
real abs(real)
real abs(complex)

vector abs(vector)
row_vector abs(row_vector)
matrix abs(matrix)

vector abs(complex_vector)
row_vector abs(complex_row_vector)
matrix abs(complex_matrix)

and then we want to have

R[] abs(T[]) 

wherever we have R abs(T)

[SteveBronder]

It seems like for many of them we implicitly cast to double for int values which is fine. One thing I worry about is the compiler getting confused between whats available to cast to when it has access to the math library and the standard library. So for instance for exp's caller struct we have a function that imports std::exp and then calls it.

struct exp_fun {
  /**
   * Return the exponential of the specified scalar argument.
   *
   * @tparam T type of argument
   * @param[in] x argument
   * @return Exponential of argument.
   */
  template <typename T>
  static inline auto fun(const T& x) {
    using std::exp;
    return exp(x);
  }
};

The using std::exp is the troublesome part. If we have an exp() in math that can also take in / convert an integer as well then I think the compiler could get confused about which one to use. I think we want to have a function in stan that directly calls the std version so there's never confusion for the compiler about which one we want.

template <typename T, require_arithmetic_t<T>* = nullptr>
inline auto exp(T x) {
  return std::exp(x);
}  // Now delete the 'using std::exp' in the struct's fun caller

[SteveBronder]

I can do that in this PR or I can do it in a seperate PR. It might be good to nip that sort of issue now

[WardBrian]

Stanc always fully qualifies function names with stan::math, so I think we avoid that issue in generated code, but it is definitely a concern in the library itself

[SteveBronder]

This is why we want to write clear issues.

My trouble with writing clear issues is that I'm not sure I would have noticed the problem (and what I think is the solution) before I dove into the code!

I think the answer to abs() and others that want int -> double or int -> int is not to do that conversion at the apply_scalar_unary level but at the calling function level

[SteveBronder]

Stanc always fully qualifies function names with stan::math, so I think we avoid that issue in generated code, but it is definitely a concern in the library itself

Yes if I remember right the compiler confusion on implicit casts for std library functions is exactly why we did that

[WardBrian]

It was actually for overloading, previously we had a list of certain functions we would prefix with std:: to avoid the confusion, but when we wanted overloading to work for stan-math functions we needed qualifiers everywhere

[bob-carpenter]

Wouldn't we just get rid of using std::exp and explicitly call stan::math::exp if it's fully overloaded. The using business was for when we were expecting to call the standard library math functions. I believe the convention is that we now overload them all in the stan::math namespace and then the stanc3 compiler includes the stan::math qualifier in calling them.

[SteveBronder]

Wouldn't we just get rid of using std::exp and explicitly call stan::math::exp if it's fully overloaded. The using business was for when we were expecting to call the standard library math functions. I believe the convention is that we now overload them all in the stan::math namespace and then the stanc3 compiler includes the stan::math qualifier in calling them.

Yes so the full change would be

template <typename T, require_arithmetic_t<T>* = nullptr>
inline auto exp(T x) {
  return std::exp(x);
}

struct exp_fun {
  /**
   * Return the exponential of the specified scalar argument.
   *
   * @tparam T type of argument
   * @param[in] x argument
   * @return Exponential of argument.
   */
  template <typename T>
  static inline auto fun(const T& x) {
    return stan::math::exp(x);
  }
};

I don't think we need to call stan::math::exp() explicitly since now std::exp() is not in the lookup namespace but it doesn't hurt to have

[rok-cesnovar]

I don't think we need to call stan::math::exp() explicitly since now std::exp() is not in the lookup namespace but it doesn't hurt to have

There is a problem with some functions also being declared in the global namespace (I think if you include <math.h> if I remember correctly). So we definitely should keep it there.

[rok-cesnovar]

Thanks @SteveBronder!