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: */
 10: #ifndef PETSCMATH_H
 11: #define PETSCMATH_H

 13: #include <math.h>
 14: #include <petscmacros.h>
 15: #include <petscsystypes.h>

 17: /* SUBMANSEC = Sys */

 19: /*

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

 25: */

 27: /*
 28:     Real number definitions
 29:  */
 30: #if defined(PETSC_USE_REAL_SINGLE)
 31:   #define PetscSqrtReal(a)        sqrtf(a)
 32:   #define PetscCbrtReal(a)        cbrtf(a)
 33:   #define PetscHypotReal(a, b)    hypotf(a, b)
 34:   #define PetscAtan2Real(a, b)    atan2f(a, b)
 35:   #define PetscPowReal(a, b)      powf(a, b)
 36:   #define PetscExpReal(a)         expf(a)
 37:   #define PetscLogReal(a)         logf(a)
 38:   #define PetscLog10Real(a)       log10f(a)
 39:   #define PetscLog2Real(a)        log2f(a)
 40:   #define PetscSinReal(a)         sinf(a)
 41:   #define PetscCosReal(a)         cosf(a)
 42:   #define PetscTanReal(a)         tanf(a)
 43:   #define PetscAsinReal(a)        asinf(a)
 44:   #define PetscAcosReal(a)        acosf(a)
 45:   #define PetscAtanReal(a)        atanf(a)
 46:   #define PetscSinhReal(a)        sinhf(a)
 47:   #define PetscCoshReal(a)        coshf(a)
 48:   #define PetscTanhReal(a)        tanhf(a)
 49:   #define PetscAsinhReal(a)       asinhf(a)
 50:   #define PetscAcoshReal(a)       acoshf(a)
 51:   #define PetscAtanhReal(a)       atanhf(a)
 52:   #define PetscErfReal(a)         erff(a)
 53:   #define PetscCeilReal(a)        ceilf(a)
 54:   #define PetscFloorReal(a)       floorf(a)
 55:   #define PetscFmodReal(a, b)     fmodf(a, b)
 56:   #define PetscCopysignReal(a, b) copysignf(a, b)
 57:   #define PetscTGamma(a)          tgammaf(a)
 58:   #if defined(PETSC_HAVE_LGAMMA_IS_GAMMA)
 59:     #define PetscLGamma(a) gammaf(a)
 60:   #else
 61:     #define PetscLGamma(a) lgammaf(a)
 62:   #endif

 64: #elif defined(PETSC_USE_REAL_DOUBLE)
 65:   #define PetscSqrtReal(a)        sqrt(a)
 66:   #define PetscCbrtReal(a)        cbrt(a)
 67:   #define PetscHypotReal(a, b)    hypot(a, b)
 68:   #define PetscAtan2Real(a, b)    atan2(a, b)
 69:   #define PetscPowReal(a, b)      pow(a, b)
 70:   #define PetscExpReal(a)         exp(a)
 71:   #define PetscLogReal(a)         log(a)
 72:   #define PetscLog10Real(a)       log10(a)
 73:   #define PetscLog2Real(a)        log2(a)
 74:   #define PetscSinReal(a)         sin(a)
 75:   #define PetscCosReal(a)         cos(a)
 76:   #define PetscTanReal(a)         tan(a)
 77:   #define PetscAsinReal(a)        asin(a)
 78:   #define PetscAcosReal(a)        acos(a)
 79:   #define PetscAtanReal(a)        atan(a)
 80:   #define PetscSinhReal(a)        sinh(a)
 81:   #define PetscCoshReal(a)        cosh(a)
 82:   #define PetscTanhReal(a)        tanh(a)
 83:   #define PetscAsinhReal(a)       asinh(a)
 84:   #define PetscAcoshReal(a)       acosh(a)
 85:   #define PetscAtanhReal(a)       atanh(a)
 86:   #define PetscErfReal(a)         erf(a)
 87:   #define PetscCeilReal(a)        ceil(a)
 88:   #define PetscFloorReal(a)       floor(a)
 89:   #define PetscFmodReal(a, b)     fmod(a, b)
 90:   #define PetscCopysignReal(a, b) copysign(a, b)
 91:   #define PetscTGamma(a)          tgamma(a)
 92:   #if defined(PETSC_HAVE_LGAMMA_IS_GAMMA)
 93:     #define PetscLGamma(a) gamma(a)
 94:   #else
 95:     #define PetscLGamma(a) lgamma(a)
 96:   #endif

 98: #elif defined(PETSC_USE_REAL___FLOAT128)
 99:   #define PetscSqrtReal(a)        sqrtq(a)
