Herezh++

Bonjour Julien,

je pense savoir ce qui se passe et ... je crois qu'il n'y a pas grand chose à faire.
Une CLL consiste a imposer une relation entre des ddl de noeuds. Pour cela, classiquement on dispose de deux méthodes : les multiplicateurs de Lagrange et la pénalisation. La méthode actuellement implantée correspond à un multiplicateur de Lagrange + une condensation statique de l'inconnue supplémentaire correspondant au multiplicateur.
La condensation du multiplicateur s'effectue via une modification de la raideur et/ou matrice masse et du second membre via une multiplication par une matrice bande.

Or lorsque l'on multiplie une matrice par une matrice bande, la largeur de bande en noeud finale vaut la somme des bandes de chaque matrice initiale -1 (si je me rappelle bien).
Donc à moins d'avoir une matrice multiplicative diagonale, il y a systématiquement augmentation de la largeur de bande à chaque fois que l'on ajoute une CLL.

Au bilan, je pense que la méthode actuelle est intéressante pour quelques CLL (1 à 5 par exemple) mais devient très rapidement abominable en temps cpu, due à l'augmentation catastrophique de la taille des matrices.

Je pense donc que dans l'avenir je vais implanter d'autres méthodes: pénalisation (à l'image du contact qui pour une partie est un cas particulier de CL non linéaire), et peut-être multiplicateur.

Donc, pour l'instant je n'ai pas de réponse satisfaisante à ton problème concernant les CLL.

Peut-être serait-il possible d'imaginer une solution en "boucle ouverte", c-a-d à chaque itération
- calcul du centre de gravité d'un set de noeuds > ses coordonnées
utilisation de ces coordonnées dans une condition limite classique
tout cela via une fonction nD externe:
> les coordonnées du set de noeuds sont définit en variables globales
> la fct nD est définie pour piloter la conditions limites classique

bon... à voir

j'ai fait un premier test. J'ai déclaré plus de 23000 variables utilisateurs. Sur ce point, Herezh n'a pas l'air de faire de difficultés pour les lire.
Ensuite, je définis une fonction nD qui utilise plus de 2000 de ces variables. Herezh++ n'a pas l'air capable de lire autant de variables sur une seule ligne.

est-ce que c'est normal ou réparable ?

quel type de fonction nD ? littérale ?
ce n'est peut-être pas la lecture qui pose pb, mais l'appel d'initialisation d'une fonction via muparser ? je pense que le nombre de variable dans muparser est limité (5 ou 7 ??)
c'est pourquoi je pensais à une fct nD externe qui n'aura pas cette limitation.
Ceci étant il y a quand même actuellement une limitation pour les fct nD externe dans le passage de paramètres.
Au niveau des pipes de dialogue, j'ai limité pour l'instant à environ 100 réels (928 caractères exactement cf. doc).
Je ne sais pas trop ce que cela ferait de * par 20, notamment au niveau efficacité.
Il faudrait peut-être faire plusieurs appels ?? à moins que le système s'en charge... bref il faudrait approfondir
Mais en tout cas cela semble possible a priori.
Par contre, ce serait bien de savoir avant si la manip 100 réels fonctionne et permet une convergence en boucle ouverte...

oui c'est une fonction littérale.

donc que ce soit externe ou interne, il y aura une limitation. Je préfère pour l'instant ne pas me lancer dans les fonctions externes par simplicité.

J'ai continué à fouiller le cas des fonctions internes.
je suis passé à un second test dont voici une archive.

Alors c'est un peu compliqué parce que plusieurs choses se téléscopent. L'idée générale était de vérifier le résultat de mes fonctions nD de calcul de centre de gravité. Je ne connais pas de moyen direct pour sortir le résultat d'une fonction nD. Donc j'ai décidé de vérifier indirectement en les utilisant pour régir la position d'un noeud via une condition en position dans les "blocages". Comme ça, j'aurais vu de mes yeux en regardant le résultat dans Gmsh (même logique que dans mon cas précédent : pb_calcul_2elts_barre.tar).

D'autre part, j'ai voulu aussi voir combien de variables utilisateur je pouvais utiliser dans une fonction nD suite à mon message précédent.

Voici en gros le cheminement :
1) j'ai voulu faire un test sur une plaque à 10000 noeuds. Je lui ajoute une barre. Je fixe un noeud de la barre et le but est de régir la position de l'autre noeud de la barre via des fonctions nD qui calculent un centre de gravité. Le centre de gravité sera celui d'un groupe de noeuds de la plaque. Je fais bouger la plaque dans l'espace pour voir si le noeud de la barre suit bien le centre de gravité.

Le fichier .info correspondant à ce calcul est : test.info
Le noeud de la barre qui est sensé suivre le centre de gravité est : N_biel_2

