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