Actual source code: petscmath.h

  1: /*
  2:    
  3:       PETSc mathematics include file. Defines certain basic mathematical 
  4:     constants and functions for working with single and double precision
  5:     floating point numbers as well as complex and integers.

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

  9: */

 13: #include <math.h>

 18: /*

 20:      Defines operations that are different for complex and real numbers;
 21:    note that one cannot really mix the use of complex and real in the same 
 22:    PETSc program. All PETSc objects in one program are built around the object
 23:    PetscScalar which is either always a double or a complex.

 25: */

 27: #define PetscExpPassiveScalar(a) PetscExpScalar()

 29: #if defined(PETSC_USE_COMPLEX)
 30: #if defined(PETSC_CLANGUAGE_CXX)
 31: /*
 32:    C++ support of complex numbers: Original support
 33: */
 34: #include <complex>

 36: #define PetscRealPart(a)        (a).real()
 37: #define PetscImaginaryPart(a)   (a).imag()
 38: #define PetscAbsScalar(a)   std::abs(a)
 39: #define PetscConj(a)        std::conj(a)
 40: #define PetscSqrtScalar(a)  std::sqrt(a)
 41: #define PetscPowScalar(a,b) std::pow(a,b)
 42: #define PetscExpScalar(a)   std::exp(a)
 43: #define PetscSinScalar(a)   std::sin(a)
 44: #define PetscCosScalar(a)   std::cos(a)

 46: typedef std::complex<double> PetscScalar;
 47: #else
 48: #include <complex.h>

 50: /* 
 51:    C support of complex numbers: Warning it needs a 
 52:    C90 compliant compiler to work...
 53:  */

 55: #define PetscRealPart(a)        creal(a)
 56: #define PetscImaginaryPart(a)   cimag(a)
 57: #define PetscAbsScalar(a)   cabs(a)
 58: #define PetscConj(a)        conj(a)
 59: #define PetscSqrtScalar(a)  csqrt(a)
 60: #define PetscPowScalar(a,b) cpow(a,b)
 61: #define PetscExpScalar(a)   cexp(a)
 62: #define PetscSinScalar(a)   csin(a)
 63: #define PetscCosScalar(a)   ccos(a)

 65: typedef double complex PetscScalar;
 66: #endif

 69: #define MPIU_SCALAR         MPIU_COMPLEX
 70: #if defined(PETSC_USE_MAT_SINGLE)
 71: #define MPIU_MATSCALAR        ??Notdone
 72: #else
 73: #define MPIU_MATSCALAR      MPIU_COMPLEX
 74: #endif

 76: /* Compiling for real numbers only */
 77: #else
 78: #  if defined(PETSC_USE_SINGLE)
 79: #    define MPIU_SCALAR           MPI_FLOAT
 80: #  elif defined(PETSC_USE_LONG_DOUBLE)
 81: #    define MPIU_SCALAR           MPI_LONG_DOUBLE
 82: #  elif defined(PETSC_INT)
 83: #    define MPIU_INT              MPI_INT
 84: #  else
 85: #    define MPIU_SCALAR           MPI_DOUBLE
 86: #  endif
 87: #  if defined(PETSC_USE_MAT_SINGLE) || defined(PETSC_USE_SINGLE)
 88: #    define MPIU_MATSCALAR        MPI_FLOAT
 89: #  elif defined(PETSC_USE_LONG_DOUBLE)
 90: #    define MPIU_MATSCALAR        MPI_LONG_DOUBLE
 91: #  elif defined(PETSC_USE_INT)
 92: #    define MPIU_MATSCALAR        MPI_INT
 93: #  else
 94: #    define MPIU_MATSCALAR        MPI_DOUBLE
 95: #  endif
 96: #  define PetscRealPart(a)      (a)
 97: #  define PetscImaginaryPart(a) (0)
 98: #  define PetscAbsScalar(a)     (((a)<0.0)   ? -(a) : (a))
 99: #  define PetscConj(a)          (a)
100: #  define PetscSqrtScalar(a)    sqrt(a)
101: #  define PetscPowScalar(a,b)   pow(a,b)
102: #  define PetscExpScalar(a)     exp(a)
103: #  define PetscSinScalar(a)     sin(a)
104: #  define PetscCosScalar(a)     cos(a)

106: #  if defined(PETSC_USE_SINGLE)
107:   typedef float PetscScalar;
108: #  elif defined(PETSC_USE_LONG_DOUBLE)
109:   typedef long double PetscScalar;
110: #  elif defined(PETSC_USE_INT)
111:   typedef int PetscScalar;
112: #  else
113:   typedef double PetscScalar;
114: #  endif
115: #endif

117: #if defined(PETSC_USE_SINGLE)
118: #  define MPIU_REAL   MPI_FLOAT
119: #elif defined(PETSC_USE_LONG_DOUBLE)
120: #  define MPIU_REAL   MPI_LONG_DOUBLE
121: #elif defined(PETSC_USE_INT)
122: #  define MPIU_REAL   MPI_INT
123: #else
124: #  define MPIU_REAL   MPI_DOUBLE
125: #endif

127: #define PetscSign(a) (((a) >= 0) ? ((a) == 0 ? 0 : 1) : -1)
128: #define PetscAbs(a)  (((a) >= 0) ? (a) : -(a))
129: /*
130:        Allows compiling PETSc so that matrix values are stored in 
131:    single precision but all other objects still use double
132:    precision. This does not work for complex numbers in that case
133:    it remains double

135:           EXPERIMENTAL! NOT YET COMPLETELY WORKING
136: */

