Implementation of FEM Application on GPU with
StarPU

Satoshi Ohshima, Satoshi Katagiri, Kengo Nakajima, Samuel Thibault,
Raymond Namyst

► To cite this version:

Satoshi Ohshima, Satoshi Katagiri, Kengo Nakajima, Samuel Thibault, Raymond Namyst. Implementation of FEM Application on GPU with StarPU. CS&E 2013 - SIAM Conference on Computational Science and Engineering, SIAM, Feb 2013, Boston, United States. hal-00926144

HAL Id: hal-00926144

https://inria.hal.science/hal-00926144v1

Submitted on 31 Jul 2025

**HAL** is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Satoshi Ohshima*, Takahiro Katagiri*, Kengo Nakajima*,
Samuel Thibault** and Raymond Namyst**

* The University of Tokyo, Japan

** University of Bordeaux, France

# Implementation of FEM Application
on GPU with StarPU
# Introduction : about our project

## JST-ANR collaboration project

* FP3C : Collaboration Project between Japan and France on Framework and Programming for Post Petascale Computing

* Japan side : U.Tsukuba, U.Tokyo, Titech, Kyoto U., RIKEN AICS

* France side : INRIA, CEA, IRIT, CNRS

## This research is based on one of the FP3C collaboration works.

* U.Tokyo : FEM application on GPU

* INRIA Bordeaux(U.Bordeaux 1) : StarPU runtime framework
  -> FEM application on GPU with StarPU

## This talk reports current status of this work.

* (To be honest, current performance is not excellent.)
# Outline

## Target application: FEM

- brief overview

- existing implementation and performance

## StarPU runtime framework

- short introduction about StarPU

## Implementation and Evaluation

- FEM on GPU with StarPU

- performance evaluation

- Conclusion and Future work
# Target application

- FEM code for 3D Solid Mechanics

- based on existing FEM application for CPU

