Fortran程序CUDA并行化总结

时间：2021-08-22 15:43:32

关键字： Fortran C CUDA Fortran程序

[导读]摘要：虽然Fortran常用来进行科学计算，但是面对计算量大的程序仍然很耗时。通常人们用MPI进行粗粒度的并行来提高程序的运行效率，近年来随着GPU计算能力的提高，将程序进行细粒度GPU并行化成为一种趋势。文章基于NVIDIA公司的CUDA框架，就Fortran程序向CUDA移植过程中的一些问题进行总结，并给出了相应的解决方案。

0 引言

Fortran 是常用的科学计算语言，其突出的特性就是能实现自然描述且描述接近数学公式，有较好的执行效率，但是由于在计算流体力学、现代医学影像、分子动力学等领域的模拟中，存在大量的程序计算，仍然耗时很多，有的计算需要几天甚至几十天才能完成。为了提高程序的计算效率，我们将 Fortran 代码并行化。通常人们用 MPI 进行粗粒度的并行来提高程序的运行效率，近年来随着 GPU 计算能力的提高，将程序进行细粒度 GPU 并行化成为一种趋势。

CUDA 是 NVIDIA 公司推出的一种用于 GPU 高性能计算的软硬件架构，它是对 C 语言的扩展。在其编程模型中，CPU 作为主机（Host）端调度整个程序，GPU 作为计算设备（device）对计算量大、数据并行性强的程序并行处理。运行在 GPU上的并行计算程序称为 kernel，其必须通过 __global__ 函数类型限定符定义，由 host 端程序调用启动。

Fortran 程序 CUDA 并行化的完成一般分为编码、编译、测试、优化几个阶段，以下从这几个方面分别进行总结阐述。

1 编码

在编码阶段，Fortran 程序 CUDA 并行化即为 Fortran → C → CUDA 的过程。

1.1 Fortran → C 的转化过程

从 Fortran 到 C 的转化过程，只需要在掌握二者语法的基础上，逐行翻译即可，但翻译工作中需要注意以下细节。

（1）数组

C 语言中数组的起始编号为 0，而 Fortran 的默认起始编号为 1，但也可以用（idx1 ：idx2）的方式自己定义，这就需要我们把程序中的每个数组变量的定义弄明白，翻译时对默认定义的数组标号减 1，非默认定义的，则用 [i-idx1] 来计算实际标号。

其次是多维阵列。虽然 C 和 Fortran 中所谓的多维阵列都是一个连续的一维存储空间，但是它们对于行列的分割却相反。如图 1（a）和图 1（b）分别给出了 C 和 Fortran 对于数组 a[3][2] 各自的数组分割方式。由此，我们在翻译过程中定义和使用多维数组时都须将数组的行列转换。如 a（3 ：2）变为 C 时应为 a[2][3]，对应图 1（b）。
Fortran程序CUDA并行化总结

（2）函数参数传递 Fortran 中函数调用时一般传递的是参数的内存地址，而 C 既可以进行值传递，也可以进行地址传递，一般需要返回多个参数值时用地址传递。在翻译中，为了方便，所有函数都采用地址传递。

（3）函数重载在 Fortran 中为了共享数据的方便一般会用 common，如下例所示，Fortran 代码为：

Integer ：：I，j，k

Integer ：：kk（3）

COMMON/test/I，j，k

用 C 表达，需要用以下的方式：

int I，j，k ；

int kk[3] ；

int*test[3]={&I，&J，&K} ；

即把变量的地址连续地存储到一个数组中。在函数参数传递时，在 Fortran 中调用 addkt 函数就可以传递数组 kk，也可以传递 test，代码如下：

subroutineaddkt（ kd，kt ）

Integer，intent（inout）：：kd（3）

但是在 C 语言中要传递以上两种参数就出现函数重载问题（一个为一维数组，另一个为二维数组），但对于过程化语言 C 则没有该功能，我们只能把函数 addt 定义为 addkt1（int *kd，intkt）和 addkt2（int *kd[]，intkt）两个函数。

1.2 C → CUDA 的转化过程

CUDA 是一种数据并行性而非功能并行性的并行计算解决方案。在 C 到 CUDA 的转化过程中，最关键的就是分析整个程序，找到最耗时的代码部分，分析整个的可并行性，在对整个物理过程理解的基础上，进行算法设计，然后并行化。

以核物理中的蒙卡输运程序为例，蒙特卡罗（MC）方法采用随机方法模拟物理过程，应用数理统计获取计算结果的计算方法。蒙卡的整个输运框架如图 2 所示，其中，n 为粒子编号，N 为粒子总数。由于每个粒子输运过程相对独立、粒子间通信量少、循环次数多，因此，可以一个粒子对应一个线程来并行。

当然，若是有对程序足够的理解，并且 Fortran 和Ｃ都很精通的情况下，则可以直接将 Fortran 程序 CUDA 并行化。

2 测试

我们借助 GDB 调试工具，将测试过程分为由下到上，和由上到下两步。首先，由下到上的对单个程序逐级测试；然后，根据程序写出多种输入参数，由上往下整体测试。

判定程序正确的方法就是此程序有相同的输入和输出，有随机数的程序会给我们的测试带来很大的困难。如上面提到的蒙卡输运程序，我们既要保证程序中输运过程的随机性，又要通过测试保证程序的正确性。一般大家会想到产生一个很大的随机数文件，分别读入到 Fortran 和 C 程序中，此方法可行，但是在粒子数很多的情况下，效率很低。文章就此问题提出了很好的解决方案，此处不再详细解述。

2.1 Fortan → C 测试

对于变量少的程序，我们只需要手动打印出需要检测的变量值进行对比，但是对于有几千个全局变量的计算程序，逐一手动输出效率会很低。因此我们首先要找到程序中用到的全局变量，然后根据这些变量书写 main 函数测试。这个问题我们可以借助 Perl、Shell 等脚本语言完成。其步骤如下：

（1）人工找到程序中所有的全局变量，其类型，维数，每一维的长度按照某种格式写到文件 A中；

（2）人工找出程序中的所有程序、子程序及函数的名字，写到 B 中；

（3）书写脚本，根据文件 A，在所要测试的 Fortran 程序的初始化部分，打印出所有全局变量的值，作为 C 程序的输入；在所要测试的 Fortran 程序结束之前，打印出所有全局变量的值，作为判断 C 程序的标准；

（4）书写脚本，实现初始化函数，即给 C 语言的变量初始化；

（5）书写脚本，实现读入 Fortran 的输出值，判断 C 程序的正确性。

在具体的脚本实现过程中，需要对 Fortran 语法详细分解，如一个程序可能有多处结束，而在每个结束前都需要输出打印。为了方便，同时打印出变量的类型、维数及每一维的长度。

2.2 Fortan → CUDA 测试

CUDA 程序的测试并不像 C 程序的测试那样简单，因为由于算法原因，在功能不变的情况下整个程序可能会被打乱，甚至对变量数组进行排序，由于原来的数组会打乱，致使不能用 GDB 调试工具进行对比，并且无法定位错误或者判断计算中的正确性。虽然程序及数组顺序乱了，但是程序的功能不会变，并且在相同功能的地方这些数组之和不会变，所以我们可以在功能相同的地方对数据之和进行对比。

3 结语

面对科学计算中串行程序的瓶颈，我们需要应用并行化方法来解决，本文就 Fortran 向基于 CUDA 架构的 GPU 移植过程中所遇到的一些问题进行总结。当然为了更高的计算效率，可以对初步的程序优化。