Actual source code: petscmath.h

  1: /*

  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.

  7:     This file is included by petscsys.h and should not be used directly.

  9: */

 11: #if !defined(PETSCMATH_H)
 12: #define PETSCMATH_H
 13: #include <math.h>
 14: #include <petscsystypes.h>

 16: /*

 18:    Defines operations that are different for complex and real numbers.
 19:    All PETSc objects in one program are built around the object
 20:    PetscScalar which is either always a real or a complex.

 22: */

 24: /*
 25:     Real number definitions
 26:  */
 27: #if defined(PETSC_USE_REAL_SINGLE)
 28: #define PetscSqrtReal(a)       sqrtf(a)
 29: #define PetscCbrtReal(a)       cbrtf(a)
 30: #define PetscHypotReal(a,b)    hypotf(a,b)
 31: #define PetscAtan2Real(a,b)    atan2f(a,b)
 32: #define PetscPowReal(a,b)      powf(a,b)
 33: #define PetscExpReal(a)        expf(a)
 34: #define PetscLogReal(a)        logf(a)
 35: #define PetscLog10Real(a)      log10f(a)
 36: #define PetscLog2Real(a)       log2f(a)
 37: #define PetscSinReal(a)        sinf(a)
 38: #define PetscCosReal(a)        cosf(a)
 39: #define PetscTanReal(a)        tanf(a)
 40: #define PetscAsinReal(a)       asinf(a)
 41: #define PetscAcosReal(a)       acosf(a)
 42: #define PetscAtanReal(a)       atanf(a)
 43: #define PetscSinhReal(a)       sinhf(a)
 44: #define PetscCoshReal(a)       coshf(a)
 45: #define PetscTanhReal(a)       tanhf(a)
 46: #define PetscAsinhReal(a)      asinhf(a)
 47: #define PetscAcoshReal(a)      acoshf(a)
 48: #define PetscAtanhReal(a)      atanhf(a)
 49: #define PetscCeilReal(a)       ceilf(a)
 50: #define PetscFloorReal(a)      floorf(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

 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 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

 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 PetscCeilReal(a)       ceilq(a)
116: #define PetscFloorReal(a)      floorq(a)
117: #define PetscFmodReal(a,b)     fmodq(a,b)
118: #define PetscCopysignReal(a,b) copysignq(a,b)
119: #define PetscTGamma(a)         tgammaq(a)
120: #if defined(PETSC_HAVE_LGAMMA_IS_GAMMA)
121: #define PetscLGamma(a)         gammaq(a)
122: #else
123: #define PetscLGamma(a)         lgammaq(a)
124: #endif

126: #elif defined(PETSC_USE_REAL___FP16)
127: #define PetscSqrtReal(a)       sqrtf(a)
128: #define PetscCbrtReal(a)       cbrtf(a)
129: #define PetscHypotReal(a,b)    hypotf(a,b)
130: #define PetscAtan2Real(a,b)    atan2f(a,b)
131: #define PetscPowReal(a,b)      powf(a,b)
132: #define PetscExpReal(a)        expf(a)
133: #define PetscLogReal(a)        logf(a)
134: #define PetscLog10Real(a)      log10f(a)
135: #define PetscLog2Real(a)       log2f(a)
136: #define PetscSinReal(a)        sinf(a)
137: #define PetscCosReal(a)        cosf(a)
138: #define PetscTanReal(a)        tanf(a)
139: #define PetscAsinReal(a)       asinf(a)
140: #define PetscAcosReal(a)       acosf(a)
141: #define PetscAtanReal(a)       atanf(a)
142: #define PetscSinhReal(a)       sinhf(a)
143: #define PetscCoshReal(a)       coshf(a)
144: #define PetscTanhReal(a)       tanhf(a)
145: #define PetscAsinhReal(a)      asinhf(a)
146: #define PetscAcoshReal(a)      acoshf(a)
147: #define PetscAtanhReal(a)      atanhf(a)
148: #define PetscCeilReal(a)       ceilf(a)
149: #define PetscFloorReal(a)      floorf(a)
150: #define PetscFmodReal(a,b)     fmodf(a,b)
151: #define PetscCopySignReal(a,b) copysignf(a,b)
152: #define PetscTGamma(a)         tgammaf(a)
153: #if defined(PETSC_HAVE_LGAMMA_IS_GAMMA)
154: #define PetscLGamma(a)         gammaf(a)
155: #else
156: #define PetscLGamma(a)         lgammaf(a)
157: #endif

159: #endif /* PETSC_USE_REAL_* */

161: PETSC_STATIC_INLINE PetscReal PetscSignReal(PetscReal a)
162: {
163:   return (PetscReal)((a < (PetscReal)0) ? -1 : ((a > (PetscReal)0) ? 1 : 0));
164: }

