5. Construction et résolution des systèmes#

Une partie importante de l’algorithme non-linéaire consiste à résoudre des systèmes linéaires, que l’on construit en quatre temps:

  • Calcul et assemblage des chargements;

  • Calcul et assemblage des seconds membres;

  • Calcul et assemblage de la matrice;

  • Résolution du système linéaire;

Dans cette partie, nous allons décrire les différentes phases et les principales routines à considérer.

5.1. Systèmes à résoudre#

Actuellement, il existe plusieurs systèmes différents qui sont résolus dans l’opérateur op0070. Ils auront toujours la forme suivante (sauf pour le post-traitement, qui utilise les opérateurs de recherche de valeurs propres):

(5.1)#\[\begin{split}\left[\begin{array}{cc}K& {B}^{T}\\ B& 0\end{array}\right].\left\lbrace \begin{array}{c}\delta r\\ \delta s\end{array}\right\rbrace =\left\lbrace \begin{array}{c}{L}_{r}\\ {L}_{s}\end{array}\right\rbrace\end{split}\]

\(\left[\begin{array}{c}K\end{array}\right]\) est la matrice de rigidité ou une combinaison linéaire de matrices et \(\left[\begin{array}{c}B\end{array}\right]\) est la matrice des Lagrange de Dirichlet.

\(\left\lbrace \begin{array}{c}\delta r\end{array}\right\rbrace\) est l’incrément de solution des valeurs nodales inconnues du système et \(\left\lbrace \begin{array}{c}\delta s\end{array}\right\rbrace\) est l’incrément des paramètres de Lagrange associés à la dualisation des conditions limites. Concrètement, les «valeurs nodales inconnues» du systèmes peuvent être des déplacements, des vitesses, des accélérations, des pressions, des températures, etc. Les deux seconds membres sont également des valeurs nodales. Le détail des systèmes est donné dans les documentations de référence [R5.03.01] et [R5.05.05].

5.2. La routine merimo#

Pour le calcul des matrices tangentes (options FULL_MECA, FULL_MECA_ELAS, RIGI_MECA_TANG, RIGI_MECA_ELAS, RIGI_MECA, RIGI_MECA_IMPLEX) et les efforts intérieurs (option RAPH_MECA), on passe systématiquement par la routine merimo. Cette routine prépare les champs d’entrées (nombreux dans ce cas), elle est appelée à la fois pour le calcul des vecteurs élémentaires pour les efforts intérieurs (voir § 5.5.2 ) mais aussi pour le calcul des matrices élémentaires tangentes de rigidité (voir § 5.3 ).

5.3. Calcul desmatrices#

Nous considérons ici des matrices élémentaires (MEELEM) et des matrices assemblées (MEASSE). Tout comme pour le chargement (§ 5.5.1 ) , ce calcul utilise un système en deux temps:

  • Création d’une liste locale des matrices à calculer ou assembler;

  • Calcul ou assemblage des matrices;

Pour la phase de création de la liste des matrices à calculer ou assembler, il y a plusieurs routines (contrairement aux chargements qui n’utilisent que la routine nmchar). Ces routines sont les suivantes:

Opérations

Routine

Construction de la liste locale pour les matrices utilisées en dynamique explicite

ndxprm

Construction de la liste locale pour les matrices utilisées en correction

nmcoma

Construction de la liste locale pour les matrices utilisées en prédiction

nmprma

Construction de la liste locale pour les matrices utilisées en post-traitement (modes vibratoires ou modes de flambement)

nmflma

Construction de la liste locale pour les matrices constantes durant tout le calcul

nminmc

Pour la phase de calcul ou assemblage des matrices, on retrouve des routines similaires au cas des chargements.

Opérations

Routine

Ajout d’une matrice à assembler et/ou calculer dans la liste locale Initialisation de la liste locale

nmcmat nmcmat

Aiguillage calcul et/ou assemblage des matrices

nmxmat

Calcul des matrices

nmcalm

Assemblage des matrices

nmassm

Calcul et assemblage des matrices de rigidité

nmrigi

