zoukankan      html  css  js  c++  java
  • 查找素数Eratosthenes筛法的mpi程序

      思路:

      只保留奇数

      (1)由输入的整数n确定存储奇数(不包括1)的数组大小:

           n=(n%2==0)?(n/2-1):((n-1)/2);//n为存储奇数的数组大小,不包括基数1

      (2)由数组大小n、进程号id和进程数p,确定每个进程负责的基数数组的第一个数、最后一个数和数组维度:

             low_value = 3 + 2*(id*(n)/p);//进程的第一个数

         high_value = 3 + 2*((id+1)*(n)/p-1);//进程的最后一个数

         size = (high_value - low_value)/2 + 1;  //进程处理的数组大小

      (3)寻找奇数的第一个倍数的下标,经过反复思考,有如下规律:

         if (prime * prime > low_value)

            first = (prime * prime - low_value)/2;

        else {

            if (!(low_value % prime))

                         first = 0;

            else

                         first=((prime-low_value%prime)%2==0)?((prime-low_value%prime)/2):((prime-low_value%prime+prime)/2);

        }

      code:

      1 #include "mpi.h"
      2 #include <math.h>
      3 #include <stdio.h>
      4 #define MIN(a,b)  ((a)<(b)?(a):(b))
      5 
      6 int main (int argc, char *argv[])
      7 {
      8    int    count;        /* Local prime count */
      9    double elapsed_time; /* Parallel execution time */
     10    int    first;        /* Index of first multiple */
     11    int    global_count; /* Global prime count */
     12    int    high_value;   /* Highest value on this proc */
     13    int    i;
     14    int    id;           /* Process ID number */
     15    int    index;        /* Index of current prime */
     16    int    low_value;    /* Lowest value on this proc */
     17    char  *marked;       /* Portion of 2,...,'n' */
     18    int    n,m;            /* Sieving from 2, ..., 'n' */
     19    int    p;            /* Number of processes */
     20    int    proc0_size;   /* Size of proc 0's subarray */
     21    int    prime;        /* Current prime */
     22    int    size;         /* Elements in 'marked' */
     23 
     24    MPI_Init (&argc, &argv);
     25 
     26    /* Start the timer */
     27 
     28    MPI_Comm_rank (MPI_COMM_WORLD, &id);
     29    MPI_Comm_size (MPI_COMM_WORLD, &p);
     30    MPI_Barrier(MPI_COMM_WORLD);
     31    elapsed_time = -MPI_Wtime();
     32 
     33    if (argc != 2) {
     34       if (!id) printf ("Command line: %s <m>\n", argv[0]);
     35       MPI_Finalize();
     36       exit (1);
     37     }
     38 
     39    n = atoi(argv[1]);
     40    m=n;//
     41    n=(n%2==0)?(n/2-1):((n-1)/2);//将输入的整数n转换为存储奇数的数组大小,不包括奇数1
     42    //if (!id) printf ("Number of odd integers:%d    Maximum value of odd integers:%d\n",n+1,3+2*(n-1));
     43    if (n==0) {//输入2时,输出1 prime,结束
     44     if (!id) printf ("There are 1 prime less than or equal to %d\n",m);
     45       MPI_Finalize();
     46       exit (1);
     47     }
     48    /* Figure out this process's share of the array, as
     49       well as the integers represented by the first and
     50       last array elements */
     51 
     52    low_value = 3 + 2*(id*(n)/p);//进程的第一个数
     53    high_value = 3 + 2*((id+1)*(n)/p-1);//进程的最后一个数
     54    size = (high_value - low_value)/2 + 1;    //进程处理的数组大小
     55 
     56 
     57    /* Bail out if all the primes used for sieving are
     58       not all held by process 0 */
     59 
     60    proc0_size = (n-1)/p;
     61 
     62    if ((3 + 2*(proc0_size-1)) < (int) sqrt((double) (3+2*(n-1)))) {//
     63       if (!id) printf ("Too many processes\n");
     64       MPI_Finalize();
     65       exit (1);
     66     }
     67 
     68    /* Allocate this process's share of the array. */
     69 
     70    marked = (char *) malloc (size);
     71 
     72    if (marked == NULL) {
     73       printf ("Cannot allocate enough memory\n");
     74       MPI_Finalize();
     75       exit (1);
     76     }
     77 
     78    for (i = 0; i < size; i++) marked[i] = 0;
     79    if (!id) index = 0;
     80    prime = 3;//从素数3开始
     81    do {
     82     //确定奇数的第一个倍数的下标
     83       if (prime * prime > low_value)
     84          first = (prime * prime - low_value)/2;
     85       else {
     86          if (!(low_value % prime)) 
     87         first = 0;
     88          else 
     89         first=((prime-low_value%prime)%2==0)?((prime-low_value%prime)/2):((prime-low_value%prime+prime)/2);
     90         }
     91 
     92       for (i = first; i < size; i += prime)  marked[i] = 1;
     93       if (!id) {
     94          while (marked[++index]);
     95          prime = 2*index + 3;//下一个未被标记的素数
     96         }
     97       if (p > 1) MPI_Bcast (&prime,  1, MPI_INT, 0, MPI_COMM_WORLD);
     98    } while (prime * prime <= 3+2*(n-1));//
     99 
    100    count = 0;
    101    for (i = 0; i < size; i++)
    102       if (!marked[i])  count++; 
    103    if (p > 1) MPI_Reduce (&count, &global_count, 1, MPI_INT, MPI_SUM,
    104       0, MPI_COMM_WORLD);
    105 
    106    /* Stop the timer */
    107 
    108    elapsed_time += MPI_Wtime();
    109 
    110 
    111    /* Print the results */
    112 
    113    if (!id) {
    114       printf ("There are %d primes less than or equal to %d\n",
    115          global_count+1, m);//前面程序是从素数3开始标记,忽略了素数2,所以素数个数要加1
    116       printf ("SIEVE (%d) %10.6f\n", p, elapsed_time);
    117    }
    118    MPI_Finalize ();
    119    return 0;
    120 }
  • 相关阅读:
    Spring IOC 容器源码分析
    OAuth协议简介
    MySQL安装步骤
    C# read and compute the code lines number of cs files based on given directory
    C# StreamWriter log batch message sync and async
    HttpClient SendAsync
    WebRequest, WebRequest.Create GetResponse() GetResponseStream()
    C# FileSystemWatcher
    Access Volumn via extern and invoke win 32 dll
    Change file readonly property File.SetAttribute and new FileInfo readonly property
  • 原文地址:https://www.cnblogs.com/LCcnblogs/p/6015096.html
Copyright © 2011-2022 走看看