[2024 Summer GeekPie HPC Lecture Notes] Week 2: Running

2024-07-19
13 min read
Notes
2024 Summer GeekPie HPC Lectures

[2024 Summer GeekPie HPC Lecture Notes] Catalog of Notes

调度器 (Scheduler)

k8s/PBS/slurm

评测 Benchmark

1. 评测设备性能

  • 工具：HPL/HPCG
  • 性能单位：FLOPS
  • 机器学习模型 工具：MLPerf Inference

可以从给的.dat 摘过来 HPL

2. 程序性能

评测程序到底优化得怎么样

收集设备的三种状况：Core Bound, Memory Bound, I/O Bound

  • Core Bound：这里体现了程序运行时的CPU占用情况，通常以计算时间、CPI (Cycles per Instruction) 等来衡量。影响计算时间的因素有很多，包括你的算法、你的并行度，以及CPI。如果你已经对体系结构有所了解，你会知道不同的指令会对CPI有一定的影响，Cache Miss率也会影响CPI。
  • Memory Bound：这里体现了程序运行时的访存效率。程序访问“距离”CPU越近的数据越频繁，这里的表现就越好。依据不同元件“距离”CPU的远近，我们可以有广义上的Cache：寄存器是CPU缓存的Cache，CPU缓存是DRAM的Cache，DRAM是硬盘的Cache。
  • I/O Bound：这里体现了程序运行时，不同进程间的通信效率，以及程序单纯做输入输出时的效率。我们“几乎”可以忽略程序输入输出的影响（除非程序需要从硬盘中读入大量数据）。而程序以多进程的方式运行时，不同进程之间会有通信。这会受到不同节点之间的带宽，当前的数据传输量影响。当短时间内有大量进程相互通信时，这些信息就会挤占带宽，从而影响整体的数据传输效率（想想DDos）。

Weak Scaling & Strong Scaling

strong scaling: 固定测试集规模，测试不同设备规模

weak scaling: 固定设备规模，测试不同测试集规模下的表现

Parallel: SIMD

数据级别并行

SSE/AVX/AVX-512

void mul(int N, float* A1, float* A2, float* B){
    int i = 0;
    for(; i < N - 7; i+=8){
        __m256 a1 = _mm256_loadu_ps(A1 + i);
        __m256 a2 = _mm256_loadu_ps(A2 + i);
        __m256 b = _mm256_mul_ps(a1, a2);
        _mm256_storeu_ps(B + i, b);
    }
    
    // tail: 尾巴处理
    for(; i < N; ++i){
        B[i] = A1[i] * A2[i];
    }
}

练习：110 lab11

Parallel: OpenMP

线程级别并行

线程 (Thread)，是能够被调度器调度，运行在CPU核心上的最小指令序列单元。一个进程 (Process) 可以包含多个线程。而线程所拥有的内存，既可以被自己单独使用，即私有 (private)，也可以在多个线程之间共享 (shared)。这样一来，每个线程就可以在私有的内存中完成自己的任务，同时在共享的内存中进行通信，或是读入数据，或是写入数据。

线程数量并不是越多越好

最好：它应当是当前CPU所能支持的同时运行的最大线程数。

