Actual source code: fp.c


  2: /*
  3:    IEEE error handler for all machines. Since each OS has
  4:    enough slight differences we have completely separate codes for each one.
  5: */

  7: /*
  8:   This feature test macro provides FE_NOMASK_ENV on GNU.  It must be defined
  9:   at the top of the file because other headers may pull in fenv.h even when
 10:   not strictly necessary.  Strictly speaking, we could include ONLY petscconf.h,
 11:   check PETSC_HAVE_FENV_H, and only define _GNU_SOURCE in that case, but such
 12:   shenanigans ought to be unnecessary.
 13: */
 14: #if !defined(_GNU_SOURCE)
 15: #define _GNU_SOURCE
 16: #endif

 18: #include <petsc/private/petscimpl.h>
 19: #include <signal.h>

 21: struct PetscFPTrapLink {
 22:   PetscFPTrap            trapmode;
 23:   struct PetscFPTrapLink *next;
 24: };
 25: static PetscFPTrap            _trapmode = PETSC_FP_TRAP_OFF; /* Current trapping mode; see PetscDetermineInitialFPTrap() */
 26: static struct PetscFPTrapLink *_trapstack;                   /* Any pushed states of _trapmode */

 28: /*@
 29:    PetscFPTrapPush - push a floating point trapping mode, restored using PetscFPTrapPop()

 31:    Not Collective

 33:    Input Parameter:
 34: .    trap - PETSC_FP_TRAP_ON or PETSC_FP_TRAP_OFF

 36:    Level: advanced

 38:    Notes:
 39:      This only changes the trapping if the new mode is different than the current mode.

 41:      This routine is called to turn off trapping for certain LAPACK routines that assume that dividing
 42:      by zero is acceptable. In particular the routine ieeeck().

 44:      Most systems by default have all trapping turned off, but certain Fortran compilers have
 45:      link flags that turn on trapping before the program begins.
 46: $       gfortran -ffpe-trap=invalid,zero,overflow,underflow,denormal
 47: $       ifort -fpe0

 49: .seealso: PetscFPTrapPop(), PetscSetFPTrap(), PetscDetermineInitialFPTrap()
 50: @*/
 51: PetscErrorCode PetscFPTrapPush(PetscFPTrap trap)
 52: {
 53:   struct PetscFPTrapLink *link;

 55:   PetscNew(&link);
 56: #if defined(PETSC_HAVE_THREADSAFETY) && defined(PETSC_HAVE_OPENMP)
 57: #pragma omp critical
 58: #endif
 59:   {
 60:     link->trapmode = _trapmode;
 61:     link->next     = _trapstack;
 62:     _trapstack     = link;
 63:   }
 64:   if (trap != _trapmode) PetscSetFPTrap(trap);
 65:   return 0;
 66: }

 68: /*@
 69:    PetscFPTrapPop - push a floating point trapping mode, to be restored using PetscFPTrapPop()

 71:    Not Collective

 73:    Level: advanced

 75: .seealso: PetscFPTrapPush(), PetscSetFPTrap(), PetscDetermineInitialFPTrap()
 76: @*/
 77: PetscErrorCode PetscFPTrapPop(void)
 78: {
 79:   struct PetscFPTrapLink *link;

 81:   if (_trapstack->trapmode != _trapmode) PetscSetFPTrap(_trapstack->trapmode);
 82: #if defined(PETSC_HAVE_THREADSAFETY) && defined(PETSC_HAVE_OPENMP)
 83: #pragma omp critical
 84: #endif
 85:   {
 86:     link       = _trapstack;
 87:     _trapstack = _trapstack->next;
 88:   }
 89:   PetscFree(link);
 90:   return 0;
 91: }

 93: /*--------------------------------------- ---------------------------------------------------*/
 94: #if defined(PETSC_HAVE_SUN4_STYLE_FPTRAP)
 95: #include <floatingpoint.h>

 97: PETSC_EXTERN PetscErrorCode ieee_flags(char*,char*,char*,char**);
 98: PETSC_EXTERN PetscErrorCode ieee_handler(char*,char*,sigfpe_handler_type(int,int,struct sigcontext*,char*));

