I my attempts to implement custom probability densities I have never encountered this error. I am trying to implement a SICHEL lpmf the reparametrisation here.

I could not find much online on this error

```
Chain 1: 0.378937 0.0977251 0.355928
-5.8373
-6.85597
-6.92309
-6.44776
-5.93319
-5.69185
-5.18989
-6.28615
-5.78904
-7.88485
-5.74056
-5.36186
-5.49455
-5.13792
Chain 1: Exception: Exception: inner_integral_grad_v: The gradient of inner_integral is nan 0 (in 'model3f1f757baa8e_poisson_GIG' at line 10)
(in 'model3f1f757baa8e_poisson_GIG' at line 36)
[1] "Error in sampler$call_sampler(args_list[[i]]) : "
[2] " Exception: Exception: inner_integral_grad_v: The gradient of inner_integral is nan 0 (in 'model3f1f757baa8e_poisson_GIG' at line 10)"
[3] " (in 'model3f1f757baa8e_poisson_GIG' at line 36)"
```

I get log probabilities as negative number as expected, and I am not even sampling with the sichel function, just printing

The model:

```
functions{
real log_besselk_frac(real v, real z);
real sichel_lpmf(int y, real gamma, real alpha, real epsilon){
// https://www.researchgate.net/publication/247109928_Parameter_Estimation_for_the_Sichel_Distribution_and_Its_Multivariate_Extension
real w = hypot(epsilon, alpha) - epsilon;
return
( (gamma * (log(w) - log(alpha))) - log_besselk_frac(gamma, w) ) +
(y * (log(epsilon) + log(w) - log(alpha))) - lgamma(y+1) +
log_besselk_frac(y + gamma, alpha);
}
}
data{
int<lower=0> N;
int y[N];
}
parameters{
real<lower=-0.5, upper=0.5> gamma;
real<lower=0, upper=0.5> alpha;
real<lower=0, upper=0.5> epsilon;
}
model {
gamma ~ normal(0,1);
alpha ~ normal(0,1);
epsilon ~ normal(0,1);
print(gamma, " " , alpha, " " , epsilon);
for(i in 1:N) print(sichel_lpmf(y[i] | gamma, alpha , epsilon));
}
```