Anomalie #330
bug contact
Description
Gérard,
comme dit à l'oral, je suis en train de regarder la partie contact pour contribuer à améliorer la gestion du contact. Je ne sais pas s'il faudra faire un ticket par problème que je vais rencontrer.
Sachant que mon objectif principal est de repérer les cas où Herezh++ perd le contact de manière anormale.
Néanmoins, voici un premier bug que j'ai trouvé qui n'a pas forcément de lien avec la pénétration. Je suis tombé sur un cas où j'obtiens le message d'erreur :
Erreur : la norme du vecteur est trop petite !
norme = nan
Droite::Droite (Coordonnee& B,Coordonnee& vec)
voir le Readme fourni pour le détail.
En gros, c'est un calcul de joint torique axisymétrique écrasé selon Y entre une chemise et un fouloir. Le maillage esclave (joint) est LINEAIRE. Les maillages maitre (chemise et fouloir) sont QUADRACOMPL.
A noter que ce bug n'intervient pas si on prend des maillages LINEAIRE pour chemise et fouloir (repertoire version_LINEAIRE/). C'est très possible que le problème ne soit pas lié à l'interpolation. J'imagine que c'est tombé sur un cas particulier.
je pense qu'à l'avenir, je vais faire une moulinette pour générer aléatoirement des mises en données et essayer de tomber sur les bugs en faisant varier quelques paramètres. D'où ma question initiale : est-ce qu'il faudra faire un ticket pour chaque bug ? J'ai mis un titre volontairement général à ce ticket au cas où.
a+
Julien
Fichiers
Mis à jour par Gérard Rio il y a presque 2 ans
- Statut changé de Nouveau à En cours
- % réalisé changé de 0 à 30
En utilisant la version non fast, on obtient une erreur au bout de qq incréments, c'est dû au fait que les noeuds de l'axe de rotation doivent être bloqués en x: donc j'ai bloqué les noeuds et là : plus d'erreur.
D'où un pb pour moi:
- version fast arrêt à l'incrément 24 du à un nan
- version non fast : pas d'erreur, le calcul va à terme !!
d'où galère !! car pas possible d'utiliser le debug.
Après mal pas de manip longues et laborieuses, j'ai localisé la cause et le lieu où le pb se construisait, maintenant il faut que je trouve un moyen d'éviter la dérive:
a priori dans une recherche itérative d'un point d'intersection, je trouve un point à l'infini ...
bon... affaire à suivre
Mis à jour par Julien Troufflard il y a presque 2 ans
- Fichier test_1_modif.tar test_1_modif.tar ajouté
"c'est dû au fait que les noeuds de l'axe de rotation doivent être bloqués en x"
=> ouhla! ah oui effectivement c'est un oubli de ma part. Mais j'ai trouvé pire dans ma mise en données. Malgré des années d'utilisation de stamm, je continue à me faire avoir. Quand on génère un maillage rectangulaire et que l'on indique un nombre d'éléments selon X moins grand que selon Y, stamm fait automatiquement une rotation des listes de référence. Par exemple, on se retrouve avec N_N sur le côté "Ouest" du maillage. Enfin, il me semble qu'il fait ça, systématiquement? je sais pas, mais en l'occurrence c'est le cas pour les maillages fouloir et chemise.
Donc les conditions limites n'ont aucune logique dans le calcul que je t'ai fourni. Tant mieux si ça tombe sur un cas à débugger, c'est le but de ce ticket.
Néanmoins, je t'ai mis ci-joint une nouvelle version de ce calcul (version avec joint linéaire et outils quadracompl + version tout linéaire, avec les fichiers .log pour les 2 calculs). Le calcul est cette fois réalisé avec la version 7.005 (fast).
Et avec cette nouvelle version, je tombe sur un des cas particuliers que je recherchais. à savoir : un calcul qui diverge quand l'un des noeuds du maillage esclave se retrouve en vis-à-vis du maillage maitre. Autrement dit, mon interprétation est que Herezh semble avoir du mal à gérer certains cas où il faut changer d'élément pour projeter le noeud esclave sur le maitre.
je t'ai mis dans le répertoire joint également 2 images png qui montre la situation où le calcul diverge. J'obtiens cette divergence exactement pour la même configuration dans le cas où les outils sont quadracompl (...QC.png) et dans le cas à maillages entièrement linéaires (...LIN.png).
je pense que ce cas de figure est également similaire à un autre cas que je n'arrive pas à retrouver sur ce cas axisymétrique => cas d'un maillage esclave qui traverse le maitre (perte du contact). De mémoire, ce cas intervenait au moment où un noeud esclave doit passer d'un élément du maitre à un autre.
Mis à jour par Gérard Rio il y a presque 2 ans
- % réalisé changé de 30 à 40
J'ai trouvé la raison du pb initial: dans le cas axi, les coordonnées sont en 3D. Dans la recherche de contact avec une droite j'ai besoin d'une normale à la droite. Or en 3D général, je génère une normale aléatoirement dans le plan normal. En axi, tous les calculs doivent se faire dans le plan xy ce qui n'était pas le cas ici ! d'où des normales inexploitable en axi.
Sans doute que l'aléa ne marche pas exactement pareil en fast et en debug d'où un résultat différent.
Bon... je regarde le nouveau cas test
à suivre
Mis à jour par Gérard Rio il y a presque 2 ans
J'ai été étonné du fait que le calcul avait du mal à converger.
En regardant l'évolution des forces de contact on s'aperçoit que lorsque le calcul a du mal à convergence, au niveau des premières itérations les forces de contact sont gigantesques.
Je pense que c'est dû à la programmation du type 4:
- quand on a une pénétration anormalement importante, l'idée est d'augmenter le facteur de pénalisation proportionnellement au dépassement de pénétration constaté. Or si la pénétration visée est faible (ici elle est de 1.e-7) et que l'on veut des pas de charge en déplacement important (typiquement très supérieure à la pénétration visée) on obtient un facteur de pénalisation très important dans les premières itérations d'où une raideur mal conditionnée et une mauvaise convergence.
Il est possible de limité la force maxi. En regardant les valeurs, au niveau des premiers incréments on a des valeurs de forces de contact de l'ordre de e3, j'ai donc plafonné à e5, et là la convergence est ... bien bien meilleurs !!
Pour cela il faut utiliser le mot clé:
FORCE_CONTACT_NOEUD_MAXI 100000.
Avec ce paramètre et un nombre maxi d'incrément de 100, j'arrive à 0.82 de la charge maxi
Et en particulier le noeud du haut passe devant 2 noeuds maître sans pb a priori !!
bon... je relance le calcul sans limite d'incrément
affaire à suivre
Mis à jour par Gérard Rio il y a presque 2 ans
- Fichier clipboard-202301192014-si7n0.png clipboard-202301192014-si7n0.png ajouté
- % réalisé changé de 40 à 80
et bien je confirme que le calcul va à terme sans pb...
durant tout le calcul, il y a seulement 2 non-convergence d'itération, qui se situent en dehors du passage devant un noeud maître (pour le noeud du haut).
du coup cela semble dire qu'effectivement on peut avoir des pb de passage d'un élément à l'autre mais que la raison est peut-être du à des pb de convergence globale qui conduisent à des déplacements incohérents ??
Le contact avec pénalisation influe sur la convergence, et en retour dépend aussi de cette convergence...
Une force maxi fixe n'est pas forcément la solution optimale car si jamais la force réel est supérieur la convergence ne peut pas être atteinte. Du coup il faut forcément choisir un sup de ce que l'on imagine en fin de calcul ce qui peut être pénalisant en début de calcul (?).
Une solution pourrait-être de piloter cette force maxi via une fct nD qui utiliserait la force de contact calculée à l'incrément précédent et en prenant par exemple 10 fois cette force.
C'est théoriquement possible (cf.doc utilisateur) avec la syntaxe:
FORCE_CONTACT_NOEUD_MAXI FCT_ND_FORCE_CONTACT_NOEUD_MAXI <nom d’une fonction nD>
J'ai bien envie d'arrêter là le ticket vue les autres qui attendent ... ?
Mis à jour par Julien Troufflard il y a presque 2 ans
merci. Je vais regarder pour voir si il y a moyen de paramétrer.
aucun problème pour mettre en suspens le ticket. Il faut que je regarde également un cas de gonflage axisymétrique où j'avais une perte de contact. Et l'étude du paramétrage par défaut serait intéressant à regarder. etc...
En fait, la logique de ce ticket est de ne jamais le fermer jusqu'à qu'il soit possible de lancer des calculs avec contact dans une boucle d'identification inverse par exemple. Je me demande s'il y aurait un endroit plus approprié sur le site d'herezh pour ce qui me parait être plus un travail d'évolution qu'un ticket classique.
Mis à jour par Gérard Rio il y a presque 2 ans
- Statut changé de En cours à Résolu
- % réalisé changé de 80 à 100
j'ai créé un élément de la roadmap, qui correspond à la mise au point, test et évolution du contact.
Les prochaines questions relatives au contact pourront donc ce mettre dans ce ticket général lié au roadmap.
Du coup, je ferme ce ticket
Mis à jour par Gérard Rio il y a presque 2 ans
NB: il faut utiliser la dernière version V 7.006 : corr axi + contact
pour tester