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 #if !defined(PETSCMATH_H) 11 #define PETSCMATH_H 12 13 #include <math.h> 14 #include <petscmacros.h> 15 #include <petscsystypes.h> 16 17 /* SUBMANSEC = Sys */ 18 19 /* 20 21 Defines operations that are different for complex and real numbers. 22 All PETSc objects in one program are built around the object 23 PetscScalar which is either always a real or a complex. 24 25 */ 26 27 /* 28 Real number definitions 29 */ 30 #if defined(PETSC_USE_REAL_SINGLE) 31 #define PetscSqrtReal(a) sqrtf(a) 32 #define PetscCbrtReal(a) cbrtf(a) 33 #define PetscHypotReal(a, b) hypotf(a, b) 34 #define PetscAtan2Real(a, b) atan2f(a, b) 35 #define PetscPowReal(a, b) powf(a, b) 36 #define PetscExpReal(a) expf(a) 37 #define PetscLogReal(a) logf(a) 38 #define PetscLog10Real(a) log10f(a) 39 #define PetscLog2Real(a) log2f(a) 40 #define PetscSinReal(a) sinf(a) 41 #define PetscCosReal(a) cosf(a) 42 #define PetscTanReal(a) tanf(a) 43 #define PetscAsinReal(a) asinf(a) 44 #define PetscAcosReal(a) acosf(a) 45 #define PetscAtanReal(a) atanf(a) 46 #define PetscSinhReal(a) sinhf(a) 47 #define PetscCoshReal(a) coshf(a) 48 #define PetscTanhReal(a) tanhf(a) 49 #define PetscAsinhReal(a) asinhf(a) 50 #define PetscAcoshReal(a) acoshf(a) 51 #define PetscAtanhReal(a) atanhf(a) 52 #define PetscErfReal(a) erff(a) 53 #define PetscCeilReal(a) ceilf(a) 54 #define PetscFloorReal(a) floorf(a) 55 #define PetscFmodReal(a, b) fmodf(a, b) 56 #define PetscCopysignReal(a, b) copysignf(a, b) 57 #define PetscTGamma(a) tgammaf(a) 58 #if defined(PETSC_HAVE_LGAMMA_IS_GAMMA) 59 #define PetscLGamma(a) gammaf(a) 60 #else 61 #define PetscLGamma(a) lgammaf(a) 62 #endif 63 64 #elif defined(PETSC_USE_REAL_DOUBLE) 65 #define PetscSqrtReal(a) sqrt(a) 66 #define PetscCbrtReal(a) cbrt(a) 67 #define PetscHypotReal(a, b) hypot(a, b) 68 #define PetscAtan2Real(a, b) atan2(a, b) 69 #define PetscPowReal(a, b) pow(a, b) 70 #define PetscExpReal(a) exp(a) 71 #define PetscLogReal(a) log(a) 72 #define PetscLog10Real(a) log10(a) 73 #define PetscLog2Real(a) log2(a) 74 #define PetscSinReal(a) sin(a) 75 #define PetscCosReal(a) cos(a) 76 #define PetscTanReal(a) tan(a) 77 #define PetscAsinReal(a) asin(a) 78 #define PetscAcosReal(a) acos(a) 79 #define PetscAtanReal(a) atan(a) 80 #define PetscSinhReal(a) sinh(a) 81 #define PetscCoshReal(a) cosh(a) 82 #define PetscTanhReal(a) tanh(a) 83 #define PetscAsinhReal(a) asinh(a) 84 #define PetscAcoshReal(a) acosh(a) 85 #define PetscAtanhReal(a) atanh(a) 86 #define PetscErfReal(a) erf(a) 87 #define PetscCeilReal(a) ceil(a) 88 #define PetscFloorReal(a) floor(a) 89 #define PetscFmodReal(a, b) fmod(a, b) 90 #define PetscCopysignReal(a, b) copysign(a, b) 91 #define PetscTGamma(a) tgamma(a) 92 #if defined(PETSC_HAVE_LGAMMA_IS_GAMMA) 93 #define PetscLGamma(a) gamma(a) 94 #else 95 #define PetscLGamma(a) lgamma(a) 96 #endif 97 98 #elif defined(PETSC_USE_REAL___FLOAT128) 99 #define PetscSqrtReal(a) sqrtq(a) 100 #define PetscCbrtReal(a) cbrtq(a) 101 #define PetscHypotReal(a, b) hypotq(a, b) 102 #define PetscAtan2Real(a, b) atan2q(a, b) 103 #define PetscPowReal(a, b) powq(a, b) 104 #define PetscExpReal(a) expq(a) 105 #define PetscLogReal(a) logq(a) 106 #define PetscLog10Real(a) log10q(a) 107 #define PetscLog2Real(a) log2q(a) 108 #define PetscSinReal(a) sinq(a) 109 #define PetscCosReal(a) cosq(a) 110 #define PetscTanReal(a) tanq(a) 111 #define PetscAsinReal(a) asinq(a) 112 #define PetscAcosReal(a) acosq(a) 113 #define PetscAtanReal(a) atanq(a) 114 #define PetscSinhReal(a) sinhq(a) 115 #define PetscCoshReal(a) coshq(a) 116 #define PetscTanhReal(a) tanhq(a) 117 #define PetscAsinhReal(a) asinhq(a) 118 #define PetscAcoshReal(a) acoshq(a) 119 #define PetscAtanhReal(a) atanhq(a) 120 #define PetscErfReal(a) erfq(a) 121 #define PetscCeilReal(a) ceilq(a) 122 #define PetscFloorReal(a) floorq(a) 123 #define PetscFmodReal(a, b) fmodq(a, b) 124 #define PetscCopysignReal(a, b) copysignq(a, b) 125 #define PetscTGamma(a) tgammaq(a) 126 #if defined(PETSC_HAVE_LGAMMA_IS_GAMMA) 127 #define PetscLGamma(a) gammaq(a) 128 #else 129 #define PetscLGamma(a) lgammaq(a) 130 #endif 131 132 #elif defined(PETSC_USE_REAL___FP16) 133 #define PetscSqrtReal(a) sqrtf(a) 134 #define PetscCbrtReal(a) cbrtf(a) 135 #define PetscHypotReal(a, b) hypotf(a, b) 136 #define PetscAtan2Real(a, b) atan2f(a, b) 137 #define PetscPowReal(a, b) powf(a, b) 138 #define PetscExpReal(a) expf(a) 139 #define PetscLogReal(a) logf(a) 140 #define PetscLog10Real(a) log10f(a) 141 #define PetscLog2Real(a) log2f(a) 142 #define PetscSinReal(a) sinf(a) 143 #define PetscCosReal(a) cosf(a) 144 #define PetscTanReal(a) tanf(a) 145 #define PetscAsinReal(a) asinf(a) 146 #define PetscAcosReal(a) acosf(a) 147 #define PetscAtanReal(a) atanf(a) 148 #define PetscSinhReal(a) sinhf(a) 149 #define PetscCoshReal(a) coshf(a) 150 #define PetscTanhReal(a) tanhf(a) 151 #define PetscAsinhReal(a) asinhf(a) 152 #define PetscAcoshReal(a) acoshf(a) 153 #define PetscAtanhReal(a) atanhf(a) 154 #define PetscErfReal(a) erff(a) 155 #define PetscCeilReal(a) ceilf(a) 156 #define PetscFloorReal(a) floorf(a) 157 #define PetscFmodReal(a, b) fmodf(a, b) 158 #define PetscCopySignReal(a, b) copysignf(a, b) 159 #define PetscTGamma(a) tgammaf(a) 160 #if defined(PETSC_HAVE_LGAMMA_IS_GAMMA) 161 #define PetscLGamma(a) gammaf(a) 162 #else 163 #define PetscLGamma(a) lgammaf(a) 164 #endif 165 166 #endif /* PETSC_USE_REAL_* */ 167 168 static inline PetscReal PetscSignReal(PetscReal a) { 169 return (PetscReal)((a < (PetscReal)0) ? -1 : ((a > (PetscReal)0) ? 1 : 0)); 170 } 171 172 #if !defined(PETSC_HAVE_LOG2) 173 #undef PetscLog2Real 174 static inline PetscReal PetscLog2Real(PetscReal a) { 175 return PetscLogReal(a) / PetscLogReal((PetscReal)2); 176 } 177 #endif 178 179 #if defined(PETSC_HAVE_REAL___FLOAT128) 180 PETSC_EXTERN MPI_Datatype MPIU___FLOAT128 PETSC_ATTRIBUTE_MPI_TYPE_TAG(__float128); 181 #endif 182 #if defined(PETSC_HAVE_REAL___FP16) 183 PETSC_EXTERN MPI_Datatype MPIU___FP16 PETSC_ATTRIBUTE_MPI_TYPE_TAG(__fp16); 184 #endif 185 186 /*MC 187 MPIU_REAL - Portable MPI datatype corresponding to `PetscReal` independent of what precision `PetscReal` is in 188 189 Notes: 190 In MPI calls that require an MPI datatype that matches a `PetscReal` or array of `PetscReal` values, pass this value. 191 192 Level: beginner 193 194 .seealso: `PetscReal`, `PetscScalar`, `PetscComplex`, `PetscInt`, `MPIU_SCALAR`, `MPIU_COMPLEX`, `MPIU_INT` 195 M*/ 196 #if defined(PETSC_USE_REAL_SINGLE) 197 #define MPIU_REAL MPI_FLOAT 198 #elif defined(PETSC_USE_REAL_DOUBLE) 199 #define MPIU_REAL MPI_DOUBLE 200 #elif defined(PETSC_USE_REAL___FLOAT128) 201 #define MPIU_REAL MPIU___FLOAT128 202 #elif defined(PETSC_USE_REAL___FP16) 203 #define MPIU_REAL MPIU___FP16 204 #endif /* PETSC_USE_REAL_* */ 205 206 /* 207 Complex number definitions 208 */ 209 #if defined(PETSC_HAVE_COMPLEX) 210 #if defined(__cplusplus) && !defined(PETSC_USE_REAL___FLOAT128) 211 /* C++ support of complex number */ 212 213 #define PetscRealPartComplex(a) (a).real() 214 #define PetscImaginaryPartComplex(a) (a).imag() 215 #define PetscAbsComplex(a) petsccomplexlib::abs(a) 216 #define PetscArgComplex(a) petsccomplexlib::arg(a) 217 #define PetscConjComplex(a) petsccomplexlib::conj(a) 218 #define PetscSqrtComplex(a) petsccomplexlib::sqrt(a) 219 #define PetscPowComplex(a, b) petsccomplexlib::pow(a, b) 220 #define PetscExpComplex(a) petsccomplexlib::exp(a) 221 #define PetscLogComplex(a) petsccomplexlib::log(a) 222 #define PetscSinComplex(a) petsccomplexlib::sin(a) 223 #define PetscCosComplex(a) petsccomplexlib::cos(a) 224 #define PetscTanComplex(a) petsccomplexlib::tan(a) 225 #define PetscAsinComplex(a) petsccomplexlib::asin(a) 226 #define PetscAcosComplex(a) petsccomplexlib::acos(a) 227 #define PetscAtanComplex(a) petsccomplexlib::atan(a) 228 #define PetscSinhComplex(a) petsccomplexlib::sinh(a) 229 #define PetscCoshComplex(a) petsccomplexlib::cosh(a) 230 #define PetscTanhComplex(a) petsccomplexlib::tanh(a) 231 #define PetscAsinhComplex(a) petsccomplexlib::asinh(a) 232 #define PetscAcoshComplex(a) petsccomplexlib::acosh(a) 233 #define PetscAtanhComplex(a) petsccomplexlib::atanh(a) 234 235 /* TODO: Add configure tests 236 237 #if !defined(PETSC_HAVE_CXX_TAN_COMPLEX) 238 #undef PetscTanComplex 239 static inline PetscComplex PetscTanComplex(PetscComplex z) 240 { 241 return PetscSinComplex(z)/PetscCosComplex(z); 242 } 243 #endif 244 245 #if !defined(PETSC_HAVE_CXX_TANH_COMPLEX) 246 #undef PetscTanhComplex 247 static inline PetscComplex PetscTanhComplex(PetscComplex z) 248 { 249 return PetscSinhComplex(z)/PetscCoshComplex(z); 250 } 251 #endif 252 253 #if !defined(PETSC_HAVE_CXX_ASIN_COMPLEX) 254 #undef PetscAsinComplex 255 static inline PetscComplex PetscAsinComplex(PetscComplex z) 256 { 257 const PetscComplex j(0,1); 258 return -j*PetscLogComplex(j*z+PetscSqrtComplex(1.0f-z*z)); 259 } 260 #endif 261 262 #if !defined(PETSC_HAVE_CXX_ACOS_COMPLEX) 263 #undef PetscAcosComplex 264 static inline PetscComplex PetscAcosComplex(PetscComplex z) 265 { 266 const PetscComplex j(0,1); 267 return j*PetscLogComplex(z-j*PetscSqrtComplex(1.0f-z*z)); 268 } 269 #endif 270 271 #if !defined(PETSC_HAVE_CXX_ATAN_COMPLEX) 272 #undef PetscAtanComplex 273 static inline PetscComplex PetscAtanComplex(PetscComplex z) 274 { 275 const PetscComplex j(0,1); 276 return 0.5f*j*PetscLogComplex((1.0f-j*z)/(1.0f+j*z)); 277 } 278 #endif 279 280 #if !defined(PETSC_HAVE_CXX_ASINH_COMPLEX) 281 #undef PetscAsinhComplex 282 static inline PetscComplex PetscAsinhComplex(PetscComplex z) 283 { 284 return PetscLogComplex(z+PetscSqrtComplex(z*z+1.0f)); 285 } 286 #endif 287 288 #if !defined(PETSC_HAVE_CXX_ACOSH_COMPLEX) 289 #undef PetscAcoshComplex 290 static inline PetscComplex PetscAcoshComplex(PetscComplex z) 291 { 292 return PetscLogComplex(z+PetscSqrtComplex(z*z-1.0f)); 293 } 294 #endif 295 296 #if !defined(PETSC_HAVE_CXX_ATANH_COMPLEX) 297 #undef PetscAtanhComplex 298 static inline PetscComplex PetscAtanhComplex(PetscComplex z) 299 { 300 return 0.5f*PetscLogComplex((1.0f+z)/(1.0f-z)); 301 } 302 #endif 303 304 */ 305 306 #else /* C99 support of complex number */ 307 308 #if defined(PETSC_USE_REAL_SINGLE) 309 #define PetscRealPartComplex(a) crealf(a) 310 #define PetscImaginaryPartComplex(a) cimagf(a) 311 #define PetscAbsComplex(a) cabsf(a) 312 #define PetscArgComplex(a) cargf(a) 313 #define PetscConjComplex(a) conjf(a) 314 #define PetscSqrtComplex(a) csqrtf(a) 315 #define PetscPowComplex(a, b) cpowf(a, b) 316 #define PetscExpComplex(a) cexpf(a) 317 #define PetscLogComplex(a) clogf(a) 318 #define PetscSinComplex(a) csinf(a) 319 #define PetscCosComplex(a) ccosf(a) 320 #define PetscTanComplex(a) ctanf(a) 321 #define PetscAsinComplex(a) casinf(a) 322 #define PetscAcosComplex(a) cacosf(a) 323 #define PetscAtanComplex(a) catanf(a) 324 #define PetscSinhComplex(a) csinhf(a) 325 #define PetscCoshComplex(a) ccoshf(a) 326 #define PetscTanhComplex(a) ctanhf(a) 327 #define PetscAsinhComplex(a) casinhf(a) 328 #define PetscAcoshComplex(a) cacoshf(a) 329 #define PetscAtanhComplex(a) catanhf(a) 330 331 #elif defined(PETSC_USE_REAL_DOUBLE) 332 #define PetscRealPartComplex(a) creal(a) 333 #define PetscImaginaryPartComplex(a) cimag(a) 334 #define PetscAbsComplex(a) cabs(a) 335 #define PetscArgComplex(a) carg(a) 336 #define PetscConjComplex(a) conj(a) 337 #define PetscSqrtComplex(a) csqrt(a) 338 #define PetscPowComplex(a, b) cpow(a, b) 339 #define PetscExpComplex(a) cexp(a) 340 #define PetscLogComplex(a) clog(a) 341 #define PetscSinComplex(a) csin(a) 342 #define PetscCosComplex(a) ccos(a) 343 #define PetscTanComplex(a) ctan(a) 344 #define PetscAsinComplex(a) casin(a) 345 #define PetscAcosComplex(a) cacos(a) 346 #define PetscAtanComplex(a) catan(a) 347 #define PetscSinhComplex(a) csinh(a) 348 #define PetscCoshComplex(a) ccosh(a) 349 #define PetscTanhComplex(a) ctanh(a) 350 #define PetscAsinhComplex(a) casinh(a) 351 #define PetscAcoshComplex(a) cacosh(a) 352 #define PetscAtanhComplex(a) catanh(a) 353 354 #elif defined(PETSC_USE_REAL___FLOAT128) 355 #define PetscRealPartComplex(a) crealq(a) 356 #define PetscImaginaryPartComplex(a) cimagq(a) 357 #define PetscAbsComplex(a) cabsq(a) 358 #define PetscArgComplex(a) cargq(a) 359 #define PetscConjComplex(a) conjq(a) 360 #define PetscSqrtComplex(a) csqrtq(a) 361 #define PetscPowComplex(a, b) cpowq(a, b) 362 #define PetscExpComplex(a) cexpq(a) 363 #define PetscLogComplex(a) clogq(a) 364 #define PetscSinComplex(a) csinq(a) 365 #define PetscCosComplex(a) ccosq(a) 366 #define PetscTanComplex(a) ctanq(a) 367 #define PetscAsinComplex(a) casinq(a) 368 #define PetscAcosComplex(a) cacosq(a) 369 #define PetscAtanComplex(a) catanq(a) 370 #define PetscSinhComplex(a) csinhq(a) 371 #define PetscCoshComplex(a) ccoshq(a) 372 #define PetscTanhComplex(a) ctanhq(a) 373 #define PetscAsinhComplex(a) casinhq(a) 374 #define PetscAcoshComplex(a) cacoshq(a) 375 #define PetscAtanhComplex(a) catanhq(a) 376 377 #endif /* PETSC_USE_REAL_* */ 378 #endif /* (__cplusplus) */ 379 380 /* 381 PETSC_i is the imaginary number, i 382 */ 383 PETSC_EXTERN PetscComplex PETSC_i; 384 385 /* 386 Try to do the right thing for complex number construction: see 387 http://www.open-std.org/jtc1/sc22/wg14/www/docs/n1464.htm 388 for details 389 */ 390 static inline PetscComplex PetscCMPLX(PetscReal x, PetscReal y) { 391 #if defined(__cplusplus) && !defined(PETSC_USE_REAL___FLOAT128) 392 return PetscComplex(x, y); 393 #elif defined(_Imaginary_I) 394 return x + y * _Imaginary_I; 395 #else 396 { /* In both C99 and C11 (ISO/IEC 9899, Section 6.2.5), 397 398 "For each floating type there is a corresponding real type, which is always a real floating 399 type. For real floating types, it is the same type. For complex types, it is the type given 400 by deleting the keyword _Complex from the type name." 401 402 So type punning should be portable. */ 403 union 404 { 405 PetscComplex z; 406 PetscReal f[2]; 407 } uz; 408 409 uz.f[0] = x; 410 uz.f[1] = y; 411 return uz.z; 412 } 413 #endif 414 } 415 416 #define MPIU_C_COMPLEX MPI_C_COMPLEX PETSC_DEPRECATED_MACRO("GCC warning \"MPIU_C_COMPLEX macro is deprecated use MPI_C_COMPLEX (since version 3.15)\"") 417 #define MPIU_C_DOUBLE_COMPLEX MPI_C_DOUBLE_COMPLEX PETSC_DEPRECATED_MACRO("GCC warning \"MPIU_C_DOUBLE_COMPLEX macro is deprecated use MPI_C_DOUBLE_COMPLEX (since version 3.15)\"") 418 419 #if defined(PETSC_HAVE_REAL___FLOAT128) 420 // if complex is not used, then quadmath.h won't be included by petscsystypes.h 421 #if defined(PETSC_USE_COMPLEX) 422 #define MPIU___COMPLEX128_ATTR_TAG PETSC_ATTRIBUTE_MPI_TYPE_TAG(__complex128) 423 #else 424 #define MPIU___COMPLEX128_ATTR_TAG 425 #endif 426 427 PETSC_EXTERN MPI_Datatype MPIU___COMPLEX128 MPIU___COMPLEX128_ATTR_TAG; 428 429 #undef MPIU___COMPLEX128_ATTR_TAG 430 #endif /* PETSC_HAVE_REAL___FLOAT128 */ 431 432 /*MC 433 MPIU_COMPLEX - Portable MPI datatype corresponding to `PetscComplex` independent of the precision of `PetscComplex` 434 435 Notes: 436 In MPI calls that require an MPI datatype that matches a `PetscComplex` or array of `PetscComplex` values, pass this value. 437 438 Level: beginner 439 440 .seealso: `PetscReal`, `PetscScalar`, `PetscComplex`, `PetscInt`, `MPIU_REAL`, `MPIU_SCALAR`, `MPIU_COMPLEX`, `MPIU_INT`, `PETSC_i` 441 M*/ 442 #if defined(PETSC_USE_REAL_SINGLE) 443 #define MPIU_COMPLEX MPI_C_COMPLEX 444 #elif defined(PETSC_USE_REAL_DOUBLE) 445 #define MPIU_COMPLEX MPI_C_DOUBLE_COMPLEX 446 #elif defined(PETSC_USE_REAL___FLOAT128) 447 #define MPIU_COMPLEX MPIU___COMPLEX128 448 #elif defined(PETSC_USE_REAL___FP16) 449 #define MPIU_COMPLEX MPI_C_COMPLEX 450 #endif /* PETSC_USE_REAL_* */ 451 452 #endif /* PETSC_HAVE_COMPLEX */ 453 454 /* 455 Scalar number definitions 456 */ 457 #if defined(PETSC_USE_COMPLEX) && defined(PETSC_HAVE_COMPLEX) 458 /*MC 459 MPIU_SCALAR - Portable MPI datatype corresponding to `PetscScalar` independent of the precision of `PetscScalar` 460 461 Notes: 462 In MPI calls that require an MPI datatype that matches a `PetscScalar` or array of `PetscScalar` values, pass this value. 463 464 Level: beginner 465 466 .seealso: `PetscReal`, `PetscScalar`, `PetscComplex`, `PetscInt`, `MPIU_REAL`, `MPIU_COMPLEX`, `MPIU_INT` 467 M*/ 468 #define MPIU_SCALAR MPIU_COMPLEX 469 470 /*MC 471 PetscRealPart - Returns the real part of a `PetscScalar` 472 473 Synopsis: 474 #include <petscmath.h> 475 PetscReal PetscRealPart(PetscScalar v) 476 477 Not Collective 478 479 Input Parameter: 480 . v - value to find the real part of 481 482 Level: beginner 483 484 .seealso: `PetscScalar`, `PetscImaginaryPart()`, `PetscMax()`, `PetscClipInterval()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscSqr()` 485 486 M*/ 487 #define PetscRealPart(a) PetscRealPartComplex(a) 488 489 /*MC 490 PetscImaginaryPart - Returns the imaginary part of a `PetscScalar` 491 492 Synopsis: 493 #include <petscmath.h> 494 PetscReal PetscImaginaryPart(PetscScalar v) 495 496 Not Collective 497 498 Input Parameter: 499 . v - value to find the imaginary part of 500 501 Level: beginner 502 503 Notes: 504 If PETSc was configured for real numbers then this always returns the value 0 505 506 .seealso: `PetscScalar`, `PetscRealPart()`, `PetscMax()`, `PetscClipInterval()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscSqr()` 507 508 M*/ 509 #define PetscImaginaryPart(a) PetscImaginaryPartComplex(a) 510 511 #define PetscAbsScalar(a) PetscAbsComplex(a) 512 #define PetscArgScalar(a) PetscArgComplex(a) 513 #define PetscConj(a) PetscConjComplex(a) 514 #define PetscSqrtScalar(a) PetscSqrtComplex(a) 515 #define PetscPowScalar(a, b) PetscPowComplex(a, b) 516 #define PetscExpScalar(a) PetscExpComplex(a) 517 #define PetscLogScalar(a) PetscLogComplex(a) 518 #define PetscSinScalar(a) PetscSinComplex(a) 519 #define PetscCosScalar(a) PetscCosComplex(a) 520 #define PetscTanScalar(a) PetscTanComplex(a) 521 #define PetscAsinScalar(a) PetscAsinComplex(a) 522 #define PetscAcosScalar(a) PetscAcosComplex(a) 523 #define PetscAtanScalar(a) PetscAtanComplex(a) 524 #define PetscSinhScalar(a) PetscSinhComplex(a) 525 #define PetscCoshScalar(a) PetscCoshComplex(a) 526 #define PetscTanhScalar(a) PetscTanhComplex(a) 527 #define PetscAsinhScalar(a) PetscAsinhComplex(a) 528 #define PetscAcoshScalar(a) PetscAcoshComplex(a) 529 #define PetscAtanhScalar(a) PetscAtanhComplex(a) 530 531 #else /* PETSC_USE_COMPLEX */ 532 #define MPIU_SCALAR MPIU_REAL 533 #define PetscRealPart(a) (a) 534 #define PetscImaginaryPart(a) ((PetscReal)0) 535 #define PetscAbsScalar(a) PetscAbsReal(a) 536 #define PetscArgScalar(a) (((a) < (PetscReal)0) ? PETSC_PI : (PetscReal)0) 537 #define PetscConj(a) (a) 538 #define PetscSqrtScalar(a) PetscSqrtReal(a) 539 #define PetscPowScalar(a, b) PetscPowReal(a, b) 540 #define PetscExpScalar(a) PetscExpReal(a) 541 #define PetscLogScalar(a) PetscLogReal(a) 542 #define PetscSinScalar(a) PetscSinReal(a) 543 #define PetscCosScalar(a) PetscCosReal(a) 544 #define PetscTanScalar(a) PetscTanReal(a) 545 #define PetscAsinScalar(a) PetscAsinReal(a) 546 #define PetscAcosScalar(a) PetscAcosReal(a) 547 #define PetscAtanScalar(a) PetscAtanReal(a) 548 #define PetscSinhScalar(a) PetscSinhReal(a) 549 #define PetscCoshScalar(a) PetscCoshReal(a) 550 #define PetscTanhScalar(a) PetscTanhReal(a) 551 #define PetscAsinhScalar(a) PetscAsinhReal(a) 552 #define PetscAcoshScalar(a) PetscAcoshReal(a) 553 #define PetscAtanhScalar(a) PetscAtanhReal(a) 554 555 #endif /* PETSC_USE_COMPLEX */ 556 557 /* 558 Certain objects may be created using either single or double precision. 559 This is currently not used. 560 */ 561 typedef enum { 562 PETSC_SCALAR_DOUBLE, 563 PETSC_SCALAR_SINGLE, 564 PETSC_SCALAR_LONG_DOUBLE, 565 PETSC_SCALAR_HALF 566 } PetscScalarPrecision; 567 568 /* --------------------------------------------------------------------------*/ 569 570 /*MC 571 PetscAbs - Returns the absolute value of a number 572 573 Synopsis: 574 #include <petscmath.h> 575 type PetscAbs(type v) 576 577 Not Collective 578 579 Input Parameter: 580 . v - the number 581 582 Note: 583 The type can be integer or real floating point value, but cannot be complex 584 585 Level: beginner 586 587 .seealso: `PetscAbsInt()`, `PetscAbsReal()`, `PetscAbsScalar()` 588 589 M*/ 590 #define PetscAbs(a) (((a) >= 0) ? (a) : (-(a))) 591 592 /*MC 593 PetscSign - Returns the sign of a number as an integer 594 595 Synopsis: 596 #include <petscmath.h> 597 int PetscSign(type v) 598 599 Not Collective 600 601 Input Parameter: 602 . v - the number 603 604 Note: 605 The type can be integer or real floating point value 606 607 Level: beginner 608 609 M*/ 610 #define PetscSign(a) (((a) >= 0) ? ((a) == 0 ? 0 : 1) : -1) 611 612 /*MC 613 PetscMin - Returns minimum of two numbers 614 615 Synopsis: 616 #include <petscmath.h> 617 type PetscMin(type v1,type v2) 618 619 Not Collective 620 621 Input Parameters: 622 + v1 - first value to find minimum of 623 - v2 - second value to find minimum of 624 625 Note: 626 The type can be integer or floating point value 627 628 Level: beginner 629 630 .seealso: `PetscMax()`, `PetscClipInterval()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscSqr()` 631 632 M*/ 633 #define PetscMin(a, b) (((a) < (b)) ? (a) : (b)) 634 635 /*MC 636 PetscMax - Returns maxium of two numbers 637 638 Synopsis: 639 #include <petscmath.h> 640 type max PetscMax(type v1,type v2) 641 642 Not Collective 643 644 Input Parameters: 645 + v1 - first value to find maximum of 646 - v2 - second value to find maximum of 647 648 Note: 649 The type can be integer or floating point value 650 651 Level: beginner 652 653 .seealso: `PetscMin()`, `PetscClipInterval()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscSqr()` 654 655 M*/ 656 #define PetscMax(a, b) (((a) < (b)) ? (b) : (a)) 657 658 /*MC 659 PetscClipInterval - Returns a number clipped to be within an interval 660 661 Synopsis: 662 #include <petscmath.h> 663 type clip PetscClipInterval(type x,type a,type b) 664 665 Not Collective 666 667 Input Parameters: 668 + x - value to use if within interval [a,b] 669 . a - lower end of interval 670 - b - upper end of interval 671 672 Note: 673 The type can be integer or floating point value 674 675 Level: beginner 676 677 .seealso: `PetscMin()`, `PetscMax()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscSqr()` 678 679 M*/ 680 #define PetscClipInterval(x, a, b) (PetscMax((a), PetscMin((x), (b)))) 681 682 /*MC 683 PetscAbsInt - Returns the absolute value of an integer 684 685 Synopsis: 686 #include <petscmath.h> 687 int abs PetscAbsInt(int v1) 688 689 Input Parameter: 690 . v1 - the integer 691 692 Level: beginner 693 694 .seealso: `PetscMax()`, `PetscMin()`, `PetscAbsReal()`, `PetscSqr()` 695 696 M*/ 697 #define PetscAbsInt(a) (((a) < 0) ? (-(a)) : (a)) 698 699 /*MC 700 PetscAbsReal - Returns the absolute value of an real number 701 702 Synopsis: 703 #include <petscmath.h> 704 Real abs PetscAbsReal(PetscReal v1) 705 706 Input Parameter: 707 . v1 - the double 708 709 Level: beginner 710 711 .seealso: `PetscMax()`, `PetscMin()`, `PetscAbsInt()`, `PetscSqr()` 712 713 M*/ 714 #if defined(PETSC_USE_REAL_SINGLE) 715 #define PetscAbsReal(a) fabsf(a) 716 #elif defined(PETSC_USE_REAL_DOUBLE) 717 #define PetscAbsReal(a) fabs(a) 718 #elif defined(PETSC_USE_REAL___FLOAT128) 719 #define PetscAbsReal(a) fabsq(a) 720 #elif defined(PETSC_USE_REAL___FP16) 721 #define PetscAbsReal(a) fabsf(a) 722 #endif 723 724 /*MC 725 PetscSqr - Returns the square of a number 726 727 Synopsis: 728 #include <petscmath.h> 729 type sqr PetscSqr(type v1) 730 731 Not Collective 732 733 Input Parameter: 734 . v1 - the value 735 736 Note: 737 The type can be integer or floating point value 738 739 Level: beginner 740 741 .seealso: `PetscMax()`, `PetscMin()`, `PetscAbsInt()`, `PetscAbsReal()` 742 743 M*/ 744 #define PetscSqr(a) ((a) * (a)) 745 746 /* ----------------------------------------------------------------------------*/ 747 748 #if defined(PETSC_USE_REAL_SINGLE) 749 #define PetscRealConstant(constant) constant##F 750 #elif defined(PETSC_USE_REAL_DOUBLE) 751 #define PetscRealConstant(constant) constant 752 #elif defined(PETSC_USE_REAL___FLOAT128) 753 #define PetscRealConstant(constant) constant##Q 754 #elif defined(PETSC_USE_REAL___FP16) 755 #define PetscRealConstant(constant) constant##F 756 #endif 757 758 /* 759 Basic constants 760 */ 761 #define PETSC_PI PetscRealConstant(3.1415926535897932384626433832795029) 762 #define PETSC_PHI PetscRealConstant(1.6180339887498948482045868343656381) 763 #define PETSC_SQRT2 PetscRealConstant(1.4142135623730950488016887242096981) 764 765 #if !defined(PETSC_USE_64BIT_INDICES) 766 #define PETSC_MAX_INT 2147483647 767 #define PETSC_MIN_INT (-PETSC_MAX_INT - 1) 768 #else 769 #define PETSC_MAX_INT 9223372036854775807L 770 #define PETSC_MIN_INT (-PETSC_MAX_INT - 1) 771 #endif 772 #define PETSC_MAX_UINT16 65535 773 774 #if defined(PETSC_USE_REAL_SINGLE) 775 #define PETSC_MAX_REAL 3.40282346638528860e+38F 776 #define PETSC_MIN_REAL (-PETSC_MAX_REAL) 777 #define PETSC_MACHINE_EPSILON 1.19209290e-07F 778 #define PETSC_SQRT_MACHINE_EPSILON 3.45266983e-04F 779 #define PETSC_SMALL 1.e-5F 780 #elif defined(PETSC_USE_REAL_DOUBLE) 781 #define PETSC_MAX_REAL 1.7976931348623157e+308 782 #define PETSC_MIN_REAL (-PETSC_MAX_REAL) 783 #define PETSC_MACHINE_EPSILON 2.2204460492503131e-16 784 #define PETSC_SQRT_MACHINE_EPSILON 1.490116119384766e-08 785 #define PETSC_SMALL 1.e-10 786 #elif defined(PETSC_USE_REAL___FLOAT128) 787 #define PETSC_MAX_REAL FLT128_MAX 788 #define PETSC_MIN_REAL (-FLT128_MAX) 789 #define PETSC_MACHINE_EPSILON FLT128_EPSILON 790 #define PETSC_SQRT_MACHINE_EPSILON 1.38777878078144567552953958511352539e-17Q 791 #define PETSC_SMALL 1.e-20Q 792 #elif defined(PETSC_USE_REAL___FP16) 793 #define PETSC_MAX_REAL 65504.0F 794 #define PETSC_MIN_REAL (-PETSC_MAX_REAL) 795 #define PETSC_MACHINE_EPSILON .0009765625F 796 #define PETSC_SQRT_MACHINE_EPSILON .03125F 797 #define PETSC_SMALL 5.e-3F 798 #endif 799 800 #define PETSC_INFINITY (PETSC_MAX_REAL / 4) 801 #define PETSC_NINFINITY (-PETSC_INFINITY) 802 803 PETSC_EXTERN PetscBool PetscIsInfReal(PetscReal); 804 PETSC_EXTERN PetscBool PetscIsNanReal(PetscReal); 805 PETSC_EXTERN PetscBool PetscIsNormalReal(PetscReal); 806 static inline PetscBool PetscIsInfOrNanReal(PetscReal v) { 807 return PetscIsInfReal(v) || PetscIsNanReal(v) ? PETSC_TRUE : PETSC_FALSE; 808 } 809 static inline PetscBool PetscIsInfScalar(PetscScalar v) { 810 return PetscIsInfReal(PetscAbsScalar(v)); 811 } 812 static inline PetscBool PetscIsNanScalar(PetscScalar v) { 813 return PetscIsNanReal(PetscAbsScalar(v)); 814 } 815 static inline PetscBool PetscIsInfOrNanScalar(PetscScalar v) { 816 return PetscIsInfOrNanReal(PetscAbsScalar(v)); 817 } 818 static inline PetscBool PetscIsNormalScalar(PetscScalar v) { 819 return PetscIsNormalReal(PetscAbsScalar(v)); 820 } 821 822 PETSC_EXTERN PetscBool PetscIsCloseAtTol(PetscReal, PetscReal, PetscReal, PetscReal); 823 PETSC_EXTERN PetscBool PetscEqualReal(PetscReal, PetscReal); 824 PETSC_EXTERN PetscBool PetscEqualScalar(PetscScalar, PetscScalar); 825 826 /* 827 These macros are currently hardwired to match the regular data types, so there is no support for a different 828 MatScalar from PetscScalar. We left the MatScalar in the source just in case we use it again. 829 */ 830 #define MPIU_MATSCALAR MPIU_SCALAR 831 typedef PetscScalar MatScalar; 832 typedef PetscReal MatReal; 833 834 struct petsc_mpiu_2scalar { 835 PetscScalar a, b; 836 }; 837 PETSC_EXTERN MPI_Datatype MPIU_2SCALAR PETSC_ATTRIBUTE_MPI_TYPE_TAG_LAYOUT_COMPATIBLE(struct petsc_mpiu_2scalar); 838 839 /* MPI Datatypes for composite reductions */ 840 struct petsc_mpiu_real_int { 841 PetscReal v; 842 PetscInt i; 843 }; 844 845 struct petsc_mpiu_scalar_int { 846 PetscScalar v; 847 PetscInt i; 848 }; 849 850 PETSC_EXTERN MPI_Datatype MPIU_REAL_INT PETSC_ATTRIBUTE_MPI_TYPE_TAG_LAYOUT_COMPATIBLE(struct petsc_mpiu_real_int); 851 PETSC_EXTERN MPI_Datatype MPIU_SCALAR_INT PETSC_ATTRIBUTE_MPI_TYPE_TAG_LAYOUT_COMPATIBLE(struct petsc_mpiu_scalar_int); 852 853 #if defined(PETSC_USE_64BIT_INDICES) 854 struct /* __attribute__((packed, aligned(alignof(PetscInt *)))) */ petsc_mpiu_2int { 855 PetscInt a; 856 PetscInt b; 857 }; 858 /* 859 static_assert(sizeof(struct petsc_mpiu_2int) == 2 * sizeof(PetscInt), ""); 860 static_assert(alignof(struct petsc_mpiu_2int) == alignof(PetscInt *), ""); 861 static_assert(alignof(struct petsc_mpiu_2int) == alignof(PetscInt[2]), ""); 862 863 clang generates warnings that petsc_mpiu_2int is not layout compatible with PetscInt[2] or 864 PetscInt *, even though (with everything else uncommented) both of the static_asserts above 865 pass! So we just comment it out... 866 */ 867 PETSC_EXTERN MPI_Datatype MPIU_2INT /* PETSC_ATTRIBUTE_MPI_TYPE_TAG_LAYOUT_COMPATIBLE(struct petsc_mpiu_2int) */; 868 #else 869 #define MPIU_2INT MPI_2INT 870 #endif 871 PETSC_EXTERN MPI_Datatype MPI_4INT; 872 PETSC_EXTERN MPI_Datatype MPIU_4INT; 873 874 static inline PetscInt PetscPowInt(PetscInt base, PetscInt power) { 875 PetscInt result = 1; 876 while (power) { 877 if (power & 1) result *= base; 878 power >>= 1; 879 base *= base; 880 } 881 return result; 882 } 883 884 static inline PetscInt64 PetscPowInt64(PetscInt base, PetscInt power) { 885 PetscInt64 result = 1; 886 while (power) { 887 if (power & 1) result *= base; 888 power >>= 1; 889 base *= base; 890 } 891 return result; 892 } 893 894 static inline PetscReal PetscPowRealInt(PetscReal base, PetscInt power) { 895 PetscReal result = 1; 896 if (power < 0) { 897 power = -power; 898 base = ((PetscReal)1) / base; 899 } 900 while (power) { 901 if (power & 1) result *= base; 902 power >>= 1; 903 base *= base; 904 } 905 return result; 906 } 907 908 static inline PetscScalar PetscPowScalarInt(PetscScalar base, PetscInt power) { 909 PetscScalar result = (PetscReal)1; 910 if (power < 0) { 911 power = -power; 912 base = ((PetscReal)1) / base; 913 } 914 while (power) { 915 if (power & 1) result *= base; 916 power >>= 1; 917 base *= base; 918 } 919 return result; 920 } 921 922 static inline PetscScalar PetscPowScalarReal(PetscScalar base, PetscReal power) { 923 PetscScalar cpower = power; 924 return PetscPowScalar(base, cpower); 925 } 926 927 /*MC 928 PetscApproximateLTE - Performs a less than or equal to on a given constant with a fudge for floating point numbers 929 930 Synopsis: 931 #include <petscmath.h> 932 bool PetscApproximateLTE(PetscReal x,constant float) 933 934 Not Collective 935 936 Input Parameters: 937 + x - the variable 938 - b - the constant float it is checking if x is less than or equal to 939 940 Notes: 941 The fudge factor is the value `PETSC_SMALL` 942 943 The constant numerical value is automatically set to the appropriate precision of PETSc so can just be provided as, for example, 3.2 944 945 This is used in several examples for setting initial conditions based on coordinate values that are computed with i*h that produces inexact 946 floating point results. 947 948 Level: advanced 949 950 .seealso: `PetscMax()`, `PetscMin()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscApproximateGTE()` 951 952 M*/ 953 #define PetscApproximateLTE(x, b) ((x) <= (PetscRealConstant(b) + PETSC_SMALL)) 954 955 /*MC 956 PetscApproximateGTE - Performs a greater than or equal to on a given constant with a fudge for floating point numbers 957 958 Synopsis: 959 #include <petscmath.h> 960 bool PetscApproximateGTE(PetscReal x,constant float) 961 962 Not Collective 963 964 Input Parameters: 965 + x - the variable 966 - b - the constant float it is checking if x is greater than or equal to 967 968 Notes: 969 The fudge factor is the value `PETSC_SMALL` 970 971 The constant numerical value is automatically set to the appropriate precision of PETSc so can just be provided as, for example, 3.2 972 973 This is used in several examples for setting initial conditions based on coordinate values that are computed with i*h that produces inexact 974 floating point results. 975 976 Level: advanced 977 978 .seealso: `PetscMax()`, `PetscMin()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscApproximateLTE()` 979 980 M*/ 981 #define PetscApproximateGTE(x, b) ((x) >= (PetscRealConstant(b) - PETSC_SMALL)) 982 983 /*MC 984 PetscCeilInt - Returns the ceiling of the quotation of two positive integers 985 986 Synopsis: 987 #include <petscmath.h> 988 PetscInt PetscCeilInt(PetscInt x,PetscInt y) 989 990 Not Collective 991 992 Input Parameters: 993 + x - the numerator 994 - y - the denominator 995 996 Level: advanced 997 998 .seealso: `PetscMax()`, `PetscMin()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscApproximateLTE()` 999 1000 M*/ 1001 #define PetscCeilInt(x, y) ((((PetscInt)(x)) / ((PetscInt)(y))) + ((((PetscInt)(x)) % ((PetscInt)(y))) ? 1 : 0)) 1002 1003 #define PetscCeilInt64(x, y) ((((PetscInt64)(x)) / ((PetscInt64)(y))) + ((((PetscInt64)(x)) % ((PetscInt64)(y))) ? 1 : 0)) 1004 1005 PETSC_EXTERN PetscErrorCode PetscLinearRegression(PetscInt, const PetscReal[], const PetscReal[], PetscReal *, PetscReal *); 1006 #endif 1007