138: #if defined(PETSC_USE_MAT_SINGLE)
139: typedef float MatScalar;
140: #else
141: typedef PetscScalar MatScalar;
142: #endif

144: #if defined(PETSC_USE_SINGLE)
145:   typedef float PetscReal;
146: #elif defined(PETSC_USE_LONG_DOUBLE)
147:   typedef long double PetscReal;
148: #elif defined(PETSC_USE_INT)
149:   typedef int PetscReal;
150: #else 
151:   typedef double PetscReal;
152: #endif

154: #if defined(PETSC_USE_COMPLEX)
155: typedef PetscReal MatReal;
156: #elif defined(PETSC_USE_MAT_SINGLE) || defined(PETSC_USE_SINGLE)
157: typedef float MatReal;
158: #else
159: typedef PetscReal MatReal;
160: #endif


163: /* --------------------------------------------------------------------------*/

165: /*
166:    Certain objects may be created using either single
167:   or double precision.
168: */
169: typedef enum { PETSC_SCALAR_DOUBLE,PETSC_SCALAR_SINGLE, PETSC_SCALAR_LONG_DOUBLE } PetscScalarPrecision;

171: /* PETSC_i is the imaginary number, i */

174: /*MC
175:    PetscMin - Returns minimum of two numbers

177:    Input Parameter:
178: +  v1 - first value to find minimum of
179: -  v2 - second value to find minimum of

181:    Synopsis:
182:    type PetscMin(type v1,type v2)

184:    Notes: type can be integer or floating point value

186:    Level: beginner


189: .seealso: PetscMin(), PetscAbsInt(), PetscAbsReal(), PetscSqr()

191: M*/
192: #define PetscMin(a,b)   (((a)<(b)) ?  (a) : (b))

194: /*MC
195:    PetscMax - Returns maxium of two numbers

197:    Input Parameter:
198: +  v1 - first value to find maximum of
199: -  v2 - second value to find maximum of

201:    Synopsis:
202:    type max PetscMax(type v1,type v2)

204:    Notes: type can be integer or floating point value

206:    Level: beginner

208: .seealso: PetscMin(), PetscAbsInt(), PetscAbsReal(), PetscSqr()

210: M*/
211: #define PetscMax(a,b)   (((a)<(b)) ?  (b) : (a))

213: /*MC
214:    PetscAbsInt - Returns the absolute value of an integer

216:    Input Parameter:
217: .   v1 - the integer

219:    Synopsis:
220:    int abs PetscAbsInt(int v1)


223:    Level: beginner

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

227: M*/
228: #define PetscAbsInt(a)  (((a)<0)   ? -(a) : (a))

230: /*MC
231:    PetscAbsReal - Returns the absolute value of an real number

233:    Input Parameter:
234: .   v1 - the double 

236:    Synopsis:
237:    int abs PetscAbsReal(PetscReal v1)


240:    Level: beginner

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

244: M*/
245: #define PetscAbsReal(a) (((a)<0)   ? -(a) : (a))

247: /*MC
248:    PetscSqr - Returns the square of a number

250:    Input Parameter:
251: .   v1 - the value

253:    Synopsis:
254:    type sqr PetscSqr(type v1)

256:    Notes: type can be integer or floating point value

258:    Level: beginner

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

262: M*/
263: #define PetscSqr(a)     ((a)*(a))

265: /* ----------------------------------------------------------------------------*/
266: /*
267:      Basic constants - These should be done much better
268: */
269: #define PETSC_PI                 3.14159265358979323846264
270: #define PETSC_DEGREES_TO_RADIANS 0.01745329251994
271: #define PETSC_MAX_INT            1000000000
272: #define PETSC_MIN_INT            -1000000000

274: #if defined(PETSC_USE_SINGLE)
275: #  define PETSC_MAX                     1.e30
276: #  define PETSC_MIN                    -1.e30
277: #  define PETSC_MACHINE_EPSILON         1.e-7
278: #  define PETSC_SQRT_MACHINE_EPSILON    3.e-4
279: #  define PETSC_SMALL                   1.e-5
280: #elif defined(PETSC_USE_INT)
281: #  define PETSC_MAX                     PETSC_MAX_INT
282: #  define PETSC_MIN                     PETSC_MIN_INT
283: #  define PETSC_MACHINE_EPSILON         1
284: #  define PETSC_SQRT_MACHINE_EPSILON    1
285: #  define PETSC_SMALL                   0
286: #else
287: #  define PETSC_MAX                     1.e300
288: #  define PETSC_MIN                    -1.e300
289: #  define PETSC_MACHINE_EPSILON         1.e-14
290: #  define PETSC_SQRT_MACHINE_EPSILON    1.e-7
291: #  define PETSC_SMALL                   1.e-10
292: #endif

294: EXTERN PetscErrorCode  PetscGlobalMax(PetscReal*,PetscReal*,MPI_Comm);
295: EXTERN PetscErrorCode  PetscGlobalMin(PetscReal*,PetscReal*,MPI_Comm);
296: EXTERN PetscErrorCode  PetscGlobalSum(PetscScalar*,PetscScalar*,MPI_Comm);


299: /* ----------------------------------------------------------------------------*/
300: /*
301:     PetscLogDouble variables are used to contain double precision numbers
302:   that are not used in the numerical computations, but rather in logging,
303:   timing etc.
304: */
305: typedef double PetscLogDouble;
306: #define MPIU_PETSCLOGDOUBLE MPI_DOUBLE

308: #define PassiveReal   PetscReal
309: #define PassiveScalar PetscScalar


313: #endif