xref: /petsc/include/petscmath.h (revision 193ac0bc5128976c4ec2c1a67f3f9cb026b77f22)
1 /*
2 
3       PETSc mathematics include file. Defines certain basic mathematical
4     constants and functions for working with single and double precision
5     floating point numbers as well as complex and integers.
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 PETSC_EXTERN_CXX_BEGIN
15 
16 extern  MPI_Datatype  MPIU_2SCALAR;
17 extern  MPI_Datatype  MPIU_2INT;
18 
19 /*
20 
21      Defines operations that are different for complex and real numbers;
22    note that one cannot really mix the use of complex and real in the same
23    PETSc program. All PETSc objects in one program are built around the object
24    PetscScalar which is either always a real or a complex.
25 
26 */
27 
28 #define PetscExpPassiveScalar(a) PetscExpScalar()
29 #if defined(PETSC_USE_REAL_SINGLE)
30 #define MPIU_REAL   MPI_FLOAT
31 typedef float PetscReal;
32 #elif defined(PETSC_USE_REAL_DOUBLE)
33 #define MPIU_REAL   MPI_DOUBLE
34 typedef double PetscReal;
35 #elif defined(PETSC_USE_REAL_LONG_DOUBLE)
36 #define MPIU_REAL   MPI_LONG_DOUBLE
37 typedef long double PetscReal;
38 #elif defined(PETSC_USE_REAL___FLOAT128)
39 #define MPIU_REAL MPIU___FLOAT128
40 typedef __float128 PetscReal;
41 #endif /* PETSC_USE_REAL_* */
42 
43 /*
44     Complex number definitions
45  */
46 #if defined(PETSC_USE_COMPLEX)
47 #if defined(PETSC_CLANGUAGE_CXX)
48 /* C++ support of complex number */
49 #include <complex>
50 
51 #define PetscRealPart(a)      (a).real()
52 #define PetscImaginaryPart(a) (a).imag()
53 #define PetscAbsScalar(a)     std::abs(a)
54 #define PetscConj(a)          std::conj(a)
55 #define PetscSqrtScalar(a)    std::sqrt(a)
56 #define PetscPowScalar(a,b)   std::pow(a,b)
57 #define PetscExpScalar(a)     std::exp(a)
58 #define PetscLogScalar(a)     std::log(a)
59 #define PetscSinScalar(a)     std::sin(a)
60 #define PetscCosScalar(a)     std::cos(a)
61 
62 #if defined(PETSC_USE_REAL_SINGLE)
63 typedef std::complex<float> PetscScalar;
64 #elif defined(PETSC_USE_REAL_DOUBLE)
65 typedef std::complex<double> PetscScalar;
66 #elif defined(PETSC_USE_REAL_LONG_DOUBLE)
67 typedef std::complex<long double> PetscScalar;
68 #endif /* PETSC_USE_REAL_* */
69 
70 #else /* PETSC_CLANGUAGE_CXX */
71 /*  C support of complex numbers: Requires C99 compliant compiler*/
72 #include <complex.h>
73 
74 #if defined(PETSC_USE_REAL_SINGLE)
75 typedef float complex PetscScalar;
76 
77 #define PetscRealPart(a)      crealf(a)
78 #define PetscImaginaryPart(a) cimagf(a)
79 #define PetscAbsScalar(a)     cabsf(a)
80 #define PetscConj(a)          conjf(a)
81 #define PetscSqrtScalar(a)    csqrtf(a)
82 #define PetscPowScalar(a,b)   cpowf(a,b)
83 #define PetscExpScalar(a)     cexpf(a)
84 #define PetscLogScalar(a)     clogf(a)
85 #define PetscSinScalar(a)     csinf(a)
86 #define PetscCosScalar(a)     ccosf(a)
87 
88 #elif defined(PETSC_USE_REAL_DOUBLE)
89 typedef double complex PetscScalar;
90 
91 #define PetscRealPart(a)      creal(a)
92 #define PetscImaginaryPart(a) cimag(a)
93 #define PetscAbsScalar(a)     cabs(a)
94 #define PetscConj(a)          conj(a)
95 #define PetscSqrtScalar(a)    csqrt(a)
96 #define PetscPowScalar(a,b)   cpow(a,b)
97 #define PetscExpScalar(a)     cexp(a)
98 #define PetscLogScalar(a)     clog(a)
99 #define PetscSinScalar(a)     csin(a)
100 #define PetscCosScalar(a)     ccos(a)
101 
102 #elif defined(PETSC_USE_REAL_LONG_DOUBLE)
103 typedef long double complex PetscScalar;
104 
105 #define PetscRealPart(a)      creall(a)
106 #define PetscImaginaryPart(a) cimagl(a)
107 #define PetscAbsScalar(a)     cabsl(a)
108 #define PetscConj(a)          conjl(a)
109 #define PetscSqrtScalar(a)    csqrtl(a)
110 #define PetscPowScalar(a,b)   cpowl(a,b)
111 #define PetscExpScalar(a)     cexpl(a)
112 #define PetscLogScalar(a)     clogl(a)
113 #define PetscSinScalar(a)     csinl(a)
114 #define PetscCosScalar(a)     ccosl(a)
115 
116 #endif /* PETSC_USE_REAL_* */
117 #endif /* PETSC_CLANGUAGE_CXX */
118 
119 #if !defined(PETSC_HAVE_MPI_C_DOUBLE_COMPLEX)
120 extern  MPI_Datatype  MPI_C_DOUBLE_COMPLEX;
121 extern  MPI_Datatype  MPI_C_COMPLEX;
122 #endif /* PETSC_HAVE_MPI_C_DOUBLE_COMPLEX */
123 
124 #if defined(PETSC_USE_REAL_SINGLE)
125 #define MPIU_SCALAR MPI_C_COMPLEX
126 #elif defined(PETSC_USE_REAL_DOUBLE)
127 #define MPIU_SCALAR MPI_C_DOUBLE_COMPLEX
128 #elif defined(PETSC_USE_REAL_LONG_DOUBLE)
129 #define MPIU_SCALAR error
130 #endif /* PETSC_USE_REAL_* */
131 
132 /*
133     real number definitions
134  */
135 #else /* PETSC_USE_COMPLEX */
136 #if defined(PETSC_USE_REAL_SINGLE)
137 #define MPIU_SCALAR           MPI_FLOAT
138 typedef float PetscScalar;
139 #elif defined(PETSC_USE_REAL_DOUBLE)
140 #define MPIU_SCALAR           MPI_DOUBLE
141 typedef double PetscScalar;
142 #elif defined(PETSC_USE_REAL_LONG_DOUBLE)
143 #define MPIU_SCALAR           MPI_LONG_DOUBLE
144 typedef long double PetscScalar;
145 #elif defined(PETSC_USE_REAL___FLOAT128)
146 extern MPI_Datatype MPIU___FLOAT128;
147 #define MPIU_SCALAR MPIU___FLOAT128
148 typedef __float128 PetscScalar;
149 #endif /* PETSC_USE_REAL_* */
150 #define PetscRealPart(a)      (a)
151 #define PetscImaginaryPart(a) ((PetscReal)0.)
152 PETSC_STATIC_INLINE PetscReal PetscAbsScalar(PetscScalar a) {return a < 0.0 ? -a : a;}
153 #define PetscConj(a)          (a)
154 #if !defined(PETSC_USE_REAL___FLOAT128)
155 #define PetscSqrtScalar(a)    sqrt(a)
156 #define PetscPowScalar(a,b)   pow(a,b)
157 #define PetscExpScalar(a)     exp(a)
158 #define PetscLogScalar(a)     log(a)
159 #define PetscSinScalar(a)     sin(a)
160 #define PetscCosScalar(a)     cos(a)
161 #else /* PETSC_USE_REAL___FLOAT128 */
162 #include <quadmath.h>
163 #define PetscSqrtScalar(a)    sqrtq(a)
164 #define PetscPowScalar(a,b)   powq(a,b)
165 #define PetscExpScalar(a)     expq(a)
166 #define PetscLogScalar(a)     logq(a)
167 #define PetscSinScalar(a)     sinq(a)
168 #define PetscCosScalar(a)     cosq(a)
169 #endif /* PETSC_USE_REAL___FLOAT128 */
170 
171 #endif /* PETSC_USE_COMPLEX */
172 
173 #define PetscSign(a) (((a) >= 0) ? ((a) == 0 ? 0 : 1) : -1)
174 #define PetscAbs(a)  (((a) >= 0) ? (a) : -(a))
175 
176 /* --------------------------------------------------------------------------*/
177 
178 /*
179    Certain objects may be created using either single or double precision.
180    This is currently not used.
181 */
182 typedef enum { PETSC_SCALAR_DOUBLE,PETSC_SCALAR_SINGLE, PETSC_SCALAR_LONG_DOUBLE } PetscScalarPrecision;
183 
184 /* PETSC_i is the imaginary number, i */
185 extern  PetscScalar  PETSC_i;
186 
187 /*MC
188    PetscMin - Returns minimum of two numbers
189 
190    Synopsis:
191    type PetscMin(type v1,type v2)
192 
193    Not Collective
194 
195    Input Parameter:
196 +  v1 - first value to find minimum of
197 -  v2 - second value to find minimum of
198 
199 
200    Notes: type can be integer or floating point value
201 
202    Level: beginner
203 
204 
205 .seealso: PetscMin(), PetscAbsInt(), PetscAbsReal(), PetscSqr()
206 
207 M*/
208 #define PetscMin(a,b)   (((a)<(b)) ?  (a) : (b))
209 
210 /*MC
211    PetscMax - Returns maxium of two numbers
212 
213    Synopsis:
214    type max PetscMax(type v1,type v2)
215 
216    Not Collective
217 
218    Input Parameter:
219 +  v1 - first value to find maximum of
220 -  v2 - second value to find maximum of
221 
222    Notes: type can be integer or floating point value
223 
224    Level: beginner
225 
226 .seealso: PetscMin(), PetscAbsInt(), PetscAbsReal(), PetscSqr()
227 
228 M*/
229 #define PetscMax(a,b)   (((a)<(b)) ?  (b) : (a))
230 
231 /*MC
232    PetscAbsInt - Returns the absolute value of an integer
233 
234    Synopsis:
235    int abs PetscAbsInt(int v1)
236 
237    Not Collective
238 
239    Input Parameter:
240 .   v1 - the integer
241 
242    Level: beginner
243 
244 .seealso: PetscMax(), PetscMin(), PetscAbsReal(), PetscSqr()
245 
246 M*/
247 #define PetscAbsInt(a)  (((a)<0)   ? -(a) : (a))
248 
249 /*MC
250    PetscAbsReal - Returns the absolute value of an real number
251 
252    Synopsis:
253    Real abs PetscAbsReal(PetscReal v1)
254 
255    Not Collective
256 
257    Input Parameter:
258 .   v1 - the double
259 
260 
261    Level: beginner
262 
263 .seealso: PetscMax(), PetscMin(), PetscAbsInt(), PetscSqr()
264 
265 M*/
266 #define PetscAbsReal(a) (((a)<0)   ? -(a) : (a))
267 
268 /*MC
269    PetscSqr - Returns the square of a number
270 
271    Synopsis:
272    type sqr PetscSqr(type v1)
273 
274    Not Collective
275 
276    Input Parameter:
277 .   v1 - the value
278 
279    Notes: type can be integer or floating point value
280 
281    Level: beginner
282 
283 .seealso: PetscMax(), PetscMin(), PetscAbsInt(), PetscAbsReal()
284 
285 M*/
286 #define PetscSqr(a)     ((a)*(a))
287 
288 /* ----------------------------------------------------------------------------*/
289 /*
290      Basic constants
291 */
292 #if defined(PETSC_USE_REAL___FLOAT128)
293 #define PETSC_PI                 M_PIq
294 #elif defined(M_PI)
295 #define PETSC_PI                 M_PI
296 #else
297 #define PETSC_PI                 3.14159265358979323846264
298 #endif
299 
300 
301 #define PETSC_MAX_INT            2147483647
302 #define PETSC_MIN_INT            -2147483647
303 
304 #if defined(PETSC_USE_REAL_SINGLE)
305 #if defined(MAXFLOAT)
306 #  define PETSC_MAX_REAL                 MAXFLOAT
307 #else
308 #  define PETSC_MAX_REAL                1.e30
309 #endif
310 #  define PETSC_MIN_REAL                -PETSC_MAX_REAL
311 #  define PETSC_MACHINE_EPSILON         1.e-7
312 #  define PETSC_SQRT_MACHINE_EPSILON    3.e-4
313 #  define PETSC_SMALL                   1.e-5
314 #elif defined(PETSC_USE_REAL_DOUBLE)
315 #  define PETSC_MAX_REAL                1.e300
316 #  define PETSC_MIN_REAL                -PETSC_MAX_REAL
317 #  define PETSC_MACHINE_EPSILON         1.e-14
318 #  define PETSC_SQRT_MACHINE_EPSILON    1.e-7
319 #  define PETSC_SMALL                   1.e-10
320 #elif defined(PETSC_USE_REAL_LONG_DOUBLE)
321 #  define PETSC_MAX_REAL                1.e4900L
322 #  define PETSC_MIN_REAL                -PETSC_MAX_REAL
323 #  define PETSC_MACHINE_EPSILON         1.e-18
324 #  define PETSC_SQRT_MACHINE_EPSILON    1.e-9
325 #  define PETSC_SMALL                   1.e-13
326 #elif defined(PETSC_USE_REAL___FLOAT128)
327 #  define PETSC_MAX_REAL                FLT128_MAX
328 #  define PETSC_MIN_REAL                -FLT128_MAX
329 #  define PETSC_MACHINE_EPSILON         FLT128_EPSILON
330 #  define PETSC_SQRT_MACHINE_EPSILON    1.38777878078e-17
331 #  define PETSC_SMALL                   1.e-20
332 #endif
333 
334 #if defined PETSC_HAVE_ADIC
335 /* Use MPI_Allreduce when ADIC is not available. */
336 extern PetscErrorCode  PetscGlobalMax(MPI_Comm, const PetscReal*,PetscReal*);
337 extern PetscErrorCode  PetscGlobalMin(MPI_Comm, const PetscReal*,PetscReal*);
338 extern PetscErrorCode  PetscGlobalSum(MPI_Comm, const PetscScalar*,PetscScalar*);
339 #endif
340 
341 /*MC
342       PetscIsInfOrNan - Returns 1 if the input double has an infinity for Not-a-number (Nan) value, otherwise 0.
343 
344     Input Parameter:
345 .     a - the double
346 
347 
348      Notes: uses the C99 standard isinf() and isnan() on systems where they exist.
349       Otherwises uses ( (a - a) != 0.0), note that some optimizing compiles compile
350       out this form, thus removing the check.
351 
352      Level: beginner
353 
354 M*/
355 #if defined(PETSC_USE_REAL___FLOAT128)
356 PETSC_STATIC_INLINE PetscErrorCode PetscIsInfOrNanScalar(PetscScalar a) {
357   return isinfq(PetscAbsScalar(a)) || isnanq(PetscAbsScalar(a));
358 }
359 PETSC_STATIC_INLINE PetscErrorCode PetscIsInfOrNanReal(PetscReal a) {
360   return isinfq(a) || isnanq(a);
361 }
362 #elif defined(PETSC_HAVE_ISINF) && defined(PETSC_HAVE_ISNAN) && !defined(_GLIBCXX_CMATH)
363 PETSC_STATIC_INLINE PetscErrorCode PetscIsInfOrNanScalar(PetscScalar a) {
364   return isinf(PetscAbsScalar(a)) || isnan(PetscAbsScalar(a));
365 }
366 PETSC_STATIC_INLINE PetscErrorCode PetscIsInfOrNanReal(PetscReal a) {
367   return isinf(a) || isnan(a);
368 }
369 #elif defined(PETSC_HAVE__FINITE) && defined(PETSC_HAVE__ISNAN)
370 #if defined(PETSC_HAVE_FLOAT_H)
371 #include "float.h"  /* Microsoft Windows defines _finite() in float.h */
372 #endif
373 #if defined(PETSC_HAVE_IEEEFP_H)
374 #include "ieeefp.h"  /* Solaris prototypes these here */
375 #endif
376 PETSC_STATIC_INLINE PetscErrorCode PetscIsInfOrNanScalar(PetscScalar a) {
377   return !_finite(PetscAbsScalar(a)) || _isnan(PetscAbsScalar(a));
378 }
379 PETSC_STATIC_INLINE PetscErrorCode PetscIsInfOrNanReal(PetscReal a) {
380   return !_finite(a) || _isnan(a);
381 }
382 #else
383 PETSC_STATIC_INLINE PetscErrorCode PetscIsInfOrNanScalar(PetscScalar a) {
384   return  ((a - a) != (PetscScalar)0);
385 }
386 PETSC_STATIC_INLINE PetscErrorCode PetscIsInfOrNanReal(PetscReal a) {
387   return ((a - a) != 0);
388 }
389 #endif
390 
391 
392 /* ----------------------------------------------------------------------------*/
393 /*
394     PetscLogDouble variables are used to contain double precision numbers
395   that are not used in the numerical computations, but rather in logging,
396   timing etc.
397 */
398 typedef double PetscLogDouble;
399 #define MPIU_PETSCLOGDOUBLE MPI_DOUBLE
400 
401 #define PassiveReal   PetscReal
402 #define PassiveScalar PetscScalar
403 
404 /*
405     These macros are currently hardwired to match the regular data types, so there is no support for a different
406     MatScalar from PetscScalar. We left the MatScalar in the source just in case we use it again.
407  */
408 #define MPIU_MATSCALAR MPIU_SCALAR
409 typedef PetscScalar MatScalar;
410 typedef PetscReal MatReal;
411 
412 
413 PETSC_EXTERN_CXX_END
414 #endif
415