《并行程序设计导论》02 openmp

变量的作用域

在oenmp中，变量的作用域涉及在parallel块中能够访问该变量的线程集合。一个能够被线程组中的所有线程访问的变量拥有共享作用域，而一个只能被单个线程访问的变量拥有私有作用域。
在trap函数中被每个线程使用的变量在线程的栈中分配，因此变量有私有作用域。在main函数中声明的变量对于所有线程组中被parallel指令启动的线程都是可以访问的。因此在parallel块之前被生命的变量的缺省作用域是共享的。
总之，在parallel指令前已经被声明的变量，拥有在线程组中所有线程间的共享作用域，而在块中声明的变量中有私有作用域。另外，在parallel块开始处的共享变量的值，与该变量在parallel块之前的值一样。在parallel块完成后，该变量的值是块结束的值。

归约子句

global_result=0.0;
#pragma omp parallel num_threads(thread_count)
{
    double my_result=0.0;
    my_result+=local_trap(a,b,n);
    //使用my_result可以使得对函数的调用在临界区之外，这样各个线程可以同时调用
 #pragma omp critical 
    global_result+=my_result; //要是把函数放这里，会强制各个线程顺序执行函数调用
}

openmp提供了一个更清晰的方法来避免local_trap的串行执行：将global_result定义为一个归约变量。归约操作符是一个二元操作，规约就是将相同的归约操作符重复地应用到操作数序列来得到一个结果的计算。另外，所有操作的中间结果存储在同一个变量中：归约变量。

global_result=0.0;
#pragma omp parallel num_threads(thread_count)\\预处理指令一行写不完就加一个\
    reduction(+: global_result)
gloabl_result+=local_trap(a,b,n);    
/*
这么来理解
把上面的程序写成
global_result=0.0;
#pragma omp parallel num_threads(thread_count)\\预处理指令一行写不完就加一个\
{
#   pragma omp critical 
    gloabl_result+=local_trap(a,b,n);//这样的话对local_trap的调用一次只能被一个线程执行，所以这就相当于强制各个线程顺序执行梯形积分法
    //所以上面是那么解决的，但是用归约子句就不用改了
}

*/   

代码明确了global_result是一个归约变量，加号（+）指示归约操作符是加法。openmp为每个线程有效地创建了一个私有变量，运行时系统在这个私有变量中存储每个线程的结果。openmp也创建了一个临界区，并且在临界区中，将存储在私有变量中的值进行相加。
reduction子句的语法

reduction(<operator>: <variable list>)

在c语言中，操作符可能是+,*,-,&,|,^,||中的任意一个，但是减法可能会有问题，因为减法不支持交换律和结合律。
注意，如果归约变量是float和double，当使用不同数量的线程时，结果可能会有些不同。这是因为浮点数不满足结合律。
当一个变量被包含在一个reduction子句中时，变量本身是共享的。然而线程组中的每个线程都创建自己的私有变量。在parallel块里，每当一个线程执行涉及这个变量的语句时，它使用的其实是私有变量。当parallel块结束后，私有变量中的值被整合到一个共享变量中。
最后注意，reduction子句创建的私有变量初始化为相同的值。

parallel for指令

h=(b-a)/n;
approx=(f(a)+f(b))/2.0;
#pragma omp parallel for num_threads(thread_count)\
    reduction (+:approx)
for(i=1;i<=n-1;i++)
    approx+=f(a+i*h);//approx作为归约变量是有必要的
approx=h*approx;

parallel for指令生成一组线程来执行后面的结构化代码块。然而，在parallel for指令之后的结构化块必须是for循环。另外，运用parallel for指令，系统通过在线程间划分循环迭代来并行化for循环。因此，parallel for指令与parallel指令非常不同，因为在parallel指令之前的块，一般来说其工作必须由线程本身在线程之间划分。
在一个已经被parallel for 指令并行化的for循环中，线程间的缺省划分方式时有系统决定的，大部分系统会粗略进行块划分。
在parallel指令中，所有变量的缺省作用于是共享的。但在parallel for指令中，如果循环变量i是共享的，那么变量更新++i也会是一个无保护的临界区。因此，在每一个被parallel for并行化的循环中，循环变量的缺省作用域都是私有的，在我们的代码中，每个线程组拥有它自己的i的副本。

for指令的警告

首先openmp只能并行化for循环，只能并行化那些可以在如下情况下确定迭代次数的for循环：
由for语句本身（for(…;…;…）)来确定
在执行之前确定。
例如

for(; ;)
{
    ...
}//无限循环不可以并行化
for(i=0;i<n;i++)
{
    if(...)break;
}//不能被并行化，因为迭代次数不能只从for语句来确定。这个for循环也不是一个结构化块，因为break添加了另一个从循环退出的出口

