Dans le cadre des simulations d'écrasement d'aspérités en présence de lubrifiant que je suis en train de mettre au point pour Arcelor, j'aimerais pouvoir effectuer des simulations 2D tout en faisant varier l'épaisseur hors plan de la géométrie de manière à prendre en compte l'allongement de la tôle et du lubrifiant dans la direction de laminage.

Pour y parvenir, je vais procéder en 2 étapes. Tout d'abord, permettre à l'utilisateur d'utiliser une épaisseur hors plan - en état plan de déformation - différente de 1. Cette première étape est l'objet de ce commit. La deuxième consistera à prendre en compte la variation de cette grandeur au cours de la simulation.

La modification principale a été apportée dans le calcul du Jacobien (cfr.SfIntegrationPointTemplate<Dimension2D>::computeJacobian(…)). La composante 33 du Jacobien ne vaut donc plus nécessairement 1 en état plan de déformation, mais cette valeur correspond désormais à l'épaisseur hors plan (hzz) définie par l'utilisateur. Pour prendre en compte ce paramètre dans le jeu de données, on l'introduit en argument de la fonction suivante:

geometry.setDimPlaneStrain(hzz)

A la demande de Luc, il n'y a pas de valeur par défaut, on indique donc explicitement la valeur de l'épaisseur utilisée dans le jeu de données. Cette valeur est ensuite transmise aux différentes interactions utilisées dans le calcul.

Pour le cas du contact faisant appel à AREAINCONTACT, le calcul de l'aire associée aux noeuds a été modifié (avec l'aide de Gaëtan): en pratique, on multiplie la valeur utilisée jusqu'ici par la moyenne des épaisseurs hors plan des éléments incidents aux noeuds de l'interaction de contact.

Sur les conseils de Luc, j'ai utilisé une valeur par défaut hzz = 0.0 dans les fonctions faisant appel au calcul du Jacobien (et pas 1.0) de manière à faire planter la simulation au cas où j'aurais oublié de passer l'épaisseur hors plan en argument qque part dans le code. C'était une bonne idée puisque ça m'a évite de faire une erreur dans le calcul des intégrales sur les cellules utilisées en ALE. L'épaisseur était jusqu'ici supposée unitaire en Etat plan de déformation, on considère que c'est toujours le cas en passant explicitement une épaisseur de 1.0 dans les fonctions. Comme l'épaisseur hors plan ne varie pas au cours de la phase de transfert de données, le fait d'utiliser une épaisseur unitaire n'a pas d'incidence sur le résultat.

Cette modification a été apportée dans la classe mtGeoSide ainsi que dans SfIntegratorTemplate (cfr. codes ci-dessous)

void
Side::computeIntegral( )
{
    ...

    for (int j=0; j<nbField; j++)
    {
        for(int numIP=0; numIP < nbOfIP; numIP++)
            unknownOnIP[numIP+1] = (*values[numIP])[j+1];
        integralValue[j+1] = integrator.integrate(xn, unknownOnIP,1.0); /* on passe une épaisseur hors plan = 1
        pour éviter les appels à la valeur définie dans le jeu de données qui pourrait occasionner une augmentation du temps CPU
        comme l'épaisseur ne varie pas lors du transfert de données ça n'a pas d'incidence sur le résultat. */
    }
    
    ...
}

/**
 * ATTENTION : LE VOLUME (AIRE) peut etre negatif.
 * Pas de test sur detJ
 */
template <class DIMSF>
double
SfIntegratorTemplate<DIMSF>::computeVolume(const std::vector<mtMath::Vect3> Co, bool isAxisym)
{
    typename std::vector<SfIntegrationPointTemplate<DIMSF> *>::iterator it;
    double volume = 0.0;
    mtMath::Matr3 jacobian;
    int i;
    for(it=integrationPoints.begin(), i=0; it != integrationPoints.end(); it++, i++)
    {
        (*it)->computeJacobian(Co, jacobian, 1.0, isAxisym); /* on passe une épaisseur hors plan = 1
        en état plan de déformation ce qui permet de calculer l'aire des éléments */
        volume +=  jacobian.determinant() * (*it)->getWeight();
    }
    return volume;
}

J'ai ajouté un répertoire appelé /gps/ (generalized plane strain) dans /apps/. Dans ce répertoire, j'ai ajouté les cas tests qui vérifient la validité des développements réalisés.

La première série de cas tests consiste à éraser un bloc carré de 1 mm de coté. Trois types de chargement différents sur le bord supérieur sont testés :

1. déplacement imposé
2. contact avec un outil rigide mobile
3. pression imposée

Pour chacun des cas, j'ai vérifié que la force de réaction sur le bord inférieur obtenue avec une épaisseur hors plan de 2 mm est bien 2 fois supérieure à celle calculée pour une épaisseur unitaire. Et ce, pour les éléments de type EAS - SRI - SRIPR.

[a]  apps/gps/eas_cont_ep1.py
[a]  apps/gps/eas_cont_ep2.py
[a]  apps/gps/eas_dep_ep1.py
[a]  apps/gps/eas_dep_ep2.py
[a]  apps/gps/eas_pImp_ep1.py
[a]  apps/gps/eas_pImp_ep2.py

[a]  apps/gps/sri_cont_ep1.py
[a]  apps/gps/sri_cont_ep2.py
[a]  apps/gps/sri_dep_ep1.py
[a]  apps/gps/sri_dep_ep2.py
[a]  apps/gps/sri_pImp_ep1.py
[a]  apps/gps/sri_pImp_ep2.py

[a]  apps/gps/sripr_cont_ep1.py
[a]  apps/gps/sripr_cont_ep2.py
[a]  apps/gps/sripr_dep_ep1.py
[a]  apps/gps/sripr_dep_ep2.py
[a]  apps/gps/sripr_pImp_ep1.py
[a]  apps/gps/sripr_pImp_ep2.py

[a]  apps/aej/cooks2Dep_ep2.py

[a]  apps/gps/tmCont_ep1.py
[a]  apps/gps/tmCont_ep2.py

[a]  apps/full/convectionEpd.py
[a]  apps/full/convectionEpd_ep2.py

[a]  apps/full/fluxEpd.py
[a]  apps/full/fluxEpd_ep2.py

[a]  apps/full/rayoEpd.py
[a]  apps/full/rayoEpd_ep2.py

[a]  apps.iso.tm2rayonnementEpd 
[a]  apps.iso.tm2rayonnementEpd_ep2 
[a]  apps.monosThermoMeca2.convection2dEpe_ep2 
[a]  apps.monosThermoMeca2.flux2dEpe_ep2

[a]  apps/XFEM/bloc_carre_ep2.epy  
variante du cas test existant "bloc_carre.py" 
-> test des éléments "XFEMTractionElement"
[a]  apps/XFEM/cont_ep2.py 
variante du cas test existant cont.py -> test des éléments "XFEMContactElement"