100:   #define PetscCbrtReal(a)        cbrtq(a)
101:   #define PetscHypotReal(a, b)    hypotq(a, b)
102:   #define PetscAtan2Real(a, b)    atan2q(a, b)
103:   #define PetscPowReal(a, b)      powq(a, b)
104:   #define PetscExpReal(a)         expq(a)
105:   #define PetscLogReal(a)         logq(a)
106:   #define PetscLog10Real(a)       log10q(a)
107:   #define PetscLog2Real(a)        log2q(a)
108:   #define PetscSinReal(a)         sinq(a)
109:   #define PetscCosReal(a)         cosq(a)
110:   #define PetscTanReal(a)         tanq(a)
111:   #define PetscAsinReal(a)        asinq(a)
112:   #define PetscAcosReal(a)        acosq(a)
113:   #define PetscAtanReal(a)        atanq(a)
114:   #define PetscSinhReal(a)        sinhq(a)
115:   #define PetscCoshReal(a)        coshq(a)
116:   #define PetscTanhReal(a)        tanhq(a)
117:   #define PetscAsinhReal(a)       asinhq(a)
118:   #define PetscAcoshReal(a)       acoshq(a)
119:   #define PetscAtanhReal(a)       atanhq(a)
120:   #define PetscErfReal(a)         erfq(a)
121:   #define PetscCeilReal(a)        ceilq(a)
122:   #define PetscFloorReal(a)       floorq(a)
123:   #define PetscFmodReal(a, b)     fmodq(a, b)
124:   #define PetscCopysignReal(a, b) copysignq(a, b)
125:   #define PetscTGamma(a)          tgammaq(a)
126:   #if defined(PETSC_HAVE_LGAMMA_IS_GAMMA)
127:     #define PetscLGamma(a) gammaq(a)
128:   #else
129:     #define PetscLGamma(a) lgammaq(a)
130:   #endif

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

166: #endif /* PETSC_USE_REAL_* */

168: static inline PetscReal PetscSignReal(PetscReal a)
169: {
170:   return (PetscReal)((a < (PetscReal)0) ? -1 : ((a > (PetscReal)0) ? 1 : 0));
171: }

173: #if !defined(PETSC_HAVE_LOG2)
174:   #undef PetscLog2Real
175: static inline PetscReal PetscLog2Real(PetscReal a)
176: {
177:   return PetscLogReal(a) / PetscLogReal((PetscReal)2);
178: }
179: #endif

181: #if defined(PETSC_HAVE_REAL___FLOAT128) && !defined(PETSC_SKIP_REAL___FLOAT128)
182: PETSC_EXTERN MPI_Datatype MPIU___FLOAT128 PETSC_ATTRIBUTE_MPI_TYPE_TAG(__float128);
183: #endif
184: #if defined(PETSC_HAVE_REAL___FP16) && !defined(PETSC_SKIP_REAL___FP16)
185: PETSC_EXTERN MPI_Datatype MPIU___FP16 PETSC_ATTRIBUTE_MPI_TYPE_TAG(__fp16);
186: #endif

188: /*MC
189:    MPIU_REAL - Portable MPI datatype corresponding to `PetscReal` independent of what precision `PetscReal` is in

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

194:    Level: beginner

196: .seealso: `PetscReal`, `PetscScalar`, `PetscComplex`, `PetscInt`, `MPIU_SCALAR`, `MPIU_COMPLEX`, `MPIU_INT`
197: M*/
198: #if defined(PETSC_USE_REAL_SINGLE)
199:   #define MPIU_REAL MPI_FLOAT
200: #elif defined(PETSC_USE_REAL_DOUBLE)
201:   #define MPIU_REAL MPI_DOUBLE
202: #elif defined(PETSC_USE_REAL___FLOAT128)
203:   #define MPIU_REAL MPIU___FLOAT128
204: #elif defined(PETSC_USE_REAL___FP16)
205:   #define MPIU_REAL MPIU___FP16
206: #endif /* PETSC_USE_REAL_* */

208: /*
209:     Complex number definitions
210:  */
211: #if defined(PETSC_HAVE_COMPLEX)
212:   #if defined(__cplusplus) && !defined(PETSC_USE_REAL___FLOAT128)
213:     /* C++ support of complex number */

