/* Copyright (C) 2014-2020, Dirk Krause SPDX-License-Identifier: BSD-3-Clause */ /* WARNING: This file was generated by the dkct program (see http://dktools.sourceforge.net/ for details). Changes you make here will be lost if dkct is run again! You should modify the original source and run dkct on it. Original source: dk3mad.ctr */ /** @file dk3mad.c The dk3mad module. */ #line 161 "dk3mad.ctr" #include "dk3conf.h" #include "dk3ma.h" #include "dk3const.h" #line 169 "dk3mad.ctr" #if (DK3_SIZEOF_DOUBLE == 8) && (DK3_HAVE_IEEE_754_DOUBLE) /** Type to construct double value from binary (hex) data. */ typedef union { unsigned char c[sizeof(double)]; /**< Bytes to specify hex. */ double d; /**< Double value to retrieve. */ } __dk3ma_max_double_t; /** Maximum double value. */ __dk3ma_max_double_t const __dk3ma_max_double = { #if DK3_WORDS_BIGENDIAN { 0x7F, 0xEF, 0xFF, 0xFF, 0xFF, 0xFF, 0xFF, 0xFF } #else { 0xFF, 0xFF, 0xFF, 0xFF, 0xFF, 0xFF, 0xEF, 0x7F } #endif }; /** Maximum double value for internal use. */ #define dk3ma_i_max_double __dk3ma_max_double.d #else #ifndef DK3_MAX_DOUBLE #define DK3_MAX_DOUBLE 1.7e308 #endif /** Maximum double value for internal use. */ #define dk3ma_i_max_double DK3_MAX_DOUBLE #endif double dk3ma_d_abs(double x) { #if DK3_HAVE_FABS return (fabs(x)); #else return ((0.0 <= x) ? x : (0.0 - x)); #endif } double dk3ma_d_rint(double x) { #if DK3_HAVE_RINT return (rint(x)); #else return (floor(x + 0.5)); #endif } double dk3ma_d_add_ok(double a, double b, int *ec) { #line 231 "dk3mad.ctr" if (NULL != ec) { if ((0.0 < a) && (0.0 < b)) { if ((dk3ma_i_max_double - a) < b) { #line 235 "dk3mad.ctr" *ec = DK3_ERROR_MATH_OVERFLOW; } } else { if ((0.0 > a) && (0.0 > b)) { if ( ( (-1.0 * dk3ma_i_max_double) - a ) > b ) { #line 241 "dk3mad.ctr" *ec = DK3_ERROR_MATH_OVERFLOW; } } } } #line 246 "dk3mad.ctr" return (a + b); } double dk3ma_d_sub_ok(double a, double b, int *ec) { #line 255 "dk3mad.ctr" if (NULL != ec) { if ((0.0 < a) && (0.0 > b)) { if( (dk3ma_i_max_double + b) < a ) { #line 259 "dk3mad.ctr" *ec = DK3_ERROR_MATH_OVERFLOW; } } else { if ((0.0 > a) && (0.0 < b)) { if ( ( (-1.0 * dk3ma_i_max_double) + b ) > a ) { #line 265 "dk3mad.ctr" *ec = DK3_ERROR_MATH_OVERFLOW; } } } } #line 270 "dk3mad.ctr" return (a - b); } double dk3ma_d_mul_ok(double a, double b, int *ec) { #line 279 "dk3mad.ctr" if (ec) { if (dk3ma_d_abs(a) > 1.0) { if ( (dk3ma_i_max_double / dk3ma_d_abs(a)) < dk3ma_d_abs(b) ) { #line 283 "dk3mad.ctr" *ec = DK3_ERROR_MATH_OVERFLOW; } } } #line 287 "dk3mad.ctr" return (a * b); } double dk3ma_d_div_ok(double a, double b, int *ec) { #if DK3_HAVE_FPCLASSIFY && defined(FP_ZERO) if (FP_ZERO == fpclassify(b)) { if (NULL != ec) { *ec = DK3_ERROR_MATH_DIVZERO; } if (FP_ZERO == fpclassify(a)) { return (NAN); } else { if (0.0 <= a) { return (INFINITY); } else { return (-INFINITY); } } } else { if (dk3ma_d_abs(b) >= 1.0) { return (a / b); } else { if ((dk3ma_i_max_double * dk3ma_d_abs(b)) < dk3ma_d_abs(a)) { if (ec) { #line 313 "dk3mad.ctr" *ec = DK3_ERROR_MATH_OVERFLOW; } if (((0.0 <= a) && (0.0 <= b)) || ((0.0 > a) && (0.0 > b))) { return (INFINITY); } else { return (-INFINITY); } } else { return (a / b); } } } #else #if DK3_ON_WINDOWS switch(_fpclass(b)) { case _FPCLASS_NZ: case _FPCLASS_PZ: { if (ec) { #line 331 "dk3mad.ctr" *ec = DK3_ERROR_MATH_DIVZERO; } switch(_fpclass(a)) { case _FPCLASS_NZ: case _FPCLASS_PZ: { #if defined(NAN) return (NAN); #else return (dk3ma_i_max_double); #endif } break; default: { if (0.0 <= a) { return (HUGE_VAL); } else { return (-HUGE_VAL); } } break; } } break; default: { if (dk3ma_d_abs(b) >= 1.0) { return (a / b); } else { if ((dk3ma_i_max_double * dk3ma_d_abs(b)) < dk3ma_d_abs(a)) { if (ec) { #line 357 "dk3mad.ctr" *ec = DK3_ERROR_MATH_OVERFLOW; } if (((0.0 <= a) && (0.0 <= b)) || ((0.0 > a) && (0.0 > b))) { return (HUGE_VAL); } else { return (-HUGE_VAL); } } else { return (a / b); } } } break; } #else if (dk3ma_d_abs(b) >= 1.0) { return (a / b); } else { if ((dk3ma_i_max_double * dk3ma_d_abs(b)) < dk3ma_d_abs(a)) { if (ec) { #line 376 "dk3mad.ctr" *ec = DK3_ERROR_MATH_OVERFLOW; } if (((0.0 <= a) && (0.0 <= b)) || ((0.0 > a) && (0.0 > b))) { #ifdef INFINITY return (INFINITY); #else #ifdef HUGE_VAL return (HUGE_VAL); #else return (dk3ma_i_max_double); #endif #endif } else { #ifdef INFINITY return (-INFINITY); #else #ifdef HUGE_VAL return (-HUGE_VAL); #else return (-1.0 * dk3ma_i_max_double); #endif #endif } } else { return (a / b); } } #endif #endif } double dk3ma_d_square_ok(double x, int *ec) { #if 0 double back = -1.0; #line 415 "dk3mad.ctr" if (0.0 <= x) { back = sqrt(x); } else { if (NULL != ec) { #line 419 "dk3mad.ctr" *ec = DK3_ERROR_MATH_OUT_OF_RANGE; } } #line 422 "dk3mad.ctr" return back; #else return (dk3ma_d_mul_ok(x, x, ec)); #endif } int dk3ma_d_equal(double a, double b, double epsilon) { int back = 0; int ec = 0; if (dk3ma_d_abs(dk3ma_d_sub_ok(a, b, &ec)) < epsilon) { if (0 == ec) { back = 1; } } return back; } double dk3ma_d_atan2(double y, double x) { #if DK3_HAVE_ATAN2 double back; #line 451 "dk3mad.ctr" back = atan2(y, x); #line 452 "dk3mad.ctr" while (back < 0.0) { back += (2.0 * M_PI); } while (back > (2.0 * M_PI)) { back -= (2.0 * M_PI); } #line 455 "dk3mad.ctr" return back; #else double back = -5.0 * M_PI; double v; #line 460 "dk3mad.ctr" v = dk3ma_d_div_ok(y, x, &mec); if (mec) { #line 462 "dk3mad.ctr" if(y < 0.0) { #line 463 "dk3mad.ctr" back = 1.5 * M_PI; } else { #line 465 "dk3mad.ctr" back = 0.5 * M_PI; } } else { #line 468 "dk3mad.ctr" back = atan(v); if (x < 0.0) { #line 470 "dk3mad.ctr" back += M_PI; } } #line 473 "dk3mad.ctr" while (back < 0.0) { back += (2.0 * M_PI); } while (back > (2.0 * M_PI)) { back -= (2.0 * M_PI); } #line 476 "dk3mad.ctr" return back; #endif } double dk3ma_d_restrict_digits(double x, size_t n) { double back; double newval; size_t i; size_t mult; int ec = 0; back = x; mult = 0; /* Multiplications */ for (i = 0; ((0 == ec) && (i < n)); i++) { newval = dk3ma_d_mul_ok(back, 10.0, &ec); if (0 == ec) { back = newval; mult++; } } /* Rounding */ back = dk3ma_d_rint(back); /* Divisions */ for (i = 0; i < mult; i++) { back = back / 10.0; } return back; } double dk3ma_d_restrict_digits_ceil(double x, size_t n) { double back; double newval; size_t i; size_t mult; int ec = 0; back = x; mult = 0; /* Multiplications */ for (i = 0; ((0 == ec) && (i < n)); i++) { newval = dk3ma_d_mul_ok(back, 10.0, &ec); if (0 == ec) { back = newval; mult++; } } /* Rounding */ back = ceil(back); /* Divisions */ for (i = 0; i < mult; i++) { back = back / 10.0; } return back; } double dk3ma_d_restrict_digits_floor(double x, size_t n) { double back; double newval; size_t i; size_t mult; int ec = 0; back = x; mult = 0; /* Multiplications */ for (i = 0; ((0 == ec) && (i < n)); i++) { newval = dk3ma_d_mul_ok(back, 10.0, &ec); if (0 == ec) { back = newval; mult++; } } /* Rounding */ back = floor(back); /* Divisions */ for (i = 0; i < mult; i++) { back = back / 10.0; } return back; }