100: static struct { int code_no; char *name; } error_codes[] = {
101:   { FPE_INTDIV_TRAP    ,"integer divide" },
102:   { FPE_FLTOPERR_TRAP  ,"IEEE operand error" },
103:   { FPE_FLTOVF_TRAP    ,"floating point overflow" },
104:   { FPE_FLTUND_TRAP    ,"floating point underflow" },
105:   { FPE_FLTDIV_TRAP    ,"floating pointing divide" },
106:   { FPE_FLTINEX_TRAP   ,"inexact floating point result" },
107:   { 0                  ,"unknown error" }
108: };
109: #define SIGPC(scp) (scp->sc_pc)

111: sigfpe_handler_type PetscDefaultFPTrap(int sig,int code,struct sigcontext *scp,char *addr)
112: {
113:   int err_ind = -1;

115:   for (int j = 0; error_codes[j].code_no; j++) {
116:     if (error_codes[j].code_no == code) err_ind = j;
117:   }

119:   if (err_ind >= 0) (*PetscErrorPrintf)("*** %s occurred at pc=%X ***\n",error_codes[err_ind].name,SIGPC(scp));
120:   else              (*PetscErrorPrintf)("*** floating point error 0x%x occurred at pc=%X ***\n",code,SIGPC(scp));

122:   (void)PetscError(PETSC_COMM_SELF,PETSC_ERR_FP,"User provided function","Unknown file",PETSC_ERR_FP,PETSC_ERROR_REPEAT,"floating point error");
123:   PETSCABORT(MPI_COMM_WORLD,PETSC_ERR_FP);
124:   return 0;
125: }

127: /*@
128:    PetscSetFPTrap - Enables traps/exceptions on common floating point errors.
129:                     This option may not work on certain systems.

131:    Not Collective

133:    Input Parameters:
134: .  flag - PETSC_FP_TRAP_ON, PETSC_FP_TRAP_OFF.

136:    Options Database Keys:
137: .  -fp_trap - Activates floating point trapping

139:    Level: advanced

141:    Description:
142:    On systems that support it, when called with PETSC_FP_TRAP_ON this routine causes floating point
143:    underflow, overflow, divide-by-zero, and invalid-operand (e.g., a NaN) to
144:    cause a message to be printed and the program to exit.

146:    Note:
147:    On many common systems, the floating
148:    point exception state is not preserved from the location where the trap
149:    occurred through to the signal handler.  In this case, the signal handler
150:    will just say that an unknown floating point exception occurred and which
151:    function it occurred in.  If you run with -fp_trap in a debugger, it will
152:    break on the line where the error occurred.  On systems that support C99
153:    floating point exception handling You can check which
154:    exception occurred using fetestexcept(FE_ALL_EXCEPT).  See fenv.h
155:    (usually at /usr/include/bits/fenv.h) for the enum values on your system.

157:    Caution:
158:    On certain machines, in particular the IBM PowerPC, floating point
159:    trapping may be VERY slow!

161: .seealso: PetscFPTrapPush(), PetscFPTrapPop(), PetscDetermineInitialFPTrap()
162: @*/
163: PetscErrorCode PetscSetFPTrap(PetscFPTrap flag)
164: {
165:   char *out;

167:   /* Clear accumulated exceptions.  Used to suppress meaningless messages from f77 programs */
168:   (void) ieee_flags("clear","exception","all",&out);
169:   if (flag == PETSC_FP_TRAP_ON) {
170:     /*
171:       To trap more fp exceptions, including underflow, change the line below to
172:       if (ieee_handler("set","all",PetscDefaultFPTrap)) {
173:     */
174:     if (ieee_handler("set","common",PetscDefaultFPTrap))        (*PetscErrorPrintf)("Can't set floatingpoint handler\n");
175:   } else if (ieee_handler("clear","common",PetscDefaultFPTrap)) (*PetscErrorPrintf)("Can't clear floatingpoint handler\n");

177:   _trapmode = flag;
178:   return 0;
179: }

