Fortran Coder

查看: 13605|回复: 2
打印 上一主题 下一主题

[求助] 取e小数点后20个10位数的质数,求debug!

[复制链接]

7

帖子

2

主题

0

精华

入门

F 币
42 元
贡献
25 点

规矩勋章

跳转到指定楼层
楼主
发表于 2014-9-26 21:21:31 | 只看该作者 回帖奖励 |倒序浏览 |阅读模式
运行后发现 第一部分求e的后10000位可以得到结果
第二部分就卡住了

我觉得可能是定义或者mod的问题,因为求的是10位数的质数,用tmp(10位数组)除数得到的余数不准确附了代码
求指明方向!!!!拜托了
[Fortran] 纯文本查看 复制代码
Program jiayi

  Implicit Real *8(A, C-H, O-S, U-Z)
  Integer, Parameter :: nk = 10000
  Dimension :: ns(nk), na(nk), nt(nk), prime(20)
  Integer :: i, j, b, it, l

  Do i = 1, nk
    na(i) = 0
    ns(i) = 0
    nt(i) = 0
  End Do

  na(1) = 1
  ns(1) = 2

  n = 1
  Do j = 1, 5000
    n = n + 1

    Do i = 2, nk
      nd = mod(na(i-1), n)*10 + nt(i)
      nt(i) = nd/n
      na(i) = mod(nd, n)
    End Do

    Do i = 1, nk
      na(i) = nt(i)
      ns(i) = ns(i) + nt(i)
    End Do

    Do i = nk, 2, -1
      ns(i-1) = ns(i-1) + ns(i)/10
      ns(i) = mod(ns(i), 10)
    End Do

  End Do

  Print *, 'e=2.'
  Write (*, '(8x,77i1)')(ns(i), i=2, 77)
  Write (*, '(80i1)')(ns(i), i=78, nk)

  b = 1

loop_utmost: Do i = 2, nk - 10
    tmp = 0
    If (ns(i)/=0) Then
loop_tmp: Do j = i, i + 9 !loop 10 times to get the 10 digit number, and save it as "tmp"
        tmp = tmp*10 + ns(j)
      End Do loop_tmp

      l = sqrt(tmp+1) !choose the squareroot of the number itself as the biggest divisor to save memory

loop_check: Do it = 2, l !To check whether tmp can be devided by any numbers in the range or not
        If (mod(tmp,real(it))==0) Then
          Cycle loop_utmost

        Else If (it<=l) Then
          Cycle loop_check
        Else
          b = b + 1
          If (b<=20) Then
            prime(b) = tmp
            Cycle loop_utmost
          Else
            Exit
          End If
        End If
      End Do loop_check
    End If
  End Do loop_utmost

  Print *, 'The ten digital primes in e:'
  Write (*, '(i5,f15.1)')(prime(i), i=1, b)

End Program jiayi

分享到:  微信微信
收藏收藏 点赞点赞 点踩点踩

954

帖子

0

主题

0

精华

大师

F 币
184 元
贡献
75 点

规矩勋章元老勋章新人勋章水王勋章热心勋章

QQ
沙发
发表于 2014-9-27 08:44:14 | 只看该作者
本帖最后由 vvt 于 2014-9-27 08:50 编辑

有多处问题,我已帮你修改并标出。
计算结果为
The ten digital primes in e:
   7427466391.0
   7413596629.0
   6059563073.0
   3490763233.0
   2988075319.0
   1573834187.0
   7021540891.0
   5408914993.0
   6480016847.0
   9920695517.0
   1838606261.0
   6062613313.0
   3845830007.0
   1692836819.0
   4425056953.0
   2505695369.0
   5490598793.0
   1782154249.0
   8215424999.0
   9229576351.0

[Fortran] 纯文本查看 复制代码
Program jiayi
  Implicit None !// 推荐用 Implicit None
  Integer, Parameter :: nk = 10000
  Integer :: ns(nk), na(nk), nt(nk) , n , nd
  Real(Kind=8) :: prime(20) , tmp !// tmp 要定义,很重要
  Integer :: i, j, b, it , l

  Do i = 1, nk
    na(i) = 0
    ns(i) = 0
    nt(i) = 0
  End Do

  na(1) = 1
  ns(1) = 2

  n = 1
  Do j = 1, 5000
    n = n + 1

    Do i = 2, nk
      nd = mod(na(i-1), n)*10 + nt(i)
      nt(i) = nd/n
      na(i) = mod(nd, n)
    End Do

    Do i = 1, nk
      na(i) = nt(i)
      ns(i) = ns(i) + nt(i)
    End Do

    Do i = nk, 2, -1
      ns(i-1) = ns(i-1) + ns(i)/10
      ns(i) = mod(ns(i), 10)
    End Do

  End Do

  Print *, 'e=2.'
  Write (*, '(8x,77i1)')(ns(i), i=2, 77)
  Write (*, '(80i1)')(ns(i), i=78, nk)

  b = 0

loop_utmost: &
  Do i = 2, nk - 10
    tmp = 0.0D0
    If (ns(i)==0) Cycle loop_utmost
    Do j = i, i + 9 !loop 10 times to get the 10 digit number, and save it as "tmp"
      tmp = tmp * 10.0D0 + ns(j)
    End Do
    l = sqrt(tmp+1) !choose the squareroot of the number itself as the biggest divisor to save memory
    Do it = 2, l !To check whether tmp can be devided by any numbers in the range or not
      If (abs(mod(tmp,real(it)))<0.01D0) Cycle loop_utmost !// real 类型不要做相等判断
    End Do !// 此处逻辑上有问题,修改了 Do 和 If 代码
    b = b + 1
    If (b>20) Exit
    prime(b) = tmp
    Cycle loop_utmost
  End Do loop_utmost

  Print *, 'The ten digital primes in e:'
  Write (*, '(f15.1)')(prime(i), i=1, min(b,20)) !// 此处改为 min(b,20)

End Program jiayi


评分

参与人数 1F 币 +2 贡献 +2 收起 理由
岸边的鱼 + 2 + 2 很给力!

查看全部评分

2033

帖子

12

主题

5

精华

论坛跑堂

臭石头雪球

F 币
1641 元
贡献
709 点

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

板凳
发表于 2014-9-27 18:05:50 | 只看该作者
结果好像不止20个,还挺多
您需要登录后才可以回帖 登录 | 极速注册

本版积分规则

捐赠本站|Archiver|关于我们 About Us|小黑屋|Fcode ( 京ICP备18005632-2号 )

GMT+8, 2024-12-23 04:28

Powered by Tencent X3.4

© 2013-2024 Tencent

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