- based on GeoFEM
    (http://geofem.tokyo.rist.or.jp/)

- already parallelized on CPU
    (MPI, OpenMP)

<figure><img src="image_2.png" /></figure>
# Original version of our FEM program (for CPU)

* simple and clear FEM program

* written by Fortran with OpenMP

* direct original version is by C with OpenMP

* Fortran with OpenMP -> C with OpenMP -> CUDA

* CRS format with 3x3 blocking

* special format, good for reducing cache miss

* CG solver with simple precondition method

* block SSOR and block diagonalization are available

* use block diagonalization in this talk

* (This is very simple method and doesn't have large effect for execution time. I don't talk about detail of this.)
Flow of our FEM application

<figure><img src="image_3.png" /></figure>
# Short summary of optimization on CUDA GPU

## important strategies on CUDA GPU

### high parallelism

* GPU has a lot of low performance processors.

* GPU can't obtain better performance than CPU in the case of low parallelism calculation.

### reduce the branch (WARP optimization)

* WARP is scheduling unit of CUDA.

* WARP $\div$ SIMD unit (32 CUDA Threads)

* per WARP calculation is well-known and important technique

### reduce data transfer between CPU and GPU

* PCIe is much slower than inner GPU memory performance

## implemented on Fermi GPU

* talked at some meetings (SIAM PP12, VECPAR12, etc.)
Existing FEM performance on CPU and GPU

<figure><img src="image_4.png" /></figure>

<figure><img src="image_5.png" /></figure>

(All CPU programs DON'T use hand-tuned SIMD optimization.)
# Current status of FEM on GPU

* GPU obtained good performance.

* especially in Sparse Matrix Solver part

* Next steps

* complex preconditioning

* optimization on latest CPU and GPU

* SandyBridge with AVX optimization vs Kepler

* because the performance is limited by memory performance, Kepler may beat SandyBridge

* utilize multiple GPUs

* Can we utilize StarPU ?
# About StarPU

http://runtime.bordeaux.inria.fr/StarPU/

* runtime framework system

* developed by runtime team at University of Bordeaux(INRIA Bordeaux)

* one of the France side group of Japan-France collaboration project

* describes tasks by using "codelet"

* supports CPU, CUDA, OpenCL, and Cell B.E.

* supports visualization tools

* NOT a auto parallelization system or something

* optimization of computation kernel is users' job

* there are few limit in kernel code

* latest interest is MPI support
Abstract image of using StarPU

<figure><img src="image_6.png" /></figure>
Expectation of StarPU

□ We want to obtain highest performance by utilizing all of processors (CPU+GPU) on target execution environment.

□ co-operation is important

□ for the compute environments in which the performance difference between CPU and GPU is not so large

□ for the applications that have various kernels with correspond performance specification

□ It is not very difficult to make good computation kernels for CPU and good computation kernels for GPU.

□ To assign the computation to best hardware considering data movement costs is more difficult.

□ (Maybe)StarPU can make it possible by realistic cost.

□ cost : labor, readability, time of implementation, and etc.
# Optimization techniques of StarPU

* make many tasks

* more than a number of processor cores

* reducing the synchronization between tasks

* use appropriate data property

* read, write, readwrite, reduction, scratch

* may reduce the number of synchronization and amount of memories

* don't access partitioned data at out of task region

* choose appropriate scheduling policy

* (default policy may be good for many cases)
# Implementation : FEM program on StarPU

## Existing FEM program on GPU

* single input data file

* some parts of FEM program has parallelism, GPU utilize the parallelism

* Matrix Assembly, Boundary Condition Setting, and Sparse Matrix Solver

## FEM program on StarPU

* implementation 1

* utilizes parallelism in single input data file

* implementation 2

* utilizes input data file level parallelism
# Flow of our FEM application

<figure><img src="image_7.png" /></figure>

result
Implementation 1

1. Input : single input data file based on OpenMP version program (same as existing FEM on GPU)

<figure><img src="image_8.png" /></figure>

2. Connectivity problem:
There are some global dependencies.

<figure><img src="image_9.png" /></figure>

(have to execute global summation)

-> performance penalty, memory amount
penalty the more problem size and
number of parallelism increase, the more
penalty increase

3. Matrix Assembly
-> performance penalty, memory amount
penalty the more problem size and
number of parallelism increase, the more
penalty increase

<figure><img src="image_10.png" /></figure>

4. Boundary Condition Setting PU utilizes these parallelisms, too.

<figure><img src="image_11.png" /></figure>

<figure><img src="image_12.png" /></figure>

5. Sparse Matrix Solver
These parts have parallelism in each part.
Existing implementation utilized it.

<figure><img src="image_13.png" /></figure>
Implementation 1 : detail

Pseudo code of original program

main(){
input_data();
connectivity();
matrix_assembly();
boundary_condition();
sparse_matrixCAMer();
output_data();
}

sparse_matrixCAMer() {
  for() {
    scale(vector);
    dot(vector);
    spmv matix, vector);
    etc.
    if(converged)break;
  }
}

vector

matrix

divide

<figure><img src="image_14.png" /></figure>

<figure><img src="image_15.png" /></figure>

-> task 1

-> task 2

-> task 3

<figure><img src="image_16.png" /></figure>

Each calculation handles large matrix and vector.

There are enough parallelism.

Original OpenMP program utilizes this parallelism.

Our StarPU implementation divides the matrix and
vector and assign them to tasks.
Implementation 2

1. Input : multiple input data file based on MPI version program

<figure><img src="image_17.png" /></figure>

<figure><img src="image_18.png" /></figure>

2. Connectivity

<figure><img src="image_19.png" /></figure>

<figure><img src="image_20.png" /></figure>

There are trivial parallelism.
StarPU utilizes this parallelism.

3. Matrix Assembly

<figure><img src="image_21.png" /></figure>

problem/question)
We have to divide target data before
computation. This division may be time
and memory penalty.

4. Boundary Condition Setting

<figure><img src="image_22.png" /></figure>

We DON'T discuss now.

5. Sparse Matrix Solver

<figure><img src="image_23.png" /></figure>
Implementation 2 : detail (example of CG method)

CG method

Compute $\mathbf{r}^{(0)} = \mathbf{b} - [\mathbf{A}]\mathbf{x}^{(0)}$

for i= 1, 2, ...

solve $[M]z^{(i-1)} = r^{(i-1)}$

$\rho_{i-1} = r^{(i-1)} z^{(i-1)}$

if i=1

$\mathbf{p}^{(1)} = \mathbf{z}^{(0)}$

else