181: /*@
182:    PetscDetermineInitialFPTrap - Attempts to determine the floating point trapping that exists when PetscInitialize() is called

184:    Not Collective

186:    Notes:
187:       Currently only supported on Linux and MacOS. Checks if divide by zero is enable and if so declares that trapping is on.

189:    Level: advanced

191: .seealso: PetscFPTrapPush(), PetscFPTrapPop(), PetscDetermineInitialFPTrap()
192: @*/
193: PetscErrorCode  PetscDetermineInitialFPTrap(void)
194: {
195:   PetscInfo(NULL,"Unable to determine initial floating point trapping. Assuming it is off\n");
196:   return 0;
197: }

199: /* -------------------------------------------------------------------------------------------*/
200: #elif defined(PETSC_HAVE_SOLARIS_STYLE_FPTRAP)
201: #include <sunmath.h>
202: #include <floatingpoint.h>
203: #include <siginfo.h>
204: #include <ucontext.h>

206: static struct { int code_no; char *name; } error_codes[] = {
207:   { FPE_FLTINV,"invalid floating point operand"},
208:   { FPE_FLTRES,"inexact floating point result"},
209:   { FPE_FLTDIV,"division-by-zero"},
210:   { FPE_FLTUND,"floating point underflow"},
211:   { FPE_FLTOVF,"floating point overflow"},
212:   { 0,         "unknown error"}
213: };
214: #define SIGPC(scp) (scp->si_addr)

216: void PetscDefaultFPTrap(int sig,siginfo_t *scp,ucontext_t *uap)
217: {
218:   int err_ind = -1,code = scp->si_code;

220:   for (int j = 0; error_codes[j].code_no; j++) {
221:     if (error_codes[j].code_no == code) err_ind = j;
222:   }

224:   if (err_ind >= 0) (*PetscErrorPrintf)("*** %s occurred at pc=%X ***\n",error_codes[err_ind].name,SIGPC(scp));
225:   else              (*PetscErrorPrintf)("*** floating point error 0x%x occurred at pc=%X ***\n",code,SIGPC(scp));

227:   (void)PetscError(PETSC_COMM_SELF,0,"User provided function","Unknown file",PETSC_ERR_FP,PETSC_ERROR_REPEAT,"floating point error");
228:   PETSCABORT(MPI_COMM_WORLD,PETSC_ERR_FP);
229: }

231: PetscErrorCode PetscSetFPTrap(PetscFPTrap flag)
232: {
233:   char *out;

235:   /* Clear accumulated exceptions.  Used to suppress meaningless messages from f77 programs */
236:   (void) ieee_flags("clear","exception","all",&out);
237:   if (flag == PETSC_FP_TRAP_ON) {
238:     if (ieee_handler("set","common",(sigfpe_handler_type)PetscDefaultFPTrap))        (*PetscErrorPrintf)("Can't set floating point handler\n");
239:   } else {
240:     if (ieee_handler("clear","common",(sigfpe_handler_type)PetscDefaultFPTrap)) (*PetscErrorPrintf)("Can't clear floatingpoint handler\n");
241:   }
242:   _trapmode = flag;
243:   return 0;
244: }

246: PetscErrorCode  PetscDetermineInitialFPTrap(void)
247: {
248:   PetscInfo(NULL,"Unable to determine initial floating point trapping. Assuming it is off\n");
249:   return 0;
250: }