Il y a néanmoins une exception notable: le calcul des matrices tangentes non-linéaires utilise un appel direct par nmrigi car cela nécessite des paramètres supplémentaires (appel à merimo, voir 5.2 ) dont n’ont pas besoin les autres matrices. Le type de matrice MERIGI utilise donc nmrigi à la place de nmcalm et nmassm.

Si on désire modifier une matrice, il faut donc:

  • Ajouter son calcul ou son assemblage au bon endroit (ndxprm, nmcoma, nmprma, nmflma ou nminmc ) selon le cas, voire ajouter une nouvelle routine;

  • Ajouter le calcul élémentaire ou l’assemblage dans nmcalm et nmassm ;

Les routines nmxmatetnmcmatde manipulation de la liste locale ne doivent pas être modifiées .

5.4. Calcul de la matrice résultante MATASS#

Selon le type de matrice résultante que l’on désire obtenir (en prédiction, en correction ou pour l’accélération initiale), on utilise trois routines globales qui vont s’occuper de cette construction. Ce sont les routines nmprac , nmcoma et nmprma. Ces trois routines gèrent également l’affichage (option de calcul de la matrice de rigidité) et les statistiques (temps passé dans les opérations). Elles ont le même schéma de principe:

  • Calcul ou assemblage de matrices élémentaires (MEELEM) dans des matrices assemblées (MEASSE), (voir § 5.3 ) ;

  • Gestion des paramètres de réactualisation des matrices (rigidité, masse, amortissement), routine nmchrm;

  • Gestion du type de la matrice de rigidité (nom de l’option), routine nmchoi;

  • Construction de la matrice résultante par combinaison linéaire des autres matrices. En dynamique, les différentes contributions dans la matrice résultante (matrices de masse, d’amortissement et de rigidité) se combinent en utilisant un coefficient multiplicateur dépendant du schéma en temps utilisé et du pas de temps, on récupère ces coefficients via la routine d’accès NDYNRE (voir § 2.4.2.17 ). C’est la routine nmmatr;

  • Factorisation (ou pas) de la matrice résultante par appel à la routine preres;

Opérations

Routine

Calcul de la matrice assemblée résultante pour le calcul de l’accélération initiale

nmprac

Calcul de la matrice assemblée résultante pour la phase de prédiction

nmprma

Calcul de la matrice assemblée résultante pour la phase de correction

nmcoma

Gestion des paramètres de réactualisation des matrices (rigidité, masse , amortissement)

nmchrm

Gestion du type de la matrice de rigidité (nom de l’option)

nmchoi

Choix de re-création de la numérotation

nmrenu

Calcul de la matrice résultante MATASS

nmmatr

Prise en compte des chargements suiveurs avec leur fonction multiplicatrice

ascoma

Prise en compte de la modification de la matrice résultante par le contact/frottement (méthode DISCRETE)

nmasfr

5.5. Calcul du second membre#

Le calcul du second membre est la partie la plus différente d’un système linéaire à l’autre, il va contenir:

  • Les chargements donnés (voir le § 5.5.1 )

  • Les forces internes provenant de l’intégration du comportement (§ 5.5.2 );

  • Les quantités liées aux conditions limites de Dirichlet comme les réactions d’appui (§ 5.5.3 );

  • Les quantités liées aux variables de commande (§ 5.5.4 );

  • Les différentes quantités liées au contact/frottement (pas de détails dans ce document);

  • Les contributions d’inertie et d’amortissement pour la dynamique (§ 5.5.6 );

Dans le cas de la dynamique, pour les schémas multi-pas (Newmark complet avec MODI_EQUI=”OUI”) on va de plus ajouter les contributions des pas de temps précédents.

5.5.1. Calcul des chargements#

Il y a trois modes d’évaluation des chargements et deux phases. Les modes sont les suivants:

  • Mode <ACCI> des chargements pour l’évaluation de l’accélération initiale en dynamique;

  • Mode <FIXE> des chargements pour l’évaluation des chargements fixes ne dépendant pas de la solution;

  • Mode <VARI>des chargements pour l’évaluation des chargements variables dépendant de la solution (chargements suiveurs);