当线程数超过这个值时，根据线程调度策略不同，有可能一个核心会在两个不同线程之间来回运行。这时候就会涉及到上下文切换 (Context Switch)了。这个过程是非常耗时的。

  • 同步与通信：线程使用的资源可能是受限的，即需要使用锁 (Lock) 等工具来控制访问。在两个线程同时访问这些资源时，一个线程需要等待另一个线程访问完毕才能使用。这在写入共享内存当中的数据非常常见。
  • 负载不均衡。每个线程有时并不能负担相同的任务量。假如你有8组数据，每组数据内部不可分割。当你用7个线程进行计算时，根据鸽笼原理，一定会有一个线程需要计算两组数据。这时候，并行计算的速度很大程度上就取决于这个需要计算两组数据的倒霉蛋了。
  • 内存带宽限制。当同时访问内存的线程数目很大时，我们就不得不考虑内存带宽限制的问题了。
  • 缓存争用**（False Sharing）。我们知道，缓存是以行 (Cache line) 为单位读入内存的，而一个缓存行会包含多个变量。因此，当不同线程操作内存地址在同一缓存行的变量时（其中至少要写入一个变量），如果是在不同线程之间依次高频操作**，为了保持缓存一致性，CPU会反复读入/驱逐这一条缓存行，造成极大的访存开销。此外，如果是不同线程同时对该缓存行中的变量进行上述操作，则会发生伪共享 (False Sharing)。此时，这条缓存行只能被一个线程操作，如同使用了一把锁来控制。

OpenMP & pThread

为了实现线程级并行，我们可以使用诸如pThread之类的库。而这里我们介绍的，是一个更为方便的库：OpenMP。

gcc: -fopenmp

intel icc: -qopenmp

神必纸质资料 (OpenMP 5.2 的官方文档)

OpenMPRefGuide-5.2-Web-2024.pdf

Cannon’s Algo (CS210 ShanghaiTech - Parallel Algorithms)

矩阵乘法

把整个矩阵做个分块，每个分块划分到线程或者进程执行

OpenMP Directives

Implicit Barriers/for/section/sections/Explicit Barriers

for/section：等待到最慢的一个结束 再进行

  • explicit barrier barrier: 写一个显式barrier让所有线程在这里等最慢的，能不用就不用。

atomic/critical/reduction

两个线程在运行之后，有可能两个都读到同一个内存，在序列化模型中，两个内存地址上的值应该随序列来更新。但是，现在读到了同一个没更新的值。<- 不可预知的计算结果

解决方法：原子操作/规约

  • atomic: 将对应的一条 对共享变量操作的语句 标记为原子操作。即同时只能被一个线程操作，例如：

    int x = 1;
#pragma omp parallel
{
    #pragma omp atomic
        x = 2 * x + 1;
}

    只能执行很simple的操作。

  • critical: 同atomic，但是能执行一个相对复杂的代码块。

  • reduction([<operator> : <variables>, ]<operator> : <variables>)：并行计算某个变量。<operator>只能是以下运算符，且不支持重载：+, *, &, /, ^, &&, ||。<variables>为需要进行并行计算的变量。如果有多个变量，以逗号,分隔。例如，对于这样一个表达式：sum = a[1] + a[2] + ... + a[8]，可以通过reduction来并行计算：

    #pragma omp parallel for reduction(+ : sum)
for(int i = 0; i < 8; ++i){
    sum += a[i];
}

    此时，假如有2个线程，以上代码就会在每个线程中分别计算i = [0, 3], i = [4, 7]时sum的值，最后将每个线程计算出的值进行汇总。事实上，以上代码的功能也可以通过atomic或critical实现。请在练习中进行实验。

single/nowait/master

  • single：将对应语句标记为单线程执行。在遇到single标记时，其他线程会等待执行single部分的线程，除非在[Clause]处标记nowait。

  • master：将对应语句标记为由主线程执行，这通常涉及到一些I/O操作等。在遇到master标记时，其他线程不会等待主线程运行结束。（没有implicit barrier）

  • nowait：以非阻塞的方式运行这一段并行块。在single中是不等待执行这一代码块的单线程，而在如sections, for的代码块，则是不在这样的代码块末尾等待其他线程结束，直接进行下一项工作。

OpenMP Function API/Env Vars

Parallel: MPI

Message Passing Interface

**进程级别并行 ** SPMD: Single Program Multiple Data

多节点上数据通信的规范，只是协议/接口。节点：集群上的机器

进程之间不共享内存，所以MPI的编程模型是分布式内存模型，让每个进程单独管理自己内存；MPI用来把这些单独内存通信连接。