252: /* ------------------------------------------------------------------------------------------*/
253: #elif defined(PETSC_HAVE_IRIX_STYLE_FPTRAP)
254: #include <sigfpe.h>
255: static struct { int code_no; char *name; } error_codes[] = {
256:   { _INVALID   ,"IEEE operand error" },
257:   { _OVERFL    ,"floating point overflow" },
258:   { _UNDERFL   ,"floating point underflow" },
259:   { _DIVZERO   ,"floating point divide" },
260:   { 0          ,"unknown error" }
261: } ;
262: void PetscDefaultFPTrap(unsigned exception[],int val[])
263: {
264:   int err_ind = -1,code = exception[0];

266:   for (int j = 0; error_codes[j].code_no; j++) {
267:     if (error_codes[j].code_no == code) err_ind = j;
268:   }
269:   if (err_ind >= 0) (*PetscErrorPrintf)("*** %s occurred ***\n",error_codes[err_ind].name);
270:   else              (*PetscErrorPrintf)("*** floating point error 0x%x occurred ***\n",code);

272:   (void)PetscError(PETSC_COMM_SELF,0,"User provided function","Unknown file",PETSC_ERR_FP,PETSC_ERROR_REPEAT,"floating point error");
273:   PETSCABORT(MPI_COMM_WORLD,PETSC_ERR_FP);
274: }

276: PetscErrorCode PetscSetFPTrap(PetscFPTrap flag)
277: {
278:   if (flag == PETSC_FP_TRAP_ON) handle_sigfpes(_ON,,_EN_UNDERFL|_EN_OVERFL|_EN_DIVZERO|_EN_INVALID,PetscDefaultFPTrap,_ABORT_ON_ERROR,0);
279:   else                          handle_sigfpes(_OFF,_EN_UNDERFL|_EN_OVERFL|_EN_DIVZERO|_EN_INVALID,0,_ABORT_ON_ERROR,0);
280:   _trapmode = flag;
281:   return 0;
282: }

284: PetscErrorCode  PetscDetermineInitialFPTrap(void)
285: {
286:   PetscInfo(NULL,"Unable to determine initial floating point trapping. Assuming it is off\n");
287:   return 0;
288: }

290: /* -------------------------------------------------------------------------------------------*/
291: #elif defined(PETSC_HAVE_SOLARIS_STYLE_FPTRAP)
292: #include <sunmath.h>
293: #include <floatingpoint.h>
294: #include <siginfo.h>
295: #include <ucontext.h>

297: static struct { int code_no; char *name; } error_codes[] = {
298:   { FPE_FLTINV,"invalid floating point operand"},
299:   { FPE_FLTRES,"inexact floating point result"},
300:   { FPE_FLTDIV,"division-by-zero"},
301:   { FPE_FLTUND,"floating point underflow"},
302:   { FPE_FLTOVF,"floating point overflow"},
303:   { 0,         "unknown error"}
304: };
305: #define SIGPC(scp) (scp->si_addr)

307: void PetscDefaultFPTrap(int sig,siginfo_t *scp,ucontext_t *uap)
308: {
309:   int err_ind = -1,code = scp->si_code;

311:   err_ind = -1;
312:   for (int j = 0; error_codes[j].code_no; j++) {
313:     if (error_codes[j].code_no == code) err_ind = j;
314:   }

316:   if (err_ind >= 0) (*PetscErrorPrintf)("*** %s occurred at pc=%X ***\n",error_codes[err_ind].name,SIGPC(scp));
317:   else              (*PetscErrorPrintf)("*** floating point error 0x%x occurred at pc=%X ***\n",code,SIGPC(scp));

319:   (void)PetscError(PETSC_COMM_SELF,0,"User provided function","Unknown file",PETSC_ERR_FP,PETSC_ERROR_REPEAT,"floating point error");
320:   PETSCABORT(MPI_COMM_WORLD,PETSC_ERR_FP);
321: }