166: #if !defined(PETSC_HAVE_LOG2)
167: #undef PetscLog2Real
168: PETSC_STATIC_INLINE PetscReal PetscLog2Real(PetscReal a)
169: {
170:   return PetscLogReal(a)/PetscLogReal((PetscReal)2);
171: }
172: #endif

174: #if defined(PETSC_USE_REAL___FLOAT128)
175: PETSC_EXTERN MPI_Datatype MPIU___FLOAT128 PetscAttrMPITypeTag(__float128);
176: #endif
177: #if defined(PETSC_USE_REAL___FP16)
178: PETSC_EXTERN MPI_Datatype MPIU___FP16 PetscAttrMPITypeTag(__fp16);
179: #endif

181: /*MC
182:    MPIU_REAL - MPI datatype corresponding to PetscReal

184:    Notes:
185:    In MPI calls that require an MPI datatype that matches a PetscReal or array of PetscReal values, pass this value.

187:    Level: beginner

189: .seealso: PetscReal, PetscScalar, PetscComplex, PetscInt, MPIU_SCALAR, MPIU_COMPLEX, MPIU_INT
190: M*/
191: #if defined(PETSC_USE_REAL_SINGLE)
192: #  define MPIU_REAL MPI_FLOAT
193: #elif defined(PETSC_USE_REAL_DOUBLE)
194: #  define MPIU_REAL MPI_DOUBLE
195: #elif defined(PETSC_USE_REAL___FLOAT128)
196: #  define MPIU_REAL MPIU___FLOAT128
197: #elif defined(PETSC_USE_REAL___FP16)
198: #  define MPIU_REAL MPIU___FP16
199: #endif /* PETSC_USE_REAL_* */

201: /*
202:     Complex number definitions
203:  */
204: #if defined(PETSC_HAVE_COMPLEX)
205: #if defined(__cplusplus) && !defined(PETSC_USE_REAL___FLOAT128)
206: /* C++ support of complex number */

208: #define PetscRealPartComplex(a)      (a).real()
209: #define PetscImaginaryPartComplex(a) (a).imag()
210: #define PetscAbsComplex(a)           petsccomplexlib::abs(a)
211: #define PetscArgComplex(a)           petsccomplexlib::arg(a)
212: #define PetscConjComplex(a)          petsccomplexlib::conj(a)
213: #define PetscSqrtComplex(a)          petsccomplexlib::sqrt(a)
214: #define PetscPowComplex(a,b)         petsccomplexlib::pow(a,b)
215: #define PetscExpComplex(a)           petsccomplexlib::exp(a)
216: #define PetscLogComplex(a)           petsccomplexlib::log(a)
217: #define PetscSinComplex(a)           petsccomplexlib::sin(a)
218: #define PetscCosComplex(a)           petsccomplexlib::cos(a)
219: #define PetscTanComplex(a)           petsccomplexlib::tan(a)
220: #define PetscAsinComplex(a)          petsccomplexlib::asin(a)
221: #define PetscAcosComplex(a)          petsccomplexlib::acos(a)
222: #define PetscAtanComplex(a)          petsccomplexlib::atan(a)
223: #define PetscSinhComplex(a)          petsccomplexlib::sinh(a)
224: #define PetscCoshComplex(a)          petsccomplexlib::cosh(a)
225: #define PetscTanhComplex(a)          petsccomplexlib::tanh(a)
226: #define PetscAsinhComplex(a)         petsccomplexlib::asinh(a)
227: #define PetscAcoshComplex(a)         petsccomplexlib::acosh(a)
228: #define PetscAtanhComplex(a)         petsccomplexlib::atanh(a)

