kanch.blog

posted at 2015-11-03 22:32:39 +0000

C语言证明哥德巴赫猜想算法以及算法优化相关问题

上周,我们C语言老师给我们布置了一道小组作业:用C语言编写程序,证明哥德巴赫猜想在4—1000,000之间成立。最后成绩以各组算法在他机器上运行时间为准。特别注意:只能使用《C语言程序设计》前五章的知识,即不能用指针,不能用数组,不能用自定义函数等等…

我们才学C语言没多久,就布置这样的题目?跪了。不过大体想了想,要证明哥德巴赫猜在4—1000,000之间成立,其实还是简单的,关键是如何优化算法,让它的时间最短。

首先我们先来熟悉下哥德巴赫猜想:哥德巴赫猜想主张每个大于等于4的偶数都是哥德巴赫数——可表示成两个素数之和的数。

基本思路

来说说基本解题思路:对于任意给定偶数a,我们要先把它拆分成x和y,并满足x+y=a。然后循环判断x与y都为质数即可。

所以就有了如下代码:

#include<stdio.h>
#include<time.h>
#include<math.h>

int main() {
  clock_t s, f;
  int x = 4, m, n, t = 2, r = 0, T = 0, R = 0;
  s = clock();
  //process start 
  //第一个循环,用来取出4——1000000之间的,所有偶数x 
  for (; x < 1000000; x += 2) {
    m = x / 2;
    n = m;
    R = m;
    //这个循环,用来对取出的偶数x的分解出的2个数做x/2次判断因为m和n变到0和x都需x/2次自减自加 
    for (r = 0; r < R; r++, m--, n++) {
      T = sqrt(n);
      //这个循环用来判断偶数x分解出的2个数m,n是否为质数,因为 n>=m 恒成立,所以我们做sqrt(n)次判断 
      for (t = 2; t < T; t++) {
        //如果n和m中有一个数能够被大于2的整数i整除,则跳出这个for循环,对m,n变换后继续判断 
        if (m % t == 0 || n % t == 0) {
          break;
        }
      }
    }
    //如果循环结束的时候,n!>=x或者m!<=0则说明,这个偶数不满祝哥德巴赫的猜想 
    if ((!(n >= x)) || (!(m <= 0))) {
      printf("%d不满足n", x);
    }
  }
  //process compltet 
  f = clock();
  printf("%.3lf secondsn", (double)(f - s) / CLOCKS_PER_SEC);
  return 0;
}

但是,这个算法在我的电脑上大概需要跑7–9秒才可以证明完毕。所以就引出了如何优化的问题。

优化

我们看到,上面代码对于分解出来的每个因子x和y,判断质数都是从2开始,一直到√x或者√y,所以,这里就有了第一个优化方法:

只判断奇数,即只用奇数去除待判断数字m或n

//第一个循环,用来取出4——1000000之间的,所有偶数x 
for (; x < 1000000; x += 2) {
  m = x / 2;
  n = m;
  R = m;
  //这个循环,用来对取出的偶数x的分解出的2个数做x/2次判断因为m和n变到0和x都需x/2次自减自加 
  for (r = 0; r < R; r++, m--, n++) {
    T = sqrt(n);
    //这个循环用来判断偶数x分解出的2个数m,n是否为质数,因为 n>=m 恒成立,所以我们做sqrt(n)次判断 
    for (t = 3; t < T; t += 2) {
      //如果n和m中有一个数能够被大于2的整数i整除,则跳出这个for循环,对m,n变换后继续判断 
      if (m % t == 0 || n % t == 0) {
        break;
      }
    }
  } //如果循环结束的时候,n!>=x或者m!<=0则说明,这个偶数不满祝哥德巴赫的猜想 
  if ((!(n >= x)) || (!(m <= 0))) {
    printf("%d不满足n", x);
  }
}

其次呢,我们可以把判断

