OpenMP C++算法:求最小值、最大值、中位数和平均数

28

我在Google上搜索一些简单的OpenMP算法页面。可能有一个示例可以从巨大的数据数组中计算最小值、最大值、中位数和平均值,但我找不到。

至少我通常会尝试将数组分成每个核心的一个块,然后进行一些边界计算,以便得到完整数组的结果。

我只是不想重复发明轮子。


附加说明: 我知道有成千上万的例子可以使用简单的规约来工作。 例如:计算π。

const int num_steps = 100000; 
double x, sum = 0.0; 
const double step = 1.0/double(num_steps); 
#pragma omp parallel for reduction(+:sum) private(x) 
for (int i=1;i<= num_steps; i++){ 
  x = double(i-0.5)*step; 
  sum += 4.0/(1.0+x*x); 
} 
const double pi = step * sum;

但是,当这些类型的算法不可用时,几乎没有其他例子可以用来缩小算法。


是的,我同意,在 OpenMP 上很难找到教程和示例... http://openmp.blogspot.com 我昨天偶然发现了这个,可能会有用,所以想在这里分享一下。 - anshu
4个回答

25

OpenMP(至少2.0版本)支持一些简单操作的归约,但不支持最大值和最小值。

在下面的例子中,使用reduction子句来进行求和,并使用critical区域使用线程本地变量更新共享变量,以避免冲突。

#include <iostream>
#include <cmath>

int main()
{
  double sum = 0;
  uint64_t ii;
  uint64_t maxii = 0;
  uint64_t maxii_shared = 0;
#pragma omp parallel shared(maxii_shared) private(ii) firstprivate(maxii)
  {
#pragma omp for reduction(+:sum) nowait
    for(ii=0; ii<10000000000; ++ii)
      {
        sum += std::pow((double)ii, 2.0);
        if(ii > maxii) maxii = ii;
      }
#pragma omp critical 
    {
      if(maxii > maxii_shared) maxii_shared = maxii;
    }
  }
  std::cerr << "Sum: " << sum << " (" << maxii_shared << ")" << std::endl;
}

编辑:更简洁的实现:

#include <cmath>
#include <limits>
#include <vector>
#include <iostream>
#include <algorithm>
#include <tr1/random>

// sum the elements of v
double sum(const std::vector<double>& v)
{
  double sum = 0.0;
#pragma omp parallel for reduction(+:sum)
  for(size_t ii=0; ii< v.size(); ++ii)
    {
      sum += v[ii];
    }
  return sum;
}

// extract the minimum of v
double min(const std::vector<double>& v)
{
  double shared_min;
#pragma omp parallel 
  {
    double min = std::numeric_limits<double>::max();
#pragma omp for nowait
    for(size_t ii=0; ii<v.size(); ++ii)
      {
        min = std::min(v[ii], min);
      }
#pragma omp critical 
    {
      shared_min = std::min(shared_min, min);
    }
  }
  return shared_min;
}

// generate a random vector and use sum and min functions.
int main()
{
  using namespace std;
  using namespace std::tr1;

  std::tr1::mt19937 engine(time(0));
  std::tr1::uniform_real<> unigen(-1000.0,1000.0);
  std::tr1::variate_generator<std::tr1::mt19937, 
    std::tr1::uniform_real<> >gen(engine, unigen);

  std::vector<double> random(1000000);
  std::generate(random.begin(), random.end(), gen);

  cout << "Sum: " << sum(random) << " Mean:" << sum(random)/random.size()
       << " Min:" << min(random) << endl;
}

3
自OpenMP 3.1版本起,它支持最小值和最大值约简操作,在gcc 4.7及以上版本中可用。 - Sameer

11

从OpenMP 3.1开始,可以通过reduction子句实现for min和max,在此链接中有详细的示例。


希望有人能够实现3.1版本,这将使生活变得更加轻松。 - Totonga
您可以从 GCC 4.7 及以上版本中找到 OpenMP 3.1。 - Mahesh

7

OpenMP不支持这些约简操作。考虑使用Intel Threading Building Blocks的parallel_reduce算法,您可以实现任意算法。

这里有一个例子。它使用部分结果的求和。您可以实现任何想要的函数。

#include <stdio.h>
#include <tbb/blocked_range.h>
#include <tbb/parallel_reduce.h>
#include <tbb/task_scheduler_init.h>


///////////////////////////////////////////////////////////////////////////////


class PiCalculation
{
private:
    long num_steps;
    double step;

public:

    // Pi partial value
    double pi;

    // Calculate partial value
    void operator () (const tbb::blocked_range<long> &r) 
    {
        double sum = 0.0;

        long end = r.end();

        for (int i = r.begin(); i != end; i++)
        {
            double x = (i + 0.5) * step;
            sum += 4.0/(1.0 + x * x);
        }

        pi += sum * step;
    }

    // Combine results. Here you can implement any functions
    void join(PiCalculation &p)
    {
        pi += p.pi;
    }

    PiCalculation(PiCalculation &p, tbb::split)
    {
        pi = 0.0;
        num_steps = p.num_steps;
        step = p.step;
    }

    PiCalculation(long steps)
    {
        pi = 0.0;
        num_steps = steps;
        step = 1./(double)num_steps;
    }
};


///////////////////////////////////////////////////////////////////////////////


int main()
{
    tbb::task_scheduler_init init;

    const long steps = 100000000;

    PiCalculation pi(steps);

    tbb::parallel_reduce(tbb::blocked_range<long>(0, steps, 1000000), pi);

    printf ("Pi is %3.20f\n", pi.pi);
}

请查看此链接以获取更多减少算法。 http://cache-www.intel.com/cd/00/00/30/11/301132_301132.pdf#page=19 请仔细阅读第3.3.1段。有一个在数组中找到最小值的示例。


1
这种归约在OpenMP中非常容易。而且有一个巨大的优势,即代码从串行到多线程不会有任何区别。但它最终只能实现简单的归约功能。const int num_steps = 100000; double x, sum = 0.0; const double step = 1.0 / double(num_steps); #pragma omp parallel for reduction(+:sum) private(x) for (int i = 1; i <= num_steps; i++) { x = double(i - 0.5) * step; sum += 4.0 / (1.0 + x * x); } const double pi = step * sum; - Totonga
1
亲爱的Totonga!OpenMP在归约函数上受到限制,只能使用少数算术运算符:+、-、*、/。而在TBB中,您可以实现任意归约函数。这是它的优势。 - Vladimir Obrizan

3

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接