230: /* TODO: Add configure tests

232: #if !defined(PETSC_HAVE_CXX_TAN_COMPLEX)
233: #undef PetscTanComplex
234: PETSC_STATIC_INLINE PetscComplex PetscTanComplex(PetscComplex z)
235: {
236:   return PetscSinComplex(z)/PetscCosComplex(z);
237: }
238: #endif

240: #if !defined(PETSC_HAVE_CXX_TANH_COMPLEX)
241: #undef PetscTanhComplex
242: PETSC_STATIC_INLINE PetscComplex PetscTanhComplex(PetscComplex z)
243: {
244:   return PetscSinhComplex(z)/PetscCoshComplex(z);
245: }
246: #endif

248: #if !defined(PETSC_HAVE_CXX_ASIN_COMPLEX)
249: #undef PetscAsinComplex
250: PETSC_STATIC_INLINE PetscComplex PetscAsinComplex(PetscComplex z)
251: {
252:   const PetscComplex j(0,1);
253:   return -j*PetscLogComplex(j*z+PetscSqrtComplex(1.0f-z*z));
254: }
255: #endif

257: #if !defined(PETSC_HAVE_CXX_ACOS_COMPLEX)
258: #undef PetscAcosComplex
259: PETSC_STATIC_INLINE PetscComplex PetscAcosComplex(PetscComplex z)
260: {
261:   const PetscComplex j(0,1);
262:   return j*PetscLogComplex(z-j*PetscSqrtComplex(1.0f-z*z));
263: }
264: #endif

266: #if !defined(PETSC_HAVE_CXX_ATAN_COMPLEX)
267: #undef PetscAtanComplex
268: PETSC_STATIC_INLINE PetscComplex PetscAtanComplex(PetscComplex z)
269: {
270:   const PetscComplex j(0,1);
271:   return 0.5f*j*PetscLogComplex((1.0f-j*z)/(1.0f+j*z));
272: }
273: #endif

275: #if !defined(PETSC_HAVE_CXX_ASINH_COMPLEX)
276: #undef PetscAsinhComplex
277: PETSC_STATIC_INLINE PetscComplex PetscAsinhComplex(PetscComplex z)
278: {
279:   return PetscLogComplex(z+PetscSqrtComplex(z*z+1.0f));
280: }
281: #endif

283: #if !defined(PETSC_HAVE_CXX_ACOSH_COMPLEX)
284: #undef PetscAcoshComplex
285: PETSC_STATIC_INLINE PetscComplex PetscAcoshComplex(PetscComplex z)
286: {
287:   return PetscLogComplex(z+PetscSqrtComplex(z*z-1.0f));
288: }
289: #endif

291: #if !defined(PETSC_HAVE_CXX_ATANH_COMPLEX)
292: #undef PetscAtanhComplex
293: PETSC_STATIC_INLINE PetscComplex PetscAtanhComplex(PetscComplex z)
294: {
295:   return 0.5f*PetscLogComplex((1.0f+z)/(1.0f-z));
296: }
297: #endif

299: */

301: #else /* C99 support of complex number */

303: #if defined(PETSC_USE_REAL_SINGLE)
304: #define PetscRealPartComplex(a)      crealf(a)
305: #define PetscImaginaryPartComplex(a) cimagf(a)
306: #define PetscAbsComplex(a)           cabsf(a)
307: #define PetscArgComplex(a)           cargf(a)
308: #define PetscConjComplex(a)          conjf(a)
309: #define PetscSqrtComplex(a)          csqrtf(a)
310: #define PetscPowComplex(a,b)         cpowf(a,b)
311: #define PetscExpComplex(a)           cexpf(a)
312: #define PetscLogComplex(a)           clogf(a)
313: #define PetscSinComplex(a)           csinf(a)
314: #define PetscCosComplex(a)           ccosf(a)
315: #define PetscTanComplex(a)           ctanf(a)
316: #define PetscAsinComplex(a)          casinf(a)
317: #define PetscAcosComplex(a)          cacosf(a)
318: #define PetscAtanComplex(a)          catanf(a)
319: #define PetscSinhComplex(a)          csinhf(a)
320: #define PetscCoshComplex(a)          ccoshf(a)
321: #define PetscTanhComplex(a)          ctanhf(a)
322: #define PetscAsinhComplex(a)         casinhf(a)
323: #define PetscAcoshComplex(a)         cacoshf(a)
324: #define PetscAtanhComplex(a)         catanhf(a)

326: #elif defined(PETSC_USE_REAL_DOUBLE)
327: #define PetscRealPartComplex(a)      creal(a)
328: #define PetscImaginaryPartComplex(a) cimag(a)
329: #define PetscAbsComplex(a)           cabs(a)
330: #define PetscArgComplex(a)           carg(a)
331: #define PetscConjComplex(a)          conj(a)
332: #define PetscSqrtComplex(a)          csqrt(a)
333: #define PetscPowComplex(a,b)         cpow(a,b)
334: #define PetscExpComplex(a)           cexp(a)
335: #define PetscLogComplex(a)           clog(a)
336: #define PetscSinComplex(a)           csin(a)
337: #define PetscCosComplex(a)           ccos(a)
338: #define PetscTanComplex(a)           ctan(a)
339: #define PetscAsinComplex(a)          casin(a)
340: #define PetscAcosComplex(a)          cacos(a)
341: #define PetscAtanComplex(a)          catan(a)
342: #define PetscSinhComplex(a)          csinh(a)
343: #define PetscCoshComplex(a)          ccosh(a)
344: #define PetscTanhComplex(a)          ctanh(a)
345: #define PetscAsinhComplex(a)         casinh(a)
346: #define PetscAcoshComplex(a)         cacosh(a)
347: #define PetscAtanhComplex(a)         catanh(a)

