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 #define __RETCAST_REAL(x) ( \ 64 __type2(__IS_FP(x) && sizeof((x)+I) == sizeof(float complex), float, \ 65 __type2(sizeof((x)+1.0+I) == sizeof(double complex), double, \ 66 long double))) 67 /* add complex to the type of x */ 68 #define __RETCAST_CX(x) (__typeof__(__RETCAST(x)0+I)) 69 #else 70 #define __RETCAST(x) 71 #define __RETCAST_2(x, y) 72 #define __RETCAST_3(x, y, z) 73 #define __RETCAST_REAL(x) 74 #define __RETCAST_CX(x) 75 #endif 76 77 /* function selection */ 78 79 #define __tg_real_nocast(fun, x) ( \ 80 __FLT(x) ? fun ## f (x) : \ 81 __LDBL(x) ? fun ## l (x) : \ 82 fun(x) ) 83 84 #define __tg_real(fun, x) (__RETCAST(x)__tg_real_nocast(fun, x)) 85 86 #define __tg_real_2_1(fun, x, y) (__RETCAST(x)( \ 87 __FLT(x) ? fun ## f (x, y) : \ 88 __LDBL(x) ? fun ## l (x, y) : \ 89 fun(x, y) )) 90 91 #define __tg_real_2(fun, x, y) (__RETCAST_2(x, y)( \ 92 __FLT(x) && __FLT(y) ? fun ## f (x, y) : \ 93 __LDBL((x)+(y)) ? fun ## l (x, y) : \ 94 fun(x, y) )) 95 96 #define __tg_complex(fun, x) (__RETCAST_CX(x)( \ 97 __FLTCX((x)+I) && __IS_FP(x) ? fun ## f (x) : \ 98 __LDBLCX((x)+I) ? fun ## l (x) : \ 99 fun(x) )) 100 101 #define __tg_complex_retreal(fun, x) (__RETCAST_REAL(x)( \ 102 __FLTCX((x)+I) && __IS_FP(x) ? fun ## f (x) : \ 103 __LDBLCX((x)+I) ? fun ## l (x) : \ 104 fun(x) )) 105 106 #define __tg_real_complex(fun, x) (__RETCAST(x)( \ 107 __FLTCX(x) ? c ## fun ## f (x) : \ 108 __DBLCX(x) ? c ## fun (x) : \ 109 __LDBLCX(x) ? c ## fun ## l (x) : \ 110 __FLT(x) ? fun ## f (x) : \ 111 __LDBL(x) ? fun ## l (x) : \ 112 fun(x) )) 113 114 /* special cases */ 115 116 #define __tg_real_remquo(x, y, z) (__RETCAST_2(x, y)( \ 117 __FLT(x) && __FLT(y) ? remquof(x, y, z) : \ 118 __LDBL((x)+(y)) ? remquol(x, y, z) : \ 119 remquo(x, y, z) )) 120 121 #define __tg_real_fma(x, y, z) (__RETCAST_3(x, y, z)( \ 122 __FLT(x) && __FLT(y) && __FLT(z) ? fmaf(x, y, z) : \ 123 __LDBL((x)+(y)+(z)) ? fmal(x, y, z) : \ 124 fma(x, y, z) )) 125 126 #define __tg_real_complex_pow(x, y) (__RETCAST_2(x, y)( \ 127 __FLTCX((x)+(y)) && __IS_FP(x) && __IS_FP(y) ? cpowf(x, y) : \ 128 __FLTCX((x)+(y)) ? cpow(x, y) : \ 129 __DBLCX((x)+(y)) ? cpow(x, y) : \ 130 __LDBLCX((x)+(y)) ? cpowl(x, y) : \ 131 __FLT(x) && __FLT(y) ? powf(x, y) : \ 132 __LDBL((x)+(y)) ? powl(x, y) : \ 133 pow(x, y) )) 134 135 #define __tg_real_complex_fabs(x) (__RETCAST_REAL(x)( \ 136 __FLTCX(x) ? cabsf(x) : \ 137 __DBLCX(x) ? cabs(x) : \ 138 __LDBLCX(x) ? cabsl(x) : \ 139 __FLT(x) ? fabsf(x) : \ 140 __LDBL(x) ? fabsl(x) : \ 141 fabs(x) )) 142 143 /* suppress any macros in math.h or complex.h */ 144 145 #undef acos 146 #undef acosh 147 #undef asin 148 #undef asinh 149 #undef atan 150 #undef atan2 151 #undef atanh 152 #undef carg 153 #undef cbrt 154 #undef ceil 155 #undef cimag 156 #undef conj 157 #undef copysign 158 #undef cos 159 #undef cosh 160 #undef cproj 161 #undef creal 162 #undef erf 163 #undef erfc 164 #undef exp 165 #undef exp2 166 #undef expm1 167 #undef fabs 168 #undef fdim 169 #undef floor 170 #undef fma 171 #undef fmax 172 #undef fmin 173 #undef fmod 174 #undef frexp 175 #undef hypot 176 #undef ilogb 177 #undef ldexp 178 #undef lgamma 179 #undef llrint 180 #undef llround 181 #undef log 182 #undef log10 183 #undef log1p 184 #undef log2 185 #undef logb 186 #undef lrint 187 #undef lround 188 #undef nearbyint 189 #undef nextafter 190 #undef nexttoward 191 #undef pow 192 #undef remainder 193 #undef remquo 194 #undef rint 195 #undef round 196 #undef scalbln 197 #undef scalbn 198 #undef sin 199 #undef sinh 200 #undef sqrt 201 #undef tan 202 #undef tanh 203 #undef tgamma 204 #undef trunc 205 206 /* tg functions */ 207 208 #define acos(x) __tg_real_complex(acos, (x)) 209 #define acosh(x) __tg_real_complex(acosh, (x)) 210 #define asin(x) __tg_real_complex(asin, (x)) 211 #define asinh(x) __tg_real_complex(asinh, (x)) 212 #define atan(x) __tg_real_complex(atan, (x)) 213 #define atan2(x,y) __tg_real_2(atan2, (x), (y)) 214 #define atanh(x) __tg_real_complex(atanh, (x)) 215 #define carg(x) __tg_complex_retreal(carg, (x)) 216 #define cbrt(x) __tg_real(cbrt, (x)) 217 #define ceil(x) __tg_real(ceil, (x)) 218 #define cimag(x) __tg_complex_retreal(cimag, (x)) 219 #define conj(x) __tg_complex(conj, (x)) 220 #define copysign(x,y) __tg_real_2(copysign, (x), (y)) 221 #define cos(x) __tg_real_complex(cos, (x)) 222 #define cosh(x) __tg_real_complex(cosh, (x)) 223 #define cproj(x) __tg_complex(cproj, (x)) 224 #define creal(x) __tg_complex_retreal(creal, (x)) 225 #define erf(x) __tg_real(erf, (x)) 226 #define erfc(x) __tg_real(erfc, (x)) 227 #define exp(x) __tg_real_complex(exp, (x)) 228 #define exp2(x) __tg_real(exp2, (x)) 229 #define expm1(x) __tg_real(expm1, (x)) 230 #define fabs(x) __tg_real_complex_fabs(x) 231 #define fdim(x,y) __tg_real_2(fdim, (x), (y)) 232 #define floor(x) __tg_real(floor, (x)) 233 #define fma(x,y,z) __tg_real_fma((x), (y), (z)) 234 #define fmax(x,y) __tg_real_2(fmax, (x), (y)) 235 #define fmin(x,y) __tg_real_2(fmin, (x), (y)) 236 #define fmod(x,y) __tg_real_2(fmod, (x), (y)) 237 #define frexp(x,y) __tg_real_2_1(frexp, (x), (y)) 238 #define hypot(x,y) __tg_real_2(hypot, (x), (y)) 239 #define ilogb(x) __tg_real_nocast(ilogb, (x)) 240 #define ldexp(x,y) __tg_real_2_1(ldexp, (x), (y)) 241 #define lgamma(x) __tg_real(lgamma, (x)) 242 #define llrint(x) __tg_real_nocast(llrint, (x)) 243 #define llround(x) __tg_real_nocast(llround, (x)) 244 #define log(x) __tg_real_complex(log, (x)) 245 #define log10(x) __tg_real(log10, (x)) 246 #define log1p(x) __tg_real(log1p, (x)) 247 #define log2(x) __tg_real(log2, (x)) 248 #define logb(x) __tg_real(logb, (x)) 249 #define lrint(x) __tg_real_nocast(lrint, (x)) 250 #define lround(x) __tg_real_nocast(lround, (x)) 251 #define nearbyint(x) __tg_real(nearbyint, (x)) 252 #define nextafter(x,y) __tg_real_2(nextafter, (x), (y)) 253 #define nexttoward(x,y) __tg_real_2(nexttoward, (x), (y)) 254 #define pow(x,y) __tg_real_complex_pow((x), (y)) 255 #define remainder(x,y) __tg_real_2(remainder, (x), (y)) 256 #define remquo(x,y,z) __tg_real_remquo((x), (y), (z)) 257 #define rint(x) __tg_real(rint, (x)) 258 #define round(x) __tg_real(round, (x)) 259 #define scalbln(x,y) __tg_real_2_1(scalbln, (x), (y)) 260 #define scalbn(x,y) __tg_real_2_1(scalbn, (x), (y)) 261 #define sin(x) __tg_real_complex(sin, (x)) 262 #define sinh(x) __tg_real_complex(sinh, (x)) 263 #define sqrt(x) __tg_real_complex(sqrt, (x)) 264 #define tan(x) __tg_real_complex(tan, (x)) 265 #define tanh(x) __tg_real_complex(tanh, (x)) 266 #define tgamma(x) __tg_real(tgamma, (x)) 267 #define trunc(x) __tg_real(trunc, (x)) 268 269 #endif 270