Actual source code: slepcmath.h

slepc-3.18.1 2022-11-02
Report Typos and Errors
  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: #define SlepcDefaultTol(tol) (((tol)==PETSC_DEFAULT)?SLEPC_DEFAULT_TOL:(tol))

 35: /*@C
 36:    SlepcAbs - Returns sqrt(x**2+y**2), taking care not to cause unnecessary
 37:    overflow. It is based on LAPACK's DLAPY2.

 39:    Not Collective

 41:    Input parameters:
 42: .  x,y - the real numbers

 44:    Output parameter:
 45: .  return - the result

 47:    Note:
 48:    This function is not available from Fortran.

 50:    Level: developer
 51: @*/
 52: static inline PetscReal SlepcAbs(PetscReal x,PetscReal y)
 53: {
 54:   PetscReal w,z,t,xabs=PetscAbs(x),yabs=PetscAbs(y);

 56:   w = PetscMax(xabs,yabs);
 57:   z = PetscMin(xabs,yabs);
 58:   if (PetscUnlikely(z == 0.0)) return w;
 59:   t = z/w;
 60:   return w*PetscSqrtReal(1.0+t*t);
 61: }

 63: /*MC
 64:    SlepcAbsEigenvalue - Returns the absolute value of a complex number given
 65:    its real and imaginary parts.

 67:    Synopsis:
 68:    PetscReal SlepcAbsEigenvalue(PetscScalar x,PetscScalar y)

 70:    Not Collective

 72:    Input parameters:
 73: +  x  - the real part of the complex number
 74: -  y  - the imaginary part of the complex number

 76:    Notes:
 77:    This function computes sqrt(x**2+y**2), taking care not to cause unnecessary
 78:    overflow. It is based on LAPACK's DLAPY2.

 80:    In complex scalars, only the first argument is used.

 82:    This function is not available from Fortran.

 84:    Level: developer
 85: M*/
 86: #if !defined(PETSC_USE_COMPLEX)
 87: #define SlepcAbsEigenvalue(x,y) SlepcAbs(x,y)
 88: #else
 89: #define SlepcAbsEigenvalue(x,y) PetscAbsScalar(x)
 90: #endif

 92: #endif

 94: /*
 95:    SlepcSetFlushToZero - Set the FTZ flag in floating-point arithmetic.
 96: */
 97: static inline PetscErrorCode SlepcSetFlushToZero(unsigned int *state)
 98: {
 99: #if defined(PETSC_HAVE_XMMINTRIN_H) && defined(_MM_FLUSH_ZERO_ON) && defined(__SSE__)
100:   *state = _MM_GET_FLUSH_ZERO_MODE();
101:   _MM_SET_FLUSH_ZERO_MODE(_MM_FLUSH_ZERO_ON);
102: #else
103:   *state = 0;
104: #endif
105:   return 0;
106: }

108: /*
109:    SlepcResetFlushToZero - Reset the FTZ flag in floating-point arithmetic.
110: */
111: static inline PetscErrorCode SlepcResetFlushToZero(unsigned int *state)
112: {
113: #if defined(PETSC_HAVE_XMMINTRIN_H) && defined(_MM_FLUSH_ZERO_MASK) && defined(__SSE__)
114:   _MM_SET_FLUSH_ZERO_MODE(*state & _MM_FLUSH_ZERO_MASK);
115: #else
116:   *state = 0;
117: #endif
118:   return 0;
119: }