国产加速器海光DCU&GPGPU深算处理器异构编程实战（下）

技术瘾君子1573

已于 2024-08-26 07:00:55 修改

阅读量2k

点赞数 67

分类专栏： Linux并行计算&HPC高性能计算文章标签：开发语言 OpenMP Fortran DCU 异构编程

于 2024-08-26 06:30:00 首次发布

本文链接：https://blog.csdn.net/qq_27815483/article/details/141326281

版权

四、使用DCU加速Fortran程序

4.1 概述

4.2 在Fortran程序中使用HIP编程接口

第一二章内容请看国产加速器海光DCU&GPGPU深算处理器异构编程实战（上）-CSDN博客

第三章内容请看国产加速器海光DCU&GPGPU深算处理器异构编程实战（中）-CSDN博客

四、使用DCU加速Fortran程序

Fortran语言自20世纪50年代诞生至今，一直被广泛应用于计算物理学、计算化学、流体力学、气候预测、有限元分析等计算密集的学科领域。同C/C++语言一样，Fortran语言也是科学/工程计算领域的主流编程语言，使用Fortran语言编写的应用程序运行在世界各地大大小小的计算集群上，为人类的各种研究与创新提供数据支撑。

本章将介绍如何让使用Fortran语言编写的程序也能利用DCU进行并行加速。

4.1 概述

在上一章，我们介绍了如何使用C/C++语言编写DCU程序，在3.1.7节中我们总结了编写DCU程序主要有两方面工作要做，其一是使用HIP编程接口实现DCU的设备管理、内存管理、DCU内存与主机内存间的数据传输等功能；另一方面是编写DCU核函数让DCU执行计算任务。我们之所以能够做到这些，是基于HIP提供了一套C/C++的编程接口，并且对C/C++语言做了扩展，使之能够编写DCU核函数。

对于Fortran编写的程序，虽然DTK中也提供了HIP的Fortran编程接口，但没有对Fortran语言提供扩展，不能使用Fortran语言编写核函数。因此如果我们想使用DCU加速Fortran编写的程序，如果是直接编写DCU核函数，只能使用C/C++来编写，再让Fortran程序调用核函数，换言之使用DCU加速Fortran程序实际上是一个Fortran/C互用性（interoperability）的问题；另一条路线则是使用OpenMP，用导语的方式指导编译器将并行计算部分卸载到DCU上。

本章主要介绍使用C/C++编写核函数并在Fortran程序中调用的方式，使用OpenMP的方式将在第5章进行介绍。

4.2 在Fortran程序中使用HIP编程接口

我们知道在C++程序中使用HIP编程接口时需要通过#include在代码中插入头文件hip_runtime.h，使用其它DTK中的库时还需要插入其他相应的头文件，如hiprand.h、hipblas.h等等（这些头文件都在DTK中的include目录下）。在Fortran程序中使用HIP编程接口和其他DTK库时也需要在程序通过USE语句使用相应的Fortran模块。

这些Fortran模块的源代码在23.04版本的DTK中的src/hipfor路径下，主要模块总结在表4-1中。

表4-1 hipfor中的主要模块

模块	主要内容
hipfort	HIP运行时编程接口
hipfort_types	dim3类型和一些枚举类型
hipfort_check	错误检查接口
hipfort_hipfft	hipFFT库的接口
hipfort_hipblas	hipBLAS库的接口
hipfort_hiprand	hipRAND库的接口
hipfort_hipsolver	hipSOLVER库的接口
hipfort_hipsparse	hipSPARSE库的接口

查看这些模块的源代码，你会发现这些模块包含着使用Fortran 2003标准中的内置模块iso_c_binding编写的大量Fortran接口，这些接口绑定到HIP的相关C语言接口，在Fortran中使用HIP接口就是以这种方式实现的。

我们结合Fortran程序来介绍，在3.5.1节中的C程序use_libs.cc中我们实现了使用DCU随机生成两个矩阵并计算矩阵乘。在程序中我们没有编写DCU核函数，只使用了HIP运行时编程接口和hipRAND、hipBLAS两个库中的函数，刚好适合作为本节例子。

我们仿照use_libs.cc的流程，用Fortran语言来实现如下(use_libs.f90)：