运行时系统必须能够在执行前确定迭代的数量，但是唯一例外的是：在循环体中可以有一个exit调用

for指令数据依赖性

在循环中，迭代中的计算依赖于一个或多个先前的迭代结果。这样使用for指令的结果是不可预计的。
两个要点
(1）openmp编译器不检查被parallel for指令并行化的循环所包含的迭代间的依赖关系，而是由程序员来识别这些依赖
(2) 一个或多个迭代结果依赖于其它迭代的循环，一般不能被openmp正确地并行化。

f[i]=f[i-1]+f[i-2],此时f之间的依赖关系称为数据依赖或循环依赖

寻找循环依赖

当我们试图使用parallel for指令时，首先就要注意循环依赖。
重点关注在一个迭代中被读或被写，而在另一个迭代中被写的变量

pai值的估计

在串行代码

double factor=1.0;
double sum=0.0;
for(k=0;k<n;k++)
{
    sum+=factor/(2*k+1);
    factor=-factor;
}
pi_approx=4.0*sum;

直观的openmp版本

double factor=1.0;
double sum=0.0;
#pragma omp parallel for num_threads(thread_count) \
    reduction(+:sum)
for(k=0;k<n;k++)
{
    sum+=factor/(2*k+1);
    factor=-factor;
}
pi_approx=4.0*sum;

发现上面的第七行和第八行存在循环依赖。
改进的openmp

double factor=1.0;
double sum=0.0;
#pragma omp parallel for num_threads(thread_count) \
    reduction(+:sum)
for(k=0;k<n;k++)
{
   factor=(k%2==0)?1.0:-1.0
   sum+=factor/(2*k+1);
}
pi_approx=4.0*sum;

但上面的代码仍有错误。
在一个已经被parallel for指令并行化的块中，缺省情况下任何在循环前声明的变量（唯一的例外是循环变量）在线程之间都是共享的。因此factor被共享。例如，线程0可能会给他赋值1，但是在他用这个值更新sum值前，线程1可能给他赋值-1了。因此，除了消除计算factor时的循环依赖外，我们还需要保证每个线程有它自己的factor副本。就是说，为了使代码正确，我们需要保证factor有私有作用域，通过添加一个private子句到parallel指令来实现

double sum=0.0;
#pragma omp parallel for num_threads(thread_count) \
    reduction(+:sum) private(factor)
for(k=0;k<n;k++)
{
   factor=(k%2==0)?1.0:-1.0
   sum+=factor/(2*k+1);
}
pi_approx=4.0*sum;

在parallel子句内列举的变量，在每个线程上都有一个私有副本被创建。因此一个线程对factor的更新不会影响另一个线程的factor值。
要更记住的是，一个私有作用域的变量的值在parallel块或者parallel for块的开始处是未指定的。他的值在parallel或者parallel for块完成之后也是未指定的。

int x=5;//这里的x跟下面的x不是一回事
#pragma omp parallel num_threads(thread_count) private(x)
{
    int my_rank=omp_get_thread_num();
    printf("Thread %d>before initialization,x=%d\n",my_rank,x);
    //这个printf语句的输出是非确定的，因为在它被显式初始化之前就打印了私有变量x
    x=2*my_rank+2;
    printf("Thread %d>after initialization,x=%d\n",my_rank,x);
}
printf("after parallel block,x=%d\n",x);
//这里的printf也是非确定的，因为他在parallel块完成之后打印x

关于作用域的更多问题

openmp提供了一个子句default，该子句显式地要求我们这样做。如果我们1添加子句

default(none)

到parallel或parallel for指令中，那么编译器将要求我们明确在这个块中使用的每个变量和已经在块之外声明的变量的作用域。（在一个块中声明的变量都是私有的，因为他们会被分配给线程的栈）
例如使用default子句，对pi的计算如下

double sum=0.0;
#pragma omp parallel for num_threads(thread_count) \
    default(none) reduction(+:sum) private(k,factor) shared(n)
for(k=0;k<n;k++)
{
   factor=(k%2==0)?1.0:-1.0
   sum+=factor/(2*k+1);
}
pi_approx=4.0*sum;

在这个例子中，for循环中使用4个变量。由于default子句，我们需要明确每个变量的作用域。
sum是一个归约变量，同时拥有私有和共享作用域的属性。
factor和k应该拥有私有作用域
从未在parallel块中更行的变量n，应该被安全地共享。共享变量在块内具有parallel块之前同样的值，在块之后的值与块内的最后一个值相同。