215:     #define PetscRealPartComplex(a)      (static_cast<PetscComplex>(a)).real()
216:     #define PetscImaginaryPartComplex(a) (static_cast<PetscComplex>(a)).imag()
217:     #define PetscAbsComplex(a)           petsccomplexlib::abs(static_cast<PetscComplex>(a))
218:     #define PetscArgComplex(a)           petsccomplexlib::arg(static_cast<PetscComplex>(a))
219:     #define PetscConjComplex(a)          petsccomplexlib::conj(static_cast<PetscComplex>(a))
220:     #define PetscSqrtComplex(a)          petsccomplexlib::sqrt(static_cast<PetscComplex>(a))
221:     #define PetscPowComplex(a, b)        petsccomplexlib::pow(static_cast<PetscComplex>(a), static_cast<PetscComplex>(b))
222:     #define PetscExpComplex(a)           petsccomplexlib::exp(static_cast<PetscComplex>(a))
223:     #define PetscLogComplex(a)           petsccomplexlib::log(static_cast<PetscComplex>(a))
224:     #define PetscSinComplex(a)           petsccomplexlib::sin(static_cast<PetscComplex>(a))
225:     #define PetscCosComplex(a)           petsccomplexlib::cos(static_cast<PetscComplex>(a))
226:     #define PetscTanComplex(a)           petsccomplexlib::tan(static_cast<PetscComplex>(a))
227:     #define PetscAsinComplex(a)          petsccomplexlib::asin(static_cast<PetscComplex>(a))
228:     #define PetscAcosComplex(a)          petsccomplexlib::acos(static_cast<PetscComplex>(a))
229:     #define PetscAtanComplex(a)          petsccomplexlib::atan(static_cast<PetscComplex>(a))
230:     #define PetscSinhComplex(a)          petsccomplexlib::sinh(static_cast<PetscComplex>(a))
231:     #define PetscCoshComplex(a)          petsccomplexlib::cosh(static_cast<PetscComplex>(a))
232:     #define PetscTanhComplex(a)          petsccomplexlib::tanh(static_cast<PetscComplex>(a))
233:     #define PetscAsinhComplex(a)         petsccomplexlib::asinh(static_cast<PetscComplex>(a))
234:     #define PetscAcoshComplex(a)         petsccomplexlib::acosh(static_cast<PetscComplex>(a))
235:     #define PetscAtanhComplex(a)         petsccomplexlib::atanh(static_cast<PetscComplex>(a))

237:     /* TODO: Add configure tests

239: #if !defined(PETSC_HAVE_CXX_TAN_COMPLEX)
240: #undef PetscTanComplex
241: static inline PetscComplex PetscTanComplex(PetscComplex z)
242: {
243:   return PetscSinComplex(z)/PetscCosComplex(z);
244: }
245: #endif

247: #if !defined(PETSC_HAVE_CXX_TANH_COMPLEX)
248: #undef PetscTanhComplex
249: static inline PetscComplex PetscTanhComplex(PetscComplex z)
250: {
251:   return PetscSinhComplex(z)/PetscCoshComplex(z);
252: }
253: #endif

255: #if !defined(PETSC_HAVE_CXX_ASIN_COMPLEX)
256: #undef PetscAsinComplex
257: static inline PetscComplex PetscAsinComplex(PetscComplex z)
258: {
259:   const PetscComplex j(0,1);
260:   return -j*PetscLogComplex(j*z+PetscSqrtComplex(1.0f-z*z));
261: }
262: #endif

264: #if !defined(PETSC_HAVE_CXX_ACOS_COMPLEX)
265: #undef PetscAcosComplex
266: static inline PetscComplex PetscAcosComplex(PetscComplex z)
267: {
268:   const PetscComplex j(0,1);
269:   return j*PetscLogComplex(z-j*PetscSqrtComplex(1.0f-z*z));
270: }
271: #endif

273: #if !defined(PETSC_HAVE_CXX_ATAN_COMPLEX)
274: #undef PetscAtanComplex
275: static inline PetscComplex PetscAtanComplex(PetscComplex z)
276: {
277:   const PetscComplex j(0,1);
278:   return 0.5f*j*PetscLogComplex((1.0f-j*z)/(1.0f+j*z));
279: }
280: #endif

282: #if !defined(PETSC_HAVE_CXX_ASINH_COMPLEX)
283: #undef PetscAsinhComplex
284: static inline PetscComplex PetscAsinhComplex(PetscComplex z)
285: {
286:   return PetscLogComplex(z+PetscSqrtComplex(z*z+1.0f));
287: }
288: #endif

290: #if !defined(PETSC_HAVE_CXX_ACOSH_COMPLEX)
291: #undef PetscAcoshComplex
292: static inline PetscComplex PetscAcoshComplex(PetscComplex z)
293: {
294:   return PetscLogComplex(z+PetscSqrtComplex(z*z-1.0f));
295: }
296: #endif

298: #if !defined(PETSC_HAVE_CXX_ATANH_COMPLEX)
299: #undef PetscAtanhComplex
300: static inline PetscComplex PetscAtanhComplex(PetscComplex z)
301: {
302:   return 0.5f*PetscLogComplex((1.0f+z)/(1.0f-z));
303: }
304: #endif

306: */

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

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

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

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

379:     #endif /* PETSC_USE_REAL_* */
380:   #endif   /* (__cplusplus) */

382: /*
383:    PETSC_i is the imaginary number, i
384: */
385: PETSC_EXTERN PetscComplex PETSC_i;

