Differences

This shows you the differences between two versions of the page.

--- commit:futur:gaetan [2016/04/23 16:16] – wautelet
+++ commit:futur:gaetan [2019/07/14 11:11] (current) – [Boundary Volume Hierarchy] wautelet
@@ Line 1: / Line 1: @@
-====== Commit 2016-04-19 ======
+====== Commit 2018-04-23 ======
 ===== Modifications principales =====
-Ce commit porte principalement sur l'amélioration du Lagrangien Augmenté.
+Ce commit porte principalement sur l'amélioration du Line Search et de la mise en rotation quasi-statique avec Lagrangien Augmenté.
-==== Algorithme du Lagrangien Augmenté ====
-Je vous met à profit les derniers développements sur l'algorithme du Lagrangien Augmenté.
+==== Algorithme du Line Search Structurel  ====
+Je vous mets à profit les derniers développements sur l'algorithme du Line Search Structurel.
 === Algorithme Général ===
-J'ai constaté à ma grande surprise qu'il faut forcer le calcul d'une itération mécanique après une augmentation.
-Après chaque augmentation (modification des lagrangiens augmentés à configuration cinématique fixée), on résout l'équilibre mécanique par la méthode de Newton Raphson. Cependant, il se peut que la modification du vecteur des forces externes ne soit pas suffisamment conséquente pour changer la norme moyenne du résidu  (alors que la norme max du résidu est nettement plus sensible à ces changements là ...) et dès lors forcer une nouvelle itération mécanique, c'est à dire une correction des positions et obtenir une nouvelle configuration cinématique. Donc, l'algorithme continue de modifier les lagrangiens augmentés à configuration cinématique donnée jusqu'à espérer avoir une variation suffisante des forces externes pour avoir une itération mécanique. Si la configuration donnée ne change pas, les gaps normal et tangentiel ne peuvent pas évoluer et nous ne pouvons pas converger vers les tolérances imposées à l'algorithme du lagrangien augmenté. En résumé, l'algorithme stagne et l'intégration temporelle risque de s'arrêter tout simplement !
-Il est vrai que les articles de référence en la matière travaillent en général en terme de norme énergétique (ce qui nécessite d'office le calcul d'une itération mécanique !).
-Dans le jeux de données Python, on a
 <code>
@@ Line 24: / Line 19: @@
 </code>