323: PetscErrorCode PetscSetFPTrap(PetscFPTrap flag)
324: {
325:   char *out;

327:   /* Clear accumulated exceptions.  Used to suppress meaningless messages from f77 programs */
328:   (void) ieee_flags("clear","exception","all",&out);
329:   if (flag == PETSC_FP_TRAP_ON) {
330:     if (ieee_handler("set","common",(sigfpe_handler_type)PetscDefaultFPTrap))        (*PetscErrorPrintf)("Can't set floating point handler\n");
331:   } else {
332:     if (ieee_handler("clear","common",(sigfpe_handler_type)PetscDefaultFPTrap)) (*PetscErrorPrintf)("Can't clear floatingpoint handler\n");
333:   }
334:   _trapmode = flag;
335:   return 0;
336: }

338: PetscErrorCode  PetscDetermineInitialFPTrap(void)
339: {
340:   PetscInfo(NULL,"Unable to determine initial floating point trapping. Assuming it is off\n");
341:   return 0;
342: }

344: /*----------------------------------------------- --------------------------------------------*/
345: #elif defined(PETSC_HAVE_RS6000_STYLE_FPTRAP)
346: /* In "fast" mode, floating point traps are imprecise and ignored.
347:    This is the reason for the fptrap(FP_TRAP_SYNC) call */
348: struct sigcontext;
349: #include <fpxcp.h>
350: #include <fptrap.h>
351: #define FPE_FLTOPERR_TRAP (fptrap_t)(0x20000000)
352: #define FPE_FLTOVF_TRAP   (fptrap_t)(0x10000000)
353: #define FPE_FLTUND_TRAP   (fptrap_t)(0x08000000)
354: #define FPE_FLTDIV_TRAP   (fptrap_t)(0x04000000)
355: #define FPE_FLTINEX_TRAP  (fptrap_t)(0x02000000)

357: static struct { int code_no; char *name; } error_codes[] = {
358:   {FPE_FLTOPERR_TRAP   ,"IEEE operand error" },
359:   { FPE_FLTOVF_TRAP    ,"floating point overflow" },
360:   { FPE_FLTUND_TRAP    ,"floating point underflow" },
361:   { FPE_FLTDIV_TRAP    ,"floating point divide" },
362:   { FPE_FLTINEX_TRAP   ,"inexact floating point result" },
363:   { 0                  ,"unknown error" }
364: } ;
365: #define SIGPC(scp) (0) /* Info MIGHT be in scp->sc_jmpbuf.jmp_context.iar */
366: /*
367:    For some reason, scp->sc_jmpbuf does not work on the RS6000, even though
368:    it looks like it should from the include definitions.  It is probably
369:    some strange interaction with the "POSIX_SOURCE" that we require.
370: */

372: void PetscDefaultFPTrap(int sig,int code,struct sigcontext *scp)
373: {
375:   int            err_ind,j;
376:   fp_ctx_t       flt_context;

378:   fp_sh_trap_info(scp,&flt_context);

380:   err_ind = -1;
381:   for (j = 0; error_codes[j].code_no; j++) {
382:     if (error_codes[j].code_no == flt_context.trap) err_ind = j;
383:   }

385:   if (err_ind >= 0) (*PetscErrorPrintf)("*** %s occurred ***\n",error_codes[err_ind].name);
386:   else              (*PetscErrorPrintf)("*** floating point error 0x%x occurred ***\n",flt_context.trap);

388:   PetscError(PETSC_COMM_SELF,0,"User provided function","Unknown file",PETSC_ERR_FP,PETSC_ERROR_REPEAT,"floating point error");
389:   PETSCABORT(MPI_COMM_WORLD,PETSC_ERR_FP);
390: }

