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