349: #elif defined(PETSC_USE_REAL___FLOAT128)
350: #define PetscRealPartComplex(a)      crealq(a)
351: #define PetscImaginaryPartComplex(a) cimagq(a)
352: #define PetscAbsComplex(a)           cabsq(a)
353: #define PetscArgComplex(a)           cargq(a)
354: #define PetscConjComplex(a)          conjq(a)
355: #define PetscSqrtComplex(a)          csqrtq(a)
356: #define PetscPowComplex(a,b)         cpowq(a,b)
357: #define PetscExpComplex(a)           cexpq(a)
358: #define PetscLogComplex(a)           clogq(a)
359: #define PetscSinComplex(a)           csinq(a)
360: #define PetscCosComplex(a)           ccosq(a)
361: #define PetscTanComplex(a)           ctanq(a)
362: #define PetscAsinComplex(a)          casinq(a)
363: #define PetscAcosComplex(a)          cacosq(a)
364: #define PetscAtanComplex(a)          catanq(a)
365: #define PetscSinhComplex(a)          csinhq(a)
366: #define PetscCoshComplex(a)          ccoshq(a)
367: #define PetscTanhComplex(a)          ctanhq(a)
368: #define PetscAsinhComplex(a)         casinhq(a)
369: #define PetscAcoshComplex(a)         cacoshq(a)
370: #define PetscAtanhComplex(a)         catanhq(a)

372: #endif /* PETSC_USE_REAL_* */
373: #endif /* (__cplusplus) */

375: /*
376:    PETSC_i is the imaginary number, i
377: */
378: PETSC_EXTERN PetscComplex PETSC_i;

380: /*
381:    Try to do the right thing for complex number construction: see
382:    http://www.open-std.org/jtc1/sc22/wg14/www/docs/n1464.htm
383:    for details
384: */
385: PETSC_STATIC_INLINE PetscComplex PetscCMPLX(PetscReal x, PetscReal y)
386: {
387: #if defined(__cplusplus) && !defined(PETSC_USE_REAL___FLOAT128)
388:   return PetscComplex(x,y);
389: #elif defined(_Imaginary_I)
390:   return x + y * _Imaginary_I;
391: #else
392:   { /* In both C99 and C11 (ISO/IEC 9899, Section 6.2.5),

394:        "For each floating type there is a corresponding real type, which is always a real floating
395:        type. For real floating types, it is the same type. For complex types, it is the type given
396:        by deleting the keyword _Complex from the type name."

398:        So type punning should be portable. */
399:     union { PetscComplex z; PetscReal f[2]; } uz;

401:     uz.f[0] = x;
402:     uz.f[1] = y;
403:     return uz.z;
404:   }
405: #endif
406: }

408: #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)\"")
409: #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)\"")

411: #if defined(PETSC_USE_REAL___FLOAT128)
412: PETSC_EXTERN MPI_Datatype MPIU___COMPLEX128 PetscAttrMPITypeTag(__complex128);
413: #endif /* PETSC_USE_REAL___FLOAT128 */

415: /*MC
416:    MPIU_COMPLEX - MPI datatype corresponding to PetscComplex

418:    Notes:
419:    In MPI calls that require an MPI datatype that matches a PetscComplex or array of PetscComplex values, pass this value.

421:    Level: beginner

423: .seealso: PetscReal, PetscScalar, PetscComplex, PetscInt, MPIU_REAL, MPIU_SCALAR, MPIU_COMPLEX, MPIU_INT, PETSC_i
424: M*/
425: #if defined(PETSC_USE_REAL_SINGLE)
426: #  define MPIU_COMPLEX MPI_C_COMPLEX
427: #elif defined(PETSC_USE_REAL_DOUBLE)
428: #  define MPIU_COMPLEX MPI_C_DOUBLE_COMPLEX
429: #elif defined(PETSC_USE_REAL___FLOAT128)
430: #  define MPIU_COMPLEX MPIU___COMPLEX128
431: #elif defined(PETSC_USE_REAL___FP16)
432: #  define MPIU_COMPLEX MPI_C_COMPLEX
433: #endif /* PETSC_USE_REAL_* */

435: #endif /* PETSC_HAVE_COMPLEX */

437: /*
438:     Scalar number definitions
439:  */
440: #if defined(PETSC_USE_COMPLEX) && defined(PETSC_HAVE_COMPLEX)
441: /*MC
442:    MPIU_SCALAR - MPI datatype corresponding to PetscScalar

444:    Notes:
445:    In MPI calls that require an MPI datatype that matches a PetscScalar or array of PetscScalar values, pass this value.

447:    Level: beginner

449: .seealso: PetscReal, PetscScalar, PetscComplex, PetscInt, MPIU_REAL, MPIU_COMPLEX, MPIU_INT
450: M*/
451: #define MPIU_SCALAR MPIU_COMPLEX

453: /*MC
454:    PetscRealPart - Returns the real part of a PetscScalar

456:    Synopsis:
457: #include <petscmath.h>
458:    PetscReal PetscRealPart(PetscScalar v)

460:    Not Collective

462:    Input Parameter:
463: .  v - value to find the real part of

465:    Level: beginner

467: .seealso: PetscScalar, PetscImaginaryPart(), PetscMax(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr()

469: M*/
470: #define PetscRealPart(a)      PetscRealPartComplex(a)