392: PetscErrorCode PetscSetFPTrap(PetscFPTrap on)
393: {
394:   if (on == PETSC_FP_TRAP_ON) {
395:     signal(SIGFPE,(void (*)(int))PetscDefaultFPTrap);
396:     fp_trap(FP_TRAP_SYNC);
397:     fp_enable(TRP_INVALID | TRP_DIV_BY_ZERO | TRP_OVERFLOW | TRP_UNDERFLOW);
398:     /* fp_enable(mask) for individual traps.  Values are:
399:        TRP_INVALID
400:        TRP_DIV_BY_ZERO
401:        TRP_OVERFLOW
402:        TRP_UNDERFLOW
403:        TRP_INEXACT
404:        Can OR then together.
405:        fp_enable_all(); for all traps.
406:     */
407:   } else {
408:     signal(SIGFPE,SIG_DFL);
409:     fp_disable(TRP_INVALID | TRP_DIV_BY_ZERO | TRP_OVERFLOW | TRP_UNDERFLOW);
410:     fp_trap(FP_TRAP_OFF);
411:   }
412:   _trapmode = on;
413:   return 0;
414: }

416: PetscErrorCode  PetscDetermineInitialFPTrap(void)
417: {
418:   PetscInfo(NULL,"Unable to determine initial floating point trapping. Assuming it is off\n");
419:   return 0;
420: }

422: /* ------------------------------------------------------------*/
423: #elif defined(PETSC_HAVE_WINDOWS_COMPILERS)
424: #include <float.h>
425: void PetscDefaultFPTrap(int sig)
426: {
427:   (*PetscErrorPrintf)("*** floating point error occurred ***\n");
428:   PetscError(PETSC_COMM_SELF,0,"User provided function","Unknown file",PETSC_ERR_FP,PETSC_ERROR_REPEAT,"floating point error");
429:   PETSCABORT(MPI_COMM_WORLD,PETSC_ERR_FP);
430: }

432: PetscErrorCode  PetscSetFPTrap(PetscFPTrap on)
433: {
434:   unsigned int cw;

436:   if (on == PETSC_FP_TRAP_ON) {
437:     cw = _EM_INVALID | _EM_ZERODIVIDE | _EM_OVERFLOW | _EM_UNDERFLOW;
439:   } else {
440:     cw = 0;
442:   }
443:   (void)_controlfp(0, cw);
444:   _trapmode = on;
445:   return 0;
446: }

448: PetscErrorCode  PetscDetermineInitialFPTrap(void)
449: {
450:   PetscInfo(NULL,"Unable to determine initial floating point trapping. Assuming it is off\n");
451:   return 0;
452: }

454: /* ------------------------------------------------------------*/
455: #elif defined(PETSC_HAVE_FENV_H) && !defined(__cplusplus)
456: /*
457:    C99 style floating point environment.

459:    Note that C99 merely specifies how to save, restore, and clear the floating
460:    point environment as well as defining an enumeration of exception codes.  In
461:    particular, C99 does not specify how to make floating point exceptions raise
462:    a signal.  Glibc offers this capability through FE_NOMASK_ENV (or with finer
463:    granularity, feenableexcept()), xmmintrin.h offers _MM_SET_EXCEPTION_MASK().
464: */
465: #include <fenv.h>
466: typedef struct {int code; const char *name;} FPNode;
467: static const FPNode error_codes[] = {
468:   {FE_DIVBYZERO,"divide by zero"},
469:   {FE_INEXACT,  "inexact floating point result"},
470:   {FE_INVALID,  "invalid floating point arguments (domain error)"},
471:   {FE_OVERFLOW, "floating point overflow"},
472:   {FE_UNDERFLOW,"floating point underflow"},
473:   {0           ,"unknown error"}
474: };