387: /*
388:    Try to do the right thing for complex number construction: see
389:    http://www.open-std.org/jtc1/sc22/wg14/www/docs/n1464.htm
390:    for details
391: */
392: static inline PetscComplex PetscCMPLX(PetscReal x, PetscReal y)
393: {
394:   #if defined(__cplusplus) && !defined(PETSC_USE_REAL___FLOAT128)
395:   return PetscComplex(x, y);
396:   #elif defined(_Imaginary_I)
397:   return x + y * _Imaginary_I;
398:   #else
399:   { /* In both C99 and C11 (ISO/IEC 9899, Section 6.2.5),

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

405:        So type punning should be portable. */
406:     union
407:     {
408:       PetscComplex z;
409:       PetscReal    f[2];
410:     } uz;

412:     uz.f[0] = x;
413:     uz.f[1] = y;
414:     return uz.z;
415:   }
416:   #endif
417: }

419:   #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)\"")
420:   #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)\"")

422:   #if defined(PETSC_HAVE_REAL___FLOAT128) && !defined(PETSC_SKIP_REAL___FLOAT128)
423:     // if complex is not used, then quadmath.h won't be included by petscsystypes.h
424:     #if defined(PETSC_USE_COMPLEX)
425:       #define MPIU___COMPLEX128_ATTR_TAG PETSC_ATTRIBUTE_MPI_TYPE_TAG(__complex128)
426:     #else
427:       #define MPIU___COMPLEX128_ATTR_TAG
428:     #endif

430: PETSC_EXTERN MPI_Datatype MPIU___COMPLEX128 MPIU___COMPLEX128_ATTR_TAG;

432:     #undef MPIU___COMPLEX128_ATTR_TAG
433:   #endif /* PETSC_HAVE_REAL___FLOAT128 */

435:   /*MC
436:    MPIU_COMPLEX - Portable MPI datatype corresponding to `PetscComplex` independent of the precision of `PetscComplex`

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

441:    Level: beginner

443: .seealso: `PetscReal`, `PetscScalar`, `PetscComplex`, `PetscInt`, `MPIU_REAL`, `MPIU_SCALAR`, `MPIU_COMPLEX`, `MPIU_INT`, `PETSC_i`
444: M*/
445:   #if defined(PETSC_USE_REAL_SINGLE)
446:     #define MPIU_COMPLEX MPI_C_COMPLEX
447:   #elif defined(PETSC_USE_REAL_DOUBLE)
448:     #define MPIU_COMPLEX MPI_C_DOUBLE_COMPLEX
449:   #elif defined(PETSC_USE_REAL___FLOAT128)
450:     #define MPIU_COMPLEX MPIU___COMPLEX128
451:   #elif defined(PETSC_USE_REAL___FP16)
452:     #define MPIU_COMPLEX MPI_C_COMPLEX
453:   #endif /* PETSC_USE_REAL_* */

455: #endif /* PETSC_HAVE_COMPLEX */

457: /*
458:     Scalar number definitions
459:  */
460: #if defined(PETSC_USE_COMPLEX) && defined(PETSC_HAVE_COMPLEX)
461:   /*MC
462:    MPIU_SCALAR - Portable MPI datatype corresponding to `PetscScalar` independent of the precision of `PetscScalar`

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

467:    Level: beginner

469: .seealso: `PetscReal`, `PetscScalar`, `PetscComplex`, `PetscInt`, `MPIU_REAL`, `MPIU_COMPLEX`, `MPIU_INT`
470: M*/
471:   #define MPIU_SCALAR MPIU_COMPLEX

473:   /*MC
474:    PetscRealPart - Returns the real part of a `PetscScalar`

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

480:    Not Collective

482:    Input Parameter:
483: .  v - value to find the real part of

485:    Level: beginner

487: .seealso: `PetscScalar`, `PetscImaginaryPart()`, `PetscMax()`, `PetscClipInterval()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscSqr()`
488: M*/
489:   #define PetscRealPart(a) PetscRealPartComplex(a)

491:   /*MC
492:    PetscImaginaryPart - Returns the imaginary part of a `PetscScalar`

494:    Synopsis:
495: #include <petscmath.h>
496:    PetscReal PetscImaginaryPart(PetscScalar v)

498:    Not Collective

500:    Input Parameter:
501: .  v - value to find the imaginary part of

503:    Level: beginner

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

508: .seealso: `PetscScalar`, `PetscRealPart()`, `PetscMax()`, `PetscClipInterval()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscSqr()`
509: M*/
510:   #define PetscImaginaryPart(a) PetscImaginaryPartComplex(a)