PROGRAM MATRIX_MULTIPLICATION
    USE hipfort
    USE hipfort_check
    USE hipfort_hiprand
    USE hipfort_hipblas
    IMPLICIT NONE
    INTEGER :: argc, N
    INTEGER(C_LONG_LONG) :: SEED
    INTEGER(C_SIZE_T) :: num_elements
    CHARACTER(LEN=80) :: argv
    TYPE(C_PTR) :: gen, start, stop_, handle
    REAL(4), TARGET :: dt
    REAL(8) :: alpha, beta
    REAL(8), ALLOCATABLE :: A(:,:), B(:,:), C(:,:)
    REAL(8), POINTER :: dA(:,:), dB(:,:), dC(:,:)
    INTEGER,PARAMETER :: max_size = 5 ! 打印的最大维度
    INTEGER :: TIME ! 对于flang需要声明

    ! 获取命令行参数个数
    argc = COMMAND_ARGUMENT_COUNT()

    IF (argc /= 1) THEN
        CALL GET_COMMAND_ARGUMENT(0,argv)
        WRITE (*,*) "Usage: ", TRIM(argv), " N"
        STOP 1
    ELSE
        CALL GET_COMMAND_ARGUMENT(1,argv)
        READ (argv, *) N
        IF (N < 2) THEN
            WRITE (*,*) "N must larger than 2."
            STOP 1
        END IF
    END IF

    ! 初始化hipRAND随机数发生器
    CALL hipCheck(hiprandCreateGenerator &
                  (gen, HIPRAND_RNG_PSEUDO_DEFAULT))
    SEED = TIME()
    CALL hipCheck(hiprandSetPseudoRandomGeneratorSeed(gen, SEED))

    ! 计时器
    CALL hipCheck(hipEventCreate(start))
    CALL hipCheck(hipEventCreate(stop_))

    ! 分配矩阵内存空间
    ALLOCATE(A(N,N), B(N,N), C(N,N))
    CALL hipCheck(hipMalloc(dA, N, N))
    CALL hipCheck(hipMalloc(dB, N, N))
    CALL hipCheck(hipMalloc(dC, N, N))

    ! 随机生成矩阵A和B
    CALL hipCheck(hipEventRecord(start,C_NULL_PTR))
    num_elements = N * N
    CALL hipCheck(hiprandGenerateUniformDouble &
                  (gen, C_LOC(dA), num_elements))
    CALL hipCheck(hiprandGenerateUniformDouble &
                  (gen, C_LOC(dB), num_elements))
    CALL hipCheck(hipEventRecord(stop_,C_NULL_PTR))
    CALL hipCheck(hipEventSynchronize(stop_))
    CALL hipCheck(hipEventElapsedTime(C_LOC(dt),start,stop_))
    WRITE(*,'("Generating Matrices took ", F8.3, " ms")') dt

    ! 打印矩阵A,B
    CALL hipCheck(hipMemcpy(A, dA, N * N, hipMemcpyDeviceToHost))
    CALL hipCheck(hipMemcpy(B, dB, N * N, hipMemcpyDeviceToHost))
    WRITE(*,*) "Matrix A:"
    CALL print_matrix(N, A)
    WRITE(*,*)
    WRITE(*,*) "Matrix B:"
    CALL print_matrix(N, B)
    WRITE(*,*)

    ! 创建hipBLAS句柄
    CALL hipCheck(hipblasCreate(handle))

    ! 矩阵乘法，将结果存储在矩阵C中
    alpha = 1.0
    beta = 0.0
    ! 预热
    CALL hipCheck(hipblasDgemm(handle, &
                  HIPBLAS_OP_N, HIPBLAS_OP_N, 2, 2, 2, &
                  alpha, dA, N, dB, N, beta, dC, N))
    ! 正式计算
    CALL hipCheck(hipEventRecord(start,C_NULL_PTR));
    CALL hipCheck(hipblasDgemm(handle, &
                  HIPBLAS_OP_N, HIPBLAS_OP_N, N, N, N, &
                  alpha, dA, N, dB, N, beta, dC, N))
    CALL hipCheck(hipEventRecord(stop_,C_NULL_PTR))
    CALL hipCheck(hipEventSynchronize(stop_))
    CALL hipCheck(hipEventElapsedTime(C_LOC(dt), start, stop_))
    WRITE(*,'("Matrix multiplication (hipblasDgemm) ", &
          "took " F8.3, " ms")') dt

    ! 打印矩阵C
    CALL hipCheck(hipMemcpy(C, dC, N * N, hipMemcpyDeviceToHost))
    WRITE(*,*) "Matrix C:"
    CALL print_matrix(N,C)

    ! 清理
    CALL hipCheck(hipFree(dA))
    CALL hipCheck(hipFree(dB))
    CALL hipCheck(hipFree(dC))
    CALL hipCheck(hipEventDestroy(start))
    CALL hipCheck(hipEventDestroy(stop_))
    CALL hipCheck(hiprandDestroyGenerator(gen))
    CALL hipCheck(hipblasDestroy(handle))
    DEALLOCATE(A,B,C)

    CONTAINS

    SUBROUTINE print_matrix(N, M)
        IMPLICIT NONE
        REAL(8) :: M(N,N)
        INTEGER :: N, i, j
        INTEGER,PARAMETER :: max_size = 5

        DO i = 1, MIN(N, max_size)
            DO j = 1, MIN(N, max_size)
                WRITE(*,'(F13.6)',ADVANCE='NO') M(i,j)
            ENDDO
            IF (N > max_size) THEN
                WRITE(*,'(" ...")',ADVANCE='NO')
            ENDIF
            WRITE(*,*)
        ENDDO
        IF (N > max_size) THEN
            WRITE(*,'("."/ "." / ".")')
        ENDIF
    END SUBROUTINE
