Sparse Matrix-Vector Multiplication on GPU:

When Is Rows Reordering Worthwhile?

Paula Prata

Instituto de Telecomunicações

Departamento de Informática

Universidade da Beira Interior

João Muranho

Instituto de Telecomunicações

IMAR – Instituto do Mar,

Departamento de Ciências da Vida

Universidade de Coimbra

Motivação

Evolução da capacidade de cálculo das placas gráficas

(GPUs – Graphics Processing Units).

Desenvolvimento de interfaces de programação para GPU, que

permitem programar a placa gráfica com linguagens de alto nível

( C/C++, Phyton, Java, …):

- CUDA1 (NVIDIA), - Brook+ (AMD/ATI), - OpenCL.

Aumento da investigação sobre como usar a GPU para aplicações

não gráficas:

GPGPU – General Purpose computation on GPU.s .

GPGPU - Áreas de Aplicação

• Problemas com paralelismo de dados , o mesmo código é executado

simultaneamente em diferentes segmentos de dados.

• Exemplos :

- Simulação de modelos moleculares,

- Previsão meteorológica,

- Processamento de sinal, finanças, ...

• Problemas de cálculo cientifico que envolvem a manipulação de matrizes

de grande dimensão .

Vários estudos mostram que para problemas que manipulam matrizes

densas, a GPU permite grandes ganhos de desempenho.

Objetivos

• Estudo da operação: produto matriz esparsa - vector

Operação dominante em problemas de resolução de sistemas de

equações lineares, e no cálculo de valores próprios

Formato de representação de matrizes esparsas condiciona o

desempenho

• Análise do impacto da ordenação das linhas pelo número de elementos

não zero

Formato CSR

Formato ELL

Arquitectura GPU da NVIDIA: GeForce GTX 295

Array de multi-processadores

Arquitectura GPU da NVIDIA: GeForce GTX 295

Cada multiprocessador com:

8 processadores (cores)

Um conjunto de registos

Área de memória partilhada

Uma unidade de operações em virgula flutuante de

precisão dupla (capability 1.3)

Modelo de Programação CUDA

Um programa em execução na CPU (host) pode:

Copiar dados da CPU para a GPU e vice-versa

Lançar a execução de funções na GPU (Kernel.s)

Executar operações de sincronização

…

Cada Kernel é executado por múltiplas threads em simultâneo sobre

diferentes conjuntos de dados

Modelo de execução:

• Em cada multiprocessador – Simple Instruction Multiple Data

• Na GPU – Simple Program Multiple Data.

Modelo de Programação CUDA

Threads agrupadas em blocos

dimensionáveis pelo utilizador.

É criada uma grelha na qual

são distribuídos os blocos de

threads.

Modelo de Programação CUDA

• A grelha é associada à placa gráfica.

• Cada bloco é associado a um multiprocessador.

• As threads de um bloco são executadas pelos núcleos do

multiprocessador associado ao bloco.

• Unidade de

escalonamento ( warp) =32

Hardware e Linguagens

Intel Core 2 Quad Q9550 a 2.83 GHz, com 4 GB de RAM

Nvidia Geforce GTX 295

(30 multiprocessadores, 8 cores cada a 1,24 GHz, 2GB memória

global)

CUDA – versão 2.3

Visual Studio, C/C++

Matrizes: sintéticas e “Williams multi-core benchmarking”

Matrizes Esparsas – Formatos de Armazenamento

1 0 0 2

0 3 0 0

0 4 5 0

0 0 0 6 COO – Coordinate Format

Linhas 0 0 1 2 2 3

Colunas 0 3 1 1 2 3

Valores 1 2 3 4 5 6

Formato CSR – Compressed Sparse Row

0 2 3 5 6

0 3 1 1 2 3

1 2 3 4 5 6

ptr =

índices das colunas=

dados =

Os elementos são armazenados por linhas

Número de

não zeros

Formatos de Armazenamento: CSR

Formatos de Armazenamento: ELL - R

Formato ELLPACK/ITPACK – ELL (ELL-R)

Os elementos são armazenados por colunas:

0 3

1 *

1 2

3 *

1 2

3 *

4 5

6 *

dados = 2

1

2

1

índices = tamanho linhas =

1 3 4 6 2 * 5 *

Algoritmos Estudados

“Thread per row”

Cada linha da matriz é atribuída a uma thread

Formato CSR: as threads de cada warp acedem a posições de

memória não contíguas

Formato ELL: as threads de cada warp acedem a posições

contíguas de memória (mais eficiente)

Algoritmos Estudados

“Warp per row”

Cada linha da matriz é atribuída às threads de um warp

Format CSR: Todas as threads acedem a elementos da

mesma linha logo a posições contíguas de memória

Eficiente se as linhas tiverem tamanho suficiente para todas

as threads terem trabalho (>=32)

Ordenar linhas, porquê?

Modelo de execução SIMD =>

o desempenho é tanto maior quanto, num mesmo warp for:

- menor a divergência no acesso à memória

