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