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