472: /*MC
473:    PetscImaginaryPart - Returns the imaginary part of a PetscScalar

475:    Synopsis:
476: #include <petscmath.h>
477:    PetscReal PetscImaginaryPart(PetscScalar v)

479:    Not Collective

481:    Input Parameter:
482: .  v - value to find the imaginary part of

484:    Level: beginner

486:    Notes:
487:        If PETSc was configured for real numbers then this always returns the value 0

489: .seealso: PetscScalar, PetscRealPart(), PetscMax(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr()

491: M*/
492: #define PetscImaginaryPart(a) PetscImaginaryPartComplex(a)

494: #define PetscAbsScalar(a)     PetscAbsComplex(a)
495: #define PetscArgScalar(a)     PetscArgComplex(a)
496: #define PetscConj(a)          PetscConjComplex(a)
497: #define PetscSqrtScalar(a)    PetscSqrtComplex(a)
498: #define PetscPowScalar(a,b)   PetscPowComplex(a,b)
499: #define PetscExpScalar(a)     PetscExpComplex(a)
500: #define PetscLogScalar(a)     PetscLogComplex(a)
501: #define PetscSinScalar(a)     PetscSinComplex(a)
502: #define PetscCosScalar(a)     PetscCosComplex(a)
503: #define PetscTanScalar(a)     PetscTanComplex(a)
504: #define PetscAsinScalar(a)    PetscAsinComplex(a)
505: #define PetscAcosScalar(a)    PetscAcosComplex(a)
506: #define PetscAtanScalar(a)    PetscAtanComplex(a)
507: #define PetscSinhScalar(a)    PetscSinhComplex(a)
508: #define PetscCoshScalar(a)    PetscCoshComplex(a)
509: #define PetscTanhScalar(a)    PetscTanhComplex(a)
510: #define PetscAsinhScalar(a)   PetscAsinhComplex(a)
511: #define PetscAcoshScalar(a)   PetscAcoshComplex(a)
512: #define PetscAtanhScalar(a)   PetscAtanhComplex(a)

514: #else /* PETSC_USE_COMPLEX */
515: #define MPIU_SCALAR MPIU_REAL
516: #define PetscRealPart(a)      (a)
517: #define PetscImaginaryPart(a) ((PetscReal)0)
518: #define PetscAbsScalar(a)     PetscAbsReal(a)
519: #define PetscArgScalar(a)     (((a) < (PetscReal)0) ? PETSC_PI : (PetscReal)0)
520: #define PetscConj(a)          (a)
521: #define PetscSqrtScalar(a)    PetscSqrtReal(a)
522: #define PetscPowScalar(a,b)   PetscPowReal(a,b)
523: #define PetscExpScalar(a)     PetscExpReal(a)
524: #define PetscLogScalar(a)     PetscLogReal(a)
525: #define PetscSinScalar(a)     PetscSinReal(a)
526: #define PetscCosScalar(a)     PetscCosReal(a)
527: #define PetscTanScalar(a)     PetscTanReal(a)
528: #define PetscAsinScalar(a)    PetscAsinReal(a)
529: #define PetscAcosScalar(a)    PetscAcosReal(a)
530: #define PetscAtanScalar(a)    PetscAtanReal(a)
531: #define PetscSinhScalar(a)    PetscSinhReal(a)
532: #define PetscCoshScalar(a)    PetscCoshReal(a)
533: #define PetscTanhScalar(a)    PetscTanhReal(a)
534: #define PetscAsinhScalar(a)   PetscAsinhReal(a)
535: #define PetscAcoshScalar(a)   PetscAcoshReal(a)
536: #define PetscAtanhScalar(a)   PetscAtanhReal(a)

538: #endif /* PETSC_USE_COMPLEX */

540: /*
541:    Certain objects may be created using either single or double precision.
542:    This is currently not used.
543: */
544: typedef enum { PETSC_SCALAR_DOUBLE, PETSC_SCALAR_SINGLE, PETSC_SCALAR_LONG_DOUBLE, PETSC_SCALAR_HALF } PetscScalarPrecision;

546: /* --------------------------------------------------------------------------*/

548: /*MC
549:    PetscAbs - Returns the absolute value of a number

551:    Synopsis:
552: #include <petscmath.h>
553:    type PetscAbs(type v)

555:    Not Collective

557:    Input Parameter:
558: .  v - the number

560:    Notes:
561:     type can be integer or real floating point value

563:    Level: beginner

565: .seealso: PetscAbsInt(), PetscAbsReal(), PetscAbsScalar()

567: M*/
568: #define PetscAbs(a)  (((a) >= 0) ? (a) : (-(a)))