- menor a divergência de execução

Se num mesmo warp houver threads a processar linhas de diferentes

comprimentos, a execução do warp só termina quando terminar o

processamento da maior das linhas, isto é, da linha que tiver maior

número de valores não zero.

Resultados – matrizes sintéticas

Melhores tempos de execução (em milissegundos) obtidos para

matrizes com 10% de não zeros gerados aleatoriamente

Matrix

order

GPU, CSR (thread per row) GPU, ELL-R (thread per row) GPU, CSR

row sorted row sorted warp per row

float double float double float double float double float double

4096 4.24 4.73 3.72 4.42 1.40 1.44 1.23 1.31 0.83 1.08

8192 20.22 22.32 18.80 21.00 4.64 5.03 3.89 4.33 2.98 3.99

16384 93.37 100.66 88.14 94.49 17.55 19.28 14.11 15.55 11.45 15.31

Resultados – matrizes sintéticas

Precisão simples mais rápido que precisão dupla

“Thread per row” - há sempre ganho com a ordenação

O ganho com a ordenação das linhas é maior para o formato ELL-R

O algoritmo “warp per row” é sempre o melhor para estas matrizes

Resultados – matrizes sintéticas

Quando a percentagem de não zeros decresce,

Algoritmo “thread per row”, formato ELL-R, com ordenação das

linhas é o melhor quando % de não zeros <= 2%

Quando a percentagem de linhas, com tamanho <= 32, cresce

Algoritmo “thread per row”, formato ELL-R, com ordenação das

linhas é o melhor quando % linhas com tamanho <=32 é >= 70%

Resultados – Williams multi-core benchmarking

Matrix Order

CPU GPU GPU texture

CSR ELL

(std) ELL-R

ELL-R

sorted

CSR

wpr ELL-R

ELL-R

sorted

Economics 206500 5.24 2.10 1.60 0.945 3.81 1.21 0.708

Accelerator 121192 5.27 0.825 0.665 0.830 2.25 0.504 0.685

Cantilever 62451 7.10 0.867 0.751 0.976 1.35 0.588 0.622

Epidemiology 525825 7.02 0.776 0.745 0.753 7.54 0.685 0.687

Protein 36417 7.76 1.658 1.08 1.36 1.26 0.787 0.857

Spheres 83334 10.12 1.41 1.02 1.45 1.98 0.827 0.915

Ship 140874 14.35 2.26 1.54 2.84 2.87 1.23 2.13

Wind Tunnel 217918 20.60 6.56 2.01 3.56 4.27 1.55 2.14

Circuit 170998 4.619 9.90 2.11 0.919 3.07 1.78 0.793

Harbor 46835 4.578 1.69 1.27 1.03 1.41 0.811 0.614

Resultados – Williams multi-core benchmarking

Ordenar as linhas pelo seu tamanho tem vantagem para 3 matrizes:

Economics, Circuit, Harbor

Nos restantes casos, não ordenar as linhas tem melhores resultados

A perda de localidade no acesso ao vector é responsável pelo pior

desempenho da ordenação

O algoritmo “warp per row” nunca é o melhor

A utilização de texturas (memória constante) para armazenar o vector,

melhora os tempos de execução mas não inverte os resultados.

Resultados – Williams multi-core benchmarking

Matrix N. of nz % of

nz

Av.of

nz

/row

Sort.

Time

(ms)

GFlops

ELL-R

Texture

GFlops

ELL-R

Sorted

Texture

Average of warp lengths

Before

Sorting After Sorting

Economics 1273389 0.003 6 12.2 2.1 3.6 22.4 6.18 27%

Accelerator 1362087 0.009 22 7.3 5.4 4.0 16.1 11.2 70%

Cantilever 2034917 0.05 65 3.9 6.9 6.5 37.6 32.6 87%

Epidemiology 2100225 0.0003 4 9.8 6.1 6.1 3.99 3.99 100%

Protein 2190591 0.165 119 3.6 5.6 5.1 90.3 60.3 67%

Spheres 3046907 0.044 72 5.0 7.4 6.7 42.2 36.6 87%

Ship 3977139 0.015 28 10.5 6.5 3.7 40.0 28.3 71%

Wind Tunnel 5926171 0.012 53 5.3 7.6 5.5 31.0 27.2 88%

Circuit 958936 0.003 6 9.0 1.1 2.4 14.1 5.66 40%

Harbor 2374001 0.11 50 3.5 5.8 7.7 76.5 50.7 66%

Conclusões

Calculando o tamanho médio dos “warps” antes e depois de ordenar

as linhas, verifica-se que quando esse valor decresce para cerca de

66% ou menos do valor inicial, há vantagem em ordenar as linhas.

Nestes casos o facto de cada warp ter uma carga de trabalho mais

equilibrada compensa a falta de localidade no acesso ao vector.

Trabalho futuro

Estudar outras matrizes

Estudar comportamento dos algoritmos na nova arquitectura da

NVIDIA, Fermi.

Estudar outra arquitecturas e algoritmos