1 /* 2 3 PETSc mathematics include file. Defines certain basic mathematical 4 constants and functions for working with single, double, and quad precision 5 floating point numbers as well as complex single and double. 6 7 This file is included by petscsys.h and should not be used directly. 8 9 */ 10 11 #if !defined(__PETSCMATH_H) 12 #define __PETSCMATH_H 13 #include <math.h> 14 15 /* 16 17 Defines operations that are different for complex and real numbers; 18 note that one cannot mix the use of complex and real in the same 19 PETSc program. All PETSc objects in one program are built around the object 20 PetscScalar which is either always a real or a complex. 21 22 */ 23 24 #if defined(PETSC_USE_REAL_SINGLE) 25 #define MPIU_REAL MPI_FLOAT 26 typedef float PetscReal; 27 #define PetscSqrtReal(a) sqrt(a) 28 #define PetscExpReal(a) exp(a) 29 #define PetscLogReal(a) log(a) 30 #define PetscLog10Real(a) log10(a) 31 #ifdef PETSC_HAVE_LOG2 32 #define PetscLog2Real(a) log2(a) 33 #endif 34 #define PetscSinReal(a) sin(a) 35 #define PetscCosReal(a) cos(a) 36 #define PetscTanReal(a) tan(a) 37 #define PetscAsinReal(a) asin(a) 38 #define PetscAcosReal(a) acos(a) 39 #define PetscAtanReal(a) atan(a) 40 #define PetscAtan2Real(a,b) atan2(a,b) 41 #define PetscSinhReal(a) sinh(a) 42 #define PetscCoshReal(a) cosh(a) 43 #define PetscTanhReal(a) tanh(a) 44 #define PetscPowReal(a,b) pow(a,b) 45 #define PetscCeilReal(a) ceil(a) 46 #define PetscFloorReal(a) floor(a) 47 #define PetscFmodReal(a,b) fmod(a,b) 48 #define PetscTGamma(a) tgammaf(a) 49 #elif defined(PETSC_USE_REAL_DOUBLE) 50 #define MPIU_REAL MPI_DOUBLE 51 typedef double PetscReal; 52 #define PetscSqrtReal(a) sqrt(a) 53 #define PetscExpReal(a) exp(a) 54 #define PetscLogReal(a) log(a) 55 #define PetscLog10Real(a) log10(a) 56 #ifdef PETSC_HAVE_LOG2 57 #define PetscLog2Real(a) log2(a) 58 #endif 59 #define PetscSinReal(a) sin(a) 60 #define PetscCosReal(a) cos(a) 61 #define PetscTanReal(a) tan(a) 62 #define PetscAsinReal(a) asin(a) 63 #define PetscAcosReal(a) acos(a) 64 #define PetscAtanReal(a) atan(a) 65 #define PetscAtan2Real(a,b) atan2(a,b) 66 #define PetscSinhReal(a) sinh(a) 67 #define PetscCoshReal(a) cosh(a) 68 #define PetscTanhReal(a) tanh(a) 69 #define PetscPowReal(a,b) pow(a,b) 70 #define PetscCeilReal(a) ceil(a) 71 #define PetscFloorReal(a) floor(a) 72 #define PetscFmodReal(a,b) fmod(a,b) 73 #define PetscTGamma(a) tgamma(a) 74 #elif defined(PETSC_USE_REAL___FLOAT128) 75 #if defined(__cplusplus) 76 extern "C" { 77 #endif 78 #include <quadmath.h> 79 #if defined(__cplusplus) 80 } 81 #endif 82 PETSC_EXTERN MPI_Datatype MPIU___FLOAT128 PetscAttrMPITypeTag(__float128); 83 #define MPIU_REAL MPIU___FLOAT128 84 typedef __float128 PetscReal; 85 #define PetscSqrtReal(a) sqrtq(a) 86 #define PetscExpReal(a) expq(a) 87 #define PetscLogReal(a) logq(a) 88 #define PetscLog10Real(a) log10q(a) 89 #ifdef PETSC_HAVE_LOG2 90 #define PetscLog2Real(a) log2q(a) 91 #endif 92 #define PetscSinReal(a) sinq(a) 93 #define PetscCosReal(a) cosq(a) 94 #define PetscTanReal(a) tanq(a) 95 #define PetscAsinReal(a) asinq(a) 96 #define PetscAcosReal(a) acosq(a) 97 #define PetscAtanReal(a) atanq(a) 98 #define PetscAtan2Real(a,b) atan2q(a,b) 99 #define PetscSinhReal(a) sinhq(a) 100 #define PetscCoshReal(a) coshq(a) 101 #define PetscTanhReal(a) tanhq(a) 102 #define PetscPowReal(a,b) powq(a,b) 103 #define PetscCeilReal(a) ceilq(a) 104 #define PetscFloorReal(a) floorq(a) 105 #define PetscFmodReal(a,b) fmodq(a,b) 106 #define PetscTGamma(a) tgammaq(a) 107 #endif /* PETSC_USE_REAL_* */ 108 109 /* 110 Complex number definitions 111 */ 112 #if defined(__cplusplus) && defined(PETSC_HAVE_CXX_COMPLEX) && !defined(PETSC_USE_REAL___FLOAT128) 113 #if !defined(PETSC_SKIP_COMPLEX) 114 #define PETSC_HAVE_COMPLEX 1 115 /* C++ support of complex number */ 116 #if defined(PETSC_HAVE_CUSP) 117 #define complexlib cusp 118 #include <cusp/complex.h> 119 #elif defined(PETSC_HAVE_VECCUDA) && __CUDACC_VER_MAJOR__ > 6 120 /* complex headers in thrust only available in CUDA 7.0 and above */ 121 #define complexlib thrust 122 #include <thrust/complex.h> 123 #else 124 #define complexlib std 125 #include <complex> 126 #endif 127 128 #define PetscRealPartComplex(a) (a).real() 129 #define PetscImaginaryPartComplex(a) (a).imag() 130 #define PetscAbsComplex(a) complexlib::abs(a) 131 #define PetscConjComplex(a) complexlib::conj(a) 132 #define PetscSqrtComplex(a) complexlib::sqrt(a) 133 #define PetscPowComplex(a,b) complexlib::pow(a,b) 134 #define PetscExpComplex(a) complexlib::exp(a) 135 #define PetscLogComplex(a) complexlib::log(a) 136 #define PetscSinComplex(a) complexlib::sin(a) 137 #define PetscCosComplex(a) complexlib::cos(a) 138 #define PetscAsinComplex(a) complexlib::asin(a) 139 #define PetscAcosComplex(a) complexlib::acos(a) 140 #if defined(PETSC_HAVE_TANCOMPLEX) 141 #define PetscTanComplex(a) complexlib::tan(a) 142 #else 143 #define PetscTanComplex(a) PetscSinComplex(a)/PetscCosComplex(a) 144 #endif 145 #define PetscSinhComplex(a) complexlib::sinh(a) 146 #define PetscCoshComplex(a) complexlib::cosh(a) 147 #if defined(PETSC_HAVE_TANHCOMPLEX) 148 #define PetscTanhComplex(a) complexlib::tanh(a) 149 #else 150 #define PetscTanhComplex(a) PetscSinhComplex(a)/PetscCoshComplex(a) 151 #endif 152 153 #if defined(PETSC_USE_REAL_SINGLE) 154 typedef complexlib::complex<float> PetscComplex; 155 #if defined(PETSC_USE_CXX_COMPLEX_FLOAT_WORKAROUND) 156 static inline PetscComplex operator+(const PetscComplex& lhs, const double& rhs) { return lhs + float(rhs); } 157 static inline PetscComplex operator+(const double& lhs, const PetscComplex& rhs) { return float(lhs) + rhs; } 158 static inline PetscComplex operator-(const PetscComplex& lhs, const double& rhs) { return lhs - float(rhs); } 159 static inline PetscComplex operator-(const double& lhs, const PetscComplex& rhs) { return float(lhs) - rhs; } 160 static inline PetscComplex operator*(const PetscComplex& lhs, const double& rhs) { return lhs * float(rhs); } 161 static inline PetscComplex operator*(const double& lhs, const PetscComplex& rhs) { return float(lhs) * rhs; } 162 static inline PetscComplex operator/(const PetscComplex& lhs, const double& rhs) { return lhs / float(rhs); } 163 static inline PetscComplex operator/(const double& lhs, const PetscComplex& rhs) { return float(lhs) / rhs; } 164 static inline bool operator==(const PetscComplex& lhs, const double& rhs) { return lhs.imag() == float(0) && lhs.real() == float(rhs); } 165 static inline bool operator==(const double& lhs, const PetscComplex& rhs) { return rhs.imag() == float(0) && rhs.real() == float(lhs); } 166 static inline bool operator!=(const PetscComplex& lhs, const double& rhs) { return lhs.imag() != float(0) || lhs.real() != float(rhs); } 167 static inline bool operator!=(const double& lhs, const PetscComplex& rhs) { return rhs.imag() != float(0) || rhs.real() != float(lhs); } 168 #endif /* PETSC_USE_CXX_COMPLEX_FLOAT_WORKAROUND */ 169 #elif defined(PETSC_USE_REAL_DOUBLE) 170 typedef complexlib::complex<double> PetscComplex; 171 #elif defined(PETSC_USE_REAL___FLOAT128) 172 typedef complexlib::complex<__float128> PetscComplex; /* Notstandard and not expected to work, use __complex128 */ 173 PETSC_EXTERN MPI_Datatype MPIU___COMPLEX128; 174 #endif /* PETSC_USE_REAL_ */ 175 #endif /* ! PETSC_SKIP_COMPLEX */ 176 177 #elif !defined(__cplusplus) && defined(PETSC_HAVE_C99_COMPLEX) 178 #if !defined(PETSC_SKIP_COMPLEX) 179 #define PETSC_HAVE_COMPLEX 1 180 #include <complex.h> 181 182 #if defined(PETSC_USE_REAL_SINGLE) 183 typedef float _Complex PetscComplex; 184 185 #define PetscRealPartComplex(a) crealf(a) 186 #define PetscImaginaryPartComplex(a) cimagf(a) 187 #define PetscAbsComplex(a) cabsf(a) 188 #define PetscConjComplex(a) conjf(a) 189 #define PetscSqrtComplex(a) csqrtf(a) 190 #define PetscPowComplex(a,b) cpowf(a,b) 191 #define PetscExpComplex(a) cexpf(a) 192 #define PetscLogComplex(a) clogf(a) 193 #define PetscSinComplex(a) csinf(a) 194 #define PetscCosComplex(a) ccosf(a) 195 #define PetscAsinComplex(a) casinf(a) 196 #define PetscAcosComplex(a) cacosf(a) 197 #define PetscTanComplex(a) ctanf(a) 198 #define PetscSinhComplex(a) csinhf(a) 199 #define PetscCoshComplex(a) ccoshf(a) 200 #define PetscTanhComplex(a) ctanhf(a) 201 202 #elif defined(PETSC_USE_REAL_DOUBLE) 203 typedef double _Complex PetscComplex; 204 205 #define PetscRealPartComplex(a) creal(a) 206 #define PetscImaginaryPartComplex(a) cimag(a) 207 #define PetscAbsComplex(a) cabs(a) 208 #define PetscConjComplex(a) conj(a) 209 #define PetscSqrtComplex(a) csqrt(a) 210 #define PetscPowComplex(a,b) cpow(a,b) 211 #define PetscExpComplex(a) cexp(a) 212 #define PetscLogComplex(a) clog(a) 213 #define PetscSinComplex(a) csin(a) 214 #define PetscCosComplex(a) ccos(a) 215 #define PetscAsinComplex(a) casin(a) 216 #define PetscAcosComplex(a) cacos(a) 217 #define PetscTanComplex(a) ctan(a) 218 #define PetscSinhComplex(a) csinh(a) 219 #define PetscCoshComplex(a) ccosh(a) 220 #define PetscTanhComplex(a) ctanh(a) 221 222 #elif defined(PETSC_USE_REAL___FLOAT128) 223 typedef __complex128 PetscComplex; 224 PETSC_EXTERN MPI_Datatype MPIU___COMPLEX128 PetscAttrMPITypeTag(__complex128); 225 226 #define PetscRealPartComplex(a) crealq(a) 227 #define PetscImaginaryPartComplex(a) cimagq(a) 228 #define PetscAbsComplex(a) cabsq(a) 229 #define PetscConjComplex(a) conjq(a) 230 #define PetscSqrtComplex(a) csqrtq(a) 231 #define PetscPowComplex(a,b) cpowq(a,b) 232 #define PetscExpComplex(a) cexpq(a) 233 #define PetscLogComplex(a) clogq(a) 234 #define PetscSinComplex(a) csinq(a) 235 #define PetscCosComplex(a) ccosq(a) 236 #define PetscAsinComplex(a) casinq(a) 237 #define PetscAcosComplex(a) cacosq(a) 238 #define PetscTanComplex(a) ctanq(a) 239 #define PetscSinhComplex(a) csinhq(a) 240 #define PetscCoshComplex(a) ccoshq(a) 241 #define PetscTanhComplex(a) ctanhq(a) 242 243 #endif /* PETSC_USE_REAL_* */ 244 #elif (defined(PETSC_USE_COMPLEX) && !defined(PETSC_SKIP_COMPLEX)) 245 #error "PETSc was configured --with-scalar-type=complex, but a language-appropriate complex library is not available" 246 #endif /* !PETSC_SKIP_COMPLEX */ 247 #endif /* (__cplusplus && PETSC_HAVE_CXX_COMPLEX) else-if (!__cplusplus && PETSC_HAVE_C99_COMPLEX) */ 248 249 #if defined(PETSC_HAVE_COMPLEX) 250 #if defined(PETSC_HAVE_MPI_C_DOUBLE_COMPLEX) 251 #define MPIU_C_DOUBLE_COMPLEX MPI_C_DOUBLE_COMPLEX 252 #define MPIU_C_COMPLEX MPI_C_COMPLEX 253 #else 254 # if defined(__cplusplus) && defined(PETSC_HAVE_CXX_COMPLEX) 255 typedef complexlib::complex<double> petsc_mpiu_c_double_complex; 256 typedef complexlib::complex<float> petsc_mpiu_c_complex; 257 # elif !defined(__cplusplus) && defined(PETSC_HAVE_C99_COMPLEX) 258 typedef double _Complex petsc_mpiu_c_double_complex; 259 typedef float _Complex petsc_mpiu_c_complex; 260 # else 261 typedef struct {double real,imag;} petsc_mpiu_c_double_complex; 262 typedef struct {float real,imag;} petsc_mpiu_c_complex; 263 # endif 264 PETSC_EXTERN MPI_Datatype MPIU_C_DOUBLE_COMPLEX PetscAttrMPITypeTagLayoutCompatible(petsc_mpiu_c_double_complex); 265 PETSC_EXTERN MPI_Datatype MPIU_C_COMPLEX PetscAttrMPITypeTagLayoutCompatible(petsc_mpiu_c_complex); 266 #endif /* PETSC_HAVE_MPI_C_DOUBLE_COMPLEX */ 267 #endif /* PETSC_HAVE_COMPLEX */ 268 269 #if defined(PETSC_HAVE_COMPLEX) 270 # if defined(PETSC_USE_REAL_SINGLE) 271 # define MPIU_COMPLEX MPIU_C_COMPLEX 272 # elif defined(PETSC_USE_REAL_DOUBLE) 273 # define MPIU_COMPLEX MPIU_C_DOUBLE_COMPLEX 274 # elif defined(PETSC_USE_REAL___FLOAT128) 275 # define MPIU_COMPLEX MPIU___COMPLEX128 276 # endif /* PETSC_USE_REAL_* */ 277 #endif 278 279 #if (defined(PETSC_USE_COMPLEX) && !defined(PETSC_SKIP_COMPLEX)) 280 typedef PetscComplex PetscScalar; 281 #define PetscRealPart(a) PetscRealPartComplex(a) 282 #define PetscImaginaryPart(a) PetscImaginaryPartComplex(a) 283 #define PetscAbsScalar(a) PetscAbsComplex(a) 284 #define PetscConj(a) PetscConjComplex(a) 285 #define PetscSqrtScalar(a) PetscSqrtComplex(a) 286 #define PetscPowScalar(a,b) PetscPowComplex(a,b) 287 #define PetscExpScalar(a) PetscExpComplex(a) 288 #define PetscLogScalar(a) PetscLogComplex(a) 289 #define PetscSinScalar(a) PetscSinComplex(a) 290 #define PetscCosScalar(a) PetscCosComplex(a) 291 #define PetscAsinScalar(a) PetscAsinComplex(a) 292 #define PetscAcosScalar(a) PetscAcosComplex(a) 293 #define PetscTanScalar(a) PetscTanComplex(a) 294 #define PetscSinhScalar(a) PetscSinhComplex(a) 295 #define PetscCoshScalar(a) PetscCoshComplex(a) 296 #define PetscTanhScalar(a) PetscTanhComplex(a) 297 #define MPIU_SCALAR MPIU_COMPLEX 298 299 /* 300 real number definitions 301 */ 302 #else /* PETSC_USE_COMPLEX */ 303 typedef PetscReal PetscScalar; 304 #define MPIU_SCALAR MPIU_REAL 305 306 #define PetscRealPart(a) (a) 307 #define PetscImaginaryPart(a) ((PetscReal)0.) 308 PETSC_STATIC_INLINE PetscReal PetscAbsScalar(PetscScalar a) {return a < 0.0 ? -a : a;} 309 #define PetscConj(a) (a) 310 #if !defined(PETSC_USE_REAL___FLOAT128) 311 #define PetscSqrtScalar(a) sqrt(a) 312 #define PetscPowScalar(a,b) pow(a,b) 313 #define PetscExpScalar(a) exp(a) 314 #define PetscLogScalar(a) log(a) 315 #define PetscSinScalar(a) sin(a) 316 #define PetscCosScalar(a) cos(a) 317 #define PetscAsinScalar(a) asin(a) 318 #define PetscAcosScalar(a) acos(a) 319 #define PetscTanScalar(a) tan(a) 320 #define PetscSinhScalar(a) sinh(a) 321 #define PetscCoshScalar(a) cosh(a) 322 #define PetscTanhScalar(a) tanh(a) 323 #else /* PETSC_USE_REAL___FLOAT128 */ 324 #define PetscSqrtScalar(a) sqrtq(a) 325 #define PetscPowScalar(a,b) powq(a,b) 326 #define PetscExpScalar(a) expq(a) 327 #define PetscLogScalar(a) logq(a) 328 #define PetscSinScalar(a) sinq(a) 329 #define PetscCosScalar(a) cosq(a) 330 #define PetscAsinScalar(a) asinq(a) 331 #define PetscAcosScalar(a) acosq(a) 332 #define PetscTanScalar(a) tanq(a) 333 #define PetscSinhScalar(a) sinhq(a) 334 #define PetscCoshScalar(a) coshq(a) 335 #define PetscTanhScalar(a) tanhq(a) 336 #endif /* PETSC_USE_REAL___FLOAT128 */ 337 338 #endif /* PETSC_USE_COMPLEX */ 339 340 #define PetscSign(a) (((a) >= 0) ? ((a) == 0 ? 0 : 1) : -1) 341 #define PetscAbs(a) (((a) >= 0) ? (a) : -(a)) 342 343 /* --------------------------------------------------------------------------*/ 344 345 /* 346 Certain objects may be created using either single or double precision. 347 This is currently not used. 348 */ 349 typedef enum { PETSC_SCALAR_DOUBLE,PETSC_SCALAR_SINGLE, PETSC_SCALAR_LONG_DOUBLE } PetscScalarPrecision; 350 351 #if defined(PETSC_HAVE_COMPLEX) 352 /* PETSC_i is the imaginary number, i */ 353 PETSC_EXTERN PetscComplex PETSC_i; 354 #endif 355 356 /*MC 357 PetscMin - Returns minimum of two numbers 358 359 Synopsis: 360 #include <petscmath.h> 361 type PetscMin(type v1,type v2) 362 363 Not Collective 364 365 Input Parameter: 366 + v1 - first value to find minimum of 367 - v2 - second value to find minimum of 368 369 Notes: type can be integer or floating point value 370 371 Level: beginner 372 373 .seealso: PetscMin(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr() 374 375 M*/ 376 #define PetscMin(a,b) (((a)<(b)) ? (a) : (b)) 377 378 /*MC 379 PetscMax - Returns maxium of two numbers 380 381 Synopsis: 382 #include <petscmath.h> 383 type max PetscMax(type v1,type v2) 384 385 Not Collective 386 387 Input Parameter: 388 + v1 - first value to find maximum of 389 - v2 - second value to find maximum of 390 391 Notes: type can be integer or floating point value 392 393 Level: beginner 394 395 .seealso: PetscMin(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr() 396 397 M*/ 398 #define PetscMax(a,b) (((a)<(b)) ? (b) : (a)) 399 400 /*MC 401 PetscClipInterval - Returns a number clipped to be within an interval 402 403 Synopsis: 404 #include <petscmath.h> 405 type clip PetscClipInterval(type x,type a,type b) 406 407 Not Collective 408 409 Input Parameter: 410 + x - value to use if within interval (a,b) 411 . a - lower end of interval 412 - b - upper end of interval 413 414 Notes: type can be integer or floating point value 415 416 Level: beginner 417 418 .seealso: PetscMin(), PetscMax(), PetscAbsInt(), PetscAbsReal(), PetscSqr() 419 420 M*/ 421 #define PetscClipInterval(x,a,b) (PetscMax((a),PetscMin((x),(b)))) 422 423 /*MC 424 PetscAbsInt - Returns the absolute value of an integer 425 426 Synopsis: 427 #include <petscmath.h> 428 int abs PetscAbsInt(int v1) 429 430 Not Collective 431 432 Input Parameter: 433 . v1 - the integer 434 435 Level: beginner 436 437 .seealso: PetscMax(), PetscMin(), PetscAbsReal(), PetscSqr() 438 439 M*/ 440 #define PetscAbsInt(a) (((a)<0) ? -(a) : (a)) 441 442 /*MC 443 PetscAbsReal - Returns the absolute value of an real number 444 445 Synopsis: 446 #include <petscmath.h> 447 Real abs PetscAbsReal(PetscReal v1) 448 449 Not Collective 450 451 Input Parameter: 452 . v1 - the double 453 454 455 Level: beginner 456 457 .seealso: PetscMax(), PetscMin(), PetscAbsInt(), PetscSqr() 458 459 M*/ 460 #define PetscAbsReal(a) (((a)<0) ? -(a) : (a)) 461 462 /*MC 463 PetscSqr - Returns the square of a number 464 465 Synopsis: 466 #include <petscmath.h> 467 type sqr PetscSqr(type v1) 468 469 Not Collective 470 471 Input Parameter: 472 . v1 - the value 473 474 Notes: type can be integer or floating point value 475 476 Level: beginner 477 478 .seealso: PetscMax(), PetscMin(), PetscAbsInt(), PetscAbsReal() 479 480 M*/ 481 #define PetscSqr(a) ((a)*(a)) 482 483 /* ----------------------------------------------------------------------------*/ 484 /* 485 Basic constants 486 */ 487 #if defined(PETSC_USE_REAL___FLOAT128) 488 #define PETSC_PI M_PIq 489 #elif defined(M_PI) 490 #define PETSC_PI M_PI 491 #else 492 #define PETSC_PI 3.14159265358979323846264338327950288419716939937510582 493 #endif 494 495 #if !defined(PETSC_USE_64BIT_INDICES) 496 #define PETSC_MAX_INT 2147483647 497 #define PETSC_MIN_INT (-PETSC_MAX_INT - 1) 498 #else 499 #define PETSC_MAX_INT 9223372036854775807L 500 #define PETSC_MIN_INT (-PETSC_MAX_INT - 1) 501 #endif 502 503 #if defined(PETSC_USE_REAL_SINGLE) 504 # define PETSC_MAX_REAL 3.40282346638528860e+38F 505 # define PETSC_MIN_REAL -PETSC_MAX_REAL 506 # define PETSC_MACHINE_EPSILON 1.19209290e-07F 507 # define PETSC_SQRT_MACHINE_EPSILON 3.45266983e-04F 508 # define PETSC_SMALL 1.e-5 509 #elif defined(PETSC_USE_REAL_DOUBLE) 510 # define PETSC_MAX_REAL 1.7976931348623157e+308 511 # define PETSC_MIN_REAL -PETSC_MAX_REAL 512 # define PETSC_MACHINE_EPSILON 2.2204460492503131e-16 513 # define PETSC_SQRT_MACHINE_EPSILON 1.490116119384766e-08 514 # define PETSC_SMALL 1.e-10 515 #elif defined(PETSC_USE_REAL___FLOAT128) 516 # define PETSC_MAX_REAL FLT128_MAX 517 # define PETSC_MIN_REAL -FLT128_MAX 518 # define PETSC_MACHINE_EPSILON FLT128_EPSILON 519 # define PETSC_SQRT_MACHINE_EPSILON 1.38777878078e-17q 520 # define PETSC_SMALL 1.e-20q 521 #endif 522 523 #define PETSC_INFINITY PETSC_MAX_REAL/4.0 524 #define PETSC_NINFINITY -PETSC_INFINITY 525 526 PETSC_EXTERN PetscErrorCode PetscIsInfOrNanReal(PetscReal); 527 PETSC_EXTERN PetscErrorCode PetscIsNanReal(PetscReal); 528 PETSC_EXTERN PetscBool PetscIsNormalReal(PetscReal); 529 PETSC_STATIC_INLINE PetscErrorCode PetscIsInfOrNanScalar(PetscScalar v) {return PetscIsInfOrNanReal(PetscAbsScalar(v));} 530 PETSC_STATIC_INLINE PetscErrorCode PetscIsNanScalar(PetscScalar v) {return PetscIsNanReal(PetscAbsScalar(v));} 531 PETSC_STATIC_INLINE PetscErrorCode PetscIsNormalScalar(PetscScalar v) {return PetscIsNormalReal(PetscAbsScalar(v));} 532 533 /* 534 These macros are currently hardwired to match the regular data types, so there is no support for a different 535 MatScalar from PetscScalar. We left the MatScalar in the source just in case we use it again. 536 */ 537 #define MPIU_MATSCALAR MPIU_SCALAR 538 typedef PetscScalar MatScalar; 539 typedef PetscReal MatReal; 540 541 struct petsc_mpiu_2scalar {PetscScalar a,b;}; 542 PETSC_EXTERN MPI_Datatype MPIU_2SCALAR PetscAttrMPITypeTagLayoutCompatible(struct petsc_mpiu_2scalar); 543 #if defined(PETSC_USE_64BIT_INDICES) || !defined(MPI_2INT) 544 struct petsc_mpiu_2int {PetscInt a,b;}; 545 PETSC_EXTERN MPI_Datatype MPIU_2INT PetscAttrMPITypeTagLayoutCompatible(struct petsc_mpiu_2int); 546 #else 547 #define MPIU_2INT MPI_2INT 548 #endif 549 550 PETSC_STATIC_INLINE PetscInt PetscPowInt(PetscInt base,PetscInt power) 551 { 552 PetscInt result = 1; 553 while (power) { 554 if (power & 1) result *= base; 555 power >>= 1; 556 base *= base; 557 } 558 return result; 559 } 560 561 PETSC_STATIC_INLINE PetscReal PetscPowRealInt(PetscReal base,PetscInt power) 562 { 563 PetscReal result = 1; 564 if (power < 0) { 565 power = -power; 566 base = ((PetscReal)1)/base; 567 } 568 while (power) { 569 if (power & 1) result *= base; 570 power >>= 1; 571 base *= base; 572 } 573 return result; 574 } 575 576 PETSC_STATIC_INLINE PetscScalar PetscPowScalarInt(PetscScalar base,PetscInt power) 577 { 578 PetscScalar result = 1; 579 if (power < 0) { 580 power = -power; 581 base = ((PetscReal)1)/base; 582 } 583 while (power) { 584 if (power & 1) result *= base; 585 power >>= 1; 586 base *= base; 587 } 588 return result; 589 } 590 591 PETSC_STATIC_INLINE PetscScalar PetscPowScalarReal(PetscScalar base,PetscReal power) 592 { 593 PetscScalar cpower = power; 594 return PetscPowScalar(base,cpower); 595 } 596 597 #ifndef PETSC_HAVE_LOG2 598 PETSC_STATIC_INLINE PetscReal PetscLog2Real(PetscReal n) 599 { 600 return PetscLogReal(n)/PetscLogReal(2.0); 601 } 602 #endif 603 #endif 604