Il y a également deux phases. Soit on est en correction, soit on est en prédiction. Ces phases ne sont utiles que pour des cas très particuliers: amortissement modal, méthode IMPLEX et impédance modale.

La routine principale qui calcule les chargements est nmchar. Dans cette routine, suivant les fonctionnalités activées (list_func_acti) et suivant le mode/phase de calcul, on va provoquer le calcul des vecteurs élémentaires (variable-chapeau VEELEM) et leur assemblage en vecteurs assemblées (variable-chapeau VEASSE ). Pour cela, on utilise un certain nombre de routines utilitaires qui gèrent des listes locales à nmchar.

Opérations

Routine

Ajout d’un chargement à assembler ou calculer, avec une option éventuelle Initialisation de la liste des chargements (initialisations listes locales de nmchar)

nmcvec nmcvec

Aiguillage calcul ou assemblage des chargements

nmxvec

Calcul des chargements

nmcalv

Assemblage des chargements

nmassv

Certains chargements n’ont pas de phase de calcul élémentaire mais seulement la création d’un vecteur assemblé directement. Si on désire ajouter un chargement, il faut donc:

  • Définir quelle phase et quel mode va l’activer;

  • Ajouter son calcul ou son assemblage dans nmchar selon la fonctionnalité activée;

  • Ajouter le calcul élémentaire et/ou l’assemblage dans nmcalv et nmassv;

Les routines nmxvec,nmcvec ne doivent pas être modifiées .

5.5.2. Calcul des quantités liées aux efforts intérieurs#

Le calcul des efforts intérieurs peut être réalisé de deux manières:

  • Par intégration de la loi de comportement, c’est l’option de calcul RAPH_MECA;

  • Par ré-utilisation des contraintes calculées par ailleurs, c’est l’option de calcul FORC_NODA;

Cette distinction est importante car elle n’engendre pas les mêmes coûts (le second cas est beaucoup plus rapide). L’intégration du comportement est une opération coûteuse qui implique de modifier variables internes et contraintes, et, souvent, de résoudre localement (à chaque point de Gauss), un système non-linéaire. Comme le calcul des contraintes est également nécessaire lors de l’évaluation des matrices tangentes, le comportement est aussi intégré lors du calcul des matrices tangentes (option FULL_MECA). En statique, la phase de prédiction calcule les efforts intérieurs par l’option FORC_NODA. En dynamique, on intègre systématiquement le comportement.

Opérations

Routine

Calcul des vecteurs élémentaires pour les efforts intérieurs (intégration du comportement)

nmfint

Assemblage des vecteurs élémentaires pour les efforts intérieurs (intégration du comportement)

nmaint

Le calcul de l’option FORC_NODA est fait par le mécanisme d’évaluation des chargements (§ 5.5.1 ).

5.5.3. Calcul des quantités liées aux conditions limites dualisées#

La dualisation des conditions limites de Dirichlet (voir [R3.01.01]) entraîne le calcul de deux quantités de type vecteur assemblée:

  • Les réactions d’appui, par le produit de la matrice des conditions limites dualisées avec les Lagrange correspondant aux degrés de liberté \(\left[B\right].\left\lbrace \lambda \right\rbrace\) , il s’agit du vecteur identifié par CNDIRI;

  • La valeur des degrés de liberté imposée par dualisation \(\left[B\right].\left\lbrace u\right\rbrace\) . A convergence, cette quantité sera égale aux conditions limites données \(\left\lbrace {u}^{d}\right\rbrace\) , le vecteur est identifié par CNBUDI;

Opérations

Routine

Calcul des vecteurs élémentaires pour les réactions d’appui

nmdiri

Assemblage des vecteurs élémentaires pour les réactions d’appui

nmadir

Calcul et assemblage des quantités imposées par dualisation

nmbudi

5.5.4. Calcul des quantités liées aux variables de commande#

Le calcul de ces quantitésse fait par le mécanisme d’évaluation des chargements (§ 5.5.1 ), type CNVCF0, CNVCF1 et CNVCPR.

5.5.5. Calcul des quantités constantes pendant le calcul#

