xref: /petsc/include/petscmath.h (revision 9371c9d470a9602b6d10a8bf50c9b2280a79e45a)
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 #if !defined(PETSCMATH_H)
11 #define PETSCMATH_H
12 
13 #include <math.h>
14 #include <petscsystypes.h>
15 
16 /* SUBMANSEC = Sys */
17 
18 /*
19 
20    Defines operations that are different for complex and real numbers.
21    All PETSc objects in one program are built around the object
22    PetscScalar which is either always a real or a complex.
23 
24 */
25 
26 /*
27     Real number definitions
28  */
29 #if defined(PETSC_USE_REAL_SINGLE)
30 #define PetscSqrtReal(a)        sqrtf(a)
31 #define PetscCbrtReal(a)        cbrtf(a)
32 #define PetscHypotReal(a, b)    hypotf(a, b)
33 #define PetscAtan2Real(a, b)    atan2f(a, b)
34 #define PetscPowReal(a, b)      powf(a, b)
35 #define PetscExpReal(a)         expf(a)
36 #define PetscLogReal(a)         logf(a)
37 #define PetscLog10Real(a)       log10f(a)
38 #define PetscLog2Real(a)        log2f(a)
39 #define PetscSinReal(a)         sinf(a)
40 #define PetscCosReal(a)         cosf(a)
41 #define PetscTanReal(a)         tanf(a)
42 #define PetscAsinReal(a)        asinf(a)
43 #define PetscAcosReal(a)        acosf(a)
44 #define PetscAtanReal(a)        atanf(a)
45 #define PetscSinhReal(a)        sinhf(a)
46 #define PetscCoshReal(a)        coshf(a)
47 #define PetscTanhReal(a)        tanhf(a)
48 #define PetscAsinhReal(a)       asinhf(a)
49 #define PetscAcoshReal(a)       acoshf(a)
50 #define PetscAtanhReal(a)       atanhf(a)
51 #define PetscErfReal(a)         erff(a)
52 #define PetscCeilReal(a)        ceilf(a)
53 #define PetscFloorReal(a)       floorf(a)
54 #define PetscFmodReal(a, b)     fmodf(a, b)
55 #define PetscCopysignReal(a, b) copysignf(a, b)
56 #define PetscTGamma(a)          tgammaf(a)
57 #if defined(PETSC_HAVE_LGAMMA_IS_GAMMA)
58 #define PetscLGamma(a) gammaf(a)
59 #else
60 #define PetscLGamma(a) lgammaf(a)
61 #endif
62 
63 #elif defined(PETSC_USE_REAL_DOUBLE)
64 #define PetscSqrtReal(a)        sqrt(a)
65 #define PetscCbrtReal(a)        cbrt(a)
66 #define PetscHypotReal(a, b)    hypot(a, b)
67 #define PetscAtan2Real(a, b)    atan2(a, b)
68 #define PetscPowReal(a, b)      pow(a, b)
69 #define PetscExpReal(a)         exp(a)
70 #define PetscLogReal(a)         log(a)
71 #define PetscLog10Real(a)       log10(a)
72 #define PetscLog2Real(a)        log2(a)
73 #define PetscSinReal(a)         sin(a)
74 #define PetscCosReal(a)         cos(a)
75 #define PetscTanReal(a)         tan(a)
76 #define PetscAsinReal(a)        asin(a)
77 #define PetscAcosReal(a)        acos(a)
78 #define PetscAtanReal(a)        atan(a)
79 #define PetscSinhReal(a)        sinh(a)
80 #define PetscCoshReal(a)        cosh(a)
81 #define PetscTanhReal(a)        tanh(a)
82 #define PetscAsinhReal(a)       asinh(a)
83 #define PetscAcoshReal(a)       acosh(a)
84 #define PetscAtanhReal(a)       atanh(a)
85 #define PetscErfReal(a)         erf(a)
86 #define PetscCeilReal(a)        ceil(a)
87 #define PetscFloorReal(a)       floor(a)
88 #define PetscFmodReal(a, b)     fmod(a, b)
89 #define PetscCopysignReal(a, b) copysign(a, b)
90 #define PetscTGamma(a)          tgamma(a)
91 #if defined(PETSC_HAVE_LGAMMA_IS_GAMMA)
92 #define PetscLGamma(a) gamma(a)
93 #else
94 #define PetscLGamma(a) lgamma(a)
95 #endif
96 
97 #elif defined(PETSC_USE_REAL___FLOAT128)
98 #define PetscSqrtReal(a)        sqrtq(a)
99 #define PetscCbrtReal(a)        cbrtq(a)
100 #define PetscHypotReal(a, b)    hypotq(a, b)
101 #define PetscAtan2Real(a, b)    atan2q(a, b)
102 #define PetscPowReal(a, b)      powq(a, b)
103 #define PetscExpReal(a)         expq(a)
104 #define PetscLogReal(a)         logq(a)
105 #define PetscLog10Real(a)       log10q(a)
106 #define PetscLog2Real(a)        log2q(a)
107 #define PetscSinReal(a)         sinq(a)
108 #define PetscCosReal(a)         cosq(a)
109 #define PetscTanReal(a)         tanq(a)
110 #define PetscAsinReal(a)        asinq(a)
111 #define PetscAcosReal(a)        acosq(a)
112 #define PetscAtanReal(a)        atanq(a)
113 #define PetscSinhReal(a)        sinhq(a)
114 #define PetscCoshReal(a)        coshq(a)
115 #define PetscTanhReal(a)        tanhq(a)
116 #define PetscAsinhReal(a)       asinhq(a)
117 #define PetscAcoshReal(a)       acoshq(a)
118 #define PetscAtanhReal(a)       atanhq(a)
119 #define PetscErfReal(a)         erfq(a)
120 #define PetscCeilReal(a)        ceilq(a)
121 #define PetscFloorReal(a)       floorq(a)
122 #define PetscFmodReal(a, b)     fmodq(a, b)
123 #define PetscCopysignReal(a, b) copysignq(a, b)
124 #define PetscTGamma(a)          tgammaq(a)
125 #if defined(PETSC_HAVE_LGAMMA_IS_GAMMA)
126 #define PetscLGamma(a) gammaq(a)
127 #else
128 #define PetscLGamma(a) lgammaq(a)
129 #endif
130 
131 #elif defined(PETSC_USE_REAL___FP16)
132 #define PetscSqrtReal(a)        sqrtf(a)
133 #define PetscCbrtReal(a)        cbrtf(a)
134 #define PetscHypotReal(a, b)    hypotf(a, b)
135 #define PetscAtan2Real(a, b)    atan2f(a, b)
136 #define PetscPowReal(a, b)      powf(a, b)
137 #define PetscExpReal(a)         expf(a)
138 #define PetscLogReal(a)         logf(a)
139 #define PetscLog10Real(a)       log10f(a)
140 #define PetscLog2Real(a)        log2f(a)
141 #define PetscSinReal(a)         sinf(a)
142 #define PetscCosReal(a)         cosf(a)
143 #define PetscTanReal(a)         tanf(a)
144 #define PetscAsinReal(a)        asinf(a)
145 #define PetscAcosReal(a)        acosf(a)
146 #define PetscAtanReal(a)        atanf(a)
147 #define PetscSinhReal(a)        sinhf(a)
148 #define PetscCoshReal(a)        coshf(a)
149 #define PetscTanhReal(a)        tanhf(a)
150 #define PetscAsinhReal(a)       asinhf(a)
151 #define PetscAcoshReal(a)       acoshf(a)
152 #define PetscAtanhReal(a)       atanhf(a)
153 #define PetscErfReal(a)         erff(a)
154 #define PetscCeilReal(a)        ceilf(a)
155 #define PetscFloorReal(a)       floorf(a)
156 #define PetscFmodReal(a, b)     fmodf(a, b)
157 #define PetscCopySignReal(a, b) copysignf(a, b)
158 #define PetscTGamma(a)          tgammaf(a)
159 #if defined(PETSC_HAVE_LGAMMA_IS_GAMMA)
160 #define PetscLGamma(a) gammaf(a)
161 #else
162 #define PetscLGamma(a) lgammaf(a)
163 #endif
164 
165 #endif /* PETSC_USE_REAL_* */
166 
167 static inline PetscReal PetscSignReal(PetscReal a) {
168   return (PetscReal)((a < (PetscReal)0) ? -1 : ((a > (PetscReal)0) ? 1 : 0));
169 }
170 
171 #if !defined(PETSC_HAVE_LOG2)
172 #undef PetscLog2Real
173 static inline PetscReal PetscLog2Real(PetscReal a) {
174   return PetscLogReal(a) / PetscLogReal((PetscReal)2);
175 }
176 #endif
177 
178 #if defined(PETSC_HAVE_REAL___FLOAT128)
179 PETSC_EXTERN MPI_Datatype MPIU___FLOAT128 PetscAttrMPITypeTag(__float128);
180 #endif
181 #if defined(PETSC_USE_REAL___FP16)
182 PETSC_EXTERN MPI_Datatype MPIU___FP16 PetscAttrMPITypeTag(__fp16);
183 #endif
184 
185 /*MC
186    MPIU_REAL - Portable MPI datatype corresponding to `PetscReal` independent of what precision `PetscReal` is in
187 
188    Notes:
189    In MPI calls that require an MPI datatype that matches a `PetscReal` or array of `PetscReal` values, pass this value.
190 
191    Level: beginner
192 
193 .seealso: `PetscReal`, `PetscScalar`, `PetscComplex`, `PetscInt`, `MPIU_SCALAR`, `MPIU_COMPLEX`, `MPIU_INT`
194 M*/
195 #if defined(PETSC_USE_REAL_SINGLE)
196 #define MPIU_REAL MPI_FLOAT
197 #elif defined(PETSC_USE_REAL_DOUBLE)
198 #define MPIU_REAL MPI_DOUBLE
199 #elif defined(PETSC_USE_REAL___FLOAT128)
200 #define MPIU_REAL MPIU___FLOAT128
201 #elif defined(PETSC_USE_REAL___FP16)
202 #define MPIU_REAL MPIU___FP16
203 #endif /* PETSC_USE_REAL_* */
204 
205 /*
206     Complex number definitions
207  */
208 #if defined(PETSC_HAVE_COMPLEX)
209 #if defined(__cplusplus) && !defined(PETSC_USE_REAL___FLOAT128)
210 /* C++ support of complex number */
211 
212 #define PetscRealPartComplex(a)      (a).real()
213 #define PetscImaginaryPartComplex(a) (a).imag()
214 #define PetscAbsComplex(a)           petsccomplexlib::abs(a)
215 #define PetscArgComplex(a)           petsccomplexlib::arg(a)
216 #define PetscConjComplex(a)          petsccomplexlib::conj(a)
217 #define PetscSqrtComplex(a)          petsccomplexlib::sqrt(a)
218 #define PetscPowComplex(a, b)        petsccomplexlib::pow(a, b)
219 #define PetscExpComplex(a)           petsccomplexlib::exp(a)
220 #define PetscLogComplex(a)           petsccomplexlib::log(a)
221 #define PetscSinComplex(a)           petsccomplexlib::sin(a)
222 #define PetscCosComplex(a)           petsccomplexlib::cos(a)
223 #define PetscTanComplex(a)           petsccomplexlib::tan(a)
224 #define PetscAsinComplex(a)          petsccomplexlib::asin(a)
225 #define PetscAcosComplex(a)          petsccomplexlib::acos(a)
226 #define PetscAtanComplex(a)          petsccomplexlib::atan(a)
227 #define PetscSinhComplex(a)          petsccomplexlib::sinh(a)
228 #define PetscCoshComplex(a)          petsccomplexlib::cosh(a)
229 #define PetscTanhComplex(a)          petsccomplexlib::tanh(a)
230 #define PetscAsinhComplex(a)         petsccomplexlib::asinh(a)
231 #define PetscAcoshComplex(a)         petsccomplexlib::acosh(a)
232 #define PetscAtanhComplex(a)         petsccomplexlib::atanh(a)
233 
234 /* TODO: Add configure tests
235 
236 #if !defined(PETSC_HAVE_CXX_TAN_COMPLEX)
237 #undef PetscTanComplex
238 static inline PetscComplex PetscTanComplex(PetscComplex z)
239 {
240   return PetscSinComplex(z)/PetscCosComplex(z);
241 }
242 #endif
243 
244 #if !defined(PETSC_HAVE_CXX_TANH_COMPLEX)
245 #undef PetscTanhComplex
246 static inline PetscComplex PetscTanhComplex(PetscComplex z)
247 {
248   return PetscSinhComplex(z)/PetscCoshComplex(z);
249 }
250 #endif
251 
252 #if !defined(PETSC_HAVE_CXX_ASIN_COMPLEX)
253 #undef PetscAsinComplex
254 static inline PetscComplex PetscAsinComplex(PetscComplex z)
255 {
256   const PetscComplex j(0,1);
257   return -j*PetscLogComplex(j*z+PetscSqrtComplex(1.0f-z*z));
258 }
259 #endif
260 
261 #if !defined(PETSC_HAVE_CXX_ACOS_COMPLEX)
262 #undef PetscAcosComplex
263 static inline PetscComplex PetscAcosComplex(PetscComplex z)
264 {
265   const PetscComplex j(0,1);
266   return j*PetscLogComplex(z-j*PetscSqrtComplex(1.0f-z*z));
267 }
268 #endif
269 
270 #if !defined(PETSC_HAVE_CXX_ATAN_COMPLEX)
271 #undef PetscAtanComplex
272 static inline PetscComplex PetscAtanComplex(PetscComplex z)
273 {
274   const PetscComplex j(0,1);
275   return 0.5f*j*PetscLogComplex((1.0f-j*z)/(1.0f+j*z));
276 }
277 #endif
278 
279 #if !defined(PETSC_HAVE_CXX_ASINH_COMPLEX)
280 #undef PetscAsinhComplex
281 static inline PetscComplex PetscAsinhComplex(PetscComplex z)
282 {
283   return PetscLogComplex(z+PetscSqrtComplex(z*z+1.0f));
284 }
285 #endif
286 
287 #if !defined(PETSC_HAVE_CXX_ACOSH_COMPLEX)
288 #undef PetscAcoshComplex
289 static inline PetscComplex PetscAcoshComplex(PetscComplex z)
290 {
291   return PetscLogComplex(z+PetscSqrtComplex(z*z-1.0f));
292 }
293 #endif
294 
295 #if !defined(PETSC_HAVE_CXX_ATANH_COMPLEX)
296 #undef PetscAtanhComplex
297 static inline PetscComplex PetscAtanhComplex(PetscComplex z)
298 {
299   return 0.5f*PetscLogComplex((1.0f+z)/(1.0f-z));
300 }
301 #endif
302 
303 */
304 
305 #else /* C99 support of complex number */
306 
307 #if defined(PETSC_USE_REAL_SINGLE)
308 #define PetscRealPartComplex(a)      crealf(a)
309 #define PetscImaginaryPartComplex(a) cimagf(a)
310 #define PetscAbsComplex(a)           cabsf(a)
311 #define PetscArgComplex(a)           cargf(a)
312 #define PetscConjComplex(a)          conjf(a)
313 #define PetscSqrtComplex(a)          csqrtf(a)
314 #define PetscPowComplex(a, b)        cpowf(a, b)
315 #define PetscExpComplex(a)           cexpf(a)
316 #define PetscLogComplex(a)           clogf(a)
317 #define PetscSinComplex(a)           csinf(a)
318 #define PetscCosComplex(a)           ccosf(a)
319 #define PetscTanComplex(a)           ctanf(a)
320 #define PetscAsinComplex(a)          casinf(a)
321 #define PetscAcosComplex(a)          cacosf(a)
322 #define PetscAtanComplex(a)          catanf(a)
323 #define PetscSinhComplex(a)          csinhf(a)
324 #define PetscCoshComplex(a)          ccoshf(a)
325 #define PetscTanhComplex(a)          ctanhf(a)
326 #define PetscAsinhComplex(a)         casinhf(a)
327 #define PetscAcoshComplex(a)         cacoshf(a)
328 #define PetscAtanhComplex(a)         catanhf(a)
329 
330 #elif defined(PETSC_USE_REAL_DOUBLE)
331 #define PetscRealPartComplex(a)      creal(a)
332 #define PetscImaginaryPartComplex(a) cimag(a)
333 #define PetscAbsComplex(a)           cabs(a)
334 #define PetscArgComplex(a)           carg(a)
335 #define PetscConjComplex(a)          conj(a)
336 #define PetscSqrtComplex(a)          csqrt(a)
337 #define PetscPowComplex(a, b)        cpow(a, b)
338 #define PetscExpComplex(a)           cexp(a)
339 #define PetscLogComplex(a)           clog(a)
340 #define PetscSinComplex(a)           csin(a)
341 #define PetscCosComplex(a)           ccos(a)
342 #define PetscTanComplex(a)           ctan(a)
343 #define PetscAsinComplex(a)          casin(a)
344 #define PetscAcosComplex(a)          cacos(a)
345 #define PetscAtanComplex(a)          catan(a)
346 #define PetscSinhComplex(a)          csinh(a)
347 #define PetscCoshComplex(a)          ccosh(a)
348 #define PetscTanhComplex(a)          ctanh(a)
349 #define PetscAsinhComplex(a)         casinh(a)
350 #define PetscAcoshComplex(a)         cacosh(a)
351 #define PetscAtanhComplex(a)         catanh(a)
352 
353 #elif defined(PETSC_USE_REAL___FLOAT128)
354 #define PetscRealPartComplex(a)      crealq(a)
355 #define PetscImaginaryPartComplex(a) cimagq(a)
356 #define PetscAbsComplex(a)           cabsq(a)
357 #define PetscArgComplex(a)           cargq(a)
358 #define PetscConjComplex(a)          conjq(a)
359 #define PetscSqrtComplex(a)          csqrtq(a)
360 #define PetscPowComplex(a, b)        cpowq(a, b)
361 #define PetscExpComplex(a)           cexpq(a)
362 #define PetscLogComplex(a)           clogq(a)
363 #define PetscSinComplex(a)           csinq(a)
364 #define PetscCosComplex(a)           ccosq(a)
365 #define PetscTanComplex(a)           ctanq(a)
366 #define PetscAsinComplex(a)          casinq(a)
367 #define PetscAcosComplex(a)          cacosq(a)
368 #define PetscAtanComplex(a)          catanq(a)
369 #define PetscSinhComplex(a)          csinhq(a)
370 #define PetscCoshComplex(a)          ccoshq(a)
371 #define PetscTanhComplex(a)          ctanhq(a)
372 #define PetscAsinhComplex(a)         casinhq(a)
373 #define PetscAcoshComplex(a)         cacoshq(a)
374 #define PetscAtanhComplex(a)         catanhq(a)
375 
376 #endif /* PETSC_USE_REAL_* */
377 #endif /* (__cplusplus) */
378 
379 /*
380    PETSC_i is the imaginary number, i
381 */
382 PETSC_EXTERN PetscComplex PETSC_i;
383 
384 /*
385    Try to do the right thing for complex number construction: see
386    http://www.open-std.org/jtc1/sc22/wg14/www/docs/n1464.htm
387    for details
388 */
389 static inline PetscComplex PetscCMPLX(PetscReal x, PetscReal y) {
390 #if defined(__cplusplus) && !defined(PETSC_USE_REAL___FLOAT128)
391   return PetscComplex(x, y);
392 #elif defined(_Imaginary_I)
393   return x + y * _Imaginary_I;
394 #else
395   { /* In both C99 and C11 (ISO/IEC 9899, Section 6.2.5),
396 
397        "For each floating type there is a corresponding real type, which is always a real floating
398        type. For real floating types, it is the same type. For complex types, it is the type given
399        by deleting the keyword _Complex from the type name."
400 
401        So type punning should be portable. */
402     union
403     {
404       PetscComplex z;
405       PetscReal    f[2];
406     } uz;
407 
408     uz.f[0] = x;
409     uz.f[1] = y;
410     return uz.z;
411   }
412 #endif
413 }
414 
415 #define MPIU_C_COMPLEX        MPI_C_COMPLEX PETSC_DEPRECATED_MACRO("GCC warning \"MPIU_C_COMPLEX macro is deprecated use MPI_C_COMPLEX (since version 3.15)\"")
416 #define MPIU_C_DOUBLE_COMPLEX MPI_C_DOUBLE_COMPLEX PETSC_DEPRECATED_MACRO("GCC warning \"MPIU_C_DOUBLE_COMPLEX macro is deprecated use MPI_C_DOUBLE_COMPLEX (since version 3.15)\"")
417 
418 #if defined(PETSC_HAVE_REAL___FLOAT128)
419 PETSC_EXTERN MPI_Datatype MPIU___COMPLEX128 PetscAttrMPITypeTag(__complex128);
420 #endif /* PETSC_HAVE_REAL___FLOAT128 */
421 
422 /*MC
423    MPIU_COMPLEX - Portable MPI datatype corresponding to `PetscComplex` independent of the precision of `PetscComplex`
424 
425    Notes:
426    In MPI calls that require an MPI datatype that matches a `PetscComplex` or array of `PetscComplex` values, pass this value.
427 
428    Level: beginner
429 
430 .seealso: `PetscReal`, `PetscScalar`, `PetscComplex`, `PetscInt`, `MPIU_REAL`, `MPIU_SCALAR`, `MPIU_COMPLEX`, `MPIU_INT`, `PETSC_i`
431 M*/
432 #if defined(PETSC_USE_REAL_SINGLE)
433 #define MPIU_COMPLEX MPI_C_COMPLEX
434 #elif defined(PETSC_USE_REAL_DOUBLE)
435 #define MPIU_COMPLEX MPI_C_DOUBLE_COMPLEX
436 #elif defined(PETSC_USE_REAL___FLOAT128)
437 #define MPIU_COMPLEX MPIU___COMPLEX128
438 #elif defined(PETSC_USE_REAL___FP16)
439 #define MPIU_COMPLEX MPI_C_COMPLEX
440 #endif /* PETSC_USE_REAL_* */
441 
442 #endif /* PETSC_HAVE_COMPLEX */
443 
444 /*
445     Scalar number definitions
446  */
447 #if defined(PETSC_USE_COMPLEX) && defined(PETSC_HAVE_COMPLEX)
448 /*MC
449    MPIU_SCALAR - Portable MPI datatype corresponding to `PetscScalar` independent of the precision of `PetscScalar`
450 
451    Notes:
452    In MPI calls that require an MPI datatype that matches a `PetscScalar` or array of `PetscScalar` values, pass this value.
453 
454    Level: beginner
455 
456 .seealso: `PetscReal`, `PetscScalar`, `PetscComplex`, `PetscInt`, `MPIU_REAL`, `MPIU_COMPLEX`, `MPIU_INT`
457 M*/
458 #define MPIU_SCALAR MPIU_COMPLEX
459 
460 /*MC
461    PetscRealPart - Returns the real part of a `PetscScalar`
462 
463    Synopsis:
464    #include <petscmath.h>
465    PetscReal PetscRealPart(PetscScalar v)
466 
467    Not Collective
468 
469    Input Parameter:
470 .  v - value to find the real part of
471 
472    Level: beginner
473 
474 .seealso: `PetscScalar`, `PetscImaginaryPart()`, `PetscMax()`, `PetscClipInterval()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscSqr()`
475 
476 M*/
477 #define PetscRealPart(a) PetscRealPartComplex(a)
478 
479 /*MC
480    PetscImaginaryPart - Returns the imaginary part of a `PetscScalar`
481 
482    Synopsis:
483    #include <petscmath.h>
484    PetscReal PetscImaginaryPart(PetscScalar v)
485 
486    Not Collective
487 
488    Input Parameter:
489 .  v - value to find the imaginary part of
490 
491    Level: beginner
492 
493    Notes:
494        If PETSc was configured for real numbers then this always returns the value 0
495 
496 .seealso: `PetscScalar`, `PetscRealPart()`, `PetscMax()`, `PetscClipInterval()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscSqr()`
497 
498 M*/
499 #define PetscImaginaryPart(a) PetscImaginaryPartComplex(a)
500 
501 #define PetscAbsScalar(a)    PetscAbsComplex(a)
502 #define PetscArgScalar(a)    PetscArgComplex(a)
503 #define PetscConj(a)         PetscConjComplex(a)
504 #define PetscSqrtScalar(a)   PetscSqrtComplex(a)
505 #define PetscPowScalar(a, b) PetscPowComplex(a, b)
506 #define PetscExpScalar(a)    PetscExpComplex(a)
507 #define PetscLogScalar(a)    PetscLogComplex(a)
508 #define PetscSinScalar(a)    PetscSinComplex(a)
509 #define PetscCosScalar(a)    PetscCosComplex(a)
510 #define PetscTanScalar(a)    PetscTanComplex(a)
511 #define PetscAsinScalar(a)   PetscAsinComplex(a)
512 #define PetscAcosScalar(a)   PetscAcosComplex(a)
513 #define PetscAtanScalar(a)   PetscAtanComplex(a)
514 #define PetscSinhScalar(a)   PetscSinhComplex(a)
515 #define PetscCoshScalar(a)   PetscCoshComplex(a)
516 #define PetscTanhScalar(a)   PetscTanhComplex(a)
517 #define PetscAsinhScalar(a)  PetscAsinhComplex(a)
518 #define PetscAcoshScalar(a)  PetscAcoshComplex(a)
519 #define PetscAtanhScalar(a)  PetscAtanhComplex(a)
520 
521 #else /* PETSC_USE_COMPLEX */
522 #define MPIU_SCALAR           MPIU_REAL
523 #define PetscRealPart(a)      (a)
524 #define PetscImaginaryPart(a) ((PetscReal)0)
525 #define PetscAbsScalar(a)     PetscAbsReal(a)
526 #define PetscArgScalar(a)     (((a) < (PetscReal)0) ? PETSC_PI : (PetscReal)0)
527 #define PetscConj(a)          (a)
528 #define PetscSqrtScalar(a)    PetscSqrtReal(a)
529 #define PetscPowScalar(a, b)  PetscPowReal(a, b)
530 #define PetscExpScalar(a)     PetscExpReal(a)
531 #define PetscLogScalar(a)     PetscLogReal(a)
532 #define PetscSinScalar(a)     PetscSinReal(a)
533 #define PetscCosScalar(a)     PetscCosReal(a)
534 #define PetscTanScalar(a)     PetscTanReal(a)
535 #define PetscAsinScalar(a)    PetscAsinReal(a)
536 #define PetscAcosScalar(a)    PetscAcosReal(a)
537 #define PetscAtanScalar(a)    PetscAtanReal(a)
538 #define PetscSinhScalar(a)    PetscSinhReal(a)
539 #define PetscCoshScalar(a)    PetscCoshReal(a)
540 #define PetscTanhScalar(a)    PetscTanhReal(a)
541 #define PetscAsinhScalar(a)   PetscAsinhReal(a)
542 #define PetscAcoshScalar(a)   PetscAcoshReal(a)
543 #define PetscAtanhScalar(a)   PetscAtanhReal(a)
544 
545 #endif /* PETSC_USE_COMPLEX */
546 
547 /*
548    Certain objects may be created using either single or double precision.
549    This is currently not used.
550 */
551 typedef enum {
552   PETSC_SCALAR_DOUBLE,
553   PETSC_SCALAR_SINGLE,
554   PETSC_SCALAR_LONG_DOUBLE,
555   PETSC_SCALAR_HALF
556 } PetscScalarPrecision;
557 
558 /* --------------------------------------------------------------------------*/
559 
560 /*MC
561    PetscAbs - Returns the absolute value of a number
562 
563    Synopsis:
564    #include <petscmath.h>
565    type PetscAbs(type v)
566 
567    Not Collective
568 
569    Input Parameter:
570 .  v - the number
571 
572    Note:
573    The type can be integer or real floating point value, but cannot be complex
574 
575    Level: beginner
576 
577 .seealso: `PetscAbsInt()`, `PetscAbsReal()`, `PetscAbsScalar()`
578 
579 M*/
580 #define PetscAbs(a) (((a) >= 0) ? (a) : (-(a)))
581 
582 /*MC
583    PetscSign - Returns the sign of a number as an integer
584 
585    Synopsis:
586    #include <petscmath.h>
587    int PetscSign(type v)
588 
589    Not Collective
590 
591    Input Parameter:
592 .  v - the number
593 
594    Note:
595    The type can be integer or real floating point value
596 
597    Level: beginner
598 
599 M*/
600 #define PetscSign(a) (((a) >= 0) ? ((a) == 0 ? 0 : 1) : -1)
601 
602 /*MC
603    PetscMin - Returns minimum of two numbers
604 
605    Synopsis:
606    #include <petscmath.h>
607    type PetscMin(type v1,type v2)
608 
609    Not Collective
610 
611    Input Parameters:
612 +  v1 - first value to find minimum of
613 -  v2 - second value to find minimum of
614 
615    Note:
616    The type can be integer or floating point value
617 
618    Level: beginner
619 
620 .seealso: `PetscMax()`, `PetscClipInterval()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscSqr()`
621 
622 M*/
623 #define PetscMin(a, b) (((a) < (b)) ? (a) : (b))
624 
625 /*MC
626    PetscMax - Returns maxium of two numbers
627 
628    Synopsis:
629    #include <petscmath.h>
630    type max PetscMax(type v1,type v2)
631 
632    Not Collective
633 
634    Input Parameters:
635 +  v1 - first value to find maximum of
636 -  v2 - second value to find maximum of
637 
638    Note:
639    The type can be integer or floating point value
640 
641    Level: beginner
642 
643 .seealso: `PetscMin()`, `PetscClipInterval()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscSqr()`
644 
645 M*/
646 #define PetscMax(a, b) (((a) < (b)) ? (b) : (a))
647 
648 /*MC
649    PetscClipInterval - Returns a number clipped to be within an interval
650 
651    Synopsis:
652    #include <petscmath.h>
653    type clip PetscClipInterval(type x,type a,type b)
654 
655    Not Collective
656 
657    Input Parameters:
658 +  x - value to use if within interval [a,b]
659 .  a - lower end of interval
660 -  b - upper end of interval
661 
662    Note:
663    The type can be integer or floating point value
664 
665    Level: beginner
666 
667 .seealso: `PetscMin()`, `PetscMax()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscSqr()`
668 
669 M*/
670 #define PetscClipInterval(x, a, b) (PetscMax((a), PetscMin((x), (b))))
671 
672 /*MC
673    PetscAbsInt - Returns the absolute value of an integer
674 
675    Synopsis:
676    #include <petscmath.h>
677    int abs PetscAbsInt(int v1)
678 
679    Input Parameter:
680 .   v1 - the integer
681 
682    Level: beginner
683 
684 .seealso: `PetscMax()`, `PetscMin()`, `PetscAbsReal()`, `PetscSqr()`
685 
686 M*/
687 #define PetscAbsInt(a) (((a) < 0) ? (-(a)) : (a))
688 
689 /*MC
690    PetscAbsReal - Returns the absolute value of an real number
691 
692    Synopsis:
693    #include <petscmath.h>
694    Real abs PetscAbsReal(PetscReal v1)
695 
696    Input Parameter:
697 .   v1 - the double
698 
699    Level: beginner
700 
701 .seealso: `PetscMax()`, `PetscMin()`, `PetscAbsInt()`, `PetscSqr()`
702 
703 M*/
704 #if defined(PETSC_USE_REAL_SINGLE)
705 #define PetscAbsReal(a) fabsf(a)
706 #elif defined(PETSC_USE_REAL_DOUBLE)
707 #define PetscAbsReal(a) fabs(a)
708 #elif defined(PETSC_USE_REAL___FLOAT128)
709 #define PetscAbsReal(a) fabsq(a)
710 #elif defined(PETSC_USE_REAL___FP16)
711 #define PetscAbsReal(a) fabsf(a)
712 #endif
713 
714 /*MC
715    PetscSqr - Returns the square of a number
716 
717    Synopsis:
718    #include <petscmath.h>
719    type sqr PetscSqr(type v1)
720 
721    Not Collective
722 
723    Input Parameter:
724 .   v1 - the value
725 
726    Note:
727    The type can be integer or floating point value
728 
729    Level: beginner
730 
731 .seealso: `PetscMax()`, `PetscMin()`, `PetscAbsInt()`, `PetscAbsReal()`
732 
733 M*/
734 #define PetscSqr(a) ((a) * (a))
735 
736 /* ----------------------------------------------------------------------------*/
737 
738 #if defined(PETSC_USE_REAL_SINGLE)
739 #define PetscRealConstant(constant) constant##F
740 #elif defined(PETSC_USE_REAL_DOUBLE)
741 #define PetscRealConstant(constant) constant
742 #elif defined(PETSC_USE_REAL___FLOAT128)
743 #define PetscRealConstant(constant) constant##Q
744 #elif defined(PETSC_USE_REAL___FP16)
745 #define PetscRealConstant(constant) constant##F
746 #endif
747 
748 /*
749      Basic constants
750 */
751 #define PETSC_PI    PetscRealConstant(3.1415926535897932384626433832795029)
752 #define PETSC_PHI   PetscRealConstant(1.6180339887498948482045868343656381)
753 #define PETSC_SQRT2 PetscRealConstant(1.4142135623730950488016887242096981)
754 
755 #if !defined(PETSC_USE_64BIT_INDICES)
756 #define PETSC_MAX_INT 2147483647
757 #define PETSC_MIN_INT (-PETSC_MAX_INT - 1)
758 #else
759 #define PETSC_MAX_INT 9223372036854775807L
760 #define PETSC_MIN_INT (-PETSC_MAX_INT - 1)
761 #endif
762 #define PETSC_MAX_UINT16 65535
763 
764 #if defined(PETSC_USE_REAL_SINGLE)
765 #define PETSC_MAX_REAL             3.40282346638528860e+38F
766 #define PETSC_MIN_REAL             (-PETSC_MAX_REAL)
767 #define PETSC_MACHINE_EPSILON      1.19209290e-07F
768 #define PETSC_SQRT_MACHINE_EPSILON 3.45266983e-04F
769 #define PETSC_SMALL                1.e-5F
770 #elif defined(PETSC_USE_REAL_DOUBLE)
771 #define PETSC_MAX_REAL             1.7976931348623157e+308
772 #define PETSC_MIN_REAL             (-PETSC_MAX_REAL)
773 #define PETSC_MACHINE_EPSILON      2.2204460492503131e-16
774 #define PETSC_SQRT_MACHINE_EPSILON 1.490116119384766e-08
775 #define PETSC_SMALL                1.e-10
776 #elif defined(PETSC_USE_REAL___FLOAT128)
777 #define PETSC_MAX_REAL             FLT128_MAX
778 #define PETSC_MIN_REAL             (-FLT128_MAX)
779 #define PETSC_MACHINE_EPSILON      FLT128_EPSILON
780 #define PETSC_SQRT_MACHINE_EPSILON 1.38777878078144567552953958511352539e-17Q
781 #define PETSC_SMALL                1.e-20Q
782 #elif defined(PETSC_USE_REAL___FP16)
783 #define PETSC_MAX_REAL             65504.0F
784 #define PETSC_MIN_REAL             (-PETSC_MAX_REAL)
785 #define PETSC_MACHINE_EPSILON      .0009765625F
786 #define PETSC_SQRT_MACHINE_EPSILON .03125F
787 #define PETSC_SMALL                5.e-3F
788 #endif
789 
790 #define PETSC_INFINITY  (PETSC_MAX_REAL / 4)
791 #define PETSC_NINFINITY (-PETSC_INFINITY)
792 
793 PETSC_EXTERN PetscBool  PetscIsInfReal(PetscReal);
794 PETSC_EXTERN PetscBool  PetscIsNanReal(PetscReal);
795 PETSC_EXTERN PetscBool  PetscIsNormalReal(PetscReal);
796 static inline PetscBool PetscIsInfOrNanReal(PetscReal v) {
797   return PetscIsInfReal(v) || PetscIsNanReal(v) ? PETSC_TRUE : PETSC_FALSE;
798 }
799 static inline PetscBool PetscIsInfScalar(PetscScalar v) {
800   return PetscIsInfReal(PetscAbsScalar(v));
801 }
802 static inline PetscBool PetscIsNanScalar(PetscScalar v) {
803   return PetscIsNanReal(PetscAbsScalar(v));
804 }
805 static inline PetscBool PetscIsInfOrNanScalar(PetscScalar v) {
806   return PetscIsInfOrNanReal(PetscAbsScalar(v));
807 }
808 static inline PetscBool PetscIsNormalScalar(PetscScalar v) {
809   return PetscIsNormalReal(PetscAbsScalar(v));
810 }
811 
812 PETSC_EXTERN PetscBool PetscIsCloseAtTol(PetscReal, PetscReal, PetscReal, PetscReal);
813 PETSC_EXTERN PetscBool PetscEqualReal(PetscReal, PetscReal);
814 PETSC_EXTERN PetscBool PetscEqualScalar(PetscScalar, PetscScalar);
815 
816 /*
817     These macros are currently hardwired to match the regular data types, so there is no support for a different
818     MatScalar from PetscScalar. We left the MatScalar in the source just in case we use it again.
819  */
820 #define MPIU_MATSCALAR MPIU_SCALAR
821 typedef PetscScalar MatScalar;
822 typedef PetscReal   MatReal;
823 
824 struct petsc_mpiu_2scalar {
825   PetscScalar a, b;
826 };
827 PETSC_EXTERN MPI_Datatype MPIU_2SCALAR PetscAttrMPITypeTagLayoutCompatible(struct petsc_mpiu_2scalar);
828 
829 /*
830    MPI Datatypes for composite reductions:
831    MPIU_REAL_INT -> struct { PetscReal; PetscInt; }
832    MPIU_SCALAR_INT -> struct { PetscScalar; PetscInt; }
833 */
834 PETSC_EXTERN MPI_Datatype MPIU_REAL_INT;
835 PETSC_EXTERN MPI_Datatype MPIU_SCALAR_INT;
836 
837 #if defined(PETSC_USE_64BIT_INDICES)
838 struct petsc_mpiu_2int {
839   PetscInt a, b;
840 };
841 PETSC_EXTERN MPI_Datatype MPIU_2INT PetscAttrMPITypeTagLayoutCompatible(struct petsc_mpiu_2int);
842 #else
843 #define MPIU_2INT MPI_2INT
844 #endif
845 PETSC_EXTERN MPI_Datatype MPI_4INT;
846 PETSC_EXTERN MPI_Datatype MPIU_4INT;
847 
848 static inline PetscInt PetscPowInt(PetscInt base, PetscInt power) {
849   PetscInt result = 1;
850   while (power) {
851     if (power & 1) result *= base;
852     power >>= 1;
853     base *= base;
854   }
855   return result;
856 }
857 
858 static inline PetscInt64 PetscPowInt64(PetscInt base, PetscInt power) {
859   PetscInt64 result = 1;
860   while (power) {
861     if (power & 1) result *= base;
862     power >>= 1;
863     base *= base;
864   }
865   return result;
866 }
867 
868 static inline PetscReal PetscPowRealInt(PetscReal base, PetscInt power) {
869   PetscReal result = 1;
870   if (power < 0) {
871     power = -power;
872     base  = ((PetscReal)1) / base;
873   }
874   while (power) {
875     if (power & 1) result *= base;
876     power >>= 1;
877     base *= base;
878   }
879   return result;
880 }
881 
882 static inline PetscScalar PetscPowScalarInt(PetscScalar base, PetscInt power) {
883   PetscScalar result = (PetscReal)1;
884   if (power < 0) {
885     power = -power;
886     base  = ((PetscReal)1) / base;
887   }
888   while (power) {
889     if (power & 1) result *= base;
890     power >>= 1;
891     base *= base;
892   }
893   return result;
894 }
895 
896 static inline PetscScalar PetscPowScalarReal(PetscScalar base, PetscReal power) {
897   PetscScalar cpower = power;
898   return PetscPowScalar(base, cpower);
899 }
900 
901 /*MC
902     PetscApproximateLTE - Performs a less than or equal to on a given constant with a fudge for floating point numbers
903 
904    Synopsis:
905    #include <petscmath.h>
906    bool PetscApproximateLTE(PetscReal x,constant float)
907 
908    Not Collective
909 
910    Input Parameters:
911 +   x - the variable
912 -   b - the constant float it is checking if x is less than or equal to
913 
914    Notes:
915      The fudge factor is the value `PETSC_SMALL`
916 
917      The constant numerical value is automatically set to the appropriate precision of PETSc so can just be provided as, for example, 3.2
918 
919      This is used in several examples for setting initial conditions based on coordinate values that are computed with i*h that produces inexact
920      floating point results.
921 
922    Level: advanced
923 
924 .seealso: `PetscMax()`, `PetscMin()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscApproximateGTE()`
925 
926 M*/
927 #define PetscApproximateLTE(x, b) ((x) <= (PetscRealConstant(b) + PETSC_SMALL))
928 
929 /*MC
930     PetscApproximateGTE - Performs a greater than or equal to on a given constant with a fudge for floating point numbers
931 
932    Synopsis:
933    #include <petscmath.h>
934    bool PetscApproximateGTE(PetscReal x,constant float)
935 
936    Not Collective
937 
938    Input Parameters:
939 +   x - the variable
940 -   b - the constant float it is checking if x is greater than or equal to
941 
942    Notes:
943      The fudge factor is the value `PETSC_SMALL`
944 
945      The constant numerical value is automatically set to the appropriate precision of PETSc so can just be provided as, for example, 3.2
946 
947      This is used in several examples for setting initial conditions based on coordinate values that are computed with i*h that produces inexact
948      floating point results.
949 
950    Level: advanced
951 
952 .seealso: `PetscMax()`, `PetscMin()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscApproximateLTE()`
953 
954 M*/
955 #define PetscApproximateGTE(x, b) ((x) >= (PetscRealConstant(b) - PETSC_SMALL))
956 
957 /*MC
958     PetscCeilInt - Returns the ceiling of the quotation of two positive integers
959 
960    Synopsis:
961    #include <petscmath.h>
962    PetscInt PetscCeilInt(PetscInt x,PetscInt y)
963 
964    Not Collective
965 
966    Input Parameters:
967 +   x - the numerator
968 -   y - the denominator
969 
970    Level: advanced
971 
972 .seealso: `PetscMax()`, `PetscMin()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscApproximateLTE()`
973 
974 M*/
975 #define PetscCeilInt(x, y) ((((PetscInt)(x)) / ((PetscInt)(y))) + ((((PetscInt)(x)) % ((PetscInt)(y))) ? 1 : 0))
976 
977 #define PetscCeilInt64(x, y) ((((PetscInt64)(x)) / ((PetscInt64)(y))) + ((((PetscInt64)(x)) % ((PetscInt64)(y))) ? 1 : 0))
978 
979 PETSC_EXTERN PetscErrorCode PetscLinearRegression(PetscInt, const PetscReal[], const PetscReal[], PetscReal *, PetscReal *);
980 #endif
981