xref: /petsc/include/petscmath.h (revision 9ae82921df069a58776bfe4da82b38e8ff7dd41c)
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 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 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 #define PetscSqrtReal(a)    sqrt(a)
33 #elif defined(PETSC_USE_REAL_DOUBLE)
34 #define MPIU_REAL   MPI_DOUBLE
35 typedef double PetscReal;
36 #define PetscSqrtReal(a)    sqrt(a)
37 #elif defined(PETSC_USE_REAL___FLOAT128)
38 #define MPIU_REAL MPIU___FLOAT128
39 typedef __float128 PetscReal;
40 #define PetscSqrtReal(a)    sqrtq(a)
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 /*
50 #include <complex>
51 
52 #define PetscRealPart(a)      (a).real()
53 #define PetscImaginaryPart(a) (a).imag()
54 #define PetscAbsScalar(a)     std::abs(a)
55 #define PetscConj(a)          std::conj(a)
56 #define PetscSqrtScalar(a)    std::sqrt(a)
57 #define PetscPowScalar(a,b)   std::pow(a,b)
58 #define PetscExpScalar(a)     std::exp(a)
59 #define PetscLogScalar(a)     std::log(a)
60 #define PetscSinScalar(a)     std::sin(a)
61 #define PetscCosScalar(a)     std::cos(a)
62 
63 #if defined(PETSC_USE_REAL_SINGLE)
64 typedef std::complex<float> PetscScalar;
65 #elif defined(PETSC_USE_REAL_DOUBLE)
66 typedef std::complex<double> PetscScalar;
67 #endif */ /* PETSC_USE_REAL_* */
68 
69 /* CUSP support of complex number */
70 
71 #include <cusp/complex.h>
72 
73 #define PetscRealPart(a)      (a).real()
74 #define PetscImaginaryPart(a) (a).imag()
75 #define PetscAbsScalar(a)     cusp::abs(a)
76 #define PetscConj(a)          cusp::conj(a)
77 #define PetscSqrtScalar(a)    cusp::sqrt(a)
78 #define PetscPowScalar(a,b)   cusp::pow(a,b)
79 #define PetscExpScalar(a)     cusp::exp(a)
80 #define PetscLogScalar(a)     cusp::log(a)
81 #define PetscSinScalar(a)     cusp::sin(a)
82 #define PetscCosScalar(a)     cusp::cos(a)
83 
84 #if defined(PETSC_USE_REAL_SINGLE)
85 typedef cusp::complex<float> PetscScalar;
86 #elif defined(PETSC_USE_REAL_DOUBLE)
87 typedef cusp::complex<double> PetscScalar;
88 #endif
89 
90 #else /* PETSC_CLANGUAGE_CXX */
91 /*  C support of complex numbers: Requires C99 compliant compiler*/
92 #include <cusp::complex.h>
93 
94 #if defined(PETSC_USE_REAL_SINGLE)
95 typedef float complex PetscScalar;
96 
97 #define PetscRealPart(a)      crealf(a)
98 #define PetscImaginaryPart(a) cimagf(a)
99 #define PetscAbsScalar(a)     cabsf(a)
100 #define PetscConj(a)          conjf(a)
101 #define PetscSqrtScalar(a)    csqrtf(a)
102 #define PetscPowScalar(a,b)   cpowf(a,b)
103 #define PetscExpScalar(a)     cexpf(a)
104 #define PetscLogScalar(a)     clogf(a)
105 #define PetscSinScalar(a)     csinf(a)
106 #define PetscCosScalar(a)     ccosf(a)
107 
108 #elif defined(PETSC_USE_REAL_DOUBLE)
109 typedef double complex PetscScalar;
110 
111 #define PetscRealPart(a)      creal(a)
112 #define PetscImaginaryPart(a) cimag(a)
113 #define PetscAbsScalar(a)     cabs(a)
114 #define PetscConj(a)          conj(a)
115 #define PetscSqrtScalar(a)    csqrt(a)
116 #define PetscPowScalar(a,b)   cpow(a,b)
117 #define PetscExpScalar(a)     cexp(a)
118 #define PetscLogScalar(a)     clog(a)
119 #define PetscSinScalar(a)     csin(a)
120 #define PetscCosScalar(a)     ccos(a)
121 
122 #endif /* PETSC_USE_REAL_* */
123 #endif /* PETSC_CLANGUAGE_CXX */
124 
125 #if defined(PETSC_HAVE_MPI_C_DOUBLE_COMPLEX)
126 #define MPIU_C_DOUBLE_COMPLEX MPI_C_DOUBLE_COMPLEX
127 #define MPIU_C_COMPLEX MPI_C_COMPLEX
128 #else
129 extern MPI_Datatype  MPIU_C_DOUBLE_COMPLEX;
130 extern MPI_Datatype  MPIU_C_COMPLEX;
131 #endif /* PETSC_HAVE_MPI_C_DOUBLE_COMPLEX */
132 
133 #if defined(PETSC_USE_REAL_SINGLE)
134 #define MPIU_SCALAR MPIU_C_COMPLEX
135 #elif defined(PETSC_USE_REAL_DOUBLE)
136 #define MPIU_SCALAR MPIU_C_DOUBLE_COMPLEX
137 #endif /* PETSC_USE_REAL_* */
138 
139 /*
140     real number definitions
141  */
142 #else /* PETSC_USE_COMPLEX */
143 #if defined(PETSC_USE_REAL_SINGLE)
144 #define MPIU_SCALAR           MPI_FLOAT
145 typedef float PetscScalar;
146 #elif defined(PETSC_USE_REAL_DOUBLE)
147 #define MPIU_SCALAR           MPI_DOUBLE
148 typedef double PetscScalar;
149 #elif defined(PETSC_USE_REAL___FLOAT128)
150 extern MPI_Datatype MPIU___FLOAT128;
151 #define MPIU_SCALAR MPIU___FLOAT128
152 typedef __float128 PetscScalar;
153 #endif /* PETSC_USE_REAL_* */
154 #define PetscRealPart(a)      (a)
155 #define PetscImaginaryPart(a) ((PetscReal)0.)
156 PETSC_STATIC_INLINE PetscReal PetscAbsScalar(PetscScalar a) {return a < 0.0 ? -a : a;}
157 #define PetscConj(a)          (a)
158 #if !defined(PETSC_USE_REAL___FLOAT128)
159 #define PetscSqrtScalar(a)    sqrt(a)
160 #define PetscPowScalar(a,b)   pow(a,b)
161 #define PetscExpScalar(a)     exp(a)
162 #define PetscLogScalar(a)     log(a)
163 #define PetscSinScalar(a)     sin(a)
164 #define PetscCosScalar(a)     cos(a)
165 #else /* PETSC_USE_REAL___FLOAT128 */
166 #include <quadmath.h>
167 #define PetscSqrtScalar(a)    sqrtq(a)
168 #define PetscPowScalar(a,b)   powq(a,b)
169 #define PetscExpScalar(a)     expq(a)
170 #define PetscLogScalar(a)     logq(a)
171 #define PetscSinScalar(a)     sinq(a)
172 #define PetscCosScalar(a)     cosq(a)
173 #endif /* PETSC_USE_REAL___FLOAT128 */
174 
175 #endif /* PETSC_USE_COMPLEX */
176 
177 #define PetscSign(a) (((a) >= 0) ? ((a) == 0 ? 0 : 1) : -1)
178 #define PetscAbs(a)  (((a) >= 0) ? (a) : -(a))
179 
180 /* --------------------------------------------------------------------------*/
181 
182 /*
183    Certain objects may be created using either single or double precision.
184    This is currently not used.
185 */
186 typedef enum { PETSC_SCALAR_DOUBLE,PETSC_SCALAR_SINGLE, PETSC_SCALAR_LONG_DOUBLE } PetscScalarPrecision;
187 
188 /* PETSC_i is the imaginary number, i */
189 extern  PetscScalar  PETSC_i;
190 
191 /*MC
192    PetscMin - Returns minimum of two numbers
193 
194    Synopsis:
195    type PetscMin(type v1,type v2)
196 
197    Not Collective
198 
199    Input Parameter:
200 +  v1 - first value to find minimum of
201 -  v2 - second value to find minimum of
202 
203 
204    Notes: type can be integer or floating point value
205 
206    Level: beginner
207 
208 
209 .seealso: PetscMin(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr()
210 
211 M*/
212 #define PetscMin(a,b)   (((a)<(b)) ?  (a) : (b))
213 
214 /*MC
215    PetscMax - Returns maxium of two numbers
216 
217    Synopsis:
218    type max PetscMax(type v1,type v2)
219 
220    Not Collective
221 
222    Input Parameter:
223 +  v1 - first value to find maximum of
224 -  v2 - second value to find maximum of
225 
226    Notes: type can be integer or floating point value
227 
228    Level: beginner
229 
230 .seealso: PetscMin(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr()
231 
232 M*/
233 #define PetscMax(a,b)   (((a)<(b)) ?  (b) : (a))
234 
235 /*MC
236    PetscClipInterval - Returns a number clipped to be within an interval
237 
238    Synopsis:
239    type clip PetscClipInterval(type x,type a,type b)
240 
241    Not Collective
242 
243    Input Parameter:
244 +  x - value to use if within interval (a,b)
245 .  a - lower end of interval
246 -  b - upper end of interval
247 
248    Notes: type can be integer or floating point value
249 
250    Level: beginner
251 
252 .seealso: PetscMin(), PetscMax(), PetscAbsInt(), PetscAbsReal(), PetscSqr()
253 
254 M*/
255 #define PetscClipInterval(x,a,b)   (PetscMax((a),PetscMin((x),(b))))
256 
257 /*MC
258    PetscAbsInt - Returns the absolute value of an integer
259 
260    Synopsis:
261    int abs PetscAbsInt(int v1)
262 
263    Not Collective
264 
265    Input Parameter:
266 .   v1 - the integer
267 
268    Level: beginner
269 
270 .seealso: PetscMax(), PetscMin(), PetscAbsReal(), PetscSqr()
271 
272 M*/
273 #define PetscAbsInt(a)  (((a)<0)   ? -(a) : (a))
274 
275 /*MC
276    PetscAbsReal - Returns the absolute value of an real number
277 
278    Synopsis:
279    Real abs PetscAbsReal(PetscReal v1)
280 
281    Not Collective
282 
283    Input Parameter:
284 .   v1 - the double
285 
286 
287    Level: beginner
288 
289 .seealso: PetscMax(), PetscMin(), PetscAbsInt(), PetscSqr()
290 
291 M*/
292 #define PetscAbsReal(a) (((a)<0)   ? -(a) : (a))
293 
294 /*MC
295    PetscSqr - Returns the square of a number
296 
297    Synopsis:
298    type sqr PetscSqr(type v1)
299 
300    Not Collective
301 
302    Input Parameter:
303 .   v1 - the value
304 
305    Notes: type can be integer or floating point value
306 
307    Level: beginner
308 
309 .seealso: PetscMax(), PetscMin(), PetscAbsInt(), PetscAbsReal()
310 
311 M*/
312 #define PetscSqr(a)     ((a)*(a))
313 
314 /* ----------------------------------------------------------------------------*/
315 /*
316      Basic constants
317 */
318 #if defined(PETSC_USE_REAL___FLOAT128)
319 #define PETSC_PI                 M_PIq
320 #elif defined(M_PI)
321 #define PETSC_PI                 M_PI
322 #else
323 #define PETSC_PI                 3.14159265358979323846264338327950288419716939937510582
324 #endif
325 
326 #if !defined(PETSC_USE_64BIT_INDICES)
327 #define PETSC_MAX_INT            2147483647
328 #define PETSC_MIN_INT            (-PETSC_MAX_INT - 1)
329 #else
330 #define PETSC_MAX_INT            9223372036854775807L
331 #define PETSC_MIN_INT            (-PETSC_MAX_INT - 1)
332 #endif
333 
334 #if defined(PETSC_USE_REAL_SINGLE)
335 #  define PETSC_MAX_REAL                3.40282346638528860e+38F
336 #  define PETSC_MIN_REAL                -PETSC_MAX_REAL
337 #  define PETSC_MACHINE_EPSILON         1.19209290e-07F
338 #  define PETSC_SQRT_MACHINE_EPSILON    3.45266983e-04F
339 #  define PETSC_SMALL                   1.e-5
340 #elif defined(PETSC_USE_REAL_DOUBLE)
341 #  define PETSC_MAX_REAL                1.7976931348623157e+308
342 #  define PETSC_MIN_REAL                -PETSC_MAX_REAL
343 #  define PETSC_MACHINE_EPSILON         2.2204460492503131e-16
344 #  define PETSC_SQRT_MACHINE_EPSILON    1.490116119384766e-08
345 #  define PETSC_SMALL                   1.e-10
346 #elif defined(PETSC_USE_REAL___FLOAT128)
347 #  define PETSC_MAX_REAL                FLT128_MAX
348 #  define PETSC_MIN_REAL                -FLT128_MAX
349 #  define PETSC_MACHINE_EPSILON         FLT128_EPSILON
350 #  define PETSC_SQRT_MACHINE_EPSILON    1.38777878078e-17
351 #  define PETSC_SMALL                   1.e-20
352 #endif
353 
354 #if defined PETSC_HAVE_ADIC
355 /* Use MPI_Allreduce when ADIC is not available. */
356 extern PetscErrorCode  PetscGlobalMax(MPI_Comm, const PetscReal*,PetscReal*);
357 extern PetscErrorCode  PetscGlobalMin(MPI_Comm, const PetscReal*,PetscReal*);
358 extern PetscErrorCode  PetscGlobalSum(MPI_Comm, const PetscScalar*,PetscScalar*);
359 #endif
360 
361 extern PetscErrorCode PetscIsInfOrNanScalar(PetscScalar);
362 extern PetscErrorCode PetscIsInfOrNanReal(PetscReal);
363 
364 /* ----------------------------------------------------------------------------*/
365 #define PassiveReal   PetscReal
366 #define PassiveScalar PetscScalar
367 
368 /*
369     These macros are currently hardwired to match the regular data types, so there is no support for a different
370     MatScalar from PetscScalar. We left the MatScalar in the source just in case we use it again.
371  */
372 #define MPIU_MATSCALAR MPIU_SCALAR
373 typedef PetscScalar MatScalar;
374 typedef PetscReal MatReal;
375 
376 
377 PETSC_EXTERN_CXX_END
378 #endif
379