2) problème : pour vérifier si ça marche, j'ai voulu appliquer une condition limite sur la position absolue du noeud N_biel_2 pour voir de visu s'il suit le nuage de points dans l'animation Gmsh. Mais il n'est pas possible actuellement de faire une condition limite sur la position (ddl : X1 X2 X3).
Si on met une condition du genre 'X1= Fonction_nD_CHARGE: etc...', Herezh++ considère qu'il s'agit d'un déplacement (donc UX UY UZ).
Pour contourner cela, j'ai voulu passer par des cll. ça fonctionne pour fixer des positions en mettant "enu= X1". Mais a priori, il n'y a pas moyen d'utiliser des fonctions nD dans les cll (dans AvecFonctionsDeCharges_). Donc je n'ai pas pu aboutir.

=> si tu veux tester cll+fonctions nD, j'ai mis un fichier test_cll_fct_nd.info à ce sujet. En l'état, Herezh renvoie une erreur avec ce .info (il s'attend à trouver des courbes 1D). Ce n'est pas ma problématique mais si tu veux modifier cela, ça sera sûrement très utile pour tout le monde.

dans l'impossibilité de régir en position, j'ai donc décidé de régir en déplacement. Le moyen de vérification sera alors de regarder via le .maple si les déplacements (UX UY UZ) du noeud N_biel_2 est bien égal à la position (X1 X2 X3) du centre de gravité des noeuds de la plaque.

Pour cela, j'ai créé un script verif_U_cdg.pl qui permet de comparer à chaque temps le déplacement de N_biel_2 au centre de gravité du nuage de points. Pour l'utiliser, il faut indiquer le nombre de noeuds du nuage de points (exemple : verif_U_cdg.pl 3 ; cf Readme).

3) dans un premier temps, le but était de voir combien de variables utilisateur je pouvais mettre dans les fonctions nD sans faire planter. Pour se faire, j'ai créé un script dont je rappelle le rôle dans le fichier Readme. Ce script (genere_fic.pl) permet de créer les fichiers de variables utilisateurs (cdg_fps.var_util), de fonctions nD (cdg_fps.fct_nD) et un fichier .lis qui contient le set des noeuds de la plaque qui serviront de nuage de points pour calculer un centre de gravité (cdg_fps.lis). On lui donne en argument un nombre de noeuds.
Tu pourras toi aussi faire varier le nombre de noeuds de la plaque utilisés pour le centre de gravité si tu veux. Les fichiers actuels cdg_fps.* ont été générés pour 3 noeuds (Les 3 premiers noeuds en bas à gauche de la plaque)

j'en ai déduit que :
- pour 200 noeuds, Herezh++ n'arrive pas à lire le fichier .info.
- pour 199 noeuds, le calcul se lance, arrive au bout, les fichiers résultats sont créés correctement, mais néanmoins j'ai un segmentation fault comme dernier affichage (ça donne ça : ============================================================= | fin HEREZH++ | =============================================================
Segmentation fault

aucune idée si ça altère quelque chose quelque part. Je n'ai pas constaté de pb dans Gmsh. Tout à l'air écrit correctement.

- pour 3 noeuds, le calcul tourne et finit sans erreur
- j'ai testé à nouveau pour 3 noeuds dans les fonctions nD mais en générant le fichier cdg_fps.var_util pour 10000 noeuds (soit une déclaration de 30000 variables utilisateurs) => le calcul tourne sans problème et l'occupation mémoire n'a augmenté que de 6 ou 7 Mo (le cas précédent utilisait déjà 215 Mo de ram, c'est donc une augmentation négligeable)

donc grosso modo en conclusion, il y a moyen d'utiliser 199 variables dans les fonctions nD et l'occupation mémoire n'est pas un problème. Cependant un pb mineur de segmentation fault serait à vérifier dans le cas à 199 noeuds.

4) dans un second temps, j'ai vérifié le résultat des fonctions nD. Donc en regardant ce qu'affiche le script verif_U_cdg.pl. J'obtiens un affichage du genre :
temps=2.500000000000e-01 :

deplacement N_biel_2 => (0 1 0.000000000000e+00)

deplacement N_biel_2 => (2.5 3.5 2.500000000000e+00)

deplacement N_biel_2 => (5 6 5.000000000000e+00)

deplacement N_biel_2 => (7.5 8.5 7.500000000000e+00)

on voit bien que le "deplacement N_biel_2" d'un certain temps est égal à "position cdg" du temps précédent. Donc, les fonctions nD fonctionnent bien. Mais je constate que les variables utilisateurs (V__X1_1, etc...) envoient la valeur de position des noeuds au début de l'incrément donc à t et non à t+dt. Et pourtant j'ai bien mis l'option TEMPS_tdt dans le fichier cdg_fps.var_util

ce point me pose un réel problème car dans mon futur cas de ballon, la position t+dt risque d'être radicalement différente de la position à t. Il est donc impératif que mes fonctions nD exploitent bien la valeur à t+dt courante à l'itération i de relaxation dynamique. Y a-t-il une solution ?

ok. merci.
oui c'est dense, je dirais même que ça a été une lutte.

