-
-
Notifications
You must be signed in to change notification settings - Fork 188
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Vectorize everything with broadcasting #202
Comments
I'm reopening this issue because the pull request that was merged in only provides the infrastructure to test these vectorize functions. However, if you think that we should close this issue and open a new one, I can do that. |
Either way is OK. You need to make the checklist of all On Mar 15, 2016, at 4:52 PM, rayleigh [email protected] wrote:
|
Okay. Copying over your response to a question of where to put files:
This should declare all the basic vectorizations and
Also, where should the test files go? Should they go in |
Thanks and yes, that's the right place for the test files, because
|
Thanks. Looking through the 2.8.0 Stan manual, I noticed that |
Ack, I think we undeprecated abs and are going to
|
Sounds good; I vectorized both. I've been trying to vectorize the function |
I think we should test derivatives. We had this discussion in the past regarding defined values and undefined On Mon, Mar 21, 2016 at 7:00 PM, rayleigh [email protected] wrote:
|
Okay. Just to clarify, should the vectorize-all function or the individual scalar functions include the error checking? Or, both? |
That was a general comment. What test am I supposed to run? Can you provide the full path? On Tue, Mar 22, 2016 at 11:34 AM, rayleigh [email protected] wrote:
|
I pushed more changes to the branch |
I pushed the change to |
Is it only exp() that's causing the issue?
|
No. This isn't an include-order issue. Here's where I'm at so far with exp:
I think if I can figure out why that specialization is not being used, it Daniel On Tue, Mar 22, 2016 at 11:09 PM, Bob Carpenter [email protected]
|
Ok. Those headers aren't included in stan/math/prim/mat.hpp or rev/mat.hpp I'll fix that now. On Tue, Mar 22, 2016 at 11:22 PM, Daniel Lee [email protected] wrote:
|
I just pushed, but it's not working yet. Bob, I might need your help apply_scalar_unary() in prim is being called first. Even when the mix If you include the rev version as the first thing in the test file, it I'm not sure how to get past this problem. I'm sure we've figured this out On Tue, Mar 22, 2016 at 11:24 PM, Daniel Lee [email protected] wrote:
|
Any luck getting past this problem? |
Copying this over from an email chain because it might have the solution (I haven't tried yet): If you edit out the very first include line in: /test/unit/math/mix/mat/vectorize/foo_test.cpp, everything passes. So the new file looks like:
This makes me think there may be a problem in the include order in mix/mat.hpp. You might want to read this: http://stackoverflow.com/questions/14626033/implicit-instantiation-of-undefined-template-class It's all about steps 1, 2, and 3 in the answer. |
Other questions:
So, is trunc.hpp designed to work on matrix?
|
Listing why certain unary functions weren't vectorized: Unvectorized exp2, log2, and log1p because of an ambiguity error when making Stan model's C++ code for untemplated prim versions. trunc and step aren't vectorized because both aren't compatible with matrices. Other (real):real functions aren't vectorized because they're missing a mode and their prim version can't be untemplated. |
Let's just get what you have working in and then create new issues
|
Sounds good though I disagree about the solution to 1:
namespace { class acosh_vari : public op_v_vari { public: acosh_vari(double val, vari* avi) : op_v_vari(val, avi) { } void chain() { avi_->adj_ += adj_ / std::sqrt(avi_->val_ * avi_->val_ - 1.0); } }; } In the rev version of digamma.hpp (which also uses code from namespace { class digamma_vari : public op_v_vari { public: explicit digamma_vari(vari* avi) : op_v_vari(digamma(avi->val_), avi) { } void chain() { avi_->adj_ += adj_ * trigamma(avi_->val_); } }; } The vectorized version of digamma handles error consistently across all versions so I think the rev and fwd mode code needs to be fixed for acosh, atanh, and asinh. |
Here's the actual function for acosh: inline var acosh(const var& a) { We could move the call to ::acosh() into acosh_vari constructor, So I'm not sure what you're suggesting or where there is
|
Thanks for pointing that out. I took another look and I think the includes might be the issue because for
So, unless it's being compiled on a Visual C++ compiler, I don't think it'll use |
That conditional include is nasty. Maybe we should just use Is there ever a reason to include both <math.h> and ? The error handling can be fixed by just adding it to the
|
@rayleigh Is the current checklist above up to date? Is there a reason the dozen or so |
Closing in favor of newer issue #347 given that it's partially done. |
Can we reopen this? The logical operators (at least) seem to remain un-vectorized
Thanks, |
@jessexknight Probably better to start a new issue specifically or logical operations, because it's not immediately obvious how they should behave. This issue was specifically for real-valued operations and used a specific vectorization assuming real valued functions with autodiff---here there are no derivatives because we get integer values out. I can imagine two alternatives for the logical operators:
Both approaches can broadcast scalars to containers. If we go with (2), then we probably want functions This is related to the issue of comparing containers to elements to return booleans, the way that R allows. That'd give us something like this int[] xs = { 1, 0, 2, 3, 0, -1};
int[] ys = (xs == 0); // after eval, ys = { 2, 5 } |
Thanks, I've evidently added an issue... Just want to add that I would strongly prefer (2) to allow more flexibility. Not sure I follow the last issue you mentioned -- I'm not familiar with this behaviour in R, and personally wouldn't find this a priority vs the vectorization of the base logical functions. |
Thanks for opening the issue. I edited it slightly to make it easier for our devs. In R, you can do this: > sex = c(0, 1, 1, 0, 0, 0, 1)
> age = c(23, 29, 30, 12, 15, 18)
> age[sex == 0]
[1] 23 12 15 18
> age[sex == 1]
[1] 29 30 NA
> age = c(23, 29, 30, 22, 25, 18, 31)
> sex = c(0, 1, 1 , 0, 0 , 0, 1)
> sex == 0
[1] TRUE FALSE FALSE TRUE TRUE TRUE FALSE
> age[sex == 0]
[1] 23 22 25 18 As you can see, it's useful for picking subgroups out of parallel sequences. But it relies on really odd behavior where if you give R a list of boolean arguments, it'll include the ones that are > age[c(1, 0, 0, 1, 1, 1, 0)]
[1] 23 23 23 23 This is a very R result in that it seems to just ignore the out of range |
Oh, I see. I think this logical indexing is available in Numpy and Matlab too. In fact, I think logical indexing can be faster in Numpy, besides the fact that I suppose this raises the idea of a logical data type in Stan, but I think this type of indexing might be one of the only use cases ... |
It's always going to be bound by having to evaluate the condition for every element of the container. You could potentially do it without constructing the intermediate
We've already committed to the C-style coding of 0 for false and everything else being true, so if we did go down the boolean route, we'd at least need the type to be promotable to |
As a bit of background for this: |
Vectorize all functions for matching shapes and allow broadcasting of lower-dimensional types.
This will require
There is an example for the function
inv()
instan/math/prim/mat/fun/inv.hpp
in the branchfeature/issue-202-vectorize-all
.This issue involves the first step above for all of the Stan functions, so it may need to be broken down into stages, starting with unary functions.
For the stan-dev/stan issue: stan-dev/stan#1683
For each of thee functions, we need (1) functor, (2) general template definition, (3) doc in manual, (4) signatures [need function for this], (5) signature tests for all instantiations.
When I did the example for
inv()
, I had to remove the templating forinv()
itself and replace it withdouble
to remove the ambiguity. For functions that only have a template version (i.e., there is not also a version forfwd
andrev
), I believe there will need to be anenable_if
on the base definition that restricts application to the base class.We'll also need a testing framework. This should be easy to do following the way
apply_scalar_unary
is defined (acrossprim
,rev
, andfwd
). It will require a simple templated testing functor for each function being tested.(int):int
abs
int_step
operator!
operator+
operator-
(int,int):int
max
min
operator*
operator+
operator-
operator/
operator<
operator<=
operator>
operator>=
operator==
operator&&
operator||
(int,real):real
bessel_first_kind
bessel_second_kind
modified_bessel_first_kind
modified_bessel_second_kind
(real):int
int_step
is_inf
is_nan
operator!
(real):real
abs
acos
acosh
asin
asinh
atan
atanh
cbrt
ceil
cos
cosh
digamma
erf
erfc
exp
exp2
expm1
fabs
floor
inv
inv_cloglog
inv_logit
inv_Phi
inv_sqrt
inv_square
lgamma
log
log10
log1m
log1m_exp
log1m_inv_logit
log1p
log1p_exp
log2
log_inv_logit
logit
Phi
Phi_approx
round
sin
sinh
sqrt
square
step
tan
tanh
tgamma
trigamma
trunc
(real,real):int
operator<
operator<=
operator>
operator>=
operator==
operator&&
operator||
(real,real):real
atan2
binomial_coefficient_log
falling_factorial
fdim
fmax
fmin
fmod
gamma_p
gamma_q
hypot
lbeta
log_diff_exp
log_falling_factorial
log_rising_factorial
log_sum_exp
max
min
multiply_log
operator*
operator+
operator-
operator/
operator^
owens_t
pow
rising_factorial
(int, real):real
binary_log_loss
lmgamma
(real,real,real):real
fma
log_mix
(int,real,real):real
if_else
The text was updated successfully, but these errors were encountered: