2. Les différents solveurs linéaires disponibles#

Ces solveurs linéaires sont en fait omniprésents dans le déroulement des opérateurs de Code_Aster car ils sont souvent enfouis au plus profond d’autres algorithmes numériques: schéma non linéaire, intégration en temps, analyse modale etc. Ils en consomment souvent la majeure partie du temps CPU et de la mémoire. Le choix et le paramétrage complet du solveur linéaire requis s’effectue via le mot-clé facteur SOLVEUR. Il est présent dans la plupart des commandes de calcul (STAT_NON_LINE, THER_LINEAIRE, CALC_MODES…).

Ce mot-clé permet de choisir entre les deux classes de solveurs: les directs et les itératifs. Concernant les directs , on dispose du classique algorithme de «Gauss» (SOLVEUR/METHODE=”LDLT”), d’une factorisation multifrontale (“MULT_FRONT”) et d’une externe (“MUMPS”). Pour les itératifs , il est possible de faire appel à un gradient conjugué (“GCPC”) ou à certains outils de la librairie publique PETSc (“PETSC”).

Seuls MULT_FRONT,MUMPS, PETSC sont parallélisés . Le premier en OpenMP, les deux autres en MPI et en OpenMP (parallélisme hybride). Mais tous les solveurs sont compatibles avec un traitement parallèle ( via MPI) des calculs élémentaires et des assemblages. Et ce, que ces traitements soient initiés juste avant l’utilisation du solveur linéaire proprement dit, ou, dans un autre opérateur (par exemple de pré/post-traitements).

Détaillons un peu le fonctionnement de chacun d’entre eux:

Solveurs directs

/”MUMPS”

Solveur direct de type multifrontale avec pivotage. Ce solveur est obtenu en appelant le produit externe MUMPS développé par CERFACS/IRIT/INRIA/CNRS (cf. Copyright §4). Le stockage matriciel hors solveur est MORSE. En entrée du solveur, on fait la conversion au format interne de MUMPS: \(i\) , \(j\) , \({K}_{ij}\) , centralisée ou distribuée. Pour Code_Aster , son intérêt principal réside dans sa capacité à pivoter lignes et/ou colonnes de la matrice lors de la factorisation en cas de pivot petit. Cette possibilité est utile (voire indispensable) pour les modèles conduisant à des matrices non définies positives (hors conditions aux limites); Par exemple, les éléments “mixtes” ayant des degrés de liberté de type “Lagrange” (éléments incompressibles…). Cette méthode est parallélisée en mémoire distribuée (MPI) et en mémoire partagée (OpenMP). Elle peut être exécutée sur plusieurs processeurs (via l’interface Astk menu Options/Optionsdelancement/ncpus&mpi_nbcpu&mpi_nbnoeud). En parallèle MPI (mpi_nbcpu&mpi_nbnoeud) , MUMPS distribue naturellement ses données (matrice, factorisée…) entre les cœurs des différents nœuds de calcul alloués. Ce qui accélère grandement les calculs et permet de réduire l’occupation mémoire nécessaire, par processus MPI, pour lancer le calcul. Cette consommation RAM peut encore plus être réduite via les paramètres GESTION_MEMOIRE ou RENUMde MUMPS. En terme de mémoire, le goulet d’étranglement peut alors se trouver au niveau de l’espace JEVEUX. Pour réduire ce dernier, on peut distribuer la matrice Aster (MATR_ASSE) via l’option MATR_DISTRIBUEE. MUMPS s’appuie aussi sur une deuxième niveau de parallélisme imbriqué au sein du parallélisme MPI (parallélisme hybride ) et basé sur OpenMP. Ce deuxième niveau de parallélisme est principalement activé dans les appels aux librairies mathématiques sous-jacentes (BLAS, LAPACK). Contrairement au MPI, c’est un parallélisme à mémoire partagée et donc limité aux seuls cœurs d’un nœuds de calcul (ncpus). Pour accélérer les grosses études (N au moins > 2.106ddls) on peut en plus activer des options d’accélération/compression de MUMPS (à partir de la v5.1.0) via les mots-clés ACCELERATION/LOW_RANK_SEUIL.

/”MULT_FRONT”

Solveur direct de type multifrontale développé en interne EDF R&D . Le stockage matriciel est MORSE (ou “CSC” pour “Compressed Sparse Column”) et donc proscrit tout pivotage. Cette méthode est parallélisée en mémoire partagée (OpenMP) et peut être exécutée sur plusieurs processeurs (via l’interface Astk menu Options/Options de lancement/ncpus). La matrice initiale est stockée dans un seul objet JEVEUXet sa factorisée est répartie sur plusieurs, donc elle peut être déchargée partiellement et automatiquement sur disque.

/”LDLT”

Solveur direct avec factorisation de Crout par blocs (sans pivotage) développé en interne EDF R&D . Le stockage matriciel hors solveur est MORSE. En entrée du solveur, on fait la conversion au format interne de LDLT: ‘ligne de ciel’ (‘SKYLINE’). On a une pagination de la mémoire complètement paramétrable (la matrice est décomposée en blocs gérés en mémoire de façon indépendante et déchargés sur disque au fur et à mesure) qui permet de passer de gros cas mais qui se paye par des accès disques coûteux. En outre, ce solveur permet de ne factoriser que partiellement la matrice. Cette possibilité est “historique”. Elle permet de factoriser la matrice en plusieurs “fois” (plusieurs travaux) voire de modifier à la volée les dernières lignes de cette factorisée. Aujourd’hui, on n’imagine pas bien l’intérêt de cette fonctionnalité hormis pour certaines méthodes (dite discrètes) de contact-frottement où l’on a, à dessein, placé dans les dernières lignes de la matrice les termes concernant les noeuds susceptibles d’être en contact. Ainsi, au fur et à mesure des itérations d’appariement, les relations entre ces noeuds changeant, on efface puis recalcule que ces dernières contributions de la factorisée. C’est un exemple typique où l’usage astucieux d’un algorithme assez frustre peu amener des gains majeurs (en temps).

Solveurs itératifs

/”GCPC”

Solveur itératif de type gradient conjugué avec préconditionnement \(\mathit{ILU}(k)\) ou basé sur une factorisée simple précision ( via MUMPS). Le stockage de la matrice est alors MORSE. Avec le Cholesky incomplet, la matrice initiale et sa factorisée incomplète sont stockées, chacune, dans un seul objet JEVEUX. Avec la factorisée simple précision, le préconditionneur est beaucoup plus coûteux (en CPU/RAM), mais il est effectué en simple précision et son calcul peut être mutualisé pendant plusieurs résolutions (problème de type multiples seconds membres, par ex. STAT_NON_LINE).

/”PETSC”

Solveurs itératifs issus de la librairie externe PETSc (Laboratoire Argonne). Le stockage matriciel hors solveur est MORSE. En entrée du solveur, on fait la conversion au format interne de PETSc: ‘CSR” pour “Compressed Sparse Row”. PETSc alloue des blocs de lignes contigus pas processeur. Cette méthode est parallélisée en mémoire distribuée (MPI) et peut être exécutée sur plusieurs processeurs (via l’interface Astk menuOptions/Optionsdelancement/mpi_nbcpu&mpi_nbnoeud). Lorsque PETSc utilise un préconditionneur basé sur MUMPS (PRECOND=”LDLT_SP”), il bénéficie aussi du parallélisme potentiellement hybride du solveur direct. Mais les meilleures accélérations de PETSc se produisent plutôt en privilégiant le premier niveau de parallélisme, MPI. Attention: les solveurs PETSCet MUMPSétant incompatibles en séquentiel, on privilégie généralement MUMPS. Pour utiliser PETSC, il faut donc souvent lancer une version parallèle de Code_Aster (quitte à ne solliciter qu’un seul processeur).

Périmètre Solveur

Robustesse

CPU

Mémoire

Détails

Direct

MULT_FRONT

Solveur universel. A déconseiller pour les modélisations nécessitant du pivotage (EF mixtes de X-FEM, incompressible…).

+++

Séq: ++ //: ++ (speed-up~2)

OOC[1]_

Sur 4/8 coeurs.

MUMPS

Solveur universel, solveur de référence

+++

Séq: ++ //: +++ (sp~30)

– IC + OOC

Jusqu’à 512 coeurs

LDLT

Solveur universel (mais très lent sur de gros cas). A déconseiller pour les modélisations nécessitant du pivotage (EF mixtes de X-FEM, incompressible…). Factorisation partielle possible (contact discret).

+++

Séq: +

OOC

Plutôt les petits cas ou les cas de taille moyenne pouvant tirer partie de factorisation partielle.

Itératif

GCPC

Problèmes symétriques réels sauf ceux requérant obligatoirement une détection de singularité (calcul modal, flambement).

  • (LDLT_INC)

  • (LDLT_SP)

Séq: +++

+++ (LDLT_INC) ++ (LDLT_SP)

Avec LDLT_INCne pas trop augmenter le niveau de préconditionnement. Parfois très efficace (thermique…) surtout en non linéaire avec LDLT_SP.

PETSC

Idem GCPCmais compatible avec le non symétrique.

  • (LDLT_INC)

  • (LDLT_SP)

Séq: +++ //: +++ (sp~4)

++ IC

Algorithmes plutôt robustes: GMRES. Souvent très efficace en non linéaire avec LDLT_SP.

Figure 2.-1: Synoptique des solveurs linéaires disponibles dans Code_Aster.

Remarque:

  • Pour être complètement exhaustif, on peut préciser que quelques (rares) opérations numériques sont opérées avec un paramétrage solveur «figé». Celui-ci est donc inaccessible aux utilisateurs (sans surcharge logiciel). Parmi celles-ci, on trouve certaines résolutions des méthodes discrètes de contact-frottement (réservées à LDLT), les recherche de modes rigides en modal, les calculs de modes d’interfaces (réservées à MUMPS ou à LDLT)… A chaque fois des raisons fonctionnelles ou de performance expliquent ce choix peu transparent.