Fortran Coder

查看: 141|回复: 3

[求助] 一个小的程序问题(误差积累)

[复制链接]

9

帖子

5

主题

0

精华

入门

F 币
57 元
贡献
35 点
发表于 2021-3-30 17:57:40 | 显示全部楼层 |阅读模式
[Fortran] 纯文本查看 复制代码
program main
implicit none


integer i,j,k,maxx,maxy,maxz,i0,jb,kc,pp,pp1,ct
real dx,dy,dz,ww,ww1,rrb1,rrb
real,allocatable ::dz_1(:)

        dz=0.05
        maxx=74;maxy=74;maxz=74
        i0=23;jb=52;kc=52
        ct=1.0

        allocate(dz_1(0:maxz))
        dz_1=0
        do k=1,maxz
        dz_1(k)=dz_1(k-1)+dz
        enddo

        DO k = 1,maxz-1
        DO j = 1,maxy-1
          DO i = 1,maxx



!            rrb=(i-0.5)*dx*st*cp+  &
!             j*dy*st*sp+    &
!             k*dz*ct               !(6-7-13)
            rrb= k*dz*ct   
!            rrb=(dx_1(i-1)+(dx_1(i)-dx_1(i-1))/2.)*st*cp+  &
!             dy_1(j)*st*sp+             &
!             dz_1(k)*ct     !+ dz/1.e+10
          rrb1=dz_1(k)*ct

              pp=int(rrb/dz)
          pp1=int(rrb1/dz)
          ww=abs(mod(rrb,dz))
          ww1=abs(mod(rrb1,dz))
!          ww=rrb/dz-pp
          if(i.eq.i0.and.j.eq.jb.and.k.eq.kc) then
          print*,rrb,pp,ww
          print*,rrb1,pp1,ww1
          pause
          endif

          if(abs(1-(pp+1)/rrb).lt.1.e-5)then
          pp=pp+1
          endif   
          if(ww.lt.1.e-6.or.abs(1-ww/dz).lt.1.e-5)then
          ww=0
          else
          ww=ww/dz
          endif

           ENDDO
          ENDDO
        ENDDO

end

    不知道rrb,pp,ww和rrb1,pp1,ww1为什么不一样,试过单双精度,对结果没有改变。



回复

使用道具 举报

1593

帖子

12

主题

5

精华

论坛跑堂

臭石头雪球

F 币
1771 元
贡献
1135 点

美女勋章热心勋章星光勋章新人勋章贡献勋章管理勋章帅哥勋章爱心勋章规矩勋章元老勋章水王勋章

发表于 2021-3-30 20:18:26 | 显示全部楼层
因为 dz_1 是累加起来的。累加会让误差积累。
你以为 dz 是 0.05 ,其实不是。它可能是 0.0500000007451
累加起来,误差就会越来越大。

而 rrb= k*dz*ct 因为用的乘法,虽然也有误差,但不会积累。

如果你难以理解为什么计算机用二进制不能精确存储 0.05
你可以联想一下,十进制不能精确存储 1/3,但三进制却可以。

9

帖子

5

主题

0

精华

入门

F 币
57 元
贡献
35 点
 楼主| 发表于 2021-3-30 20:37:00 | 显示全部楼层
fcode 发表于 2021-3-30 20:18
因为 dz_1 是累加起来的。累加会让误差积累。
你以为 dz 是 0.05 ,其实不是。它可能是 0.0500000007451
累 ...

感谢女神指导

118

帖子

2

主题

0

精华

大师

F 币
914 元
贡献
463 点

规矩勋章

发表于 2021-3-30 21:25:17 | 显示全部楼层

是大神,不能以貌识人。
您需要登录后才可以回帖 登录 | 极速注册

本版积分规则

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

GMT+8, 2021-4-19 08:49

Powered by Discuz! X3.2

© 2001-2017 Comsenz Inc.

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