Fortran程序CUDA并行化总结

段红英

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

关键词：Fortran ；C；CUDA；Fortran程序

中图分类号：TP311 文献标识码：A 文章编号：2095-1302（2015）11-00-02

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）。

（a） C语言数组分割图 （b）Fortran数组分割图

图1 不同数组的分割图

（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为粒子总数。由于每个粒子输运过程相对独立、粒子间通信量少、循环次数多，因此，可以一个粒子对应一个线程来并行。

图2 蒙卡输运框架

当然，若是有对程序足够的理解，并且Fortran和C都很精通的情况下，则可以直接将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移植过程中所遇到的一些问题进行总结。当然为了更高的计算效率，可以对初步的程序优化。

参考文献

[1]龚春叶.面向异构体系结构的粒子输运并行算法研究[D].长沙：国防科学技术大学，2011.

[2] NVIDIA.CUDA C PROGRAMMING GUIDE[Z].2012.

[3] Stephen J.Chapman. Fortran 95/2003程序设计[M].刘瑾，庞岩梅，赵越，等，译.北京：中国电力出版社，2009.

[4] Jason Sanders，EdwardKandrot.GPU高性能编程CUDA实战[M]. 聂雪军，等，译.北京：机械工业出版社，2011：43-68.

[5] Randal L.Scbwartz，TomPboenix brain d foy. Perl语言入门[M].盛春，蒋永清，王晖，译.南京：东南大学出版社，2007.