开启fortran浮点数异常处理

fortran对于浮点数异常默认是忽略的,如果要开启,有两种方法,一种可以通过调用c函数去开启这部分异常处理,并且需要将函数打包成库,链接到fortran中,另一个是通过编译选项设置。

第一种方法的代码:

trapfpe.c

/*
* ref:
* https://gcc.gnu.org/onlinedocs/gcc-3.3.6/g77/Floating_002dpoint-Exception-Handling.html
* https://riptutorial.com/fortran/example/7149/calling-c-from-fortran
* https://stackoverflow.com/questions/17845931/calling-c-function-subroutine-in-fortran-code
* https://stackoverflow.com/questions/1202494/why-doesnt-attribute-constructor-work-in-a-static-library
*/

#define _GNU_SOURCE 1
#include <fenv.h>
#include <stdio.h>

// 在main函数开始前调用
// fortran对于浮点数异常默认是忽略的,这里只能通过调用c函数去开启这部分异常处理
// 并且需要将函数打包成库,链接到fortran中
static void __attribute__ ((constructor))
trapfpe ()
{
  /* Enable some exceptions.  At startup all exceptions are masked.  */
  printf("trapfpe: enabling exceptions\n");
  feenableexcept (FE_INVALID|FE_DIVBYZERO|FE_OVERFLOW);
}

void fun()
{
    printf("func\n");
}

// static void before(void) __attribute__((constructor));
// static void after(void) __attribute__((destructor));

// static void before()
// {
//     printf("before main\n");
// }

// static void after(void)
// {

//     printf("after main\n");

// }

main.f90

program main
    implicit none

! 必须显式调用静态库中一个的函数,才能让trapefpe生效
! 因此这里定义一个接口,让fortran调用c函数
INTERFACE
  SUBROUTINE fun() BIND(C)
    USE, INTRINSIC :: ISO_C_BINDING, ONLY: C_INT
    IMPLICIT NONE
  END SUBROUTINE fun
END INTERFACE

    double precision :: a
    double precision :: b
    double precision :: c

    a = 1.0
    b = 0.0
    c = a / b
    print *, c
    call fun()
end program

编译方法:

# 生成trapfe.o
gcc -c trapfpe.c
# 生成静态库文件libtrapfe.a
ar -r libtrapfpe.a trapfpe.o
# 将.目录下的trapfe.a作为静态库,链接到main中
gfortran -g main.f90 -L. -ltrapfpe

运行结果:

➜  test ./a.out 
trapfpe: enabling exceptions

Program received signal SIGFPE: Floating-point exception - erroneous arithmetic operation.

Backtrace for this error:
#0  0x7fd098838d21 in ???
#1  0x7fd098837ef5 in ???
#2  0x7fd09851a0bf in ???
#3  0x5619e64951bf in MAIN__
        at /mnt/f/programing/test/main.f90:19
#4  0x5619e649526b in main
        at /mnt/f/programing/test/main.f90:22
[1]    990 floating point exception  ./a.out

可以看到,第19行出现了浮点数异常,因为这里做了除0操作。

第二种方法的代码

fortran

! https://gcc.gnu.org/onlinedocs/gfortran/Debugging-Options.html
program main 
    implicit none 
    real a,b,c
    b = 1.0
    c = 0.0
    a = b / c 
    write(*,*) a
end program main

编译:

gfortran -g main2.f90 -ffpe-trap=invalid,zero,overflow

输出:

➜  fortran-exception-handle ./a.out 

Program received signal SIGFPE: Floating-point exception - erroneous arithmetic operation.

Backtrace for this error:
#0  0x7fa4d0733d21 in ???
#1  0x7fa4d0732ef5 in ???
#2  0x7fa4d05640bf in ???
#3  0x5620d91911af in MAIN__
        at /mnt/f/programing/tools-code/fortran-exception-handle/main2.f90:6
#4  0x5620d9191260 in main
        at /mnt/f/programing/tools-code/fortran-exception-handle/main2.f90:8
[1]    1161 floating point exception  ./a.out
posted @ 2022-04-08 11:55  JayYin  阅读(192)  评论(0编辑  收藏  举报