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