OpenMP的局限

OpenMP：线程并行，逻辑核 <= 400

物理核：实际有的；

逻辑核：不断切换造成的假象

更大的数据集，OMP+MPI可以一起使用

安装

MPICH (mpicc, mpicxx)

spack install mpich
spack load mpich

Intel MPI (mpiicc, mpiicxx …)

spack install intel-oneapi-compilers
spack load intel-oneapi-mpi@2021.11.0 ## 2021.13.0?
mpicc -show

mpirun 和run的时候可修改的环境变量

讲师也不知道。看教程

<mpi.h> APIs

MPI API的类型大致可以分为四类。历史遗留问题->传参需求特别大

本身只是通信上的东西，没有优化。只是负责了信号的传输，对于MPI的优化：就是如何让传输的延迟更低；如何掩盖传输 计算

docs.open-mpi.org

1. 基本函数

#include <mpi.h>
// init MPI env, 传入main的
int MPI_Init(int *argc, char **argv[]);

// end MPI env
int MPI_Finalize(void);

// get rank id // 获得进程编号，放到rank
int MPI_Comm_rank(MPI_Comm comm, int *rank);

// get size of ranks // 获得进程有几个，放到size
int MPI_Comm_size(MPI_Comm comm, int *size);

2. 点对点传递和相关

阻塞性：send和recv(receive)操作做完了才下一步

a->b，需要先send后recv。缺了也会一直卡着

流程可以是这样的：process a 先执行 send; process b recv; process a 下一步；process b 下一步

也可以先process b 先执行到recv，这个时候等a执行send。这之后，process a/b 下一步。

// 将发送缓冲区buf中count个datatype数据类型的数据
// 发送到comm中标识号为dst的目标进程，本次发送标识符为tag
// input: buf, count, datatype, dest, tag, comm

int MPI_Send(const void *buf, int count, MPI_Datatype datatype, 
        int dest,int tag, MPI_Comm comm)

// 从comm中标识符为source的目标进程
// 接收count个datatype数据类型的数据
// 接收数据入缓冲区buf中，本次接收标识符为tag
// input: count, datatype, source, tag, comm
// output: buf, status

int MPI_Recv(void *buf, int count, MPI_Datatype datatype,
        int source, int tag, MPI_Comm comm, MPI_Status *status)

非阻塞性：就没死锁了

img

// 相比MPI_Send多一个参数request用于检查操作完成
int MPI_Isend(const void *buf, int count, MPI_Datatype datatype, 
        int dest,int tag, MPI_Comm comm, MPI_Request *request)
    
// 相比MPI_Recv多一个参数request用于检查操作完成
int MPI_Irecv(void *buf, int count, MPI_Datatype datatype,
        int source, int tag, MPI_Comm comm, MPI_Request *request)
    
// 等待MPI发送或接收结束，然后返回
int MPI_Wait(MPI_Request *request, MPI_Status *status)

// 若flag为true，代表操作完成。反之代表操作未完成
// 该函数会立即返回
// input: request
// output: flag, status
int MPI_Test(MPI_Request *request, int *flag, MPI_Status *status)

3. 集合传递

4. 数据类型函数

用于**自定义MPI消息传递的数据类型**，从而减少消息传递次数，增大通信粒度，同时可以减少消息传递时数据在内存中的拷贝。

MPI只有一点数据类型，我们有时候要传些vector之类的就要在这里写

样例代码之后贴

常用函数及其参数如下。

MPI HelloWorld

#include <mpi.h>
#include <string.h>

