# HPC et calcul des structures : de la difficulté de rendre ultra-scalable un code généraliste

Thomas de Soza, Nicolas Sellenet

► To cite this version:

Thomas de Soza, Nicolas Sellenet. HPC et calcul des structures : de la difficulté de rendre ultra-scalable un code généraliste. 11e colloque national en calcul des structures, CSMA, May 2013, Giens, France. hal-01717096

HAL Id: hal-01717096

https://hal.science/hal-01717096

Submitted on 25 Feb 2018

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.

Public Domain
# HPC et calcul des structures : de la difficulté de rendre *ultra-scalable* un code généraliste

Thomas DE SOZA <fnref n="1" />, Nicolas SELLENET <fnref n="2" />

**Résumé** — Les codes implicites de calcul des structures ont jusqu'à maintenant utilisé de manière contrastée le calcul parallèle. Les codes généralistes n'utilisent souvent pas plus de quelques dizaines de processeurs alors que des codes de recherche sont capables de s'exécuter efficacement sur plusieurs milliers de cœurs. *Code_Aster*, logiciel libre de calcul des structures développé par EDF-R&D n'échappe pas à cette règle : cette communication illustre la difficulté à le rendre *ultra-scalable* qui a conduit dans un premier temps à abandonner la généricité. Un calcul d'élasticité linéaire sur un cube à 500 millions de degrés de liberté a alors pu être résolu sur 700 processeurs tandis que l'étude d'une pompe de géométrie complexe comprenant 40 millions de degrés de liberté a nécessité moins de 5 minutes.

**Mots clés** — calcul parallèle, mécanique des structures implicite, *Code_Aster*, PETSc, multigrille, solveur linéaire

## 1 Introduction

The R&D department of French utility Electricité de France (EDF) has been developing for more than 20 years an all-purpose numerical simulation software for implicit structural mechanics released under the GNU GPL license [1]. Widely used in-house for verification and maintenance of assets and power plants, this Finite Element Analysis (FEA) solution covers a large range of applications such as 3D thermal and mechanical analysis in linear and non-linear statics and dynamics. It can be applied to machines, pressure vessels, pipes, civil engineering structures and so on.

Beyond the standard functionalities of an all-purpose FEA software for structural mechanics, *Code_Aster* is also the placeholder for knowledge sharing of numerical models related to the strong research involvement of EDF in various fields : fatigue, damage, fracture, contact, earth materials and porous media, multi-physics coupling, vibration analysis and rotating machinery.

## 2 Background and motivation

Nowadays numerical studies in implicit mechanics involve larger and larger models to simulate more and more complex phenomenon related to the lifespan of a structure especially in the nuclear sector. Time to solution for these problems may be up to a few days which hinders any effort towards stochastic analysis. Sometimes due to the huge size of the numerical model they may even exceed memory requirements preventing them to be easily run in the engineering units.

Meanwhile High Performance Computing (HPC) infrastructures get more and more cores, with fewer memory per core. Moreover, the performance per core may not see tremendous improvements in the near future and may even decrease in favour of many cores on one socket.

Traditionally, the field of implicit structural mechanics has not been very implicated in HPC and while most all-purpose FEA software take advantage of parallelism, they aren't massively parallel contrary to computational fluid dynamics codes or explicit dynamics software. That is to say calculations hardly ever use more than 100 cores and are usually efficient only up to 32.

There are examples of highly parallel FEA software in mechanics that are able to scale up to a thousand cores such as Elmer [2]. However they often cover neither the wide range of features offered by an

<footnote n="1">¹ EDF R&D, Analyses Mécaniques et Acoustique, thomas.de-soza at edf.fr</footnote>

<footnote n="2">² EDF R&D, Analyses Mécaniques et Acoustique, nicolas.sellenet at edf.fr</footnote>
<figure><img src="image_1.png" /><figcaption>Fig. 1: A simple refinement may lead to a rapid increase of degrees of freedom (20 000 dofs in the case of one screw)</figcaption></figure>

all-purpose software nor the specificities of specialized fields such as civil engineering or porous media. As such, EDF R&D faces a real challenge: on the one hand ensure that in the years to come *Code_Aster* engages on the road towards massively scalable software. On the other hand, make sure that current needs in studies' performance are met by perfecting the current state of parallelism in the code. This talk presents how these objectives are being tackled in an ever-evolving package such as *Code_Aster*.

## 3 Present capability of *Code_Aster*