512:   #define PetscAbsScalar(a)    PetscAbsComplex(a)
513:   #define PetscArgScalar(a)    PetscArgComplex(a)
514:   #define PetscConj(a)         PetscConjComplex(a)
515:   #define PetscSqrtScalar(a)   PetscSqrtComplex(a)
516:   #define PetscPowScalar(a, b) PetscPowComplex(a, b)
517:   #define PetscExpScalar(a)    PetscExpComplex(a)
518:   #define PetscLogScalar(a)    PetscLogComplex(a)
519:   #define PetscSinScalar(a)    PetscSinComplex(a)
520:   #define PetscCosScalar(a)    PetscCosComplex(a)
521:   #define PetscTanScalar(a)    PetscTanComplex(a)
522:   #define PetscAsinScalar(a)   PetscAsinComplex(a)
523:   #define PetscAcosScalar(a)   PetscAcosComplex(a)
524:   #define PetscAtanScalar(a)   PetscAtanComplex(a)
525:   #define PetscSinhScalar(a)   PetscSinhComplex(a)
526:   #define PetscCoshScalar(a)   PetscCoshComplex(a)
527:   #define PetscTanhScalar(a)   PetscTanhComplex(a)
528:   #define PetscAsinhScalar(a)  PetscAsinhComplex(a)
529:   #define PetscAcoshScalar(a)  PetscAcoshComplex(a)
530:   #define PetscAtanhScalar(a)  PetscAtanhComplex(a)

532: #else /* PETSC_USE_COMPLEX */
533:   #define MPIU_SCALAR           MPIU_REAL
534:   #define PetscRealPart(a)      (a)
535:   #define PetscImaginaryPart(a) ((PetscReal)0)
536:   #define PetscAbsScalar(a)     PetscAbsReal(a)
537:   #define PetscArgScalar(a)     (((a) < (PetscReal)0) ? PETSC_PI : (PetscReal)0)
538:   #define PetscConj(a)          (a)
539:   #define PetscSqrtScalar(a)    PetscSqrtReal(a)
540:   #define PetscPowScalar(a, b)  PetscPowReal(a, b)
541:   #define PetscExpScalar(a)     PetscExpReal(a)
542:   #define PetscLogScalar(a)     PetscLogReal(a)
543:   #define PetscSinScalar(a)     PetscSinReal(a)
544:   #define PetscCosScalar(a)     PetscCosReal(a)
545:   #define PetscTanScalar(a)     PetscTanReal(a)
546:   #define PetscAsinScalar(a)    PetscAsinReal(a)
547:   #define PetscAcosScalar(a)    PetscAcosReal(a)
548:   #define PetscAtanScalar(a)    PetscAtanReal(a)
549:   #define PetscSinhScalar(a)    PetscSinhReal(a)
550:   #define PetscCoshScalar(a)    PetscCoshReal(a)
551:   #define PetscTanhScalar(a)    PetscTanhReal(a)
552:   #define PetscAsinhScalar(a)   PetscAsinhReal(a)
553:   #define PetscAcoshScalar(a)   PetscAcoshReal(a)
554:   #define PetscAtanhScalar(a)   PetscAtanhReal(a)

556: #endif /* PETSC_USE_COMPLEX */

558: /*
559:    Certain objects may be created using either single or double precision.
560:    This is currently not used.
561: */
562: typedef enum {
563:   PETSC_SCALAR_DOUBLE,
564:   PETSC_SCALAR_SINGLE,
565:   PETSC_SCALAR_LONG_DOUBLE,
566:   PETSC_SCALAR_HALF
567: } PetscScalarPrecision;

569: /*MC
570:    PetscAbs - Returns the absolute value of a number

572:    Synopsis:
573: #include <petscmath.h>
574:    type PetscAbs(type v)

576:    Not Collective

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

581:    Level: beginner

583:    Note:
584:    The type can be integer or real floating point value, but cannot be complex

586: .seealso: `PetscAbsInt()`, `PetscAbsReal()`, `PetscAbsScalar()`, `PetscSign()`
587: M*/
588: #define PetscAbs(a) (((a) >= 0) ? (a) : (-(a)))

590: /*MC
591:    PetscSign - Returns the sign of a number as an integer

593:    Synopsis:
594: #include <petscmath.h>
595:    int PetscSign(type v)

597:    Not Collective

599:    Input Parameter:
600: .  v - the number

602:    Level: beginner

604:    Note:
605:    The type can be integer or real floating point value

607: .seealso: `PetscAbsInt()`, `PetscAbsReal()`, `PetscAbsScalar()`
608: M*/
609: #define PetscSign(a) (((a) >= 0) ? ((a) == 0 ? 0 : 1) : -1)

611: /*MC
612:    PetscMin - Returns minimum of two numbers

614:    Synopsis:
615: #include <petscmath.h>
616:    type PetscMin(type v1,type v2)

618:    Not Collective

620:    Input Parameters:
621: +  v1 - first value to find minimum of
622: -  v2 - second value to find minimum of

624:    Level: beginner

626:    Note:
627:    The type can be integer or floating point value

629: .seealso: `PetscMax()`, `PetscClipInterval()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscSqr()`
630: M*/
631: #define PetscMin(a, b) (((a) < (b)) ? (a) : (b))

