Actual source code: petscmath.h
petsc-3.9.4 2018-09-11
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: */
13: #include <math.h>
15: /*
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.
22: */
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_* */
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
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
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 */
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>
223: #if defined(PETSC_USE_REAL_SINGLE) || defined(PETSC_USE_REAL___FP16)
224: typedef float _Complex PetscComplex;
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)
243: #elif defined(PETSC_USE_REAL_DOUBLE)
244: typedef double _Complex PetscComplex;
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)
263: #elif defined(PETSC_USE_REAL___FLOAT128)
264: typedef __complex128 PetscComplex;
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: #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) */
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 */
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
320: #if (defined(PETSC_USE_COMPLEX) && !defined(PETSC_SKIP_COMPLEX))
321: typedef PetscComplex PetscScalar;
323: /*MC
324: PetscRealPart - Returns the real part of a PetscScalar
326: Synopsis:
327: #include <petscmath.h>
328: PetscScalar PetscRealPart(PetscScalar v)
330: Not Collective
332: Input Parameter:
333: . v - value to find the real part of
335: Level: beginner
337: .seealso: PetscScalar, PetscImaginaryPart(), PetscMax(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr()
339: M*/
340: #define PetscRealPart(a) PetscRealPartComplex(a)
342: /*MC
343: PetscImaginaryPart - Returns the imaginary part of a PetscScalar
345: Synopsis:
346: #include <petscmath.h>
347: PetscScalar PetscImaginaryPart(PetscScalar v)
349: Not Collective
351: Input Parameter:
352: . v - value to find the imaginary part of
354: Level: beginner
356: Notes:
357: If PETSc was configured for real numbers then this always returns the value 0
359: .seealso: PetscScalar, PetscRealPart(), PetscMax(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr()
361: M*/
362: #define PetscImaginaryPart(a) PetscImaginaryPartComplex(a)
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
380: /*
381: real number definitions
382: */
383: #else /* PETSC_USE_COMPLEX */
384: typedef PetscReal PetscScalar;
385: #define MPIU_SCALAR MPIU_REAL
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 */
432: #endif /* PETSC_USE_COMPLEX */
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)))
438: /* --------------------------------------------------------------------------*/
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;
446: #if defined(PETSC_HAVE_COMPLEX)
447: /* PETSC_i is the imaginary number, i */
448: PETSC_EXTERN PetscComplex PETSC_i;
450: /* Try to do the right thing for complex number construction: see
452: http://www.open-std.org/jtc1/sc22/wg14/www/docs/n1464.htm
454: for details
455: */
456: PETSC_STATIC_INLINE PetscComplex PetscCMPLX(PetscReal x, PetscReal y)
457: {
458: #if defined(__cplusplus) && !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),
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."
469: So type punning should be portable. */
470: union { PetscComplex z; PetscReal f[2]; } uz;
472: uz.f[0] = x;
473: uz.f[1] = y;
474: return uz.z;
475: }
476: #endif
477: }
478: #endif
481: /*MC
482: PetscMin - Returns minimum of two numbers
484: Synopsis:
485: #include <petscmath.h>
486: type PetscMin(type v1,type v2)
488: Not Collective
490: Input Parameter:
491: + v1 - first value to find minimum of
492: - v2 - second value to find minimum of
494: Notes: type can be integer or floating point value
496: Level: beginner
498: .seealso: PetscMax(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr()
500: M*/
501: #define PetscMin(a,b) (((a)<(b)) ? (a) : (b))
503: /*MC
504: PetscMax - Returns maxium of two numbers
506: Synopsis:
507: #include <petscmath.h>
508: type max PetscMax(type v1,type v2)
510: Not Collective
512: Input Parameter:
513: + v1 - first value to find maximum of
514: - v2 - second value to find maximum of
516: Notes: type can be integer or floating point value
518: Level: beginner
520: .seealso: PetscMin(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr()
522: M*/
523: #define PetscMax(a,b) (((a)<(b)) ? (b) : (a))
525: /*MC
526: PetscClipInterval - Returns a number clipped to be within an interval
528: Synopsis:
529: #include <petscmath.h>
530: type clip PetscClipInterval(type x,type a,type b)
532: Not Collective
534: Input Parameter:
535: + x - value to use if within interval (a,b)
536: . a - lower end of interval
537: - b - upper end of interval
539: Notes: type can be integer or floating point value
541: Level: beginner
543: .seealso: PetscMin(), PetscMax(), PetscAbsInt(), PetscAbsReal(), PetscSqr()
545: M*/
546: #define PetscClipInterval(x,a,b) (PetscMax((a),PetscMin((x),(b))))
548: /*MC
549: PetscAbsInt - Returns the absolute value of an integer
551: Synopsis:
552: #include <petscmath.h>
553: int abs PetscAbsInt(int v1)
555: Not Collective
557: Input Parameter:
558: . v1 - the integer
560: Level: beginner
562: .seealso: PetscMax(), PetscMin(), PetscAbsReal(), PetscSqr()
564: M*/
565: #define PetscAbsInt(a) (((a)<0) ? (-(a)) : (a))
567: /*MC
568: PetscAbsReal - Returns the absolute value of an real number
570: Synopsis:
571: #include <petscmath.h>
572: Real abs PetscAbsReal(PetscReal v1)
574: Not Collective
576: Input Parameter:
577: . v1 - the double
580: Level: beginner
582: .seealso: PetscMax(), PetscMin(), PetscAbsInt(), PetscSqr()
584: M*/
585: #if defined(PETSC_USE_REAL_SINGLE)
586: #define PetscAbsReal(a) fabsf(a)
587: #elif defined(PETSC_USE_REAL_DOUBLE)
588: #define PetscAbsReal(a) fabs(a)
589: #elif defined(PETSC_USE_REAL___FLOAT128)
590: #define PetscAbsReal(a) fabsq(a)
591: #elif defined(PETSC_USE_REAL___FP16)
592: #define PetscAbsReal(a) fabsf(a)
593: #endif
595: /*MC
596: PetscSqr - Returns the square of a number
598: Synopsis:
599: #include <petscmath.h>
600: type sqr PetscSqr(type v1)
602: Not Collective
604: Input Parameter:
605: . v1 - the value
607: Notes: type can be integer or floating point value
609: Level: beginner
611: .seealso: PetscMax(), PetscMin(), PetscAbsInt(), PetscAbsReal()
613: M*/
614: #define PetscSqr(a) ((a)*(a))
616: /* ----------------------------------------------------------------------------*/
618: #if defined(PETSC_USE_REAL_SINGLE)
619: #define PetscRealConstant(constant) constant##F
620: #elif defined(PETSC_USE_REAL___FLOAT128)
621: #define PetscRealConstant(constant) constant##Q
622: #else
623: #define PetscRealConstant(constant) constant
624: #endif
626: /*
627: Basic constants
628: */
629: #define PETSC_PI PetscRealConstant(3.1415926535897932384626433832795029)
630: #define PETSC_PHI PetscRealConstant(1.6180339887498948482045868343656381)
632: #if !defined(PETSC_USE_64BIT_INDICES)
633: #define PETSC_MAX_INT 2147483647
634: #define PETSC_MIN_INT (-PETSC_MAX_INT - 1)
635: #else
636: #define PETSC_MAX_INT 9223372036854775807L
637: #define PETSC_MIN_INT (-PETSC_MAX_INT - 1)
638: #endif
640: #if defined(PETSC_USE_REAL_SINGLE)
641: # define PETSC_MAX_REAL 3.40282346638528860e+38F
642: # define PETSC_MIN_REAL (-PETSC_MAX_REAL)
643: # define PETSC_MACHINE_EPSILON 1.19209290e-07F
644: # define PETSC_SQRT_MACHINE_EPSILON 3.45266983e-04F
645: # define PETSC_SMALL 1.e-5F
646: #elif defined(PETSC_USE_REAL_DOUBLE)
647: # define PETSC_MAX_REAL 1.7976931348623157e+308
648: # define PETSC_MIN_REAL (-PETSC_MAX_REAL)
649: # define PETSC_MACHINE_EPSILON 2.2204460492503131e-16
650: # define PETSC_SQRT_MACHINE_EPSILON 1.490116119384766e-08
651: # define PETSC_SMALL 1.e-10
652: #elif defined(PETSC_USE_REAL___FLOAT128)
653: # define PETSC_MAX_REAL FLT128_MAX
654: # define PETSC_MIN_REAL (-FLT128_MAX)
655: # define PETSC_MACHINE_EPSILON FLT128_EPSILON
656: # define PETSC_SQRT_MACHINE_EPSILON 1.38777878078144567552953958511352539e-17Q
657: # define PETSC_SMALL 1.e-20Q
658: #elif defined(PETSC_USE_REAL___FP16) /* maybe should use single precision values for these? */
659: # define PETSC_MAX_REAL 65504.
660: # define PETSC_MIN_REAL (-PETSC_MAX_REAL)
661: # define PETSC_MACHINE_EPSILON .00097656
662: # define PETSC_SQRT_MACHINE_EPSILON .0312
663: # define PETSC_SMALL 5.e-3
664: #endif
666: #define PETSC_INFINITY (PETSC_MAX_REAL/4)
667: #define PETSC_NINFINITY (-PETSC_INFINITY)
669: PETSC_EXTERN PetscBool PetscIsInfReal(PetscReal);
670: PETSC_EXTERN PetscBool PetscIsNanReal(PetscReal);
671: PETSC_EXTERN PetscBool PetscIsNormalReal(PetscReal);
672: PETSC_STATIC_INLINE PetscBool PetscIsInfOrNanReal(PetscReal v) {return PetscIsInfReal(v) || PetscIsNanReal(v) ? PETSC_TRUE : PETSC_FALSE;}
673: PETSC_STATIC_INLINE PetscBool PetscIsInfScalar(PetscScalar v) {return PetscIsInfReal(PetscAbsScalar(v));}
674: PETSC_STATIC_INLINE PetscBool PetscIsNanScalar(PetscScalar v) {return PetscIsNanReal(PetscAbsScalar(v));}
675: PETSC_STATIC_INLINE PetscBool PetscIsInfOrNanScalar(PetscScalar v) {return PetscIsInfOrNanReal(PetscAbsScalar(v));}
676: PETSC_STATIC_INLINE PetscBool PetscIsNormalScalar(PetscScalar v) {return PetscIsNormalReal(PetscAbsScalar(v));}
678: PETSC_EXTERN PetscBool PetscEqualReal(PetscReal,PetscReal);
679: PETSC_EXTERN PetscBool PetscEqualScalar(PetscScalar,PetscScalar);
681: /*
682: These macros are currently hardwired to match the regular data types, so there is no support for a different
683: MatScalar from PetscScalar. We left the MatScalar in the source just in case we use it again.
684: */
685: #define MPIU_MATSCALAR MPIU_SCALAR
686: typedef PetscScalar MatScalar;
687: typedef PetscReal MatReal;
689: struct petsc_mpiu_2scalar {PetscScalar a,b;};
690: PETSC_EXTERN MPI_Datatype MPIU_2SCALAR PetscAttrMPITypeTagLayoutCompatible(struct petsc_mpiu_2scalar);
691: #if defined(PETSC_USE_64BIT_INDICES) || !defined(MPI_2INT)
692: struct petsc_mpiu_2int {PetscInt a,b;};
693: PETSC_EXTERN MPI_Datatype MPIU_2INT PetscAttrMPITypeTagLayoutCompatible(struct petsc_mpiu_2int);
694: #else
695: #define MPIU_2INT MPI_2INT
696: #endif
698: PETSC_STATIC_INLINE PetscInt PetscPowInt(PetscInt base,PetscInt power)
699: {
700: PetscInt result = 1;
701: while (power) {
702: if (power & 1) result *= base;
703: power >>= 1;
704: base *= base;
705: }
706: return result;
707: }
709: PETSC_STATIC_INLINE PetscReal PetscPowRealInt(PetscReal base,PetscInt power)
710: {
711: PetscReal result = 1;
712: if (power < 0) {
713: power = -power;
714: base = ((PetscReal)1)/base;
715: }
716: while (power) {
717: if (power & 1) result *= base;
718: power >>= 1;
719: base *= base;
720: }
721: return result;
722: }
724: PETSC_STATIC_INLINE PetscScalar PetscPowScalarInt(PetscScalar base,PetscInt power)
725: {
726: PetscScalar result = 1;
727: if (power < 0) {
728: power = -power;
729: base = ((PetscReal)1)/base;
730: }
731: while (power) {
732: if (power & 1) result *= base;
733: power >>= 1;
734: base *= base;
735: }
736: return result;
737: }
739: PETSC_STATIC_INLINE PetscScalar PetscPowScalarReal(PetscScalar base,PetscReal power)
740: {
741: PetscScalar cpower = power;
742: return PetscPowScalar(base,cpower);
743: }
745: #ifndef PETSC_HAVE_LOG2
746: PETSC_STATIC_INLINE PetscReal PetscLog2Real(PetscReal n)
747: {
748: return PetscLogReal(n)/PetscLogReal(2.0);
749: }
750: #endif
751: #endif