x和y是否为质数中,大的那个数字,放到最内层循环,来减少CPU跨切循环层次数

 T = sqrt(m);
 //首先判断m是否为质数,如果是质数,就进入到下个if块,判断n是否为质数 
 //如果m不是质数,跳出#2,开始新一次#1,并对m--,n++,因为一个不满足,这一对就不满足 
 for (t = 3; t <= T; t += 2) //#2 
 {
   if (m % t == 0) {
     break;
   }
 } //得到了第一个质数m,我们再判断n是否为质数,是的话就跳出#1,不是的话,回到#2继续寻找下一个质数m 
 if (t > T) {
   n = x - m;
   T = sqrt(n); //判断n是否为质数。如果是,就执行下个if块,并跳出#1,回到#0,判断下一个偶数是否满足
   for (t = 3, T += 1; t <= T; t += 2) //#3
   {
     if (n % t == 0) {
       break;
     }
   }
   if (t > T) {
     TT = 0;
     //printf(" %d+%d=%d",n,m,x); 
     continue;
   }
 }
 ++m, ++m;
 goto label;

再者,我们每次循环的过程中都要判断x%i==0,即x能否被i整除,在以上代码中,我们可以看到:每次判断x能被i整除后,我们用的是break语句,跳出这个循环,继续执行循环后的语句(即if块)。但是由于

我们已经知道这个数不是质数。所以我们可以用一个goto语句,直接跳到顶部,改变x和y的值,继续判断。这样便去除了冗余的if判断

做完以上优化后,最终版本如下:

//第一组作品,by Kanch 算法优化参考了:http://m.blog.csdn.net/blog/sunjiajiang/7887724
#include < stdio.h > #include < time.h > #include < math.h >

  int main() {
    clock_t s, f;
    int x = 6, m, n, t = 2, T = 0, TT = 0;
    s = clock();
    //第一个循环,用来取出4——1000000之间的,所有偶数x 
    for (; x < 1000000; ++x, ++x) //#0 { 
      //这个循环,用来对取出的偶数x的分解出的2个数做x/2次判断因为m和n变到0和x都需x/2次自减自加 
      //由于每次判断后,都要对m++,n--,所以,m+n=x 
      m = 3; 
      label_dm:
      T = sqrt(m);
    //首先判断m是否为质数,如果是质数,就进入到下个if块,判断n是否为质数
    //如果m不是质数,跳出#2,开始新一次#1,并对m--,n++,因为一个不满足,这一对就不满足
    for (t = 3; t <= T; t += 2) //#2
    {
      if (m % t == 0) {
        //break;
        ++m, ++m;
        goto label_dm;
      }
    }
    //得到了第一个质数m,我们再判断n是否为质数,是的话就跳出#1,不是的话,回到#2继续寻找下一个质数m
    if (t > T) {
      n = x - m;
      T = sqrt(n);
      //判断n是否为质数。如果是,就执行下个if块,并跳出#1,回到#0,判断下一个偶数是否满足
      for (t = 3, ++T; t <= T; t += 2) //#3
      {
        if (n % t == 0) {
          break;
        }
      }
      if (t > T) {
        TT = 0;
        //printf(" %d+%d=%d",n,m,x);
        continue;
      }
    }
    ++m, ++m;
    goto label_dm;

    //因为每次最外层循环开始前我都会把TT设置为1,这样,只要无法执行到51行,就说明这个偶数不满足哥德巴赫猜想啊
    if (TT) {
      printf("%d不满足n", x);
      TT = 1;
    }
  }
f = clock();
printf("%.6lf secondsn", (double)(f - s) / CLOCKS_PER_SEC);
return 0;

}

<2015-11-5日更新,这里再补上一张算法解释图>

gdbh

<kanch@20151103>


最好的科学上网服务->西游:https://xiyou4you.com/r/?s=2219145



© kanchzl AT kanchz DOT com

last updated on 2019-09-19 10:00:09 +0000