pour compléter, j'ai été revoir un ticket précédent à propos de l'utilisation des variables utilisateurs au cours des itérations. Dans le ticket #265 ( https://herezh.irdl.fr/issues/265 ), tu m'y indiquais que :
"""""
À condition que tu indiques que le temps considéré pour la variable est "TEMPS_tdt", pendant toute l'itération, c-a-d durant les calculs intermédiaires entre l'itération i et i+1, la variable globale correspondante au ddl aura comme valeur celle obtenue à l'issue de l'itération i.
Donc par exemple si durant l'itération, le ddl change, la variable globale n'est pas modifiée.
Plus précisément, supposons que le ddl soit == une coordonnée d'un noeud noe, et qu'il y a une condition de contact sur le noeud noe.
Durant l'itération, le noeud peut éventuellement bouger du à la méthode de contact (dans herezh via la méthode de contact 1), mais cela ne changera pas la valeur de la variable globale correspondante. Celle-ci sera mise à jour uniquement après la résolution de l'équation d'équilibre global.
"""""

je me suis donc dit que mon cas précédent était peut-être spécial puisque tous les ddl ont une condition, donc il n'y a pas vraiment d'itérations, le résidu est toujours nul => ça converge tout de suite. Bon, en réalité, dans le cas de la relaxation dynamique, Herezh++ fait quand même 3 ou 4 itérations par incrément malgré 0 dans la norme de convergence, ça je ne sais pas pourquoi. bref.

et donc, supposons qu'il n'y ait qu'une seule itération par incrément. A une itération i+1, les variables utilisateurs sont sensées avoir pour valeur les données à l'itération i. En l'occurrence, l'itération i serait la situation au début de l'incrément donc avec les données à t. Et donc cela expliquerait pourquoi j'ai un incrément de retard comme précisé dans mon point 4) de mon message précédent.

J'ai donc refait un nouveau test qui forcerait Herezh à faire plusieurs itérations d'équilibre à chaque incrément. J'ai repris le test précédent, donc plaque + barre. J'ai :
- réduit la taille du maillage pour que ça calcule plus vite
- passé en mode non_dynamique
- mis des conditions limites sur la plaque qui lui font subir une sorte de cisaillement avec un côté fixe et côté mobile

Comme prévu, il y a plusieurs itérations à chaque incrément. Mais malgré tout, je retrouve la même chose que précédemment, à savoir que tout se passe comme si les valeurs des variables utilisateurs restent égales à leurs valeurs à t.

Dans le répertoire joint, le calcul est déjà fait. Comme montré dans le fichier test.log, il y a eu plusieurs itérations par incrément. Si on lance "verif_U_cdg.pl 3", on obtient l'affichage suivant :
temps=2.500000000000e-01 :

deplacement N_biel_2 => (0 0 0.000000000000e+00)

deplacement N_biel_2 => (-2.2776631516 2.8064943535 2.500000000000e+00)

deplacement N_biel_2 => (-4.5562160268 0.614214886900001 5.000000000000e+00)

deplacement N_biel_2 => (-6.83494044989999 -1.57484878243 7.500000000000e+00)

encore une fois, le déplacement "deplacement N_biel_2" à un temps est égal à la position "position cdg" du temps précédent. Les positions intermédiaires prises au cours des itérations n'ont pas l'air de modifier les variables utilisateurs.

oui, j'ai trouvé ce qui se passe.
En fait en implicite, les conditions limites sont imposées avant les itérations et ne sont pas modifiées pendant les itérations.
A priori c'est normal car les itérations ont pour objectif de trouver une solution pour les ddl libres, tandis-que les ddl bloqués restent fixes. Du coup, normalement il n'y a pas à les modifier donc a regarder la fct nD qui régie la condition.
D'où effectivement il n'y a pas de mise à jour de la CL pendant les itérations.

Ce fonctionnement me paraît normal mais ne convient pas dans ton cas.
Si je mets à jour les CL bloqués à chaque itération, cela ne pose pas de pb a priori, par contre cela va entraîner un temps de calcul supplémentaire inutile dans la majorité des cas. Du coup je ne sais pas trop quoi faire.
Peut-être une solution serait d'introduire un indicateur supplémentaire dans la mise en données qui indiquerait que pour une CL donnée, il faut en tenir compte à chaque itération....
Cela introduirait un peu de temps cpu ... car il faudra à chaque itération reconsulter toutes les CL... cela ne me parait pas terrible,

Je pourrais aussi définir un mot clé au niveau de l'algo lui-même qui globalement, imposerait de reconsulter les CL à chaque itération ??
C'est sans doute ce dernier cas qui serait le moins impactant pour des calculs classiques (il serait même quasiment inodore) et cela simplifierait la mise en données...
Qu'en penses-tu ?

suite à l'évolution #274 (cf. roadmap) la version 6.977 intègre
- la possibilité de relire les CLL à chaque itération
- la possibilité de faire une statistique sur une ref de noeud (cf. doc) ce qui génère une variable globale. Ensuite utilisation d'une variable relais pour utiliser les infos de la statistique dans une fonction nD

Me dire si cela permet de répondre à la demande.
Merci

j'ai ajouté l'option dans le fichier pièce jointe. ça ne change rien (j'espère l'avoir mis au bon endroit)