int main(int argc, char **argv)
{
    // preprocessing
    char message[20];
    int myrank;
    int size;
    MPI_Status status;

    // MPI initialization
    MPI_Init( &argc, &argv );

    // Compute and Message passing
    MPI_Comm_rank( MPI_COMM_WORLD, &myrank );
    MPI_Comm_size( MPI_COMM_WORLD, &size);
    if (myrank == 0) /* code for process zero */
    {
        strcpy(message,"Hello");
	for(int i = 1; i < size; i++)
	{
		printf("Send to rank %d from rank %d\n", i, myrank);
		MPI_Send(message, strlen(message)+1, MPI_CHAR, i, 99, MPI_COMM_WORLD);
	}
    }
    else /* code for process one */
    {
        MPI_Recv(message, 20, MPI_CHAR, 0, 99, MPI_COMM_WORLD, &status);
        printf("received :%s: by rank %d\n", message, myrank);
    }

    // MPI finalization
    MPI_Finalize();

    // postprocessing
    return 0;
}

MPI Deadlock

CA3课程里

练习：并行计算 BBP计算pi

Parallel: CUDA

编译器：nvcc

CUDA环境变量

  • CUDA_VISIBLE_DEVICES=1,2，表示序号1，2的设备对程序可见。
  • CUDA_VISIBLE_DEVICES=^1，表示除了序号1设备，其他设备对程序可见。
  • CUDA_VISIBLE_DEVICES=-1，表示禁用所有设备

CUDA程序结构

典型的CUDA程序的结构如下

  1. 设置device
  2. 在device/host上分配内存
  3. 从host上拷贝原始数据至device
  4. 在host上调用核函数
  5. 从device上拷贝处理后数据至host
  6. 释放device/host上的内存

使用CUDA编写的VectorAdd程序

// 这是驱动库
// #include <cuda.h>
// 这是运行时库
#include <cuda_runtime.h>
#include <iostream>

__device__ int add(int a, int b)
{
    return a + b;
}

__global__ void VectorAdd(int *a, int *b, int *c, int n)
{
    int i = threadIdx.x;
    if (i < n)
    {
        c[i] = add(a[i], b[i]);
    }
}

__host__ void print_ans(int *a, int n)
{
    for (int i = 0; i < n; i++)
    {
        std::cout << a[i] << " ";
    }
    std::cout << std::endl;
}

int main(int argc, char const *argv[])
{
    // 1. 设置device
    // 不调用则默认为0
    cudaSetDevice(0); // device id! `nvidia-smi` 可以看到

    int n = 10;
    int *a, *b, *c;
    int *d_a, *d_b, *d_c;

    // 2. 在device上分配内存
    a = (int *)malloc(n * sizeof(int)); // host
    b = (int *)malloc(n * sizeof(int)); // host
    c = (int *)malloc(n * sizeof(int)); // host

    cudaMalloc(&d_a, n * sizeof(int)); // gpu device
    cudaMalloc(&d_b, n * sizeof(int));
    cudaMalloc(&d_c, n * sizeof(int));

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

    // 3. 从host上拷贝原始数据至device
    cudaMemcpy(d_a, a, n * sizeof(int), cudaMemcpyHostToDevice);
    cudaMemcpy(d_b, b, n * sizeof(int), cudaMemcpyHostToDevice);

    // 4. 在host上调用核函数
    VectorAdd<<<1, n>>>(d_a, d_b, d_c, n);

    // 5. 从device上拷贝处理后数据至host
    cudaMemcpy(c, d_c, n * sizeof(int), cudaMemcpyDeviceToHost);

    print_ans(c, n);

    // 6. 释放device上的内存
    free(a);
    free(b);
    free(c);

    cudaFree(d_a);
    cudaFree(d_b);
    cudaFree(d_c);

    return 0;
}

上述代码可以用nvcc $filename进行编译。运行后能够在命令行中看到如下输出

0 2 4 6 8 10 12 14 16 18。

APIs

API: cudaSetDevice(id); cudaMalloc(&device_var, size); cudaMemcpy; cudaFree

CUDA对于C代码的扩展