633: /*MC
634:    PetscMax - Returns maximum of two numbers

636:    Synopsis:
637: #include <petscmath.h>
638:    type max PetscMax(type v1,type v2)

640:    Not Collective

642:    Input Parameters:
643: +  v1 - first value to find maximum of
644: -  v2 - second value to find maximum of

646:    Level: beginner

648:    Note:
649:    The type can be integer or floating point value

651: .seealso: `PetscMin()`, `PetscClipInterval()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscSqr()`
652: M*/
653: #define PetscMax(a, b) (((a) < (b)) ? (b) : (a))

655: /*MC
656:    PetscClipInterval - Returns a number clipped to be within an interval

658:    Synopsis:
659: #include <petscmath.h>
660:    type clip PetscClipInterval(type x,type a,type b)

662:    Not Collective

664:    Input Parameters:
665: +  x - value to use if within interval [a,b]
666: .  a - lower end of interval
667: -  b - upper end of interval

669:    Level: beginner

671:    Note:
672:    The type can be integer or floating point value

674: .seealso: `PetscMin()`, `PetscMax()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscSqr()`
675: M*/
676: #define PetscClipInterval(x, a, b) (PetscMax((a), PetscMin((x), (b))))

678: /*MC
679:    PetscAbsInt - Returns the absolute value of an integer

681:    Synopsis:
682: #include <petscmath.h>
683:    int abs PetscAbsInt(int v1)

685:    Input Parameter:
686: .   v1 - the integer

688:    Level: beginner

690: .seealso: `PetscMax()`, `PetscMin()`, `PetscAbsReal()`, `PetscSqr()`
691: M*/
692: #define PetscAbsInt(a) (((a) < 0) ? (-(a)) : (a))

694: /*MC
695:    PetscAbsReal - Returns the absolute value of an real number

697:    Synopsis:
698: #include <petscmath.h>
699:    Real abs PetscAbsReal(PetscReal v1)

701:    Input Parameter:
702: .   v1 - the `PetscReal` value

704:    Level: beginner

706: .seealso: `PetscMax()`, `PetscMin()`, `PetscAbsInt()`, `PetscSqr()`
707: M*/
708: #if defined(PETSC_USE_REAL_SINGLE)
709:   #define PetscAbsReal(a) fabsf(a)
710: #elif defined(PETSC_USE_REAL_DOUBLE)
711:   #define PetscAbsReal(a) fabs(a)
712: #elif defined(PETSC_USE_REAL___FLOAT128)
713:   #define PetscAbsReal(a) fabsq(a)
714: #elif defined(PETSC_USE_REAL___FP16)
715:   #define PetscAbsReal(a) fabsf(a)
716: #endif

718: /*MC
719:    PetscSqr - Returns the square of a number

721:    Synopsis:
722: #include <petscmath.h>
723:    type sqr PetscSqr(type v1)

725:    Not Collective

727:    Input Parameter:
728: .   v1 - the value

730:    Level: beginner

732:    Note:
733:    The type can be integer or floating point value

735: .seealso: `PetscMax()`, `PetscMin()`, `PetscAbsInt()`, `PetscAbsReal()`
736: M*/
737: #define PetscSqr(a) ((a) * (a))

739: #if defined(PETSC_USE_REAL_SINGLE)
740:   #define PetscRealConstant(constant) constant##F
741: #elif defined(PETSC_USE_REAL_DOUBLE)
742:   #define PetscRealConstant(constant) constant
743: #elif defined(PETSC_USE_REAL___FLOAT128)
744:   #define PetscRealConstant(constant) constant##Q
745: #elif defined(PETSC_USE_REAL___FP16)
746:   #define PetscRealConstant(constant) constant##F
747: #endif