For a long time, FEA software in structural analysis used to be sequential from mesh input to linear solving as schematized in Figure 2(a). When it comes to distributed memory parallelism, several approaches exist. The most common is to adopt parallelism at the solving phase (matrices building distribution and linear solver parallelism) with a master process still reading a sequential mesh (see Figure 2(b)). This approach enables to greatly decrease wall-clock time since the solving phase is usually the most CPU-consuming part of a computation. Moreover this approach guarantees most of the features of the software will remain available in parallel.

Nowadays, *Code_Aster* enables parallel processing at the following steps: computation and assembly of matrices and vectors, solving the system of equations. Every solving operator is available including non-linear statics and dynamics. Parallelism therefore also benefits to advanced features such as contact (though geometrical pairing of contact zones is not parallel).

<figure><img src="image_2.png" /><figcaption>Fig. 2: Different models of parallelism. (a) sequential. (b) parallel element computations and solve. (c) fully distributed approach</figcaption></figure>
On a simple yet large and non-linear contact problem as shown in Figure 3 the computation times presented in Table 1 are obtained. These were done on a cluster composed of Intel Xeon X5570 nodes (2 x 2.93Ghz quad-cores) running Debian Lenny 64bits. The peculiarity of this case is that the number of contacts represents more than 10% of the total number of degrees of freedom.

<figure><img src="image_3.png" /><figcaption>Fig. 3: Contact between two beams and a rigid foundation</figcaption></figure>

<table><caption>Table 1: Wall-clock time (s) on cases with increasing size</caption><thead><tr><th>Case</th><th>1 core</th><th>4 cores</th><th>8 cores</th></tr></thead><tbody><tr><td>46 000 dofs ( 4 300 contacts)</td><td>90</td><td>-</td><td>-</td></tr><tr><td>205 000 dofs (20 500 contacts)</td><td>-</td><td>390</td><td>-</td></tr><tr><td>320 000 dofs (32 000 contacts)</td><td>770</td><td>580</td><td>430</td></tr></tbody></table>

In this example the elementary computations and assemblies display almost ideal speed-up and the drop in efficiency is mainly due to the linear solver which is known to be the bottleneck in structural mechanics. Indeed, even when relying on iterative linear solvers, preconditioning usually requires working with direct solvers especially for ill-conditioned problems arising in structural mechanics.

The non-linear solver for this problem was a generalized Newton method for contact relying on inexact Krylov resolutions [3]. The linear solver was GMRES preconditioned by a single precision factorization of the tangent matrix.

# 4 Paving the way for a massively parallel FEA software

## 4.1 Description of the work

The common type of parallelism described in the previous section is already beneficial to the end user. However it cannot address very big models where memory is the limiting factor.

In order for a FEA software to be fully scalable, task distribution must occur at the very beginning of the calculation and not be delayed till the solving step : that is a partitioned mesh must be input (see Figure 2(c)). Instead of using a full-fledge Domain Decomposition (DD) approach which glues unknowns at the interface with an iterative algorithm, the strategy adopted in *Code\_Aster* is to assemble distributed matrices and vectors, as it would be done for a DD approach, but to solve the same system as in a sequential run, e.g. the global problem.

The first step is to decompose the finite element mesh into subdomains (see Figure 4). This is done using
dedicated mesh partitioners (and algorithms based on the METIS or SCOTCH libraries for example). These will create joints between subdomains that provide matching between subsets of the subdomains' interfaces (see Figure 5).

<figure><img src="image_4.png" /><figcaption>Fig. 4: Partition of a mesh into four subdomains which shall be distributed on four processes</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>Fig. 5: Matching of nodes on the interface between subdomains 0 and 1</figcaption></figure>

Apart from a mesh partitioner, two ingredients are needed to achieve the fixed goal: a global ordering for the unknowns must be devised out of the partitioned domains then a parallel linear solver shall be selected. This solver must handle distributed input to fully leverage the advantages of the decomposition. In order to build a global ordering and make sure communications are done in an efficient way, couplings between sub-domains may be determined using the joints. These couplings will enable point-to-point communications between processes as well as overlapping, saving time in the operation. An example of couplings for the partition presented above is given in Figure 6.

Finally the linear solver must be massively scalable that's why an iterative Krylov algorithm (Gener-
alized Conjugate Residual, GCR) was used in this work. This solver was provided as part of the PETSc
library [4] which is known to scale very well. When used in conjunction with an algebraic multi-grid
acting as a preconditioner [5], the iterations count in the Krylov solver tends to remain constant as prob-
lem size or number of cores increase, at least for well-conditioned cases. For these reasons, good scaling
can be expected.