570: /*MC
571:    PetscSign - Returns the sign of a number as an integer

573:    Synopsis:
574: #include <petscmath.h>
575:    int PetscSign(type v)

577:    Not Collective

579:    Input Parameter:
580: .  v - the number

582:    Notes:
583:     type can be integer or real floating point value

585:    Level: beginner

587: M*/
588: #define PetscSign(a) (((a) >= 0) ? ((a) == 0 ? 0 : 1) : -1)

590: /*MC
591:    PetscMin - Returns minimum of two numbers

593:    Synopsis:
594: #include <petscmath.h>
595:    type PetscMin(type v1,type v2)

597:    Not Collective

599:    Input Parameters:
600: +  v1 - first value to find minimum of
601: -  v2 - second value to find minimum of

603:    Notes:
604:     type can be integer or floating point value

606:    Level: beginner

608: .seealso: PetscMax(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr()

610: M*/
611: #define PetscMin(a,b)   (((a)<(b)) ?  (a) : (b))

613: /*MC
614:    PetscMax - Returns maxium of two numbers

616:    Synopsis:
617: #include <petscmath.h>
618:    type max PetscMax(type v1,type v2)

620:    Not Collective

622:    Input Parameters:
623: +  v1 - first value to find maximum of
624: -  v2 - second value to find maximum of

626:    Notes:
627:     type can be integer or floating point value

629:    Level: beginner

631: .seealso: PetscMin(), PetscClipInterval(), PetscAbsInt(), PetscAbsReal(), PetscSqr()

633: M*/
634: #define PetscMax(a,b)   (((a)<(b)) ?  (b) : (a))

636: /*MC
637:    PetscClipInterval - Returns a number clipped to be within an interval

639:    Synopsis:
640: #include <petscmath.h>
641:    type clip PetscClipInterval(type x,type a,type b)

643:    Not Collective

645:    Input Parameters:
646: +  x - value to use if within interval [a,b]
647: .  a - lower end of interval
648: -  b - upper end of interval

650:    Notes:
651:     type can be integer or floating point value

653:    Level: beginner

655: .seealso: PetscMin(), PetscMax(), PetscAbsInt(), PetscAbsReal(), PetscSqr()

657: M*/
658: #define PetscClipInterval(x,a,b)   (PetscMax((a),PetscMin((x),(b))))

660: /*MC
661:    PetscAbsInt - Returns the absolute value of an integer

663:    Synopsis:
664: #include <petscmath.h>
665:    int abs PetscAbsInt(int v1)

667:    Not Collective

669:    Input Parameter:
670: .   v1 - the integer

672:    Level: beginner

674: .seealso: PetscMax(), PetscMin(), PetscAbsReal(), PetscSqr()

676: M*/
677: #define PetscAbsInt(a)  (((a)<0)   ? (-(a)) : (a))

679: /*MC
680:    PetscAbsReal - Returns the absolute value of an real number

682:    Synopsis:
683: #include <petscmath.h>
684:    Real abs PetscAbsReal(PetscReal v1)

686:    Not Collective

688:    Input Parameter:
689: .   v1 - the double

691:    Level: beginner

693: .seealso: PetscMax(), PetscMin(), PetscAbsInt(), PetscSqr()

695: M*/
696: #if defined(PETSC_USE_REAL_SINGLE)
697: #define PetscAbsReal(a) fabsf(a)
698: #elif defined(PETSC_USE_REAL_DOUBLE)
699: #define PetscAbsReal(a) fabs(a)
700: #elif defined(PETSC_USE_REAL___FLOAT128)
701: #define PetscAbsReal(a) fabsq(a)
702: #elif defined(PETSC_USE_REAL___FP16)
703: #define PetscAbsReal(a) fabsf(a)
704: #endif

706: /*MC
707:    PetscSqr - Returns the square of a number

709:    Synopsis:
710: #include <petscmath.h>
711:    type sqr PetscSqr(type v1)

713:    Not Collective

715:    Input Parameter:
716: .   v1 - the value

718:    Notes:
719:     type can be integer or floating point value

721:    Level: beginner

723: .seealso: PetscMax(), PetscMin(), PetscAbsInt(), PetscAbsReal()

725: M*/
726: #define PetscSqr(a)     ((a)*(a))

728: /* ----------------------------------------------------------------------------*/

730: #if defined(PETSC_USE_REAL_SINGLE)
731: #define PetscRealConstant(constant) constant##F
732: #elif defined(PETSC_USE_REAL_DOUBLE)
733: #define PetscRealConstant(constant) constant
734: #elif defined(PETSC_USE_REAL___FLOAT128)
735: #define PetscRealConstant(constant) constant##Q
736: #elif defined(PETSC_USE_REAL___FP16)
737: #define PetscRealConstant(constant) constant##F
738: #endif