749: /*
750:      Basic constants
751: */
752: #define PETSC_PI    PetscRealConstant(3.1415926535897932384626433832795029)
753: #define PETSC_PHI   PetscRealConstant(1.6180339887498948482045868343656381)
754: #define PETSC_SQRT2 PetscRealConstant(1.4142135623730950488016887242096981)

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_REAL_MIN             1.1754944e-38F
760:   #define PETSC_MACHINE_EPSILON      1.19209290e-07F
761:   #define PETSC_SQRT_MACHINE_EPSILON 3.45266983e-04F
762:   #define PETSC_SMALL                1.e-5F
763: #elif defined(PETSC_USE_REAL_DOUBLE)
764:   #define PETSC_MAX_REAL             1.7976931348623157e+308
765:   #define PETSC_MIN_REAL             (-PETSC_MAX_REAL)
766:   #define PETSC_REAL_MIN             2.225073858507201e-308
767:   #define PETSC_MACHINE_EPSILON      2.2204460492503131e-16
768:   #define PETSC_SQRT_MACHINE_EPSILON 1.490116119384766e-08
769:   #define PETSC_SMALL                1.e-10
770: #elif defined(PETSC_USE_REAL___FLOAT128)
771:   #define PETSC_MAX_REAL             FLT128_MAX
772:   #define PETSC_MIN_REAL             (-FLT128_MAX)
773:   #define PETSC_REAL_MIN             FLT128_MIN
774:   #define PETSC_MACHINE_EPSILON      FLT128_EPSILON
775:   #define PETSC_SQRT_MACHINE_EPSILON 1.38777878078144567552953958511352539e-17Q
776:   #define PETSC_SMALL                1.e-20Q
777: #elif defined(PETSC_USE_REAL___FP16)
778:   #define PETSC_MAX_REAL             65504.0F
779:   #define PETSC_MIN_REAL             (-PETSC_MAX_REAL)
780:   #define PETSC_REAL_MIN             .00006103515625F
781:   #define PETSC_MACHINE_EPSILON      .0009765625F
782:   #define PETSC_SQRT_MACHINE_EPSILON .03125F
783:   #define PETSC_SMALL                5.e-3F
784: #endif

786: #define PETSC_INFINITY  (PETSC_MAX_REAL / 4)
787: #define PETSC_NINFINITY (-PETSC_INFINITY)

789: PETSC_EXTERN PetscBool  PetscIsInfReal(PetscReal);
790: PETSC_EXTERN PetscBool  PetscIsNanReal(PetscReal);
791: PETSC_EXTERN PetscBool  PetscIsNormalReal(PetscReal);
792: static inline PetscBool PetscIsInfOrNanReal(PetscReal v)
793: {
794:   return PetscIsInfReal(v) || PetscIsNanReal(v) ? PETSC_TRUE : PETSC_FALSE;
795: }
796: static inline PetscBool PetscIsInfScalar(PetscScalar v)
797: {
798:   return PetscIsInfReal(PetscAbsScalar(v));
799: }
800: static inline PetscBool PetscIsNanScalar(PetscScalar v)
801: {
802:   return PetscIsNanReal(PetscAbsScalar(v));
803: }
804: static inline PetscBool PetscIsInfOrNanScalar(PetscScalar v)
805: {
806:   return PetscIsInfOrNanReal(PetscAbsScalar(v));
807: }
808: static inline PetscBool PetscIsNormalScalar(PetscScalar v)
809: {
810:   return PetscIsNormalReal(PetscAbsScalar(v));
811: }

813: PETSC_EXTERN PetscBool PetscIsCloseAtTol(PetscReal, PetscReal, PetscReal, PetscReal);
814: PETSC_EXTERN PetscBool PetscEqualReal(PetscReal, PetscReal);
815: PETSC_EXTERN PetscBool PetscEqualScalar(PetscScalar, PetscScalar);

817: /*
818:     These macros are currently hardwired to match the regular data types, so there is no support for a different
819:     MatScalar from PetscScalar. We left the MatScalar in the source just in case we use it again.
820:  */
821: #define MPIU_MATSCALAR MPIU_SCALAR
822: typedef PetscScalar MatScalar;
823: typedef PetscReal   MatReal;

825: struct petsc_mpiu_2scalar {
826:   PetscScalar a, b;
827: };
828: PETSC_EXTERN MPI_Datatype MPIU_2SCALAR PETSC_ATTRIBUTE_MPI_TYPE_TAG_LAYOUT_COMPATIBLE(struct petsc_mpiu_2scalar);

830: /* MPI Datatypes for composite reductions */
831: struct petsc_mpiu_real_int {
832:   PetscReal v;
833:   PetscInt  i;
834: };

836: struct petsc_mpiu_scalar_int {
837:   PetscScalar v;
838:   PetscInt    i;
839: };

841: PETSC_EXTERN MPI_Datatype MPIU_REAL_INT PETSC_ATTRIBUTE_MPI_TYPE_TAG_LAYOUT_COMPATIBLE(struct petsc_mpiu_real_int);
842: PETSC_EXTERN MPI_Datatype MPIU_SCALAR_INT PETSC_ATTRIBUTE_MPI_TYPE_TAG_LAYOUT_COMPATIBLE(struct petsc_mpiu_scalar_int);

