Herezh++

Les figures ne sont pas passées dans le message initial

FORCE_CONTACT :

NORMALE_CONTACT :

La raison est:
- la sortie de la normale est de norme 1 si le noeud fait partie d'un seul élément de contact et de n si le noeud fait partie de n éléments de contact.
Même moi je ne m'en rappelais plus, mais c'est marqué dans la doc.
En fait j'avais mis cela pour justement voir si le noeud était en contact multiple.

Bref ici il s'agit d'un noeud qui appartient à 2 éléments de contact. Mais ici ce n'est pas normal car il se trouve que les deux éléments de contact sont identiques ! et normalement il y a un filtre qui interdit cela ... il semble qu'il y ait quelque chose de pas normale !
Il faut que je trouve le pourquoi ...
à suivre

oui, tout à fait, mais c'était fait un peu à dessein.
Donc effectivement deux fois le même élément est comptabilisé ce qui double la force de contact.
Une des raisons de l'existence de ces 2 contacts simultanés est que lorsque l'on glisse d'une frontière à l'autre on peut avoir un contact sur 2 frontières, car le calcul prend en compte un léger débordement de la frontière pour le calcul de la projection et actuellement j'avais permis d'avoir un noeud en contact avec plusieurs frontières en même temps.
Du coup après réflexion voici ce que je propose:
- un noeud ne peut avoir qu'un seul contact par zone de contact.
Je pense que le fait de pouvoir définir plusieurs zones de contact permet quand d'avoir plusieurs contacts, mais ceux-ci sont alors bien identifiés. Par exemple on peut avoir sur un même noeud esclave un contact vertical + un contact horizontal

Et donc lorsque je détecte 2 contacts sur une même zone, je choisis celui dont la distance du noeud à tdt à sa projection est la plus faible.

a priori cela semble assez bien fonctionner sur l'exemple de l'écrasement du joint.
- il n'y a plus de pb de normale
- la convergence est ok

Par contre au cours de mes tests j'ai relevé les points suivants:
1) comme la loi proposée est hyperélastique avec un module de compressibilité très grand par rapport au module de cisaillement, à partir d'un certain taux de compression, la convergence globale devient difficile. Cela n'a, apriori rien à voir avec le contact, c'est tout simplement que l'on observe un effet de blocage numérique venant du fait que la partie sphérique écrase complètement la partie déviatorique au niveau de la raideur, et cette dernière manque pour l'unicité de la solution.
Je fais le même essai avec un K 10 fois plus petit et là aucun pb pour écraser fortement: hauteur initiale 7.8 , déplacement imposé -5.6
donc il faut se méfier des non-convergences supposées dues au contact. Malgré le fait que le contact peut entraîner pas mal de pb sur la convergence.
2) toujours concernant la convergence due au contact, j'ai noté que lorsque l'on atteignait une forte compression, la géométrie finale avait tendance à osciller, ceci due à ce que j'ai expliqué au 1). Du coup certains noeuds en contact avaient tendance momentanément à sortir du maître. Or dans la mise en données on indique:
PENETRATION_BORNE_REGULARISATION 1e-5
ce qui veut dire que si le noeud sort de plus de 1.e-5, on ne tient plus compte de sa réaction. Et du coup sur la visualisation on voit le noeud qui s'écarte violemment du maître avec lequel il n'est donc plus accroché.
Pour éviter ce phénomène, il faut introduire un autre paramètre:
PENETRATION_BORNE_REGULARISATION_PLUS 0.1
Dans ce cas tant que le noeud ne s'écarte pas plus de 0.1, le noeud reste accroché.

Donc en utilisant ce paramètre, cela permet d'éviter des décollements intempestifs lorsque globalement la solution n'est pas parfaitement stable.

Néanmoins pour la loi considérée dans l'exemple, il arrive un moment où le calcul patine due au 1) et là il n'y a pas de solution en implicite. En tout cas je n'en ai pas trouvé. Peut-être en RD ?

La nouvelle version est la : V 7.061

NB: On peut également y tester les sorties individuelles pour chaque maillage

ok pour les pb de compressibilité.

J'ai construit un cas test en AXI assez intéressant. Au départ, c'était pour valider cette nouvelle version Herezh. Mais finalement, il montre des pb pour les 2 versions : celle que j'utilise (7.056) et la dernière 7.061.

ci-joint les fichiers de calcul dans test.tar, le résultat Gmsh par version et le log des 2 calcul.

Pour ce calcul, j'ai pris une loi Mooney-Rivlin assez compressible, équivalente en petite déf à une loi élastique log E=100, NU=0.35

Le principe est un carré (maillage plot_AXI.her) écrasé dans un logement en U (logement_AXI.her). La situation est la suivante :

Comme indiqué, le but est d'avoir une config qui teste :
- le cas d'un noeud esclave qui franchit un noeud maitre (ici 2 cas : un cas horizontal et cas vertical).
- un noeud dans le coin qui va subir un contact vertical + horizontal
- il y a un petit jeu initial selon Y de 1e-4 (pas de contact initial) et un gros jeu initial selon X pour décomposer le contact vertical/horizontal au niveau du coin