740: /*
741:      Basic constants
742: */
743: #define PETSC_PI    PetscRealConstant(3.1415926535897932384626433832795029)
744: #define PETSC_PHI   PetscRealConstant(1.6180339887498948482045868343656381)
745: #define PETSC_SQRT2 PetscRealConstant(1.4142135623730950488016887242096981)

747: #if !defined(PETSC_USE_64BIT_INDICES)
748: #define PETSC_MAX_INT            2147483647
749: #define PETSC_MIN_INT            (-PETSC_MAX_INT - 1)
750: #else
751: #define PETSC_MAX_INT            9223372036854775807L
752: #define PETSC_MIN_INT            (-PETSC_MAX_INT - 1)
753: #endif
754: #define PETSC_MAX_UINT16         65535

756: #if defined(PETSC_USE_REAL_SINGLE)
757: #  define PETSC_MAX_REAL                3.40282346638528860e+38F
758: #  define PETSC_MIN_REAL                (-PETSC_MAX_REAL)
759: #  define PETSC_MACHINE_EPSILON         1.19209290e-07F
760: #  define PETSC_SQRT_MACHINE_EPSILON    3.45266983e-04F
761: #  define PETSC_SMALL                   1.e-5F
762: #elif defined(PETSC_USE_REAL_DOUBLE)
763: #  define PETSC_MAX_REAL                1.7976931348623157e+308
764: #  define PETSC_MIN_REAL                (-PETSC_MAX_REAL)
765: #  define PETSC_MACHINE_EPSILON         2.2204460492503131e-16
766: #  define PETSC_SQRT_MACHINE_EPSILON    1.490116119384766e-08
767: #  define PETSC_SMALL                   1.e-10
768: #elif defined(PETSC_USE_REAL___FLOAT128)
769: #  define PETSC_MAX_REAL                FLT128_MAX
770: #  define PETSC_MIN_REAL                (-FLT128_MAX)
771: #  define PETSC_MACHINE_EPSILON         FLT128_EPSILON
772: #  define PETSC_SQRT_MACHINE_EPSILON    1.38777878078144567552953958511352539e-17Q
773: #  define PETSC_SMALL                   1.e-20Q
774: #elif defined(PETSC_USE_REAL___FP16)
775: #  define PETSC_MAX_REAL                65504.0F
776: #  define PETSC_MIN_REAL                (-PETSC_MAX_REAL)
777: #  define PETSC_MACHINE_EPSILON         .0009765625F
778: #  define PETSC_SQRT_MACHINE_EPSILON    .03125F
779: #  define PETSC_SMALL                   5.e-3F
780: #endif

782: #define PETSC_INFINITY               (PETSC_MAX_REAL/4)
783: #define PETSC_NINFINITY              (-PETSC_INFINITY)

785: PETSC_EXTERN PetscBool PetscIsInfReal(PetscReal);
786: PETSC_EXTERN PetscBool PetscIsNanReal(PetscReal);
787: PETSC_EXTERN PetscBool PetscIsNormalReal(PetscReal);
788: PETSC_STATIC_INLINE PetscBool PetscIsInfOrNanReal(PetscReal v) {return PetscIsInfReal(v) || PetscIsNanReal(v) ? PETSC_TRUE : PETSC_FALSE;}
789: PETSC_STATIC_INLINE PetscBool PetscIsInfScalar(PetscScalar v) {return PetscIsInfReal(PetscAbsScalar(v));}
790: PETSC_STATIC_INLINE PetscBool PetscIsNanScalar(PetscScalar v) {return PetscIsNanReal(PetscAbsScalar(v));}
791: PETSC_STATIC_INLINE PetscBool PetscIsInfOrNanScalar(PetscScalar v) {return PetscIsInfOrNanReal(PetscAbsScalar(v));}
792: PETSC_STATIC_INLINE PetscBool PetscIsNormalScalar(PetscScalar v) {return PetscIsNormalReal(PetscAbsScalar(v));}

794: PETSC_EXTERN PetscBool PetscIsCloseAtTol(PetscReal,PetscReal,PetscReal,PetscReal);
795: PETSC_EXTERN PetscBool PetscEqualReal(PetscReal,PetscReal);
796: PETSC_EXTERN PetscBool PetscEqualScalar(PetscScalar,PetscScalar);

798: /*
799:     These macros are currently hardwired to match the regular data types, so there is no support for a different
800:     MatScalar from PetscScalar. We left the MatScalar in the source just in case we use it again.
801:  */
802: #define MPIU_MATSCALAR MPIU_SCALAR
803: typedef PetscScalar MatScalar;
804: typedef PetscReal MatReal;

806: struct petsc_mpiu_2scalar {PetscScalar a,b;};
807: PETSC_EXTERN MPI_Datatype MPIU_2SCALAR PetscAttrMPITypeTagLayoutCompatible(struct petsc_mpiu_2scalar);

