xref: /petsc/include/petscmath.h (revision adcec1e53b1fc19aa0d23dc8f01b9147301e4f11)
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 PETSC_EXTERN MPI_Datatype MPIU_2SCALAR;
16 PETSC_EXTERN MPI_Datatype MPIU_2INT;
17 
18 /*
19 
20      Defines operations that are different for complex and real numbers;
21    note that one cannot mix the use of complex and real in the same
22    PETSc program. All PETSc objects in one program are built around the object
23    PetscScalar which is either always a real or a complex.
24 
25 */
26 
27 #define PetscExpPassiveScalar(a) PetscExpScalar()
28 #if defined(PETSC_USE_REAL_SINGLE)
29 #define MPIU_REAL   MPI_FLOAT
30 typedef float PetscReal;
31 #define PetscSqrtReal(a)    sqrt(a)
32 #define PetscExpReal(a)     exp(a)
33 #define PetscLogReal(a)     log(a)
34 #define PetscSinReal(a)     sin(a)
35 #define PetscCosReal(a)     cos(a)
36 #elif defined(PETSC_USE_REAL_DOUBLE)
37 #define MPIU_REAL   MPI_DOUBLE
38 typedef double PetscReal;
39 #define PetscSqrtReal(a)    sqrt(a)
40 #define PetscExpReal(a)     exp(a)
41 #define PetscLogReal(a)     log(a)
42 #define PetscSinReal(a)     sin(a)
43 #define PetscCosReal(a)     cos(a)
44 #elif defined(PETSC_USE_REAL___FLOAT128)
45 #if defined(__cplusplus)
46 extern "C" {
47 #endif
48 #include <quadmath.h>
49 #if defined(__cplusplus)
50 }
51 #endif
52 PETSC_EXTERN MPI_Datatype MPIU___FLOAT128;
53 #define MPIU_REAL MPIU___FLOAT128
54 typedef __float128 PetscReal;
55 #define PetscSqrtReal(a)    sqrtq(a)
56 #define PetscExpReal(a)     expq(a)
57 #define PetscLogReal(a)     logq(a)
58 #define PetscSinReal(a)     sinq(a)
59 #define PetscCosReal(a)     cosq(a)
60 #endif /* PETSC_USE_REAL_* */
61 
62 /*
63     Complex number definitions
64  */
65 #if defined(PETSC_CLANGUAGE_CXX) && defined(PETSC_HAVE_CXX_COMPLEX)
66 #if defined(PETSC_USE_COMPLEX) || defined(PETSC_DESIRE_COMPLEX)
67 #define PETSC_HAVE_COMPLEX 1
68 /* C++ support of complex number */
69 #if defined(PETSC_HAVE_CUSP)
70 #define complexlib cusp
71 #include <cusp/complex.h>
72 #else
73 #define complexlib std
74 #include <complex>
75 #endif
76 
77 #define PetscRealPartComplex(a)      (a).real()
78 #define PetscImaginaryPartComplex(a) (a).imag()
79 #define PetscAbsComplex(a)           complexlib::abs(a)
80 #define PetscConjComplex(a)          complexlib::conj(a)
81 #define PetscSqrtComplex(a)          complexlib::sqrt(a)
82 #define PetscPowComplex(a,b)         complexlib::pow(a,b)
83 #define PetscExpComplex(a)           complexlib::exp(a)
84 #define PetscLogComplex(a)           complexlib::log(a)
85 #define PetscSinComplex(a)           complexlib::sin(a)
86 #define PetscCosComplex(a)           complexlib::cos(a)
87 
88 #if defined(PETSC_USE_REAL_SINGLE)
89 typedef complexlib::complex<float> PetscComplex;
90 #elif defined(PETSC_USE_REAL_DOUBLE)
91 typedef complexlib::complex<double> PetscComplex;
92 #elif defined(PETSC_USE_REAL___FLOAT128)
93 typedef complexlib::complex<__float128> PetscComplex; /* Notstandard and not expected to work, use __complex128 */
94 #endif  /* PETSC_USE_REAL_ */
95 #endif  /* PETSC_USE_COMPLEX && PETSC_DESIRE_COMPLEX */
96 
97 #elif defined(PETSC_CLANGUAGE_C) && defined(PETSC_HAVE_C99_COMPLEX)
98 /* Use C99 _Complex for the type. Do not include complex.h by default to define "complex" because of symbol conflicts in Hypre. */
99 /* Compilation units that can safely use complex should define PETSC_DESIRE_COMPLEX before including any headers */
100 #if defined(PETSC_USE_COMPLEX) || defined(PETSC_DESIRE_COMPLEX)
101 #define PETSC_HAVE_COMPLEX 1
102 #include <complex.h>
103 
104 #if defined(PETSC_USE_REAL_SINGLE)
105 typedef float _Complex PetscComplex;
106 
107 #define PetscRealPartComplex(a)      crealf(a)
108 #define PetscImaginaryPartComplex(a) cimagf(a)
109 #define PetscAbsComplex(a)           cabsf(a)
110 #define PetscConjComplex(a)          conjf(a)
111 #define PetscSqrtComplex(a)          csqrtf(a)
112 #define PetscPowComplex(a,b)         cpowf(a,b)
113 #define PetscExpComplex(a)           cexpf(a)
114 #define PetscLogComplex(a)           clogf(a)
115 #define PetscSinComplex(a)           csinf(a)
116 #define PetscCosComplex(a)           ccosf(a)
117 
118 #elif defined(PETSC_USE_REAL_DOUBLE)
119 typedef double _Complex PetscComplex;
120 
121 #define PetscRealPartComplex(a)      creal(a)
122 #define PetscImaginaryPartComplex(a) cimag(a)
123 #define PetscAbsComplex(a)           cabs(a)
124 #define PetscConjComplex(a)          conj(a)
125 #define PetscSqrtComplex(a)          csqrt(a)
126 #define PetscPowComplex(a,b)         cpow(a,b)
127 #define PetscExpComplex(a)           cexp(a)
128 #define PetscLogComplex(a)           clog(a)
129 #define PetscSinComplex(a)           csin(a)
130 #define PetscCosComplex(a)           ccos(a)
131 
132 #elif defined(PETSC_USE_REAL___FLOAT128)
133 typedef __complex128 PetscComplex;
134 PETSC_EXTERN MPI_Datatype MPIU___COMPLEX128;
135 
136 #define PetscRealPartComplex(a)      crealq(a)
137 #define PetscImaginaryPartComplex(a) cimagq(a)
138 #define PetscAbsComplex(a)           cabsq(a)
139 #define PetscConjComplex(a)          conjq(a)
140 #define PetscSqrtComplex(a)          csqrtq(a)
141 #define PetscPowComplex(a,b)         cpowq(a,b)
142 #define PetscExpComplex(a)           cexpq(a)
143 #define PetscLogComplex(a)           clogq(a)
144 #define PetscSinComplex(a)           csinq(a)
145 #define PetscCosComplex(a)           ccosq(a)
146 #endif /* PETSC_USE_REAL_* */
147 #elif defined(PETSC_USE_COMPLEX)
148 #error "PETSc was configured --with-scalar-type=complex, but a language-appropriate complex library is not available"
149 #endif /* PETSC_USE_COMPLEX || PETSC_DESIRE_COMPLEX */
150 #endif /* (PETSC_CLANGUAGE_CXX && PETSC_HAVE_CXX_COMPLEX) else-if (PETSC_CLANGUAGE_C && PETSC_HAVE_C99_COMPLEX) */
151 
152 #if defined(PETSC_HAVE_MPI_C_DOUBLE_COMPLEX)
153 #define MPIU_C_DOUBLE_COMPLEX MPI_C_DOUBLE_COMPLEX
154 #define MPIU_C_COMPLEX MPI_C_COMPLEX
155 #else
156 PETSC_EXTERN MPI_Datatype MPIU_C_DOUBLE_COMPLEX;
157 PETSC_EXTERN MPI_Datatype MPIU_C_COMPLEX;
158 #endif /* PETSC_HAVE_MPI_C_DOUBLE_COMPLEX */
159 
160 #if defined(PETSC_HAVE_COMPLEX)
161 #  if defined(PETSC_USE_REAL_SINGLE)
162 #    define MPIU_COMPLEX MPIU_C_COMPLEX
163 #  elif defined(PETSC_USE_REAL_DOUBLE)
164 #    define MPIU_COMPLEX MPIU_C_DOUBLE_COMPLEX
165 #  elif defined(PETSC_USE_REAL___FLOAT128)
166 #    define MPIU_COMPLEX MPIU___COMPLEX128
167 #  endif /* PETSC_USE_REAL_* */
168 #endif
169 
170 #if defined(PETSC_USE_COMPLEX)
171 typedef PetscComplex PetscScalar;
172 #define PetscRealPart(a)      PetscRealPartComplex(a)
173 #define PetscImaginaryPart(a) PetscImaginaryPartComplex(a)
174 #define PetscAbsScalar(a)     PetscAbsComplex(a)
175 #define PetscConj(a)          PetscConjComplex(a)
176 #define PetscSqrtScalar(a)    PetscSqrtComplex(a)
177 #define PetscPowScalar(a,b)   PetscPowComplex(a,b)
178 #define PetscExpScalar(a)     PetscExpComplex(a)
179 #define PetscLogScalar(a)     PetscLogComplex(a)
180 #define PetscSinScalar(a)     PetscSinComplex(a)
181 #define PetscCosScalar(a)     PetscCosComplex(a)
182 
183 #define MPIU_SCALAR MPIU_COMPLEX
184 
185 /*
186     real number definitions
187  */
188 #else /* PETSC_USE_COMPLEX */
189 typedef PetscReal PetscScalar;
190 #define MPIU_SCALAR MPIU_REAL
191 
192 #define PetscRealPart(a)      (a)
193 #define PetscImaginaryPart(a) ((PetscReal)0.)
194 PETSC_STATIC_INLINE PetscReal PetscAbsScalar(PetscScalar a) {return a < 0.0 ? -a : a;}
195 #define PetscConj(a)          (a)
196 #if !defined(PETSC_USE_REAL___FLOAT128)
197 #define PetscSqrtScalar(a)    sqrt(a)
198 #define PetscPowScalar(a,b)   pow(a,b)
199 #define PetscExpScalar(a)     exp(a)
200 #define PetscLogScalar(a)     log(a)
201 #define PetscSinScalar(a)     sin(a)
202 #define PetscCosScalar(a)     cos(a)
203 #else /* PETSC_USE_REAL___FLOAT128 */
204 #define PetscSqrtScalar(a)    sqrtq(a)
205 #define PetscPowScalar(a,b)   powq(a,b)
206 #define PetscExpScalar(a)     expq(a)
207 #define PetscLogScalar(a)     logq(a)
208 #define PetscSinScalar(a)     sinq(a)
209 #define PetscCosScalar(a)     cosq(a)
210 #endif /* PETSC_USE_REAL___FLOAT128 */
211 
212 #endif /* PETSC_USE_COMPLEX */
213 
214 #define PetscSign(a) (((a) >= 0) ? ((a) == 0 ? 0 : 1) : -1)
215 #define PetscAbs(a)  (((a) >= 0) ? (a) : -(a))
216 
217 /* --------------------------------------------------------------------------*/
218 
219 /*
220    Certain objects may be created using either single or double precision.
221    This is currently not used.
222 */
223 typedef enum { PETSC_SCALAR_DOUBLE,PETSC_SCALAR_SINGLE, PETSC_SCALAR_LONG_DOUBLE } PetscScalarPrecision;
224 
225 #if defined(PETSC_HAVE_COMPLEX)
226 /* PETSC_i is the imaginary number, i */
227 PETSC_EXTERN PetscComplex PETSC_i;
228 #endif
229 
230 /*MC
231    PetscMin - Returns minimum of two numbers
232 
233    Synopsis:
234    type PetscMin(type v1,type v2)
235 
236    Not Collective
237 
238    Input Parameter:
239 +  v1 - first value to find minimum of
240 -  v2 - second value to find minimum of
241 
242 
243    Notes: type can be integer or floating point value
244 
245    Level: beginner
246 
247 
248 .seealso: PetscMin(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr()
249 
250 M*/
251 #define PetscMin(a,b)   (((a)<(b)) ?  (a) : (b))
252 
253 /*MC
254    PetscMax - Returns maxium of two numbers
255 
256    Synopsis:
257    type max PetscMax(type v1,type v2)
258 
259    Not Collective
260 
261    Input Parameter:
262 +  v1 - first value to find maximum of
263 -  v2 - second value to find maximum of
264 
265    Notes: type can be integer or floating point value
266 
267    Level: beginner
268 
269 .seealso: PetscMin(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr()
270 
271 M*/
272 #define PetscMax(a,b)   (((a)<(b)) ?  (b) : (a))
273 
274 /*MC
275    PetscClipInterval - Returns a number clipped to be within an interval
276 
277    Synopsis:
278    type clip PetscClipInterval(type x,type a,type b)
279 
280    Not Collective
281 
282    Input Parameter:
283 +  x - value to use if within interval (a,b)
284 .  a - lower end of interval
285 -  b - upper end of interval
286 
287    Notes: type can be integer or floating point value
288 
289    Level: beginner
290 
291 .seealso: PetscMin(), PetscMax(), PetscAbsInt(), PetscAbsReal(), PetscSqr()
292 
293 M*/
294 #define PetscClipInterval(x,a,b)   (PetscMax((a),PetscMin((x),(b))))
295 
296 /*MC
297    PetscAbsInt - Returns the absolute value of an integer
298 
299    Synopsis:
300    int abs PetscAbsInt(int v1)
301 
302    Not Collective
303 
304    Input Parameter:
305 .   v1 - the integer
306 
307    Level: beginner
308 
309 .seealso: PetscMax(), PetscMin(), PetscAbsReal(), PetscSqr()
310 
311 M*/
312 #define PetscAbsInt(a)  (((a)<0)   ? -(a) : (a))
313 
314 /*MC
315    PetscAbsReal - Returns the absolute value of an real number
316 
317    Synopsis:
318    Real abs PetscAbsReal(PetscReal v1)
319 
320    Not Collective
321 
322    Input Parameter:
323 .   v1 - the double
324 
325 
326    Level: beginner
327 
328 .seealso: PetscMax(), PetscMin(), PetscAbsInt(), PetscSqr()
329 
330 M*/
331 #define PetscAbsReal(a) (((a)<0)   ? -(a) : (a))
332 
333 /*MC
334    PetscSqr - Returns the square of a number
335 
336    Synopsis:
337    type sqr PetscSqr(type v1)
338 
339    Not Collective
340 
341    Input Parameter:
342 .   v1 - the value
343 
344    Notes: type can be integer or floating point value
345 
346    Level: beginner
347 
348 .seealso: PetscMax(), PetscMin(), PetscAbsInt(), PetscAbsReal()
349 
350 M*/
351 #define PetscSqr(a)     ((a)*(a))
352 
353 /* ----------------------------------------------------------------------------*/
354 /*
355      Basic constants
356 */
357 #if defined(PETSC_USE_REAL___FLOAT128)
358 #define PETSC_PI                 M_PIq
359 #elif defined(M_PI)
360 #define PETSC_PI                 M_PI
361 #else
362 #define PETSC_PI                 3.14159265358979323846264338327950288419716939937510582
363 #endif
364 
365 #if !defined(PETSC_USE_64BIT_INDICES)
366 #define PETSC_MAX_INT            2147483647
367 #define PETSC_MIN_INT            (-PETSC_MAX_INT - 1)
368 #else
369 #define PETSC_MAX_INT            9223372036854775807L
370 #define PETSC_MIN_INT            (-PETSC_MAX_INT - 1)
371 #endif
372 
373 #if defined(PETSC_USE_REAL_SINGLE)
374 #  define PETSC_MAX_REAL                3.40282346638528860e+38F
375 #  define PETSC_MIN_REAL                -PETSC_MAX_REAL
376 #  define PETSC_MACHINE_EPSILON         1.19209290e-07F
377 #  define PETSC_SQRT_MACHINE_EPSILON    3.45266983e-04F
378 #  define PETSC_SMALL                   1.e-5
379 #elif defined(PETSC_USE_REAL_DOUBLE)
380 #  define PETSC_MAX_REAL                1.7976931348623157e+308
381 #  define PETSC_MIN_REAL                -PETSC_MAX_REAL
382 #  define PETSC_MACHINE_EPSILON         2.2204460492503131e-16
383 #  define PETSC_SQRT_MACHINE_EPSILON    1.490116119384766e-08
384 #  define PETSC_SMALL                   1.e-10
385 #elif defined(PETSC_USE_REAL___FLOAT128)
386 #  define PETSC_MAX_REAL                FLT128_MAX
387 #  define PETSC_MIN_REAL                -FLT128_MAX
388 #  define PETSC_MACHINE_EPSILON         FLT128_EPSILON
389 #  define PETSC_SQRT_MACHINE_EPSILON    1.38777878078e-17
390 #  define PETSC_SMALL                   1.e-20
391 #endif
392 
393 PETSC_EXTERN PetscErrorCode PetscIsInfOrNanScalar(PetscScalar);
394 PETSC_EXTERN PetscErrorCode PetscIsInfOrNanReal(PetscReal);
395 
396 /* ----------------------------------------------------------------------------*/
397 #define PassiveReal   PetscReal
398 #define PassiveScalar PetscScalar
399 
400 /*
401     These macros are currently hardwired to match the regular data types, so there is no support for a different
402     MatScalar from PetscScalar. We left the MatScalar in the source just in case we use it again.
403  */
404 #define MPIU_MATSCALAR MPIU_SCALAR
405 typedef PetscScalar MatScalar;
406 typedef PetscReal MatReal;
407 
408 
409 #endif
410