Actual source code: slepcmath.h
slepc-3.19.2 2023-09-05
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
10: /*
11: SLEPc mathematics include file. Defines basic operations and functions.
12: This file is included by slepcsys.h and should not be used directly.
13: */
15: #if !defined(SLEPCMATH_H)
16: #define SLEPCMATH_H
18: /* SUBMANSEC = sys */
20: /*
21: Default tolerance for the different solvers, depending on the precision
22: */
23: #if defined(PETSC_USE_REAL_SINGLE)
24: # define SLEPC_DEFAULT_TOL 1e-5
25: #elif defined(PETSC_USE_REAL_DOUBLE)
26: # define SLEPC_DEFAULT_TOL 1e-8
27: #elif defined(PETSC_USE_REAL___FLOAT128)
28: # define SLEPC_DEFAULT_TOL 1e-16
29: #elif defined(PETSC_USE_REAL___FP16)
30: # define SLEPC_DEFAULT_TOL 1e-2
31: #endif
33: static inline PetscReal SlepcDefaultTol(PetscReal tol)
34: {
35: return tol == (PetscReal)PETSC_DEFAULT ? SLEPC_DEFAULT_TOL : tol;
36: }
38: /*@C
39: SlepcAbs - Returns sqrt(x**2+y**2), taking care not to cause unnecessary
40: overflow. It is based on LAPACK's DLAPY2.
42: Not Collective
44: Input parameters:
45: . x,y - the real numbers
47: Output parameter:
48: . return - the result
50: Note:
51: This function is not available from Fortran.
53: Level: developer
54: @*/
55: static inline PetscReal SlepcAbs(PetscReal x,PetscReal y)
56: {
57: PetscReal w,z,t,xabs=PetscAbs(x),yabs=PetscAbs(y);
59: w = PetscMax(xabs,yabs);
60: z = PetscMin(xabs,yabs);
61: if (PetscUnlikely(z == 0.0)) return w;
62: t = z/w;
63: return w*PetscSqrtReal(1.0+t*t);
64: }
66: /*MC
67: SlepcAbsEigenvalue - Returns the absolute value of a complex number given
68: its real and imaginary parts.
70: Synopsis:
71: PetscReal SlepcAbsEigenvalue(PetscScalar x,PetscScalar y)
73: Not Collective
75: Input parameters:
76: + x - the real part of the complex number
77: - y - the imaginary part of the complex number
79: Notes:
80: This function computes sqrt(x**2+y**2), taking care not to cause unnecessary
81: overflow. It is based on LAPACK's DLAPY2.
83: In complex scalars, only the first argument is used.
85: This function is not available from Fortran.
87: Level: developer
88: M*/
89: #if !defined(PETSC_USE_COMPLEX)
90: #define SlepcAbsEigenvalue(x,y) SlepcAbs(x,y)
91: #else
92: #define SlepcAbsEigenvalue(x,y) PetscAbsScalar(x)
93: #endif
95: #endif
97: /*
98: SlepcSetFlushToZero - Set the FTZ flag in floating-point arithmetic.
99: */
100: static inline PetscErrorCode SlepcSetFlushToZero(unsigned int *state)
101: {
102: PetscFunctionBegin;
103: #if defined(PETSC_HAVE_XMMINTRIN_H) && defined(_MM_FLUSH_ZERO_ON) && defined(__SSE__)
104: *state = _MM_GET_FLUSH_ZERO_MODE();
105: _MM_SET_FLUSH_ZERO_MODE(_MM_FLUSH_ZERO_ON);
106: #else
107: *state = 0;
108: #endif
109: PetscFunctionReturn(PETSC_SUCCESS);
110: }
112: /*
113: SlepcResetFlushToZero - Reset the FTZ flag in floating-point arithmetic.
114: */
115: static inline PetscErrorCode SlepcResetFlushToZero(unsigned int *state)
116: {
117: PetscFunctionBegin;
118: #if defined(PETSC_HAVE_XMMINTRIN_H) && defined(_MM_FLUSH_ZERO_MASK) && defined(__SSE__)
119: _MM_SET_FLUSH_ZERO_MODE(*state & _MM_FLUSH_ZERO_MASK);
120: #else
121: *state = 0;
122: #endif
123: PetscFunctionReturn(PETSC_SUCCESS);
124: }