$\beta_{i-1} = \rho_{i-1} / \rho_{i-2}$

$\mathbf{p}^{(i)} = \mathbf{z}^{(i-1)} + \beta_{i-1} \mathbf{z}^{(i)}$

endif

$\mathbf{q}^{(i)} = [\mathbf{A}]\mathbf{p}^{(i)}$ (SpMV)

$\alpha_i = \rho_{i-1} / \rho^{(i)}$

$\mathbf{x}^{(i)} = \mathbf{x}^{(i-1)} + \alpha_i \mathbf{p}^{(i)}$

$\mathbf{r}^{(i)} = \mathbf{r}^{(i-1)} - \alpha_i \mathbf{q}^{(i)}$

check convergence $|r|$

end

Original MPI program (on each process)

endif

`dataǒxchangebyusingMPI`

q (i) = [A]p (i) @ local data

$\alpha_i = \rho_{i-1} / \rho^{(i)}$ @ local data

<figure><img src="image_24.png" /></figure>

StarPU program (on main thread)

endif

`data_exchange_by_using_memcopy`

exec“q (i) = [A]p (i)”as StarPU task

exec“$\alpha_i = \rho_{i-1} / \rho^{(i)}$”as StarPU task
Performance 1 : use 1 socket CPU

■ Matrix Assembly

■ Sparse Matrix Solver

■ initialize GPU (only Existing CPU+GPU)

■ etc. (Includes I/O)

Westmere(X5650)*2

+ Fermi(C2050)*2

gcc 4.4.6, cuda 4.2

Starpu 1.0.4

multiple input data file

<figure><img src="image_25.png" /></figure>

fast
Performance 1: use 1 socket CPU

Westmere(X5650)*2
+ Fermi(C2050)*2
gcc 4.4.6, cuda 4.2
Starpu 1.0.4

poor performance
some overheads?

single input data
file

good performance
vs hand-tuned?

Multiple input data file

■ Matrix Assembly

■ Sparse Matrix G

■ initialize GPU

■ etc. (Includes

<figure><img src="image_26.png" /></figure>

fast
Performance 2 : use 2 sockets CPUs

Westmere(X5650)*2
+ Fermi(C2050)*2
gcc 4.4.6, cuda 4.2
Starpu 1.0.4
multiple input data file

- Matrix Assembly
- Sparse Matrix Solver
- initialize GPU (only Existing CPU+GPU)
- etc. (Includes I/O)

<figure><img src="image_27.png" /></figure>
Performance 2 : use 2 sockets CPUs

speed up
but there is a slight effectiveness
hand-tuned is more faster

speed down
maybe, data transfer between CPU sockets
is one of the important overhead
(effect of NUMA architecture)

<figure><img src="image_28.png" /></figure>

fast
Compute, Idle, Data Fetching

ViTE http://vite.gforge.inria.fr/

StarPU execution trace (1/5)

CPU 4cores+2GPU, single input data file, whole application

<figure><img src="image_29.png" /></figure>
## Compute, Idle, Data Fetching

# StarPU execution trace (2/5)

CPU 4cores+2GPU, single input data file

<figure><img src="image_30.png" /></figure>
## Compute, Idle, Data Fetching

# StarPU execution trace (3/5)

## CPU 10cores+2GPU, single input data file

<figure><img src="image_31.png" /></figure>
Compute, Idle, Data Fetching

<figure><img src="image_32.png" /><figcaption># StarPU execution trace (4/5)</figcaption></figure>

CPU 4cores+2GPU, multi input data file
Compute, Idle, Data Fetching

StarPU execution trace (5/5)

CPU 10cores+2GPU, multi input data file

<figure><img src="image_33.png" /></figure>
# Conclusion & Future Work

We implemented our FEM application on GPU with StarPU.

* We obtained limited performance improvement.

* Excellent improvement is not obtained.

* Visualization (trace) makes some performance issues clear.

* Future work: improve the performance

* reduce the overhead

* single input data file version:

* reduce the idle time and fetching time

* multiple input data file version:

* compare the performance with hand-tuned 2GPUs program

* How to improve performance on 2 sockets CPUs ?

* How to resolve the NUMA issue ?
Questions?

contact:
Satoshi OHSHIMA
ohshima@cc.u-tokyo.ac.jp