OpenMP

基本概念

OpenMP是一种用于共享内存并行系统的多线程程序设计方案，支持的编程语言包括C、C++和Fortran。OpenMP提供了对并行算法的高层抽象描述，特别适合在多核CPU机器上的并行程序设计。编译器根据程序中添加的pragma指令，自动将程序并行处理，使用OpenMP降低了并行编程的难度和复杂度。当编译器不支持OpenMP时，程序会退化成普通（串行）程序。程序中已有的OpenMP指令不会影响程序的正常编译运行。在VS中启用OpenMP很简单，很多主流的编译环境都内置了OpenMP。在项目上右键->属性->配置属性->C/C+±>语言->OpenMP支持，选择“是”即可。

OpenMP编程模型

内存共享模型：OpenMP是专为多处理器/核，共享内存机器所设计的。底层架构可以是UMA和NUMA。即(Uniform Memory Access和Non-Uniform Memory Access)

基于线程的并行性

• OpenMP仅通过线程来完成并行
• 一个线程的运行是可由操作系统调用的最小处理单
• 线程们存在于单个进程的资源中，没有了这个进程，线程也不存在了
• 通常，线程数与机器的处理器/核数相匹配，然而，实际使用取决与应用程序

明确的并行

OpenMP是一种显式（非自动）编程模型，为程序员提供对并行化的完全控制。 一方面，并行化可像执行串行程序和插入编译指令那样简单。 另一方面，像插入子程序来设置多级并行、锁、甚至嵌套锁一样复杂

数据范围

由于OpenMP时是共享内存模型，默认情况下，在共享区域的大部分数据是被共享的, 并行区域中的所有线程可以同时访问这个共享的数据. 如果不需要默认的共享作用域，OpenMP为程序员提供一种“显示”指定数据作用域的方法

嵌套并行

API提供在其它并行区域放置并行区域，实际实现也可能不支持

动态线程

API为运行环境提供动态的改变用于执行并行区域的线程数， 实际实现也可能不支持

执行模式(Fork join)

OpenMP采用fork-join的执行模式。开始的时候只存在一个主线程，当需要进行并行计算的时候，派生出若干个分支线程来执行并行任务。当并行代码执行完成之后，分支线程会合，并把控制流程交给单独的主线程。
一个典型的fork-join执行模型的示意图如下：

OpenMP就是采用Fork-Join模型
所有的OpenML程序都以一个单个进程——master thread开始，master threads按顺序执行知道遇到第一个并行区域
Fork：主线程创造一个并行线程组
Join：当线程组完成并行区域的语句时，它们同步、终止，仅留下主线程

1.Types of variables

在并行节中，变量可以是私有的（每个线程拥有该变量的副本）或在所有线程之间共享。共享变量必须小心使用，因为它们会导致竞争条件

1. shared：共享并行区域内的数据，这意味着所有线程同时可见和可访问。默认情况下，工作共享区域中的所有变量都被共享，除了循环迭代计数器。
2. private:并行区域内的数据对每个线程都是专用的，这意味着每个线程将具有本地副本并将其用作临时变量。私有变量不会被初始化，并且该值不会保留供在并行区域之外使用。默认情况下，OpenMP循环结构中的循环迭代计数器是私有的。

Variables scope in OpenMP

1. shared(x) – all threads access the same memory location
2. private(x)
   (1) each thread has its own private copy of x
   (2) all local instances of x are not initialized
   (3) local updates to x are lost when exiting the parallel region
   (4) the original value of x is retained at the end of the block (OpenMP ≥ 3.0 only)
3. firstprivate(x)
   (1) each thread has its own private copy of x
   (2) all copies of x are initialized with the current value of x
   (3) local updates to x are lost when exiting the parallel region the original value of x is retained at the end of the block (OpenMP ≥ 3.0 only)
4. default(shared) or default(none)
   (1) affects all the variables not specified in other clauses
   (2) default(none)ensures that you must specify the scope of each variable used in the parallel block that the compiler can not figure out by itself (highly recommended!!)

2. Synchronization

OpenMP使您可以指定如何同步线程。选项：

1. critical: 包含的代码块一次只能由一个线程执行，而不能由多个线程同时执行。它通常用于保护共享数据免受竞争条件的影响
2. atomic:the memory update (write, or read-modify-write) in the next instruction will be performed atomically. It does not make the entire statement atomic; only the memory update is atomic. A compiler might use special hardware instructions for better performance than when using critical.
3. ordered:the structured block is executed in the order in which iterations would be executed in a sequential loop
4. barrier: 每个线程都等待，直到团队的所有其他线程都达到了这一点。A work-sharing construct has an implicit barrier synchronization at the end.
5. nowait: 指定可以在不等待团队中所有线程完成的情况下继续执行完成分配工作的线程。 In the absence of this clause, threads encounter a barrier synchronization at the end of the work sharing construct.
6. single

3. Task

1. OpenMP允许应用程序创建特定任务.A task is composed of:

• Code that will be executed
• Data environment (inputs/outputs)
• A location where the task will be executed (thread)
  When a thread encounters a task construct, 它可以选择立即执行任务或将其执行推迟到以后. 如果推迟，则将任务放置在与当前并行区域关联的概念库中

2. The task associated with a task construct will be executed only once. A task is tied to the code if the code is executed by the same thread from beginning to end. A task is untied the code can be executed by more than one thread.

4.编程要素

OpenMP编程模型以线程为基础，通过编译制导指令制导并行化，有三种编程要素可以实现并行化控制，他们分别是编译制导、API函数集和环境变量。

1. 编译器指令

OpenMP通过在串行程序中插入编译制导指令, 来实现并行化, 支持OpenMP的编译器可以识别, 处理这些指令并实现对应的功能. 所有的编译制导指令都是以#pragma omp开始, 后面跟具体的功能指令(directive)或者命令. 一般格式如下所示:

