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