476: void PetscDefaultFPTrap(int sig)
477: {
478:   const FPNode *node;
479:   int          code;
480:   PetscBool    matched = PETSC_FALSE;

482:   /* Note: While it is possible for the exception state to be preserved by the
483:    * kernel, this seems to be rare which makes the following flag testing almost
484:    * useless.  But on a system where the flags can be preserved, it would provide
485:    * more detail.
486:    */
487:   code = fetestexcept(FE_ALL_EXCEPT);
488:   for (node=&error_codes[0]; node->code; node++) {
489:     if (code & node->code) {
490:       matched = PETSC_TRUE;
491:       (*PetscErrorPrintf)("*** floating point error \"%s\" occurred ***\n",node->name);
492:       code &= ~node->code; /* Unset this flag since it has been processed */
493:     }
494:   }
495:   if (!matched || code) { /* If any remaining flags are set, or we didn't process any flags */
496:     (*PetscErrorPrintf)("*** unknown floating point error occurred ***\n");
497:     (*PetscErrorPrintf)("The specific exception can be determined by running in a debugger.  When the\n");
498:     (*PetscErrorPrintf)("debugger traps the signal, the exception can be found with fetestexcept(0x%x)\n",FE_ALL_EXCEPT);
499:     (*PetscErrorPrintf)("where the result is a bitwise OR of the following flags:\n");
500:     (*PetscErrorPrintf)("FE_INVALID=0x%x FE_DIVBYZERO=0x%x FE_OVERFLOW=0x%x FE_UNDERFLOW=0x%x FE_INEXACT=0x%x\n",FE_INVALID,FE_DIVBYZERO,FE_OVERFLOW,FE_UNDERFLOW,FE_INEXACT);
501:   }

503:   (*PetscErrorPrintf)("Try option -start_in_debugger\n");
504: #if PetscDefined(USE_DEBUG)
505:   (*PetscErrorPrintf)("likely location of problem given in stack below\n");
506:   (*PetscErrorPrintf)("---------------------  Stack Frames ------------------------------------\n");
507:   PetscStackView(PETSC_STDOUT);
508: #else
509:   (*PetscErrorPrintf)("configure using --with-debugging=yes, recompile, link, and run \n");
510:   (*PetscErrorPrintf)("with -start_in_debugger to get more information on the crash.\n");
511: #endif
512:   PetscError(PETSC_COMM_SELF,0,"User provided function","Unknown file",PETSC_ERR_FP,PETSC_ERROR_INITIAL,"trapped floating point error");
513:   PETSCABORT(MPI_COMM_WORLD,PETSC_ERR_FP);
514: }

516: PetscErrorCode  PetscSetFPTrap(PetscFPTrap on)
517: {
518:   if (on == PETSC_FP_TRAP_ON) {
519:     /* Clear any flags that are currently set so that activating trapping will not immediately call the signal handler. */
521: #if defined(FE_NOMASK_ENV)
522:     /* Could use fesetenv(FE_NOMASK_ENV), but that causes spurious exceptions (like gettimeofday() -> PetscLogDouble). */
524: #elif defined PETSC_HAVE_XMMINTRIN_H
525:    _MM_SET_EXCEPTION_MASK(_MM_GET_EXCEPTION_MASK() & ~_MM_MASK_DIV_ZERO);
526:    _MM_SET_EXCEPTION_MASK(_MM_GET_EXCEPTION_MASK() & ~_MM_MASK_UNDERFLOW);
527:    _MM_SET_EXCEPTION_MASK(_MM_GET_EXCEPTION_MASK() & ~_MM_MASK_OVERFLOW);
528:    _MM_SET_EXCEPTION_MASK(_MM_GET_EXCEPTION_MASK() & ~_MM_MASK_INVALID);
529: #else
530:     /* C99 does not provide a way to modify the environment so there is no portable way to activate trapping. */
531: #endif
533:   } else {
535:     /* can use _MM_SET_EXCEPTION_MASK(_MM_GET_EXCEPTION_MASK() | _MM_MASK_UNDERFLOW); if PETSC_HAVE_XMMINTRIN_H exists */
537:   }
538:   _trapmode = on;
539:   return 0;
540: }