#pragma omp directive [clause [[,] clause]...]
    structured bloc

Parallel Construct(并行域结构)

为了使程序可以并行执行, 我们首先要构造一个并行域(parallel region), 在这里我们使用parallel指令来实现并行域的构造, 其语法形式如下

#pragma omp parallel  [clause [[,] clause]...]
     structured block

parallel ：用在一个结构块之前，表示这段代码将被多个线程并行执行；该指令只保证代码以并行的方式执行, 但是并不负责线程之间的任务分发. 在并行域执行结束之后, 会有一个隐式的屏障(barrier), 来同步所有的该区域内的所有线程. 下面是一个使用示例:

#include <iostream>
#include <omp.h>

using namespace std;

void fun1()

{
#pragma omp parallel num_threads(5)  //定义5个线程，每个线程都将运行{}内代码，运行结果：输出6次Test
    {
        cout << "Test" << endl;
    }

    system("pause");

}
int main()
{
    fun1();
}
/*TestTestTestTest
Test*/

parallel指令后面可以跟一些子句(clause), 如下所示

if(scalar-expression)
num_threads(integer-expression)
private(list)
firstprivate(list)
shared(list)
default(none | shared)
copyin(list)
reduction(operator:list)

后面介绍这些从句用法

Work-sharing Construct(任务分担结构)

任务分担指令主要用于为线程分配不同的任务, 一个任务分担域(work-sharing region)必须要和一个活跃(active)的并行域(parellel region)关联, 如果任务分担指令处于一个不活跃的并行域或者处于一个串行域中, 那么该指令就会被忽略. 在C/C++有3个任务分担指令: for、sections、single, 严格意义上讲只有for和sections是任务分担指令, 而single只是协助任务分担的指令.

Syntax in C++	功能
#pragma omp for	Distribute iterations over the threads
#pragma omp sections	Distribute independent work units
#pragma omp single	Only one thread executes the code block

for

用于for循环语句之前，表示将循环计算任务分配到多个线程中并行执行，以实现任务分担，必须由编程人员自己保证每次循环之间无数据相关性；
for指令一般可以和parallel指令合起来形成parallel for指令使用，也可以单独用在parallel语句的并行块中。parallel for用于生成一个并行域，并将计算任务在多个线程之间分配，用于分担任务。程序示例如下：

void parallel_for() {
    int n = 9;
    int i = 0;
    #pragma omp parallel shared(n) private(i) 
    {
        #pragma omp for
        for(i = 0; i < n; i++) {
            printf("Thread %d executes loop iteration %d\n", omp_get_thread_num(),i);
        }
    }
}
/*Thread 2 executes loop iteration 5
Thread 2 executes loop iteration 6
Thread 3 executes loop iteration 7
Thread 3 executes loop iteration 8
Thread 0 executes loop iteration 0
Thread 0 executes loop iteration 1
Thread 0 executes loop iteration 2
Thread 1 executes loop iteration 3
Thread 1 executes loop iteration 4
*/

在上面的程序中共有4个线程执行9次循环, 线程0分到了3次, 剩余的线程分到了2次, 这是一种常用的调度方式, 即假设有n次循环迭代, t个线程, 那么每个线程分配到n/t 或者 n/t + 1 次连续的迭代计算, 但是某些情况下使用这种方式并不是最好的选择, 我们可以使用schedule 来指定调度方式, 在后面会具体介绍. 下面是for 指令后面可以跟的一些子句:

private(list)
fistprivate(list)
lastprivate(list)
reduction(operator:list)
ordered
schedule(kind[,chunk_size])
nowait

下面还有一些示例

#include <iostream>
#include <omp.h>

using namespace std;

int main()
{
#pragma omp parallel  
    for (int i = 0; i < 10; i++)
    {
        cout << i;
    }

    return 0;
}
/*00102132435671203010102430182132456756479849342516278956788345956934576677889899*/

咋一看程序好像抽风了，产生了4个线程同时执行了for循环。通常这不是我们想要的，我们想要的是把for中的任务等分成4份，分别由4个线程各执行其中的一份。这样做其实很简单，只要在parallel后面加上for就可以了。

#include <iostream>
#include <omp.h>

using namespace std;


int main()
{
#pragma omp parallel for
    for (int i = 0; i < 10; i++)
    {
        cout << i;
    }

    return 0;
}

//result: 2890156734

结果为对的, 因为在同一段并行代码中，我们并不能保证各线程执行的先后顺序。扩散开来，在并行代码中，各线程的工作不能有依赖性，比如如果一个线程的输入和另一个线程的输出相依赖，那么此程序不适合并行计算。
下面再看下面的例子

#include <iostream>
#include <omp.h>

using namespace std;


int main()
{
    //for循环并行化声明形式1  
#pragma omp parallel   
    {
        cout << "OK" << endl;
#pragma omp for   
        for (int i = 0; i < 4; ++i)
        {
            cout << i << endl;
        }
    }

    cout << "形式2" << endl;

    //for循环并行化声明形式2  
#pragma omp parallel for  
    //cout << "ERROR" << endl; 
    for (int j = 0; j < 4; ++j)
    {
        cout << j << endl;
    }
    return 0;
}

这个结果为

形式1与形式2相比, 其优点是在for循环体前可以有其他执行代码，当然在一个#pragma omp parallel块内，可以有多个#pragma omp parallel for循环体。

for循环并行化的约束条件

尽管OpenMP可以方便地对for循环进行并行化，但并不是所有的for循环都可以进行并行化。以下几种情况不能进行并行化：

1. for循环中的循环变量必须是有符号整形。例如，for (unsigned int i = 0; i < 10; ++i){}会编译不通过；

2. for循环中比较操作符必须是<, <=, >, >=。例如for (int i = 0; i ！= 10; ++i){}会编译不通过；

