使用 pthreads 和 OpenMp 在 C 中并行化 Simpson 的方法

Parallelizing Simpson's Method in C using pthreads and OpenMp

我对整个多处理和并行化世界还很陌生。我目前正在研究一种算法来实现 Simpson 的方法和 Simpson 3/8。我已经在 C 中以串行形式实现了该算法。

我需要一些帮助来使用 OpenMp 和 Pthreads 并行化这两种算法。欢迎和赞赏任何建议。

//Simpson.c
#include <stdio.h>
#include <stdlib.h>

double function(double x) //función a la cual se le aplicará el método
{
  return 0.5*x;    
}

int main(int argc, char*argv[])
{
  double resultado = 0.0; //resultado de la integral tras aplicar el método de simpson
  double a = 0.0;  //límite inferior de la integral
  double b = 0.0;  //límite superior de la integral
  int n = 0; //cantidad de particiones entre los intervalos, tiene que ser un número par
  int i = 0; //variable para poder iterar sobre las diferentes particiones
  double h = 0.0; //variable para guardar el incremento que habrá entre cada valor de "x"
  double calc = 0.0; //variable intermedia para ir acumulando el resultado de evaluar las diferentes "x" en la función
  
  //variables para almacenar los diferentes valores de "x", es decir los límites inferiores y superiores y el valor intermedio de cada partición  
  double x0 = 0.0;
  double x1 = 0.0;
  double x2 = 0.0;  
    
  printf("Introduce el límite inferior: ");
  scanf("%lf", &a);
  printf("Introduce el límite superior: ");
  scanf("%lf", &b);
  printf("Introduce la cantidad de particiones (número par): ");
  scanf("%d", &n);  
    
  h = (b-a)/n; //se calcula el incremento para poder iterar sobre todas las particiones
    
      //Para el cálculo de la integral se utilizan xj, x(j+1) y x(j+2) en el subintervalo [xj, x(j+2)] para cada i = 1, 3, 5, ..., n
    
  for (i=0; i<=n; i+=2)
  {   
    x0 = a + h*i; //límite inferior
    calc += function(x0); //se evalua la x en la función
    
    x1 = a + h*(i+1); //valor intermedio
    calc += 4*function(x1);
      
    x2 = a + h*(i+2); //límite superior
    calc += function(x2);
      
    calc *= ((x2-x0)/6) ; //se vuelve a calcular el incremento entre nuestros límites, y se divide entre 6
      
    resultado += calc; //variable para ir acumulando los cálculos intermedios de cada "x"
  }
    
  printf("La aproximación a la integral es: %f \n", resultado);

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

// función a integrar
#define f(x) (0.5*x)

int main()
{
  float b; //Limite superior
  float a; //Limite inferior
  float resultado=0.0; 
  float dx;
  float k;
  int i; 
  int N;  //numero de intervalos o iteraciones

  
 // Entradas
 printf("Dame el limite inferior: ");
 scanf("%f", &a);
 printf("Dame el limite superior: ");
 scanf("%f", &b);
 printf("Número de iteraciones/intervalos (número par): ");
 scanf("%d", &N);

  //Calculando delta de x
  dx = (b - a)/N;
  
  //Sumas de la integracion
  resultado = f(a) + f(b); //primera suma
  
  for(i=1; i<= N-1; i++)
  {
    k = a + i*dx; //Suma del intervalo
    
    if(i%3 == 0)
    {
      resultado = resultado + 2 * f(k);
    }
    else
    {
      resultado = resultado + 3 * f(k);
    }
  }
  
  resultado = resultado * dx*3/8;
  
  printf("\nEl valor de la integral es: %f \n", resultado);
  printf("\n");
  
  return 0;
}

我的第一个建议是首先阅读大量有关共享内存范例及其问题的文章。然后阅读 Pthreads 和 OpenMP,然后查看具体的 OpenMP 和 Pthread 并行化,以及它们如何处理共享内存范例的一些问题。如果您必须同时学习 Pthread 和 OpenMP,在我看来,如果我是您,我会先研究 Pthread,然后再研究 OpenMP。后者有助于解决前者编码的一些繁琐细节和陷阱。

我会尽量给你一个非常笼统的(不是很详细公式关于如何处理代码的并行化,我将使用您的第二个代码和 OpenMP 作为示例。

如果您选择并行化算法,首先您需要查看值得并行化的代码部分(即,并行化的开销被它带来的加速)。因为这种分析是基本的,但是,对于像您这样的较小代码,很容易弄清楚,它通常是循环,例如:

for(i=1; i<= N-1; i++)
  {
    k = a + i*dx; //Suma del intervalo
    
    if(i%3 == 0)
    {
      resultado = resultado + 2 * f(k);
    }
    else
    {
      resultado = resultado + 3 * f(k);
    }
  }

要并行化此代码,我应该使用哪个构造函数 tasksparallel for ?您可以在 SO thread 上查看该主题的精彩答案。在您的代码中,很明显我们应该使用 parallel for,即 #pragma omp parallel for)。因此,我们告诉 OpenMP 在线程之间划分该循环的迭代。

现在你应该考虑并行区域内使用的变量在线程之间应该是private还是shared。为此,openMP 提供了 privatefirstprivatelastprivateshared 等构造函数。可以找到关于该主题的另一个很棒的 SO 线程 here。要对此进行评估,您 确实 需要了解 shared memory 范例以及 OpenMP 如何处理它。

然后查看潜在的竞争条件、变量之间的相互依赖性等。 OpenMP 提供了 criticalatomic operationsreductions 等构造函数来解决这些问题(其他出色的 SO 线程 atomic vs critical and )。 TL;DR :您应该选择能够在不影响正确性的情况下提供最佳性能的解决方案。

在您的情况下,变量“resultado”的 reduction 显然是最佳选择:

 #pragma omp parallel reduction ( +: resultado) 
 for(i=1; i<= N-1; i++)
  {
    k = a + i*dx; //Suma del intervalo
    
    ... 
  }

保证代码的正确性后,您可以查看负载不平衡问题(每个线程执行的工作之间的差异)。为了解决这个问题,openMP 提供了并行 for 调度策略,例如 dynamicguided(另一个 SO thread)以及调整块大小的可能性那些发行版。

在您的代码中,线程将执行相同数量的并行工作,因此无需使用 dynamicguided 计划,您可以选择 static 计划. static 优于其他方法的好处是它不会引入必须在运行时协调线程之间的任务分配的开销。

要从您的代码 运行 所在的体系结构中提取最大性能,还有很多东西需要研究。但就目前而言,我认为这是一个足够好的开始。