Actual source code: petscmath.h

  1: /* $Id: petscmath.h,v 1.21 2001/04/10 19:37:48 bsmith Exp $ */
  2: /*
  3:    
  4:       PETSc mathematics include file. Defines certain basic mathematical 
  5:     constants and functions for working with single and double precision
  6:     floating point numbers as well as complex and integers.

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

 10: */

 14: #include <math.h>

 16: /*

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

 23: */
 24: #if defined(PETSC_USE_COMPLEX)

 26: #if defined (PETSC_HAVE_STD_COMPLEX)
 27: #include <complex>
 28: #elif defined(PETSC_HAVE_NONSTANDARD_COMPLEX_H)
 29: #include PETSC_HAVE_NONSTANDARD_COMPLEX_H
 30: #else
 31: #include <complex.h>
 32: #endif

 34: extern  MPI_Datatype        MPIU_COMPLEX;
 35: #define MPIU_SCALAR         MPIU_COMPLEX
 36: #if defined(PETSC_USE_MAT_SINGLE)
 37: #define MPIU_MATSCALAR        ??Notdone
 38: #else
 39: #define MPIU_MATSCALAR      MPIU_COMPLEX
 40: #endif

 42: #if defined (PETSC_HAVE_STD_COMPLEX)
 43: #define PetscRealPart(a)        (a).real()
 44: #define PetscImaginaryPart(a)   (a).imag()
 45: #define PetscAbsScalar(a)   std::abs(a)
 46: #define PetscConj(a)        std::conj(a)
 47: #define PetscSqrtScalar(a)  std::sqrt(a)
 48: #define PetscPowScalar(a,b) std::pow(a,b)
 49: #define PetscExpScalar(a)   std::exp(a)
 50: #define PetscSinScalar(a)   std::sin(a)
 51: #define PetscCosScalar(a)   std::cos(a)
 52: #else
 53: #define PetscRealPart(a)        real(a)
 54: #define PetscImaginaryPart(a)   imag(a)
 55: #define PetscAbsScalar(a)   abs(a)
 56: #define PetscConj(a)        conj(a)
 57: #define PetscSqrtScalar(a)  sqrt(a)
 58: #define PetscPowScalar(a,b) pow(a,b)
 59: #define PetscExpScalar(a)   exp(a)
 60: #define PetscSinScalar(a)   sin(a)
 61: #define PetscCosScalar(a)   cos(a)
 62: #endif
 63: /*
 64:   The new complex class for GNU C++ is based on templates and is not backward
 65:   compatible with all previous complex class libraries.
 66: */
 67: #if defined(PETSC_HAVE_STD_COMPLEX)
 68: #define Scalar             std::complex<double>
 69: #elif defined(PETSC_HAVE_TEMPLATED_COMPLEX)
 70: #define Scalar             complex<double>
 71: #else
 72: #define Scalar             complex
 73: #endif

 75: /* Compiling for real numbers only */
 76: #else
 77: #  define MPIU_SCALAR           MPI_DOUBLE
 78: #  if defined(PETSC_USE_MAT_SINGLE)
 79: #    define MPIU_MATSCALAR        MPI_FLOAT
 80: #  else
 81: #    define MPIU_MATSCALAR        MPI_DOUBLE
 82: #  endif
 83: #  define PetscRealPart(a)      (a)
 84: #  define PetscImaginaryPart(a) (a)
 85: #  define PetscAbsScalar(a)     (((a)<0.0)   ? -(a) : (a))
 86: #  define PetscConj(a)          (a)
 87: #  define PetscSqrtScalar(a)    sqrt(a)
 88: #  define PetscPowScalar(a,b)   pow(a,b)
 89: #  define PetscExpScalar(a)     exp(a)
 90: #  define PetscSinScalar(a)     sin(a)
 91: #  define PetscCosScalar(a)     cos(a)

 93: #  if defined(PETSC_USE_SINGLE)
 94: #    define Scalar            float
 95: #  else
 96: #    define Scalar            double
 97: #  endif
 98: #endif

100: /*
101:        Allows compiling PETSc so that matrix values are stored in 
102:    single precision but all other objects still use double
103:    precision. This does not work for complex numbers in that case
104:    it remains double

106:           EXPERIMENTAL! NOT YET COMPLETELY WORKING
107: */
108: #if defined(PETSC_USE_COMPLEX)

110: #define MatScalar Scalar 
111: #define MatReal   double

113: #elif defined(PETSC_USE_MAT_SINGLE)

115: #define MatScalar float
116: #define MatReal   float

118: #else

120: #define MatScalar Scalar
121: #define MatReal   double

123: #endif

125: #if defined(PETSC_USE_SINGLE)
126: #define PetscReal float
127: #else 
128: #define PetscReal double
129: #endif

131: /* --------------------------------------------------------------------------*/

133: /*
134:    Certain objects may be created using either single
135:   or double precision.
136: */
137: typedef enum { SCALAR_DOUBLE,SCALAR_SINGLE } ScalarPrecision;

139: /* PETSC_i is the imaginary number, i */
140: extern  Scalar            PETSC_i;

142: #define PetscMin(a,b)      (((a)<(b)) ? (a) : (b))
143: #define PetscMax(a,b)      (((a)<(b)) ? (b) : (a))
144: #define PetscAbsInt(a)     (((a)<0)   ? -(a) : (a))
145: #define PetscAbsDouble(a)  (((a)<0)   ? -(a) : (a))

147: /* ----------------------------------------------------------------------------*/
148: /*
149:      Basic constants
150: */
151: #define PETSC_PI                 3.14159265358979323846264
152: #define PETSC_DEGREES_TO_RADIANS 0.01745329251994
153: #define PETSC_MAX                1.e300
154: #define PETSC_MIN                -1.e300
155: #define PETSC_MAX_INT            1000000000;
156: #define PETSC_MIN_INT            -1000000000;

158: /* ----------------------------------------------------------------------------*/
159: /*
160:     PetscLogDouble variables are used to contain double precision numbers
161:   that are not used in the numerical computations, but rather in logging,
162:   timing etc.
163: */
164: typedef double PetscLogDouble;
165: /*
166:       Once PETSc is compiling with a ADIC enhanced version of MPI
167:    we will create a new MPI_Datatype for the inactive double variables.
168: */
169: #if defined(AD_DERIV_H)
170: /* extern  MPI_Datatype  MPIU_PETSCLOGDOUBLE; */
171: #else
172: #if !defined(USING_MPIUNI)
173: #define MPIU_PETSCLOGDOUBLE MPI_DOUBLE
174: #endif
175: #endif

177: #define PETSCMAP1_a(a,b)  a ## _ ## b
178: #define PETSCMAP1_b(a,b)  PETSCMAP1_a(a,b)
179: #define PETSCMAP1(a)      PETSCMAP1_b(a,Scalar)
180: #endif