Actual source code: slepcmath.h

slepc-3.13.3 2020-06-14
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2020, 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: /*
 19:     Default tolerance for the different solvers, depending on the precision
 20: */
 21: #if defined(PETSC_USE_REAL_SINGLE)
 22: #  define SLEPC_DEFAULT_TOL   1e-6
 23: #elif defined(PETSC_USE_REAL_DOUBLE)
 24: #  define SLEPC_DEFAULT_TOL   1e-8
 25: #elif defined(PETSC_USE_REAL___FLOAT128)
 26: #  define SLEPC_DEFAULT_TOL   1e-16
 27: #else
 28: #  define SLEPC_DEFAULT_TOL   1e-7
 29: #endif

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

 35:    Not Collective

 37:    Input parameters:
 38: .  x,y - the real numbers

 40:    Output parameter:
 41: .  return - the result

 43:    Note:
 44:    This function is not available from Fortran.

 46:    Level: developer
 47: @*/
 48: PETSC_STATIC_INLINE PetscReal SlepcAbs(PetscReal x,PetscReal y)
 49: {
 50:   PetscReal w,z,t,xabs=PetscAbs(x),yabs=PetscAbs(y);

 52:   w = PetscMax(xabs,yabs);
 53:   z = PetscMin(xabs,yabs);
 54:   if (z == 0.0) return w;
 55:   t = z/w;
 56:   return w*PetscSqrtReal(1.0+t*t);
 57: }

 59: /*MC
 60:    SlepcAbsEigenvalue - Returns the absolute value of a complex number given
 61:    its real and imaginary parts.

 63:    Synopsis:
 64:    PetscReal SlepcAbsEigenvalue(PetscScalar x,PetscScalar y)

 66:    Not Collective

 68:    Input parameters:
 69: +  x  - the real part of the complex number
 70: -  y  - the imaginary part of the complex number

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

 76:    This function is not available from Fortran.

 78:    Level: developer
 79: M*/
 80: #if !defined(PETSC_USE_COMPLEX)
 81: #define SlepcAbsEigenvalue(x,y) SlepcAbs(x,y)
 82: #else
 83: #define SlepcAbsEigenvalue(x,y) PetscAbsScalar(x)
 84: #endif

 86: #endif


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

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