Actual source code: petscmath.h
1: /* $Id: petscmath.h,v 1.32 2001/08/30 20:37:06 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: extern MPI_Datatype MPIU_2SCALAR;
17: /*
19: Defines operations that are different for complex and real numbers;
20: note that one cannot really mix the use of complex and real in the same
21: PETSc program. All PETSc objects in one program are built around the object
22: PetscScalar which is either always a double or a complex.
24: */
26: #if defined(PETSC_USE_COMPLEX)
28: /*
29: PETSc now only supports std::complex
30: */
31: #include <complex>
33: extern MPI_Datatype MPIU_COMPLEX;
34: #define MPIU_SCALAR MPIU_COMPLEX
35: #if defined(PETSC_USE_MAT_SINGLE)
36: #define MPIU_MATSCALAR ??Notdone
37: #else
38: #define MPIU_MATSCALAR MPIU_COMPLEX
39: #endif
41: #define PetscRealPart(a) (a).real()
42: #define PetscImaginaryPart(a) (a).imag()
43: #define PetscAbsScalar(a) std::abs(a)
44: #define PetscConj(a) std::conj(a)
45: #define PetscSqrtScalar(a) std::sqrt(a)
46: #define PetscPowScalar(a,b) std::pow(a,b)
47: #define PetscExpScalar(a) std::exp(a)
48: #define PetscSinScalar(a) std::sin(a)
49: #define PetscCosScalar(a) std::cos(a)
51: typedef std::complex<double> PetscScalar;
53: /* Compiling for real numbers only */
54: #else
55: # if defined(PETSC_USE_SINGLE)
56: # define MPIU_SCALAR MPI_FLOAT
57: # else
58: # define MPIU_SCALAR MPI_DOUBLE
59: # endif
60: # if defined(PETSC_USE_MAT_SINGLE) || defined(PETSC_USE_SINGLE)
61: # define MPIU_MATSCALAR MPI_FLOAT
62: # else
63: # define MPIU_MATSCALAR MPI_DOUBLE
64: # endif
65: # define PetscRealPart(a) (a)
66: # define PetscImaginaryPart(a) (a)
67: # define PetscAbsScalar(a) (((a)<0.0) ? -(a) : (a))
68: # define PetscConj(a) (a)
69: # define PetscSqrtScalar(a) sqrt(a)
70: # define PetscPowScalar(a,b) pow(a,b)
71: # define PetscExpScalar(a) exp(a)
72: # define PetscSinScalar(a) sin(a)
73: # define PetscCosScalar(a) cos(a)
75: # if defined(PETSC_USE_SINGLE)
76: typedef float PetscScalar;
77: # else
78: typedef double PetscScalar;
79: # endif
80: #endif
82: #if defined(PETSC_USE_SINGLE)
83: # define MPIU_REAL MPI_FLOAT
84: #else
85: # define MPIU_REAL MPI_DOUBLE
86: #endif
88: #define PetscSign(a) (((a) >= 0) ? ((a) == 0 ? 0 : 1) : -1)
89: #define PetscAbs(a) (((a) >= 0) ? a : -a)
90: /*
91: Allows compiling PETSc so that matrix values are stored in
92: single precision but all other objects still use double
93: precision. This does not work for complex numbers in that case
94: it remains double
96: EXPERIMENTAL! NOT YET COMPLETELY WORKING
97: */
99: #if defined(PETSC_USE_MAT_SINGLE)
100: typedef float MatScalar;
101: #else
102: typedef PetscScalar MatScalar;
103: #endif
105: #if defined(PETSC_USE_COMPLEX)
106: typedef double MatReal;
107: #elif defined(PETSC_USE_MAT_SINGLE) || defined(PETSC_USE_SINGLE)
108: typedef float MatReal;
109: #else
110: typedef double MatReal;
111: #endif
113: #if defined(PETSC_USE_SINGLE)
114: typedef float PetscReal;
115: #else
116: typedef double PetscReal;
117: #endif
119: /* --------------------------------------------------------------------------*/
121: /*
122: Certain objects may be created using either single
123: or double precision.
124: */
125: typedef enum { PETSC_SCALAR_DOUBLE,PETSC_SCALAR_SINGLE } PetscScalarPrecision;
127: /* PETSC_i is the imaginary number, i */
128: extern PetscScalar PETSC_i;
130: #define PetscMin(a,b) (((a)<(b)) ? (a) : (b))
131: #define PetscMax(a,b) (((a)<(b)) ? (b) : (a))
132: #define PetscAbsInt(a) (((a)<0) ? -(a) : (a))
133: #define PetscAbsReal(a) (((a)<0) ? -(a) : (a))
134: #define PetscSqr(a) ((a)*(a))
136: /* ----------------------------------------------------------------------------*/
137: /*
138: Basic constants
139: */
140: #define PETSC_PI 3.14159265358979323846264
141: #define PETSC_DEGREES_TO_RADIANS 0.01745329251994
142: #define PETSC_MAX_INT 1000000000
143: #define PETSC_MIN_INT -1000000000
145: #if defined(PETSC_USE_SINGLE)
146: # define PETSC_MAX 1.e30
147: # define PETSC_MIN -1.e30
148: # define PETSC_MACHINE_EPSILON 1.e-7
149: # define PETSC_SQRT_MACHINE_EPSILON 3.e-4
150: # define PETSC_SMALL 1.e-5
151: #else
152: # define PETSC_MAX 1.e300
153: # define PETSC_MIN -1.e300
154: # define PETSC_MACHINE_EPSILON 1.e-14
155: # define PETSC_SQRT_MACHINE_EPSILON 1.e-7
156: # define PETSC_SMALL 1.e-10
157: #endif
159: extern int PetscGlobalMax(PetscReal*,PetscReal*,MPI_Comm);
160: extern int PetscGlobalMin(PetscReal*,PetscReal*,MPI_Comm);
161: extern int PetscGlobalSum(PetscScalar*,PetscScalar*,MPI_Comm);
164: /* ----------------------------------------------------------------------------*/
165: /*
166: PetscLogDouble variables are used to contain double precision numbers
167: that are not used in the numerical computations, but rather in logging,
168: timing etc.
169: */
170: typedef double PetscLogDouble;
171: /*
172: Once PETSc is compiling with a ADIC enhanced version of MPI
173: we will create a new MPI_Datatype for the inactive double variables.
174: */
175: #if defined(AD_DERIV_H)
176: /* extern MPI_Datatype MPIU_PETSCLOGDOUBLE; */
177: #else
178: #if !defined(_petsc_mpi_uni)
179: #define MPIU_PETSCLOGDOUBLE MPI_DOUBLE
180: #endif
181: #endif
183: #define PassiveReal PetscReal
184: #define PassiveScalar PetscScalar
186: #define PETSCMAP1_a(a,b) a ## _ ## b
187: #define PETSCMAP1_b(a,b) PETSCMAP1_a(a,b)
188: #define PETSCMAP1(a) PETSCMAP1_b(a,PetscScalar)
189: #endif