3. for循环中的第三个表达式，必须是整数的加减，并且加减的值必须是一个循环不变量。例如for (int i = 0; i != 10; i = i + 1){}会编译不通过；感觉只能++i; i++; --i; 或i–；

4. 如果for循环中的比较操作为<或<=，那么循环变量只能增加；反之亦然。例如for (int i = 0; i != 10; --i)会编译不通过；

5. 循环必须是单入口、单出口，也就是说循环内部不允许能够达到循环以外的跳转语句，exit除外。异常的处理也必须在循环体内处理。例如：若循环体内的break或goto会跳转到循环体外，那么会编译不通过。

section

sections指令可以为不同的线程分配不同的任务, 语法如下所示:

#pragma omp sections [clause[[,] clause]...] 
    {
        [#pragma omp section]
            structured block
        [#pragma omp section]
            structured block
        ...
    }

从上面的代码中我们可以看到, sections将代码分为多个section, 每个线程处理一个section, 注意区分sections和section）下面是一个使用示例:

/**
 * 使用#pragma omp sections 和 #pragma omp section, 来使不同的线程执行不同的任务
 * 如果线程数量大于section数量, 那么多余的线程会处于空闲状态(idle)
 * 如果线程数量少于section数量, 那么一个线程会执行多个section代码
 */

void funcA() {
    printf("In funcA: this section is executed by thread %d\n",
            omp_get_thread_num());
}

void funcB() {
    printf("In funcB: this section is executed by thread %d\n",
            omp_get_thread_num());
}

void parallel_section() {
    #pragma omp parallel
    {
        #pragma omp sections
        {
            #pragma omp section 
            {
                (void)funcA();
            }

            #pragma omp section 
            {
                (void)funcB();
            }
        }
    } 
}
/*In funcA: this section is executed by thread 0
In funcB: this section is executed by thread 1
*/

下面是sections后面可以跟的一些子句

private(list)
firstprivate(list)
lastprivate(list)
reduction(operator:list)
nowait

single

single 指令用来指定某段代码块只能被一个线程来执行, 如果没有nowait字句, 所有线程在 single 指令结束处隐式同步点同步, 如果single指令有nowait从句, 则别的线程直接往下执行. 不过single指令并不指定哪个线程来执行. 语法如下所示:

#pragma omp single [clause[[,] clause]...]
    structured block

下面是一个示例

// Author: Mario Talevski
#include <iostream>
#include <cstdlib>
#include <fstream>
#include <string>
#include <sstream>
#include <omp.h>

void parallel_single() {
    int a = 0;
    int n = 10;
    int i;
    int b[10];
#pragma omp parallel shared(a, b) private(i)
    {
        // 只有一个线程会执行这段代码, 其他线程会等待该线程执行完毕
#pragma omp single 
        {
            a = 10;
            printf("Single construct executed by thread %d\n", omp_get_thread_num());
        }

        // A barrier is automatically inserted here

#pragma omp for
        for (i = 0; i < n; i++) {
            b[i] = a;
        }
    }

    printf("After the parallel region:\n");
    for (i = 0; i < n; i++)
        printf("b[%d] = %d\n", i, b[i]);
}


int main() {
    parallel_single();
}

/*Single construct executed by thread 0
After the parallel region:
b[0] = 10
b[1] = 10
b[2] = 10
b[3] = 10
b[4] = 10
b[5] = 10
b[6] = 10
b[7] = 10
b[8] = 10
b[9] = 10
*/

下面是single指令后面可以跟随的子句:

private(list)
firstprivate(list)
copyprivate(list)
nowait

Combined Parallel Work-Sharing Constructs

将parallel指令和work-sharing指令结合起来, 使代码更加简洁. 如下面的代码

#pragma omp parallel
{
    #pragma omp for
    for(.....)
}

也可以写为

#pragma omp parallel
{
    #pragma omp for
    for(.....)
}

具体参见下图

使用这些组合结构体(combined constructs)不仅增加程序的可读性, 而且对程序的性能有一定的帮助. 当使用这些组合结构体的时候, 编译器可以知道下一步要做什么, 从而可能会生成更高效的代码.

Clauses to Control Parallel and Work-Sharing Constructs

OpenMP指令后面可以跟一些子句, 用来控制构造器的行为. 下面介绍一些常用的子句.

shared

shared子句可以用于声明一个或多个变量为共享变量。所谓的共享变量，是值在一个并行区域的team内的所有线程只拥有变量的一个内存地址，所有线程访问同一地址。所以，对于并行区域内的共享变量，需要考虑数据竞争条件，要防止竞争，需要增加对应的保护, 下面是其使用方法:

#pragma omp parallel for shared(a)
    for(i = 0; i < n; i++)
    {
        a[i] += i;
    }

在并行域中使用共享变量时, 如果存在写操作, 需要对共享变量加以保存, 因为可能存在多个线程同时修改共享变量或者在一个线程读取共享变量时另外一个变量在更新共享变量的情况, 而这些情况都可能会引起程序错误. 程序示例如下：

#include <iostream>
#include <omp.h>
using namespace std;

#define COUNT  10000
int main()
{
    int sum = 0;
#pragma omp parallel for shared(sum)
    for (int i = 0; i < COUNT; i++)
    {
        sum = sum + i;
    }
    printf("%d\n", sum);
    return 0;
}
//result=46367978

多次运行，结果可能不一样. 需要注意的是：循环迭代变量 i 在循环构造区域里是私有的。即使使用shared来修饰循环迭代变量，也不会改变循环迭代变量在循环构造区域中是私有的这一特点。程序示例如下：

#include <iostream>
#include <omp.h>

using namespace std;

#define COUNT  10

int main()

{
    int sum = 0;
    int i = 0;
#pragma omp parallel for shared(sum, i)
    for (i = 0; i < COUNT; i++)
    {
        sum = sum + i;
    }
    printf("%d\n", i);
    printf("%d\n", sum);
    return 0;
}
/*0
 45*/

上述程序中，循环迭代变量i的输出值为0，尽管这里使用shared修饰变量i。注意，这里的规则只是针对循环并行区域，对于其他的并行区域没有这样的要求。同时在循环并行区域内，循环迭代变量是不可修改的。即在上述程序中，不能在for循环体内对循环迭代变量i进行修改。

#pragma omp parallel for shared(a)
    for(i = 0; i < n; i++)
    {
        a[i] += i;
    }

private

private子句用来指定哪些数据是线程私有的, 即每个线程具有变量的私有副本, 线程之间互不影响. 其语法形式为private(list), 使用方法如下:

// Author: Mario Talevski
#include <omp.h>
#include <iostream>


using namespace std;

void test_private() {
    int n = 8;
    int i = 2, a = 3;
    // i,a 定义为private之后不改变原先的值
#pragma omp parallel for private(i, a)
    for (i = 0; i < n; i++)
    {
        a = i + 1;
        printf("In for: thread %d has a value of a = %d for i = %d\n", omp_get_thread_num(), a, i);
    }

    printf("\n");
    printf("Out for: thread %d has a value of a = %d for i = %d\n", omp_get_thread_num(), a, i);
}

int main(int argc, char* argv[])
{
    test_private();
	return 0;
}

/*
In for: thread 3 has a value of a = 4 for i = 3
In for: thread 0 has a value of a = 1 for i = 0
In for: thread 5 has a value of a = 6 for i = 5
In for: thread 2 has a value of a = 3 for i = 2
In for: thread 4 has a value of a = 5 for i = 4
In for: thread 1 has a value of a = 2 for i = 1
In for: thread 7 has a value of a = 8 for i = 7
In for: thread 6 has a value of a = 7 for i = 6

Out for: thread 0 has a value of a = 3 for i = 2*/

对于private子句中的变量, 需要注意一下两点:

• 不论该变量之前有没有初始值, 在进入并行域之后都是未初始化的.
• 并行域中对变量的修改只在该域中起作用, 当离开并行域后, 变量值仍然是未进入并行域之前的值

lastprivate

lastprivate会在退出并行域时, 将其修饰变量的最后取值(last value)保存下来, 可以作用于 for 和 sections, 语法格式为lastprivate(list). 关于last value的定义: 如果是作用于for指令, 那么last value就是指串行执行的最后一次循环的值；如果是作用于sections指令, 那么last value就是执行完最后一个包含该变量的section之后的值. 使用方法如下:

#include <omp.h>
#include <iostream>


using namespace std;

void test_last_private() {
    int n = 8;
    int i = 2, a = 3;
    // lastprivate 将for中最后一次循环(i == n-1) a 的值赋给a    
#pragma omp parallel for private(i) lastprivate(a)
    for (i = 0; i < n; i++)
    {
        a = i + 1;
        printf("In for: thread %d has a value of a = %d for i = %d\n", omp_get_thread_num(), a, i);
    }

    printf("\n");
    printf("Out for: thread %d has a value of a = %d for i = %d\n", omp_get_thread_num(), a, i);
}

int main(int argc, char* argv[])
{
    test_last_private();
	return 0;
}
/*
In for: thread 2 has a value of a = 3 for i = 2
In for: thread 3 has a value of a = 4 for i = 3
In for: thread 5 has a value of a = 6 for i = 5
In for: thread 0 has a value of a = 1 for i = 0
In for: thread 6 has a value of a = 7 for i = 6
In for: thread 1 has a value of a = 2 for i = 1
In for: thread 4 has a value of a = 5 for i = 4
In for: thread 7 has a value of a = 8 for i = 7

Out for: thread 0 has a value of a = 8 for i = 2*/

firstprivate

firstprivate 子句用于为private变量提供初始值. 使用firstprivate修饰的变量会使用在前面定义的同名变量的值作为其初始值. 语法形式为firstprivate(list), 使用方法如下:

#include <omp.h>
#include <iostream>


using namespace std;


void test_first_private() {
    int n = 8;
    int i = 0, a[8];

    for (i = 0; i < n; i++) {
        a[i] = i + 1;
    }
#pragma omp parallel for private(i) firstprivate(a)
    for (i = 0; i < n; i++)
    {
        printf("thread %d: a[%d] is %d\n", omp_get_thread_num(), i, a[i]);
    }
}

int main(int argc, char* argv[])
{
    test_first_private();
	return 0;
}

/*
thread 3: a[3] is 4
thread 4: a[4] is 5
thread 0: a[0] is 1
thread 5: a[5] is 6
thread 2: a[2] is 3
thread 1: a[1] is 2
thread 6: a[6] is 7
thread 7: a[7] is 8
*/

下面是合并的例子

```cpp
#include <iostream>
#include <omp.h>

using namespace std;

int main()
{
    int k = 100;

#pragma omp parallel for firstprivate(k),lastprivate(k)

    for (int i = 0; i < 3; i++)
    {
        k += i;
        cout<<"kin"<<k<<endl;
    }
    cout << "kout " << k << endl;
}
/*
kinkin102100
kin101

kout 102*/

打印结果看出, 退出for循环的并行区域后，共享变量k的值变成了103，而不是保持原来的100不变。OpenMP规范中指出，如果是循环迭代，那么是将最后一次循环迭代中的值赋给对应的共享变量；如果是section构造，那么是最后一个section语句中的值赋给对应的共享变量。注意这里说的最后一个section是指程序语法上的最后一个，而不是实际运行时的最后一个运行完的。

default

default指定并行区域内变量的属性，C++的OpenMP中default的参数只能为shared或none。default(shared)：表示并行区域内的共享变量在不指定的情况下都是shared属性

default(none)：表示必须显式指定所有共享变量的数据属性，否则会报错，除非变量有明确的属性定义（比如循环并行区域的循环迭代变量只能是私有的）如果一个并行区域，没有使用default子句，那么其默认行为为default(shared)

nowait

用于取消任务分担结构(work-sharing constructs)中的隐式屏障(implicit barrier)(即并行区域中所有线程执行完毕之后，主线程才继续执行), 下面是一个使用示例:

#include <iostream>
#include <omp.h>
using namespace std;

void test_nowait() {
    int i, n = 3;
    //#pragma omp parallel 
    {
        //#pragma omp for nowait 
        for (i = 0; i < n; i++) {
            printf("thread %d: ++++\n", omp_get_thread_num());
        }

        #pragma omp for
        for (i = 0; i < n; i++) {
            printf("thread %d: ----\n", omp_get_thread_num());
        }
    }
}
int main()
{
    test_nowait();
}
/*
thread 3: ++++
thread 0: ++++
thread 0: ++++
thread 2: ++++
thread 1: ++++
thread 1: ++++
thread 0: ----
thread 0: ----
thread 3: ----
thread 1: ----
thread 1: ----
thread 2: ----

*/

因为for指令有一个隐式的屏障, 会同步所有的线程直到第一个for循环执行完, 再继续往下执行. 加上 nowait 之后就消除了这个屏障, 使线程执行完第一个for循环之后无需再等待其他线程就可以去执行第二个for循环的内容, 下面是加上nowait之后的输出:

#include <iostream>
#include <omp.h>


using namespace std;


void test_nowait() {
    int i, n = 3;
    #pragma omp parallel 
    {
        #pragma omp for nowait 
        for (i = 0; i < n; i++) {
            printf("thread %d: ++++\n", omp_get_thread_num());
        }

        #pragma omp for
        for (i = 0; i < n; i++) {
            printf("thread %d: ----\n", omp_get_thread_num());
        }
    }
}
int main()
{
    test_nowait();
}
/*
thread 0: ++++
thread 0: ----
thread 1: ++++
thread 1: ----
thread 2: ++++
thread 2: ----
*/

使用nowait时需要注意前后for之间有没有依赖关系, 如果第二个for循环需要用到第一个for循环的结果, 那么使用nowait就可能会造成程序错误.

schedule

schedule子句只作用于循环结构(loop construct), 它用来设置循环任务的调度方式. 语法形式为schedule(kind[,chunk_size]), 其中kind的取值有static, dynamic, guided,auto ,runtime chunk_size是可选项,可以指定也可以不指定. 下面是使用方法:

#include <iostream>
#include <omp.h>
using namespace std;

void test_schedule() {
    int i, n = 10;

#pragma omp parallel for default(none) schedule(static, 2) private(i) shared(n)
    for (i = 0; i < n; i++) {
        printf("Iteration %d executed by thread %d\n", i, omp_get_thread_num());
    }
}
int main()
{
    test_schedule();
}
/*
Iteration 8 executed by thread 4
Iteration 9 executed by thread 4
Iteration 0 executed by thread 0
Iteration 1 executed by thread 0
Iteration 2 executed by thread 1
Iteration 3 executed by thread 1
Iteration 4 executed by thread 2
Iteration 5 executed by thread 2
Iteration 6 executed by thread 3
Iteration 7 executed by thread 3
*/

下面介绍一下各个取值的含义, 假设有n次循环, t个线程

static

静态调度, 如果不指定chunk_size , 那么会为每个线程分配 n/t 或者 n/t+1(不能除尽)次连续的迭代计算, 如果指定了 chunk_size, 那么每次为线程分配chunk_size次迭代计算, 如果第一轮没有分配完, 则循环进行下一轮分配, 假设n=8, t=4, 下表给出了chunk_size未指定、等于1、等于3时的分配情况.

dynamic

动态调度, 动态为线程分配迭代计算, 只要线程空闲就为其分配任务, 计算快的线程分配到更多的迭代. 如果不指定chunk_size参数, 则每次为一个线程分配一次迭代循环(相当于chunk_size=1), 若指定chunk_size, 则每次为一个线程分配chunk_size次迭代循环. 在动态调度下, 分配结果是不固定的, 重复执行同一个程序, 每次的分配结果一般来说是不同的, 下面给出n=12, t=4时, chunk_size未指定、等于2时的分配情况(运行两次)

使用动态动态可以一定程度减少负载不均衡的问题, 但是需要注意任务动态申请时也会有一定的开销.

guided

guided调度是一种指定性的启发式自调度方法. 开始时每个线程会分配到较大的迭代块, 之后分配到的迭代块的大小会逐渐递减. 如果指定chunk_size, 则迭代块会按指数级下降到指定的chunk_size大小, 如果没有指定size参数, 那么迭代块大小最小会降到1(相当于chunk_size=1). 和动态调度一样, 执行块的线程会分到更多的任务, 不同的是这里迭代块的大小是变化的. 同样使用guided调度的分配结果也不是固定的, 重复执行会得到不同的分配结果. 下面给出n=20, t=4, chunk_size未指定、chunk_size=3时的分配情况(执行两次)

当设置chunk_size=3时, 因为最后只剩下18、19两次循环, 所以最后执行的那个线程只分配到2次循环.

下面的图展示了当循环次数为200次, 线程数量为4时, static 、 (dynamic,7) 、(guided, 7) 3种调度方式的分配情况

负载不均衡

在for循环中, 如果每次循环之间花费的时间是不同的, 那么就可能出现负载不均衡问题,比如,

static   : use time 8.67s
static,2 : use time 6.42s
dynamic  : use time 5.62s
dynamic,2: use time 6.43s
guided   : use time 5.92s
guided,2 : use time 6.43s

对于static调度, 如果不指定chunk_size的值, 则会将最后几次循环分给最后一个线程, 而最后几次循环是最耗时的, 其他线程执行完各自的工作需要等待这个线程执行完毕, 浪费了系统资源, 这就造成了负载不均衡. dynamic和guided可以在一定程度上减轻负载不均衡, 但是也不是绝对的, 最终选用哪种方式还是要根据具体的问题.

Synchronization Constructs(同步)

同步指令主要用来控制多个线程之间对于共享变量的访问. 它可以保证线程以一定的顺序更新共享变量, 或者保证两个或多个线程不同时修改共享变量.

barrier

同步路障(barrier), 当线程遇到路障时必须要停下等待, 直到并行区域中的所有线程都到达路障点, 线程才继续往下执行. 在每一个并行域和任务分担域的结束处都会有一个隐式的同步路障, 即在parallel、for、sections、single构造的区域之后会有一个隐式的路障, 因此在很多时候我们无需显示的插入路障. 下面示例:

#include <iostream>
#include <omp.h>
#include <chrono>
#include <ctime>
#include <vector>
using namespace std;



void test_barrier() {
#pragma omp parallel num_threads(6)
{
    #pragma omp critical
    cout << omp_get_thread_num() << " ";
     #pragma omp critical
    cout << omp_get_thread_num() + 10 << " ";
    }
}

int main() {
    test_barrier();
}
/*
0 10 3 13 4 14 2 12 5 15 1 11*/

我们加上barrier

#include <iostream>
#include <omp.h>
#include <chrono>
#include <ctime>
#include <vector>
using namespace std;



void test_barrier() {
#pragma omp parallel num_threads(6)
{
    #pragma omp critical
    cout << omp_get_thread_num() << " ";
    #pragma omp barrier
     #pragma omp critical
    cout << omp_get_thread_num() + 10 << " ";
    }
}

int main() {
    test_barrier();
}
/*
4 1 2 0 3 5 15 11 13 12 14 10*/

可以看到，这时一位数数字打印完了才开始打印两位数数字，因为，所有线程执行到第5行代码时，都要等待所有线程都执行到第5行，这时所有线程再都继续执行第7行及以后的代码，即所谓同步.
再来说说for, sections, single directives的隐含barrier，以及nowait clause如下示例：

#include <iostream>
#include <omp.h>
#include <chrono>
#include <ctime>
#include <vector>
using namespace std;



void test_barrier() {
#pragma omp parallel num_threads(6)
{   
    #pragma omp for
        for (int i = 0; i < 10; ++i) {
        #pragma omp critical
            cout << omp_get_thread_num() << " ";
    }
        //this is an implicit barrier here
     #pragma omp critical
    cout << omp_get_thread_num() + 10 << " ";
    }
}

int main() {
    test_barrier();
}
/*
0 0 3 3 1 1 5 4 2 2 14 12 13 15 11 10*/

加上nowait将会disable implicit barrier

#include <iostream>
#include <omp.h>
#include <chrono>
#include <ctime>
#include <vector>
using namespace std;



void test_barrier() {
#pragma omp parallel num_threads(6)
{   
    #pragma omp for nowait
        for (int i = 0; i < 10; ++i) {
        #pragma omp critical
            cout << omp_get_thread_num() << " ";
    }
        //The implicit barrier here is disabled by nowait.
     #pragma omp critical
    cout << omp_get_thread_num() + 10 << " ";
    }
}

int main() {
    test_barrier();
}

/*
0 0 10 2 2 12 1 1 11 5 15 3 3 13 4 14*/

sections, single directives是类似的.

ordered

ordered结构允许在并行域中以串行的顺序执行一段代码, 如果我们在并行域中想按照顺序打印被不同的线程计算的数据, 就可以使用这个子句. 在使用时注意以下两点

• ordered 只作用于循环结构(loop construct)
• 使用ordered时需要在构造并行域的时候加上ordered子句, 如#pragma omp parallel for ordered

#include <iostream>
#include <omp.h>
#include <chrono>
#include <ctime>
#include <vector>
using namespace std;

void test_order() {
    int i, tid, n = 5;
    int a[5];
    for (i = 0; i < n; i++) {
        a[i] = 0;
    }

    #pragma omp parallel for default(none) ordered  schedule(dynamic) private (i, tid) shared(n, a)
    for (i = 0; i < n; i++) {
        tid = omp_get_thread_num();
        printf("Thread %d updates a[%d]\n", tid, i);

        a[i] += i;

         #pragma omp ordered
        {
            printf("Thread %d printf value of a[%d] = %d\n", tid, i, a[i]);
        }
    }
}

int main() {
    test_order();
}

/*
Thread 0 updates a[0]
Thread 0 printf value of a[0] = 0
Thread 4 updates a[2]
Thread 3 updates a[4]
Thread 1 updates a[3]
Thread 2 updates a[1]
Thread 2 printf value of a[1] = 1
Thread 4 printf value of a[2] = 2
Thread 1 printf value of a[3] = 3
Thread 3 printf value of a[4] = 4
*/

从输出结果我们可以看到, 在update时是以乱序的顺序更新, 但是在打印时是以串行顺序的形式打印.

critical

临界区(critical), 临界区保证在任意一个时间段内只有一个线程执行该区域中的代码, 一个线程要进入临界区必须要等待临界区处于空闲状态, 下面是语法形式

#pragma omp critical [(name)]
    structured block

其中name是为临界区指定的一个名字. 下面是一个求和的使用示例, 注意这里只是用来说明临界区的作用, 对于求和操作我们可以使用reduction指令

#include <iostream>
#include <omp.h>
#include <chrono>
#include <ctime>
#include <vector>
using namespace std;



void test_critical() {
    int n = 100, sum = 0, sumLocal, i, tid;
    int a[100];
    for (i = 0; i < n; i++) {
        a[i] = i;
    }

#pragma omp parallel shared(n, a, sum) private (tid, sumLocal)
    {
        tid = omp_get_thread_num();
        sumLocal = 0;
        #pragma omp for
        for (i = 0; i < n; i++) {
            sumLocal += a[i];
        }

        #pragma omp critical(update_sum) 
        {
            sum += sumLocal;
            printf("Thread %d: sumLocal = %d sum =%d\n", tid, sumLocal, sum);
        }
    }

    printf("Value of sum after parallel region: %d\n", sum);
}

int main() {
    test_critical();
}
/*
Thread 6: sumLocal = 978 sum =978
Thread 0: sumLocal = 78 sum =1056
Thread 2: sumLocal = 416 sum =1472
Thread 1: sumLocal = 247 sum =1719
Thread 4: sumLocal = 690 sum =2409
Thread 5: sumLocal = 834 sum =3243
Thread 3: sumLocal = 585 sum =3828
Thread 7: sumLocal = 1122 sum =4950
Value of sum after parallel region: 4950
下面是将临界区去掉的运行结果(运行结果不是固定的, 这里只是其中一种情况):
Thread 5: sumLocal = 834 sum =834
Thread 2: sumLocal = 416 sum =3531
Thread 6: sumLocal = 978 sum =3115
Thread 1: sumLocal = 247 sum =247
Thread 7: sumLocal = 1122 sum =1447
Thread 3: sumLocal = 585 sum =4116
Thread 0: sumLocal = 78 sum =325
Thread 4: sumLocal = 690 sum =2137
Value of sum after parallel region: 4116
*/

通过对比我们可以看到临界区保证了程序的正确性.

atomic

原子操作, 可以锁定一个特殊的存储单元(可以是一个单独的变量,也可以是数组元素), 使得该存储单元只能原子的更新, 而不允许让多个线程同时去写. atomic只能作用于单条赋值语句, 而不能作用于代码块. 语法形式为:

#pragma omp atomic
    statement

在C/C++中, statement必须是下列形式之一

• x++, x--, ++x, --x,
• x binop= expr其中binop是二元操作符:+, -, *, /, &, ^, |, <<, >>之一
  atomic的可以有效的利用的硬件的原子操作机制来控制多个线程对共享变量的写操作, 效率较高, 下面是一个使用示例

#include <iostream>
#include <omp.h>
#include <chrono>
#include <ctime>
#include <vector>
using namespace std;


void test_atomic() {
    int counter = 0, n = 1000000, i;

#pragma omp parallel for shared(counter, n)
    for (i = 0; i < n; i++) {
        #pragma omp atomic
        counter += 1;
    }

    printf("counter is %d\n", counter);
}

int main() {
    test_atomic();
}
/*
counter is 1000000
*/

对于下面的情况

#pragma omp atomic
ic += func();

atomic只保证ic的更新是原子的, 即不会被多个线程同时更新, 但是不会保证func函数的执行是原子的, 即多个线程可以同时执行func函数, 如果要使func的执行也是原子的, 可以使用临界区.

master

用于指定一段代码只由主线程执行. master指令和single指令的区别如下:

• master指令包含的代码段只有主线程执行, 而single指令包含的代码可以由任意一个线程执行.
• master指令在结束处没有隐式同步, 也不可以使用nowait从句
  下面使用一个示例

#include <iostream>
#include <omp.h>
#include <chrono>
#include <ctime>
#include <vector>
using namespace std;


void test_master() {
    int a, i, n = 5;
    int b[5];
#pragma omp parallel shared(a, b) private(i)
    {
         #pragma omp master
        {
            a = 10;
            printf("Master construct is executed by thread %d\n", omp_get_thread_num());
        }
         #pragma omp barrier

          #pragma omp for
        for (i = 0; i < n; i++)
            b[i] = a;
    }

    printf("After the parallel region:\n");
    for (i = 0; i < n; i++)
        printf("b[%d] = %d\n", i, b[i]);
}

int main() {
    test_master();
}
/*
Master construct is executed by thread 0
After the parallel region:
b[0] = 10
b[1] = 10
b[2] = 10
b[3] = 10
b[4] = 10
*/

flush：保证各个OpenMP线程的数据影像的一致性；

threadprivate：用于指定一个或多个变量是线程专用，后面会解释线程专有和私有的区别。

相应的OpenMP子句为：

reduction：用来指定一个或多个变量是私有的，并且在并行处理结束后这些变量要执行指定的归约运算，并将结果返回给主线程同名变量；
num_threads：指定并行域内的线程的数目；

copyprivate：配合single指令，将指定线程的专有变量广播到并行域内其他线程的同名变量中；
copyin n：用来指定一个threadprivate类型的变量需要用主线程同名变量进行初始化；

2. API函数

除上述编译制导指令之外，OpenMP还提供了一组API函数用于控制并发线程的某些行为，下面是一些常用的OpenMP API函数以及说明：

// 设置并行线程数
_OMPIMP void _OMPAPI omp_set_num_threads(int _Num_threads);
 
// 获取当前并行线程数
_OMPIMP int  _OMPAPI omp_get_num_threads(void);
 
// 获取当前系统最大可并行运行的线程数
_OMPIMP int  _OMPAPI omp_get_max_threads(void);
 
// 获取当前运行线程的ID，注意和操作系统中的线程ID不同
_OMPIMP int  _OMPAPI omp_get_thread_num(void);
 
// 获取当前系统中处理器数目
_OMPIMP int  _OMPAPI omp_get_num_procs(void);
 
_OMPIMP void _OMPAPI omp_set_dynamic(int _Dynamic_threads);
 
_OMPIMP int  _OMPAPI omp_get_dynamic(void);
 
_OMPIMP int  _OMPAPI omp_in_parallel(void);
 
_OMPIMP void _OMPAPI omp_set_nested(int _Nested);
 
_OMPIMP int  _OMPAPI omp_get_nested(void);
 
_OMPIMP void _OMPAPI omp_init_lock(omp_lock_t * _Lock);
 
_OMPIMP void _OMPAPI omp_destroy_lock(omp_lock_t * _Lock);
 
_OMPIMP void _OMPAPI omp_set_lock(omp_lock_t * _Lock);
 
_OMPIMP void _OMPAPI omp_unset_lock(omp_lock_t * _Lock);
 
_OMPIMP int  _OMPAPI omp_test_lock(omp_lock_t * _Lock);
 
_OMPIMP void _OMPAPI omp_init_nest_lock(omp_nest_lock_t * _Lock);
 
_OMPIMP void _OMPAPI omp_destroy_nest_lock(omp_nest_lock_t * _Lock);
 
_OMPIMP void _OMPAPI omp_set_nest_lock(omp_nest_lock_t * _Lock);
 
_OMPIMP void _OMPAPI omp_unset_nest_lock(omp_nest_lock_t * _Lock);
 
_OMPIMP int  _OMPAPI omp_test_nest_lock(omp_nest_lock_t * _Lock);
 
_OMPIMP double _OMPAPI omp_get_wtime(void);
 
_OMPIMP double _OMPAPI omp_get_wtick(void);

使用上述函数的一些例子

#include <iostream>
#include <omp.h>

using namespace std;


void TestAPIs()
{   // 获取当前并行线程数
    cout << "Num of Procs: " << omp_get_num_procs() << endl;
    // 获取当前系统最大可并行运行的线程数
    cout << "Max Threads: " << omp_get_max_threads() << endl;
    cout << "Set Num of Threads = 2 " << endl;
    // 设置并行线程数
    omp_set_num_threads(2);

#pragma omp parallel
    cout << "Get Thread Num: " << omp_get_thread_num() << endl;

  // 设置并行线程数
    omp_set_num_threads(omp_get_num_procs() - 1);
#pragma omp parallel 
    {
        cout << "OPENMP\n";
    }
}

int main()
{
    TestAPIs();
}

3.环境变量

OpenMP提供了一些环境变量，用来在运行时对并行代码的执行进行控制。这些环境变量可以控制：1）设置线程数；2）指定循环如何划分；3）将线程绑定到处理器；4）启用/禁用嵌套并行，设置最大的嵌套并行级别；5）启用/禁用动态线程；6）设置线程堆栈大小；7）设置线程等待策略。常用的环境变量：

OMP_SCHEDULE ：用于for循环并行化后的调度，它的值就是循环调度的类型；
OMP_NUM_THREADS：用于设置并行域中的线程数；
OMP_DYNAMIC：通过设定变量值，来确定是否允许动态设定并行域内的线程数；
OMP_NESTED：指出是否可以并行嵌套。

OpenMP指令及其用法

4.threadprivate

7. Copyin

copyin子句用于将主线程中threadprivate变量的值拷贝到执行并行区域的各个线程的threadprivate变量中，从而使得team内的子线程都拥有和主线程同样的初始值。程序示例如下：

8. Copyprivate

OpenMP中的任务调度

OpenMP中，任务调度主要用于并行的for循环中，当循环中每次迭代的计算量不相等时，如果简单地给各个线程分配相同次数的迭代的话，会造成各个线程计算负载不均衡，这会使得有些线程先执行完，有些后执行完，造成某些CPU核空闲，影响程序性能。OpenMP提供了schedule子句来实现任务的调度。
schedule(type,[size])
参数type是指调度的类型，可以取值为static，dynamic，guided，runtime四种值。其中runtime允许在运行时确定调度类型，因此实际调度策略只有前面三种。

参数size表示每次调度的迭代数量，必须是整数。该参数是可选的。当type的值是runtime时，不能够使用该参数。

1. 静态调度static

大部分编译器在没有使用schedule子句的时候，默认是static调度。static在编译的时候就已经确定了，那些循环由哪些线程执行。假设有n次循环迭代，t个线程，那么给每个线程静态分配大约n/t次迭代计算。n/t不一定是整数，因此实际分配的迭代次数可能存在差1的情况。

在不使用size参数时，分配给每个线程的是n/t次连续的迭代，若循环次数为10，线程数为2，则线程0得到了0～4次连续迭代，线程1得到5～9次连续迭代。

当使用size时，将每次给线程分配size次迭代。若循环次数为10，线程数为2，指定size为2则0、1次迭代分配给线程0，2、3次迭代分配给线程1，以此类推。

2.动态调度

动态调度依赖于运行时的状态动态确定线程所执行的迭代，也就是线程执行完已经分配的任务后，会去领取还有的任务（与静态调度最大的不同，每个线程完成的任务数量可能不一样）。由于线程启动和执行完的时间不确定，所以迭代被分配到哪个线程是无法事先知道的。

当不使用size 时，是将迭代逐个地分配到各个线程。当使用size 时，逐个分配size个迭代给各个线程，这个用法类似静态调度。

3. 启发式调度

采用启发式调度方法进行调度，每次分配给线程迭代次数不同，开始比较大，以后逐渐减小。开始时每个线程会分配到较大的迭代块，之后分配到的迭代块会逐渐递减。迭代块的大小会按指数级下降到指定的size大小，如果没有指定size参数，那么迭代块大小最小会降到1。

size表示每次分配的迭代次数的最小值，由于每次分配的迭代次数会逐渐减少，少到size时，将不再减少。具体采用哪一种启发式算法，需要参考具体的编译器和相关手册的信息。

OpenMP程序设计技巧总结

1.当循环次数较少时，如果分成过多的线程来执行的话，可能会使得总的运行时间高于较少线程或一个线程的执行情况，并且会增加能耗；

2.如果设置的线程数量远大于CPU的核数的话，那么存在着大量的任务切换和调度的开销，也会降低整体的效率。

3.在嵌套循环中，如果外层循环迭代次数较少时，如果将来CPU核数增加到一定程度时，创建的线程数将可能小于CPU核数。另外如果内层循环存在负载平衡的情况下，很难调度外层循环使之达到负载平衡。