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