Autre point : le calcul utilise un suite_point_info pour décomposer en 2 phases
1) Il ne se passe rien entre t=0 et t=1 => le maillage du plot ne doit subir aucune déformation. Il ne doit y avoir aucune force de contact (c'est en référence au ticket https://herezh.irdl.fr/issues/391 car je testerai ensuite en 2D et ainsi qu'en 3D en faisant une extrusion Z de ces maillages)
2) L'écrasement a lieu entre t=1 et t=2 par déplacement imposé selon Y sur le bord du haut du plot => le contact vertical doit être détecté dès le premier pas de temps, puis interviendra le contact horizontal

Pour des raisons de nombre de fichiers pièce jointe, je vais montrer quelques résultats dans les messages suivants.

résultats pour l'ancienne version 7.056

pas pb pour phase 1) du calcul entre t=0 et t=1 => pas de déformation, pas de force de contact (le contact est quand même détecté car il y a une normale de contact en certains noeuds). Donc ok pour cette partie.

Ensuite, dans la phase d'écrasement, voici quelques exemples :

au temps 1.14 : le noeud du coin se prépare au contact horizontal. Il a bien une normale à 45°.
Problème : il a par contre une norme de normale égale à sqrt(2) au lieu de 1

au temps 1.21 : le plot est au contact des 2 parois avec des forces de contact non nulles selon X et Y
Problème : la norme de contact au coin est maintenant inférieure à 1

au temps 1.23 : il y a eu franchissement d'un noeud esclave en bas.
Problème : on retrouve le problème initial de ce ticket. Ce noeud a désormais une norme de contact de norme 2

au temps 1.9 : malgré les défauts de norme, le contact se passe correctement avant ce temps.
Problème : à ce moment du calcul, le coin perd le contact et pénètre dans le maitre. La normale de contact n'a plus ni la bonne norme, ni la bonne direction

au temps final : le calcul converge malgré la perte totale de contact du noeud coin.

ps : oubli de préciser pour la version 7.056 => il n'y a pas de pb de norme de normale pour le noeud esclave en haut à droite qui franchit un noeud maitre sur la paroi verticale. Elle reste à 1 contrairement au cas du noeud du noeud en bas à gauche sur la paroi horizontale. bizarre ??

résultats pour la nouvelle version 7.061

ok pour la phase 1) => pas de déformation ni forces de contact

pour la partie écrasement : en général, il n'y a plus de pb de norme de normale pour tous les noeuds. Le souci est désormais uniquement sur le noeud coin qui oscille entre une normale purement horizontale ou purement verticale avec une impossibilité de trouver le contact correctement.

Voici 2 exemples parmi tant d'autres :

au temps 1.21 : normale purement verticale et noeud a pénétré dans le maitre

au temps 1.26 : idem avec un normale purement horizontale

au temps final : on aboutit à peu près à la même situation qu'avec version 7.056 => calcul converge malgré un coin qui a perdu totalement le contact.

oui tu as sûrement raison. Je crois même que ça fait parti des bonnes pratiques EF d'ajouter un petit rayon aux angles vifs, même fictifs (i.e pas issu d'une pièce réelle).
Donc, il faudra bel et bien penser à décomposer les zones de contact. J'y étais réticent, mais finalement, ça parait d'une part logique et d'autre part largement faisable dans la grande majorité des cas que l'on a traiter.

Si je décompose le logement en 2 zones de contact, ça conduit les 2 versions Herezh au bon résultat final (pas de pénétration).
avec néanmoins une normale étrange dans le cas de l'ancienne version 7.056

version 7.056 :

1) avant que le contact ne devienne effectif sur la partie verticale, la normale au coin est dans la bonne direction

2) ensuite, la normale n'est plus dans le bon sens. L'exemple ci-dessous est représentatif du reste du calcul avec une normale au coin rentrante dans la matière et une normale de norme 2 au noeud en bas à gauche du plot. Mais néanmoins, le calcul converge avec un contact correct. Bizarre pour que ça fonctionne vue la normale inversée??? (Mais bon inutile d'investiguer puisque ce ne sera pas le cas avec la version 7.061)

version 7.061 :

1) l'exemple ci-dessous est représentatif de la quasi intégralité des incréments du calcul => normale au coin dans la bonne direction

2) dans les derniers incréments, la normale au coin devient purement verticale. Mais ça n'empêche pas le caclul de converger et avec un contact correct sans pénétration.

Au final, on a une norme au coin qui peut être différente de 1 (somme des normales de contact). Je ne sais pas si c'est un réel problème. A priori non (chaque contact est calculé avec une norme 1 multiplié par un coef de pénalisation => la somme des 2 donne une résultante pour l'affichage Gmsh).

Par curiosité, je vais tenter d'écraser un peu plus (pour voir si la normale finale pose problème, car c'est quand même étrange de ne plus avoir une résultante horizontale+verticale).
Puis ensuite je vais faire l'ajout d'un rayon de courbure au coin et voir ce que ça donne avec une zone de contact unique. Et mettre le résultat ici.

Ensuite je vais donc passer au cas 2D, et donc sans doute basculer sur le ticket https://herezh.irdl.fr/issues/391

j'ai oublié de te donner les nouvelles versions de maillage logement_AXI.her et fichier .info (avec la déclaration de 2 zones de contact et le déplacement -0.75)