diff options
Diffstat (limited to 'src/base/rng/poisson.c')
| -rw-r--r-- | src/base/rng/poisson.c | 21 |
1 files changed, 10 insertions, 11 deletions
diff --git a/src/base/rng/poisson.c b/src/base/rng/poisson.c index 3ec15c9..776070c 100644 --- a/src/base/rng/poisson.c +++ b/src/base/rng/poisson.c @@ -28,7 +28,7 @@ log1pmx(double x, double off) return x*x*r; } - return log(1+x) - off; + return math·log(1+x) - off; } static inline @@ -49,14 +49,14 @@ procf(double mu, double s, int64 K, double *px, double *py, double *fx, double * if(K < 10) { *px = -mu; - *py = pow(mu,K) / factorial[K]; + *py = math·pow(mu,K) / factorial[K]; }else{ d = 0.08333333333333333/K; d = d - 4.8*d*d*d; V = (mu - K) / K; *px = K*log1pmx(V,mu-K) - d; - *py = 0.3989422804014327/sqrt(K); + *py = 0.3989422804014327/math·sqrt(K); } X = (K - mu + 0.5)/s; @@ -66,21 +66,20 @@ procf(double mu, double s, int64 K, double *px, double *py, double *fx, double * return c; } -static inline -uint64 +static inline uint64 bigpoisson(double mu) { int64 L,K; double G,s,d,U,E,T; double px,py,fx,fy,c; - s = sqrt(mu); + s = math·sqrt(mu); d = 6*mu*mu; - L = floor(mu - 1.1484); + L = math·floor(mu - 1.1484); stepN: G = mu + s*rng·normal(); - K = floor(G); + K = math·floor(G); if(K<0) goto stepP; stepI: @@ -99,13 +98,13 @@ stepE: E = rng·exponential(1.0); U = rng·random(); U = U + U - 1; - T = 1.8 + copysign(E,U); + T = 1.8 + math·copysign(E,U); if(T < 0.6744) goto stepE; - K = floor(mu + s*T); + K = math·floor(mu + s*T); c = procf(mu, s, K, &px, &py, &fx, &fy); stepH: - if(c*fabs(U) > (py*exp(px + E) - fy*exp(fx + E))) + if(c*math·fabs(U) > (py*math·exp(px + E) - fy*math·exp(fx + E))) goto stepE; return K; } |