1 #ifndef _TGMATH_H 2 #define _TGMATH_H 3 4 /* 5 the return types are only correct with gcc (__GNUC__) 6 otherwise they are long double or long double complex 7 8 the long double version of a function is never chosen when 9 sizeof(double) == sizeof(long double) 10 (but the return type is set correctly with gcc) 11 */ 12 13 #include <math.h> 14 #include <complex.h> 15 16 #define __IS_FP(x) (sizeof((x)+1ULL) == sizeof((x)+1.0f)) 17 #define __IS_CX(x) (__IS_FP(x) && sizeof(x) == sizeof((x)+I)) 18 #define __IS_REAL(x) (__IS_FP(x) && 2*sizeof(x) == sizeof((x)+I)) 19 20 #define __FLT(x) (__IS_REAL(x) && sizeof(x) == sizeof(float)) 21 #define __LDBL(x) (__IS_REAL(x) && sizeof(x) == sizeof(long double) && sizeof(long double) != sizeof(double)) 22 23 #define __FLTCX(x) (__IS_CX(x) && sizeof(x) == sizeof(float complex)) 24 #define __DBLCX(x) (__IS_CX(x) && sizeof(x) == sizeof(double complex)) 25 #define __LDBLCX(x) (__IS_CX(x) && sizeof(x) == sizeof(long double complex) && sizeof(long double) != sizeof(double)) 26 27 /* return type */ 28 29 #ifdef __GNUC__ 30 /* 31 the result must be casted to the right type 32 (otherwise the result type is determined by the conversion 33 rules applied to all the function return types so it is long 34 double or long double complex except for integral functions) 35 36 this cannot be done in c99, so the typeof gcc extension is 37 used and that the type of ?: depends on wether an operand is 38 a null pointer constant or not 39 (in c11 _Generic can be used) 40 41 the c arguments below must be integer constant expressions 42 so they can be in null pointer constants 43 (__IS_FP above was carefully chosen this way) 44 */ 45 /* if c then t else void */ 46 #define __type1(c,t) __typeof__(*(0?(t*)0:(void*)!(c))) 47 /* if c then t1 else t2 */ 48 #define __type2(c,t1,t2) __typeof__(*(0?(__type1(c,t1)*)0:(__type1(!(c),t2)*)0)) 49 /* cast to double when x is integral, otherwise use typeof(x) */ 50 #define __RETCAST(x) ( \ 51 __type2(__IS_FP(x), __typeof__(x), double)) 52 /* 2 args case, should work for complex types (cpow) */ 53 #define __RETCAST_2(x, y) ( \ 54 __type2(__IS_FP(x) && __IS_FP(y), \ 55 __typeof__((x)+(y)), \ 56 __typeof__((x)+(y)+1.0))) 57 /* 3 args case (fma only) */ 58 #define __RETCAST_3(x, y, z) ( \ 59 __type2(__IS_FP(x) && __IS_FP(y) && __IS_FP(z), \ 60 __typeof__((x)+(y)+(z)), \ 61 __typeof__((x)+(y)+(z)+1.0))) 62 /* drop complex from the type of x */ 63 /* TODO: wrong when sizeof(long double)==sizeof(double) */ 64 #define __RETCAST_REAL(x) ( \ 65 __type2(__IS_FP(x) && sizeof((x)+I) == sizeof(float complex), float, \ 66 __type2(sizeof((x)+1.0+I) == sizeof(double complex), double, \ 67 long double))) 68 /* add complex to the type of x */ 69 #define __RETCAST_CX(x) (__typeof__(__RETCAST(x)0+I)) 70 #else 71 #define __RETCAST(x) 72 #define __RETCAST_2(x, y) 73 #define __RETCAST_3(x, y, z) 74 #define __RETCAST_REAL(x) 75 #define __RETCAST_CX(x) 76 #endif 77 78 /* function selection */ 79 80 #define __tg_real_nocast(fun, x) ( \ 81 __FLT(x) ? fun ## f (x) : \ 82 __LDBL(x) ? fun ## l (x) : \ 83 fun(x) ) 84 85 #define __tg_real(fun, x) (__RETCAST(x)__tg_real_nocast(fun, x)) 86 87 #define __tg_real_2_1(fun, x, y) (__RETCAST(x)( \ 88 __FLT(x) ? fun ## f (x, y) : \ 89 __LDBL(x) ? fun ## l (x, y) : \ 90 fun(x, y) )) 91 92 #define __tg_real_2(fun, x, y) (__RETCAST_2(x, y)( \ 93 __FLT(x) && __FLT(y) ? fun ## f (x, y) : \ 94 __LDBL((x)+(y)) ? fun ## l (x, y) : \ 95 fun(x, y) )) 96 97 #define __tg_complex(fun, x) (__RETCAST_CX(x)( \ 98 __FLTCX((x)+I) && __IS_FP(x) ? fun ## f (x) : \ 99 __LDBLCX((x)+I) ? fun ## l (x) : \ 100 fun(x) )) 101 102 #define __tg_complex_retreal(fun, x) (__RETCAST_REAL(x)( \ 103 __FLTCX((x)+I) && __IS_FP(x) ? fun ## f (x) : \ 104 __LDBLCX((x)+I) ? fun ## l (x) : \ 105 fun(x) )) 106 107 #define __tg_real_complex(fun, x) (__RETCAST(x)( \ 108 __FLTCX(x) ? c ## fun ## f (x) : \ 109 __DBLCX(x) ? c ## fun (x) : \ 110 __LDBLCX(x) ? c ## fun ## l (x) : \ 111 __FLT(x) ? fun ## f (x) : \ 112 __LDBL(x) ? fun ## l (x) : \ 113 fun(x) )) 114 115 /* special cases */ 116 117 #define __tg_real_remquo(x, y, z) (__RETCAST_2(x, y)( \ 118 __FLT(x) && __FLT(y) ? remquof(x, y, z) : \ 119 __LDBL((x)+(y)) ? remquol(x, y, z) : \ 120 remquo(x, y, z) )) 121 122 #define __tg_real_fma(x, y, z) (__RETCAST_3(x, y, z)( \ 123 __FLT(x) && __FLT(y) && __FLT(z) ? fmaf(x, y, z) : \ 124 __LDBL((x)+(y)+(z)) ? fmal(x, y, z) : \ 125 fma(x, y, z) )) 126 127 #define __tg_real_complex_pow(x, y) (__RETCAST_2(x, y)( \ 128 __FLTCX((x)+(y)) && __IS_FP(x) && __IS_FP(y) ? cpowf(x, y) : \ 129 __FLTCX((x)+(y)) ? cpow(x, y) : \ 130 __DBLCX((x)+(y)) ? cpow(x, y) : \ 131 __LDBLCX((x)+(y)) ? cpowl(x, y) : \ 132 __FLT(x) && __FLT(y) ? powf(x, y) : \ 133 __LDBL((x)+(y)) ? powl(x, y) : \ 134 pow(x, y) )) 135 136 #define __tg_real_complex_fabs(x) (__RETCAST_REAL(x)( \ 137 __FLTCX(x) ? cabsf(x) : \ 138 __DBLCX(x) ? cabs(x) : \ 139 __LDBLCX(x) ? cabsl(x) : \ 140 __FLT(x) ? fabsf(x) : \ 141 __LDBL(x) ? fabsl(x) : \ 142 fabs(x) )) 143 144 /* suppress any macros in math.h or complex.h */ 145 146 #undef acos 147 #undef acosh 148 #undef asin 149 #undef asinh 150 #undef atan 151 #undef atan2 152 #undef atanh 153 #undef carg 154 #undef cbrt 155 #undef ceil 156 #undef cimag 157 #undef conj 158 #undef copysign 159 #undef cos 160 #undef cosh 161 #undef cproj 162 #undef creal 163 #undef erf 164 #undef erfc 165 #undef exp 166 #undef exp2 167 #undef expm1 168 #undef fabs 169 #undef fdim 170 #undef floor 171 #undef fma 172 #undef fmax 173 #undef fmin 174 #undef fmod 175 #undef frexp 176 #undef hypot 177 #undef ilogb 178 #undef ldexp 179 #undef lgamma 180 #undef llrint 181 #undef llround 182 #undef log 183 #undef log10 184 #undef log1p 185 #undef log2 186 #undef logb 187 #undef lrint 188 #undef lround 189 #undef nearbyint 190 #undef nextafter 191 #undef nexttoward 192 #undef pow 193 #undef remainder 194 #undef remquo 195 #undef rint 196 #undef round 197 #undef scalbln 198 #undef scalbn 199 #undef sin 200 #undef sinh 201 #undef sqrt 202 #undef tan 203 #undef tanh 204 #undef tgamma 205 #undef trunc 206 207 /* tg functions */ 208 209 #define acos(x) __tg_real_complex(acos, (x)) 210 #define acosh(x) __tg_real_complex(acosh, (x)) 211 #define asin(x) __tg_real_complex(asin, (x)) 212 #define asinh(x) __tg_real_complex(asinh, (x)) 213 #define atan(x) __tg_real_complex(atan, (x)) 214 #define atan2(x,y) __tg_real_2(atan2, (x), (y)) 215 #define atanh(x) __tg_real_complex(atanh, (x)) 216 #define carg(x) __tg_complex_retreal(carg, (x)) 217 #define cbrt(x) __tg_real(cbrt, (x)) 218 #define ceil(x) __tg_real(ceil, (x)) 219 #define cimag(x) __tg_complex_retreal(cimag, (x)) 220 #define conj(x) __tg_complex(conj, (x)) 221 #define copysign(x,y) __tg_real_2(copysign, (x), (y)) 222 #define cos(x) __tg_real_complex(cos, (x)) 223 #define cosh(x) __tg_real_complex(cosh, (x)) 224 #define cproj(x) __tg_complex(cproj, (x)) 225 #define creal(x) __tg_complex_retreal(creal, (x)) 226 #define erf(x) __tg_real(erf, (x)) 227 #define erfc(x) __tg_real(erfc, (x)) 228 #define exp(x) __tg_real_complex(exp, (x)) 229 #define exp2(x) __tg_real(exp2, (x)) 230 #define expm1(x) __tg_real(expm1, (x)) 231 #define fabs(x) __tg_real_complex_fabs(x) 232 #define fdim(x,y) __tg_real_2(fdim, (x), (y)) 233 #define floor(x) __tg_real(floor, (x)) 234 #define fma(x,y,z) __tg_real_fma((x), (y), (z)) 235 #define fmax(x,y) __tg_real_2(fmax, (x), (y)) 236 #define fmin(x,y) __tg_real_2(fmin, (x), (y)) 237 #define fmod(x,y) __tg_real_2(fmod, (x), (y)) 238 #define frexp(x,y) __tg_real_2_1(frexp, (x), (y)) 239 #define hypot(x,y) __tg_real_2(hypot, (x), (y)) 240 #define ilogb(x) __tg_real_nocast(ilogb, (x)) 241 #define ldexp(x,y) __tg_real_2_1(ldexp, (x), (y)) 242 #define lgamma(x) __tg_real(lgamma, (x)) 243 #define llrint(x) __tg_real_nocast(llrint, (x)) 244 #define llround(x) __tg_real_nocast(llround, (x)) 245 #define log(x) __tg_real_complex(log, (x)) 246 #define log10(x) __tg_real(log10, (x)) 247 #define log1p(x) __tg_real(log1p, (x)) 248 #define log2(x) __tg_real(log2, (x)) 249 #define logb(x) __tg_real(logb, (x)) 250 #define lrint(x) __tg_real_nocast(lrint, (x)) 251 #define lround(x) __tg_real_nocast(lround, (x)) 252 #define nearbyint(x) __tg_real(nearbyint, (x)) 253 #define nextafter(x,y) __tg_real_2(nextafter, (x), (y)) 254 #define nexttoward(x,y) __tg_real_2(nexttoward, (x), (y)) 255 #define pow(x,y) __tg_real_complex_pow((x), (y)) 256 #define remainder(x,y) __tg_real_2(remainder, (x), (y)) 257 #define remquo(x,y,z) __tg_real_remquo((x), (y), (z)) 258 #define rint(x) __tg_real(rint, (x)) 259 #define round(x) __tg_real(round, (x)) 260 #define scalbln(x,y) __tg_real_2_1(scalbln, (x), (y)) 261 #define scalbn(x,y) __tg_real_2_1(scalbn, (x), (y)) 262 #define sin(x) __tg_real_complex(sin, (x)) 263 #define sinh(x) __tg_real_complex(sinh, (x)) 264 #define sqrt(x) __tg_real_complex(sqrt, (x)) 265 #define tan(x) __tg_real_complex(tan, (x)) 266 #define tanh(x) __tg_real_complex(tanh, (x)) 267 #define tgamma(x) __tg_real(tgamma, (x)) 268 #define trunc(x) __tg_real(trunc, (x)) 269 270 #endif 271