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 #define PetscExpPassiveScalar(a) PetscExpScalar() 25 #if defined(PETSC_USE_REAL_SINGLE) 26 #define MPIU_REAL MPI_FLOAT 27 typedef float PetscReal; 28 #define PetscSqrtReal(a) sqrt(a) 29 #define PetscExpReal(a) exp(a) 30 #define PetscLogReal(a) log(a) 31 #define PetscSinReal(a) sin(a) 32 #define PetscCosReal(a) cos(a) 33 #elif defined(PETSC_USE_REAL_DOUBLE) 34 #define MPIU_REAL MPI_DOUBLE 35 typedef double PetscReal; 36 #define PetscSqrtReal(a) sqrt(a) 37 #define PetscExpReal(a) exp(a) 38 #define PetscLogReal(a) log(a) 39 #define PetscSinReal(a) sin(a) 40 #define PetscCosReal(a) cos(a) 41 #elif defined(PETSC_USE_REAL___FLOAT128) 42 #if defined(__cplusplus) 43 extern "C" { 44 #endif 45 #include <quadmath.h> 46 #if defined(__cplusplus) 47 } 48 #endif 49 PETSC_EXTERN MPI_Datatype MPIU___FLOAT128 PetscAttrMPITypeTag(__float128); 50 #define MPIU_REAL MPIU___FLOAT128 51 typedef __float128 PetscReal; 52 #define PetscSqrtReal(a) sqrtq(a) 53 #define PetscExpReal(a) expq(a) 54 #define PetscLogReal(a) logq(a) 55 #define PetscSinReal(a) sinq(a) 56 #define PetscCosReal(a) cosq(a) 57 #endif /* PETSC_USE_REAL_* */ 58 59 /* 60 Complex number definitions 61 */ 62 #if defined(PETSC_CLANGUAGE_CXX) && defined(PETSC_HAVE_CXX_COMPLEX) 63 #if defined(PETSC_USE_COMPLEX) || defined(PETSC_DESIRE_COMPLEX) 64 #define PETSC_HAVE_COMPLEX 1 65 /* C++ support of complex number */ 66 #if defined(PETSC_HAVE_CUSP) 67 #define complexlib cusp 68 #include <cusp/complex.h> 69 #else 70 #define complexlib std 71 #include <complex> 72 #endif 73 74 #define PetscRealPartComplex(a) (a).real() 75 #define PetscImaginaryPartComplex(a) (a).imag() 76 #define PetscAbsComplex(a) complexlib::abs(a) 77 #define PetscConjComplex(a) complexlib::conj(a) 78 #define PetscSqrtComplex(a) complexlib::sqrt(a) 79 #define PetscPowComplex(a,b) complexlib::pow(a,b) 80 #define PetscExpComplex(a) complexlib::exp(a) 81 #define PetscLogComplex(a) complexlib::log(a) 82 #define PetscSinComplex(a) complexlib::sin(a) 83 #define PetscCosComplex(a) complexlib::cos(a) 84 85 #if defined(PETSC_USE_REAL_SINGLE) 86 typedef complexlib::complex<float> PetscComplex; 87 #elif defined(PETSC_USE_REAL_DOUBLE) 88 typedef complexlib::complex<double> PetscComplex; 89 #elif defined(PETSC_USE_REAL___FLOAT128) 90 typedef complexlib::complex<__float128> PetscComplex; /* Notstandard and not expected to work, use __complex128 */ 91 #endif /* PETSC_USE_REAL_ */ 92 #endif /* PETSC_USE_COMPLEX && PETSC_DESIRE_COMPLEX */ 93 94 #elif defined(PETSC_CLANGUAGE_C) && defined(PETSC_HAVE_C99_COMPLEX) 95 /* Use C99 _Complex for the type. Do not include complex.h by default to define "complex" because of symbol conflicts in Hypre. */ 96 /* Compilation units that can safely use complex should define PETSC_DESIRE_COMPLEX before including any headers */ 97 #if defined(PETSC_USE_COMPLEX) || defined(PETSC_DESIRE_COMPLEX) 98 #define PETSC_HAVE_COMPLEX 1 99 #include <complex.h> 100 101 #if defined(PETSC_USE_REAL_SINGLE) 102 typedef float _Complex PetscComplex; 103 104 #define PetscRealPartComplex(a) crealf(a) 105 #define PetscImaginaryPartComplex(a) cimagf(a) 106 #define PetscAbsComplex(a) cabsf(a) 107 #define PetscConjComplex(a) conjf(a) 108 #define PetscSqrtComplex(a) csqrtf(a) 109 #define PetscPowComplex(a,b) cpowf(a,b) 110 #define PetscExpComplex(a) cexpf(a) 111 #define PetscLogComplex(a) clogf(a) 112 #define PetscSinComplex(a) csinf(a) 113 #define PetscCosComplex(a) ccosf(a) 114 115 #elif defined(PETSC_USE_REAL_DOUBLE) 116 typedef double _Complex PetscComplex; 117 118 #define PetscRealPartComplex(a) creal(a) 119 #define PetscImaginaryPartComplex(a) cimag(a) 120 #define PetscAbsComplex(a) cabs(a) 121 #define PetscConjComplex(a) conj(a) 122 #define PetscSqrtComplex(a) csqrt(a) 123 #define PetscPowComplex(a,b) cpow(a,b) 124 #define PetscExpComplex(a) cexp(a) 125 #define PetscLogComplex(a) clog(a) 126 #define PetscSinComplex(a) csin(a) 127 #define PetscCosComplex(a) ccos(a) 128 129 #elif defined(PETSC_USE_REAL___FLOAT128) 130 typedef __complex128 PetscComplex; 131 PETSC_EXTERN MPI_Datatype MPIU___COMPLEX128 PetscAttrMPITypeTag(__complex128); 132 133 #define PetscRealPartComplex(a) crealq(a) 134 #define PetscImaginaryPartComplex(a) cimagq(a) 135 #define PetscAbsComplex(a) cabsq(a) 136 #define PetscConjComplex(a) conjq(a) 137 #define PetscSqrtComplex(a) csqrtq(a) 138 #define PetscPowComplex(a,b) cpowq(a,b) 139 #define PetscExpComplex(a) cexpq(a) 140 #define PetscLogComplex(a) clogq(a) 141 #define PetscSinComplex(a) csinq(a) 142 #define PetscCosComplex(a) ccosq(a) 143 #endif /* PETSC_USE_REAL_* */ 144 #elif defined(PETSC_USE_COMPLEX) 145 #error "PETSc was configured --with-scalar-type=complex, but a language-appropriate complex library is not available" 146 #endif /* PETSC_USE_COMPLEX || PETSC_DESIRE_COMPLEX */ 147 #endif /* (PETSC_CLANGUAGE_CXX && PETSC_HAVE_CXX_COMPLEX) else-if (PETSC_CLANGUAGE_C && PETSC_HAVE_C99_COMPLEX) */ 148 149 #if defined(PETSC_HAVE_MPI_C_DOUBLE_COMPLEX) 150 #define MPIU_C_DOUBLE_COMPLEX MPI_C_DOUBLE_COMPLEX 151 #define MPIU_C_COMPLEX MPI_C_COMPLEX 152 #else 153 # if defined(PETSC_CLANGUAGE_CXX) && defined(PETSC_HAVE_CXX_COMPLEX) 154 typedef complexlib::complex<double> petsc_mpiu_c_double_complex; 155 typedef complexlib::complex<float> petsc_mpiu_c_complex; 156 # elif defined(PETSC_CLANGUAGE_C) && defined(PETSC_HAVE_C99_COMPLEX) 157 typedef double _Complex petsc_mpiu_c_double_complex; 158 typedef float _Complex petsc_mpiu_c_complex; 159 # else 160 typedef struct {double real,imag;} petsc_mpiu_c_double_complex; 161 typedef struct {float real,imag;} petsc_mpiu_c_complex; 162 # endif 163 PETSC_EXTERN MPI_Datatype MPIU_C_DOUBLE_COMPLEX PetscAttrMPITypeTagLayoutCompatible(petsc_mpiu_c_double_complex); 164 PETSC_EXTERN MPI_Datatype MPIU_C_COMPLEX PetscAttrMPITypeTagLayoutCompatible(petsc_mpiu_c_complex); 165 #endif /* PETSC_HAVE_MPI_C_DOUBLE_COMPLEX */ 166 167 #if defined(PETSC_HAVE_COMPLEX) 168 # if defined(PETSC_USE_REAL_SINGLE) 169 # define MPIU_COMPLEX MPIU_C_COMPLEX 170 # elif defined(PETSC_USE_REAL_DOUBLE) 171 # define MPIU_COMPLEX MPIU_C_DOUBLE_COMPLEX 172 # elif defined(PETSC_USE_REAL___FLOAT128) 173 # define MPIU_COMPLEX MPIU___COMPLEX128 174 # endif /* PETSC_USE_REAL_* */ 175 #endif 176 177 #if defined(PETSC_USE_COMPLEX) 178 typedef PetscComplex PetscScalar; 179 #define PetscRealPart(a) PetscRealPartComplex(a) 180 #define PetscImaginaryPart(a) PetscImaginaryPartComplex(a) 181 #define PetscAbsScalar(a) PetscAbsComplex(a) 182 #define PetscConj(a) PetscConjComplex(a) 183 #define PetscSqrtScalar(a) PetscSqrtComplex(a) 184 #define PetscPowScalar(a,b) PetscPowComplex(a,b) 185 #define PetscExpScalar(a) PetscExpComplex(a) 186 #define PetscLogScalar(a) PetscLogComplex(a) 187 #define PetscSinScalar(a) PetscSinComplex(a) 188 #define PetscCosScalar(a) PetscCosComplex(a) 189 190 #define MPIU_SCALAR MPIU_COMPLEX 191 192 /* 193 real number definitions 194 */ 195 #else /* PETSC_USE_COMPLEX */ 196 typedef PetscReal PetscScalar; 197 #define MPIU_SCALAR MPIU_REAL 198 199 #define PetscRealPart(a) (a) 200 #define PetscImaginaryPart(a) ((PetscReal)0.) 201 PETSC_STATIC_INLINE PetscReal PetscAbsScalar(PetscScalar a) {return a < 0.0 ? -a : a;} 202 #define PetscConj(a) (a) 203 #if !defined(PETSC_USE_REAL___FLOAT128) 204 #define PetscSqrtScalar(a) sqrt(a) 205 #define PetscPowScalar(a,b) pow(a,b) 206 #define PetscExpScalar(a) exp(a) 207 #define PetscLogScalar(a) log(a) 208 #define PetscSinScalar(a) sin(a) 209 #define PetscCosScalar(a) cos(a) 210 #else /* PETSC_USE_REAL___FLOAT128 */ 211 #define PetscSqrtScalar(a) sqrtq(a) 212 #define PetscPowScalar(a,b) powq(a,b) 213 #define PetscExpScalar(a) expq(a) 214 #define PetscLogScalar(a) logq(a) 215 #define PetscSinScalar(a) sinq(a) 216 #define PetscCosScalar(a) cosq(a) 217 #endif /* PETSC_USE_REAL___FLOAT128 */ 218 219 #endif /* PETSC_USE_COMPLEX */ 220 221 #define PetscSign(a) (((a) >= 0) ? ((a) == 0 ? 0 : 1) : -1) 222 #define PetscAbs(a) (((a) >= 0) ? (a) : -(a)) 223 224 /* --------------------------------------------------------------------------*/ 225 226 /* 227 Certain objects may be created using either single or double precision. 228 This is currently not used. 229 */ 230 typedef enum { PETSC_SCALAR_DOUBLE,PETSC_SCALAR_SINGLE, PETSC_SCALAR_LONG_DOUBLE } PetscScalarPrecision; 231 232 #if defined(PETSC_HAVE_COMPLEX) 233 /* PETSC_i is the imaginary number, i */ 234 PETSC_EXTERN PetscComplex PETSC_i; 235 #endif 236 237 /*MC 238 PetscMin - Returns minimum of two numbers 239 240 Synopsis: 241 type PetscMin(type v1,type v2) 242 243 Not Collective 244 245 Input Parameter: 246 + v1 - first value to find minimum of 247 - v2 - second value to find minimum of 248 249 250 Notes: type can be integer or floating point value 251 252 Level: beginner 253 254 255 .seealso: PetscMin(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr() 256 257 M*/ 258 #define PetscMin(a,b) (((a)<(b)) ? (a) : (b)) 259 260 /*MC 261 PetscMax - Returns maxium of two numbers 262 263 Synopsis: 264 type max PetscMax(type v1,type v2) 265 266 Not Collective 267 268 Input Parameter: 269 + v1 - first value to find maximum of 270 - v2 - second value to find maximum of 271 272 Notes: type can be integer or floating point value 273 274 Level: beginner 275 276 .seealso: PetscMin(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr() 277 278 M*/ 279 #define PetscMax(a,b) (((a)<(b)) ? (b) : (a)) 280 281 /*MC 282 PetscClipInterval - Returns a number clipped to be within an interval 283 284 Synopsis: 285 type clip PetscClipInterval(type x,type a,type b) 286 287 Not Collective 288 289 Input Parameter: 290 + x - value to use if within interval (a,b) 291 . a - lower end of interval 292 - b - upper end of interval 293 294 Notes: type can be integer or floating point value 295 296 Level: beginner 297 298 .seealso: PetscMin(), PetscMax(), PetscAbsInt(), PetscAbsReal(), PetscSqr() 299 300 M*/ 301 #define PetscClipInterval(x,a,b) (PetscMax((a),PetscMin((x),(b)))) 302 303 /*MC 304 PetscAbsInt - Returns the absolute value of an integer 305 306 Synopsis: 307 int abs PetscAbsInt(int v1) 308 309 Not Collective 310 311 Input Parameter: 312 . v1 - the integer 313 314 Level: beginner 315 316 .seealso: PetscMax(), PetscMin(), PetscAbsReal(), PetscSqr() 317 318 M*/ 319 #define PetscAbsInt(a) (((a)<0) ? -(a) : (a)) 320 321 /*MC 322 PetscAbsReal - Returns the absolute value of an real number 323 324 Synopsis: 325 Real abs PetscAbsReal(PetscReal v1) 326 327 Not Collective 328 329 Input Parameter: 330 . v1 - the double 331 332 333 Level: beginner 334 335 .seealso: PetscMax(), PetscMin(), PetscAbsInt(), PetscSqr() 336 337 M*/ 338 #define PetscAbsReal(a) (((a)<0) ? -(a) : (a)) 339 340 /*MC 341 PetscSqr - Returns the square of a number 342 343 Synopsis: 344 type sqr PetscSqr(type v1) 345 346 Not Collective 347 348 Input Parameter: 349 . v1 - the value 350 351 Notes: type can be integer or floating point value 352 353 Level: beginner 354 355 .seealso: PetscMax(), PetscMin(), PetscAbsInt(), PetscAbsReal() 356 357 M*/ 358 #define PetscSqr(a) ((a)*(a)) 359 360 /* ----------------------------------------------------------------------------*/ 361 /* 362 Basic constants 363 */ 364 #if defined(PETSC_USE_REAL___FLOAT128) 365 #define PETSC_PI M_PIq 366 #elif defined(M_PI) 367 #define PETSC_PI M_PI 368 #else 369 #define PETSC_PI 3.14159265358979323846264338327950288419716939937510582 370 #endif 371 372 #if !defined(PETSC_USE_64BIT_INDICES) 373 #define PETSC_MAX_INT 2147483647 374 #define PETSC_MIN_INT (-PETSC_MAX_INT - 1) 375 #else 376 #define PETSC_MAX_INT 9223372036854775807L 377 #define PETSC_MIN_INT (-PETSC_MAX_INT - 1) 378 #endif 379 380 #if defined(PETSC_USE_REAL_SINGLE) 381 # define PETSC_MAX_REAL 3.40282346638528860e+38F 382 # define PETSC_MIN_REAL -PETSC_MAX_REAL 383 # define PETSC_MACHINE_EPSILON 1.19209290e-07F 384 # define PETSC_SQRT_MACHINE_EPSILON 3.45266983e-04F 385 # define PETSC_SMALL 1.e-5 386 #elif defined(PETSC_USE_REAL_DOUBLE) 387 # define PETSC_MAX_REAL 1.7976931348623157e+308 388 # define PETSC_MIN_REAL -PETSC_MAX_REAL 389 # define PETSC_MACHINE_EPSILON 2.2204460492503131e-16 390 # define PETSC_SQRT_MACHINE_EPSILON 1.490116119384766e-08 391 # define PETSC_SMALL 1.e-10 392 #elif defined(PETSC_USE_REAL___FLOAT128) 393 # define PETSC_MAX_REAL FLT128_MAX 394 # define PETSC_MIN_REAL -FLT128_MAX 395 # define PETSC_MACHINE_EPSILON FLT128_EPSILON 396 # define PETSC_SQRT_MACHINE_EPSILON 1.38777878078e-17 397 # define PETSC_SMALL 1.e-20 398 #endif 399 400 PETSC_EXTERN PetscErrorCode PetscIsInfOrNanScalar(PetscScalar); 401 PETSC_EXTERN PetscErrorCode PetscIsInfOrNanReal(PetscReal); 402 403 /* ----------------------------------------------------------------------------*/ 404 #define PassiveReal PetscReal 405 #define PassiveScalar PetscScalar 406 407 /* 408 These macros are currently hardwired to match the regular data types, so there is no support for a different 409 MatScalar from PetscScalar. We left the MatScalar in the source just in case we use it again. 410 */ 411 #define MPIU_MATSCALAR MPIU_SCALAR 412 typedef PetscScalar MatScalar; 413 typedef PetscReal MatReal; 414 415 struct petsc_mpiu_2scalar {PetscScalar a,b;}; 416 PETSC_EXTERN MPI_Datatype MPIU_2SCALAR PetscAttrMPITypeTagLayoutCompatible(struct petsc_mpiu_2scalar); 417 #if defined(PETSC_USE_64BIT_INDICES) || !defined(MPI_2INT) 418 struct petsc_mpiu_2int {PetscInt a,b;}; 419 PETSC_EXTERN MPI_Datatype MPIU_2INT PetscAttrMPITypeTagLayoutCompatible(struct petsc_mpiu_2int); 420 #else 421 #define MPIU_2INT MPI_2INT 422 #endif 423 424 #endif 425