-Par défaut, ce paramètre est mis à False !
-=== Mise à jour des lagrangiens augmentés ===
+==== Algorithme du Lagrangien Augmenté dans le contexte de la mise en charge avec vitesse de rotation initiale ====
-J'ai développé une autre approche pour mettre à jour les lagrangiens augmentés sur base de l'histoire de l'évolution de la pression de contact en fonction du gap normal et le cisaillement de contact en fonction du gap tangentiel. Ces lagrangiens augmentés "extrapolés" sont mis à jour si besoin via le critère de frottement de Coulomb ou de Tresca. Au lieu de mettre à jour les coefficients de pénalité pour chaque nœud de contact à chaque augmentation ou au cours de chaque itération mécanique, on utilise  cette autre approche développée par Wriggers and Zavarise en contact sans frottement qui possède un ordre de convergence super-linéaire, mais qui n'est plus spécialement monotone en convergence sans garde fou !
+==== Divers ====
-Avec ces derniers gardes fous, l'algorithme semble être encore plus efficace que l'approche standard car on tient de l'historique au cours des augmentations !
-<code>
-augLagAugmentation = AugLagExtrapolationAugmentation(alm)
-augLagAugmentation.setUseStabilityCriterion(True)
-augLagAugmentation.setUseExtrapolationCriterion(True)
-augLagAugmentation.setExtrapolationFactor(_factor)
-augLagAugmentation.setResetPreviousAugmentedLagrangians(True)
-</code>
-La dernière option permet d'activer la procédure dès la première augmentation en indiquant False, sinon il faut au moins une augmentation selon la procédure d'Uzawa avant de lancer celle par extrapolation.
+=== VectorToScalarOperator ===
-Par défaut, tous les critères sont activés et le facteur d'extrapolation est 0.1 (c'est à dire on autorise une variation de maximum 10 % de la valeur courante de la pression de contact et du cisaillement de contact.).
+=== MiscValueExtractor ===
-Cette procédure est activée si on a le nœud de contact est en contact au cours de deux résolutions successives de l'équilibre mécanique (Le lagrangien normal est mis à jour en fonction de la valeur à l'augmentation k du couple (contact pressure, normal gap) et l'augmentation k-1 du même couple), sinon les lagrangiens augmentés sont mis à jour selon la méthode d'Uzawa.
-=== Critère d'arrêt géométrique ===
+==== Boundary Volume Hierarchy ====
-J'ai ajouté une nouvelle méthode de calcul de longueur de contact utilisée lors des critères géométriques normalisés.
-Dans les estimateurs de pas de temps critique dans les schémas explicites, la notion de dimension équivalente d'un élément apparaît naturellement dans la définition du pas de temps critique. Cependant, ce choix reste vague (distance minimale entre deux nœuds de la maille ou autre) et dès lors, il propose une analyse des valeurs propres à partir de laquelle il dérive une longueur équivalente (plus conservative !). Leur raisonnement est indépendant à priori du type de maille et de la dimension de l'espace.  C'est pour ça, que j'ai repris leur développement et je l'ai codé dans Metafor pour le cas 2D, 3D et axisymmétrique.
-En résumé, on a les méthodes suivantes :
-  * _method = ALM_AIC_LENGTH (Par défaut) : on calcule une longueur équivalente basée sur la surface de contact (en 2D, c'est la moyenne des demi-longueurs des courbes adjacentes au nœud et en 3D, c'est la moyenne de la racine carrée  des quart aires des faces adjacentes au noeud.). On normalise par la même valeur pour les gaps normal et tangentiel.
+=== Cas Tests Modifiés ===
-  * _method = ALM_CRITICAL_LENGTH : on calcule une longueur équivalente basée sur l'estimateur du pas de temps critique en explicit des éléments volumiques adjacents au nœud de contact.  On normalise par la même valeur pour les gaps normal et tangentiel.
-  * _method = ALM_GEO_LENGTH : Pour le gap normal, on calcule une distance normale des mailles volumiques par rapport au plan tangent de contact ce qui donne une épaisseur des éléments volumiques sous-jacents. On divise cette distance par deux pour obtenir la longueur caractéristique. Pour le gap tangentiel, on procède de la même manière que dans la première méthode mais on travaille dans le plan tangent de contact ! Dès lors, ce critère est par définition anisotrope.
-Par exemple, pour un élément quadrangulaire rectangulaire $a \times b $(avec$ b > a > 0$), on a
+En activant l'utilisation d'une recherche sous forme d'arbre pour la détection du contact globale dans le cas test du pied milieu que Luc a ajouté dans la batterie, j'ai remarqué un temps CPU excessif lors de la construction de l'arbre. Après investigation, j'ai trouvé un bug suite à un objet passé par copie au lieu de par référence. Depuis lors vu la quantité de triangle sur les corps rigides maîtres, j'ai paramétrisé le cas test pour activer l'utilisation des arbres lors de la détection du contact globale.
-  * _method = ALM_AIC_LENGTH :  $Lc_N = Lc_T = \frac{b}{2}$
-  * _method = ALM_GEO_LENGTH :  $Lc_N = \frac{a}{2}$  et  $Lc_T = \frac{b}{2}$
-  * _method = ALM_CRITICAL_LENGTH :  $Lc_N = Lc_T =  \frac{a}{\sqrt{1+\frac{b^2}{a^2}}}$  (= Moitié de la diagonale si on a un carré  $b = a$ .)
-Dans le jeux de données Python, on choisit les options comme ceci :
-<code>
+Nastran test
-augLagGeoCriterion = AugLagNormalisedGeoCriterion(alm)
-augLagGeoCriterion.setCharacteristicLengthMeth(_method)
-</code>
-Sur un pas de temps donné, la maille peut subir de grande déformation d'où le fait que l'on recalcule à chaque fois la longueur caractéristique au cours des augmentations.
-=== Mise à jour du pas de temps ===
+J'ai activé l'utilisation du boundary volume hierarchy dans ce cas test pour pouvoir
-L'algorithme du Lagrangien augmenté permet à priori de résoudre un problème de contact sur un pas de temps plus élevé que l'équivalent en Méthode de Pénalité. Pour mettre à jour le pas de temps, j'ai proposé plusieurs méthodes :
-  * _equivalentIteMethod = ALM_TSCM_FIRST : Nombre d'itérations mécaniques avant la première augmentation (DEFAULT)
+apps.toolbox.createContactTests.py
-  * _equivalentIteMethod = ALM_TSCM_MAX : Nombre d'itérations mécaniques maximales au cours des augmentations
-  * _equivalentIteMethod = ALM_TSCM_MEAN : Nombre d'itérations mécaniques moyennes au cours des augmentations
-  * _equivalentIteMethod = ALM_TSCM_WEIGHTEDMEAN : Nombre d'itérations mécaniques moyennes pondérées au cours des augmentations (Les poids sont le rapport entre le nombre d'itération mécanique divisée par le nombre d'itération mécanique maximale.)
-Cette valeur d'itération "équivalente" est insérée dans l'algorithme de mise à jour du pas de temps proposé par JPP dans sa thèse.
-Dans le jeux de données Python, on a
+La fonction defineMim() passe de cinq arguments d'entrée à trois arguments d'entrée pour configurer le Newton Raphson mécanique. L'algorithme du line search est activé avec une fonction defineMls() qui prend trois arguments d'entrée.
-<code>
-tscm = NbOfAugLagMechNRIterationsTimeStepComputationMethod(metafor)
-tscm.setTimeStepDivisionFactor(_timeStepDivisionFactor)
-tscm.setEquivalentIteMethod(_equivalentIteMethod)
-tscm.setMaxNbOfIteInAdaptOfTimeStep(_maxNbOfIteInAdaptOfTimeStep)
-tscm.setNbOptiIte(_nbOptiIte)
-tsm.setTimeStepComputationMethod(tscm)
-</code>
-Par exemple, si on a le profil suivant de convergence sur un pas de temps de donné (12 augmentations pour 30 itérations mécaniques) :
+SurfaceSelector() Nouveau selecteur...
-^ ALM ITE NO     ^ MECH ITE NUMBER  ^
-| 0     | 2     |
-| 1     | 3     |
-| 2     | 4     |
-| 3     | 3     |
-| 4     | 3     |
-| 5     | 3     |
-| 6     | 2     |
-| 7     | 2     |
-| 8     | 2     |
-| 9     | 2     |
-| 10     | 2     |
-| 11     | 1     |
-| 12     | 1     |
-on obtient l'itération mécanique équivalente suivante :
-^ ALMTSCM     ^ EQUIVALENT MECH ITE NUMBER  ^
-| ALM_TSCM_FIRST    | 2     |
-| ALM_TSCM_MAX    | 4     |
-| ALM_TSCM_MEAN    | 2     |
-| ALM_TSCM_WEIGHTEDMEAN    | 3     |
-=== Benchmark - Cont2 ===
-==== Informations sur Interaction de Contact ====
+=== Vector To Scalar Operator ===
-Il est possible de demander d'afficher des informations sur le status des nœuds de contact d'une interaction.
+J'ai ajouté plusieurs nouveaux opérateurs de po
-<code>
-ci = RdContactInteraction(1)
-ci.setShowInformation(True)
-</code>
-Voici un exemple de sortie.
-<code>
-mechanical iteration 0  : Rmean = 1.000000e+000  Rmax = 2.244123e+002
-RdContactInteraction #2
-Actual/Total Number of Nodes in Contact = 4 / 11
-Number of Nodes in Sticking/Sliding Contact = 4 / 0
-mechanical iteration 1  : Rmean = 2.192696e-001  Rmax = 4.097902e+001
-RdContactInteraction #2
-Actual/Total Number of Nodes in Contact = 4 / 11
-Number of Nodes in Sticking/Sliding Contact = 4 / 0
-mechanical iteration 2  : Rmean = 2.016662e-001  Rmax = 2.821996e+001
-RdContactInteraction #2
-Actual/Total Number of Nodes in Contact = 3 / 11
-Number of Nodes in Sticking/Sliding Contact = 3 / 0
-mechanical iteration 3  : Rmean = 1.084675e-001  Rmax = 9.532517e+000
-RdContactInteraction #2
-Actual/Total Number of Nodes in Contact = 3 / 11
-Number of Nodes in Sticking/Sliding Contact = 3 / 0
-mechanical iteration 4  : Rmean = 2.831632e-002  Rmax = 2.978541e+000
-RdContactInteraction #2
-Actual/Total Number of Nodes in Contact = 3 / 11
-Number of Nodes in Sticking/Sliding Contact = 3 / 0
-mechanical iteration 5  : Rmean = 1.060746e-003  Rmax = 1.218047e-001
-RdContactInteraction #2
-Actual/Total Number of Nodes in Contact = 3 / 11
-Number of Nodes in Sticking/Sliding Contact = 3 / 0
-mechanical iteration 6  : Rmean = 1.744038e-005  Rmax = 1.792596e-003
-RdContactInteraction #2
-Actual/Total Number of Nodes in Contact = 3 / 11
-Number of Nodes in Sticking/Sliding Contact = 3 / 0
-</code>
-Cela peut permettre de comprendre la raison pour laquelle un pas de temps ne converge pas (Status collant glissant qui change ou status en contact ou hors contact qui change).
-==== Nombre de conditionnement réciproque ====
-Il est désormais possible d'estimer le nombre de conditionnement réciproque d'une matrice carrée (Matrix de Metafor) quelconque via les routines LAPACK. Les matrices utilisées pour la résolution du système linéaire sous-jacente au calcul de l'équilibre implicite se trouvent sous un forme creux liée au solveur utilisé. Cependant, il n'existe pas (sauf erreur de ma part) d'estimateur du nombre de conditionnement réciproque dans les routines LAPACK, donc nous devons transformer une matrice creuse sous un format plein pour pouvoir effectuer le calcul proprement dit (ce n'est pas top du tout mais bon faute de mieux ...). J'ai réuni un peu de documentation sur ces estimateurs, si qqn veut avoir les infos.
-On peut demander de calculer à chaque résolution du système linéaire, le nombre de conditionnement réciproque. Pour se faire, il suffit de faire :
-<code>
-solvermanager = metafor.getSolverManager()
-solver = DSSolver()
-solver.setShowRCond(True)
-solvermanager.setSolver(solver)
-</code>
-Également, il est possible de faire un extracteur de cette grandeur (Nombre de conditionnement réciproque de la matrice FreeJacobian du schéma d'intégration temporelle) via les commandes suivantes :
-<code>
-valuesmanager.add(1, MiscValueExtractor(metafor, EXT_RCOND), 'RCond')
-</code>
-Dans Matlab, on compare ce nombre de conditionnement réciproque avec l'epsilon machine multipliée par la taille du système à résoudre. Pour rappel, ce nombre est proche de 1 lorsque le système est bien conditionné et ce nombre est proche de 0 lorsque le système est mal conditionné.
-Finalement, cet extracteur peut aider à quantifier l'influence du coefficient de pénalité sur la résolution du système linéaire.
-==== Nouveaux extracteurs de contact ====
-Par la même occasion, j'ai ajouté plein de nouveau extracteurs de variables de contact. J'ai ajouté directement les infos dans la partie sur les extracteurs de contact. Finalement, la documentation est à jour à ce niveau là.
+MinOfNonZeroOperator()
+MaxOfNonZeroOperator()
+AbsMinOfNonZeroOperator()
+AbsMaxOfNonZeroOperator()
+MinAbsOfNonZeroOperator()
+MaxAbsOfNonZeroOperator()
@@ Line 199: / Line 77: @@
 <code>
 [r]:
 [a]:
 </code>
- --- //[[gwautelet@ulg.ac.be|Gaëtan WAUTELET]] 2016/04/22//
+ --- //[[gwautelet@ulg.ac.be|Gaëtan WAUTELET]] 2016/04/23//