844: #if defined(PETSC_USE_64BIT_INDICES)
845: struct /* __attribute__((packed, aligned(alignof(PetscInt *)))) */ petsc_mpiu_2int {
846:   PetscInt a;
847:   PetscInt b;
848: };
849: /*
850:  static_assert(sizeof(struct petsc_mpiu_2int) == 2 * sizeof(PetscInt), "");
851:  static_assert(alignof(struct petsc_mpiu_2int) == alignof(PetscInt *), "");
852:  static_assert(alignof(struct petsc_mpiu_2int) == alignof(PetscInt[2]), "");

854:  clang generates warnings that petsc_mpiu_2int is not layout compatible with PetscInt[2] or
855:  PetscInt *, even though (with everything else uncommented) both of the static_asserts above
856:  pass! So we just comment it out...
857: */
858: PETSC_EXTERN MPI_Datatype MPIU_2INT /* PETSC_ATTRIBUTE_MPI_TYPE_TAG_LAYOUT_COMPATIBLE(struct petsc_mpiu_2int) */;
859: #else
860:   #define MPIU_2INT MPI_2INT
861: #endif
862: PETSC_EXTERN MPI_Datatype MPI_4INT;
863: PETSC_EXTERN MPI_Datatype MPIU_4INT;

865: static inline PetscInt PetscPowInt(PetscInt base, PetscInt power)
866: {
867:   PetscInt result = 1;
868:   while (power) {
869:     if (power & 1) result *= base;
870:     power >>= 1;
871:     if (power) base *= base;
872:   }
873:   return result;
874: }

876: static inline PetscInt64 PetscPowInt64(PetscInt base, PetscInt power)
877: {
878:   PetscInt64 result = 1;
879:   while (power) {
880:     if (power & 1) result *= base;
881:     power >>= 1;
882:     if (power) base *= base;
883:   }
884:   return result;
885: }

887: static inline PetscReal PetscPowRealInt(PetscReal base, PetscInt power)
888: {
889:   PetscReal result = 1;
890:   if (power < 0) {
891:     power = -power;
892:     base  = ((PetscReal)1) / base;
893:   }
894:   while (power) {
895:     if (power & 1) result *= base;
896:     power >>= 1;
897:     if (power) base *= base;
898:   }
899:   return result;
900: }

902: static inline PetscScalar PetscPowScalarInt(PetscScalar base, PetscInt power)
903: {
904:   PetscScalar result = (PetscReal)1;
905:   if (power < 0) {
906:     power = -power;
907:     base  = ((PetscReal)1) / base;
908:   }
909:   while (power) {
910:     if (power & 1) result *= base;
911:     power >>= 1;
912:     if (power) base *= base;
913:   }
914:   return result;
915: }

917: static inline PetscScalar PetscPowScalarReal(PetscScalar base, PetscReal power)
918: {
919:   PetscScalar cpower = power;
920:   return PetscPowScalar(base, cpower);
921: }

923: /*MC
924:     PetscApproximateLTE - Performs a less than or equal to on a given constant with a fudge for floating point numbers

926:    Synopsis:
927: #include <petscmath.h>
928:    bool PetscApproximateLTE(PetscReal x,constant float)

930:    Not Collective

932:    Input Parameters:
933: +   x - the variable
934: -   b - the constant float it is checking if `x` is less than or equal to

936:    Level: advanced

938:    Notes:
939:      The fudge factor is the value `PETSC_SMALL`

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

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

946: .seealso: `PetscMax()`, `PetscMin()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscApproximateGTE()`
947: M*/
948: #define PetscApproximateLTE(x, b) ((x) <= (PetscRealConstant(b) + PETSC_SMALL))

950: /*MC
951:     PetscApproximateGTE - Performs a greater than or equal to on a given constant with a fudge for floating point numbers

953:    Synopsis:
954: #include <petscmath.h>
955:    bool PetscApproximateGTE(PetscReal x,constant float)

957:    Not Collective

959:    Input Parameters:
960: +   x - the variable
961: -   b - the constant float it is checking if `x` is greater than or equal to

963:    Level: advanced

965:    Notes:
966:      The fudge factor is the value `PETSC_SMALL`

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

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

973: .seealso: `PetscMax()`, `PetscMin()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscApproximateLTE()`
974: M*/
975: #define PetscApproximateGTE(x, b) ((x) >= (PetscRealConstant(b) - PETSC_SMALL))

977: /*MC
978:     PetscCeilInt - Returns the ceiling of the quotation of two positive integers

980:    Synopsis:
981: #include <petscmath.h>
982:    PetscInt PetscCeilInt(PetscInt x,PetscInt y)

984:    Not Collective

986:    Input Parameters:
987: +   x - the numerator
988: -   y - the denominator

990:    Level: advanced

992: .seealso: `PetscMax()`, `PetscMin()`, `PetscAbsInt()`, `PetscAbsReal()`, `PetscApproximateLTE()`
993: M*/
994: #define PetscCeilInt(x, y) ((((PetscInt)(x)) / ((PetscInt)(y))) + ((((PetscInt)(x)) % ((PetscInt)(y))) ? 1 : 0))

996: #define PetscCeilInt64(x, y) ((((PetscInt64)(x)) / ((PetscInt64)(y))) + ((((PetscInt64)(x)) % ((PetscInt64)(y))) ? 1 : 0))

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