2 #include "../dpdefs/csprogsdefs.qh"
6 #include "../dpdefs/progsdefs.qh"
7 #include "../dpdefs/dpextensions.qh"
11 int fpclassify(float x)
21 bool isfinite(float x)
23 return !(isnan(x) || isinf(x));
27 return (x != 0) && (x + x == x);
35 bool isnormal(float x)
46 return log(x + sqrt(x*x - 1));
50 return log(x + sqrt(x*x + 1));
54 return 0.5 * log((1+x) / (1-x));
58 return 0.5 * (exp(x) + exp(-x));
62 return 0.5 * (exp(x) - exp(-x));
66 return sinh(x) / cosh(x);
92 return floor(log2(fabs(x)));
94 float ldexp(float x, int e)
100 return log(x) * M_LOG10E;
108 return log(x) * M_LOG2E;
112 return floor(log2(fabs(x)));
116 return '1 0 0' * (f - trunc(f)) + '0 1 0' * trunc(f);
119 float scalbn(float x, int n)
121 return x * pow(2, n);
126 return copysign(pow(fabs(x), 1.0/3.0), x);
128 float hypot(float x, float y)
130 return sqrt(x*x + y*y);
135 // approximation taken from wikipedia
138 return copysign(sqrt(1 - exp(-y * (1.273239544735163 + 0.14001228868667 * y) / (1 + 0.14001228868667 * y))), x);
144 vector lgamma(float x)
146 // TODO improve accuracy
148 return fabs(x) * '1 0 0' + copysign(1, x) * '0 1 0';
149 if(x < 1 && x == floor(x))
150 return nan("gamma") * '1 1 1';
155 // reflection formula:
156 // gamma(1-z) * gamma(z) = pi / sin(pi*z)
157 // lgamma(1-z) + lgamma(z) = log(pi) - log(sin(pi*z))
158 // sign of gamma(1-z) = sign of gamma(z) * sign of sin(pi*z)
160 v.x = log(M_PI) - log(fabs(v.z)) - v.x;
167 return lgamma(x + 1) - log(x) * '1 0 0';
169 return (0.5 * log(2 * M_PI * x) + x * (log(x) - 1)) * '1 0 0' + '0 1 0';
171 float tgamma(float x)
175 return exp(v.x) * v.y;
187 float pymod(float x, float y)
189 return x - y * floor(x / y);
192 float nearbyint(float x)
198 return (x>=0) ? floor(x) : ceil(x);
201 float fmod(float x, float y)
203 return x - y * trunc(x / y);
205 float remainder(float x, float y)
207 return x - y * rint(x / y);
209 vector remquo(float x, float y)
218 float copysign(float x, float y)
220 return fabs(x) * ((y>0) ? 1 : -1);
222 float nan(string tag)
226 float nextafter(float x, float y)
230 return nan("nextafter");
232 return -nextafter(-x, -y);
233 // now we know that x < y
234 // so we need the next number > x
236 d = max(fabs(x), 0.00000000000000000000001);
247 float nexttoward(float x, float y)
249 return nextafter(x, y);
252 float fdim(float x, float y)
256 float fmax(float x, float y)
260 float fmin(float x, float y)
264 float fma(float x, float y, float z)
269 int isgreater(float x, float y)
273 int isgreaterequal(float x, float y)
277 int isless(float x, float y)
281 int islessequal(float x, float y)
285 int islessgreater(float x, float y)
287 return x < y || x > y;
289 int isunordered(float x, float y)
291 return !(x < y || x == y || x > y);
294 vector cross(vector a, vector b)
297 '1 0 0' * (a.y * b.z - a.z * b.y)
298 + '0 1 0' * (a.z * b.x - a.x * b.z)
299 + '0 0 1' * (a.x * b.y - a.y * b.x);