接着我们来解释下CUDA代码中的特别之处

  • 函数类型 <- Specifiers

    在CUDA编程中，每个函数都会通过Function Execution Space Specifiers标记函数的执行方与调用方。有如下：

    • __device__: 函数由device执行，由device调用
    • __global__：函数由device执行，由host调用。而这样定义的函数一般被称为核函数（kernel）
    • __host__: 函数由host执行，由host调用（当函数没有declaration specifier时默认为该项）

  • 核函数调用 <—- __global__ 函数

    从上述代码我们可以发现核函数在host上的调用是很特殊的。需要将<<< GridDim, BlockDim >>>放在函数名与传参之间，代表了核函数的线程执行配置，其含义与接下来的GPU线程模型相关。

    核函数的调用在host上是立即返回的。而cudaMemcpy等用于内存拷贝的函数会**隐式等待之前提交的所有核函数完成**再执行，并且该函数会等待操作执行完成之后再返回。

    • （host==CPU可以继续执行接下来的任务，和之前MPI的非阻塞执行一样。）
    • 其他的CUDA API都是阻塞执行

  • 设备内存分配

    device上的内存需要在host代码中显式分配。上述代码中通过cudaMalloc分配。所分配的内存被称为设备内存。此处是设备内存中的全局内存。除了cudaMalloc之外，CUDA API还提供了其他类型的分配在设备内存上的全局内存，如cudaHostAlloc用于分配页锁定内存。除了分配设备内存还可以分配共享内存等其他类型。具体会在内存模型中展开。

CUDA的线程模型：三维结构

dim3 BlockDim, GridDim;
dim3 BlockDim(2,2); // 也可以缺省维数

SM & Warp

示意图中的一行（一黄一紫一堆绿的一行）代表一个SM(streaming multiprocessor)

CPU-GPU

分配之后，SM会将Block内的Thread映射至CUDA Cores上执行，并且以Warp(文档第一段)的形式管理、调度、执行任务。

Warp先得到thread的一维索引：threadIdx.z * blockDim.x * blockDim.y + threadIdx.y * blockDim.x + threadIdx.x

一个warp是由32个并行执行的线程组成的组。这些线程在同一个SIMT（单指令多线程）单元中同步执行同一条指令，但操作不同的数据。

Thread会根据BlockDim组织为一个Block。一个Block内Thread可以通过threadIdx(.x/y/z)唯一确定。一个Block大小可以从blockDim(.x/y/z)获得。

Warp Divergence

单个warp中的threads只会同时做一种计算，所以要是threads进入了不同分支，其中一个分支的threads会被暂停，成为序列化执行，降低计算效率。

GPU的内存分割

Global Memory 全局内存

  • 可以在全局作用域用__device__(线性内存)或__constant__(常量内存)修饰符声明变量。

  • 若需要动态分配可以在host代码中使用先声明需要的指针变量，在host上通过cudaMalloc，cudaFree管理。

    • 静态分配：全局作用域（最外面）__device__ deviceVar; __constant__ constantVar;
    • 动态分配：在main()/__host__代码中声明，如float *devArray，然后在host上通过cudaMalloc分配、通过cudaFree管理
    #include <cuda_runtime.h>
#include <iostream>

// 声明全局内存中的变量  **静态分配**
__device__ float globalVar = 10.0f;

// 声明常量内存中的变量  **静态分配**
__constant__ float constVar = 20.0f;

// 核函数，使用全局变量和常量变量
__global__ void useGlobalAndConstVar(float *output, int dataSize) 
{
    int idx = threadIdx.x;
    if (idx < dataSize) 
    {
            output[idx] = globalVar + constVar;
    }
}