542: PetscErrorCode  PetscDetermineInitialFPTrap(void)
543: {
544: #if defined(FE_NOMASK_ENV) || defined PETSC_HAVE_XMMINTRIN_H
545:   unsigned int flags;
546: #endif

548: #if defined(FE_NOMASK_ENV)
549:   flags = fegetexcept();
550:   if (flags & FE_DIVBYZERO) {
551: #elif defined PETSC_HAVE_XMMINTRIN_H
552:   flags = _MM_GET_EXCEPTION_MASK();
553:   if (!(flags & _MM_MASK_DIV_ZERO)) {
554: #else
555:   PetscInfo(NULL,"Floating point trapping unknown, assuming off\n");
556:   return 0;
557: #endif
558: #if defined(FE_NOMASK_ENV) || defined PETSC_HAVE_XMMINTRIN_H
559:     _trapmode = PETSC_FP_TRAP_ON;
560:     PetscInfo(NULL,"Floating point trapping is on by default %d\n",flags);
561:   } else {
562:     _trapmode = PETSC_FP_TRAP_OFF;
563:     PetscInfo(NULL,"Floating point trapping is off by default %d\n",flags);
564:   }
565:   return 0;
566: #endif
567: }

569: /* ------------------------------------------------------------*/
570: #elif defined(PETSC_HAVE_IEEEFP_H)
571: #include <ieeefp.h>
572: void PetscDefaultFPTrap(int sig)
573: {
574:   (*PetscErrorPrintf)("*** floating point error occurred ***\n");
575:   PetscError(PETSC_COMM_SELF,0,"User provided function","Unknown file",PETSC_ERR_FP,PETSC_ERROR_REPEAT,"floating point error");
576:   PETSCABORT(MPI_COMM_WORLD,PETSC_ERR_FP);
577: }

579: PetscErrorCode  PetscSetFPTrap(PetscFPTrap on)
580: {
581:   if (on == PETSC_FP_TRAP_ON) {
582: #if defined(PETSC_HAVE_FPPRESETSTICKY)
583:     fpresetsticky(fpgetsticky());
584: #elif defined(PETSC_HAVE_FPSETSTICKY)
585:     fpsetsticky(fpgetsticky());
586: #endif
587:     fpsetmask(FP_X_INV | FP_X_DZ | FP_X_OFL |  FP_X_OFL);
589:   } else {
590: #if defined(PETSC_HAVE_FPPRESETSTICKY)
591:     fpresetsticky(fpgetsticky());
592: #elif defined(PETSC_HAVE_FPSETSTICKY)
593:     fpsetsticky(fpgetsticky());
594: #endif
595:     fpsetmask(0);
597:   }
598:   _trapmode = on;
599:   return 0;
600: }

602: PetscErrorCode  PetscDetermineInitialFPTrap(void)
603: {
604:   PetscInfo(NULL,"Unable to determine initial floating point trapping. Assuming it is off\n");
605:   return 0;
606: }

608: /* -------------------------Default -----------------------------------*/
609: #else

611: void PetscDefaultFPTrap(int sig)
612: {
613:   (*PetscErrorPrintf)("*** floating point error occurred ***\n");
614:   PetscError(PETSC_COMM_SELF,0,"User provided function","Unknown file",PETSC_ERR_FP,PETSC_ERROR_REPEAT,"floating point error");
615:   PETSCABORT(MPI_COMM_WORLD,PETSC_ERR_FP);
616: }

618: PetscErrorCode  PetscSetFPTrap(PetscFPTrap on)
619: {
620:   if (on == PETSC_FP_TRAP_ON) {
621:     if (SIG_ERR == signal(SIGFPE,PetscDefaultFPTrap)) (*PetscErrorPrintf)("Can't set floatingpoint handler\n");
622:   } else if (SIG_ERR == signal(SIGFPE,SIG_DFL))       (*PetscErrorPrintf)("Can't clear floatingpoint handler\n");

624:   _trapmode = on;
625:   return 0;
626: }

628: PetscErrorCode  PetscDetermineInitialFPTrap(void)
629: {
630:   PetscInfo(NULL,"Unable to determine initial floating point trapping. Assuming it is off\n");
631:   return 0;
632: }
633: #endif