1 /* 2 3 PETSc mathematics include file. Defines certain basic mathematical 4 constants and functions for working with single, double, and quad precision 5 floating point numbers as well as complex single and double. 6 7 This file is included by petscsys.h and should not be used directly. 8 9 */ 10 11 #if !defined(__PETSCMATH_H) 12 #define __PETSCMATH_H 13 #include <math.h> 14 15 /* 16 17 Defines operations that are different for complex and real numbers; 18 note that one cannot mix the use of complex and real in the same 19 PETSc program. All PETSc objects in one program are built around the object 20 PetscScalar which is either always a real or a complex. 21 22 */ 23 24 #if defined(PETSC_USE_REAL_SINGLE) 25 #define MPIU_REAL MPI_FLOAT 26 typedef float PetscReal; 27 #define PetscRoundReal(a) roundf(a) 28 #define PetscSqrtReal(a) sqrtf(a) 29 #define PetscExpReal(a) expf(a) 30 #define PetscLogReal(a) logf(a) 31 #define PetscLog10Real(a) log10f(a) 32 #ifdef PETSC_HAVE_LOG2 33 #define PetscLog2Real(a) log2f(a) 34 #endif 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 PetscAtan2Real(a,b) atan2f(a,b) 42 #define PetscSinhReal(a) sinhf(a) 43 #define PetscCoshReal(a) coshf(a) 44 #define PetscTanhReal(a) tanhf(a) 45 #define PetscPowReal(a,b) powf(a,b) 46 #define PetscCeilReal(a) ceilf(a) 47 #define PetscFloorReal(a) floorf(a) 48 #define PetscFmodReal(a,b) fmodf(a,b) 49 #define PetscTGamma(a) tgammaf(a) 50 #elif defined(PETSC_USE_REAL_DOUBLE) 51 #define MPIU_REAL MPI_DOUBLE 52 typedef double PetscReal; 53 #define PetscRoundReal(a) round(a) 54 #define PetscSqrtReal(a) sqrt(a) 55 #define PetscExpReal(a) exp(a) 56 #define PetscLogReal(a) log(a) 57 #define PetscLog10Real(a) log10(a) 58 #ifdef PETSC_HAVE_LOG2 59 #define PetscLog2Real(a) log2(a) 60 #endif 61 #define PetscSinReal(a) sin(a) 62 #define PetscCosReal(a) cos(a) 63 #define PetscTanReal(a) tan(a) 64 #define PetscAsinReal(a) asin(a) 65 #define PetscAcosReal(a) acos(a) 66 #define PetscAtanReal(a) atan(a) 67 #define PetscAtan2Real(a,b) atan2(a,b) 68 #define PetscSinhReal(a) sinh(a) 69 #define PetscCoshReal(a) cosh(a) 70 #define PetscTanhReal(a) tanh(a) 71 #define PetscPowReal(a,b) pow(a,b) 72 #define PetscCeilReal(a) ceil(a) 73 #define PetscFloorReal(a) floor(a) 74 #define PetscFmodReal(a,b) fmod(a,b) 75 #define PetscTGamma(a) tgamma(a) 76 #elif defined(PETSC_USE_REAL___FLOAT128) 77 #if defined(__cplusplus) 78 extern "C" { 79 #endif 80 #include <quadmath.h> 81 #if defined(__cplusplus) 82 } 83 #endif 84 PETSC_EXTERN MPI_Datatype MPIU___FLOAT128 PetscAttrMPITypeTag(__float128); 85 #define MPIU_REAL MPIU___FLOAT128 86 typedef __float128 PetscReal; 87 #define PetscRoundReal(a) roundq(a) 88 #define PetscSqrtReal(a) sqrtq(a) 89 #define PetscExpReal(a) expq(a) 90 #define PetscLogReal(a) logq(a) 91 #define PetscLog10Real(a) log10q(a) 92 #ifdef PETSC_HAVE_LOG2 93 #define PetscLog2Real(a) log2q(a) 94 #endif 95 #define PetscSinReal(a) sinq(a) 96 #define PetscCosReal(a) cosq(a) 97 #define PetscTanReal(a) tanq(a) 98 #define PetscAsinReal(a) asinq(a) 99 #define PetscAcosReal(a) acosq(a) 100 #define PetscAtanReal(a) atanq(a) 101 #define PetscAtan2Real(a,b) atan2q(a,b) 102 #define PetscSinhReal(a) sinhq(a) 103 #define PetscCoshReal(a) coshq(a) 104 #define PetscTanhReal(a) tanhq(a) 105 #define PetscPowReal(a,b) powq(a,b) 106 #define PetscCeilReal(a) ceilq(a) 107 #define PetscFloorReal(a) floorq(a) 108 #define PetscFmodReal(a,b) fmodq(a,b) 109 #define PetscTGamma(a) tgammaq(a) 110 #elif defined(PETSC_USE_REAL___FP16) 111 PETSC_EXTERN MPI_Datatype MPIU___FP16 PetscAttrMPITypeTag(__fp16); 112 #define MPIU_REAL MPIU___FP16 113 typedef __fp16 PetscReal; 114 #define PetscRoundReal(a) roundf(a) 115 #define PetscSqrtReal(a) sqrtf(a) 116 #define PetscExpReal(a) expf(a) 117 #define PetscLogReal(a) logf(a) 118 #define PetscLog10Real(a) log10f(a) 119 #ifdef PETSC_HAVE_LOG2 120 #define PetscLog2Real(a) log2f(a) 121 #endif 122 #define PetscSinReal(a) sinf(a) 123 #define PetscCosReal(a) cosf(a) 124 #define PetscTanReal(a) tanf(a) 125 #define PetscAsinReal(a) asinf(a) 126 #define PetscAcosReal(a) acosf(a) 127 #define PetscAtanReal(a) atanf(a) 128 #define PetscAtan2Real(a,b) atan2f(a,b) 129 #define PetscSinhReal(a) sinhf(a) 130 #define PetscCoshReal(a) coshf(a) 131 #define PetscTanhReal(a) tanhf(a) 132 #define PetscPowReal(a,b) powf(a,b) 133 #define PetscCeilReal(a) ceilf(a) 134 #define PetscFloorReal(a) floorf(a) 135 #define PetscFmodReal(a,b) fmodf(a,b) 136 #define PetscTGamma(a) tgammaf(a) 137 #endif /* PETSC_USE_REAL_* */ 138 139 /* 140 Complex number definitions 141 */ 142 #if defined(__cplusplus) && defined(PETSC_HAVE_CXX_COMPLEX) && !defined(PETSC_USE_REAL___FLOAT128) 143 #if !defined(PETSC_SKIP_COMPLEX) 144 #define PETSC_HAVE_COMPLEX 1 145 /* C++ support of complex number */ 146 #if defined(PETSC_HAVE_VECCUDA) && __CUDACC_VER_MAJOR__ > 6 147 /* complex headers in thrust only available in CUDA 7.0 and above */ 148 #define complexlib thrust 149 #include <thrust/complex.h> 150 #else 151 #define complexlib std 152 #include <complex> 153 #endif 154 155 #define PetscRealPartComplex(a) (a).real() 156 #define PetscImaginaryPartComplex(a) (a).imag() 157 #define PetscAbsComplex(a) complexlib::abs(a) 158 #define PetscConjComplex(a) complexlib::conj(a) 159 #define PetscSqrtComplex(a) complexlib::sqrt(a) 160 #define PetscPowComplex(a,b) complexlib::pow(a,b) 161 #define PetscExpComplex(a) complexlib::exp(a) 162 #define PetscLogComplex(a) complexlib::log(a) 163 #define PetscSinComplex(a) complexlib::sin(a) 164 #define PetscCosComplex(a) complexlib::cos(a) 165 #define PetscAsinComplex(a) complexlib::asin(a) 166 #define PetscAcosComplex(a) complexlib::acos(a) 167 #if defined(PETSC_HAVE_TANCOMPLEX) 168 #define PetscTanComplex(a) complexlib::tan(a) 169 #else 170 #define PetscTanComplex(a) PetscSinComplex(a)/PetscCosComplex(a) 171 #endif 172 #define PetscSinhComplex(a) complexlib::sinh(a) 173 #define PetscCoshComplex(a) complexlib::cosh(a) 174 #if defined(PETSC_HAVE_TANHCOMPLEX) 175 #define PetscTanhComplex(a) complexlib::tanh(a) 176 #else 177 #define PetscTanhComplex(a) PetscSinhComplex(a)/PetscCoshComplex(a) 178 #endif 179 180 #if defined(PETSC_USE_REAL_SINGLE) 181 typedef complexlib::complex<float> PetscComplex; 182 #if defined(PETSC_USE_CXX_COMPLEX_FLOAT_WORKAROUND) 183 static inline PetscComplex operator+(const PetscComplex& lhs, const double& rhs) { return lhs + float(rhs); } 184 static inline PetscComplex operator+(const double& lhs, const PetscComplex& rhs) { return float(lhs) + rhs; } 185 static inline PetscComplex operator-(const PetscComplex& lhs, const double& rhs) { return lhs - float(rhs); } 186 static inline PetscComplex operator-(const double& lhs, const PetscComplex& rhs) { return float(lhs) - rhs; } 187 static inline PetscComplex operator*(const PetscComplex& lhs, const double& rhs) { return lhs * float(rhs); } 188 static inline PetscComplex operator*(const double& lhs, const PetscComplex& rhs) { return float(lhs) * rhs; } 189 static inline PetscComplex operator/(const PetscComplex& lhs, const double& rhs) { return lhs / float(rhs); } 190 static inline PetscComplex operator/(const double& lhs, const PetscComplex& rhs) { return float(lhs) / rhs; } 191 static inline bool operator==(const PetscComplex& lhs, const double& rhs) { return lhs.imag() == float(0) && lhs.real() == float(rhs); } 192 static inline bool operator==(const double& lhs, const PetscComplex& rhs) { return rhs.imag() == float(0) && rhs.real() == float(lhs); } 193 static inline bool operator!=(const PetscComplex& lhs, const double& rhs) { return lhs.imag() != float(0) || lhs.real() != float(rhs); } 194 static inline bool operator!=(const double& lhs, const PetscComplex& rhs) { return rhs.imag() != float(0) || rhs.real() != float(lhs); } 195 #endif /* PETSC_USE_CXX_COMPLEX_FLOAT_WORKAROUND */ 196 #elif defined(PETSC_USE_REAL_DOUBLE) 197 typedef complexlib::complex<double> PetscComplex; 198 #if defined(PETSC_USE_CXX_COMPLEX_FLOAT_WORKAROUND) 199 static inline PetscComplex operator+(const PetscComplex& lhs, const PetscInt& rhs) { return lhs + double(rhs); } 200 static inline PetscComplex operator+(const PetscInt& lhs, const PetscComplex& rhs) { return double(lhs) + rhs; } 201 static inline PetscComplex operator-(const PetscComplex& lhs, const PetscInt& rhs) { return lhs - double(rhs); } 202 static inline PetscComplex operator-(const PetscInt& lhs, const PetscComplex& rhs) { return double(lhs) - rhs; } 203 static inline PetscComplex operator*(const PetscComplex& lhs, const PetscInt& rhs) { return lhs * double(rhs); } 204 static inline PetscComplex operator*(const PetscInt& lhs, const PetscComplex& rhs) { return double(lhs) * rhs; } 205 static inline PetscComplex operator/(const PetscComplex& lhs, const PetscInt& rhs) { return lhs / double(rhs); } 206 static inline PetscComplex operator/(const PetscInt& lhs, const PetscComplex& rhs) { return double(lhs) / rhs; } 207 static inline bool operator==(const PetscComplex& lhs, const PetscInt& rhs) { return lhs.imag() == double(0) && lhs.real() == double(rhs); } 208 static inline bool operator==(const PetscInt& lhs, const PetscComplex& rhs) { return rhs.imag() == double(0) && rhs.real() == double(lhs); } 209 static inline bool operator!=(const PetscComplex& lhs, const PetscInt& rhs) { return lhs.imag() != double(0) || lhs.real() != double(rhs); } 210 static inline bool operator!=(const PetscInt& lhs, const PetscComplex& rhs) { return rhs.imag() != double(0) || rhs.real() != double(lhs); } 211 #endif /* PETSC_USE_CXX_COMPLEX_FLOAT_WORKAROUND */ 212 #elif defined(PETSC_USE_REAL___FLOAT128) 213 typedef complexlib::complex<__float128> PetscComplex; /* Notstandard and not expected to work, use __complex128 */ 214 PETSC_EXTERN MPI_Datatype MPIU___COMPLEX128; 215 #endif /* PETSC_USE_REAL_ */ 216 #endif /* ! PETSC_SKIP_COMPLEX */ 217 218 #elif defined(PETSC_HAVE_C99_COMPLEX) && !defined(PETSC_USE_REAL___FP16) 219 #if !defined(PETSC_SKIP_COMPLEX) 220 #define PETSC_HAVE_COMPLEX 1 221 #include <complex.h> 222 223 #if defined(PETSC_USE_REAL_SINGLE) || defined(PETSC_USE_REAL___FP16) 224 typedef float _Complex PetscComplex; 225 226 #define PetscRealPartComplex(a) crealf(a) 227 #define PetscImaginaryPartComplex(a) cimagf(a) 228 #define PetscAbsComplex(a) cabsf(a) 229 #define PetscConjComplex(a) conjf(a) 230 #define PetscSqrtComplex(a) csqrtf(a) 231 #define PetscPowComplex(a,b) cpowf(a,b) 232 #define PetscExpComplex(a) cexpf(a) 233 #define PetscLogComplex(a) clogf(a) 234 #define PetscSinComplex(a) csinf(a) 235 #define PetscCosComplex(a) ccosf(a) 236 #define PetscAsinComplex(a) casinf(a) 237 #define PetscAcosComplex(a) cacosf(a) 238 #define PetscTanComplex(a) ctanf(a) 239 #define PetscSinhComplex(a) csinhf(a) 240 #define PetscCoshComplex(a) ccoshf(a) 241 #define PetscTanhComplex(a) ctanhf(a) 242 243 #elif defined(PETSC_USE_REAL_DOUBLE) 244 typedef double _Complex PetscComplex; 245 246 #define PetscRealPartComplex(a) creal(a) 247 #define PetscImaginaryPartComplex(a) cimag(a) 248 #define PetscAbsComplex(a) cabs(a) 249 #define PetscConjComplex(a) conj(a) 250 #define PetscSqrtComplex(a) csqrt(a) 251 #define PetscPowComplex(a,b) cpow(a,b) 252 #define PetscExpComplex(a) cexp(a) 253 #define PetscLogComplex(a) clog(a) 254 #define PetscSinComplex(a) csin(a) 255 #define PetscCosComplex(a) ccos(a) 256 #define PetscAsinComplex(a) casin(a) 257 #define PetscAcosComplex(a) cacos(a) 258 #define PetscTanComplex(a) ctan(a) 259 #define PetscSinhComplex(a) csinh(a) 260 #define PetscCoshComplex(a) ccosh(a) 261 #define PetscTanhComplex(a) ctanh(a) 262 263 #elif defined(PETSC_USE_REAL___FLOAT128) 264 typedef __complex128 PetscComplex; 265 PETSC_EXTERN MPI_Datatype MPIU___COMPLEX128 PetscAttrMPITypeTag(__complex128); 266 267 #define PetscRealPartComplex(a) crealq(a) 268 #define PetscImaginaryPartComplex(a) cimagq(a) 269 #define PetscAbsComplex(a) cabsq(a) 270 #define PetscConjComplex(a) conjq(a) 271 #define PetscSqrtComplex(a) csqrtq(a) 272 #define PetscPowComplex(a,b) cpowq(a,b) 273 #define PetscExpComplex(a) cexpq(a) 274 #define PetscLogComplex(a) clogq(a) 275 #define PetscSinComplex(a) csinq(a) 276 #define PetscCosComplex(a) ccosq(a) 277 #define PetscAsinComplex(a) casinq(a) 278 #define PetscAcosComplex(a) cacosq(a) 279 #define PetscTanComplex(a) ctanq(a) 280 #define PetscSinhComplex(a) csinhq(a) 281 #define PetscCoshComplex(a) ccoshq(a) 282 #define PetscTanhComplex(a) ctanhq(a) 283 284 #endif /* PETSC_USE_REAL_* */ 285 #elif (defined(PETSC_USE_COMPLEX) && !defined(PETSC_SKIP_COMPLEX)) 286 #error "PETSc was configured --with-scalar-type=complex, but a language-appropriate complex library is not available" 287 #endif /* !PETSC_SKIP_COMPLEX */ 288 #endif /* (__cplusplus && PETSC_HAVE_CXX_COMPLEX) else-if (!__cplusplus && PETSC_HAVE_C99_COMPLEX) */ 289 290 #if defined(PETSC_HAVE_COMPLEX) 291 #if defined(PETSC_HAVE_MPI_C_DOUBLE_COMPLEX) 292 #define MPIU_C_DOUBLE_COMPLEX MPI_C_DOUBLE_COMPLEX 293 #define MPIU_C_COMPLEX MPI_C_COMPLEX 294 #else 295 # if defined(__cplusplus) && defined(PETSC_HAVE_CXX_COMPLEX) && !defined(PETSC_USE_REAL___FLOAT128) 296 typedef complexlib::complex<double> petsc_mpiu_c_double_complex; 297 typedef complexlib::complex<float> petsc_mpiu_c_complex; 298 # elif !defined(__cplusplus) && defined(PETSC_HAVE_C99_COMPLEX) 299 typedef double _Complex petsc_mpiu_c_double_complex; 300 typedef float _Complex petsc_mpiu_c_complex; 301 # else 302 typedef struct {double real,imag;} petsc_mpiu_c_double_complex; 303 typedef struct {float real,imag;} petsc_mpiu_c_complex; 304 # endif 305 PETSC_EXTERN MPI_Datatype MPIU_C_DOUBLE_COMPLEX PetscAttrMPITypeTagLayoutCompatible(petsc_mpiu_c_double_complex); 306 PETSC_EXTERN MPI_Datatype MPIU_C_COMPLEX PetscAttrMPITypeTagLayoutCompatible(petsc_mpiu_c_complex); 307 #endif /* PETSC_HAVE_MPI_C_DOUBLE_COMPLEX */ 308 #endif /* PETSC_HAVE_COMPLEX */ 309 310 #if defined(PETSC_HAVE_COMPLEX) 311 # if defined(PETSC_USE_REAL_SINGLE) || defined(PETSC_USE_REAL___FP16) 312 # define MPIU_COMPLEX MPIU_C_COMPLEX 313 # elif defined(PETSC_USE_REAL_DOUBLE) 314 # define MPIU_COMPLEX MPIU_C_DOUBLE_COMPLEX 315 # elif defined(PETSC_USE_REAL___FLOAT128) 316 # define MPIU_COMPLEX MPIU___COMPLEX128 317 # endif /* PETSC_USE_REAL_* */ 318 #endif 319 320 #if (defined(PETSC_USE_COMPLEX) && !defined(PETSC_SKIP_COMPLEX)) 321 typedef PetscComplex PetscScalar; 322 323 /*MC 324 PetscRealPart - Returns the real part of a PetscScalar 325 326 Synopsis: 327 #include <petscmath.h> 328 PetscScalar PetscRealPart(PetscScalar v) 329 330 Not Collective 331 332 Input Parameter: 333 . v - value to find the real part of 334 335 Level: beginner 336 337 .seealso: PetscScalar, PetscImaginaryPart(), PetscMax(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr() 338 339 M*/ 340 #define PetscRealPart(a) PetscRealPartComplex(a) 341 342 /*MC 343 PetscImaginaryPart - Returns the imaginary part of a PetscScalar 344 345 Synopsis: 346 #include <petscmath.h> 347 PetscScalar PetscImaginaryPart(PetscScalar v) 348 349 Not Collective 350 351 Input Parameter: 352 . v - value to find the imaginary part of 353 354 Level: beginner 355 356 Notes: 357 If PETSc was configured for real numbers then this always returns the value 0 358 359 .seealso: PetscScalar, PetscRealPart(), PetscMax(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr() 360 361 M*/ 362 #define PetscImaginaryPart(a) PetscImaginaryPartComplex(a) 363 364 #define PetscAbsScalar(a) PetscAbsComplex(a) 365 #define PetscConj(a) PetscConjComplex(a) 366 #define PetscSqrtScalar(a) PetscSqrtComplex(a) 367 #define PetscPowScalar(a,b) PetscPowComplex(a,b) 368 #define PetscExpScalar(a) PetscExpComplex(a) 369 #define PetscLogScalar(a) PetscLogComplex(a) 370 #define PetscSinScalar(a) PetscSinComplex(a) 371 #define PetscCosScalar(a) PetscCosComplex(a) 372 #define PetscAsinScalar(a) PetscAsinComplex(a) 373 #define PetscAcosScalar(a) PetscAcosComplex(a) 374 #define PetscTanScalar(a) PetscTanComplex(a) 375 #define PetscSinhScalar(a) PetscSinhComplex(a) 376 #define PetscCoshScalar(a) PetscCoshComplex(a) 377 #define PetscTanhScalar(a) PetscTanhComplex(a) 378 #define MPIU_SCALAR MPIU_COMPLEX 379 380 /* 381 real number definitions 382 */ 383 #else /* PETSC_USE_COMPLEX */ 384 typedef PetscReal PetscScalar; 385 #define MPIU_SCALAR MPIU_REAL 386 387 #define PetscRealPart(a) (a) 388 #define PetscImaginaryPart(a) ((PetscReal)0.) 389 PETSC_STATIC_INLINE PetscReal PetscAbsScalar(PetscScalar a) {return a < 0.0 ? -a : a;} 390 #define PetscConj(a) (a) 391 #if !defined(PETSC_USE_REAL___FLOAT128) && !defined(PETSC_USE_REAL___FP16) 392 #define PetscSqrtScalar(a) sqrt(a) 393 #define PetscPowScalar(a,b) pow(a,b) 394 #define PetscExpScalar(a) exp(a) 395 #define PetscLogScalar(a) log(a) 396 #define PetscSinScalar(a) sin(a) 397 #define PetscCosScalar(a) cos(a) 398 #define PetscAsinScalar(a) asin(a) 399 #define PetscAcosScalar(a) acos(a) 400 #define PetscTanScalar(a) tan(a) 401 #define PetscSinhScalar(a) sinh(a) 402 #define PetscCoshScalar(a) cosh(a) 403 #define PetscTanhScalar(a) tanh(a) 404 #elif defined(PETSC_USE_REAL___FP16) 405 #define PetscSqrtScalar(a) sqrtf(a) 406 #define PetscPowScalar(a,b) powf(a,b) 407 #define PetscExpScalar(a) expf(a) 408 #define PetscLogScalar(a) logf(a) 409 #define PetscSinScalar(a) sinf(a) 410 #define PetscCosScalar(a) cosf(a) 411 #define PetscAsinScalar(a) asinf(a) 412 #define PetscAcosScalar(a) acosf(a) 413 #define PetscTanScalar(a) tanf(a) 414 #define PetscSinhScalar(a) sinhf(a) 415 #define PetscCoshScalar(a) coshf(a) 416 #define PetscTanhScalar(a) tanhf(a) 417 #else /* PETSC_USE_REAL___FLOAT128 */ 418 #define PetscSqrtScalar(a) sqrtq(a) 419 #define PetscPowScalar(a,b) powq(a,b) 420 #define PetscExpScalar(a) expq(a) 421 #define PetscLogScalar(a) logq(a) 422 #define PetscSinScalar(a) sinq(a) 423 #define PetscCosScalar(a) cosq(a) 424 #define PetscAsinScalar(a) asinq(a) 425 #define PetscAcosScalar(a) acosq(a) 426 #define PetscTanScalar(a) tanq(a) 427 #define PetscSinhScalar(a) sinhq(a) 428 #define PetscCoshScalar(a) coshq(a) 429 #define PetscTanhScalar(a) tanhq(a) 430 #endif /* PETSC_USE_REAL___FLOAT128 */ 431 432 #endif /* PETSC_USE_COMPLEX */ 433 434 #define PetscSign(a) (((a) >= 0) ? ((a) == 0 ? 0 : 1) : -1) 435 #define PetscSignReal(a) (((a) >= 0.0) ? ((a) == 0.0 ? 0.0 : 1.0) : -1.0) 436 #define PetscAbs(a) (((a) >= 0) ? (a) : (-(a))) 437 438 /* --------------------------------------------------------------------------*/ 439 440 /* 441 Certain objects may be created using either single or double precision. 442 This is currently not used. 443 */ 444 typedef enum { PETSC_SCALAR_DOUBLE,PETSC_SCALAR_SINGLE, PETSC_SCALAR_LONG_DOUBLE, PETSC_SCALAR_HALF } PetscScalarPrecision; 445 446 #if defined(PETSC_HAVE_COMPLEX) 447 /* PETSC_i is the imaginary number, i */ 448 PETSC_EXTERN PetscComplex PETSC_i; 449 450 /* Try to do the right thing for complex number construction: see 451 452 http://www.open-std.org/jtc1/sc22/wg14/www/docs/n1464.htm 453 454 for details 455 */ 456 PETSC_STATIC_INLINE PetscComplex PetscCMPLX(PetscReal x, PetscReal y) 457 { 458 #if defined(__cplusplus) && defined(PETSC_HAVE_CXX_COMPLEX) && !defined(PETSC_USE_REAL___FLOAT128) 459 return PetscComplex(x,y); 460 #elif defined(_Imaginary_I) 461 return x + y * _Imaginary_I; 462 #else 463 { /* In both C99 and C11 (ISO/IEC 9899, Section 6.2.5), 464 465 "For each floating type there is a corresponding real type, which is always a real floating 466 type. For real floating types, it is the same type. For complex types, it is the type given 467 by deleting the keyword _Complex from the type name." 468 469 So type punning should be portable. */ 470 union { PetscComplex z; PetscReal f[2]; } uz; 471 472 uz.f[0] = x; 473 uz.f[1] = y; 474 return uz.z; 475 } 476 #endif 477 } 478 #endif 479 480 481 /*MC 482 PetscMin - Returns minimum of two numbers 483 484 Synopsis: 485 #include <petscmath.h> 486 type PetscMin(type v1,type v2) 487 488 Not Collective 489 490 Input Parameter: 491 + v1 - first value to find minimum of 492 - v2 - second value to find minimum of 493 494 Notes: 495 type can be integer or floating point value 496 497 Level: beginner 498 499 .seealso: PetscMax(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr() 500 501 M*/ 502 #define PetscMin(a,b) (((a)<(b)) ? (a) : (b)) 503 504 /*MC 505 PetscMax - Returns maxium of two numbers 506 507 Synopsis: 508 #include <petscmath.h> 509 type max PetscMax(type v1,type v2) 510 511 Not Collective 512 513 Input Parameter: 514 + v1 - first value to find maximum of 515 - v2 - second value to find maximum of 516 517 Notes: 518 type can be integer or floating point value 519 520 Level: beginner 521 522 .seealso: PetscMin(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr() 523 524 M*/ 525 #define PetscMax(a,b) (((a)<(b)) ? (b) : (a)) 526 527 /*MC 528 PetscClipInterval - Returns a number clipped to be within an interval 529 530 Synopsis: 531 #include <petscmath.h> 532 type clip PetscClipInterval(type x,type a,type b) 533 534 Not Collective 535 536 Input Parameter: 537 + x - value to use if within interval (a,b) 538 . a - lower end of interval 539 - b - upper end of interval 540 541 Notes: 542 type can be integer or floating point value 543 544 Level: beginner 545 546 .seealso: PetscMin(), PetscMax(), PetscAbsInt(), PetscAbsReal(), PetscSqr() 547 548 M*/ 549 #define PetscClipInterval(x,a,b) (PetscMax((a),PetscMin((x),(b)))) 550 551 /*MC 552 PetscAbsInt - Returns the absolute value of an integer 553 554 Synopsis: 555 #include <petscmath.h> 556 int abs PetscAbsInt(int v1) 557 558 Not Collective 559 560 Input Parameter: 561 . v1 - the integer 562 563 Level: beginner 564 565 .seealso: PetscMax(), PetscMin(), PetscAbsReal(), PetscSqr() 566 567 M*/ 568 #define PetscAbsInt(a) (((a)<0) ? (-(a)) : (a)) 569 570 /*MC 571 PetscAbsReal - Returns the absolute value of an real number 572 573 Synopsis: 574 #include <petscmath.h> 575 Real abs PetscAbsReal(PetscReal v1) 576 577 Not Collective 578 579 Input Parameter: 580 . v1 - the double 581 582 583 Level: beginner 584 585 .seealso: PetscMax(), PetscMin(), PetscAbsInt(), PetscSqr() 586 587 M*/ 588 #if defined(PETSC_USE_REAL_SINGLE) 589 #define PetscAbsReal(a) fabsf(a) 590 #elif defined(PETSC_USE_REAL_DOUBLE) 591 #define PetscAbsReal(a) fabs(a) 592 #elif defined(PETSC_USE_REAL___FLOAT128) 593 #define PetscAbsReal(a) fabsq(a) 594 #elif defined(PETSC_USE_REAL___FP16) 595 #define PetscAbsReal(a) fabsf(a) 596 #endif 597 598 /*MC 599 PetscSqr - Returns the square of a number 600 601 Synopsis: 602 #include <petscmath.h> 603 type sqr PetscSqr(type v1) 604 605 Not Collective 606 607 Input Parameter: 608 . v1 - the value 609 610 Notes: 611 type can be integer or floating point value 612 613 Level: beginner 614 615 .seealso: PetscMax(), PetscMin(), PetscAbsInt(), PetscAbsReal() 616 617 M*/ 618 #define PetscSqr(a) ((a)*(a)) 619 620 /* ----------------------------------------------------------------------------*/ 621 622 #if defined(PETSC_USE_REAL_SINGLE) 623 #define PetscRealConstant(constant) constant##F 624 #elif defined(PETSC_USE_REAL___FLOAT128) 625 #define PetscRealConstant(constant) constant##Q 626 #else 627 #define PetscRealConstant(constant) constant 628 #endif 629 630 /* 631 Basic constants 632 */ 633 #define PETSC_PI PetscRealConstant(3.1415926535897932384626433832795029) 634 #define PETSC_PHI PetscRealConstant(1.6180339887498948482045868343656381) 635 636 #if !defined(PETSC_USE_64BIT_INDICES) 637 #define PETSC_MAX_INT 2147483647 638 #define PETSC_MIN_INT (-PETSC_MAX_INT - 1) 639 #else 640 #define PETSC_MAX_INT 9223372036854775807L 641 #define PETSC_MIN_INT (-PETSC_MAX_INT - 1) 642 #endif 643 644 #if defined(PETSC_USE_REAL_SINGLE) 645 # define PETSC_MAX_REAL 3.40282346638528860e+38F 646 # define PETSC_MIN_REAL (-PETSC_MAX_REAL) 647 # define PETSC_MACHINE_EPSILON 1.19209290e-07F 648 # define PETSC_SQRT_MACHINE_EPSILON 3.45266983e-04F 649 # define PETSC_SMALL 1.e-5F 650 #elif defined(PETSC_USE_REAL_DOUBLE) 651 # define PETSC_MAX_REAL 1.7976931348623157e+308 652 # define PETSC_MIN_REAL (-PETSC_MAX_REAL) 653 # define PETSC_MACHINE_EPSILON 2.2204460492503131e-16 654 # define PETSC_SQRT_MACHINE_EPSILON 1.490116119384766e-08 655 # define PETSC_SMALL 1.e-10 656 #elif defined(PETSC_USE_REAL___FLOAT128) 657 # define PETSC_MAX_REAL FLT128_MAX 658 # define PETSC_MIN_REAL (-FLT128_MAX) 659 # define PETSC_MACHINE_EPSILON FLT128_EPSILON 660 # define PETSC_SQRT_MACHINE_EPSILON 1.38777878078144567552953958511352539e-17Q 661 # define PETSC_SMALL 1.e-20Q 662 #elif defined(PETSC_USE_REAL___FP16) /* maybe should use single precision values for these? */ 663 # define PETSC_MAX_REAL 65504. 664 # define PETSC_MIN_REAL (-PETSC_MAX_REAL) 665 # define PETSC_MACHINE_EPSILON .00097656 666 # define PETSC_SQRT_MACHINE_EPSILON .0312 667 # define PETSC_SMALL 5.e-3 668 #endif 669 670 #define PETSC_INFINITY (PETSC_MAX_REAL/4) 671 #define PETSC_NINFINITY (-PETSC_INFINITY) 672 673 PETSC_EXTERN PetscBool PetscIsInfReal(PetscReal); 674 PETSC_EXTERN PetscBool PetscIsNanReal(PetscReal); 675 PETSC_EXTERN PetscBool PetscIsNormalReal(PetscReal); 676 PETSC_STATIC_INLINE PetscBool PetscIsInfOrNanReal(PetscReal v) {return PetscIsInfReal(v) || PetscIsNanReal(v) ? PETSC_TRUE : PETSC_FALSE;} 677 PETSC_STATIC_INLINE PetscBool PetscIsInfScalar(PetscScalar v) {return PetscIsInfReal(PetscAbsScalar(v));} 678 PETSC_STATIC_INLINE PetscBool PetscIsNanScalar(PetscScalar v) {return PetscIsNanReal(PetscAbsScalar(v));} 679 PETSC_STATIC_INLINE PetscBool PetscIsInfOrNanScalar(PetscScalar v) {return PetscIsInfOrNanReal(PetscAbsScalar(v));} 680 PETSC_STATIC_INLINE PetscBool PetscIsNormalScalar(PetscScalar v) {return PetscIsNormalReal(PetscAbsScalar(v));} 681 682 PETSC_EXTERN PetscBool PetscEqualReal(PetscReal,PetscReal); 683 PETSC_EXTERN PetscBool PetscEqualScalar(PetscScalar,PetscScalar); 684 685 /* 686 These macros are currently hardwired to match the regular data types, so there is no support for a different 687 MatScalar from PetscScalar. We left the MatScalar in the source just in case we use it again. 688 */ 689 #define MPIU_MATSCALAR MPIU_SCALAR 690 typedef PetscScalar MatScalar; 691 typedef PetscReal MatReal; 692 693 struct petsc_mpiu_2scalar {PetscScalar a,b;}; 694 PETSC_EXTERN MPI_Datatype MPIU_2SCALAR PetscAttrMPITypeTagLayoutCompatible(struct petsc_mpiu_2scalar); 695 #if defined(PETSC_USE_64BIT_INDICES) || !defined(MPI_2INT) 696 struct petsc_mpiu_2int {PetscInt a,b;}; 697 PETSC_EXTERN MPI_Datatype MPIU_2INT PetscAttrMPITypeTagLayoutCompatible(struct petsc_mpiu_2int); 698 #else 699 #define MPIU_2INT MPI_2INT 700 #endif 701 702 PETSC_STATIC_INLINE PetscInt PetscPowInt(PetscInt base,PetscInt power) 703 { 704 PetscInt result = 1; 705 while (power) { 706 if (power & 1) result *= base; 707 power >>= 1; 708 base *= base; 709 } 710 return result; 711 } 712 713 PETSC_STATIC_INLINE PetscReal PetscPowRealInt(PetscReal base,PetscInt power) 714 { 715 PetscReal result = 1; 716 if (power < 0) { 717 power = -power; 718 base = ((PetscReal)1)/base; 719 } 720 while (power) { 721 if (power & 1) result *= base; 722 power >>= 1; 723 base *= base; 724 } 725 return result; 726 } 727 728 PETSC_STATIC_INLINE PetscScalar PetscPowScalarInt(PetscScalar base,PetscInt power) 729 { 730 PetscScalar result = 1; 731 if (power < 0) { 732 power = -power; 733 base = ((PetscReal)1)/base; 734 } 735 while (power) { 736 if (power & 1) result *= base; 737 power >>= 1; 738 base *= base; 739 } 740 return result; 741 } 742 743 PETSC_STATIC_INLINE PetscScalar PetscPowScalarReal(PetscScalar base,PetscReal power) 744 { 745 PetscScalar cpower = power; 746 return PetscPowScalar(base,cpower); 747 } 748 749 #ifndef PETSC_HAVE_LOG2 750 PETSC_STATIC_INLINE PetscReal PetscLog2Real(PetscReal n) 751 { 752 return PetscLogReal(n)/PetscLogReal(2.0); 753 } 754 #endif 755 #endif 756