source: libcfa/src/math.hfa@ 5ad6f0d

Last change on this file since 5ad6f0d was 658f3179, checked in by Andrew Beach <ajbeach@…>, 9 months ago

Moved massive function log2_u32_32 out of header.
  • Property mode set to 100644
File size: 21.7 KB
Line 
1//
2// Cforall Version 1.0.0 Copyright (C) 2016 University of Waterloo
3//
4// The contents of this file are covered under the licence agreement in the
5// file "LICENCE" distributed with Cforall.
6//
7// math.hfa --
8//
9// Author : Peter A. Buhr
10// Created On : Mon Apr 18 23:37:04 2016
11// Last Modified By : Peter A. Buhr
12// Last Modified On : Sun Jun 18 08:13:53 2023
13// Update Count : 202
14//
15
16#pragma once
17
18#include <math.h>
19#include <complex.h>
20
21//---------------------------------------
22
23#include "common.hfa"
24#include "bits/debug.hfa"
25
26//---------------------- General ----------------------
27
28static inline __attribute__((always_inline)) {
29 float ?%?( float x, float y ) { return fmodf( x, y ); }
30 float fmod( float x, float y ) { return fmodf( x, y ); }
31 double ?%?( double x, double y ) { return fmod( x, y ); }
32 // extern "C" { double fmod( double, double ); }
33 long double ?%?( long double x, long double y ) { return fmodl( x, y ); }
34 long double fmod( long double x, long double y ) { return fmodl( x, y ); }
35
36 float remainder( float x, float y ) { return remainderf( x, y ); }
37 // extern "C" { double remainder( double, double ); }
38 long double remainder( long double x, long double y ) { return remainderl( x, y ); }
39
40 float remquo( float x, float y, int * quo ) { return remquof( x, y, quo ); }
41 // extern "C" { double remquo( double x, double y, int * quo ); }
42 long double remquo( long double x, long double y, int * quo ) { return remquol( x, y, quo ); }
43 [ int, float ] remquo( float x, float y ) { int quo; x = remquof( x, y, &quo ); return [ quo, x ]; }
44 [ int, double ] remquo( double x, double y ) { int quo; x = remquo( x, y, &quo ); return [ quo, x ]; }
45 [ int, long double ] remquo( long double x, long double y ) { int quo; x = remquol( x, y, &quo ); return [ quo, x ]; }
46
47 [ float, float ] div( float x, float y ) { y = modff( x / y, &x ); return [ x, y ]; }
48 [ double, double ] div( double x, double y ) { y = modf( x / y, &x ); return [ x, y ]; }
49 [ long double, long double ] div( long double x, long double y ) { y = modfl( x / y, &x ); return [ x, y ]; }
50
51 float fma( float x, float y, float z ) { return fmaf( x, y, z ); }
52 // extern "C" { double fma( double, double, double ); }
53 long double fma( long double x, long double y, long double z ) { return fmal( x, y, z ); }
54
55 float fdim( float x, float y ) { return fdimf( x, y ); }
56 // extern "C" { double fdim( double, double ); }
57 long double fdim( long double x, long double y ) { return fdiml( x, y ); }
58
59 float nan( const char tag[] ) { return nanf( tag ); }
60 // extern "C" { double nan( const char [] ); }
61 long double nan( const char tag[] ) { return nanl( tag ); }
62} // distribution
63
64//---------------------- Exponential ----------------------
65
66static inline __attribute__((always_inline)) {
67 float exp( float x ) { return expf( x ); }
68 // extern "C" { double exp( double ); }
69 long double exp( long double x ) { return expl( x ); }
70 float _Complex exp( float _Complex x ) { return cexpf( x ); }
71 double _Complex exp( double _Complex x ) { return cexp( x ); }
72 long double _Complex exp( long double _Complex x ) { return cexpl( x ); }
73
74 float exp2( float x ) { return exp2f( x ); }
75 // extern "C" { double exp2( double ); }
76 long double exp2( long double x ) { return exp2l( x ); }
77 //float _Complex exp2( float _Complex x ) { return cexp2f( x ); }
78 //double _Complex exp2( double _Complex x ) { return cexp2( x ); }
79 //long double _Complex exp2( long double _Complex x ) { return cexp2l( x ); }
80
81 float expm1( float x ) { return expm1f( x ); }
82 // extern "C" { double expm1( double ); }
83 long double expm1( long double x ) { return expm1l( x ); }
84
85 float pow( float x, float y ) { return powf( x, y ); }
86 // extern "C" { double pow( double, double ); }
87 long double pow( long double x, long double y ) { return powl( x, y ); }
88 float _Complex pow( float _Complex x, float _Complex y ) { return cpowf( x, y ); }
89 double _Complex pow( double _Complex x, double _Complex y ) { return cpow( x, y ); }
90 long double _Complex pow( long double _Complex x, long double _Complex y ) { return cpowl( x, y ); }
91} // distribution
92
93//---------------------- Logarithm ----------------------
94
95static inline __attribute__((always_inline)) {
96 float log( float x ) { return logf( x ); }
97 // extern "C" { double log( double ); }
98 long double log( long double x ) { return logl( x ); }
99 float _Complex log( float _Complex x ) { return clogf( x ); }
100 double _Complex log( double _Complex x ) { return clog( x ); }
101 long double _Complex log( long double _Complex x ) { return clogl( x ); }
102
103 // O(1) polymorphic integer log2, using clz, which returns the number of leading 0-bits, starting at the most
104 // significant bit (single instruction on x86)
105 int log2( unsigned int n ) { return n == 0 ? -1 : sizeof(n) * __CHAR_BIT__ - 1 - __builtin_clz( n ); }
106 long int log2( unsigned long int n ) { return n == 0 ? -1 : sizeof(n) * __CHAR_BIT__ - 1 - __builtin_clzl( n ); }
107 long long int log2( unsigned long long int n ) { return n == 0 ? -1 : sizeof(n) * __CHAR_BIT__ - 1 - __builtin_clzll( n ); }
108 float log2( float x ) { return log2f( x ); }
109 // extern "C" { double log2( double ); }
110 long double log2( long double x ) { return log2l( x ); }
111 // float _Complex log2( float _Complex x ) { return clog2f( x ); }
112 // double _Complex log2( double _Complex x ) { return clog2( x ); }
113 // long double _Complex log2( long double _Complex x ) { return clog2l( x ); }
114
115 float log10( float x ) { return log10f( x ); }
116 // extern "C" { double log10( double ); }
117 long double log10( long double x ) { return log10l( x ); }
118 // float _Complex log10( float _Complex x ) { return clog10f( x ); }
119 // double _Complex log10( double _Complex x ) { return clog10( x ); }
120 // long double _Complex log10( long double _Complex x ) { return clog10l( x ); }
121
122 float log1p( float x ) { return log1pf( x ); }
123 // extern "C" { double log1p( double ); }
124 long double log1p( long double x ) { return log1pl( x ); }
125
126 int ilogb( float x ) { return ilogbf( x ); }
127 // extern "C" { int ilogb( double ); }
128 int ilogb( long double x ) { return ilogbl( x ); }
129
130 float logb( float x ) { return logbf( x ); }
131 // extern "C" { double logb( double ); }
132 long double logb( long double x ) { return logbl( x ); }
133
134 float sqrt( float x ) { return sqrtf( x ); }
135 // extern "C" { double sqrt( double ); }
136 long double sqrt( long double x ) { return sqrtl( x ); }
137 float _Complex sqrt( float _Complex x ) { return csqrtf( x ); }
138 double _Complex sqrt( double _Complex x ) { return csqrt( x ); }
139 long double _Complex sqrt( long double _Complex x ) { return csqrtl( x ); }
140
141 float cbrt( float x ) { return cbrtf( x ); }
142 // extern "C" { double cbrt( double ); }
143 long double cbrt( long double x ) { return cbrtl( x ); }
144
145 float hypot( float x, float y ) { return hypotf( x, y ); }
146 // extern "C" { double hypot( double, double ); }
147 long double hypot( long double x, long double y ) { return hypotl( x, y ); }
148} // distribution
149
150unsigned long long log2_u32_32( unsigned long long val );
151
152//---------------------- Trigonometric ----------------------
153
154static inline __attribute__((always_inline)) {
155 float sin( float x ) { return sinf( x ); }
156 // extern "C" { double sin( double ); }
157 long double sin( long double x ) { return sinl( x ); }
158 float _Complex sin( float _Complex x ) { return csinf( x ); }
159 double _Complex sin( double _Complex x ) { return csin( x ); }
160 long double _Complex sin( long double _Complex x ) { return csinl( x ); }
161
162 float cos( float x ) { return cosf( x ); }
163 // extern "C" { double cos( double ); }
164 long double cos( long double x ) { return cosl( x ); }
165 float _Complex cos( float _Complex x ) { return ccosf( x ); }
166 double _Complex cos( double _Complex x ) { return ccos( x ); }
167 long double _Complex cos( long double _Complex x ) { return ccosl( x ); }
168
169 float tan( float x ) { return tanf( x ); }
170 // extern "C" { double tan( double ); }
171 long double tan( long double x ) { return tanl( x ); }
172 float _Complex tan( float _Complex x ) { return ctanf( x ); }
173 double _Complex tan( double _Complex x ) { return ctan( x ); }
174 long double _Complex tan( long double _Complex x ) { return ctanl( x ); }
175
176 float asin( float x ) { return asinf( x ); }
177 // extern "C" { double asin( double ); }
178 long double asin( long double x ) { return asinl( x ); }
179 float _Complex asin( float _Complex x ) { return casinf( x ); }
180 double _Complex asin( double _Complex x ) { return casin( x ); }
181 long double _Complex asin( long double _Complex x ) { return casinl( x ); }
182
183 float acos( float x ) { return acosf( x ); }
184 // extern "C" { double acos( double ); }
185 long double acos( long double x ) { return acosl( x ); }
186 float _Complex acos( float _Complex x ) { return cacosf( x ); }
187 double _Complex acos( double _Complex x ) { return cacos( x ); }
188 long double _Complex acos( long double _Complex x ) { return cacosl( x ); }
189
190 float atan( float x ) { return atanf( x ); }
191 // extern "C" { double atan( double ); }
192 long double atan( long double x ) { return atanl( x ); }
193 float _Complex atan( float _Complex x ) { return catanf( x ); }
194 double _Complex atan( double _Complex x ) { return catan( x ); }
195 long double _Complex atan( long double _Complex x ) { return catanl( x ); }
196
197 float atan2( float x, float y ) { return atan2f( x, y ); }
198 // extern "C" { double atan2( double, double ); }
199 long double atan2( long double x, long double y ) { return atan2l( x, y ); }
200
201 // alternative name for atan2
202 float atan( float x, float y ) { return atan2f( x, y ); }
203 double atan( double x, double y ) { return atan2( x, y ); }
204 long double atan( long double x, long double y ) { return atan2l( x, y ); }
205} // distribution
206
207//---------------------- Hyperbolic ----------------------
208
209static inline __attribute__((always_inline)) {
210 float sinh( float x ) { return sinhf( x ); }
211 // extern "C" { double sinh( double ); }
212 long double sinh( long double x ) { return sinhl( x ); }
213 float _Complex sinh( float _Complex x ) { return csinhf( x ); }
214 double _Complex sinh( double _Complex x ) { return csinh( x ); }
215 long double _Complex sinh( long double _Complex x ) { return csinhl( x ); }
216
217 float cosh( float x ) { return coshf( x ); }
218 // extern "C" { double cosh( double ); }
219 long double cosh( long double x ) { return coshl( x ); }
220 float _Complex cosh( float _Complex x ) { return ccoshf( x ); }
221 double _Complex cosh( double _Complex x ) { return ccosh( x ); }
222 long double _Complex cosh( long double _Complex x ) { return ccoshl( x ); }
223
224 float tanh( float x ) { return tanhf( x ); }
225 // extern "C" { double tanh( double ); }
226 long double tanh( long double x ) { return tanhl( x ); }
227 float _Complex tanh( float _Complex x ) { return ctanhf( x ); }
228 double _Complex tanh( double _Complex x ) { return ctanh( x ); }
229 long double _Complex tanh( long double _Complex x ) { return ctanhl( x ); }
230
231 float asinh( float x ) { return asinhf( x ); }
232 // extern "C" { double asinh( double ); }
233 long double asinh( long double x ) { return asinhl( x ); }
234 float _Complex asinh( float _Complex x ) { return casinhf( x ); }
235 double _Complex asinh( double _Complex x ) { return casinh( x ); }
236 long double _Complex asinh( long double _Complex x ) { return casinhl( x ); }
237
238 float acosh( float x ) { return acoshf( x ); }
239 // extern "C" { double acosh( double ); }
240 long double acosh( long double x ) { return acoshl( x ); }
241 float _Complex acosh( float _Complex x ) { return cacoshf( x ); }
242 double _Complex acosh( double _Complex x ) { return cacosh( x ); }
243 long double _Complex acosh( long double _Complex x ) { return cacoshl( x ); }
244
245 float atanh( float x ) { return atanhf( x ); }
246 // extern "C" { double atanh( double ); }
247 long double atanh( long double x ) { return atanhl( x ); }
248 float _Complex atanh( float _Complex x ) { return catanhf( x ); }
249 double _Complex atanh( double _Complex x ) { return catanh( x ); }
250 long double _Complex atanh( long double _Complex x ) { return catanhl( x ); }
251} // distribution
252
253//---------------------- Error / Gamma ----------------------
254
255static inline __attribute__((always_inline)) {
256 float erf( float x ) { return erff( x ); }
257 // extern "C" { double erf( double ); }
258 long double erf( long double x ) { return erfl( x ); }
259 // float _Complex erf( float _Complex );
260 // double _Complex erf( double _Complex );
261 // long double _Complex erf( long double _Complex );
262
263 float erfc( float x ) { return erfcf( x ); }
264 // extern "C" { double erfc( double ); }
265 long double erfc( long double x ) { return erfcl( x ); }
266 // float _Complex erfc( float _Complex );
267 // double _Complex erfc( double _Complex );
268 // long double _Complex erfc( long double _Complex );
269
270 float lgamma( float x ) { return lgammaf( x ); }
271 // extern "C" { double lgamma( double ); }
272 long double lgamma( long double x ) { return lgammal( x ); }
273 float lgamma( float x, int * sign ) { return lgammaf_r( x, sign ); }
274 double lgamma( double x, int * sign ) { return lgamma_r( x, sign ); }
275 long double lgamma( long double x, int * sign ) { return lgammal_r( x, sign ); }
276
277 float tgamma( float x ) { return tgammaf( x ); }
278 // extern "C" { double tgamma( double ); }
279 long double tgamma( long double x ) { return tgammal( x ); }
280} // distribution
281
282//---------------------- Nearest Integer ----------------------
283
284inline __attribute__((always_inline)) static {
285 // force divide before multiply
286 signed char floor( signed char n, signed char align ) { return (n / align) * align; }
287 unsigned char floor( unsigned char n, unsigned char align ) { return (n / align) * align; }
288 short int floor( short int n, short int align ) { return (n / align) * align; }
289 unsigned short int floor( unsigned short int n, unsigned short int align ) { return (n / align) * align; }
290 int floor( int n, int align ) { return (n / align) * align; }
291 unsigned int floor( unsigned int n, unsigned int align ) { return (n / align) * align; }
292 long int floor( long int n, long int align ) { return (n / align) * align; }
293 unsigned long int floor( unsigned long int n, unsigned long int align ) { return (n / align) * align; }
294 long long int floor( long long int n, long long int align ) { return (n / align) * align; }
295 unsigned long long int floor( unsigned long long int n, unsigned long long int align ) { return (n / align) * align; }
296
297 // forall( T | { T ?/?( T, T ); T ?*?( T, T ); } )
298 // T floor( T n, T align ) { return (n / align) * align; }
299
300 signed char ceiling_div( signed char n, char align ) { return (n + (align - 1hh)) / align; }
301 unsigned char ceiling_div( unsigned char n, unsigned char align ) { return (n + (align - 1hhu)) / align; }
302 short int ceiling_div( short int n, short int align ) { return (n + (align - 1h)) / align; }
303 unsigned short int ceiling_div( unsigned short int n, unsigned short int align ) { return (n + (align - 1hu)) / align; }
304 int ceiling_div( int n, int align ) { return (n + (align - 1n)) / align; }
305 unsigned int ceiling_div( unsigned int n, unsigned int align ) { return (n + (align - 1nu)) / align; }
306 long int ceiling_div( long int n, long int align ) { return (n + (align - 1l)) / align; }
307 unsigned long int ceiling_div( unsigned long int n, unsigned long int align ) { return (n + (align - 1lu)) / align; }
308 long long int ceiling_div( long long int n, long long int align ) { return (n + (align - 1ll)) / align; }
309 unsigned long long int ceiling_div( unsigned long long int n, unsigned long long int align ) { return (n + (align - 1llu)) / align; }
310
311 signed char ceiling( signed char n, char align ) {
312 typeof(n) trunc = floor( n, align );
313 return n < 0 || n == trunc ? trunc : trunc + align;
314 }
315 unsigned char ceiling( unsigned char n, unsigned char align ) {
316 typeof(n) trunc = floor( n, align );
317 return n == trunc ? trunc : trunc + align;
318 }
319 short int ceiling( short int n, short int align ) {
320 typeof(n) trunc = floor( n, align );
321 return n < 0 || n == trunc ? trunc : trunc + align;
322 }
323 unsigned short int ceiling( unsigned short int n, unsigned short int align ) {
324 typeof(n) trunc = floor( n, align );
325 return n == trunc ? trunc : trunc + align;
326 }
327 int ceiling( int n, int align ) {
328 typeof(n) trunc = floor( n, align );
329 return n < 0 || n == trunc ? trunc : trunc + align;
330 }
331 unsigned int ceiling( unsigned int n, unsigned int align ) {
332 typeof(n) trunc = floor( n, align );
333 return n == trunc ? trunc : trunc + align;
334 }
335 long int ceiling( long int n, long int align ) {
336 typeof(n) trunc = floor( n, align );
337 return n < 0 || n == trunc ? trunc : trunc + align;
338 }
339 unsigned long int ceiling( unsigned long int n, unsigned long int align ) {
340 typeof(n) trunc = floor( n, align );
341 return n == trunc ? trunc : trunc + align;
342 }
343 long long int ceiling( long long int n, signed long long int align ) {
344 typeof(n) trunc = floor( n, align );
345 return n < 0 || n == trunc ? trunc : trunc + align;
346 }
347 unsigned long long int ceiling( unsigned long long int n, unsigned long long int align ) {
348 typeof(n) trunc = floor( n, align );
349 return n == trunc ? trunc : trunc + align;
350 }
351
352 float floor( float x ) { return floorf( x ); }
353 // extern "C" { double floor( double ); }
354 long double floor( long double x ) { return floorl( x ); }
355
356 float ceil( float x ) { return ceilf( x ); }
357 // extern "C" { double ceil( double ); }
358 long double ceil( long double x ) { return ceill( x ); }
359
360 float trunc( float x ) { return truncf( x ); }
361 // extern "C" { double trunc( double ); }
362 long double trunc( long double x ) { return truncl( x ); }
363
364 float rint( float x ) { return rintf( x ); }
365 // extern "C" { double rint( double x ); }
366 long double rint( long double x ) { return rintl( x ); }
367 long int rint( float x ) { return lrintf( x ); }
368 long int rint( double x ) { return lrint( x ); }
369 long int rint( long double x ) { return lrintl( x ); }
370 long long int rint( float x ) { return llrintf( x ); }
371 long long int rint( double x ) { return llrint( x ); }
372 long long int rint( long double x ) { return llrintl( x ); }
373
374 long int lrint( float x ) { return lrintf( x ); }
375 // extern "C" { long int lrint( double ); }
376 long int lrint( long double x ) { return lrintl( x ); }
377 long long int llrint( float x ) { return llrintf( x ); }
378 // extern "C" { long long int llrint( double ); }
379 long long int llrint( long double x ) { return llrintl( x ); }
380
381 float nearbyint( float x ) { return nearbyintf( x ); }
382 // extern "C" { double nearbyint( double ); }
383 long double nearbyint( long double x ) { return nearbyintl( x ); }
384
385 float round( float x ) { return roundf( x ); }
386 // extern "C" { double round( double x ); }
387 long double round( long double x ) { return roundl( x ); }
388 long int round( float x ) { return lroundf( x ); }
389 long int round( double x ) { return lround( x ); }
390 long int round( long double x ) { return lroundl( x ); }
391 long long int round( float x ) { return llroundf( x ); }
392 long long int round( double x ) { return llround( x ); }
393 long long int round( long double x ) { return llroundl( x ); }
394
395 long int lround( float x ) { return lroundf( x ); }
396 // extern "C" { long int lround( double ); }
397 long int lround( long double x ) { return lroundl( x ); }
398 long long int llround( float x ) { return llroundf( x ); }
399 // extern "C" { long long int llround( double ); }
400 long long int llround( long double x ) { return llroundl( x ); }
401} // distribution
402
403//---------------------- Manipulation ----------------------
404
405static inline __attribute__((always_inline)) {
406 float copysign( float x, float y ) { return copysignf( x, y ); }
407 // extern "C" { double copysign( double, double ); }
408 long double copysign( long double x, long double y ) { return copysignl( x, y ); }
409
410 float frexp( float x, int * ip ) { return frexpf( x, ip ); }
411 // extern "C" { double frexp( double, int * ); }
412 long double frexp( long double x, int * ip ) { return frexpl( x, ip ); }
413
414 float ldexp( float x, int exp2 ) { return ldexpf( x, exp2 ); }
415 // extern "C" { double ldexp( double, int ); }
416 long double ldexp( long double x, int exp2 ) { return ldexpl( x, exp2 ); }
417
418 [ float, float ] modf( float x ) { float i; x = modff( x, &i ); return [ i, x ]; }
419 float modf( float x, float * i ) { return modff( x, i ); }
420 [ double, double ] modf( double x ) { double i; x = modf( x, &i ); return [ i, x ]; }
421 // extern "C" { double modf( double, double * ); }
422 [ long double, long double ] modf( long double x ) { long double i; x = modfl( x, &i ); return [ i, x ]; }
423 long double modf( long double x, long double * i ) { return modfl( x, i ); }
424
425 float nextafter( float x, float y ) { return nextafterf( x, y ); }
426 // extern "C" { double nextafter( double, double ); }
427 long double nextafter( long double x, long double y ) { return nextafterl( x, y ); }
428
429 float nexttoward( float x, long double y ) { return nexttowardf( x, y ); }
430 // extern "C" { double nexttoward( double, long double ); }
431 long double nexttoward( long double x, long double y ) { return nexttowardl( x, y ); }
432
433 float scalbn( float x, int exp ) { return scalbnf( x, exp ); }
434 // extern "C" { double scalbn( double, int ); }
435 long double scalbn( long double x, int exp ) { return scalbnl( x, exp ); }
436 float scalbn( float x, long int exp ) { return scalblnf( x, exp ); }
437 double scalbn( double x, long int exp ) { return scalbln( x, exp ); }
438 long double scalbn( long double x, long int exp ) { return scalblnl( x, exp ); }
439
440 float scalbln( float x, long int exp ) { return scalblnf( x, exp ); }
441 // extern "C" { double scalbln( double, long int ); }
442 long double scalbln( long double x, long int exp ) { return scalblnl( x, exp ); }
443} // distribution
444
445//---------------------------------------
446
447static inline __attribute__((always_inline)) {
448 forall( T | { void ?{}( T &, one_t ); T ?+?( T, T ); T ?-?( T, T );T ?*?( T, T ); } )
449 T lerp( T x, T y, T a ) { return x * ((T){1} - a) + y * a; }
450
451 forall( T | { void ?{}( T &, zero_t ); void ?{}( T &, one_t ); int ?<?( T, T ); } )
452 T step( T edge, T x ) { return x < edge ? (T){0} : (T){1}; }
453
454 forall( T | { void ?{}( T &, int ); T clamp( T, T, T ); T ?-?( T, T ); T ?*?( T, T ); T ?/?( T, T ); } )
455 T smoothstep( T edge0, T edge1, T x ) { T t = clamp( (x - edge0) / (edge1 - edge0), (T){0}, (T){1} ); return t * t * ((T){3} - (T){2} * t); }
456} // distribution
457
458// Local Variables: //
459// mode: c //
460// tab-width: 4 //
461// End: //
Note: See TracBrowser for help on using the repository browser.