Fortran如何定义大整数类型以支持25位数字素数判断计算

你当前使用的integer(kind=8)是64位有符号整数类型，最大值为2^63-1即9223372036854775807，最多只能存储19位十进制整数，无法满足25位数字的存储、计算需求，可通过以下两种方案实现目标。

方案1：使用编译器原生支持的大整数类型（推荐，实现成本最低） #

主流Fortran编译器（gfortran、Intel Fortran等）基本都支持128位扩展整数类型，128位有符号整数的最大值为2^127-1，对应38位左右的十进制数，完全覆盖25位数字的使用需求。

不要硬编码kind数值（比如直接写integer(kind=16)），不同编译器的kind编号规则不统一，可移植性差。应当使用Fortran标准内置函数selected_int_kind(r)获取对应精度的整数kind值，这个函数会返回能表示所有绝对值小于10^r的整数的最小kind编号，如果编译器不支持对应精度，会返回-1。

注意：如果调用selected_int_kind(25)返回-1，说明当前编译器不支持原生128位整数，需要切换到方案2。

针对你的素数判断需求，只需要替换变量的整数类型定义，微调类型转换逻辑即可，优化后的可运行参考代码如下：

program cp
implicit none
! 获取可存储25位十进制整数的整数类型kind
integer, parameter :: long_int = selected_int_kind(25)
integer(kind=long_int) :: nmb, ind
double precision :: sqnmb
integer :: io_stat

! 先检查编译器是否支持所需精度
if (long_int == -1) then
    print *, "当前编译器不支持25位原生整数类型，请使用大整数算术方案。"
    stop
end if

print*, "当前整数类型可存储的最大值为 ", huge(nmb)

DO
    print *,"Enter the number larger than 1:"
    read (*,*, iostat=io_stat) nmb
    ! 处理输入非法、输入值超范围的情况
    if (io_stat /=0) then
        print *, "输入无效，请输入合法的正整数。"
        cycle
    end if

    IF ( nmb <= 1) THEN 
        print *, "Please enter a number larger than 1."
        cycle
    END IF
    EXIT
END DO 

! 处理2的特殊情况
if (nmb == 2) then
    print *, "The number", nmb ,"is a prime"
    stop
end if
! 先判断是否为偶数，减少一半循环量
if (mod(nmb, 2_long_int) == 0) then
    print *, "The number", nmb ,"is not a prime, it is divisible by 2"
    stop
end if

sqnmb = sqrt(real(nmb, kind=kind(1.0d0)))
! 从3开始步长为2，只判断奇数因子提升效率
do ind=3_long_int, int(sqnmb, kind=long_int), 2_long_int
    IF (mod(nmb,ind) == 0) THEN 
        print *, "The number", nmb ,"is not a prime, it is divisible by",ind
        STOP
    END IF
end do

print *, "The number", nmb ,"is a prime"
end program cp

方案2：使用大整数算术实现 #

如果你的运行环境编译器不支持128位原生整数，可以通过大整数算术方案实现：

• 轻量场景可自行实现大整数存储：25位整数长度很短，只需要把数字拆分存储到2个64位整数元素的数组中，自行实现大整数取模、整数开方、比较的基础逻辑即可，开发成本很低。
• 通用场景可直接使用成熟的Fortran大整数算术库，这类库支持任意长度整数的各类运算，直接调用库提供的取模、开方接口替换代码中的内置mod、sqrt函数即可，不需要自行实现底层运算逻辑。

内容的提问来源于stack exchange，提问作者brownser