The fact that an algebraic multigrid solver was used has implications on the type of studies allowed within this work. Indeed a constant number of degrees of freedom per node is required so mixing of models is forbidden. Moreover Lagrange multipliers cannot be used out of the box so applying complex boundary conditions may prove difficult. These limitations may be worked around for example by switching to a direct solver but scaling will be greatly impacted.
<figure><img src="image_6.png" /><figcaption>Fig. 6: Couplings between processes for the mesh partition in Figure 4</figcaption></figure>

## 4.2 Numerical examples

In order to illustrate the performance of this work, numerical examples will now be provided. Both academic and industrial applications are presented. The machine selected for these runs is a 1400-nodes cluster composed of Intel Xeon X5670 processors (2 x 2.93Ghz 6-cores).

<figure><img src="image_7.png" /><figcaption>Fig. 7: Cube in linear elasticity with prescribed displacement on opposite sides</figcaption></figure>

The example shown in Figure 7 is an academic analysis of cube in linear elasticity with displacements restrained on the bottom side and vertical displacement prescribed on the top side.

Using the previous work, a calculation on this cube totalling 500 millions degrees of freedom was successfully run on 700 cores in less than 2 minutes. Each process required a total amount of 4Gb of memory and treated around 230 000 elements. The Krylov iterations count was 24.

A smaller example of this cube (100 millions dofs) was run in different configurations. The speed-up and memory scaling for these runs are presented in Figure 8 and Figure 9. We can notice that speed-up is beginning to fall beyond 256 cores, two factors may explain this behaviour: elementary computations on
512 cores are becoming negligible so that timings may not be accurate, balancing is not perfect anymore at 512 cores. The memory measurements illustrate decent scaling up to 512 with a slight drop in the linear solver likely linked to balancing.

<figure><img src="image_8.png" /><figcaption>Fig. 8: Speed-up for a 100 millions dofs cube in linear elasticity</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>Fig. 9: Memory scaling for a 100 millions dofs cube in linear elasticity</figcaption></figure>

To finish with, an application on a pump equipment with pressure load is presented. Figure 10 shows the mesh partition of a coarse model of the pump and post-processing of the displacement solution.

Two runs were conducted, timings are given in Table 2. The important difference in performance compared to the academic example is mainly due to the "real life" nature of the problem. The complex geometry of the model has a deep impact on the distribution of the eigenvalues of the assembled matrix, to which iterative solvers are very sensitive. Furthermore, meshing this geometry is a difficult task for all meshing algorithms who often generate badly shaped tetrahedrons that deteriorate the conditioning of the stiffness matrix.
<figure><img src="image_10.png" /><figcaption>Fig. 10: Analysis of a pump equipment</figcaption></figure>

<table><caption>Table 2: Wall-clock time (s) for the pump equipment</caption><thead><tr><td>Case</td><td>64 cores</td><td>512 cores</td></tr></thead><tbody><tr><td>6 millions dofs</td><td>110</td><td>-</td></tr><tr><td>40 millions dofs</td><td>-</td><td>295</td></tr></tbody></table>

## 5 Conclusion

In this paper, we recalled the challenges being faced by mechanical FEA software users, namely the increasing size of models, and by software developers with computer architectures evolving towards more cores per node.

EDF's in-house developed FEA software was used to illustrate how these challenges are being treated in order to support current studies as well as prepare for the future.

Current capability of *Code_Aster* in non-linear mechanics was demonstrated, then a new paradigm for parallelism was presented and details of the implementation given.

Finally, results obtained with a fully distributed *Code_Aster* were discussed. A 40 million dofs analysis on a complex geometry was successfully run on 512 cores in linear elasticity.

The scope of usage of the fully distributed version is quite limited : only linear elasticity is supported and mixing of models is not allowed. Ongoing work is being carried to extend analyses to non-linear simulations.

## References

[1] *Code_Aster*. Official web site. [online] Available at http://www.code-aster.org/. [Accessed February 19th 2013].

[2] CSC. Elmer official website. [online] Available at http://www.csc.fi/english/pages/elmer/. [Accessed February 19th 2013].

[3] N. Tardieu. *A Newton-Krylov method for solid mechanics*, European Journal of Computational Mechanics, 21 (3-6), pp.374-384, 2012.

[4] S. Balay, J. Brown, K. Buschelman, V. Eijkhout, W. D. Gropp, D. Kaushik, M. G. Knepley, L. Curfman McInnes, B. F. Smith and H. Zhang. *PETSc Users manual*, ANL-95/11 - Revision 3.3, Argonne National Laboratory, 2012.
[5] M.W. Gee, C.M. Siefert, J.J. Hu, R.S. Tuminaro and M.G. Sala. *ML 5.0 Smoothed Aggregation User's Guide*, Sandia National Laboratories, 2006.