Fortran Coder

查看: 68|回复: 2

[Module] 调用module程序出错

[复制链接]

3

帖子

2

主题

0

精华

入门

F 币
39 元
贡献
19 点
发表于 2020-11-14 21:25:14 | 显示全部楼层 |阅读模式
本帖最后由 ocw^。 于 2020-11-14 21:30 编辑

   编写一个进行50位有效数字进行计算的程序,出来很多bug,但不知道这些bug的原因和如何修改,求助
[Fortran] 纯文本查看 复制代码
module MyPrec
    implicit none

    integer, parameter :: rp = selected_real_kind(p = 50)!保存50位浮点数

end module
    
module MyHmlt
    
    use MyPrec
    implicit none
    
    real(rp) :: t1, t2, p0, er1, er2, er3, er4!两个精度为50的浮点数  
    integer :: n0 !有效数字
    
end module
    
    
program play

    use MyPrec
    use MyHmlt
    implicit none
    character :: a
    
    write(*,*)"输入运算符号"
    read(*,*) a
    t1 = 1.454578787_rp
    t2 = 5.9897854685_rp
    if(a == '+') then
        p0 = 0._rp
        p0 = t1 + t2
        call jiaxiangduiwucha(t1, t2)
        call youxiaoshuzi(p0, er1)
        call jiafa(p0)
    end if
    
    
end program
    
!subroutine guiling
    
!    use MyPrec
!    use MyHmlt 
!    implicit none

!    t1 = 0._rp
!    t2 = 0._rp
!    er1 = 0._rp
!    er2 = 0._rp
!    er3 = 0._rp
!    er4 = 0._rp
!    p0 = 0._rp
!    rp2 = 0
!    n0 = 0

!end subroutine   
    
subroutine jiafa   
  
    use MyPrec
    use MyHmlt
    implicit none
    real(rp) :: x1
    integer :: i, m, n
    character(0:99) :: string, jieguo
    
    write(x1,"(E100.50)") p0
    write(string,"(I100)") x1
    
    m = 0
    n = 0
    do i = 0, 99, 1
        if(n == n0) then
            exit
        end if
        if(string(i) /= ' ') then
            jieguo(m) = string(i)
            m = m + 1
            if(string(i) == '.') then
                cycle
            end if
            n = n + 1
        end if
    end do
    
    open(30, file = "jiafa.txt", access = "append")
    write(30, *) jieguo
    close(30)
     
end subroutine
    
subroutine jiaxiangduiwucha

    use MyPrec
    use MyHmlt
    implicit none
    integer :: m1, m2
    
    m1 = floor(log10(t1))
    m2 = floor(log10(t2))
    er1 = (0.5*10**(m1+1-50) + 0.5*10**(m2+1-50))/(abs(t1+t2))
    
    write(*,*) "t1+t2的相对误差为 ", er1

end subroutine
    
subroutine jianxiangduiwucha

    use MyPrec
    use MyHmlt
    implicit none
    integer :: m1, m2
    
    m1 = floor(log10(t1))
    m2 = floor(log10(t2))
    er2 = (0.5*10**(m1+1-50) + 0.5*10**(m2+1-50))/(abs(t1-t2))
    
    write(*,*) "t1-t2的相对误差为 ", er2
    
end subroutine
    
subroutine chengxiangduiwucha
    
    use MyPrec
    use MyHmlt
    implicit none
    integer :: a1, a2
    
    a1 = FLOOR(abs(t1*10**(-LOG10(t1))))
    a2 = FLOOR(abs(t1*10**(-LOG10(t2))))
    er3 = 10**(-49)/(2*a1) + 10**(-49)/(2*a2)
    
    write(*,*) "t1*t2的相对误差为 ", er3
    
end subroutine    

subroutine chuxiangduiwucha
    
    use MyPrec
    use MyHmlt
    implicit none
    integer :: a1, a2
    
    a1 = FLOOR(abs(t1*10**(-LOG10(t1))))
    a2 = FLOOR(abs(t2*10**(-LOG10(t2))))
    er4 = abs(10**(-49)/(2*a1) - 10**(-49)/(2*a2))
    
    write(*,*) "t1/t2的相对误差为 ", er4
 
end subroutine    

subroutine youxiaoshuzi(er0)

    use MyPrec
    use MyHmlt
    implicit none
    real(rp) :: er0
    integer :: a1
    
    a1 = FLOOR(LOG10(p0))
    a1 = FLOOR(abs(p0)*10**(-a1))
    n0 = -FLOOR(LOG10(2*er0*(a1+1)))
    
    write(*,*) "有效数字为", n0
    
end subroutine




6b68c46fb7499ab3f742617dad474a4.png
回复

使用道具 举报

21

帖子

0

主题

0

精华

熟手

Vim

F 币
244 元
贡献
121 点
发表于 2020-11-15 12:40:56 | 显示全部楼层
第一步不对
[Fortran] 纯文本查看 复制代码
module MyPrec
    implicit none

    integer, parameter :: rp = selected_real_kind(p = 50)!保存50位浮点数

end module

精度最大是33(real(kind=16),没有50精度的

3

帖子

2

主题

0

精华

入门

F 币
39 元
贡献
19 点
 楼主| 发表于 2020-11-15 19:06:41 来自移动端 | 显示全部楼层
谢谢明白了
您需要登录后才可以回帖 登录 | 极速注册

本版积分规则

QQ|捐赠本站|Archiver|关于我们 About Us|群聊|Fcode

GMT+8, 2020-12-5 09:41

Powered by Discuz! X3.2

© 2001-2017 Comsenz Inc.

快速回复 返回顶部 返回列表