809: /*
810:    MPI Datatypes for composite reductions:
811:    MPIU_REAL_INT -> struct { PetscReal; PetscInt; }
812:    MPIU_SCALAR_INT -> struct { PetscScalar; PetscInt; }
813: */
814: PETSC_EXTERN MPI_Datatype MPIU_REAL_INT;
815: PETSC_EXTERN MPI_Datatype MPIU_SCALAR_INT;

817: #if defined(PETSC_USE_64BIT_INDICES)
818: struct petsc_mpiu_2int {PetscInt a,b;};
819: PETSC_EXTERN MPI_Datatype MPIU_2INT PetscAttrMPITypeTagLayoutCompatible(struct petsc_mpiu_2int);
820: #else
821: #define MPIU_2INT MPI_2INT
822: #endif
823: PETSC_EXTERN MPI_Datatype MPI_4INT;
824: PETSC_EXTERN MPI_Datatype MPIU_4INT;

826: PETSC_STATIC_INLINE PetscInt PetscPowInt(PetscInt base,PetscInt power)
827: {
828:   PetscInt result = 1;
829:   while (power) {
830:     if (power & 1) result *= base;
831:     power >>= 1;
832:     base *= base;
833:   }
834:   return result;
835: }

837: PETSC_STATIC_INLINE PetscInt64 PetscPowInt64(PetscInt base,PetscInt power)
838: {
839:   PetscInt64 result = 1;
840:   while (power) {
841:     if (power & 1) result *= base;
842:     power >>= 1;
843:     base *= base;
844:   }
845:   return result;
846: }

848: PETSC_STATIC_INLINE PetscReal PetscPowRealInt(PetscReal base,PetscInt power)
849: {
850:   PetscReal result = 1;
851:   if (power < 0) {
852:     power = -power;
853:     base  = ((PetscReal)1)/base;
854:   }
855:   while (power) {
856:     if (power & 1) result *= base;
857:     power >>= 1;
858:     base *= base;
859:   }
860:   return result;
861: }

863: PETSC_STATIC_INLINE PetscScalar PetscPowScalarInt(PetscScalar base,PetscInt power)
864: {
865:   PetscScalar result = (PetscReal)1;
866:   if (power < 0) {
867:     power = -power;
868:     base  = ((PetscReal)1)/base;
869:   }
870:   while (power) {
871:     if (power & 1) result *= base;
872:     power >>= 1;
873:     base *= base;
874:   }
875:   return result;
876: }

878: PETSC_STATIC_INLINE PetscScalar PetscPowScalarReal(PetscScalar base,PetscReal power)
879: {
880:   PetscScalar cpower = power;
881:   return PetscPowScalar(base,cpower);
882: }

884: /*MC
885:     PetscLTE - Performs a less than or equal to on a given constant with a fudge for floating point numbers

887:    Synopsis:
888: #include <petscmath.h>
889:    bool PetscLTE(PetscReal x,constant float)

891:    Not Collective

893:    Input Parameters:
894: +   x - the variable
895: -   b - the constant float it is checking if x is less than or equal to

897:    Notes:
898:      The fudge factor is the value PETSC_SMALL

900:      The constant numerical value is automatically set to the appropriate precision of PETSc so can just be provided as, for example, 3.2

902:      This is used in several examples for setting initial conditions based on coordinate values that are computed with i*h that produces inexact
903:      floating point results.

905:    Level: advanced

907: .seealso: PetscMax(), PetscMin(), PetscAbsInt(), PetscAbsReal(), PetscGTE()

909: M*/
910: #define PetscLTE(x,b)  ((x) <= (PetscRealConstant(b)+PETSC_SMALL))

912: /*MC
913:     PetscGTE - Performs a greater than or equal to on a given constant with a fudge for floating point numbers

915:    Synopsis:
916: #include <petscmath.h>
917:    bool PetscGTE(PetscReal x,constant float)

919:    Not Collective

921:    Input Parameters:
922: +   x - the variable
923: -   b - the constant float it is checking if x is greater than or equal to

925:    Notes:
926:      The fudge factor is the value PETSC_SMALL

928:      The constant numerical value is automatically set to the appropriate precision of PETSc so can just be provided as, for example, 3.2

930:      This is used in several examples for setting initial conditions based on coordinate values that are computed with i*h that produces inexact
931:      floating point results.

933:    Level: advanced

935: .seealso: PetscMax(), PetscMin(), PetscAbsInt(), PetscAbsReal(), PetscLTE()

937: M*/
938: #define PetscGTE(x,b)  ((x) >= (PetscRealConstant(b)-PETSC_SMALL))

940: PETSC_EXTERN PetscErrorCode PetscLinearRegression(PetscInt,const PetscReal[],const PetscReal[],PetscReal*,PetscReal*);
941: #endif