END PROGRAM

在上面的代码中，我们使用了hipfort、hipfort_check、hipfort_hiprand和hipfort_hipblas四个Fortran模块，此后就可以像在C语言中一样方便的使用HIP运行时API和hipRAND、hipBLAS两个库中的函数了。

4.3 hipfor的使用方法

我们先编译上节中编写的Fortran程序use_libs.f90，在程序中我们通过USE语句使用了hipfor中的4个模块，但我们需要首先编译出这些模块。此时我们有两种方式，其一是在代码的第一行添加INCLUDE语句来插入hipfor的源代码：

INCLUDE `hipfort.f`

因为我们的程序以及hipfor的代码都是Fortran的标准代码，因此不需要像C代码那样只能使用hipcc进行编译，可以使用GCC中的gfortran、Intel的ifort或DTK中包含的flang。

DTK中之所以没有提供编译好的hipfor模块，原因在于Fortran模块在不同编译器甚至不同编译器版本间不通用。

我们使用gfortran编译如下(需要加载DTK-23.04环境，如使用其他DTK版本，需要将$ROCM_PATH替换为23.04版本根目录的绝对路径)：

gfortran use_libs.f90 -o use_libs -I$ROCM_PATH/src/hipfor -L$ROCM_PATH/lib -lamdhip64 -lhiprand -lhipfft -lhipblas -lhipsolver -lhipsparse

这里的-I$ROCM_PATH/src/hipfor是让编译器能找到hipfort.f文件，此外hipfor模块中绑定了大量的HIP编程接口，因此还需要链接一系列DTK中的库。

我们的程序只包含一个文件，按这种方式使用hipfor就足够了，如果是包含多个源文件的项目，则至少需要保证编译的第一个使用hipfor的文件中插入hipfort.f，否则编译器会报找不到模块的错误。

因此，更推荐按第二种方式使用hipfor，即首先将hipfor构建好。将DTK-23.04中src路径下的hipfor目录整个拷贝到工作路径下，比如说use_lib.f90所在的路径。最好将hipfort.f重命名为hipfort.f90，这样告诉编译器hipfor中的代码是使用自由格式编写的。

用下面的命令将hipfor编译出静态库libhipfort.a以及一系列.mod文件（以使用gfortran为例，换成ifort或flang同样可以）：

gfortran hipfort.f90 -c -o libhipfort.a

再进一步编译我们的程序，并与libhipfort.a链接在一起，注意到这次的-I./hipfor是为了让编译器找到.mod后缀的模块文件：

gfortran use_libs.f90 ./hipfor/libhipfort.a -o use_libs -I./hipfor -lamdhip64 -lhiprand -lhipfft -lhipblas -lhipsolver -lhipsparse

我们运行程序，结果与C语言编写的版本是一样的，都用上了DCU：

./use_libs 4096
Generating Matrices took    5.505 ms
 Matrix A:
     0.008903     0.629933     0.709649     0.638452     0.644310 ...
     0.342506     0.346941     0.301005     0.960350     0.240107 ...