int main() 
{
    const int dataSize = 256;
    float hostArray[dataSize];

    // 分配设备内存  **动态分配全局内存**
    float *devArray;
    cudaMalloc((void **)&devArray, dataSize * sizeof(float));

    // 启动核函数
    useGlobalAndConstVar<<<1, dataSize>>>(devArray, dataSize);

    // 将结果从设备内存复制回主机内存
    cudaMemcpy(hostArray, devArray, dataSize * sizeof(float),cudaMemcpyDeviceToHost);

    // 打印结果
    for (int i = 0; i < 10; ++i) 
    { // 仅打印前10个结果
        std::cout << "Result[" << i << "]: " << hostArray[i] << std::endl;
    }

    // 清理设备内存
    cudaFree(devArray);

    return 0;
}

  • 如果需要在host控制在全局作用域中声明的变量则需要cudaGetSymbolAddress() / cudaGetSymbolSize() / cudaMemcpyToSymbol() / cudaMemcpyFromSymbol()

  • 全局内存还包括了纹理内存(Texture Memory)和表面内存(Surface Memory)。他们具有特殊优化的缓存结构。还有统一虚拟内存(UVA)，用于控制多设备与主机内存访问。Device Memory还有很多种，详细了解可以参考官网。

Shared Memory 共享内存

CUDA程序的重要优化点：共享内存的读写比全局内存快多了，所以可以让共享内存中介全局内存的直接读写！

静态声明：在

__global__ void exampleKernel() {
    // 定义一个共享内存变量
    __shared__ float sharedVariable[256];

    // 使用 sharedVariable 进行计算...
}

动态声明：

__global__ void dynamicSharedMemoryKernel(int* data, int dataSize) {
    // 动态分配共享内存
    extern __shared__ int sharedData[];

    // 使用 sharedData 进行计算...
    // 例如，将每个线程的数据复制到共享内存中
    int tid = threadIdx.x;
    if (tid < dataSize) {
        sharedData[tid] = data[tid];
    }

    // 确保所有线程都完成了共享内存的写入
    __syncthreads();

    // 其他使用共享内存的计算...
}   

// 在主函数中调用 kernel，指定动态共享内存的大小
int main() {
    int* data;
    int dataSize = 256; // 假设数据大小为 256
    int sharedMemSize = dataSize * sizeof(int);

    // 分配 data 空间、初始化 data、拷贝数据到设备等操作...

    // 启动 kernel，最后一个参数指定动态共享内存的大小
    dynamicSharedMemoryKernel<<<1, dataSize, sharedMemSize>>>(data, dataSize);

    // 清理资源等操作...

    return 0;
}
  • 作用域是Block内的所有线程。
  • 生命周期直到Block的销毁，一般是核函数的结束
  • 由于共享内存物理设计上离核更近（有些结构上与L1共用片上缓存，其配比由核函数间接决定，也可以通过API手动控制，不过不能完全控制），其大小有限但访问速度十分快。一般是优化CUDA程序访存的重点。

Local Memory 本地内存

对于每个线程私有，定义方式：在核函数内直接定义。

__global__ void exampleKernel()
{
    // 定义一个 local memory 变量
    float localVariable = 1.0f;    // 如果加 __shared__ 就是 shared memory

    // 使用 localVariable 进行计算
}

OpenACC

没讲

#HPC #GeekPie #MPI #OpenMP #CUDA #Parallel Computing #Process-Level Parallelism #GPU Parallelism & CUDA #Thread-Level Parallelism
Next [Parallel Algorithm] Parallel Prefix Sum Algorithm, mainly by MPI

Portals

[2024 Summer GeekPie HPC Lecture Notes] Week 3 Exercise
[MPI] MPI Matrix Multiplication
[MPI] MPI Scan explained
[2024 Summer GeekPie HPC Lecture Notes] Week 1 Exercise
[2024 Summer GeekPie HPC Lecture Notes] Catalog of Notes
[2024 Summer GeekPie HPC Lecture Notes] Week 2: Running
[Parallel Algorithm] Parallel Prefix Sum Algorithm, mainly by MPI

See Also

[Parallel Algorithm] Parallel Prefix Sum Algorithm, mainly by MPI
[2024 Summer GeekPie HPC Lecture Notes] Catalog of Notes