Un certain nombre de vecteurs sont constants durant le calcul, ils utilisent alors un mécanisme similaire aux chargements (§ 5.5.1 ), à l’aide des routines nmxvec, nmcvec, nmcalv et nmassv, sauf qu’on passe par la routine nminvc au lieu de nmchar.

Opérations

Routine

Calcul et assemblage des vecteurs constants durant le calcul

nminvc

5.5.6. Calcul des quantités liées à la dynamique – Forces d’inertie et d’amortissement#

Le calcul de ces quantitésse fait par le mécanisme d’évaluation des chargements (§ 5.5.1 ), type CNDYNA. Ces quantités n’existent pas à l’état élémentaire puisqu’ils sont le résultat d’un produit matrice/vecteur.

Opérations

Routine

Calcul des quantités dynamiques pour le second membre

ndfdyn

Calcul des forces d’inertie

nminer

Calcul des forces d’amortissement

nmhyst

Selon les schémas en temps, il est nécessaire d’utiliser des coefficients multiplicateurs spécifiques devant chacune de ces quantités, on récupère ces coefficients via la routine d’accès NDYNRE (voir § 2.4.2.17 ).

5.5.7. Seconds membres résultants#

Il en reste plus qu’à construire les différents seconds membres par combinaison linéaire de toutes les quantités détaillées précédemment. Les routines principales de construction de ces seconds membres sont au nombre de sept:

Opérations

Routine

Calcul du second membre pour la correction

nmassc

Calcul du second membre pour la prédiction – Option DEPL_CALCULE ou EXTRAPOLE

nmassd

Calcul du second membre pour l’accélération initiale

nmassi

Calcul du second membre pour la prédiction

nmassp

Calcul du second membre pour la prédiction – Cas statique

nsassp

Calcul du second membre pour la prédiction – Cas dynamique implicite

ndassp

Calcul du second membre pour la prédiction – Cas dynamique explicite

nmassx

Toutes ces routines reposent sur les mêmes principes:

  • Récupération des quantités nécessaires, déjà assemblées ou calculées localement dans ces routines, via la variable-chapeau VEASSE (voir les § 5.5.1 à § 5.5.6 );

  • Récupération des différents coefficients multiplicateurs. En dynamique, selon les schémas en temps, il est nécessaire d’utiliser des coefficients multiplicateurs spécifiques devant chacune de ces quantités, on récupère ces coefficients via la routine d’accès NDYNRE (voir § 2.4.2.17 );

  • Construction d’un ou de deux seconds membres par combinaison linéaire (utilisation des routines standards vtaxpy). Il y a deux seconds membres dans le cas du pilotage;

Pour améliorer la lisibilité, on utilise des routines utilitaires standardisées pour regrouper les cas les plus fréquents.

Opérations

Routine

Calcul des composantes du vecteur second membre pour les chargements de type Dirichlet, fixes au cours du pas de temps, pour le cas normal et le cas piloté

nmasdi

Calcul des composantes du vecteur second membre pour les chargements de type Neumann, fixes au cours du pas de temps, pour le cas normal et le cas piloté

nmasfi

Calcul des composantes du vecteur second membre pour les chargements de type Neumann, variables (suiveurs) au cours du pas de temps, pour le cas normal (pas de cas piloté)

nmasva

Calcul des composantes du vecteur second membre pour les chargements spécifiques de la dynamique, variables (suiveurs) au cours du pas de temps, pour le cas normal (pas de cas piloté)

ndasva

Calcul des composantes du vecteur second membre pour les chargements spécifiques de la dynamique, calcul de l’accélération initiale, variables (suiveurs) au cours du pas de temps, pour le cas normal (pas de cas piloté)

nmacva

5.6. Résolution du système#

La résolution se fait via l’intermédiaire de la routine nmresd, qui prend en entrée la matrice résultante assemblée (voir § 5.4 ) et les seconds membres idoines (deux seconds membres si on a du pilotage).

Opérations

Routine

Résolution du système

nmresd

Résolution du système sur les degrés de liberté physiques (appel à resoud)

nmreso

Résolution du système sur les degrés de liberté généralisés

nmresg

La routine de résolution va donc retourner un ou deux champs solutions DEPSO1 et DEPSO2, stockés dans la variable-chapeau SOLALG.