Download g_theta.procedur

Back to the list

   1 : * G_THETA   PROCEDUR  BP208322  13/12/20    21:15:08     7890           
   2 : DEBPROC G_THETA SUPTAB*'TABLE';
   3 : *|=====================================================================|
   4 : *|                                                                     |
   5 : *|    OBJECTIF :                                                       |
   6 : *|    ==========                                                       |
   7 : *|                                                                     |
   8 : *| 1) calculer l'integrale caracteristique de mecanique de la rupture  |
   9 : *|    a) J en elasto-plasticite ou en elasto-dynamique pour un         |
  10 : *|       materiau isotrope.                                            |
  11 : *|    b) dJ/da en elasto-plasticite, utilisable uniquement dans le cas |
  12 : *|       de materiau isotrope et homogene pour les elements massifs.   |
  13 : *|    c) C* dans le cas de fluage secondaire stationnaire pour un      |
  14 : *|       materiau isotrope.                                            |
  15 : *|    d) C*H dans le cas de fluage primaire sous un chargement radial  |
  16 : *|       pour un materiau isotrope.                                    |
  17 : *|                                                                     |
  18 : *| 2) separer les modes K1, K2 et K3 en elasticite, utilisable         |
  19 : *|    uniquement dans le cas de materiau homogene et isotrope          |
  20 : *|    pour les elements massifs.                                       |
  21 : *|                                                                     |
  22 : *|                                                                     |
  23 : *|    ENTREE :                                                         |
  24 : *|    ========                                                         |
  25 : *|                                                                     |
  26 : *| SUPTAB  objet de type TABLE. En entree, SUPTAB sert a definir les   |
  27 : *|         options et les parametres du calcul. Ses indices sont des   |
  28 : *|         objets de type MOTS (a ecrire en toutes lettres) dont voici |
  29 : *|         la liste :                                                  |
  30 : *|                                                                     |
  31 : *|                                                                     |
  32 : *|    Arguments obligatoires dans tous les cas                         |
  33 : *|    ----------------------------------------                         |
  34 : *|                                                                     |
  35 : *| SUPTAB.'OBJECTIF' = MOT pour preciser le but du calcul, vaut        |
  36 : *|                     1) 'J' pour calculer l'integrale J (ou G),      |
  37 : *|                         caracteristique en elasto-plastique.        |
  38 : *|                     2) 'J_DYNA' pour calculer l'integrale J (ou G), |
  39 : *|                         caracteristique en elasto-dynamique.        |
  40 : *|                     3) 'C*' pour calculer l'integrale C*,           |
  41 : *|                         caracteristique en fluage secondaire        |
  42 : *|                         stationnaire.                               |
  43 : *|                     4) 'C*H' pour calculer l'integrale C*(h),       |
  44 : *|                         caracteristique en fluage primaire ou       |
  45 : *|                         tertiaire.                                  |
  46 : *|                     5) 'DJ/DA' pour calculer l'integrale de la      |
  47 : *|                         derivation dJ/da, caracteristique pour      |
  48 : *|                         analyser la stabilite de propagation d'une  |
  49 : *|                         fissure ou des fissures interagissantes.    |
  50 : *|                     6) 'DECOUPLAGE' pour decouper les modes mixtes, |
  51 : *|                         c'est a dire la separation des facteurs K1, |
  52 : *|                         K2 (et K3 et 3D).                           |
  53 : *|                                                                     |
  54 : *| SUPTAB.'COUCHE' = ENTIER representant le nombre de couches          |
  55 : *|                   d'elements autour du front de la fissure          |
  56 : *|                   qui se deplacent pour simuler la propagtion       |
  57 : *|                   de la fissure. Il vaut 0 si seul la pointe de     |
  58 : *|                   la fissure se deplace, 1 si c'est la premiere     |
  59 : *|                   couche d'elements entourant la fissure, 2 si      |
  60 : *|                   c'est l'ensemble des premiere et deuxieme couches |
  61 : *|                   d'elements etc. Il convient veiller a ce que      |
  62 : *|                   l'ensemble des elements a deplacer n'atteint pas  |
  63 : *|                   le bord de la structure fissuree.                 |
  64 : *|                   Cet argument doit etre absent si l'on souhaite    |
  65 : *|                   preciser soi-meme le CHAMP_THETA (cf.8.)          |
  66 : *|                                                                     |
  67 : *| SUPTAB.'FRONT_FISSURE' = POINT en 2D ou MAILLAGE en 3D massif       |
  68 : *|                          representant le front de la fissure.       |
  69 : *|                                                                     |
  70 : *|                                                                     |
  71 : *|    Arguments obligatoires avec des elements standards               |
  72 : *|    --------------------------------------------------               |
  73 : *|                                                                     |
  74 : *| SUPTAB.'LEVRE_SUPERIEURE' = Selon la convention de definition, cet  |
  75 : *|                             objet (type MAILLAGE) representant la   |
  76 : *|                             levre superieure de la fissure.         |
  77 : *|                                                                     |
  78 : *| SUPTAB.'LEVRE_INFERIEURE' = Selon la convention de definition, cet  |
  79 : *|                             objet (type MAILLAGE) representant la   |
  80 : *|                             la levre inferieure de la fissure. Si   |
  81 : *|                             une seule levre est modelisee, un des   |
  82 : *|                             des deux mots ici (LEVRE_SUPERIEURE ou  |
  83 : *|                             LEVRE_INFERIEURE) sera suffisant pour   |
  84 : *|                             decrire la fissure.                     |
  85 : *|                                                                     |
  86 : *|                                                                     |
  87 : *|    Arguments obligatoires avec des elements enrichis (XFEM)         |
  88 : *|    --------------------------------------------------------         |
  89 : *|                                                                     |
  90 : *| SUPTAB.'PSI' =   1ere level set (CHPOINT) decrivant la fissure dans |
  91 : *|                  le cas ou l'on utilise des elements XFEM .         |
  92 : *| SUPTAB.'PHI' =   2eme level set.                                    |
  93 : *|                                                                     |
  94 : *|                                                                     |
  95 : *|                                                                     |
  96 : *|    Solution obligatoire issus de la procedure PASAPAS               |
  97 : *|    --------------------------------------------------               |
  98 : *|                                                                     |
  99 : *| SUPTAB.'SOLUTION_PASAPAS' = TABLE sortant de la procedure PASAPAS.  |
 100 : *|                                                                     |
 101 : *|                                                                     |
 102 : *|    Solution obligatoire issus de l'operateur RESO                   |
 103 : *|    ----------------------------------------------                   |
 104 : *|                                                                     |
 105 : *| SUPTAB.'SOLUTION_RESO' = CHPOINT de deplacement issus de RESO.      |
 106 : *| SUPTAB.'CARACTERISTIQUES' = Champ de caractristiques matrielles     |
 107 : *|                             et eventuellement geometriques          |
 108 : *|                             si necessaire.                          |
 109 : *| SUPTAB.'MODELE' = Objet modele (type MMODEL) englobant toute la     |
 110 : *|                   structure.                                        |
 111 : *| SUPTAB.'TEMPERATURES' = CHPOINT de temperature creant une contrainte|
 112 : *|                         thermique non nulle si elle existe.         |
 113 : *| SUPTAB.'CHARGEMENTS_MECANIQUES' = CHPOINT representant l'ensemble   |
 114 : *|                                   des forces exterieures            |
 115 : *|                                   (surfaciques, volumiques ou       |
 116 : *|                                   ponctuelles ....) appliquees sur  |
 117 : *|                                   le systeme si elles existent.     |
 118 : *| SUPTAB.'BLOCAGES_MECANIQUES' = RIGIDITE representant le blocages    |
 119 : *|                                mecanique du probleme, a fournir     |
 120 : *|                                uniquement dans le cas de calcul     |
 121 : *|                                de la derivation dJ/da.              |
 122 : *|                                                                     |
 123 : *|                                                                     |
 124 : *|    Arguments optionnels                                             |
 125 : *|    --------------------                                             |
 126 : *|                                                                     |
 127 : *|                                                                     |
 128 : *|    1 : Materiaux composites (2D massif ou 3D coque seulement)       |
 129 : *|                                                                     |
 130 : *| SUPTAB.'MODELES_COMPOSITES' = TABLE indicee par des entiers (1 2... |
 131 : *|                               M, M = nombre de Materiaux composites)|
 132 : *|                               pour donner les modeles des materiaux |
 133 : *|                               ayant des discontinutes de proprietes |
 134 : *|                               materielles.                          |
 135 : *|                                                                     |
 136 : *|    2 : Pour un front de fissure tridimensionnel massif              |
 137 : *|                                                                     |
 138 : *| SUPTAB.'NOEUDS_AVANCES' = MAILLAGE de type POI1 pour donner les     |
 139 : *|                           points du front pour lesquels le calcul   |
 140 : *|                           sera effectue. Si cet argument est        |
 141 : *|                           obsent, le calcul sera fait pour tous     |
 142 : *|                           les noeuds sur le front de la fissure.    |
 143 : *|                                                                     |
 144 : *|    3 : Calcul des termes croises de la matrice dJi/daj              |
 145 : *|        (i non egal a j) dans le cas des fisures interagissantes.    |
 146 : *|                                                                     |
 147 : *| SUPTAB.'FISSURE_2' = Objet de type MAILLAGE representant une autre  |
 148 : *|                      fissure (levres superieure + inferieure si     |
 149 : *|                      toutes les deux levres sont presentes).        |
 150 : *| SUPTAB.'FRONT_FISSURE_2' = POINT ou MAILLAGE reprsentant le front   |
 151 : *|                            de la fissure 2 decrite ci-dessus.       |
 152 : *|                                                                     |
 153 : *|                                                                     |
 154 : *|    4 : Cas d'une fissure circulaire dans une geometrie plane        |
 155 : *|                                                                     |
 156 : *| SUPTAB.'POINT_CENTRE'  = centre de la fissure circulaire            |
 157 : *|                                                                     |
 158 : *|    5 : Cas ou l'extension de la fissure correspond a une simple     |
 159 : *|        translation dans un tuyauterie droite (3D). Dans ce cas      |
 160 : *|        on effectue dans la procedure CH_THETA une transformation    |
 161 : *|        de tuyau en plaque en passant au systeme de coordonnees      |
 162 : *|        cylindriques. Il est alors necessaire de fournir :           |
 163 : *|                                                                     |
 164 : *| SUPTAB.'POINT_1' = centre du systeme de coordonnees                 |
 165 : *| SUPTAB.'POINT_2' = POINT tel que l'axe defini par POINT_1           |
 166 : *|                    vers POINT_2 soit l'axe Z poisitif               |
 167 : *| SUPTAB.'POINT_3' = POINT tel que le plan defini par les 3 points    |
 168 : *|                    POINT_1 POINT_2 POINT_3 donne l'angle theta nul  |
 169 : *|                                                                     |
 170 : *|    6 : Cas ou l'extension de la fissure ne correspond               |
 171 : *|        pas a une simple translation (3D)                            |
 172 : *|                                                                     |
 173 : *|      6.1 Fissure dans un tuyauterie droite (3D, Rotation)           |
 174 : *|                                                                     |
 175 : *| SUPTAB.'POINT_1' = Objet de type POINT                              |
 176 : *| SUPTAB.'POINT_2' = Objet de type POINT qui, avec le point POINT_1,  |
 177 : *|                    constitue l'axe perpendiculaire a la section     |
 178 : *|                    fissuree.                                        |
 179 : *|                                                                     |
 180 : *|      6.2 Fissure dans un coude (3D, rotation + transformation)      |
 181 : *|          Outre les deux points SUPTAB.'POINT_1' et SUPTAB.'POINT_2' |
 182 : *|          definis en haut on donne encore :                          |
 183 : *|                                                                     |
 184 : *| SUPTAB.'CHPOINT_TRANSFORMATION' = Objet de type CHPOINT utilise     |
 185 : *|                                   pour transformer une coude en un  |
 186 : *|                                   tuyauterie droite.                |
 187 : *| SUPTAB.'OPERATEUR' = Objet de type MOT valant 'PLUS' ou 'MOIN' pour |
 188 : *|                      indiquer l'operateur PLUS ou MOIN a utiliser   |
 189 : *|                      si l'on veut transformer la coude en un        |
 190 : *|                      tuyauterie droite.                             |
 191 : *|                                                                     |
 192 : *|    7 : Rotation rigidifiante imposee dans le calcul par PASAPAS     |
 193 : *|                                                                     |
 194 : *| SUPTAB.'ROTATION_RIGIDIFIANTE' = table indicee par entiers 0,1,2... |
 195 : *|                                  donnant les champs de deplacements |
 196 : *|                                  due a une rotation rigidifiante de |
 197 : *|                                  la piece autour d'un point. Cette  |
 198 : *|                                  rotation rigidifiante est imposee  |
 199 : *|                                  dans le calcul par PASAPAS en tant |
 200 : *|                                  d'un calcul en grand deplacement.  |
 201 : *|                                                                     |
 202 : *|    8 : Cas ou on souhaite donner soi-meme le champ THETA            |
 203 : *|                                                                     |
 204 : *| SUPTAB.'CHAMP_THETA' = Objet de type CHPOINT caracterisant la       |
 205 : *|                        propagation de la fissure. Dans ce cas,      |
 206 : *|                        ne pas fournir l'indice 'COUCHE' de SUPTAB,  |
 207 : *|                        mais fournir 'CHAMP_THETA' a chaque appel.   |
 208 : *|                                                                     |
 209 : *|    9 : Cas ou on souhaite calculer une integrale dans l epaisseur   |
 210 : *|        d une structure en coque (rapport DMT/96-317)                |
 211 : *|                                                                     |
 212 : *|        On utilise pour cela la technique de multicouche, qui        |
 213 : *|        consiste, avant d'appeler la proceduer G_THETA, a :          |
 214 : *|        1) Etablir un modele multicouches (cf MODE CONS) sur un ou   |
 215 : *|           des element(s) proche(s) de la fissure sachant qu'il faut |
 216 : *|           au moins une couche en peau inferieure, une couche en     |
 217 : *|           peau superieure, une couche en ligne moyenne {ces couches |
 218 : *|           doivent avoir une epaisseur inferieure a 1e-4*(epaisseur  |
 219 : *|           totale de la coque) et donc 2 couches intermediaires.     |
 220 : *|        2) Penser a donner un excentrement et un nom constituant     |
 221 : *|           different a ces couches.                                  |
 222 : *|        3) Assembler le modele multicouches avec le modele du reste  |
 223 : *|           de la structure.                                          |
 224 : *|        4) Effectuer le calcul des contraintes et des deplacements   |
 225 : *|           avec le modele total et le materiau qui en decoule.       |
 226 : *|        Le calcul de l'integrale avec la procedure G_THETA sera      |
 227 : *|        realise sur un seul element en multicouche et pour toutes les|
 228 : *|        couches dans cet element qui ont une epaisseur inferieure a  |
 229 : *|        1e-4*(epaisseur totale de la coque). Un tel element doit     |
 230 : *|        etre designe par l'argument suivant :                        |
 231 : *|                                                                     |
 232 : *| SUPTAB.'ELEMENT_MULTICOUCHE' = Objet MAILLAGE comportant UN SEUL    |
 233 : *|                                element modelise en multicouche. Il  |
 234 : *|                                doit etre a l'interieur de la zone   |
 235 : *|                                THETA, c'est a dire dans la zone     |
 236 : *|                                definie par le nombre SUPTAB.'COUCHE'.
 237 : *|                                Il ne doit pas etre trop loin, ni    |
 238 : *|                                trop proche de la pointe de la       |
 239 : *|                                fissure. Theoriquement, l'integrale  |
 240 : *|                                a calculer est independant du choix  |
 241 : *|                                de l'element pres de la fissure, ce  |
 242 : *|                                qui est numeriquement verifiable en  |
 243 : *|                                la determinant sur des elemens en    |
 244 : *|                                multicouche differents. NOTA : Cette |
 245 : *|                                technique necessite un maillage tres |
 246 : *|                                fin dans la zone de la pointe de la  |
 247 : *|                                fissure.                             |
 248 : *|                                                                     |
 249 : *|                                                                     |
 250 : *|    SORTIE :                                                         |
 251 : *|    ========                                                         |
 252 : *|                                                                     |
 253 : *| Les resultats du calcul correspondant a un champ THETA specifie par |
 254 : *| l'objet SUPTAB.'COUCHE' (ou SUPTAB.'CHAMP_THETA' dans le cas ou on  |
 255 : *| souhaite donner soi-meme un champ de type Theta) sont sauves de la  |
 256 : *| maniere suivante :                                                  |
 257 : *|                                                                     |
 258 : *|                                                                     |
 259 : *|    Dans tous les cas de calcul                                      |
 260 : *|    ---------------------------                                      |
 261 : *|                                                                     |
 262 : *| SUPTAB.'RESULTATS' = Objet contenant la valeur numerique du calcul. |
 263 : *|                      Son type est variable selon qu'on est en 2D ou |
 264 : *|                      3D et selon la solution du probleme traite :   |
 265 : *|                                                                     |
 266 : *|                      1) valeur de l'integrale de contour dans le cas|
 267 : *|                         d'une solution provenant de l'operateur RESO|
 268 : *|                         2D        => FLOTTANT                       |
 269 : *|                         3D massif => TABLE indicee par              |
 270 : *|                            .(points au front de fissure)            |
 271 : *|                            .'GLOBAL' pour une estimation globale    |
 272 : *|                         3D coque  => TABLE indicee par mots         |
 273 : *|                            .'SUPERI' en peau superieure             |
 274 : *|                            .'INFERI' en peau inferieure             |
 275 : *|                            .'MEDIAN' au plan median et              |
 276 : *|                            .'GLOBAL' pour une estimation globale    |
 277 : *|                                                                     |
 278 : *|                      2) valeur de l'integrale de contour a un       |
 279 : *|                         certain pas du calcul dans le cas d'une     |
 280 : *|                         solution provenant de la procedure PASAPAS  |
 281 : *|                         2D        => TABLE indicee par              |
 282 : *|                            .(numero du pas de calcul)               |
 283 : *|                         3D massif => TABLE indicees par             |
 284 : *|                            .(numero du pas de calcul).(points au    |
 285 : *|                              front de fissure)                      |
 286 : *|                         3D coque  => TABLE indicees                 |
 287 : *|                            .(numero du pas de calcul).'SUPERI'      |
 288 : *|                            .(numero du pas de calcul).'INFERI'      |
 289 : *|                            .(numero du pas de calcul).'MEDIAN' et   |
 290 : *|                            .(numero du pas de calcul).'GLOBAL'      |
 291 : *|                                                                     |
 292 : *|                      3) valeur des F.I.C. (facteurs d'intensite des |
 293 : *|                         contraintes) dans le cas de decouplage des  |
 294 : *|                         modes avec une solution provenant de        |
 295 : *|                         l'operateur RESO                            |
 296 : *|                         2D        => TABLE indicee par mots         |
 297 : *|                            .'I'  pour KI                            |
 298 : *|                            .'II' pour KII                           |
 299 : *|                         3D massif => TABLE indicees par             |
 300 : *|                            .'I'  .(points au front de fissure)      |
 301 : *|                              pour KI                                |
 302 : *|                            .'II' .(points au front de fissure)      |
 303 : *|                              pour KII                               |
 304 : *|                            .'III'.(points au front de fissure)      |
 305 : *|                              pour KIII et                           |
 306 : *|                            .'GLOBAL'.(points au front de fissure)   |
 307 : *|                                                                     |
 308 : *|                      4) valeur des F.I.C. (facteurs d'intensite des |
 309 : *|                         contraintes) a un certain pas du calcul     |
 310 : *|                         dans le cas de decouplage des modes avec    |
 311 : *|                         une solution provenant de la procedure      |
 312 : *|                         PASAPAS                                     |
 313 : *|                         2D        => TABLE indicees                 |
 314 : *|                            .'I' .(numero du pas de calcul) pour KI  |
 315 : *|                            .'II'.(numero du pas de calcul) pour KII |
 316 : *|                         3D massif => TABLE indicees par             |
 317 : *|                            .'I'  .(numero du pas de calcul).(point  |
 318 : *|                              au front de fissure) pour KI           |
 319 : *|                            .'II' .(numero du pas de calcul).(points |
 320 : *|                              au front de fissure) pour KII          |
 321 : *|                            .'III'.(numero du pas de calcul).(points |
 322 : *|                              au front de fissure) pour KIII         |
 323 : *|                                                                     |
 324 : *|                                                                     |
 325 : *|    Dans le cas de calcul effectue pas a pas                         |
 326 : *|    ----------------------------------------                         |
 327 : *|                                                                     |
 328 : *| SUPTAB.'EVOLUTION_RESULTATS' = Objet contenant l'evolution des      |
 329 : *|                                resultats en fonction du temps.      |
 330 : *|                                Son type est variable selon la       |
 331 : *|                                configuration du probleme traite :   |
 332 : *|                                                                     |
 333 : *|                             1) Evolution de l'integrale de contour  |
 334 : *|                                2D        => EVOLUTION               |
 335 : *|                                3D massif => TABLE indicee par       |
 336 : *|                                   .(points au front de fissure)     |
 337 : *|                                   .'GLOBAL' evolution pour une      |
 338 : *|                                     estimation globale              |
 339 : *|                                3D coque  => TABLE indicee par mots  |
 340 : *|                                   .'SUPERI' en peau superieure      |
 341 : *|                                   .'INFERI' en peau inferieure      |
 342 : *|                                   .'MEDIAN' au plan median et       |
 343 : *|                                   .'GLOBAL' evolution pour une      |
 344 : *|                                     estimation globale              |
 345 : *|                                                                     |
 346 : *|                             2) Evolution des F.I.C. (facteurs       |
 347 : *|                                d'intensite de contrainte)           |
 348 : *|                                2D        => TABLE indicee par       |
 349 : *|                                   .'I'  pour KI                     |
 350 : *|                                   .'II' pour KII                    |
 351 : *|                                3D massif => TABLE indicee par       |
 352 : *|                                   .'I'.  (points au front de fissure)
 353 : *|                                   .'II'. (points au front de fissure)
 354 : *|                                   .'III'.(points au front de fissure)
 355 : *|                                   .'GLOBAL' evolution pour une      |
 356 : *|                                     estimation globale              |
 357 : *|                                                                     |
 358 : *|                                                                     |
 359 : *|    Dans le cas des elements de coque                                |
 360 : *|    ---------------------------------                                |
 361 : *|                                                                     |
 362 : *| SUPTAB.'EPAISSEUR_RESULTATS' = representant l'evolution de la valeur|
 363 : *|                                des integrales dans l'epaisseur de la|
 364 : *|                                coque. Son type est variable selon la|
 365 : *|                                solution du probleme traite :        |
 366 : *|                             1) EVOLUTION dans le cas d'une solution |
 367 : *|                                provenant de l'operateur RESO        |
 368 : *|                             2) TABLE indicee par .(numero du pas de |
 369 : *|                                calcul) dans le cas d'une solution   |
 370 : *|                                provenant de la procedure PASAPAS    |
 371 : *|                                                                     |
 372 : *|                                                                     |
 373 : *|    Dans le cas de calcul elasto-plastique                           |
 374 : *|    --------------------------------------                           |
 375 : *|                                                                     |
 376 : *| SUPTAB.'CRITERE_DECHARGE' = En cas de calcul élasto-plastique       |
 377 : *|         isotrope ou cinématique, éventuellement thermique, on       |
 378 : *|         calcul un critère de décharge des contraintes défini par    |
 379 : *|         (si, F = courbe de traction ): crit = F(EPSeq)/ SIGeq.      |
 380 : *|         crit = 1. si non-décharge et crit > 1. si décharge.         |
 381 : *|         SUPTAB.'CRITERE_DECHARGE' est une table indicée par les     |
 382 : *|         temps de calcul.                                            |
 383 : *|                                                                     |
 384 : *|                                                                     |
 385 : *|    Dans le cas du contact frottement sur les levres                 |
 386 : *|    --------------------------------------                           |
 387 : *|                                                                     |
 388 : *| A ce jour, cela est traité pour le cas xfem (thèse de B.Trolle).    |
 389 : *|                                                                     |
 390 : *|=====================================================================|
 391 : fltrac = faux ;
 392 : * fltrac = VRAI ;
 393 : flmess = VRAI ;
 394 : si(flmess);
 395 :   SAUT LIGN;  mess '------------------'
 396 :   'DEBUT DE LA PROCEDURE G_THETA'   '--------------------';
 397 : finsi;
 398 : **************************************************
 399 : ************* INFORMATIONS GENERALES *************
 400 : **************************************************
 401 : SAUT 1 'LIGNE'; VALPI = 3.14159261626;
 402 : &DIME = VALE DIME; &MODE = VALE 'MODE';
 403 : &ELEM = VALE 'ELEM'; MOTAX = MOT 'AXIS' ;
 404 : CONFIG0 = FORM;
 405 : 
 406 : 
 407 : **************************************************
 408 : ***  QUELQUES MOTS POUR SIMPLIFIER L'ECRITURE  ***
 409 : **************************************************
 410 : MTS1 = MOTS 'SCAL';
 411 : MU1 = MOT 'UX'; MU2 = MOT 'UY'; MU3 = MOT 'UZ';
 412 : MF1 = MOT 'FX'; MF2 = MOT 'FY'; MF3 = MOT 'FZ';
 413 : SI (EGA MOTAX &MODE) ;
 414 :    MU1 = MOT 'UR'; MU2 = MOT 'UZ'; MU3 = MOT 'UT';
 415 :    MF1 = MOT 'FR'; MF2 = MOT 'FZ';
 416 : FINSI;
 417 : *listmots
 418 : SI (EGA &DIME 2); 
 419 : * cas axis : faut-il ajouter mu3 ?
 420 :   MU123 = mots MU1 MU2;
 421 :   MF123 = mots MF1 MF2;
 422 :   MV123 = mots 'VX' 'VY';
 423 : SINO; 
 424 :   MU123 = mots MU1 MU2 MU3;
 425 :   MF123 = mots MF1 MF2 MF3;
 426 :   MV123 = mots 'VX' 'VY' 'VZ';
 427 : FINS;
 428 : 
 429 : **************************************************
 430 : ***** DONNEES OBLIGATOIRES DANS TOUS LES CAS *****
 431 : **************************************************
 432 : SI (NON (EXIS SUPTAB 'OBJECTIF'));
 433 :    MESS 'ERREUR :IL FAUT SPECIFIER L INTEGRALE';
 434 :    MESS '        A CALCULER PAR UN MOT';
 435 :    erre 21; QUIT G_THETA;
 436 : SINON;
 437 :     IINTE = 0;
 438 :    SI (EGA SUPTAB.'OBJECTIF' 'J');
 439 :       IINTE = 1;
 440 :    FINSI;
 441 :    SI (EGA SUPTAB.'OBJECTIF' 'C*');
 442 :       IINTE = 2;
 443 :    FINSI;
 444 :    SI (EGA SUPTAB.'OBJECTIF' 'C*H');
 445 :       IINTE = 3;
 446 :    FINSI;
 447 :    SI (EGA SUPTAB.'OBJECTIF' 'DJ/DA');
 448 :       IINTE = 4;
 449 :    FINSI;
 450 :    SI (EGA SUPTAB.'OBJECTIF' 'J_DYNA');
 451 :       IINTE = 5;
 452 :    FINSI;
 453 :    SI (EGA SUPTAB.'OBJECTIF' 'DECOUPLAGE');
 454 :       IINTE = 99;
 455 :    FINSI;
 456 :    SI (EGA IINTE 0);
 457 :       MESS 'ERREUR : ON NE CONNAIT PAS L INTEGRALE SPECIFIEE';
 458 :       MESS '         A CALCULER';
 459 :       erre 21; QUIT G_THETA;
 460 :    FINSI;
 461 : FINSI;
 462 : SI (NON (EXIS SUPTAB 'FRONT_FISSURE'));
 463 :    MESS 'ERREUR : ON VEUT LE FRONT DE LA FISSURE';
 464 :    erre 21; QUIT G_THETA;
 465 : FINSI;
 466 : MESHFR1 = SUPTAB . 'FRONT_FISSURE';
 467 : SI(EGA (TYPE MESHFR1) 'POINT');
 468 :    MESHFR1 = MANU 'POI1' MESHFR1;
 469 : FINSI;
 470 : MESHFR1 = MESHFR1 COUL 'OLIV';
 471 : 
 472 : **************************************************
 473 : ****** TERMES CROISES DE LA MATRICE dJi/daj ******
 474 : **************************************************
 475 : SI (EGA IINTE 4);
 476 :   SI ((EXIS SUPTAB 'FISSURE_2') OU
 477 :         (EXIS SUPTAB 'FRONT_FISSURE_2'));
 478 :     SI (NON (EXIS SUPTAB 'FISSURE_2'));
 479 :       MESS 'ERREUR : ON VEUT AUSSI LA FISSURE 2 POUR CALCULER';
 480 :       MESS '         LES TERMES CROISES DE LA MATRICE';
 481 :       erre 21; QUIT G_THETA;
 482 :     FINSI;
 483 :     SI (NON (EXIS SUPTAB 'FRONT_FISSURE_2'));
 484 :       MESS 'ERREUR : ON VEUT AUSSI LE FROND DE LA FISSURE 2 POUR';
 485 :       MESS '         CALCULER LES TERMES CROISES DE LA MATRICE';
 486 :       erre 21; QUIT G_THETA;
 487 :     FINSI;
 488 :   SINON;
 489 :     SI (EGA SUPTAB.'COUCHE' 0);
 490 :       MESS 'ERREUR : LE NOMBRE DE COUCHES DOIT ETRE SUPERIEUR A';
 491 :       MESS '         0 POUR LE CALCUL DU TERME PRINCIPAL DJi/DAi';
 492 :       erre 21; QUIT G_THETA;
 493 :     FINSI;
 494 :   FINSI;
 495 : FINSI;
 496 : 
 497 : **************************************************
 498 : ****** DONNEES EN CAS DE CALCUL NONLINEAIRE ******
 499 : **************************************************
 500 : *
 501 : *initialisation des valeurs par defaut************
 502 : IPAP  = mot 'NONDEFINI';
 503 : IGDEP = FAUX;
 504 : IGDER = FAUX;
 505 : IPERSO1 = FAUX;
 506 : 
 507 : SI (EXIS SUPTAB 'SOLUTION_PASAPAS');
 508 :   IPAP = VRAI;
 509 : 
 510 : * recup du model et du materiau depuis WTABLE ****
 511 :   SI (EXIS SUPTAB.'SOLUTION_PASAPAS' 'WTABLE');
 512 :       WTAB=SUPTAB.'SOLUTION_PASAPAS'.'WTABLE';
 513 :       OBJMOD=WTAB.'MOD_MEC';
 514 :       OBJMAT=WTAB.'MAT_MEC';
 515 :     SI (EGA IINTE 5);
 516 :        SI (NON WTAB.'DYNAMIQUE');
 517 :           MESS 'ERREUR : IL FAUT UNE SOLUTION ELASTO-DYNAMIQUE';
 518 :           MESS '         POUR CALCULER LE J DYNAMIQUE.';
 519 :           erre 21; QUIT G_THETA;
 520 :        FINSI;
 521 :     FINSI;
 522 : 
 523 : * recup du model et du materiau depuis SOLUTION_PASAPAS ****
 524 : * rem BP: on ne devrait jamais passer par ici ...
 525 :   SINON;
 526 :     MESS 'Absence de WTABLE !  l execution continue ...';
 527 : *   on reduit le modele et le materiau au seul comportement mecanique
 528 :     OBJMOD = EXTR (SUPTAB.'SOLUTION_PASAPAS'.'MODELE')
 529 :              'FORM' 'MECANIQUE';
 530 :     OBJMAT = REDU (SUPTAB.'SOLUTION_PASAPAS'.'CARACTERISTIQUES')
 531 :                OBJMOD ;
 532 :     SI (EGA IINTE 5);
 533 :        SI (NON SUPTAB.'SOLUTION_PASAPAS'.'DYNAMIQUE');
 534 :           MESS 'ERREUR : IL FAUT UNE SOLUTION ELASTO-DYNAMIQUE';
 535 :           MESS '         POUR CALCULER LE J DYNAMIQUE.';
 536 :           erre 21; QUIT G_THETA;
 537 :        FINSI;
 538 :     FINSI;
 539 :     WTAB= SUPTAB.'SOLUTION_PASAPAS';
 540 :   FINSI;
 541 : 
 542 :   SI WTAB.'GRANDS_DEPLACEMENTS';
 543 :      IGDEP = VRAI;
 544 :   FINSI;
 545 :   SI (EXIS SUPTAB 'ROTATION_RIGIDIFIANTE');
 546 :      IGDER = VRAI;
 547 :   FINSI;
 548 : 
 549 : * cas particulier de perso1 ou l on ne calcule que     ****
 550 : * le dernier pas de temps (contenu dans la table estim)
 551 :   SI (EXIS SUPTAB 'PERSO1');
 552 :      IPERSO1 = SUPTAB . 'PERSO1';
 553 :      si(flmess); mess 'utilisation de PERSO1 en cours de dvlpt'; finsi;
 554 :      SI (EXIS SUPTAB.'SOLUTION_PASAPAS' 'ESTIMATION');
 555 :        ESTIM = SUPTAB . 'SOLUTION_PASAPAS' . 'ESTIMATION';
 556 :      SINON;
 557 :        MESS 'ERREUR : il faut une ESTIMATION dans la SOLUTION_PASAPAS';
 558 :        erre 21; QUIT G_THETA;
 559 :      FINSI;
 560 :      si(non (exis SUPTAB 'MAILLAGE_REDUIT'));
 561 :        mess 'Attention! utilisation de PERSO1 sans MAILLAGE_REDUIT';
 562 :        mess 'uniquement valable dans le cas de fissure stationnaire';
 563 :      finsi;
 564 :   FINSI;
 565 : 
 566 : FINSI;
 567 : 
 568 : 
 569 : **************************************************
 570 : ******** DONNEES EN CAS DE CALCUL LINEAIRE *******
 571 : **************************************************
 572 : SI (EXIS SUPTAB 'SOLUTION_RESO');
 573 :   SI (EGA (TYPE SUPTAB.'SOLUTION_RESO') 'CHPOINT ');
 574 :     IPAP = FAUX;
 575 :     SI ((EGA IINTE 2) OU (EGA IINTE 3));
 576 :        MESS 'ERREUR : C* OU C*H EST UNE INTEGRALE';
 577 :        MESS '         CARACTERISTIQUE EN FLUAGE';
 578 :        erre 21; QUIT G_THETA;
 579 :     FINSI;
 580 :     SI (EGA IINTE 5);
 581 :        MESS 'ERREUR : IL FAUT UNE SOLUTION DE LA PROCEDURE PASAPAS';
 582 :        MESS '         POUR CALCULER LE J EN ELASTO-DYNAMIQUE.';
 583 :        erre 21; QUIT G_THETA;
 584 :     FINSI;
 585 :     SI (NON (EXIS SUPTAB 'CARACTERISTIQUES'));
 586 :        MESS 'ERREUR : IL FAUT DONNER LE CHAMP CARACTERISTIQUE';
 587 :        erre 21; QUIT G_THETA;
 588 :     FINSI;
 589 :     SI (NON (EXIS SUPTAB 'MODELE'));
 590 :        MESS 'ERREUR : IL FAUT DONNER LE MODELE DE CALCUL';
 591 :        erre 21; QUIT G_THETA;
 592 :     FINSI;
 593 :     SI ((NON (EXIS SUPTAB 'TEMPERATURES')) ET
 594 :           (NON (EXIS SUPTAB 'CHARGEMENTS_MECANIQUES')));
 595 :        MESS 'ERREUR : IL FAUT LES CHARGEMENTS APPLIQUES :';
 596 :        MESS '         MECANIQUES, THERMIQUES OU LES DEUX';
 597 :        erre 21; QUIT G_THETA;
 598 :     FINSI;
 599 :     SI ((EGA IINTE 4) ET
 600 :           (NON (EXIS SUPTAB 'BLOCAGES_MECANIQUES')));
 601 :        MESS 'ERREUR : IL FAUT DONNER LE BLOCAGES MECANIQUES';
 602 :        erre 21; QUIT G_THETA;
 603 :     FINSI;
 604 :      OBJMOD = SUPTAB.'MODELE';
 605 :      OBJMAT = SUPTAB.'CARACTERISTIQUES';
 606 :   FINSI;
 607 : FINSI;
 608 : SI (EGA (TYPE IPAP) 'MOT     ');
 609 :    MESS 'ERREUR : IL FAUT UNE SOLUTION PROVENANT DE PASAPAS';
 610 :    MESS '         OU DE RESO POUR DETERMINER L INTEGRALE';
 611 :    erre 641; QUIT G_THETA;
 612 : FINSI;
 613 : 
 614 : **************************************************
 615 : ***** DONNEES EN CAS DE CHARGEMENT THERMIQUE *****
 616 : **************************************************
 617 : SI IPAP;
 618 :    CHAR1 = SUPTAB.'SOLUTION_PASAPAS'.'CHARGEMENT';
 619 :    ITHER = EXIS CHAR1 'T   ';
 620 :   SI ITHER;
 621 :      TALPH1 = WTAB.'TALPHA_REFERENCE';
 622 :   FINSI;
 623 : SINON;
 624 :    ITHER = EXIS SUPTAB 'TEMPERATURES';
 625 : FINSI;
 626 : 
 627 : **************************************************************
 628 : ***** DONNEES EN CAS DE CHARGEMENT DEFORMATIONS IMPOSEES *****
 629 : **************************************************************
 630 : SI IPAP;
 631 :   IDEFI = EXIS CHAR1 'DEFI';
 632 : SINON;
 633 :   IDEFI = EXIS SUPTAB 'DEFORMATIONS_IMPOSEES'; 
 634 : FINSI;
 635 : 
 636 : **************************************************
 637 : ***** DONNEES EN CAS DE CONTACT (ajout BP BT) ****
 638 : **************************************************
 639 : IFROT=faux;
 640 : SI IPAP;
 641 : * todo : pas developpe pour l'instant
 642 : SINON;
 643 :   SI (exis SUPTAB 'MODELE_FISSURE');
 644 :     IFROT = vrai;
 645 :     OBJCON = SUPTAB  . 'MODELE_FISSURE';
 646 :     MAICON = extr OBJCON 'MAILLAGE';
 647 :   FINS;
 648 : FINSI;
 649 :   
 650 : **************************************************
 651 : ******* TYPE DES ELEMENTS : COQUE OU MASSIF ******
 652 : **************************************************
 653 : * IPLAN = (EGA &ELEM 'TRI3') OU (EGA &ELEM 'QUA4') OU
 654 : *         (EGA &ELEM 'TRI6') OU (EGA &ELEM 'QUA8');
 655 : *bp: pas tres robuste => on remplace par :
 656 : MAILLAGE = EXTR OBJMOD 'MAIL' ;
 657 : LELEM1 = MAILLAGE ELEM 'TYPE' ;
 658 : IPLAN = (EXIS LELEM1 'TRI3') OU (EXIS LELEM1 'QUA4') OU
 659 :         (EXIS LELEM1 'TRI6') OU (EXIS LELEM1 'QUA8');
 660 : ICOQU = (&DIME EGA 3) ET IPLAN;
 661 : 
 662 : **************************************************
 663 : **** MODELE MULTICOUCHES DANS LE CAS DE COQUE ****
 664 : **************************************************
 665 : SI ICOQU;
 666 :     M_DETA = EXTR OBJMOD 'ZONE';
 667 :    SI (EXIS SUPTAB 'ELEMENT_MULTICOUCHE');
 668 :       ELMULT = SUPTAB.'ELEMENT_MULTICOUCHE';
 669 :      SI (NEG (TYPE ELMULT) 'MAILLAGE');
 670 :         MESS 'ERREUR : L ELEMENT EN MULTICOUCHE DOIT';
 671 :         MESS '         ETRE UN OBJET DE TYPE MAILLAGE';
 672 :         erre 21; QUIT G_THETA;
 673 :      FINSI;
 674 :      SI (NEG (NBEL ELMULT) 1);
 675 :         MESS 'ERREUR : ON VEUT UN SEUL ELEMENT EN MULTICOUCHE';
 676 :         erre 21; QUIT G_THETA;
 677 :      FINSI;
 678 :    SINON;
 679 :      MESS 'ERREUR : IL FAUT DESIGNER UN ELEMENT EN MULTICOUCHE';
 680 :      erre 21; QUIT G_THETA;
 681 :    FINSI;
 682 :     M_ELMU = EXTR (REDU OBJMOD ELMULT) 'ZONE';
 683 :    SI ((DIME M_ELMU) '<' 10);
 684 :      MESS 'ERREUR : IL FAUT AU MOINS 3 COUCHES (peau inf, peau';
 685 :      MESS '         sup, ligne moyenne) D EPAISSEUR INFERIEURE';
 686 :      MESS '         A 1E-4*(EPAISSEUR TOTALE) + 2 COUCHES';
 687 :      MESS '         INTERMEDIAIRES POUR L ELEMENT DESIGNE EN';
 688 :      MESS '         MULTICOUCHE PROCHE DE LA FISSURE.';
 689 :      erre 21; QUIT G_THETA;
 690 :    FINSI;
 691 :     PEX1 = PROG; LMO1 = LECT; MODCOU = TABLE; EPAITO = 0.;
 692 :    REPETER NBJ5 ((DIME M_ELMU)/2);
 693 :       I1 = (2 * &NBJ5) - 1;
 694 :       MODCOU.&NBJ5 = M_ELMU.I1;
 695 :       EX1 = EXTR (REDU MODCOU.&NBJ5 OBJMAT) 'EXCE' 1 1 1;
 696 :       EP1 = EXTR (REDU MODCOU.&NBJ5 OBJMAT) 'EPAI' 1 1 1;
 697 :       EPAITO = EPAITO + EP1;
 698 :       PEX1 = PEX1 ET (PROG EX1);
 699 :       LMO1 = LMO1 ET (LECT &NBJ5);
 700 :    FIN NBJ5;
 701 :     NSUPE = 0; NMOYE = 0; NINFE = 0;
 702 :    REPETER NBJ6 (DIME MODCOU);
 703 :       EX1 = EXTR PEX1 &NBJ6;
 704 :       LM1 = EXTR LMO1 &NBJ6;
 705 :      SI (EGA EX1 (EPAITO/2.) 1.E-4); NSUPE = LM1; FINSI;
 706 :      SI (EGA EX1 0. 1.E-10); NMOYE = LM1; FINSI;
 707 :      SI (EGA EX1 (EPAITO/(-2.)) 1.E-4); NINFE = LM1; FINSI;
 708 :    FIN NBJ6;
 709 :    SI (EGA NSUPE 0);
 710 :      MESS 'ERREUR : IL FAUT UNE COUCHE EN PEAU SUPERIEURE';
 711 :      MESS '         D EPAISSEUR INFERIEURE A';
 712 :      MESS '         1E-4*(EPAISSEUR TOTALE) ';
 713 :      erre 21; QUIT G_THETA;
 714 :    FINSI;
 715 :    SI (EGA NMOYE 0);
 716 :      MESS 'ERREUR : IL FAUT UNE COUCHE AU PLAN MEDIAN';
 717 :      MESS '         AYANT UN EXCENTREMENT NUL';
 718 :      erre 21; QUIT G_THETA;
 719 :    FINSI;
 720 :    SI (EGA NINFE 0);
 721 :      MESS 'ERREUR : IL FAUT UNE COUCHE EN PEAU INFERIEURE';
 722 :      MESS '         D EPAISSEUR INFERIEURE A';
 723 :      MESS '         1E-4*(EPAISSEUR TOTALE) ';
 724 :      erre 21; QUIT G_THETA;
 725 :    FINSI;
 726 :     SUPTAB.'EPAISSEUR' = EPAITO;
 727 :     M_SUPE = MODCOU.NSUPE;
 728 :     M_MOYE = MODCOU.NMOYE;
 729 :     M_INFE = MODCOU.NINFE;
 730 : FINSI;
 731 : 
 732 : **************************************************
 733 : ** MAILLAGE UTILISE DANS LA RESOLUTION PAR E.F. **
 734 : **************************************************
 735 : * MAILLAGE = EXTR OBJMOD 'MAIL' ;
 736 : * -> bp:fait + haut
 737 : SI ICOQU;
 738 :    TMULT = TABLE;
 739 :   REPETER NBJ8 ((DIME M_DETA)/2);
 740 :      M1 = M_DETA.(2*&NBJ8);
 741 :     SI (EXIS TMULT M1);
 742 :        ITER NBJ8;
 743 :     FINSI;
 744 :      M2 = EXTR (REDU OBJMOD M1) 'ZONE';
 745 :     SI ('>' (DIME M2) 2);
 746 :        TMULT.M1 = VRAI;
 747 :       REPETER NBJ9 (((DIME M2)/2) - 1);
 748 :          MAILLAGE = DIFF MAILLAGE M1;
 749 :       FIN NBJ9;
 750 :     FINSI;
 751 :   FIN NBJ8;
 752 : FINSI;
 753 : SUPTAB.'MAILLAGE' = MAILLAGE;
 754 : 
 755 : **************************************************
 756 : *************** TYPES D ELEMENTS  ****************
 757 : **************************************************
 758 : 
 759 : ******* ELEMENTS LINEAIRES OU NONLINEAIRES *******
 760 : NBNO1 = NBNO (ELEM (CHAN 'LIGNE' MAILLAGE) 1);
 761 : ILIN = EGA NBNO1 2; IQUA = EGA NBNO1 3;
 762 : 
 763 : ******* ELEMENTS XFEM OU STANDARD ****************
 764 : IXFEM = EXIS OBJMOD 'ELEM' 'XQ4R' 'XC8R';
 765 : 
 766 : **************************************************
 767 : ************* DEFINITON DE LA FISSURE ************
 768 : **************************************************
 769 : 
 770 : ******* LEVELSET PSI ET PHI POUR XFEM ************
 771 : SI(IXFEM);
 772 :   SI ((EXIS SUPTAB 'PSI') et (EXIS SUPTAB 'PHI'));
 773 :     PSI0 = SUPTAB . 'PSI';
 774 :     PHI0 = SUPTAB . 'PHI';
 775 :   SINO;
 776 :     MESS 'ERREUR : ON VEUT PSI et PHI LEVELSET DE LA FISSURE';
 777 :     erre 641; QUIT G_THETA;
 778 :   FINSI;
 779 : 
 780 : ******* LEVRE_SUPERIEURE ET INFERIEURE POUR STD ***
 781 : SINO;
 782 :   SI (NON (EXIS SUPTAB 'LEVRE_SUPERIEURE'));
 783 :      SI (EGA IINTE 99);
 784 :         MESS 'ERREUR : ON VEUT LA LEVRE SUPERIEURE DE LA FISSURE';
 785 :         erre 641; QUIT G_THETA;
 786 :      FINSI;
 787 :      SI (NON (EXIS SUPTAB 'LEVRE_INFERIEURE'));
 788 :         MESS 'ERREUR : IL FAUT DONNER LA FISSURE';
 789 :         MESS '(LEVRE_SUPERIEURE ou LEVRE_INFERIEURE ou les 2)'
 790 :         erre 641; QUIT G_THETA;
 791 :      SINON;
 792 :          SUPTAB.'FISSURE' = SUPTAB.'LEVRE_INFERIEURE';
 793 :      FINSI;
 794 :   SINON;
 795 :      SI (NON (EXIS SUPTAB 'LEVRE_INFERIEURE'));
 796 :         SI (EGA IINTE 99);
 797 :            MESS 'ERREUR : ON VEUT LA LEVRE INFERIEURE DE LA FISSURE';
 798 :            erre 641; QUIT G_THETA;
 799 :         FINSI;
 800 :          SUPTAB.'FISSURE' = SUPTAB.'LEVRE_SUPERIEURE';
 801 :      SINON;
 802 :          SUPTAB.'FISSURE' = (SUPTAB.'LEVRE_SUPERIEURE') ET
 803 :                             (SUPTAB.'LEVRE_INFERIEURE');
 804 :      FINSI;
 805 :   FINSI;
 806 :   si (IFROT);
 807 :     MESS 'ERREUR : CONTACT via MODELE_FISSURE avec XFEM seulement';
 808 :     erre 641; QUIT G_THETA;  
 809 :   fins;
 810 : FINSI;
 811 : 
 812 : ****************************************************
 813 : ******* DETERMINATION DES CHAMPS THETA ET PI *******
 814 : ******* ET DE LA ZONE DE TRAVAIL ELTETA      *******
 815 : ****************************************************
 816 : si(flmess);   mess 'DETERMINATION DES CHAMPS THETA';  fins;
 817 : 
 818 : *** nombre de COUCHE donné => on calcule tout le reste
 819 : * CHAMP_THETA + DIRTETA
 820 : SI (EXIS SUPTAB 'COUCHE');
 821 :   SI(IXFEM);
 822 :     SUPTAB.'CHAMP_THETA' UTILTETA = CH_THETX SUPTAB;
 823 :     SUPTAB.'UTILTET1' = UTILTETA;
 824 :   SINO;
 825 :     SUPTAB.'CHAMP_THETA' UTILTETA = CH_THETA SUPTAB;
 826 :     SUPTAB.'UTILTET1' = UTILTETA;
 827 :     SI (EGA IINTE 4);
 828 :       SI (NON (EXIS SUPTAB 'FRONT_FISSURE_2'));
 829 :          SUPTAB.'COUCHE' = (SUPTAB.'COUCHE') - 1;
 830 :          SUPTAB.'CHAMP_PI' UTILPI = CH_THETA SUPTAB;
 831 :          SUPTAB.'COUCHE' = (SUPTAB.'COUCHE') + 1;
 832 :       SINON;
 833 :          P1 = SUPTAB.'FRONT_FISSURE';
 834 :          SUPTAB.'FRONT_FISSURE' = SUPTAB.'FRONT_FISSURE_2';
 835 :          SUPTAB.'FISSURE' = SUPTAB.'FISSURE_2';
 836 :          SUPTAB.'CHAMP_PI' UTILPI = CH_THETA SUPTAB;
 837 :          SUPTAB.'FRONT_FISSURE' = P1;
 838 :       FINSI;
 839 :     FINSI;
 840 :   FINSI;
 841 : * ELTETA = ...
 842 :   si(exis SUPTAB 'MAILLAGE_REDUIT');
 843 : *   ELTETA = MAILLAGE fourni par l utilisateur (attention pas de test de compati
 844 :     ELTETA = SUPTAB . 'MAILLAGE_REDUIT';
 845 :   sino;
 846 : *   ELTETA = MAILLAGE OU TETA N EST PAS NUL + 1 couche
 847 :     SI (EGA (TYPE ( SUPTAB.'CHAMP_THETA')) 'CHPOINT  ');
 848 :       uu = EXTR  SUPTAB.'CHAMP_THETA' 'MAILLAGE';
 849 :     SINON;
 850 :       uu= EXTR SUPTAB.'CHAMP_THETA'.'GLOBAL'  'MAILLAGE';
 851 :     FINSI;
 852 :     ELTETA = ELEM MAILLAGE 'APPU' 'LARG' UU;
 853 :   fins;
 854 : 
 855 : *** CHAMP_THETA donné, on calcule DIRTETA sur le front de fissure
 856 : *(pour faire simple, on appelle CH_THETA pour cela, mais ce n'est pas economique
 857 : SINON;
 858 :   SI(EXIS SUPTAB 'CHAMP_THETA');
 859 :     MESS 'CHAMP_THETA FOURNI PAR L UTILISATEUR';
 860 :     q7 = SUPTAB.'CHAMP_THETA';
 861 :     SI (NEG (TYPE q7) 'CHPOINT  ');
 862 :       q7= q7 . 'GLOBAL' ;
 863 :     FINSI;
 864 : *   ELTETA = ...
 865 :     si(exis SUPTAB 'MAILLAGE_REDUIT');
 866 : *     ELTETA = MAILLAGE fourni par l utilisateur (attention pas de test de compa
 867 :       ELTETA = SUPTAB . 'MAILLAGE_REDUIT';
 868 :     sino;
 869 : *     ELTETA = MAILLAGE OU TETA N EST PAS NUL + 1 couche
 870 :       uu = EXTR  q7 'MAILLAGE';
 871 :       ELTETA = ELEM MAILLAGE 'APPU' 'LARG' UU;
 872 :     fins;
 873 : *   UTILTETA = ...
 874 :     SI(EXIS  SUPTAB 'UTILTETA');
 875 :       MESS 'UTILTETA FOURNI PAR L UTILISATEUR';
 876 :       UTILTETA = SUPTAB . 'UTILTETA';
 877 :       VECTEUR1 = UTILTETA . 'DIRECTION1';
 878 :       VECTEUR2 = UTILTETA . 'DIRECTION2';
 879 :       SI(EGA &DIME 3); VECTEUR3 = UTILTETA . 'DIRECTION3'; FINSI;
 880 :     SINON;
 881 :       UTILTETA = TABL;
 882 : *     DIRECTIONS dans TABUTIL
 883 :       VECTEUR1 = INT_COMP ELTETA q7 MESHFR1;
 884 :       NV1 = PSCA VECTEUR1 VECTEUR1 MU123 MU123;
 885 :       VECTEUR1 = (VECTEUR1 / (NV1**0.5)) 
 886 :       CHAN 'ATTRIBUT' 'NATURE' 'DIFFUS';
 887 :       UTILTETA . 'DIRECTION1' = VECTEUR1;
 888 :       SI(EGA &DIME 2);
 889 :         VECTEUR2 =  (-1.*(EXCO VECTEUR1 'UY' 'UX'))
 890 :                       et (EXCO VECTEUR1 'UX' 'UY');
 891 :         UTILTETA . 'DIRECTION2' =  VECTEUR2;
 892 :       SINO;
 893 : *       CHT CHN CHB  = FRENET  SUPTAB.'FRONT_FISSURE';
 894 : *       VECTEUR1 = -1.*CHN;
 895 : *       VECTEUR2 = -1.*CHB;
 896 : *       VECTEUR3 = -1.*CHT;
 897 :        MESS ' bp: !!! option non testée, mais on est joueur !!!';
 898 :        modfro = MODE MESHFR1 MECANIQUE ELASTIQUE 'POUT';
 899 :        VECTEUR3 = (VSUR modfro 'NORM') EXCO MV123 MU123 ;
 900 :        VECTEUR3 = CHAN 'CHPO' VECTEUR3 'MOYE';
 901 :        NV3 = PSCA VECTEUR3 VECTEUR3 MU123 MU123;
 902 :        VECTEUR3 = (VECTEUR3 / (NV3**0.5)) 
 903 :        CHAN 'ATTRIBUT' 'NATURE' 'DIFFUS';
 904 :        VECTEUR2 = (PVEC VECTEUR3 MU123 VECTEUR1 MU123 MU123)
 905 :        CHAN 'ATTRIBUT' 'NATURE' 'DIFFUS';
 906 :        UTILTETA . 'DIRECTION2' =  VECTEUR2;
 907 :        UTILTETA . 'DIRECTION3' =  VECTEUR3;
 908 :       FINS;
 909 :     FINSI;
 910 : 
 911 :   SINO;
 912 : *** ni COUCHE ni CHAMP_THETA donné, ERREUR !
 913 :     MESS 'ERREUR : ON VEUT LE NOMBRE DE COUCHEs D ELEMENTS';
 914 :     MESS '         AUTOUR DE LA FISSURE QUI SE DEPLACE';
 915 :     MESS '         ou LE CHAMP_THETA';
 916 :     MESS '         POUR SIMULER LA PROPAGATION DE LA FISSURE';
 917 :     erre 641; QUIT G_THETA;
 918 :   FINSI;
 919 : FINSI;
 920 : *
 921 : si(fltrac);
 922 :   q7 = SUPTAB.'CHAMP_THETA';
 923 :   si(neg (type q7) 'CHPOINT'); q7 = q7 . 'GLOBAL'; fins;
 924 :   vq7 = VECT q7 'DEPL' 'BLEU' ;
 925 :   trac vq7 (MAILLAGE et MESHFR1) 'TITR' 'CHAMP_THETA';
 926 : finsi;
 927 : 
 928 : 
 929 : **************************************************
 930 : ************** DIRECTIONS UTILES *****************
 931 : **************************************************
 932 : 
 933 : *** DIRECTION DE PROPAGATION DE LA FISSURE = DIRTETA
 934 : * SI ((EGA &DIME 2) OU ICOQU);
 935 : *    DIRTETA = UTILTETA . 'DIRECTION';
 936 :    DIRTETA = UTILTETA . 'DIRECTION1';
 937 : * FINSI;
 938 : * SI ((EGA &DIME 3) ET (NON ICOQU));
 939 : *   IND1 = INDE (UTILTETA.'DIRECTION');
 940 : *   DIRTETA = 0. 0. 0.;
 941 : *   REPETER BC1 ((DIME IND1) - 1);
 942 : *      DIRTETA = DIRTETA 'PLUS' (UTILTETA.'DIRECTION'.(IND1.&BC1));
 943 : *   FIN BC1;
 944 : * FINSI;
 945 : *   DIRTETA = DIRTETA / (NORM DIRTETA);
 946 : si(non ICOQU);
 947 :    DIRNORM = UTILTETA . 'DIRECTION2';
 948 : fins;
 949 : 
 950 : *** DIRECTION DE CISAILLEMENT SI SEPARATION DE MODES =DIRCISA
 951 : *SI ((EGA IINTE 99) ET (EGA &DIME 3) ET (NON ICOQU));
 952 : SI ((EGA &DIME 3) ET (NON ICOQU));
 953 : *   SI(IXFEM);
 954 :      DIRCISA = UTILTETA . 'DIRECTION3';
 955 : *   SINON;
 956 : *      F1 = PRES 'MASS' OBJMOD SUPTAB.'LEVRE_SUPERIEURE' 1.;
 957 : *      N1 = NBNO SUPTAB.'FRONT_FISSURE';
 958 : *      P1 = POIN SUPTAB.'FRONT_FISSURE' ((N1 + 1)/2);
 959 : *      V1 = EXTR F1 MF1 P1;
 960 : *      V2 = EXTR F1 MF2 P1;
 961 : *      V3 = EXTR F1 MF3 P1;
 962 : *      DIRCISA = PVEC DIRTETA (V1 V2 V3);
 963 : *   FINSI;
 964 : *   DIRCISA = DIRCISA / (NORM DIRCISA);
 965 : FINSI;
 966 : 
 967 : * si(fltrac);
 968 : * *si(vrai);
 969 : *   dx1 = coor (&DIME + 1) MESHFR1;
 970 : * *  dx1 =  maxi (prog 1. ((maxi (resu dx1)) / (nbno MESHFR1)));
 971 : *   dx1 =  (maxi (resu dx1)) / (nbno MESHFR1);
 972 : *   vdir7 = (VECT dx1 DIRTETA 'DEPL' 'BLEU')
 973 : *        et (VECT dx1 DIRNORM 'DEPL' 'ROUG');
 974 : *   si((EGA &DIME 3) ET (NON ICOQU));
 975 : *     vdir7 = vdir7 et (VECT dx1 DIRCISA 'DEPL' 'VERT');
 976 : *     trac vdir7 (MESHFR1 et (aret MAILLAGE)) TITR 'DIRECTIONS LOCALES';
 977 : *   sino;
 978 : *     trac vdir7 (MESHFR1 et (cont MAILLAGE)) TITR 'DIRECTIONS LOCALES';
 979 : *   fins;
 980 : * fins;
 981 : 
 982 : 
 983 : **************************************************
 984 : **************  ON COMPLETE ELTETA  **************
 985 : **************************************************
 986 : 
 987 : *** ajout eventuel de ELPIa ELTETA ******
 988 :  SI ((EXIS SUPTAB 'FRONT_FISSURE_2') ET (EGA IINTE 4));
 989 :     ELPI = SUPTAB.'FRONT_FISSURE_2';
 990 :     REPETER MAIL2 ((SUPTAB.'COUCHE') + 1);
 991 :        ELPI = MAILLAGE ELEM 'APPU' 'LARG' ELPI ;
 992 :     FIN MAIL2 ;
 993 :     ELTETA = ELTETA ET ELPI;
 994 :  FINSI;
 995 : 
 996 : *** AJOUT DU NOEUD SUPPORT EN DEF.PL.GENERALISEES
 997 : SI (EGA &MODE 'PLANGENE');
 998 :    ELTETA = ELTETA ET (VALE 'MODE' 'PLANGENE');
 999 : FINSI;
1000 :   ELPOI1 = CHAN ELTETA 'POI1';
1001 : 
1002 : *** L ELEMENT SUPPORTANT LE MODELE MULTICOUCHE
1003 : *** DOIT ETRE DANS LA ZONE THETA
1004 : SI ICOQU;
1005 :    N1 = NBNO ELTETA;
1006 :    N2 = NBNO (ELTETA ET (EXTR M_MOYE 'MAIL'));
1007 :   SI (NEG N1 N2);
1008 :      MESS 'ERREUR : L ELEMENT EN MULTICOUCHE DESIGNE POUR CALCULER';
1009 :      MESS '         L INTEGRALE SE TROUVE EN DEHORS DE LA ZONE';
1010 :      MESS '         DEFINIE PAR LE NOMBRE DE COUCHES DONNE.';
1011 :      erre 21; QUIT G_THETA;
1012 :   FINSI;
1013 : FINSI;
1014 : 
1015 : 
1016 : **************************************************
1017 : *********** TESTER SI REPRISE DE CALCUL **********
1018 : **************************************************
1019 : 
1020 : *** REPRISE DE CALCUL ? **************************
1021 : IREPRI = FAUX;
1022 : SI (IPAP ET (NON IPERSO1));
1023 :    N1 = DIME (SUPTAB.'SOLUTION_PASAPAS'.'TEMPS');
1024 :   SI ((EXIS SUPTAB 'IABC') ET
1025 : *bp        (EXIS SUPTAB 'COU1') ET
1026 : *bp        (EXIS SUPTAB 'CHAMP_THETA') ET
1027 :         (EXIS SUPTAB 'ELTET1') ET
1028 :         (EXIS SUPTAB 'RESULTATS') ET
1029 :         (EXIS SUPTAB 'EVOLUTION_RESULTATS'));
1030 :      IREPRI = '>' (N1 - 1) SUPTAB.'IABC';
1031 :      mess 'on tente une reprise...';
1032 :   FINSI;
1033 : FINSI;
1034 : 
1035 : *** TESTS DE COMPATIBILITE SI REPRISE DE CALCUL ***
1036 : SI IREPRI;
1037 : * on verifie que l objectif reste le meme
1038 :   SI (NEG SUPTAB.'OBJ1' SUPTAB.'OBJECTIF');
1039 :      MESS 'ERREUR : REPRISE IMPOSSIBLE CAR L OBJECTIF DU';
1040 :      MESS '         CALCUL ACTUEL N EST PAS LE MEME QUE';
1041 :      MESS '         CELUI DU CALCUL PRECEDENT';
1042 :      erre 21; QUIT G_THETA;
1043 :   FINSI;
1044 : * on doit avoir le meme nombre de couche (on suppose la fissure fixe)
1045 :   SI ((EXIS SUPTAB 'COUCHE') et (EXIS SUPTAB 'COU1'));
1046 :     SI (NEG SUPTAB.'COU1' SUPTAB.'COUCHE');
1047 :        MESS 'ERREUR : REPRISE IMPOSSIBLE CAR LE NOMBRE DE';
1048 :        MESS '         COUCHE ACTUEL N EST PAS LE MEME QUE';
1049 :        MESS '         CELUI UTILISE POUR LE CALCUL PRECEDENT';
1050 :        erre 21; QUIT G_THETA;
1051 :     FINSI;
1052 :   FINSI;
1053 : * reste a verifier la compatibilite des support de champ teta via elteta
1054 : * ELTETA doit etre inclus dans ELTET1
1055 :   ELTET1 = SUPTAB.'ELTET1';
1056 :   si(neg (nbno ELTETA)  (nbno (ELTET1 inte ELTETA)));
1057 :       MESS 'ERREUR : REPRISE IMPOSSIBLE CAR LE SUPPORT DU ';
1058 :       MESS '         CHAMP_THETA FOURNI N EST PAS INCLUS DANS';
1059 :       MESS '         CELUI UTILISE POUR LE CALCUL PRECEDENT';
1060 :       erre 21; QUIT G_THETA;
1061 :   fins;
1062 :   MESS 'REPRISE DU CALCUL AUTORISE !';
1063 : FINSI;
1064 : 
1065 : *** TESTS DE COMPATIBILITE SI UTILISATION DE PERSO1 ***
1066 : SI (IPERSO1 et (EXIS SUPTAB 'ELTET1'));
1067 :   ELTET1 = SUPTAB.'ELTET1';
1068 :   si(neg (nbno ELTETA)  (nbno (ELTET1 inte ELTETA)));
1069 :       MESS 'ERREUR : REPRISE IMPOSSIBLE CAR LE SUPPORT DU ';
1070 :       MESS '         CHAMP_THETA FOURNI N EST PAS INCLUS DANS';
1071 :       MESS '         CELUI UTILISE POUR LE CALCUL PRECEDENT';
1072 :       erre 21; QUIT G_THETA;
1073 :   fins;
1074 :   si(flmess); MESS 'POURSUITE DU CALCUL via PERSO1 AUTORISE !'; fins;
1075 : FINS;
1076 : 
1077 : 
1078 : 
1079 : **************************************************
1080 : ** MODELES ET MATERIAUX DANS LA ZONE DE TRAVAIL **
1081 : **************************************************
1082 : 
1083 : *** VERIFICATION DES DONNEES D ENTREE POUR MODELES_COMPOSITES
1084 : SI (EXIS SUPTAB 'MODELES_COMPOSITES');
1085 :   SI ((&DIME EGA 3) ET (NON ICOQU));
1086 :     MESS 'ERREUR : ON NE PEUT ENCORE TRAITER LES PROBLEMES';
1087 :     MESS '         DE MATERIAUX COMPOSITES EN 3D MASSIF';
1088 :     erre 21; QUIT G_THETA;
1089 :   FINSI;
1090 :    N1 = DIME SUPTAB.'MODELES_COMPOSITES';
1091 :   SI ('<' N1 2);
1092 :     MESS 'ERREUR : IL FAUT AU MOINS DEUX MODELES POUR';
1093 :     MESS '         DETERMINER LA LIGNE COMMUNE D INTERFACE';
1094 :     erre 21; QUIT G_THETA;
1095 :   FINSI;
1096 :    M1 = EXTR OBJMOD 'MAIL';
1097 :   REPETER BIN4 N1;
1098 :      T1 = TYPE SUPTAB.'MODELES_COMPOSITES'.&BIN4;
1099 :     SI (NEG T1 'MMODEL  ');
1100 :        MESS 'ERREUR : LE TYPE DE L OBJET No' &BIN4 'DANS LA';
1101 :        MESS '         TABLE MODELES_COMPOSITES EST INCORRECTE';
1102 :        erre 21; QUIT G_THETA;
1103 :     FINSI;
1104 :     SI (EGA &BIN4 1);
1105 :        M2 = EXTR SUPTAB.'MODELES_COMPOSITES'.&BIN4 'MAIL';
1106 :     SINON;
1107 :        M2 = M2 ET
1108 :            (EXTR SUPTAB.'MODELES_COMPOSITES'.&BIN4 'MAIL');
1109 :     FINSI;
1110 :   FIN BIN4;
1111 :   SI (NEG (NBNO M1) (NBNO M2));
1112 :     MESS 'ERREUR : TOUS LES MODELES DE MATERIAUX';
1113 :     MESS '         COMPOSITES NE SONT PAS DONNES';
1114 :     erre 21; QUIT G_THETA;
1115 :   FINSI;
1116 : FINSI;
1117 : 
1118 : *** CREATION DE OBJMOD ET TABMOD *****************
1119 : TABMOD = TABL;
1120 : SI (EXIS SUPTAB 'MODELES_COMPOSITES');
1121 : * CAS DE MODELES COMPOSITES (AVEC DISCONTINUITE) : ON A DU TRAVAIL
1122 :   REPETER BIN1 (DIME SUPTAB.'MODELES_COMPOSITES');
1123 :     M1 = SUPTAB.'MODELES_COMPOSITES'.&BIN1;
1124 :     M2 = EXTR M1 'MAIL';
1125 :     N1 = NBNO M2;
1126 :     N2 = NBNO ELTETA;
1127 :     N3 = NBNO (ELTETA ET M2);
1128 : *   si on a des noeuds en commun, ...
1129 :     SI (NEG (N1 + N2) N3);
1130 :        M2 = CHAN M2 'POI1';
1131 : *      ... on les recupere
1132 :        E1 = (DIFF ELPOI1 M2) DIFF (ELPOI1 ET M2);
1133 :        N1 = NBNO (CONT ELTETA);
1134 :        N2 = NBNO (E1 ET (CONT ELTETA));
1135 : *      si tous les noeuds en commun sont sur le contour
1136 : *      => pas d elements a recuperer => on passe au modele suivant
1137 :        SI (EGA N1 N2);
1138 :           ITER BIN1;
1139 :        FINSI;
1140 : *      sinon, on recupere les elements concernes et le modele reduit
1141 :        E1 = MAILLAGE ELEM 'APPU' 'STRI' E1;
1142 :        N1 = (DIME TABMOD) + 1;
1143 :        TABMOD.N1 = REDU M1 E1;
1144 :        SI (EGA N1 1);
1145 :          OBJMOD = TABMOD.N1;
1146 :        SINON;
1147 :          OBJMOD = OBJMOD ET TABMOD.N1;
1148 :        FINSI;
1149 :     FINSI;
1150 :   FIN BIN1;
1151 : SINON;
1152 : * CAS DE MODELES SANS DISCONTINUITE : ON A MOINS DE TRAVAIL
1153 :    OBJMOD = REDU OBJMOD ELTETA;
1154 :    TABMOD.1 = OBJMOD;
1155 : FINSI;
1156 : * list OBJMOD;
1157 : NBOBJ = DIME TABMOD;
1158 : OBJMAT = REDU OBJMAT OBJMOD;
1159 : 
1160 : *** CHAMP EPAISSEUR DANS LA ZONE DE TRAVAIL ******
1161 : SI ICOQU;
1162 :    EPAICH = (CHAN (EXCO OBJMAT 'EPAI' 'SCAL')
1163 :             'STRESSES' OBJMOD) CHAN 'TYPE' 'SCALAIRE';
1164 : FINSI;
1165 : 
1166 : 
1167 : **************************************************
1168 : ** CALCUL DE C* PAR DEUX TYPES DE MODELE FLUAGE **
1169 : **************************************************
1170 : *** ITYPEF = 1  MODELE FLUAGE POUR LEQUEL ON A UNE EXPRESSION
1171 : ***             EXPLICITE DE L'INTEGRATION DE LA VITESSE DE
1172 : ***             DEFORMATION DE FLUAGE SUR LE TEMPS
1173 : *** ITYPEF = 2  MODELE FLUAGE POUR LEQUEL ON N'OBTIENT PAS
1174 : ***             FACILEMENT CETTE EXPRESSION EXPLICITE
1175 : *** ITYPEF = 99 SI EN ELASTO OU THERMO-ELASTO-PLASTICITE
1176 : 
1177 : *** VERIF DES DONNEES ****************************
1178 : SI ((EGA IINTE 2) OU (EGA IINTE 3));
1179 :   SI (NON (EXIS OBJMOD 'MATE' 'FLUAGE'));
1180 :     MESS 'ERREUR : LA FORMULATION DU PROBLEME NE PERMET';
1181 :     MESS '         PAR DE CALCULER L INTEGRALE SPECIFIEE';
1182 :     erre 21; QUIT G_THETA;
1183 :   FINSI;
1184 : FINSI;
1185 : SI (EGA IINTE 3);
1186 :   SI ((EXIS OBJMOD 'MATE' 'FLUAGE' 'BLACKBURN') OU
1187 :       (EXIS OBJMOD 'MATE' 'FLUAGE' 'RCCMR_316') OU
1188 :       (EXIS OBJMOD 'MATE' 'FLUAGE' 'RCCMR_304') OU
1189 :       (EXIS OBJMOD 'MATE' 'FLUAGE' 'POLYNOMIAL') OU
1190 :       (EXIS OBJMOD 'MATE' 'FLUAGE' 'LEMAITRE'));
1191 :     MESS 'ERREUR : IL FAUT UN MODELE DE FLUAGE NORTON';
1192 :     MESS '         SEUL POUR CALCULER L INTEGRALE C*(H)';
1193 :     erre 21; QUIT G_THETA;
1194 :   FINSI;
1195 : FINSI;
1196 : 
1197 : *** DETERMINATION DE ITYPEF *********************
1198 : ITYPEF = 99;
1199 : SI ((EGA IINTE 2) OU (EGA IINTE 3));
1200 :   SI ((EXIS OBJMOD 'MATE' 'FLUAGE' 'NORTON') OU
1201 :       (EXIS OBJMOD 'MATE' 'FLUAGE' 'POLYNOMIAL'));
1202 :      ITYPEF = 1;
1203 :   FINSI;
1204 :   SI ((EXIS OBJMOD 'MATE' 'FLUAGE' 'BLACKBURN') OU
1205 :       (EXIS OBJMOD 'MATE' 'FLUAGE' 'RCCMR_316') OU
1206 :       (EXIS OBJMOD 'MATE' 'FLUAGE' 'RCCMR_304') OU
1207 :       (EXIS OBJMOD 'MATE' 'FLUAGE' 'LEMAITRE'));
1208 :      ITYPEF = 2;
1209 :   FINSI;
1210 : FINSI;
1211 : 
1212 : 
1213 : **************************************************
1214 : ******* INTERFACES DANS LA ZONE DE TRAVAIL *******
1215 : **************************************************
1216 : 
1217 : *** CREATION DES INTERFACES INTER-MODELE *********
1218 : IPARAL = VRAI; LINTER = TABL;
1219 : SI ((EXIS SUPTAB 'MODELES_COMPOSITES')
1220 :        ET ('>' (DIME TABMOD) 1));
1221 : * on boucle sur les modeles qui appartiennent a ELTETA
1222 :   REPETER BIN2 ((DIME TABMOD) - 1);
1223 :      M1 = EXTR (TABMOD.&BIN2) 'MAIL';
1224 :      IIN3 = &BIN2;
1225 :      NIN3 = (DIME TABMOD) - &BIN2;
1226 :      REPETER BIN3 NIN3;
1227 :        IIN3 = IIN3 + 1;
1228 :        LE1 = LECT &BIN2 IIN3;
1229 :        M2 = EXTR (TABMOD . IIN3) 'MAIL';
1230 : *      On itere si (M1 inclut dans M2) ou (M2 inclut dans M1)
1231 :        SI (EGA (NBNO (M1 DIFF M2)) 0);
1232 :          ITER BIN3;
1233 :        FINSI;
1234 : *      On itere si M1 et M2 n ont pas de noeuds communs
1235 :        SI (EGA ((NBNO M1) + (NBNO M2)) (NBNO (M1 ET M2)));
1236 :          ITER BIN3;
1237 :        FINSI;
1238 : *      On recupere l interface M1-M2
1239 :        L1 = (CONT M1) ELEM 'APPU' (CONT M2);
1240 :        N1 = NBNO M1; N2 = NBNO M2;
1241 : *      LO1=vrai <=> il existe des noeuds communs a M1 et M2
1242 :        LO1 = NEG (N1 + N2) (NBNO (M1 ET M2));
1243 : *      LO2=vrai <=> il n'y a pas 1 noeud commun a M1 et M2
1244 :        LO2 = NEG ('ABS' ((N1 + N2) - (NBNO (M1 ET M2)))) 1;
1245 : *      LO4 = M1 et M2 forment bien une interface et ne se chevauchent pas
1246 :        LO4 = NEG (NBEL L1) 0;
1247 : *       SI (LO1 ET LO2 ET LO3);
1248 :        SI (LO1 ET LO2 ET LO4);
1249 : *        on ajoute l interface car on a >1 noeuds en commun a M1 et M2
1250 :          LINTER.LE1 = L1;
1251 : *        IPARAL=vrai <=> toutes les interfaces sont // a la fissure
1252 : *        rem: si IPARAL=faux, alors il faut ajouter des termes d interfaces au c
1253 :          P1 = (POIN L1 1) 'MOIN' (POIN L1 2);
1254 :          P1 = P1 / (NORM P1);
1255 : *petite modif car DIRTETA doit etre un chpoint desormais (a verifier)...
1256 :          PDIRTETA = resu DIRTETA;
1257 :          Presu = (extr PDIRTETA 'MAIL') poin 1;
1258 :          xDIRTETA = extr PDIRTETA Presu 'UX';
1259 :          yDIRTETA = extr PDIRTETA Presu 'UY';
1260 :          si(&DIME ega 2);
1261 :            PDIRTETA = xDIRTETA yDIRTETA;
1262 :          sino;
1263 :            zDIRTETA = extr PDIRTETA Presu 'UZ';
1264 :            PDIRTETA = xDIRTETA yDIRTETA zDIRTETA;
1265 :          fins;
1266 :          PDIRTETA = PDIRTETA / (norm PDIRTETA);
1267 :          LO1 = ((EGA P1 PDIRTETA 1.E-6) OU
1268 :                 (EGA P1 (-1.*PDIRTETA) 1.E-6));
1269 :          IPARAL = IPARAL ET LO1;
1270 :        FINSI;
1271 :     FIN BIN3;
1272 :   FIN BIN2;
1273 : FINSI;
1274 : 
1275 : **** dans le cas decouplage seulement :
1276 : **** TEST SI FRONT_FISSURE EST DANS UNE INTERFACE ****
1277 : IDANS = FAUX;
1278 : SI ((NEG (dime LINTER) 0) ET (EGA IINTE 99));
1279 :   IND1 = INDE LINTER;
1280 :   REPETER BIN4 (DIME IND1);
1281 :      LE1 = IND1.&BIN4;
1282 :      M1 = CHAN LINTER.LE1 'POI1';
1283 :      N1 = NBNO M1;
1284 :      N2 = NBNO (M1 ET SUPTAB.'FRONT_FISSURE');
1285 :      SI (EGA N1 N2);
1286 :        IDANS = VRAI;
1287 :        QUIT BIN4;
1288 :      SINON;
1289 :        ITER BIN4;
1290 :      FINSI;
1291 :   FIN BIN4;
1292 : FINSI;
1293 : *** SI OUI (IDANS),ON DETERMINE LES MODELES SUP ET INF *******
1294 : MODINF = 0;  MODSUP = 0;
1295 : SI IDANS;
1296 : *  on redéfinit : IPARAL=vrai <=> l'interface a laquelle appartient la fissure e
1297 :    IPARAL=FAUX;
1298 :    M1 = EXTR TABMOD.(EXTR LE1 1) 'MAIL';
1299 :    M2 = EXTR TABMOD.(EXTR LE1 2) 'MAIL';
1300 :    LSUP = SUPTAB.'LEVRE_SUPERIEURE';
1301 :    LINF = SUPTAB.'LEVRE_INFERIEURE';
1302 :    N1 = NBNO M1;
1303 :    N2 = NBNO M2;
1304 :    NLSUP = NBNO LSUP;
1305 :    NLINF = NBNO LINF;
1306 :    N1SUP = NBNO (M1 ET LSUP);
1307 :    N2INF = NBNO (M2 ET LINF);
1308 :    SI ( ((N1 + NLSUP - N1SUP) > 1) ET ((N2 + NLINF - N2INF) > 1));
1309 : *  LSUP et Mod1 ont plus d'1 point commun  ET  idem pour LINF et Mod2
1310 :      MODSUP = TABMOD.(EXTR LE1 1);
1311 :      MODINF = TABMOD.(EXTR LE1 2);
1312 :    SINON;
1313 :       N1INF = NBNO (M1 ET LINF);
1314 :       N2SUP = NBNO (M2 ET LSUP);
1315 :       SI ( ((N1 + NLINF - N1INF) > 1) ET ((N2 + NLSUP - N2SUP) > 1));
1316 :          MODSUP = TABMOD.(EXTR LE1 2);
1317 :          MODINF = TABMOD.(EXTR LE1 1);
1318 :       SINON;
1319 :          MESS 'ERREUR : INCOMPATIBILITE ENTRE LE MODELES_COMPOSITES';
1320 :          MESS '         ET LES LEVRE_SUPERIEURE ET _INFERIEURE';
1321 :          erre 21; QUIT G_THETA;
1322 :       FINSI;
1323 :    FINSI;
1324 : *  LA FISSURE EST BIEN DANS LE PROLONGEMENT DE L' INTERFACE
1325 :    IPARAL=VRAI;
1326 : FINSI;
1327 : * REM: il faudrait egalement verifier que MODSUP et MODINF suffisent a decrire E
1328 : 
1329 : 
1330 : **************************************************
1331 : *** MODPLA : table indicée par entier pour stocker les modèles
1332 : ***          mécaniques de chaque objet MMODEL.
1333 : **************************************************
1334 : *** Elle est vide si le modèle est élastique ou élastoplastique
1335 : *** avec une courbe de traction independante de la température.
1336 : *** Dans le cas contraire la table vaut :
1337 : ***    1 si le modèle est plastique isotrope. Alors une
1338 : ***      nouvelle courbe de traction EPSE-SIGMA est faite.
1339 : ***    2 si le modèle est plastique cinématique
1340 : ***    3 si le modèle est plastique parfaite
1341 : 
1342 : YOUVARI = FAUX; NUVARI = FAUX;
1343 : ALFVARI = FAUX; MODPLA = TABLE; TABTRA = TABLE;
1344 : REPETER BCMOD1 NBOBJ;
1345 :    MODI = TABMOD.&BCMOD1;
1346 :    MATI = REDU OBJMAT MODI;
1347 : *
1348 : *** YOUVARI **************************************
1349 :    YO1 = EXCO MATI 'YOUN';
1350 :    TYPYO = TYPE (EXTR YO1 'YOUN' 1 1 1);
1351 :    SI (EGA TYPYO 'EVOLUTIO');
1352 :      YOUVARI = VRAI;
1353 :    SINON;
1354 :      TEST1 = ((MAXI YO1) - (MINI YO1))/(MINI YO1);
1355 :     SI (TEST1 '>' 1.E-10);
1356 :        YOUVARI = VRAI;
1357 :     FINSI;
1358 :    FINSI;
1359 : *
1360 : *** NUVARI **************************************
1361 :    NU1 = EXCO MATI 'NU';
1362 :    TYPNU = TYPE (EXTR NU1 'NU' 1 1 1);
1363 :    SI (EGA TYPNU 'EVOLUTIO');
1364 :      NUVARI = VRAI;
1365 :    SINON;
1366 :      TEST1 = ((MAXI NU1) - (MINI NU1))/(MINI NU1);
1367 :     SI (TEST1 '>' 1.E-10);
1368 :        NUVARI = VRAI;
1369 :     FINSI;
1370 :    FINSI;
1371 : *
1372 : *** ALFVARI **************************************
1373 :    SI ITHER;
1374 :      AL1 = EXCO MATI 'ALPH';
1375 :      TYPAL = TYPE (EXTR AL1 'ALPH' 1 1 1);
1376 :     SI (EGA TYPAL 'EVOLUTIO');
1377 :        ALFVARI = VRAI;
1378 :     SINON;
1379 :        TEST1 = ((MAXI AL1) - (MINI AL1))/(MINI AL1);
1380 :       SI (TEST1 '>' 1.E-10);
1381 :          ALFVARI = VRAI;
1382 :       FINSI;
1383 :     FINSI;
1384 :    FINSI;
1385 : *
1386 : *** courbe de TRACtion ***************************
1387 :    SI (EXIS MATI 'TRAC');
1388 :      TR1 = EXCO MATI 'TRAC';
1389 :      TYPTR = TYPE (EXTR TR1 'TRAC' 1 1 1);
1390 :      SI (EGA TYPTR 'NUAGE   ');
1391 :        MODPLA.&BCMOD1 = 1;
1392 :        TRA1 = EXTR TR1 'TRAC' 1 1 1; COM1 = EXTR TRA1 'COMP';
1393 :        NOMEVO1 = MOT 'TRAC'; NOMFLO1 = MOT 'T';
1394 :       REPETER BNUA1 (DIME TRA1 'UPLE');
1395 :         SI (EGA &BNUA1 1);
1396 :            NUA1 = EXTR TRA1 'MINI' NOMFLO1;
1397 :         SINON;
1398 :            NUA1 = EXTR TRA1 'SUPE' NOMFLO1 (T1 + 1.E-10);
1399 :         FINSI;
1400 :          T1 =  EXTR NUA1 NOMFLO1;
1401 :          EV1 = EXTR NUA1 NOMEVO1;
1402 :          PSIG1 = EXTR EV1 ORDO; PEPS1 = EXTR EV1 'ABSC';
1403 :          VYOU1 = (EXTR 2 PSIG1) / (EXTR 2 PEPS1);
1404 :          PEPS2 = PROG;
1405 :         REPETER BSIG1 ((DIME PSIG1) - 1);
1406 :            VA1 = (EXTR (&BSIG1 + 1) PEPS1) -
1407 :                 ((EXTR (&BSIG1 + 1) PSIG1) / VYOU1);
1408 :            PEPS2 = PEPS2 ET (PROG VA1);
1409 :         FIN BSIG1;
1410 :          EV1 = EVOL 'MANU' 'EPSE' PEPS2 SIGM ('ENLE' PSIG1 1);
1411 :         SI (&BNUA1 EGA 1);
1412 :            TRA2 = 'NUAG' 'COMP' NOMFLO1 T1 'COMP' NOMEVO1 EV1;
1413 :         SINON;
1414 :            TRA2 = TRA2 ET ('NUAG' 'COMP' NOMFLO1
1415 :                             T1 'COMP' NOMEVO1 EV1);
1416 :         FINSI;
1417 :       FIN BNUA1;
1418 :        TABTRA.&BCMOD1 = TRA2;
1419 : *** On enlève la courbe de traction si elle depend de
1420 : *** la temperature (operation trop couteuse pour VARI)
1421 :        MAT0 = MATI; LCOMP1 = EXTR MAT0 'COMP';
1422 :       REPETER BCOM1 (DIME LCOMP1);
1423 :          C1 = EXTR LCOMP1 &BCOM1;
1424 :         SI (NEG C1 'TRAC');
1425 :           SI (EGA &BCOM1 1);
1426 :              MATI = 'MATE' MODI C1 (EXCO C1 MAT0);
1427 :           SINON;
1428 :              MATI = MATI ET ('MATE' MODI C1 (EXCO C1 MAT0));
1429 :           FINSI;
1430 :         FINSI;
1431 :       FIN BCOM1;
1432 :     FINSI;
1433 : 
1434 : *** SIGY (et pas TRAC) *************************
1435 :   SINON;
1436 :     SI (EXIS MATI 'SIGY');
1437 :        SI1 = EXCO MATI 'SIGY';
1438 :        TYPSI = TYPE (EXTR SI1 'SIGY' 1 1 1);
1439 :       SI (EXIS MATI 'H');
1440 :          H1 = EXCO MATI 'H';
1441 :          TYPH = TYPE (EXTR H1 'H' 1 1 1);
1442 :         SI ((EGA TYPH 'EVOLUTIO') OU
1443 :               (EGA TYPSI 'EVOLUTIO'));
1444 :            MODPLA.&BCMOD1 = 2;
1445 :         FINSI;
1446 :       SINON;
1447 :         SI (EGA TYPSI 'EVOLUTIO');
1448 :            MODPLA.&BCMOD1 = 3;
1449 :         FINSI;
1450 :       FINSI;
1451 :     FINSI;
1452 :   FINSI;
1453 : FIN BCMOD1;
1454 : MATVARI = YOUVARI OU NUVARI OU ALFVARI  OU ((DIME MODPLA) '>' 0);
1455 : 
1456 : 
1457 : **************************************************
1458 : ************ CAS IMPOSSIBLE A TRAITER ************
1459 : **************************************************
1460 : 
1461 : SI (EXIS SUPTAB 'MODELES_COMPOSITES');
1462 : *  SI ((EGA IINTE 99) ET (NON IPARAL));
1463 :   SI ((EGA IINTE 99) ET (IDANS ET (NON IPARAL)));
1464 :      MESS 'ERREUR : ON NE PEUT ENCORE DECOUPLER LES MODES';
1465 :      MESS '         DANS LE CAS DES MATERIAUX COMPOSITES';
1466 :      MESS '         SI LA FISSURE N APPARTIENT PAS A L INTERFACE';
1467 :      erre 21; QUIT G_THETA;
1468 :   FINSI;
1469 : *  SI ((EGA IINTE 4) ET (NON IPARAL));
1470 :   SI ((EGA IINTE 4) ET (NEG (dime LINTER) 0));
1471 :      MESS 'ERREUR : ON NE PEUT ENCORE CALCULER DJ/DA';
1472 :      MESS '         DANS LE CAS DES MATERIAUX COMPOSITES';
1473 :      erre 21; QUIT G_THETA;
1474 :   FINSI;
1475 : SINON;
1476 :   SI ((EGA IINTE 99) ET MATVARI);
1477 :      MESS 'ERREUR : ON NE PEUT DECOUPLER LES MODES';
1478 :      MESS '         DANS LE CAS DES CARACTERISTIQUES';
1479 :      MESS '         MATERIELS VARIABLES DANS L ESPACE';
1480 :      erre 21; QUIT G_THETA;
1481 :   FINSI;
1482 :   SI ((EGA IINTE 4) ET MATVARI);
1483 :      MESS 'ERREUR : ON NE PEUT ENCORE CALCULER DJ/DA DANS';
1484 :      MESS '         LE CAS DES CARACTERISTIQUES MATERIELS';
1485 :      MESS '         VARIABLES DANS L ESPACE';
1486 :      erre 21; QUIT G_THETA;
1487 :   FINSI;
1488 : FINSI;
1489 : *
1490 : SI ((EGA IINTE 99) ET ICOQU);
1491 :    MESS 'ERREUR : ON NE PEUT ENCORE DECOUPER LES MODES';
1492 :    MESS '         DANS LE CAS DES ELEMENTS DE COQUE';
1493 :    erre 21; QUIT G_THETA;
1494 : FINSI;
1495 : SI ((EGA IINTE 4) ET ICOQU);
1496 :    MESS 'ERREUR : ON NE PEUT ENCORE CALCULER DJ/DA';
1497 :    MESS '         DANS LE CAS DES ELEMENTS DE COQUE';
1498 :    erre 21; QUIT G_THETA;
1499 : FINSI;
1500 : SI ((EGA IINTE 3) ET ITHER);
1501 :    MESS 'ERREUR : ON NE PEUT ENCORE CALCULER L INTEGRALE';
1502 :    MESS '         C*(H) DANS LE CAS DE CHARGEMENT THERMIQUE';
1503 :    erre 21; QUIT G_THETA;
1504 : FINSI;
1505 : SI ((EGA IINTE 2) ET ITHER);
1506 :    MESS 'ERREUR : ON NE PEUT ENCORE CALCULER L INTEGRALE';
1507 :    MESS '         C* DANS LE CAS DE CHARGEMENT THERMIQUE';
1508 :    erre 21; QUIT G_THETA;
1509 : FINSI;
1510 : *
1511 : * SI (IGDEP ET (NEG IINTE 1));
1512 : *BP : on propose de calculer malgré tout...
1513 : SI (IGDEP ET ((NEG IINTE 1) et (NEG IINTE 5) et (NEG IINTE 99)));
1514 :    MESS 'ERREUR : ON NE PEUT ENCORE CALCULER L INTEGRALE';
1515 :    MESS '         SPECIFIEE EN GRANDS-DEPLACEMENTS';
1516 :    erre 21; QUIT G_THETA;
1517 : FINSI;
1518 : 
1519 : 
1520 : **************************************************
1521 : *** TITRES A AFFICHER SELON LE PROBLEME TRAITER **
1522 : **************************************************
1523 : SI IPAP;
1524 :    TXMECANI= MOT '  Mecanique';
1525 :    TXTERMI = MOT '    Thermique';
1526 :    TXPRESS = MOT '    Volumique';
1527 : SINON;
1528 :    TXMECANI= MOT ' Mecanique';
1529 :    TXTERMI = MOT '    Thermique';
1530 :    TXPRESS = MOT '    Volumique';
1531 : FINSI;
1532 : MC1 = MOT ' EN VISCO-THERMO-PLASTIQUE';
1533 : MC2 = MOT ' EN VISCO-PLASTICITE';
1534 : MC3 = MOT ' EN THERMO-PLASTICITE';
1535 : MC4 = MOT ' EN ELASTO-PLASTICITE';
1536 : si(exis SUPTAB 'COUCHE');
1537 :    MC10 = CHAI ' (Theta ' SUPTAB.'COUCHE' ')';
1538 : sino;
1539 :    MC10 = CHAI ' (Theta utilisateur)';
1540 : fins;
1541 : SI (EGA IINTE 1);
1542 :    CHA1 = CHAI 'INTEGRALE J EN FONCTION DU TEMPS';
1543 :    MOTTI = MOT 'J';
1544 :    MOTCO = MOT '         J';
1545 :   SI ITHER;
1546 :      TX1 = CHAI '        INTEGRALE J' MC3;
1547 :   SINON;
1548 :      TX1 = CHAI '        INTEGRALE J' MC4;
1549 :   FINSI;
1550 : FINSI;
1551 : SI (EGA IINTE 2);
1552 :    CHA1 = CHAI 'INTEGRALE C* EN FONCTION DU TEMPS';
1553 :    MOTTI = MOT 'C*';
1554 :    MOTCO = MOT '        C*';
1555 :   SI ITHER;
1556 :      TX1 = CHAI '       INTEGRALE C*' MC1;
1557 :   SINON;
1558 :      TX1 = CHAI '        INTEGRALE C*' MC2;
1559 :   FINSI;
1560 : FINSI;
1561 : SI (EGA IINTE 3);
1562 :    CHA1 = CHAI 'INTEGRALE C*H EN FONCTION DU TEMPS';
1563 :    MOTTI = MOT 'C*H';
1564 : 
1565 :    MOTCO = MOT '       C*(H)';
1566 :   SI ITHER;
1567 :      TX1 = CHAI '          INTEGRALE C*H' MC1;
1568 :   SINON;
1569 :      TX1 = CHAI '       INTEGRALE C*H' MC2;
1570 :   FINSI;
1571 : FINSI;
1572 : SI (EGA IINTE 4);
1573 :   SI (EXIS SUPTAB 'FISSURE_2');
1574 :      CHA1 = CHAI 'INTEGRALE CROISEE DJi/DAj EN FONCTION DU TEMPS';
1575 :      MOTTI = MOT 'DJi/DAj';
1576 :      MOTCO = MOT '      dJi/dAj';
1577 :     SI ITHER;
1578 :        TX1 = CHAI '      INTEGRALE CROISEE DJi/DAj' MC3;
1579 :     SINON;
1580 :        TX1 = CHAI '  INTEGRALE CROISEE DJi/DAj' MC4;
1581 :     FINSI;
1582 :   SINON;
1583 :      CHA1 = CHAI 'INTEGRALE DJ/DA EN FONCTION DU TEMPS';
1584 :      MOTTI = MOT 'DJ/DA';
1585 :      MOTCO = MOT '       dJ/dA';
1586 :     SI ITHER;
1587 :        TX1 = CHAI '      INTEGRALE DJ/DA' MC3;
1588 :     SINON;
1589 :        TX1 = CHAI '      INTEGRALE DJ/DA' MC4;
1590 :     FINSI;
1591 :   FINSI;
1592 : FINSI;
1593 : SI (EGA IINTE 5);
1594 :    CHA1 = CHAI 'INTEGRALE J DYNAMIQUE EN FONCTION DU TEMPS';
1595 :    MOTTI = MOT 'J_DYNA';
1596 :    MOTCO = MOT '      J_DYNA';
1597 :   SI ITHER;
1598 :      TX1 = CHAI '     INTEGRALE J EN THERMO-ELASTO-DYNAMIQUE';
1599 :   SINON;
1600 :      TX1 = CHAI '         INTEGRALE J EN ELASTO-DYNAMIQUE';
1601 :   FINSI;
1602 : FINSI;
1603 : SI (EGA IINTE 99);
1604 :    CHA1=CHAI 'F.I.C. Ki EN FONCTION DU TEMPS';
1605 :    MOTTI = MOT 'Ki';
1606 :    MOTCO = MOT '         K';
1607 :   SI ITHER;
1608 :      TX1 = CHAI ' SEPARATION DES F.I.C.' MC3;
1609 :   SINON;
1610 :      TX1 = CHAI ' SEPARATION DES F.I.C.' MC4;
1611 :   FINSI;
1612 : FINSI;
1613 : *
1614 : TX2 = CHAI ' Contribution due au chargement' MC10;
1615 : TX3 = CHAI ' °°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°';
1616 : si(exis SUPTAB 'COUCHE');
1617 :    SI ('>' SUPTAB.'COUCHE' 9);
1618 :    TX3 = CHAI ' °°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°';
1619 :    FINSI;
1620 :    SI ('>' SUPTAB.'COUCHE' 99);
1621 :    TX3 = CHAI ' °°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°';
1622 :    FINSI;
1623 : sino;
1624 :   TX3 = CHAI ' °°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°';
1625 : fins;
1626 : 
1627 : ***************************************************
1628 : ****************  AFFICHAGE DU TITRE **************
1629 : ***************************************************
1630 : SI IPAP;
1631 :   SI (EGA &DIME 2);
1632 :     SI (EGA IINTE 99);
1633 :       MESS '          ' TX1;
1634 :       MESS '           ' TX2; MESS '           ' TX3;
1635 :       MESS 'Mode  No.Pas ' TXMECANI TXTERMI TXPRESS MOTCO;
1636 :     SINON;
1637 :       MESS '       ' TX1;
1638 :       MESS '          ' TX2; MESS '          ' TX3;
1639 :       MESS ' No.Pas ' TXMECANI TXTERMI TXPRESS MOTCO;
1640 :     FINSI;
1641 :   FINSI;
1642 :   SI ((EGA &DIME 3) ET (NON ICOQU));
1643 :     SI (EGA IINTE 99);
1644 :       MESS '              ' TX1;
1645 :       MESS '               ' TX2; MESS '               ' TX3;
1646 :       MESS 'Mode  Noeud  No.Pas ' TXMECANI TXTERMI TXPRESS MOTCO;
1647 :     SINON;
1648 :       MESS '           ' TX1;
1649 :       MESS '              ' TX2; MESS '              ' TX3;
1650 :       MESS '  Noeud  No.Pas ' TXMECANI TXTERMI TXPRESS MOTCO;
1651 :     FINSI;
1652 :   FINSI;
1653 :   SI ICOQU;
1654 :     MESS '           ' TX1;
1655 :     MESS;
1656 :     MESS '              ' TX2; MESS '              ' TX3;
1657 :     MESS '  Plan   No.Pas ' TXMECANI TXTERMI TXPRESS MOTCO;
1658 :   FINSI;
1659 : SINON;
1660 :   SI (EGA &DIME 2);
1661 :     SI (EGA IINTE 99);
1662 :        MESS '        ' TX2; MESS '        ' TX3;
1663 :        MESS 'Mode  ' TXMECANI TXTERMI TXPRESS MOTCO;
1664 :     SINON;
1665 :        MESS '    ' TX2; MESS '    ' TX3;
1666 :        MESS ' ' TXMECANI TXTERMI TXPRESS MOTCO;
1667 :     FINSI;
1668 :   FINSI;
1669 :   SI ((EGA &DIME 3) ET (NON ICOQU));
1670 :     SI (EGA IINTE 99);
1671 :        MESS '            ' TX2; MESS '            ' TX3;
1672 :        MESS ' Mode   Noeud  ' TXMECANI TXTERMI TXPRESS MOTCO;
1673 :     SINON;
1674 :        MESS '          ' TX2; MESS '          ' TX3;
1675 :        MESS '  Noeud  ' TXMECANI TXTERMI TXPRESS MOTCO;
1676 :     FINSI;
1677 :   FINSI;
1678 :   SI ICOQU;
1679 :     MESS '          ' TX2; MESS '          ' TX3;
1680 :     MESS '  Plan   ' TXMECANI TXTERMI TXPRESS MOTCO;
1681 :   FINSI;
1682 : FINSI;
1683 : 
1684 : 
1685 : *******************************************************
1686 : * CONDITIONS AUX LIMITES POUR LE DECOUPLAGE DES MODES *
1687 : * utile dans le cas de l utilisation d une METHODE MECANIQUE
1688 : * pour la creation des champs auxiliaires
1689 : *******************************************************
1690 : SI (EGA IINTE 99);
1691 :  PM = SUPTAB.'FRONT_FISSURE';
1692 :  SI ((EGA &DIME 3) ET (NON IDANS));
1693 : *   SI (EGA &DIME 2);
1694 : *      X1 Y1 = COOR MAILLAGE; X0 Y0 = COOR PM;
1695 : *      DIS1 = (((X1 - X0)**2) + ((Y1 - Y0)**2))**0.5;
1696 : *   SINON;
1697 : *      X1 Y1 Z1 = COOR MAILLAGE;
1698 :       SI (NON ICOQU); PM = POIN PM 'INIT'; FINSI;
1699 : *      X0 Y0 Z0 = COOR PM;
1700 : *      DIS1 = (((X1 - X0)**2) + ((Y1 - Y0)**2) + ((Z1 - Z0)**2))**0.5;
1701 : *   FINSI;
1702 : *   PLOIN1 = POIN 1 (POIN 'MAXI' DIS1);
1703 : *   SI (EGA &DIME 2);
1704 : *      X1 Y1 = COOR (MAILLAGE DIFF (SUPTAB.'LEVRE_INFERIEURE'
1705 : *              ET SUPTAB.'LEVRE_SUPERIEURE'));
1706 : *      X0 Y0 = COOR PLOIN1;
1707 : *      DIS1 = (((X1 - X0)**2) + ((Y1 - Y0)**2))**0.5;
1708 : *   SINON;
1709 : *      X1 Y1 Z1 = COOR (MAILLAGE DIFF (SUPTAB.'LEVRE_INFERIEURE'
1710 : *                 ET SUPTAB.'LEVRE_SUPERIEURE'));
1711 : *      X0 Y0 Z0 = COOR PLOIN1;
1712 : *      DIS1 = (((X1 - X0)**2) + ((Y1 - Y0)**2) + ((Z1 - Z0)**2))**0.5;
1713 : *   FINSI;
1714 : *   PLOIN2 = POIN 1 (POIN 'MAXI' DIS1);
1715 : *   BLOQ0 = (BLOQ 'DEPL' 'ROTA' PLOIN2) ET (BLOQ MU2 PLOIN1);
1716 : 
1717 : *BP : on fait + simple pour les CL utiles au calcul des champs aux.
1718 : *     + tard on fera probablement du tout analytique comme pour les xfem
1719 : *     car la methode mecanique n est valable que pour fissure plane et
1720 : *     front rectiligne a cause de la direction du chargement difficile a
1721 : *     definir sinon
1722 :     si(NON ICOQU);
1723 :       BLOQ1 = BLOQ 'DEPL' MESHFR1;
1724 :     sino;
1725 :       BLOQ1 = (BLOQ 'DEPL' MESHFR1) et (BLOQ 'ROTA' MESHFR1);
1726 :     fins;
1727 : 
1728 :  FINSI;
1729 : FINSI;
1730 : 
1731 : 
1732 : ******************************************
1733 : * FISSURE DANS LE REPERE GLOBAL ET LOCAL *
1734 : ******************************************
1735 : 
1736 : * ELEMENT FINI STANDARD ******************
1737 : SI ((EGA IINTE 99) et (non IXFEM));
1738 : 
1739 : *** CAS 2D *******************************
1740 :   SI (EGA &DIME 2);
1741 : *   Inclinaison de la fissure par rapport à l'axe global
1742 :     XG0 YG0 = COOR PM;
1743 :     SEG1 = ORDO (SUPTAB.'LEVRE_SUPERIEURE' ELEM 'APPU' 'LARG' PM) ;
1744 :     P_SUP = POIN SEG1 'INIT';
1745 :     SI (EGA P_SUP PM);
1746 :        P_SUP = POIN SEG1 'FINA';
1747 :     FINSI;
1748 :     SEG1 = ORDO (SUPTAB.'LEVRE_INFERIEURE' ELEM 'APPU' 'LARG' PM) ;
1749 :     P_INF = POIN SEG1 'INIT';
1750 :     SI (EGA P_INF PM);
1751 :        P_INF = POIN SEG1 'FINA';
1752 :     FINSI;
1753 :     XP1 = COOR 1 P_SUP; XP2 = COOR 1 P_INF;
1754 :     YP1 = COOR 2 P_SUP; YP2 = COOR 2 P_INF;
1755 :     ALPHA1 = ATG (YG0 - ((YP1 + YP2)/2.)) (XG0 - ((XP1 + XP2)/2.));
1756 : *   Coordonnées dans le repère Global et Local
1757 :     XG1 YG1 = COOR ELTETA;
1758 :     XL1 = ((XG1 - XG0)*(COS ALPHA1)) + ((YG1 - YG0)*(SIN ALPHA1));
1759 :     YL1 = ((YG1 - YG0)*(COS ALPHA1)) - ((XG1 - XG0)*(SIN ALPHA1));
1760 :     SI ('<' ('MESU'  ('DROI' 1 P_SUP P_INF)) 1.E-10);
1761 :        L1 = SUPTAB.'LEVRE_SUPERIEURE' ELEM 'APPU' ELTETA;
1762 :        C1 = MANU 'CHPO' L1 1 'SCAL' 1.E-10;
1763 :        L2 = SUPTAB.'LEVRE_INFERIEURE' ELEM 'APPU' ELTETA;
1764 :        C2 = MANU 'CHPO' L2 1 'SCAL' -1.E-10;
1765 :        YL1 = YL1 + C1 + C2;
1766 :     FINSI;
1767 : *   Coordonnées cylindriques RAY1 TETA1 (1.E-10 pour eviter erreur atg 0 0
1768 :     TETA1 = ATG YL1 (XL1 + 1.E-10);
1769 :     RAY1 = (((XL1*XL1) + (YL1*YL1))**0.5) + 1.E-10;
1770 :     M1 = ELTETA ELEM 'APPU' 'LARG' P_SUP;
1771 :     M2 = ELTETA ELEM 'APPU' 'LARG' P_INF;
1772 :     VA1 = XTY (MANU 'CHPO' M1 1 'SCAL' 1.)
1773 :               (REDU YL1 M1) MTS1 MTS1;
1774 :     VA2 = XTY (MANU 'CHPO' M2 1 'SCAL' 1.)
1775 :               (REDU YL1 M2) MTS1 MTS1;
1776 : *   On inverse afin d'avoir YL1 > 0 pour modsup et <0 pour modinf
1777 :     SI (('<' VA1 0.) ET ('>' VA2 0.));
1778 :        PPPP = P_SUP;  P_SUP = P_INF;   P_INF = PPPP;
1779 :        MMDD = MODSUP; MODSUP = MODINF; MODINF = MMDD;
1780 :     FINSI;
1781 :     SI ((EGA XP1 XP2 1.E-10) ET (EGA YP1 YP2 1.E-10));
1782 :        TETA_S = REDU TETA1 SUPTAB.'LEVRE_SUPERIEURE';
1783 :        TETA_F = REDU TETA1 SUPTAB.'LEVRE_INFERIEURE';
1784 :        TETA1 = TETA1 - TETA_S - TETA_F;
1785 : *      On se debrouille pour avoir exactement +/-180 sur levre sup/inf
1786 :        SI (('>' VA1 0.) ET ('<' VA2 0.));
1787 :          TETA1 = TETA1 + ((TETA_S*0.) + 180.) + ((TETA_F*0.) - 180.);
1788 :        SINON;
1789 :          TETA1 = TETA1 + ((TETA_F*0.) + 180.) + ((TETA_S*0.) - 180.);
1790 :        FINSI;
1791 :     FINSI;
1792 : *   valeur en radian
1793 :     TETA1rad = TETA1*VALPI/180.;
1794 : *   cas d une interface: on construit 3 zones :
1795 : *   PM1 = points du modele 1 / interface
1796 : *   PM2 = points du modele 2 / interface
1797 : *   L1  = elem de l interface
1798 :     SI IDANS;
1799 :        M1 = EXTR MODSUP 'MAIL';
1800 :        M2 = EXTR MODINF 'MAIL';
1801 :        L1 = (CONT M1) ELEM 'APPU' (CONT M2);
1802 :        M1 = M1 ELEM 'APPU' 'STRI' ELTETA;
1803 :        M2 = M2 ELEM 'APPU' 'STRI' ELTETA;
1804 :        PM1 = (CHAN M1 'POI1') DIFF (CHAN L1 'POI1');
1805 :        PM2 = (CHAN M2 'POI1') DIFF (CHAN L1 'POI1');
1806 :     FINSI;
1807 : *     si(fltrac);
1808 : *       trac RAY1 ELTETA 'TITR' 'RAY1';
1809 : *       trac TETA1 ELTETA 'TITR' 'TETA1' (prog -180. PAS 20. 180.);
1810 : *     fins;
1811 : 
1812 : *** CAS 3D *******************************
1813 :   SINO;
1814 :     mess 'on ne fait pas ici le passage local global en 3D ...?';
1815 :     mess 'et on ne calule pas RAY1 TETA1 non plus?';
1816 :   FINSI;
1817 : 
1818 : 
1819 : FINSI;
1820 : 
1821 : * XFEM************************************
1822 : SI ((EGA IINTE 99) et (IXFEM));
1823 : 
1824 : *** CAS 2D *******************************
1825 :   SI (EGA &DIME 2);
1826 : *   On recupere les level set
1827 :     PSI1 = REDU PSI0 ELTETA;
1828 :     PHI1 = REDU PHI0 ELTETA;
1829 :     LV7  = (PSI1 NOMC 'UX') ET (PHI1 NOMC 'UY');
1830 :     GLV7 = CHAN (GRAD LV7 OBJMOD) 'TYPE' 'SCALAIRE';
1831 : *trac GLV7 objmod 'TITR' 'GLV7';
1832 : *   ce repere est il direct? (si oui/non, SDIR1=+/-1)
1833 :     GLV7po = CHAN 'CHPO' OBJMOD GLV7  'MOYE';
1834 :     XDIR1 = ((EXCO GLV7po 'UX,X' 'SCAL') * (EXCO GLV7po 'UY,Y'))
1835 :           - ((EXCO GLV7po 'UX,Y' 'SCAL') * (EXCO GLV7po 'UY,X'));
1836 :     XDIR1 = MAXI (RESU XDIR1);
1837 :     SDIR1 = SIGN 'FLOTTANT' XDIR1;
1838 : *   Angle ALPHA1 de passage local -> global
1839 :     NGPSI1 =   ( ((EXCO GLV7 'UX,X' 'SCAL')**2)
1840 :                + ((EXCO GLV7 'UX,Y' 'SCAL')**2) )**(-0.5) ;
1841 :     NGPHI1 =   ( ((EXCO GLV7 'UY,X' 'SCAL')**2)
1842 :                + ((EXCO GLV7 'UY,Y' 'SCAL')**2) )**(-0.5) ;
1843 :     COS1A = 0.5 * ( ( (EXCO GLV7 'UX,X' 'SCAL') * NGPSI1)
1844 :          + (SDIR1 * ( (EXCO GLV7 'UY,Y' 'SCAL') * NGPHI1)) );
1845 :     SIN1A = 0.5 * ( ( (EXCO GLV7 'UX,Y' 'SCAL') * NGPSI1)
1846 :          - (SDIR1 * ( (EXCO GLV7 'UY,X' 'SCAL') * NGPHI1)) );
1847 :     ALPHA1 =  (MASQ SIN1A 'EGSUPE' 0.) * (ACOS COS1A);
1848 : *trac ALPHA1 objmod 'TITR' 'ALPHA1';
1849 :     COSA2 = COS1A ** 2 ;
1850 :     SINA2 = SIN1A ** 2 ;
1851 :     SINCOSA = SIN1A * COS1A ;
1852 : *   repere local de la fissure
1853 :     XL1 = (NOMC PSI1 'SCAL') CHAN 'CHAM' OBJMOD 'STRESSES';
1854 :     YL1 = (SDIR1 * (NOMC PHI1 'SCAL')) CHAN 'CHAM' OBJMOD 'STRESSES';
1855 : *   Coordonnées cylindriques RAY1 TETA1
1856 :     UN1 = MANU 'CHML' OBJMOD 'SCAL' 1. 'TYPE' 'SCALAIRE' 'STRESSES';
1857 :     ZER1= MANU 'CHML' OBJMOD 'SCAL' 0. 'TYPE' 'SCALAIRE' 'STRESSES';
1858 :     TETA1 = (90.*UN1) - (ATG (XL1*(abs(YL1**(-1)))));
1859 :     TETA1 = ((MASQ YL1 'SUPERIEUR' 0.) * TETA1)
1860 :           - ((MASQ YL1 'INFERIEUR' 0.) * TETA1);
1861 :     TETA1 = CHAN TETA1 'TYPE' 'SCALAIRE';
1862 :     RAY1 = ((XL1**2) + (YL1**2))**0.5;
1863 :     RAY1 = CHAN RAY1 'TYPE' 'SCALAIRE';
1864 : 
1865 : *** CAS 3D *******************************
1866 :   SINO;
1867 :     mess 'DECOUPLAGE 3D XFEM en cours de dvpt...';
1868 : *   On recupere les level set
1869 :     PSI1 = REDU PSI0 ELTETA;
1870 :     PHI1 = REDU PHI0 ELTETA;
1871 :     LV7  = (PSI1 NOMC 'UX') ET (PHI1 NOMC 'UY')
1872 :         et (MANU  'CHPO' ELTETA 1 'UZ  ' 0. 'NATURE' 'DIFFUS');
1873 :     GLV7 = CHAN (GRAD LV7 OBJMOD) 'TYPE' 'SCALAIRE';
1874 :     SDIR1 = 1.;
1875 : *   creation de la matrice de rotation
1876 :     V1 = UTILTETA . 'V1';
1877 :     V2 = UTILTETA . 'V2';
1878 :     V3 = UTILTETA . 'V3';
1879 :     ROT1= (EXCO V1 (mots 'UX' 'UY' 'UZ') (mots 'UX,X' 'UY,X' 'UZ,X'))
1880 :        ET (EXCO V2 (mots 'UX' 'UY' 'UZ') (mots 'UX,Y' 'UY,Y' 'UZ,Y'))
1881 :        ET (EXCO V3 (mots 'UX' 'UY' 'UZ') (mots 'UX,Z' 'UY,Z' 'UZ,Z'));
1882 :      ROT1 = EXCO ROT1
1883 :      (mots UX,X UX,Y UX,Z UY,X UY,Y UY,Z UZ,X UZ,Y UZ,Z);
1884 :     ROT1 = CHAN 'CHAM' ROT1 OBJMOD 'STRESSES' 'GRADIENT';
1885 : *     SUPTAB . 'ROT1' = ROT1;
1886 : *   repere local de la fissure
1887 :     XL1 = (NOMC PSI1 'SCAL') CHAN 'CHAM' OBJMOD 'STRESSES';
1888 :     YL1 = (NOMC PHI1 'SCAL') CHAN 'CHAM' OBJMOD 'STRESSES';
1889 : *   Coordonnées cylindriques RAY1 TETA1
1890 :     UN1 = MANU 'CHML' OBJMOD 'SCAL' 1. 'TYPE' 'SCALAIRE' 'STRESSES';
1891 :     ZER1= MANU 'CHML' OBJMOD 'SCAL' 0. 'TYPE' 'SCALAIRE' 'STRESSES';
1892 :     TETA1 = (90.*UN1) - (ATG (XL1*(abs(YL1**(-1)))));
1893 :     TETA1 = ((MASQ YL1 'SUPERIEUR' 0.) * TETA1)
1894 :           - ((MASQ YL1 'INFERIEUR' 0.) * TETA1);
1895 :     TETA1 = CHAN TETA1 'TYPE' 'SCALAIRE';
1896 :     RAY1 = ((XL1**2) + (YL1**2))**0.5;
1897 :     RAY1 = CHAN RAY1 'TYPE' 'SCALAIRE';
1898 :   FINSI;
1899 : 
1900 : *   si(fltrac);
1901 : *     trac ROT1 OBJMOD 'TITR' 'ROT1';
1902 : *     trac RAY1 objmod 'TITR' 'RAY1';
1903 : *     trac TETA1 objmod 'TITR' 'TETA1' (prog -180. PAS 20. 180.);
1904 : *   finsi;
1905 : 
1906 : FINSI;
1907 : 
1908 : 
1909 : 
1910 : ***************************************************
1911 : ****** FORCE, DEPLACEMENT ET GRADIENT NULS   ******
1912 : ***************************************************
1913 : 
1914 : FOR000 = CHAN 'CHPO' OBJMOD ('ZERO' OBJMOD 'FORCES  ');
1915 : DEP000 = CHAN 'CHPO' OBJMOD ('ZERO' OBJMOD 'DEPLACEM');
1916 : SI (EGA &MODE 'PLANGENE');
1917 :    FOR000 = MANU 'CHPO' (EXTR OBJMOD 'MAIL') 2 'FX' 0.
1918 :             'FY' 0. 'TITR' 'FORCES  ' 'NATURE' 'DIFFUS';
1919 :    FOR000 = FOR000 ET (MANU 'CHPO' (VALE 'MODE' 'PLANGENE')
1920 :             3 'FZ' 0. 'MX' 0. 'MY' 0.
1921 :            'TITR' 'FORCES  ' 'NATURE' 'DIFFUS');
1922 :    DEP000 = MANU 'CHPO' (EXTR OBJMOD 'MAIL') 2 'UX' 0.
1923 :             'UY' 0. 'TITR' 'DEPLACEM' 'NATURE' 'DIFFUS';
1924 :    DEP000 = DEP000 ET (MANU 'CHPO' (VALE 'MODE' 'PLANGENE')
1925 :             3 'UZ' 0. 'RX' 0. 'RY' 0.
1926 :            'TITR' 'DEPLACEM' 'NATURE' 'DIFFUS');
1927 : FINSI;
1928 : CMD000 = CHAN 'NOEUD' OBJMOD ('ZERO' OBJMOD 'DEPLACEM');
1929 : CMD001 = CHAN 'STRESSES' OBJMOD ('ZERO' OBJMOD 'DEPLACEM');
1930 : GRA000 = 'ZERO' OBJMOD 'GRADIENT';
1931 : VAR000 = 'ZERO' OBJMOD 'VARINTER';
1932 : 
1933 : 
1934 : **************************************************
1935 : * NOMBRE DE BOUCLE POUR LE CALCUL DES INTEGRALES *
1936 : **************************************************
1937 : SI IPAP;
1938 :    NBG = -1;
1939 :    NBDEP = DIME (SUPTAB.'SOLUTION_PASAPAS'.'TEMPS');
1940 :    SI IREPRI ;
1941 :      NBG = SUPTAB.'IABC';
1942 :      NBDEP = NBDEP - 1 - NBG;
1943 :    FINSI;
1944 :    SI IPERSO1;
1945 :      NBG = (WTAB . 'PAS') - 1;
1946 :      NBDEP = 1;
1947 :    FINSI;
1948 : SINON;
1949 :    NBG = -1;
1950 :    NBDEP = 1;
1951 : FINSI;
1952 : 
1953 : 
1954 : ***************************************************
1955 : ** SOLUTION DU PAS PRECEDENT SI REPRISE DE CALCUL *  (ou si perso1)
1956 : ***************************************************
1957 : SI (IREPRI ou (IPERSO1 et (NBG >eg 0))) ;
1958 :    SIG1 = (SUPTAB.'SOLUTION_PASAPAS'.'CONTRAINTES'.NBG)
1959 :            REDU OBJMOD;
1960 :    SI (EXIS (SUPTAB.'SOLUTION_PASAPAS') 'VARIABLES_INTERNES');
1961 :       VAR1 = (SUPTAB.'SOLUTION_PASAPAS'.'VARIABLES_INTERNES'.NBG)
1962 :              REDU OBJMOD;
1963 :    SINON;
1964 :       VAR1 = VAR000;
1965 :    FINSI;
1966 :    MAT1 = SUPTAB.'MAT1';
1967 :    WELAS = 0.5*('ENER' OBJMOD SIG1 ('ELAS' OBJMOD SIG1 MAT1));
1968 : *    WPLAS = SUPTAB.'END1' - WELAS;
1969 :    WPLAS = (redu OBJMOD SUPTAB.'END1') - WELAS;
1970 :    SI (EGA ITYPEF 2);
1971 :       VDI1 = SUPTAB.'VDI1';
1972 :    FINSI;
1973 :    SI (((DIME MODPLA) '>' 0) ET ITHER);
1974 :       WVMIS = SUPTAB.'ENV1';
1975 :    FINSI;
1976 :    si(flmess); mess 'RECUP DU PAS PRECEDENT OK'; fins;
1977 : FINSI;
1978 : 
1979 : 
1980 : ***************************************************
1981 : ******** TABLE INFORMATION COMPLEMENTAIRE *********
1982 : ***************************************************
1983 : INFTAB = TABL;
1984 : INFTAB.'MOTTI' = MOTTI;
1985 : INFTAB.'MODCOU' = MODCOU;
1986 : INFTAB.'TABMOD' = TABMOD;
1987 : INFTAB.'LINTER' = LINTER;
1988 : INFTAB.'ICOQU' = ICOQU;
1989 : INFTAB.'IGDEP' = IGDEP;
1990 : INFTAB.'IGDER' = IGDER;
1991 : INFTAB.'IREPRI' = IREPRI;
1992 : INFTAB.'IPAP' = IPAP;
1993 : INFTAB.'ILIN' = ILIN;
1994 : INFTAB.'ITHER' = ITHER;
1995 : INFTAB.'IPARAL' = IPARAL;
1996 : INFTAB.'MATVARI' = MATVARI;
1997 : INFTAB.'YOUVARI' = YOUVARI;
1998 : INFTAB.'ALFVARI' = ALFVARI;
1999 : INFTAB.'IINTE' = IINTE;
2000 : INFTAB.'IQUA' = IQUA;
2001 : INFTAB.'ELTETA' = ELTETA;
2002 : INFTAB.'OBJMOD' = OBJMOD;
2003 : INFTAB.'MODPLA' = MODPLA;
2004 : INFTAB.'FOR000' = FOR000;
2005 : INFTAB.'DEP000' = DEP000;
2006 : INFTAB.'CMD000' = CMD000;
2007 : INFTAB.'CMD001' = CMD001;
2008 : INFTAB.'GRA000' = GRA000;
2009 : INFTAB.'IXFEM' = IXFEM;
2010 : INFTAB.'ITYPEF' = ITYPEF;
2011 : INFTAB.'IPERSO1' = IPERSO1;
2012 : * ajout sm
2013 : INFTAB.'IDEFI' = IDEFI;
2014 : *ajout BP BT pour le contact frottant
2015 : INFTAB . 'IFROT'  = IFROT  ;
2016 : si (IFROT);
2017 :   INFTAB . 'OBJCON' = OBJCON ;
2018 : fins;
2019 : *fin ajout BP BT
2020 : si(IPERSO1); INFTAB . 'ESTIMATION' = ESTIM; fins;
2021 : 
2022 : 
2023 : ***********************************************
2024 : ***********************************************
2025 : ********* BOUCLE SUR LE PAS DE CALCUL *********
2026 : ***********************************************
2027 : ***********************************************
2028 : 
2029 : REPETER BOUCEXT NBDEP ;
2030 :    IABC = NBG + &BOUCEXT;
2031 : 
2032 : ***************************************************
2033 : ** DEPLACEMENTS,CONTRAINTES ... A L INSTANT INST **
2034 : ***************************************************
2035 : 
2036 : *** SOLUTION_PASAPAS ******************************
2037 :   SI IPAP;
2038 : *** Cas PERSO1
2039 :     SI IPERSO1;
2040 :       INST   =  ESTIM . 'TEMPS' ;
2041 :       DEPINT = (ESTIM . 'DEPLACEMENTS') REDU ELTETA;
2042 :       SIGF   = (ESTIM . 'CONTRAINTES' ) REDU OBJMOD;
2043 :       SI (IGDEP ET (NON IGDER));
2044 :          SIGF = 'CAPI' SIGF DEPINT OBJMOD;
2045 :       FINSI;
2046 :       SI IGDER;
2047 :          SI (NON (EXIS (SUPTAB.'ROTATION_RIGIDIFIANTE') IABC));
2048 :           MESS 'ERREUR : Le deplacement du a une rotation';
2049 :           MESS '         rigidifiante au pas ' IABC ' n est pas donne';
2050 :           erre 21; QUIT G_THETA;
2051 :          FINSI;
2052 :          DEPINT = DEPINT -
2053 :          (REDU SUPTAB.'ROTATION_RIGIDIFIANTE'.IABC ELTETA) ;
2054 :       FINSI;
2055 :       SI (EXIS ESTIM 'VARIABLES_INTERNES');
2056 :          VARF = (ESTIM . 'VARIABLES_INTERNES') REDU OBJMOD;
2057 :       SINON;
2058 :          VARF = VAR000;
2059 :       FINSI;
2060 :       SI (EGA IINTE 2);
2061 :         SI (EGA IABC 0) ;
2062 :            DELTAT = INST + 1.E+30;
2063 :            DEPINT = ESTIM . 'DEPLACEMENTS';
2064 :            VITDFI = ESTIM . 'DEFORMATIONS_INELASTIQUES';
2065 :            SIG1 = SIGF * 1.;
2066 :         FINSI;
2067 :         SI (IABC '>' 0);
2068 :            DELTAT= INST - (ESTIM . 'TEMPS');
2069 :            DEPINT= (ESTIM . 'DEPLACEMENTS') - (ESTIM . 'DEPLACEMENTS');
2070 :            VITDFI= (ESTIM . 'DEFORMATIONS_INELASTIQUES')
2071 :                  - (ESTIM . 'DEFORMATIONS_INELASTIQUES');
2072 :         FINSI;
2073 :          DEPINT = (REDU ELTETA DEPINT) / DELTAT;
2074 :          VITDFI = (REDU ELTETA VITDFI) / DELTAT;
2075 :       FINSI;
2076 :       SI (EGA IINTE 5);
2077 :          VITF = (ESTIM .'VITESSES')      REDU ELTETA;
2078 :          ACCF = (ESTIM .'ACCELERATIONS') REDU ELTETA;
2079 :       FINSI;
2080 : *** Cas ou on appelle g_theta apres pasapas
2081 :     SINO;
2082 :       INST = SUPTAB.'SOLUTION_PASAPAS'.'TEMPS'.IABC ;
2083 :       DEPINT = (SUPTAB.'SOLUTION_PASAPAS'.'DEPLACEMENTS'.IABC)
2084 :                REDU ELTETA;
2085 :       SIGF = (SUPTAB.'SOLUTION_PASAPAS'.'CONTRAINTES'.IABC)
2086 :              REDU OBJMOD;
2087 :       SI (IGDEP ET (NON IGDER));
2088 :          SIGF = 'CAPI' SIGF DEPINT OBJMOD;
2089 :       FINSI;
2090 :       SI IGDER;
2091 :          SI (NON (EXIS (SUPTAB.'ROTATION_RIGIDIFIANTE') IABC));
2092 :           MESS 'ERREUR : Le deplacement du a une rotation';
2093 :           MESS '         rigidifiante au pas ' IABC ' n est pas donne';
2094 :           erre 21; QUIT G_THETA;
2095 :          FINSI;
2096 :          DEPINT = DEPINT -
2097 :                  (REDU SUPTAB.'ROTATION_RIGIDIFIANTE'.IABC ELTETA) ;
2098 :       FINSI;
2099 :       SI (EXIS (SUPTAB.'SOLUTION_PASAPAS') 'VARIABLES_INTERNES');
2100 :          VARF = (SUPTAB.'SOLUTION_PASAPAS'.'VARIABLES_INTERNES'.IABC)
2101 :                 REDU OBJMOD;
2102 :       SINON;
2103 :          VARF = VAR000;
2104 :       FINSI;
2105 :       SI (EGA IINTE 2);
2106 :         SI (EGA IABC 0) ;
2107 :            DELTAT = INST + 1.E+30;
2108 :            DEPINT = SUPTAB.'SOLUTION_PASAPAS'.'DEPLACEMENTS'.IABC;
2109 :            VITDFI = SUPTAB.'SOLUTION_PASAPAS'.
2110 :                    'DEFORMATIONS_INELASTIQUES'.IABC;
2111 :            SIG1 = SIGF * 1.;
2112 :         FINSI;
2113 :         SI (IABC '>' 0);
2114 :         DELTAT= INST - (SUPTAB.'SOLUTION_PASAPAS'.'TEMPS'.(IABC - 1));
2115 :            DEPINT= (SUPTAB.'SOLUTION_PASAPAS'.'DEPLACEMENTS'.IABC) -
2116 :                  (SUPTAB.'SOLUTION_PASAPAS'.'DEPLACEMENTS'.(IABC - 1));
2117 :            VITDFI=
2118 :          (SUPTAB.'SOLUTION_PASAPAS'.'DEFORMATIONS_INELASTIQUES' .IABC)
2119 :  - (SUPTAB.'SOLUTION_PASAPAS'. 'DEFORMATIONS_INELASTIQUES'.(IABC - 1));
2120 :         FINSI;
2121 :          DEPINT = (REDU ELTETA DEPINT) / DELTAT;
2122 :          VITDFI = (REDU ELTETA VITDFI) / DELTAT;
2123 :       FINSI;
2124 :       SI (EGA IINTE 5);
2125 :          VITF = (SUPTAB.'SOLUTION_PASAPAS'.'VITESSES'.IABC)
2126 :                 REDU ELTETA;
2127 :          ACCF = (SUPTAB.'SOLUTION_PASAPAS'.'ACCELERATIONS'.IABC)
2128 :                 REDU ELTETA;
2129 :       FINSI;
2130 :     FINSI;
2131 : *** SOLUTION_RESO **********************************
2132 :   SINON;
2133 :      DEPINT = REDU (SUPTAB.'SOLUTION_RESO') ELTETA;
2134 :      SIGF = SIGM DEPINT OBJMOD OBJMAT;
2135 :   FINSI;
2136 : 
2137 : *** ON CHANGE LE DEPLACEMENT DEPINT EN MCHAML AU NOEUD ******
2138 : *bp: utilite?
2139 :    DEPINT = CHAN 'CHAM' DEPINT OBJMOD 'NOEUD' 'DEPLACEMENTS';
2140 : 
2141 : ****************************************************
2142 : * MODIFICATION DES CHAMPS THETA ET PI SI GRANDE ROT
2143 : ****************************************************
2144 :   SI IPAP;
2145 :     SI IGDER;
2146 :       FORM SUPTAB.'ROTATION_RIGIDIFIANTE'.IABC;
2147 :       SUPTAB.'CHAMP_THETA' UTILTETA = CH_THETA SUPTAB;
2148 :       SI (EGA IINTE 4);
2149 :         SI (NON (EXIS SUPTAB 'FRONT_FISSURE_2'));
2150 :            SUPTAB.'COUCHE' = (SUPTAB.'COUCHE') - 1;
2151 :            SUPTAB.'PI' UTILPI = CH_THETA SUPTAB;
2152 :            SUPTAB.'COUCHE' = (SUPTAB.'COUCHE') + 1;
2153 :         SINON;
2154 :            P1 = SUPTAB.'FRONT_FISSURE';
2155 :            SUPTAB.'FRONT_FISSURE' = SUPTAB.'FRONT_FISSURE_2';
2156 :            SUPTAB.'FISSURE' = SUPTAB.'FISSURE_2';
2157 :            SUPTAB.'PI' UTILPI = CH_THETA SUPTAB;
2158 :            SUPTAB.'FRONT_FISSURE' = P1;
2159 :         FINSI;
2160 :       FINSI;
2161 :     FINSI;
2162 :   FINSI;
2163 : 
2164 : ***************************************************
2165 : ********** TEMPERATURES A L INSTANT INST **********
2166 : ***************************************************
2167 :   SI ITHER ;
2168 :     SI IPAP;
2169 :        TEPINT = 'TIRE' CHAR1 INST 'T' ;
2170 :        TEPINT = TEPINT - TALPH1 ;
2171 :     SINON;
2172 :        TEPINT = SUPTAB.'TEMPERATURES';
2173 :     FINSI;
2174 :   FINSI;
2175 : 
2176 :   
2177 : ***************************************************                     
2178 : ********** DEF IMPOSEE A L INSTANT INST **********                     
2179 : ***************************************************                     
2180 :   SI IDEFI;                                                          
2181 :     SI IPAP; 
2182 :        DEFINT = 'TIRE' CHAR1 INST 'DEFI' ;      
2183 :     SINON;
2184 :        DEFINT = SUPTAB . 'DEFORMATIONS_IMPOSEES';      
2185 :     FINSI;
2186 :   FINSI;
2187 :    
2188 : ***************************************************                     
2189 : ********** CONTACT FROTTANT **********                     
2190 : ***************************************************                     
2191 :   si (IFROT);
2192 :     SI IPAP; 
2193 : *...todo    
2194 :     SINON;
2195 : *     a priori DEPLACEMENT_FISSURE pas tres utile ...
2196 :       si (exis SUPTAB 'DEPLACEMENT_FISSURE');
2197 :          WDEP  = SUPTAB . 'DEPLACEMENT_FISSURE';
2198 :       sino;
2199 :          WDEP  = REDU (SUPTAB.'SOLUTION_RESO') MAICON;    
2200 :       fins;
2201 :       toto = EXTR WDEP 'MAIL';
2202 :       si (ega (nbel toto) 0);
2203 :         MESS 'ERREUR : IL FAUT DEPLACEMENT_FISSURE si MODELE_FISSURE';
2204 :         erre 21; QUIT G_THETA;
2205 :       finsi;
2206 :       si (exis SUPTAB 'PRESSION_FISSURE');
2207 :          SIGCON = SUPTAB . 'PRESSION_FISSURE';
2208 :       sino;
2209 :          SIGCON = REDU SIGF OBJCON;
2210 :       fins;
2211 : *     peut etre faire un test sur SIGCON ...
2212 :     FINSI;
2213 :   fins;
2214 : 
2215 : *****************************************************
2216 : * CONTRAINTE RECALCULEE SI ITHER = VRAI et NON IPAP *
2217 : *****************************************************
2218 :   SI (ITHER ET (NON IPAP));
2219 :      SIGF = SIGF - ('THET' OBJMOD OBJMAT TEPINT);
2220 :   FINSI;
2221 : 
2222 : ***************************************************
2223 : ************ MATERIAU A L INSTANT INST ************
2224 : ***************************************************
2225 :   SI (MATVARI ET ITHER ET IPAP);
2226 :      TEPABS = TEPINT + TALPH1 ;
2227 :      MAT1 = VARI 'NUAG' OBJMOD OBJMAT (EXCO 'T' TEPABS 'T');
2228 :   SINON;
2229 :      MAT1 = OBJMAT;
2230 : *bp: ne devrait on pas ecrire ci dessous..? inutile?
2231 : *      MAT1 = REDU OBJMOD OBJMAT
2232 :   FINSI;
2233 : 
2234 : ***************************************************
2235 : ********* RIGIDITE TOTALE A L INSTANT INST ********
2236 : ***************************************************
2237 : *bp: utilite?
2238 :   SI ( (EGA IINTE 4) OU ((EGA IINTE 99)
2239 :         et (non IXFEM) et (EGA &DIME 3) ET (NON IDANS)) );
2240 :     SI IPAP;
2241 :       SI (MATVARI ET ITHER);
2242 :          M1 = VARI 'NUAG' OBJMOD (EXCO 'T' TEPABS 'T')
2243 :               (SUPTAB.'SOLUTION_PASAPAS'.'CARACTERISTIQUES');
2244 :       SINON;
2245 :          M1 = SUPTAB.'SOLUTION_PASAPAS'.'CARACTERISTIQUES';
2246 :       FINSI;
2247 :       RIGTOT = RIGI M1 (SUPTAB.'SOLUTION_PASAPAS'.'MODELE');
2248 :     SINON;
2249 :       RIGTOT = RIGI (SUPTAB.'CARACTERISTIQUES') (SUPTAB.'MODELE');
2250 :     FINSI;
2251 :   FINSI;
2252 : 
2253 : ***************************************************
2254 : ****** CHARGEMENT MECANIQUE A L INSTANT INST ******
2255 : ***************************************************
2256 :   PREINT = FOR000;
2257 : *** SOLUTION_PASAPAS ******************************
2258 :   SI IPAP;
2259 :     SI (EXIS CHAR1 'MECA');
2260 :        PREINT = PREINT + ((TIRE CHAR1 INST 'MECA') REDU ELTETA);
2261 :     FINSI;
2262 :     SI (EXIS CHAR1 'PSUI');
2263 :        PR = (TIRE CHAR1 'PSUI' INST) + 1.E-10;
2264 :       SI (NON ICOQU);
2265 :          M1 = EXTR (SUPTAB.'SOLUTION_PASAPAS'.'MODELE') 'ELEM'
2266 :               'TRI3' 'QUA4' 'TRI6' 'QUA8' 'ICT3' 'ICT6' 'ICQ4' 'ICQ8'
2267 :               'CUB8' 'TET4' 'PRI6' 'PYR5' 'CU20' 'TE10' 'PR15' 'PY13';
2268 :          M2 = EXTR M1 'MAIL';
2269 :          P1 = PRES 'MASS' M1 (REDU PR M2);
2270 :       FINSI;
2271 :       SI ICOQU;
2272 :         SI ILIN;
2273 :            M1 = EXTR (SUPTAB.'SOLUTION_PASAPAS'.'MODELE')
2274 :                     'ELEM' 'COQ2' 'COQ3' 'COQ4' 'DKT' 'DST' ;
2275 :            M3 = TEXT '        ';
2276 :         FINSI;
2277 :         SI IQUA;
2278 :            M1 = EXTR (SUPTAB.'SOLUTION_PASAPAS'.'MODELE')
2279 :                     'ELEM' 'COQ6' 'COQ8';
2280 :           SI (MATVARI ET ITHER);
2281 :              M3 = VARI 'NUAG' OBJMOD (EXCO 'T' TEPABS 'T')
2282 :                   (SUPTAB.'SOLUTION_PASAPAS'.'CARACTERISTIQUES');
2283 :           SINON;
2284 :              M3 = SUPTAB.'SOLUTION_PASAPAS'.'CARACTERISTIQUES';
2285 :           FINSI;
2286 :         FINSI;
2287 :          M2 = EXTR M1 'MAIL';
2288 :          P1 = PRES 'COQU' M1 (REDU PR M2) 'NORM' M3;
2289 :       FINSI;
2290 :        PREINT = (REDU P1 ELTETA) + PREINT;
2291 :     FINSI;
2292 : *** SOLUTION_RESO **********************************
2293 :   SINON ;
2294 :     SI (EXIS SUPTAB 'CHARGEMENTS_MECANIQUES');
2295 :        PREINT = PREINT +
2296 :                (SUPTAB.'CHARGEMENTS_MECANIQUES' REDU ELTETA);
2297 :     FINSI;
2298 :   FINSI;
2299 : 
2300 : 
2301 : ****************************************************
2302 : *** SOLUTIONS AUXILAIRES SI DECOUPLAGE DES MODES ***
2303 : ****************************************************
2304 : * NBMIXT = nbre d integrale a calculer (=1 si J, =2 si K1 K2, =3 si K1 K2 K3)
2305 :   NBMIXT = 1;
2306 :   IM = 0;
2307 :   SI (EGA IINTE 99);
2308 :      NBMIXT = 2;
2309 :      SI (EGA &DIME 3);
2310 :        NBMIXT = 3;
2311 :      FINSI;
2312 : 
2313 : **** CONSTANTES MATERIAUX **************************
2314 :      CHAM1 = (EXCO 'YOUN' MAT1) ET (EXCO 'NU  ' MAT1);
2315 : *   -CAS D UN MATERIAU HOMOGENE
2316 :      SI (NON IDANS);
2317 : *      on construit les champ auxilaires d'une solution d'un materiau
2318 : *      homogene => on prend les valeurs de E et nu en pointe de fissure
2319 :        CHPO1 = CHAN 'CHPO' OBJMOD CHAM1;
2320 :        si(IXFEM);
2321 :          CHPO1 = INT_COMP ELTETA CHPO1 (MANU 'POI1' PM);
2322 :        fins;
2323 :        SI (EGA (TYPE PM) 'MAILLAGE');
2324 :          VYO_1 = (maxi (resu (exco CHPO1 'YOUN'))) / (nbno PM);
2325 :          VNU_1 = (maxi (resu (exco CHPO1 'NU  '))) / (nbno PM);
2326 :        SINON;
2327 :          VYO_1 = EXTR CHPO1 'YOUN' PM;
2328 :          VNU_1 = EXTR CHPO1 'NU  ' PM;
2329 :        FINSI;
2330 : *      Constante de Kolosov
2331 :        SI (EGA &MODE 'PLANCONT');
2332 :          KAP_1 = (3. - VNU_1) / (1. + VNU_1);
2333 :        SINON;
2334 :          KAP_1 = (3. - (4. * VNU_1));
2335 :        FINSI;
2336 : *      Module de cisaillement
2337 :        MU_1 = VYO_1 / (2.*(1. + VNU_1));
2338 : *      Constante C_MATE (= 1 / E^etoile)
2339 :        SI (EGA &MODE 'PLANCONT');
2340 :          C_MATE = 1. / VYO_1;
2341 :        SINON;
2342 :          C_MATE = (1. - (VNU_1*VNU_1)) / VYO_1;
2343 :        FINSI;
2344 : *   -CAS D UN BI-MATERIAU
2345 :      SINON;
2346 :        CHPO1 = CHAN 'CHPO' MODSUP (REDU CHAM1 MODSUP);
2347 :        VYO_1 = EXTR CHPO1 'YOUN' PM;
2348 :        VNU_1 = EXTR CHPO1 'NU  ' PM;
2349 :        CHPO1 = CHAN 'CHPO' MODINF (REDU CHAM1 MODINF);
2350 :        VYO_2 = EXTR CHPO1 'YOUN' PM;
2351 :        VNU_2 = EXTR CHPO1 'NU  ' PM;
2352 :        si(flmess);
2353 :          mess ' mat1  mat2';
2354 :          mess  VYO_1 VYO_2;
2355 :          mess  VNU_1 VNU_2;
2356 :        fins;
2357 : *      Constante de Kolosov
2358 :        SI (EGA &MODE 'PLANCONT');
2359 :          KAP_1 = (3. - VNU_1) / (1. + VNU_1);
2360 :          KAP_2 = (3. - VNU_2) / (1. + VNU_2);
2361 :        SINON;
2362 :          KAP_1 = (3. - (4. * VNU_1));
2363 :          KAP_2 = (3. - (4. * VNU_2));
2364 :        FINSI;
2365 : *      Module de cisaillement
2366 :        MU_1 = VYO_1 / (2.*(1. + VNU_1));
2367 :        MU_2 = VYO_2 / (2.*(1. + VNU_2));
2368 : *      Constante bi-metallique EPS1
2369 :        VA1 = (KAP_1/MU_1) + (1./MU_2);
2370 :        VA2 = (KAP_2/MU_2) + (1./MU_1);
2371 :        EPS1 = (1./(2.*VALPI)) * (LOG (VA1/VA2));
2372 : *      Constante C_MATE (= 1 / E^etoile)
2373 :        COSH1 = VALPI * EPS1;
2374 :        COSH1 = ((EXP COSH1) + (EXP (COSH1*(-1.)))) / 2.;
2375 :        VA1 = (MU_1 + (KAP_1*MU_2)) * (MU_2 + (KAP_2*MU_1));
2376 :        VA2 = MU_1 * MU_2 * ((MU_1*(1. + KAP_2)) + (MU_2*(1. + KAP_1)));
2377 :        C_MATE = (COSH1*COSH1*VA1) / (4.*VA2);
2378 :        mess  'EPS1=' EPS1 '  C_MATE=' C_MATE;
2379 :     FINSI;
2380 : *   rappel: G = C_MATE * (K1^2 + K2^2)  +  1/MU * K3^2
2381 : *           M = 2*C_MATE * (K1*K1^aux + K2*K2^aux)  +  2/MU * K3*K3^aux
2382 : 
2383 : *   on evite le cas 3d + EF standard
2384 :     si (non ((EGA &DIME 3) et (non IXFEM)));
2385 : 
2386 : **** CHAMPS POUR SIMPLIFIER L ECRITURE ************
2387 : *    (CHPOINT pour EF std   //   CHAMELEM si XFEM)
2388 : *   -CAS D UN MATERIAU HOMOGENE
2389 :      SI (NON IDANS);
2390 :         SIN1T  = sin TETA1;       COS1T  = cos TETA1;
2391 :         SIN05T = sin (0.5*TETA1); COS05T = cos (0.5*TETA1);
2392 :         si(IXFEM);
2393 :           RM05 = (2.*VALPI*RAY1) ** -0.5 ;
2394 :           COE_GU = RM05 /  (4. * MU_1);
2395 :           SIN15T = sin (1.5*TETA1); COS15T = cos (1.5*TETA1);
2396 : *          KAP_1 = KAP_1 * UN1;
2397 :         sino;
2398 :           COE_1 = ((RAY1/(2.*VALPI))**0.5) / (2.*MU_1);
2399 :         fins;
2400 : *   -CAS D UN BI-MATERIAU
2401 :      SINON;
2402 :        EPSLGR = EPS1 * (LOG RAY1);
2403 :        VA1 = COS (EPSLGR*180./VALPI);
2404 :        VA2 = SIN (EPSLGR*180./VALPI);
2405 :        BTA1   = ((0.5*VA1) + (EPS1*VA2)) / (0.25 + (EPS1*EPS1));
2406 :        BTAPM1 = ((0.5*VA2) - (EPS1*VA1)) / (0.25 + (EPS1*EPS1));
2407 :        DTA_1 = EXP (0. - ((VALPI - TETA1rad)*EPS1));
2408 :        DTA_2 = EXP ((VALPI + TETA1rad)*EPS1);
2409 :        GAM_1 = (KAP_1*DTA_1) - (DTA_1**(-1.));
2410 :        GAM_2 = (KAP_2*DTA_2) - (DTA_2**(-1.));
2411 :        GAMPM_1 = (KAP_1*DTA_1) + (DTA_1**(-1.));
2412 :        GAMPM_2 = (KAP_2*DTA_2) + (DTA_2**(-1.));
2413 :        COS05T = COS (TETA1/2.);  SIN05T = SIN (TETA1/2.);
2414 :        D_1 = (BTA1*GAM_1*COS05T) + (BTAPM1*GAMPM_1*SIN05T);
2415 :        D_2 = (BTA1*GAM_2*COS05T) + (BTAPM1*GAMPM_2*SIN05T);
2416 :        DPM_1 = (BTAPM1*GAM_1*COS05T) - (BTA1*GAMPM_1*SIN05T);
2417 :        DPM_2 = (BTAPM1*GAM_2*COS05T) - (BTA1*GAMPM_2*SIN05T);
2418 :        GTAR1 = EPSLGR + (0.5*TETA1rad);
2419 :        CVA_1 = (SIN TETA1) * (SIN (GTAR1*180./VALPI));
2420 :        CVA_2 = (SIN TETA1) * (COS (GTAR1*180./VALPI));
2421 :        COE_1 = ((RAY1/(2.*VALPI))**0.5) / (4.*MU_1);
2422 :        COE_2 = ((RAY1/(2.*VALPI))**0.5) / (4.*MU_2);
2423 :      FINSI;
2424 : 
2425 :     fins;
2426 : 
2427 :   FINSI;
2428 : 
2429 : 
2430 : 
2431 : 
2432 : * BOUCLE SUR LES INTEGRALES A CALCULER ==============================*
2433 :   REPETER BOUCMIX NBMIXT;
2434 :   IM = IM + 1;
2435 : 
2436 : ****************************************************
2437 : **** CHAMPS AUXILIAIRES SI DECOUPLAGE **************
2438 : * Il existe plusieurs manieres de creer les champs aux :
2439 : * 1. en utilisant l expression analytique de grad(U) et sigma
2440 : * 2. en utilisant l expression analytique de U et en calculant
2441 : *    sigma, Fint=Bsigma, grad(U) ...
2442 : * 3. en appliquant une pression/cisaillement sur les faces de la fissure
2443 : *    et en resolvant le pb associe
2444 : * On utilise la 1ere lorsqu'on peut, la 2eme sinon
2445 :   SI (EGA IINTE 99);
2446 : 
2447 : 
2448 : **** MOTMIX et MOTMIA **************************
2449 :     SI (IM EGA 1); MOTMIX = MOT 'I';  MOTMIA = MOT '  I'; FINSI;
2450 :     SI (IM EGA 2); MOTMIX = MOT 'II'; MOTMIA = MOT ' II'; FINSI;
2451 :     SI (IM EGA 3); MOTMIX = MOT 'III';MOTMIA = MOT 'III'; FINSI;
2452 :   mess 'CHAMPS AUXILIAIRES mode' MOTMIA;
2453 : 
2454 : 
2455 : **** METHODE ANALYTIQUE PURE *******************
2456 : 
2457 : *  -CAS D UN MATERIAU HOMOGENE 2D et 3D --------------------
2458 : *    SI (IXFEM et (EGA &DIME 2) ET (NON IDANS));
2459 :     SI (IXFEM ET (NON IDANS));
2460 :       si(flmess);
2461 :         mess 'MATERIAU HOMOGENE XFEM: METHODE ANALYTIQUE';
2462 :       fins;
2463 : 
2464 : *     DERIVEES REPERE CYLINDRIQUE / COORDONNEES GLOBALES
2465 : *        R,X = DRDX         R,Y = DRDY
2466 : *        T,X = (1/R)*DTDX   T,Y = (1/R)*DTDY
2467 :          DRDX = (      COS1T*(EXCO GLV7 'UX,X' 'SCAL'))
2468 :               + (SDIR1*SIN1T*(EXCO GLV7 'UY,X' 'SCAL')) ;
2469 :          DRDY = (      COS1T*(EXCO GLV7 'UX,Y' 'SCAL'))
2470 :               + (SDIR1*SIN1T*(EXCO GLV7 'UY,Y' 'SCAL')) ;
2471 :          DTDX = (  -1.*SIN1T*(EXCO GLV7 'UX,X' 'SCAL'))
2472 :               + (SDIR1*COS1T*(EXCO GLV7 'UY,X' 'SCAL'));
2473 :          DTDY =  ( -1.*SIN1T*(EXCO GLV7 'UX,Y' 'SCAL'))
2474 :               + (SDIR1*COS1T*(EXCO GLV7 'UY,Y' 'SCAL'));
2475 :          si(EGA &DIME 3);
2476 :            DRDZ = (      COS1T*(EXCO GLV7 'UX,Z' 'SCAL'))
2477 :                 + (SDIR1*SIN1T*(EXCO GLV7 'UY,Z' 'SCAL')) ;
2478 :            DTDZ = (  -1.*SIN1T*(EXCO GLV7 'UX,Z' 'SCAL'))
2479 :                 + (SDIR1*COS1T*(EXCO GLV7 'UY,Z' 'SCAL'));
2480 :          fins;
2481 : 
2482 : *    -debut du cas contact frottant IFROT (btrolle 19/02/2013)
2483 : *     ajout des solutions analytiques du saut sur les lèvres de la fissure
2484 : *     projection de psi1 sur la fissure et calcul de psi,x
2485 :      'SI' IFROT;
2486 : 
2487 : *        PSI en CHAML aux PG sur la fissure utilisé pour calcul du gradient
2488 : *        on sélectionne la partie de psi2 non nulle (= dans le champ theta)
2489 :          SI (EGA &DIME 2); 
2490 : *        contour du domaine theta
2491 :          con1  = 'CONTOUR' ELTETA;
2492 : *        maillage support pour l'intégration
2493 :          mai1 = 'INCLUSION' MAICON con1 'STRI';        
2494 : *          OBJCON2 = MODE mai1 'MECANIQUE' 'ZCO2';
2495 :          OBJCON2 = REDU  OBJCON mai1;
2496 :          'SINON';
2497 : *        maillage support pour l'intégration
2498 :          mai1 = 'INCLUSION' ('EXTRAIRE' OBJCON 'MAIL') ELTETA 'VOLU'
2499 :             'STRI';        
2500 : *          OBJCON2 = 'MODELISER' mai1 'MECANIQUE' 'ZCO3';
2501 :          OBJCON2 = REDU  OBJCON mai1;
2502 :          'FINSI';
2503 :          PSI1e = CHAN 'CHAM' PSI1 OBJMOD 'NOEUDS';
2504 :          PSI2 = PROI OBJCON2 PSI1e 'STRESSES';
2505 : *        PSI en CHPOINT sur la fissure utilisé pour calcul du repère local
2506 :          PSI3 = PROI MAICON PSI1e;        
2507 :           TesPsi = PSI2 'MASQUE' 'INFERIEUR'(-1E-15) ;
2508 : *          PSI2B =  'CHANGER' 'CHAM' (TesPsi * (-1.*PSI2) ('MOTS' 'PSI') ('MOTS'
2509 : *          ('MOTS' 'PSI')) OBJCON 'STRESSES';
2510 :           PSI2B = (TesPsi * (-1.*PSI2) ('MOTS' 'PSI') ('MOTS' 'PSI')
2511 :           ('MOTS' 'PSI'));
2512 : *         'MESSAGE' '  PSI2B =  '; 'LISTE'  PSI2B;
2513 : *        terme sqrt(r/2pI) et sa dérivée,r 
2514 :          RM05B = ((PSI2B/(2.*VALPI)) **0.5);
2515 :          RM05BR =0.5*((2.*VALPI*PSI2B)**-0.5);
2516 : *        Change le nom pour pouvoir calculer grad avec ZCO  
2517 :          LV72  = (PSI3 NOMC 'AX')
2518 :          et (MANU  'CHPO' ('EXTRAIRE' OBJCON2 'MAIL') 1 'AY  ' 0.
2519 :          'NATURE' 'DIFFUS')
2520 :          et (MANU  'CHPO' ('EXTRAIRE' OBJCON2 'MAIL') 1 'AZ  ' 0.
2521 :          'NATURE' 'DIFFUS');
2522 :          GLV72 = 'CHANGER' (GRAD LV72 OBJCON2) 'TYPE' 'SCALAIRE';
2523 : *         'MESSAGE' ' GLV72 =';'LISTE'  GLV72 ;
2524 :          'SI' (&DIME 'EGA' 2);
2525 : *   Angle ALPHA1 de passage local -> global
2526 :          NGPSI2 =   ( ((EXCO GLV72 'AX,X' 'SCAL')**2)
2527 :                + ((EXCO GLV72 'AX,Y' 'SCAL')**2) )**(0.5) ;
2528 :          COS1AB =  ( (EXCO GLV72 'AX,Y' 'SCAL') / NGPSI2) ;
2529 :          SIN1AB =  ( (EXCO GLV72 'AX,X' 'SCAL') / NGPSI2) ;         
2530 :          'SINON';
2531 : *   Angle ALPHA1 de passage local -> global
2532 :          NGPSI2 =   ( ((EXCO GLV72 'AX,X' 'SCAL')**2)
2533 :                + ((EXCO GLV72 'AX,Y' 'SCAL')**2)
2534 :                + ((EXCO GLV72 'AX,Z' 'SCAL')**2) )**(0.5) ;               
2535 : *        CHAM première tangente = grad de psi
2536 :          V12X = ('NOMC' 'UX' ((EXCO GLV72 'AX,X' 'SCAL') '/' NGPSI2));
2537 :          V12Y = ('NOMC' 'UY' ((EXCO GLV72 'AX,Y' 'SCAL') '/' NGPSI2));
2538 :          V12Z = ('NOMC' 'UZ' ((EXCO GLV72 'AX,Z' 'SCAL') '/' NGPSI2));
2539 :          V12 = V12X 'ET' V12Y 'ET' V12Z;
2540 : *        CHAM des normales         
2541 :          V22 =  vsur OBJCON2 'NORM';
2542 :          V22X = 'EXCO' 'VX' V22;
2543 :          V22Y = 'EXCO' 'VY' V22;
2544 :          V22Z = 'EXCO' 'VZ' V22;
2545 : *        CHAM de la dernière tangente produit vect des 2 autres
2546 :          V32X1 = 'NOMC' 'UX' (V12Y * V22Z ('MOTS' 'UY') ('MOTS' 'VZ')
2547 :          ('MOTS' 'UX'));
2548 :          V32X2 =  'NOMC' 'UX'(V12Z * V22Y ('MOTS' 'UZ') ('MOTS' 'VY')
2549 :          ('MOTS' 'UX'));    
2550 :          V32X = V32X1 '-' V32X2;
2551 :          V32Y1 = 'NOMC' 'UY' (V12Z * V22X ('MOTS' 'UZ') ('MOTS' 'VX')
2552 :          ('MOTS' 'UY'));
2553 :          V32Y2 =  'NOMC' 'UY'(V12X * V22Z ('MOTS' 'UX') ('MOTS' 'VZ')
2554 :          ('MOTS' 'UY'));
2555 :          V32Y = V32Y1 '-' V32Y2;
2556 :          V32Z1 = 'NOMC' 'UZ' (V12X * V22Y ('MOTS' 'UX') ('MOTS' 'VY')
2557 :          ('MOTS' 'UZ'));
2558 :          V32Z2 =  'NOMC' 'UZ'(V12Y * V22X ('MOTS' 'UY') ('MOTS' 'VX')
2559 :          ('MOTS' 'UZ'));        
2560 :          V32Z = V32Z1 '-' V32Z2;
2561 :          V32 = V32X 'ET' V32Y 'ET' V32Z;         
2562 :     ROT2= (EXCO V12 (mots 'UX' 'UY' 'UZ') (mots 'AX,X' 'AY,X' 'AZ,X'))
2563 :        ET (EXCO V22 (mots 'VX' 'VY' 'VZ') (mots 'AX,Y' 'AY,Y' 'AZ,Y'))
2564 :        ET (EXCO V32 (mots 'UX' 'UY' 'UZ') (mots 'AX,Z' 'AY,Z' 'AZ,Z'));
2565 : *         ROT2 = EXCO ROT2
2566 : *         (mots AX,X AX,Y AX,Z AY,X AY,Y AY,Z AZ,X AZ,Y AZ,Z);
2567 : *         ROT2 = CHAN 'CHAM' ROT2 OBJCON 'STRESSES' 'GRADIENT';
2568 : *         'MESSAGE' 'ROT2 = ';
2569 : *        'LISTE' ROT2;
2570 :          'FINSI';
2571 :          
2572 : *        r,x (r = -PSI sur la fissure)
2573 :          dRdX2 = (EXCO GLV72 'AX,X' 'SCAL');
2574 :          dRdY2 = (EXCO GLV72 'AX,Y' 'SCAL');
2575 :          si(EGA &DIME 3);
2576 :             dRdZ2 = (EXCO GLV72 'AX,Z' 'SCAL');
2577 :          'FINSI';
2578 : 
2579 :       'FINSI';
2580 : *    -fin du cas contact frottant IFROT (btrolle 19/02/2013)
2581 :          
2582 : *     CHAMP AUX. MODE 1
2583 :       SI (IM EGA 1);
2584 : *        derivees elementaires / repere local
2585 : *        Ui,R = (1/4µ)*RM05*UiR              i=1,2
2586 : *        Ui,T = (1/4µ)*RM05*UiT              i=1,2
2587 :          U1R = ((KAP_1 - 0.5) * COS05T) - (0.5 * COS15T);
2588 :          U1T = ((0.5 - KAP_1) * SIN05T) + (1.5 * SIN15T);
2589 :          U2R = ((KAP_1 + 0.5) * SIN05T) - (0.5 * SIN15T);
2590 :          U2T = ((KAP_1 + 0.5) * COS05T) - (1.5 * COS15T);
2591 : *        contraintes dans le repere local
2592 : *        SIGij        ij={11,22,33,12}
2593 :          SIG11 = RM05 * ( COS05T * (UN1 - (SIN05T*SIN15T)) );
2594 :          SIG22 = RM05 * ( COS05T * (UN1 + (SIN05T*SIN15T)) );
2595 :          SIG33 = ZER1;
2596 :          SIG12 = RM05 * ( COS05T * (SIN05T * COS15T) ) ;
2597 : *        gradient de deplacement dans le repere local
2598 : *        COE_GU = (1/4µ)*RM05
2599 : *        U1,X  U1,Y  U2,X  U2,Y
2600 :          GU1X = COE_GU * ( ( U1R * DRDX )  + ( U1T * DTDX ) );
2601 :          GU1Y = COE_GU * ( ( U1R * DRDY )  + ( U1T * DTDY ) );
2602 :          GU2X = COE_GU * ( ( U2R * DRDX )  + ( U2T * DTDX ) );
2603 :          GU2Y = COE_GU * ( ( U2R * DRDY )  + ( U2T * DTDY ) );
2604 : *       -debut du cas contact frottant IFROT (btrolle 19/02/2013)
2605 : *        saut en ouverture
2606 :          'SI' IFROT;
2607 :          W2 = 8.*RM05B *C_MATE;
2608 : *        w2,r
2609 : *         W2R = ('CHANGER' 'CHAM' (8*C_MATE * RM05BR) OBJCON 'STRESSES');
2610 :          W2R =  (8*C_MATE * RM05BR);
2611 :          W2R = 'NOMC' W2R 'SCAL';
2612 : *        gradient saut ouverture mis à 0 car on ne veut pas modifier KI en prése
2613 :          W2X = 0.*W2R*dRdX2;
2614 :          W2Y = 0.*W2R*dRdY2;
2615 :          W1X = 0.*W2X; W1Y = W1X;
2616 :          W3X = 0.*W2X; W3Y = W3X;
2617 :          'FINSI';
2618 : *       -fin du cas contact frottant IFROT (btrolle 19/02/2013)
2619 :  
2620 :          
2621 :          si(EGA &DIME 3);
2622 :            GU1Z = COE_GU * ( ( U1R * DRDZ )  + ( U1T * DTDZ ) );
2623 :            GU2Z = COE_GU * ( ( U2R * DRDZ )  + ( U2T * DTDZ ) );
2624 :            SIG13 = ZER1; SIG23 = ZER1;
2625 :            GU3X = ZER1;   GU3Y = ZER1;    GU3Z = ZER1;
2626 :            'SI' IFROT;
2627 :              W2Z = 0.*W2R*dRdZ2;
2628 :              W1Z = W1X;W3Z = W3X;
2629 :            'FINSI';
2630 :          fins;
2631 :       FINSI;
2632 : 
2633 : *     CHAMP AUX. MODE 2
2634 :       SI (IM EGA 2);
2635 : *        derivees elementaires / repere local
2636 : *        Ui,R = (1/4µ)*RM05*UiR              i=1,2
2637 : *        Ui,T = (1/4µ)*RM05*UiT              i=1,2
2638 :          U1R = ((KAP_1 + 1.5) * SIN05T) + (0.5 * SIN15T);
2639 :          U1T = ((KAP_1 + 1.5) * COS05T) + (1.5 * COS15T);
2640 :          U2R = ((1.5 - KAP_1) * COS05T) - (0.5 * COS15T);
2641 :          U2T = ((KAP_1 - 1.5) * SIN05T) + (1.5 * SIN15T);
2642 : *        contraintes dans le repere local
2643 : *        SIGij        ij={11,22,33,12}
2644 :          SIG11 = -1.*RM05 * ( SIN05T * ((2.*UN1) + (COS05T*COS15T)) );
2645 :          SIG22 =     RM05 * ( SIN05T * (COS05T*COS15T) ) ;
2646 :          SIG33 = ZER1;
2647 :          SIG12 = RM05 * ( COS05T * (UN1 - (SIN05T * SIN15T)) ) ;
2648 : *        gradient de deplacement dans le repere local
2649 : *        COE_GU = (1/4µ)*RM05
2650 : *        U1,X  U1,Y  U2,X  U2,Y
2651 :          GU1X = COE_GU * ( ( U1R * DRDX )  + ( U1T * DTDX ) );
2652 :          GU1Y = COE_GU * ( ( U1R * DRDY )  + ( U1T * DTDY ) );
2653 :          GU2X = COE_GU * ( ( U2R * DRDX )  + ( U2T * DTDX ) );
2654 :          GU2Y = COE_GU * ( ( U2R * DRDY )  + ( U2T * DTDY ) );
2655 : 
2656 :          
2657 : *       -debut du cas contact frottant IFROT (btrolle 19/02/2013)
2658 : *        saut cisaillement
2659 :          'SI' IFROT;
2660 :          W1 = 8.*RM05B *C_MATE;
2661 : *        w1,r
2662 : *         W1R = ('CHANGER' 'CHAM' (8*C_MATE * RM05BR) OBJCON'STRESSES');
2663 : *        On divise par 2 car on saut défini comme plus - moyenne et pas (lèvre +
2664 :          W1R = (8.*C_MATE * RM05BR) '/' 2.;
2665 :          W1R = 'NOMC' W1R 'SCAL'; 
2666 : *        gradient saut cisaillement = W1,r r,X
2667 :          W1X = W1R*dRdX2; 
2668 :          W1Y = W1R*dRdY2;
2669 :          W2X = 0.*W1X; W2Y = 0.*W2X;
2670 :          W3X = 0.*W1X; W3Y = 0.*W3X;
2671 :          'FINSI';
2672 : *       -fin du cas contact frottant IFROT (btrolle 19/02/2013)
2673 :           
2674 :          si(EGA &DIME 3);
2675 :            GU1Z = COE_GU * ( ( U1R * DRDZ )  + ( U1T * DTDZ ) );
2676 :            GU2Z = COE_GU * ( ( U2R * DRDZ )  + ( U2T * DTDZ ) );
2677 :            SIG13 = ZER1; SIG23 = ZER1;
2678 :            GU3X = ZER1;   GU3Y = ZER1;    GU3Z = ZER1;
2679 :            'SI' IFROT;
2680 :                W1Z = W1R*dRdZ2;
2681 :                W2Z = 0.*W2X;W3Z = 0.*W3X;
2682 :            'FINSI';
2683 :          fins;
2684 :       FINSI;
2685 : 
2686 : *     CHAMP AUX. MODE 3  (automatiquement on a : EGA &DIME 3)
2687 :       SI (IM EGA 3);
2688 : *        derivees elementaires / repere local
2689 : *        Ui,R = (1/4µ)*RM05*UiR              i=3
2690 : *        Ui,T = (1/4µ)*RM05*UiT              i=3
2691 :          U3R = SIN05T;
2692 :          U3T = COS05T;
2693 : *        contraintes dans le repere local
2694 : *        SIGij        ij={11,22,33,12}
2695 :          SIG11 = ZER1;
2696 :          SIG22 = ZER1;
2697 :          SIG33 = ZER1;
2698 :          SIG12 = ZER1;
2699 :          SIG13 = -1.*RM05 * SIN05T;
2700 :          SIG23 =     RM05 * COS05T;
2701 : *        gradient de deplacement dans le repere local
2702 : *        COE_GU = (1/4µ)*RM05
2703 : *        U1,X  U1,Y  U2,X  U2,Y
2704 :          GU1X = ZER1; GU1Y = ZER1; GU1Z = ZER1;
2705 :          GU2X = ZER1; GU2Y = ZER1; GU2Z = ZER1;
2706 :          GU3X = COE_GU * ( ( U3R * DRDX )  + ( U3T * DTDX ) );
2707 :          GU3Y = COE_GU * ( ( U3R * DRDY )  + ( U3T * DTDY ) );
2708 :          GU3Z = COE_GU * ( ( U3R * DRDZ )  + ( U3T * DTDZ ) );
2709 : 
2710 : *    -fin du cas contact frottant IFROT (btrolle 19/02/2013)
2711 : *        saut cisaillement antiplan
2712 :          'SI' IFROT;
2713 :          W3 = 8.*RM05B *C_MATE/ (1. - VNU_1) ;
2714 : *        w3,r
2715 : *         W3R = ('CHANGER' 'CHAM' (8*C_MATE * RM05BR'/'(1. - VNU_1))
2716 : *                     OBJCON 'STRESSES');
2717 : *        On divise par 2 car on saut défini comme plus - moyenne et pas (lèvre +
2718 :          W3R =  (8*C_MATE * RM05BR'/'(1. - VNU_1)) '/' 2.;
2719 :          W3R = 'NOMC' W3R 'SCAL';
2720 : *        gradient saut cisaillement antiplan = W3,r r,X
2721 :          W3X = W3R*dRdX2;
2722 :          W3Y = W3R*dRdY2;
2723 :          W3Z = W3R*dRdZ2;
2724 :          W1X = 0.*W3X; W1Y = 0.*W1X; W1Z = 0.*W1X;
2725 :          W2X = 0.*W3X; W2Y = 0.*W2X; W2Z = 0.*W2X;
2726 :          'FINSI';
2727 : *    -fin du cas contact frottant IFROT (btrolle 19/02/2013)
2728 : 
2729 :        FINSI;
2730 : 
2731 : *     PASSAGE DANS LE REPERE GLOBAL
2732 :       si(EGA &DIME 2);
2733 : *     ... des contraintes
2734 : *     SIGIJ = [P] * SIGij * [P]**-1   local: ij={1,2..}  global: IJ={X,Y..}
2735 :       SIGXX= (COSA2*SIG11) - (2.*SINCOSA*SIG12) + (SINA2*SIG22);
2736 :       SIGYY= (SINA2*SIG11) + (2.*SINCOSA*SIG12) + (COSA2*SIG22);
2737 :       SIGZZ= SIG33;
2738 :       SIGXY= (SINCOSA*SIG11) + ((COSA2-SINA2)*SIG12) - (SINCOSA*SIG22);
2739 :       SIGAUX = (NOMC 'SMXX' SIGXX) ET  (NOMC 'SMYY' SIGYY)
2740 :             ET (NOMC 'SMZZ' SIGZZ) ET  (NOMC 'SMXY' SIGXY);
2741 :       SIGAUX = CHAN SIGAUX 'TYPE' 'CONTRAINTES';
2742 : *      TRAC SIGAUX OBJMOD 'TITR' ' SIGAUX';
2743 : *     ... des gradient de deplacement
2744 : *     UI,J = [P] * Ui,j       local: ij={1,2..}  global: IJ={X,Y..}
2745 :       GRUXX  = (GU1X*COS1A) - (GU2X*SIN1A);
2746 :       GRUXY  = (GU1Y*COS1A) - (GU2Y*SIN1A);
2747 :       GRUYX  = (GU1X*SIN1A) + (GU2X*COS1A);
2748 :       GRUYY  = (GU1Y*SIN1A) + (GU2Y*COS1A);
2749 :       GRUAUX = (NOMC 'UX,X' GRUXX) ET (NOMC 'UX,Y' GRUXY)
2750 :             ET (NOMC 'UX,Z' ZER1)
2751 :             ET (NOMC 'UY,X' GRUYX) ET (NOMC 'UY,Y' GRUYY)
2752 :             ET (NOMC 'UY,Z' ZER1)
2753 :             ET (NOMC 'UZ,X' ZER1) ET (NOMC 'UZ,Y' ZER1)
2754 :             ET (NOMC 'UZ,Z' ZER1)       ;
2755 :       GRUAUX = CHAN GRUAUX 'TYPE' 'GRADIENT';
2756 : *      TRAC GRUAUX OBJMOD 'TITR' ' GRUAUX';
2757 : *    -debut du cas contact frottant IFROT (btrolle 19/02/2013)
2758 :       'SI' IFROT;
2759 : *     ... des gradient de deplacement sur la fissure
2760 : *     UI,J = [P] * Ui,j       local: ij={1,2..}  global: IJ={X,Y..}
2761 :       GRWXX  = (W1X*COS1AB) - (W2X*SIN1AB);
2762 :       GRWXY  = (W1Y*COS1AB) - (W2Y*SIN1AB);
2763 :       GRWYX  = (W1X*SIN1AB) + (W2X*COS1AB);
2764 :       GRWYY  = (W1Y*SIN1AB) + (W2Y*COS1AB);
2765 :       GRWAUX = (NOMC 'AX,X' GRWXX) ET (NOMC 'AX,Y' GRWXY)
2766 :             ET (NOMC 'AX,Z' ZER1)
2767 :             ET (NOMC 'AY,X' GRWYX) ET (NOMC 'AY,Y' GRWYY)
2768 :             ET (NOMC 'AY,Z' ZER1)
2769 :             ET (NOMC 'AZ,X' ZER1) ET (NOMC 'AZ,Y' ZER1)
2770 :             ET (NOMC 'AZ,Z' ZER1)       ;
2771 :       GRWAUX = CHAN GRWAUX 'TYPE' 'GRADIENT';
2772 :       'FINSI';
2773 : *    -fin du cas contact frottant IFROT (btrolle 19/02/2013)
2774 : 
2775 :       sino;
2776 :       SIGAUX = (NOMC 'SMXX' SIG11) ET  (NOMC 'SMYY' SIG22)
2777 :             ET (NOMC 'SMZZ' SIG33) ET  (NOMC 'SMXY' SIG12)
2778 :             ET (NOMC 'SMXZ' SIG13) ET  (NOMC 'SMYZ' SIG23);
2779 :       SIGAUX = CHAN SIGAUX 'TYPE' 'CONTRAINTES';
2780 : *       SUPTAB . (chai 'SIGAUX_LOCAL' IM) = SIGAUX;
2781 :       SIGAUX = RTEN SIGAUX OBJMOD ROT1 'RART';
2782 : *       SUPTAB . (chai 'SIGAUX_GLOBAL' IM) = SIGAUX;
2783 : GUiJ = (nomc GU1X 'U1,X') et (nomc GU1Y 'U1,Y') et (nomc GU1Z 'U1,Z')
2784 :     et (nomc GU2X 'U2,X') et (nomc GU2Y 'U2,Y') et (nomc GU2Z 'U2,Z')
2785 :     et (nomc GU3X 'U3,X') et (nomc GU3Y 'U3,Y') et (nomc GU3Z 'U3,Z');
2786 : 
2787 :       GRUXX  = PSCA ROT1 GUiJ
2788 :        (mots 'UX,X' 'UX,Y' 'UX,Z') (mots 'U1,X' 'U2,X' 'U3,X');
2789 :       GRUYX  = PSCA ROT1 GUiJ
2790 :        (mots 'UY,X' 'UY,Y' 'UY,Z') (mots 'U1,X' 'U2,X' 'U3,X');
2791 :       GRUZX  = PSCA ROT1 GUiJ
2792 :        (mots 'UZ,X' 'UZ,Y' 'UZ,Z') (mots 'U1,X' 'U2,X' 'U3,X');
2793 :       GRUXY  = PSCA ROT1 GUiJ
2794 :        (mots 'UX,X' 'UX,Y' 'UX,Z') (mots 'U1,Y' 'U2,Y' 'U3,Y');
2795 :       GRUYY  = PSCA ROT1 GUiJ
2796 :        (mots 'UY,X' 'UY,Y' 'UY,Z') (mots 'U1,Y' 'U2,Y' 'U3,Y');
2797 :       GRUZY  = PSCA ROT1 GUiJ
2798 :        (mots 'UZ,X' 'UZ,Y' 'UZ,Z') (mots 'U1,Y' 'U2,Y' 'U3,Y');
2799 :       GRUXZ  = PSCA ROT1 GUiJ
2800 :        (mots 'UX,X' 'UX,Y' 'UX,Z') (mots 'U1,Z' 'U2,Z' 'U3,Z');
2801 :       GRUYZ  = PSCA ROT1 GUiJ
2802 :        (mots 'UY,X' 'UY,Y' 'UY,Z') (mots 'U1,Z' 'U2,Z' 'U3,Z');
2803 :       GRUZZ  = PSCA ROT1 GUiJ
2804 :        (mots 'UZ,X' 'UZ,Y' 'UZ,Z') (mots 'U1,Z' 'U2,Z' 'U3,Z');
2805 : 
2806 :       GRUAUX = (NOMC 'UX,X' GRUXX) ET (NOMC 'UX,Y' GRUXY)
2807 :             ET (NOMC 'UX,Z' GRUXZ)
2808 :             ET (NOMC 'UY,X' GRUYX) ET (NOMC 'UY,Y' GRUYY)
2809 :             ET (NOMC 'UY,Z' GRUYZ)
2810 :             ET (NOMC 'UZ,X' GRUZX) ET (NOMC 'UZ,Y' GRUZY)
2811 :             ET (NOMC 'UZ,Z' GRUZZ)       ;
2812 :       GRUAUX = CHAN GRUAUX 'TYPE' 'GRADIENT';
2813 : *       TRAC GRUAUX OBJMOD 'TITR' ' GRUAUX' (100. 0. 0.);
2814 : *       TRAC GRUAUX OBJMOD 'TITR' ' GRUAUX' (-100. 0. 0.);
2815 : *    -debut du cas contact frottant IFROT (btrolle 19/02/2013)
2816 :       'SI' IFROT;
2817 : GWiJ = (nomc W1X 'A1,X') et (nomc W1Y 'A1,Y') et (nomc W1Z 'A1,Z')
2818 :     et (nomc W2X 'A2,X') et (nomc W2Y 'A2,Y') et (nomc W2Z 'A2,Z')
2819 :     et (nomc W3X 'A3,X') et (nomc W3Y 'A3,Y') et (nomc W3Z 'A3,Z');
2820 :       GRWXX  = PSCA ROT2 GWiJ
2821 :        (mots 'AX,X' 'AX,Y' 'AX,Z') (mots 'A1,X' 'A2,X' 'A3,X');
2822 :       GRWYX  = PSCA ROT2 GWiJ
2823 :        (mots 'AY,X' 'AY,Y' 'AY,Z') (mots 'A1,X' 'A2,X' 'A3,X');
2824 :       GRWZX  = PSCA ROT2 GWiJ
2825 :        (mots 'AZ,X' 'AZ,Y' 'AZ,Z') (mots 'A1,X' 'A2,X' 'A3,X');
2826 :       GRWXY  = PSCA ROT2 GWiJ
2827 :        (mots 'AX,X' 'AX,Y' 'AX,Z') (mots 'A1,Y' 'A2,Y' 'A3,Y');
2828 :       GRWYY  = PSCA ROT2 GWiJ
2829 :        (mots 'AY,X' 'AY,Y' 'AY,Z') (mots 'A1,Y' 'A2,Y' 'A3,Y');
2830 :       GRWZY  = PSCA ROT2 GWiJ
2831 :        (mots 'AZ,X' 'AZ,Y' 'AZ,Z') (mots 'A1,Y' 'A2,Y' 'A3,Y');
2832 :       GRWXZ  = PSCA ROT2 GWiJ
2833 :        (mots 'AX,X' 'AX,Y' 'AX,Z') (mots 'A1,Z' 'A2,Z' 'A3,Z');
2834 :       GRWYZ  = PSCA ROT2 GWiJ
2835 :        (mots 'AY,X' 'AY,Y' 'AY,Z') (mots 'A1,Z' 'A2,Z' 'A3,Z');
2836 :       GRWZZ  = PSCA ROT2 GWiJ
2837 :        (mots 'AZ,X' 'AZ,Y' 'AZ,Z') (mots 'A1,Z' 'A2,Z' 'A3,Z');
2838 :       GRWAUX = (NOMC 'AX,X' GRWXX) ET (NOMC 'AX,Y' GRWXY)
2839 :             ET (NOMC 'AX,Z' GRWXZ)
2840 :             ET (NOMC 'AY,X' GRWYX) ET (NOMC 'AY,Y' GRWYY)
2841 :             ET (NOMC 'AY,Z' GRWYZ)
2842 :             ET (NOMC 'AZ,X' GRWZX) ET (NOMC 'AZ,Y' GRWZY)
2843 :             ET (NOMC 'AZ,Z' GRWZZ)       ;
2844 :     GRWAUX = CHAN GRWAUX 'TYPE' 'GRADIENT';
2845 :       'FINSI';
2846 : *    -fin du cas contact frottant IFROT (btrolle 19/02/2013)
2847 :       fins;
2848 : *     champs auxiliaires
2849 :       A_DEPI = DEP000;
2850 :       A_DEPGR= GRUAUX;
2851 :       A_SIGF = SIGAUX;
2852 :       A_PREI = FOR000;
2853 :       'SI' IFROT; B_DEPGR = GRWAUX; 'FINSI';
2854 :     FINSI;
2855 : 
2856 : 
2857 : **** METHODE U-ANALYTIQUE **************************
2858 : 
2859 : *  -CAS D UN MATERIAU HOMOGENE 2D --------------------
2860 :     SI ((non IXFEM) et (EGA &DIME 2) ET (NON IDANS));
2861 :       si(flmess);
2862 :         mess 'MATERIAU HOMOGENE 2D: METHODE U-ANALYTIQUE';
2863 :       fins;
2864 : *     champ aux. mode 1
2865 :       SI (IM EGA 1);
2866 :          UX_1 = COE_1 * COS05T * (KAP_1 - COS1T);
2867 :          UY_1 = COE_1 * SIN05T * (KAP_1 - COS1T);
2868 :       FINSI;
2869 : *     champ aux. mode 2
2870 :       SI (IM EGA 2);
2871 :          UX_1 = COE_1 * SIN05T * (KAP_1 + 2. + COS1T);
2872 :          UY_1 = (-1.)*COE_1 * COS05T * (KAP_1 - 2. + COS1T);
2873 :       FINSI;
2874 : *     deplacement dans le repère local
2875 :       UL1 = CHAN 'ATTRIBUT' UX_1 'NATURE' 'DIFFUS';
2876 :       UL2 = CHAN 'ATTRIBUT' UY_1 'NATURE' 'DIFFUS';
2877 : *     retour du repère local vers le global
2878 :       UG1 = (UL1*(COS ALPHA1)) - (UL2*(SIN ALPHA1));
2879 :       UG2 = (UL1*(SIN ALPHA1)) + (UL2*(COS ALPHA1));
2880 :       UG1 = CHAN 'ATTRIBUT' UG1 'NATURE' 'DIFFUS';
2881 :       UG2 = CHAN 'ATTRIBUT' UG2 'NATURE' 'DIFFUS';
2882 : *     champs auxiliaires
2883 :       A_DEPI = (('NOMC' UG1 MU1) ET ('NOMC' UG2 MU2)) + DEP000;
2884 :       A_SIGF = SIGM MAT1 OBJMOD A_DEPI;
2885 :       A_PREI = BSIG A_SIGF OBJMOD;
2886 :     FINSI;
2887 : 
2888 : *  -CAS D UN BI-MATERIAU 2D --------------------
2889 :     SI ((EGA &DIME 2) ET (IDANS));
2890 :       si(flmess);
2891 :         mess 'BI-MATERIAU 2D: METHODE ANALYTIQUE';
2892 :       fins;
2893 : *     champ aux. mode 1
2894 :       SI (IM EGA 1);
2895 :          UX_1 = COE_1*(D_1 + (2.*DTA_1*CVA_1));
2896 :          UX_2 = COE_2*(D_2 + (2.*DTA_2*CVA_1));
2897 :          UY_1 = (-1.)*COE_1*(DPM_1 + (2.*DTA_1*CVA_2));
2898 :          UY_2 = (-1.)*COE_2*(DPM_2 + (2.*DTA_2*CVA_2));
2899 :       FINSI;
2900 : *     champ aux. mode 2
2901 :       SI (IM EGA 2);
2902 :          UX_1 = (-1.)*COE_1*(DPM_1 - (2.*DTA_1*CVA_2));
2903 :          UX_2 = (-1.)*COE_2*(DPM_2 - (2.*DTA_2*CVA_2));
2904 :          UY_1 = (-1.)*COE_1*(D_1 - (2.*DTA_1*CVA_1));
2905 :          UY_2 = (-1.)*COE_2*(D_2 - (2.*DTA_2*CVA_1));
2906 :       FINSI;
2907 : *     deplacement dans le repère local : on moyenne sur l interface
2908 :       UX_1 = CHAN 'ATTRIBUT' UX_1 'NATURE' 'DIFFUS';
2909 :       UX_2 = CHAN 'ATTRIBUT' UX_2 'NATURE' 'DIFFUS';
2910 :       UY_1 = CHAN 'ATTRIBUT' UY_1 'NATURE' 'DIFFUS';
2911 :       UY_2 = CHAN 'ATTRIBUT' UY_2 'NATURE' 'DIFFUS';
2912 :       UL1 = UX_1 ET UX_2;
2913 :       UL2 = UY_1 ET UY_2;
2914 : *     retour du repère local vers le global
2915 :       UG1 = (UL1*(COS ALPHA1)) - (UL2*(SIN ALPHA1));
2916 :       UG2 = (UL1*(SIN ALPHA1)) + (UL2*(COS ALPHA1));
2917 :       UG1 = CHAN 'ATTRIBUT' UG1 'NATURE' 'DIFFUS';
2918 :       UG2 = CHAN 'ATTRIBUT' UG2 'NATURE' 'DIFFUS';
2919 : *     champs auxiliaires
2920 :       A_DEPI = (('NOMC' UG1 MU1) ET ('NOMC' UG2 MU2)) + DEP000;
2921 :       A_SIGF = SIGM MAT1 OBJMOD A_DEPI;
2922 :       A_PREI = BSIG A_SIGF OBJMOD;
2923 :     FINSI;
2924 : 
2925 : 
2926 : **** METHODE MECANIQUE **************************
2927 : *** (EFFORTS APPLIQUES AUX LEVRES DE LA FISSURE)
2928 : * BP : cette methode est tres couteuse car appel a resou pour chaque
2929 : * noeud du front de fissure -> passer a une methode analytique + tard
2930 : 
2931 : *  -CAS D UN MATERIAU HOMOGENE 3D  --------------------
2932 :     SI ((non IXFEM) et (EGA &DIME 3) ET (NON IDANS));
2933 :       mess 'MATERIAU HOMOGENE 3D: METHODE MECANIQUE';
2934 :       LSUP = SUPTAB.'LEVRE_SUPERIEURE';
2935 :       LINF = SUPTAB.'LEVRE_INFERIEURE';
2936 :       VCISA = RESU DIRCISA;
2937 :       VCISA = (maxi (exco VCISA 'UX')) (maxi (exco VCISA 'UY'))
2938 :               (maxi (exco VCISA 'UZ'));
2939 :       VTETA = RESU DIRTETA;
2940 :       VTETA = (maxi (exco VTETA 'UX')) (maxi (exco VTETA 'UY'))
2941 :               (maxi (exco VTETA 'UZ'));
2942 :       VNORM = RESU DIRNORM;
2943 :       VNORM = (maxi (exco VNORM 'UX')) (maxi (exco VNORM 'UY'))
2944 :               (maxi (exco VNORM 'UZ'));
2945 :       SI (EGA &DIME 3);
2946 :          LSUP = LSUP DIFF (LSUP ELEM 'APPU' 'STRI' ELTETA);
2947 :          LINF = LINF DIFF (LINF ELEM 'APPU' 'STRI' ELTETA);
2948 :       FINSI;
2949 : *bp:  pas bien compris sur quelle partie de LSUP on applique les efforts ?
2950 :       SI IPAP;
2951 :          F11=PRES 'MASS' (SUPTAB.'SOLUTION_PASAPAS'.'MODELE') LSUP 1;
2952 :          F21=PRES 'MASS' (SUPTAB.'SOLUTION_PASAPAS'.'MODELE') LINF 1;
2953 :       SINON;
2954 :          F11 = PRES 'MASS' (SUPTAB.'MODELE') LSUP 1;
2955 :          F21 = PRES 'MASS' (SUPTAB.'MODELE') LINF 1;
2956 :       FINSI;
2957 : *       SI ((IM EGA 1) OU (IM EGA 2));
2958 : *          BLOQ1 = BLOQ0 ET (BLOQ 'DEPL' 'DIRECTION' DIRCISA MAILLAGE);
2959 : *       FINSI;
2960 : *       SI (IM EGA 3);
2961 : *         BLOQ1 = BLOQ0 ET (BLOQ 'DEPL' 'DIRECTION' DIRTETA MAILLAGE);
2962 : *       FINSI;
2963 : *       BP : on BLOQUE toutes les directions + haut
2964 : *     champ de pression aux. mode 1
2965 :       SI (IM EGA 1);
2966 :          A_PREI = F11 + F21;
2967 :       FINSI;
2968 : *     champ de pression aux. mode 2
2969 :       SI (IM EGA 2);
2970 :          V1 = EXCO (F11 + F21) MF1 'SCAL';
2971 :          V2 = EXCO (F11 + F21) MF2 'SCAL';
2972 :          V3 = EXCO (F11 + F21) MF3 'SCAL';
2973 :          N1 = ((V1 * V1) + (V2 * V2) + (V3 * V3)) ** 0.5;
2974 : *          V1 = COOR 1 DIRTETA;
2975 : *          V2 = COOR 2 DIRTETA;
2976 : *          V3 = COOR 3 DIRTETA;
2977 : *        adaptation car DIRTETA est désormais un CHPOINT
2978 :          V1 = COOR 1 VTETA;
2979 :          V2 = COOR 2 VTETA;
2980 :          V3 = COOR 3 VTETA;
2981 :          A_PREI = (MANU 'CHPO' LSUP 3 MF1 V1 MF2 V2 MF3 V3) -
2982 :                   (MANU 'CHPO' LINF 3 MF1 V1 MF2 V2 MF3 V3);
2983 : * *        adaptation car DIRTETA est désormais un CHPOINT
2984 : *          A_PREI = EXCO DIRTETA (mots MU1 MU2 MU3) (mots MF1 MF2 MF3);
2985 : *          A_PREI = (REDU A_PREI LSUP) - (REDU A_PREI LINF);
2986 :          A_PREI = A_PREI * N1;
2987 :       FINSI;
2988 : *     champ de pression aux. mode 3
2989 :       SI (IM EGA 3);
2990 :          C_MATE = (1. + VNU_1) / VYO_1;
2991 :          V1 = EXCO (F11 + F21) MF1 'SCAL';
2992 :          V2 = EXCO (F11 + F21) MF2 'SCAL';
2993 :          V3 = EXCO (F11 + F21) MF3 'SCAL';
2994 :          N1 = ((V1 * V1) + (V2 * V2) + (V3 * V3)) ** 0.5;
2995 : *          V1 = COOR 1 DIRCISA;
2996 : *          V2 = COOR 2 DIRCISA;
2997 : *          V3 = COOR 3 DIRCISA;
2998 : *        adaptation car DIRCISA est désormais un CHPOINT
2999 :          V1 = COOR 1 VCISA;
3000 :          V2 = COOR 2 VCISA;
3001 :          V3 = COOR 3 VCISA;
3002 :          A_PREI = (MANU 'CHPO' LSUP 3 MF1 V1 MF2 V2 MF3 V3) -
3003 :                   (MANU 'CHPO' LINF 3 MF1 V1 MF2 V2 MF3 V3);
3004 : * *        adaptation car DIRCISA est désormais un CHPOINT
3005 : *          A_PREI = EXCO DIRCISA (mots MU1 MU2 MU3) (mots MF1 MF2 MF3);
3006 : *          A_PREI = (REDU A_PREI LSUP) - (REDU A_PREI LINF);
3007 :          A_PREI = A_PREI * N1;
3008 :       FINSI;
3009 : *     le deplacement et les contraintes sont deduites de la pression
3010 : *      trac (vect A_PREI 'FORC' 'ROUG') elteta;
3011 :       A_DEPI = RESO (RIGTOT ET BLOQ1) A_PREI;
3012 :       A_SIGF = SIGM MAT1 OBJMOD A_DEPI;
3013 :     FINSI;
3014 : 
3015 : *  -CAS D UN BI-MATERIAU 3D  --------------------
3016 :     SI ((EGA &DIME 3) ET (IDANS));
3017 :      MESS 'ERREUR : ON NE PEUT ENCORE DECOUPLER LES MODES';
3018 :      MESS '         DANS LE CAS DES MATERIAUX COMPOSITES';
3019 :      MESS '         EN 3D';
3020 :      QUIT G_THETA;
3021 :     FINSI;
3022 : 
3023 : *     A_PREI = REDU A_PREI ELTETA;
3024 : *     A_DEPI = REDU A_DEPI ELTETA;
3025 : *     A_SIGF = REDU A_SIGF OBJMOD;
3026 : 
3027 : 
3028 : *   si(fltrac);
3029 : * *    trac A_SIGF OBJMOD (defo A_DEPI ELTETA)
3030 : *     trac A_SIGF OBJMOD
3031 : *      'TITR' (chai 'sigma^aux_' IM ' (analytique)');
3032 : * *     trac (vect A_PREI 'FORC' 'ROUG')  ELTETA
3033 : * *      'TITR' (chai 'sigma^aux_' IM ' (analytique)');
3034 : *   fins;
3035 : 
3036 :   FINSI;
3037 : **FIN DU CALCUL DES CHAMP AUX.**********************
3038 : 
3039 : 
3040 : ****************************************************
3041 : ******* EN CAS DE CALCUL EN VISCO_PLASTICITE *******
3042 : ****************************************************
3043 :   SI (((EGA IINTE 2) OU (EGA IINTE 3)) ET (EGA ITYPEF 1));
3044 :      CHAF1 = (CHAN 'STRESSES' OBJMOD (EXCO 'AF1 ' MAT1 'SCAL'))
3045 :                 CHAN 'TYPE' 'SCALAIRE';
3046 :      CHAF2 = (CHAN 'STRESSES' OBJMOD (EXCO 'AF2 ' MAT1 'SCAL'))
3047 :                 CHAN 'TYPE' 'SCALAIRE';
3048 :      CHAF3 = (CHAN 'STRESSES' OBJMOD (EXCO 'AF3 ' MAT1 'SCAL'))
3049 :                 CHAN 'TYPE' 'SCALAIRE';
3050 :     SI (EXIS MAT1 'AF0 ');
3051 :        CHAF4 =( CHAN 'STRESSES' OBJMOD (EXCO 'AF4 ' MAT1 'SCAL'))
3052 :                 CHAN 'TYPE' 'SCALAIRE';
3053 :        CHAF5 =( CHAN 'STRESSES' OBJMOD (EXCO 'AF5 ' MAT1 'SCAL'))
3054 :                 CHAN 'TYPE' 'SCALAIRE';
3055 :        CHAF6 = (CHAN 'STRESSES' OBJMOD (EXCO 'AF6 ' MAT1 'SCAL'))
3056 :                 CHAN 'TYPE' 'SCALAIRE';
3057 :     FINSI;
3058 :   FINSI;
3059 : 
3060 : ****************************************************
3061 : ***** FACTEUR AFFECTE A J SI ON CALCULE C*(H) ******
3062 : ****************************************************
3063 :   FACT1 = 1.;
3064 :   SI (EGA IINTE 3);
3065 :      CHAR2 = EXTR CHAR1 'MECA';
3066 :      COE1 = ((MINI CHAF1) + (MAXI CHAF1)) / 2.;
3067 :      COE2 = ((MINI CHAF2) + (MAXI CHAF2)) / 2.;
3068 :      COE3 = ((MINI CHAF3) + (MAXI CHAF3)) / 2.;
3069 :      N1 = COE2 / COE3;
3070 :     SI (EGA INST 0. 1.E-10); INST = 1.E-10; FINSI;
3071 :      P1 = PROG 0. 'PAS' (INST / 100.) INST;
3072 :      P2 = PROG;
3073 :     REPETER BF1 (DIME P1);
3074 :        T1 = EXTR P1 &BF1;
3075 :        V1 = 0.;
3076 :       REPETER BF2 (DIME CHAR2);
3077 :          E1 = EXTR CHAR2 'EVOL' &BF2;
3078 :          V1 = V1 + ('IPOL' T1 (EXTR E1 'ABSC') (EXTR E1 ORDO));
3079 :       FIN BF2;
3080 :        P2 = P2 ET (PROG V1);
3081 :     FIN BF1;
3082 :      E1 = EVOL 'MANU' 'TEMPS' P1 'FORCE' (P2 ** N1);
3083 :      FACT1 = EXTR (SOMM E1) 1;
3084 :     SI (EGA FACT1 0. 1.E-10); FACT1 = 1.E+30; FINSI;
3085 :      FACT1 = ((V1 ** N1) / FACT1) ** COE3;
3086 :   FINSI;
3087 : 
3088 : ********************************************************
3089 : * SI LA COURBE DE TRACTION DEPEND DE LA TEMPERATURE ON *
3090 : * CALCULE LA VARIATION DE CONTRAINTES DE VON-MISES LORS*
3091 : * D'UNE AUGMENTATION (DETATE) DE LA TEMPERATURE A INST *
3092 : ********************************************************
3093 : *bp 11/08/2011 : on calcule plutot la variation de la limite d elasticite
3094 : *                ce qui est moins faux si presence de decharge...
3095 : *            => on construit VM1 comm VM2 !
3096 :   SI (((DIME MODPLA) '>' 0) ET ITHER);
3097 :     DETATE = 1.;
3098 :     TEP1 = TEPABS + (MANU 'CHPO' ELTETA 1 'T' DETATE);
3099 :     EPS1 = (EXCO VARF 'EPSE') CHAN 'TYPE' 'SCALAIRE';
3100 :     MSQ1 = 'MASQ' 'SUPERIEUR' EPS1 1.E-10;
3101 :     EPS1 = MSQ1 * EPS1;
3102 : *rem : on pourrait utiliser BORN dans le futur
3103 : *     VMI1 = (CHAN ('VMIS' OBJMOD SIGF MAT1)
3104 : *                     'TYPE' 'SCALAIRE')*MSQ1;
3105 : *     DETAVM = VMI1 * 0.;
3106 :     DETAVM = 0.;
3107 :     REPETER BCMOD2 NBOBJ;
3108 :       MODI = TABMOD.&BCMOD2;
3109 : *       VM1 = REDU VMI1 MODI;
3110 :       EPS2 = redu EPS1 MODI;
3111 :       SI (EXIS MODPLA &BCMOD2);
3112 :         SI (EGA MODPLA.&BCMOD2 1);
3113 :            MA1 = VARI 'NUAG' MODI (MATE MODI
3114 :                  'TRAC' TABTRA.&BCMOD2) TEPABS;
3115 :            VM1 = VARI 'NUAG' MODI (EXCO 'TRAC' MA1 'SIGM')
3116 :                  EPS2 'STRESSES' 'SCALAIRE';
3117 :            VM1 = (EXCO VM1 'SIGM' 'SCAL') CHAN 'TYPE' 'SCALAIRE';
3118 :            MA2 = 'VARI' 'NUAG' MODI (MATE MODI
3119 :                  'TRAC' TABTRA.&BCMOD2) TEP1;
3120 : *            VM2 = 'VARI' 'NUAG' MODI (EXCO 'TRAC' MA2 SIGM)
3121 : *                  (REDU EPS1 MODI) 'STRESSES' 'SCALAIRE';
3122 : *            VM2 = ((EXCO VM2 SIGM 'SCAL')CHAN 'TYPE' 'SCALAIRE')
3123 : *                         * (REDU MSQ1 MODI);
3124 :            VM2 = VARI 'NUAG' MODI (EXCO 'TRAC' MA2 'SIGM')
3125 :                  EPS2 'STRESSES' 'SCALAIRE';
3126 :            VM2 = (EXCO VM2 SIGM 'SCAL') CHAN 'TYPE' 'SCALAIRE';
3127 :         FINSI;
3128 :         SI ((EGA MODPLA.&BCMOD2 2) OU
3129 :             (EGA MODPLA.&BCMOD2 3));
3130 :           MA1 = VARI 'NUAG' MODI (REDU OBJMAT MODI) TEPABS;
3131 :           VM1 = CHAN 'STRESSES' MODI (EXCO MA1 'SIGY' 'SCAL');
3132 :           MA2 = VARI 'NUAG' MODI (REDU OBJMAT MODI) TEP1;
3133 : *           EPS2 = REDU EPS1 MODI;
3134 :           VM2 = CHAN 'STRESSES' MODI (EXCO MA2 'SIGY' 'SCAL');
3135 :           SI (EGA MODPLA.&BCMOD2 2);
3136 :              HSCAL1 = (EXCO MA1 'H' 'SCAL') CHAN 'TYPE' 'SCALAIRE';
3137 :              VM1 = VM1 + ((CHAN 'STRESSES' MODI HSCAL1)*EPS2);
3138 :              HSCAL2 = (EXCO MA2 'H' 'SCAL') CHAN 'TYPE' 'SCALAIRE';
3139 :              VM2 = VM2 + ((CHAN 'STRESSES' MODI HSCAL2)*EPS2);
3140 :           FINSI;
3141 : *           VM2 = VM2 * (REDU MSQ1 MODI);
3142 :         FINSI;
3143 :         si(ega (type DETAVM) 'FLOTTANT');
3144 :           DETAVM = ((VM2 - VM1) / DETATE);
3145 :         sinon;
3146 :           DETAVM = DETAVM + ((VM2 - VM1) / DETATE);
3147 :         finsi;
3148 :       FINSI;
3149 :     FIN BCMOD2;
3150 :   FINSI;
3151 : 
3152 : *******************************************************
3153 : **** ENERGIE DE DEFORMATION ELASTIQUE ET PLASTIQUE ****
3154 : *******************************************************
3155 : ***
3156 : *** DENSITE D'ENERGIE EN ELASTO OU THERMO-ELASTO-PLASTICITE ET
3157 : *** DENSITE D'ENERGIE LIEE A LA VARIATION DE COURBE DE TRACTION
3158 : ***
3159 :   SI (NEG IINTE 2);
3160 :     WELAS = 0.5*('ENER' OBJMOD SIGF ('ELAS' OBJMOD SIGF MAT1));
3161 : *     SI (IPAP ET (NEG IINTE 99));
3162 :     SI (IPAP);
3163 :       SI (EGA IABC 0);
3164 :         VMI1 = CHAN ('VMIS' OBJMOD SIGF MAT1) 'TYPE' 'SCALAIRE';
3165 :         SI (ICOQU ET ILIN); VMI1 = VMI1*OBJMOD EPAICH; FINSI;
3166 :         WPLAS=0.5*VMI1*((EXCO VARF 'EPSE')CHAN 'TYPE' 'SCALAIRE');
3167 :         SI (((DIME MODPLA) '>' 0) ET ITHER);
3168 :            WVMIS = 0.5*DETAVM*((EXCO VARF 'EPSE')
3169 :                     CHAN 'TYPE' 'SCALAIRE') ;
3170 :         FINSI;
3171 :       SINON ;
3172 :         VMI1 = CHAN (0.5*(('VMIS' OBJMOD SIG1 MAT1) +
3173 :                  ('VMIS' OBJMOD SIGF MAT1))) 'TYPE' 'SCALAIRE';
3174 :         SI (ICOQU ET ILIN); VMI1 = VMI1*OBJMOD EPAICH; FINSI;
3175 :         WPLAS = WPLAS + (VMI1*((EXCO (VARF - VAR1) 'EPSE')
3176 :                CHAN 'TYPE' 'SCALAIRE')) ;
3177 :         SI (((DIME MODPLA) '>' 0) ET ITHER);
3178 :            WVMIS = WVMIS + ((0.5*(DETAV1 + DETAVM))*
3179 :                    ((EXCO (VARF - VAR1) 'EPSE')
3180 :                  CHAN 'TYPE' 'SCALAIRE') );
3181 :         FINSI;
3182 :       FINSI ;
3183 :       ENERM = WELAS + WPLAS;
3184 :       SIG11 = SIG1 ; SIG1 = SIGF*1.; VAR11 = VAR1 ; VAR1 = VARF*1.;
3185 :       SI (((DIME MODPLA) '>' 0) ET ITHER);
3186 :          DETAV1 = DETAVM;
3187 :       FINSI;
3188 :     SINON;
3189 :       ENERM = WELAS;
3190 :     FINSI;
3191 :   FINSI;
3192 : ***
3193 : *** DENSITE D'ENERGIE POUR LES FLUAGES DONT ON A UNE
3194 : *** EXPRESSION EXPLICITE DE L'INTEGRATION SUR LE TEMPS
3195 : ***
3196 :   SI ((EGA ITYPEF 1) ET (IINTE EGA 2));
3197 :      UN1 = CHAF1 * (CHAF1 ** (-1.));
3198 :     SI (EXIS MAT1 'AF0 ');
3199 :        COE1 = ((MINI (CHAF2 + UN1)) + (MAXI (CHAF2 + UN1)))/2.;
3200 :        COE2 = ((MINI (CHAF4 + UN1)) + (MAXI (CHAF4 + UN1)))/2.;
3201 :        COE3 = ((MINI (CHAF6 + UN1)) + (MAXI (CHAF6 + UN1)))/2.;
3202 :        VMI1 = (EXCO ('VMIS' OBJMOD SIGF MAT1) 'SCAL')
3203 :                        CHAN 'TYPE' 'SCALAIRE'         ;
3204 :        ENERM1 = (CHAF2*((CHAF2 + UN1)**(-1.)))*CHAF1*(VMI1**COE1);
3205 :        ENERM2 = (CHAF4*((CHAF4 + UN1)**(-1.)))*CHAF3*(VMI1**COE2);
3206 :        ENERM3 = (CHAF6*((CHAF6 + UN1)**(-1.)))*CHAF5*(VMI1**COE3);
3207 :        ENERM = ENERM1 + ENERM2 + ENERM3;
3208 :     SINON;
3209 :        COE1 = ((MINI (CHAF2 + UN1)) + (MAXI (CHAF2 + UN1)))/2.;
3210 :        COE2 = ((MINI CHAF3) + (MAXI CHAF3))/2.;
3211 :        VMI1 = (EXCO ('VMIS' OBJMOD SIGF MAT1) 'SCAL')
3212 :               CHAN 'TYPE' 'SCALAIRE'     ;
3213 :       SI ((EGA INST 0. 1.E-10) ET ('<'(COE2 - 1) 0.));
3214 :          V1 = 0.;
3215 :       SINON;
3216 :          V1 = INST**(COE2 - 1);
3217 :       FINSI;
3218 :        ENERM = (CHAF2*((CHAF2 + UN1)**(-1.)))*
3219 :                 CHAF1*(VMI1**COE1)*CHAF3*V1;
3220 :     FINSI;
3221 :     SI (ICOQU ET ILIN); ENERM = ENERM*OBJMOD EPAICH; FINSI;
3222 :   FINSI;
3223 : ***
3224 : *** ON N'A PAS UNE EXPRESSION EXPLICITE DE
3225 : *** L'INTEGRATION DU FLUAGE SUR LE TEMPS
3226 : ***
3227 :   SI ((EGA ITYPEF 2) ET (IINTE EGA 2));
3228 :      SIGMOY = 0.5*(SIG1 + SIGF);
3229 :     SI ((EGA IABC 0) ET (NON IREPRI));
3230 :        ENERM = 'ENER' OBJMOD VITDFI SIGMOY;
3231 :     SINON ;
3232 :        ENERM = ENERM + ('ENER' OBJMOD (VITDFI - VDI1) SIGMOY);
3233 :     FINSI ;
3234 :      SIG11 = SIG1 ; SIG1 = SIGF; VDI1 = VITDFI;
3235 :   FINSI;
3236 : 
3237 : ***************************************************
3238 : *********** APPEL A LA PROCEDURE G_CALCUL *********
3239 : ***************************************************
3240 :      INFTAB.'IABC' = IABC;
3241 :      INFTAB.'&BOUCEXT' = &BOUCEXT;
3242 :      INFTAB.'&BOUCMIX' = &BOUCMIX;
3243 :      INFTAB.'INST' = INST;
3244 :      INFTAB.'FACT1' = FACT1;
3245 :      INFTAB.'C_MATE' = C_MATE;
3246 :      INFTAB.'MAT1' = MAT1;
3247 :      INFTAB.'RIGTOT' = RIGTOT;
3248 :      INFTAB.'ENERM' = ENERM;
3249 :      INFTAB.'WVMIS' = WVMIS;
3250 :      INFTAB.'TALPH1' = TALPH1;
3251 :      INFTAB.'TEPINT' = TEPINT;
3252 :      INFTAB.'TEPABS' = TEPABS;
3253 :      INFTAB.'PREINT' = PREINT;
3254 :      INFTAB.'DEPINT' = DEPINT;
3255 :      INFTAB.'SIGF' = SIGF;
3256 :      INFTAB.'SIG1' = SIG11;
3257 :      INFTAB.'VARF' = VARF;
3258 :      INFTAB.'VITF' = VITF;
3259 :      INFTAB.'ACCF' = ACCF;
3260 :      INFTAB.'MOTMIX' = MOTMIX;
3261 :      INFTAB.'MOTMIA' = MOTMIA;
3262 :      INFTAB.'A_PREI' = A_PREI;
3263 :      INFTAB.'A_DEPI' = A_DEPI;
3264 :      INFTAB.'A_SIGF' = A_SIGF;
3265 :      INFTAB.'A_DEPGR'= A_DEPGR;
3266 : *    ajout sm
3267 :      INFTAB.'DEFINT' = DEFINT;
3268 : *    ajout BP BT
3269 :      si(IFROT);
3270 :        INFTAB . 'WDEP'  = WDEP  ;
3271 :        INFTAB . 'SIGCON' = SIGCON ;
3272 :        INFTAB . 'B_DEPGR' =  B_DEPGR;
3273 :        INFTAB . 'OBJCON2' =  OBJCON2;
3274 :      fins;
3275 :      
3276 : 
3277 : *si(flmess);   mess 'appel a  G_CALCUL';  fins;
3278 : *      mess 'avant g_calcul';
3279 : *      si(exis SUPTAB 'RESULTATS'); list SUPTAB.'RESULTATS' ; fins;
3280 :      G_CALCUL SUPTAB INFTAB;
3281 : *      mess 'apres g_calcul';
3282 : *si(flmess); list SUPTAB.'RESULTATS' ;   fins;
3283 : 
3284 :   FIN BOUCMIX;
3285 : * FIN DE BOUCLE SUR LES INTEGRALES A CALCULER ========================*
3286 : 
3287 : 
3288 : ****************************************************
3289 : ****** ON RECUPERE LA CONFIGURATION INITIALE *******
3290 : ****************************************************
3291 :   SI IPAP;
3292 :     SI IGDER;
3293 :        FORM CONFIG0;
3294 :     FINSI;
3295 :   FINSI;
3296 : 
3297 : ****************************************************
3298 : ******* FIN DE BOUCLE SUR LES PAS DE CALCUL ********
3299 : ****************************************************
3300 : FIN BOUCEXT ;
3301 : 
3302 : 
3303 : ****************************************************
3304 : ** STOCKAGE DES RESULTATS DANS L OBJET EVOLUTIONS **
3305 : ****************************************************
3306 : * en plus de ce qui existe dans SUPTAB.'RESULTATS'
3307 : * et dans CHPO_RESULTATS (fait dans G_calcul)
3308 : SI (IPAP et (non IPERSO1));
3309 :   IND1 = INDE (SUPTAB.'RESULTATS');
3310 : 
3311 : * Cas 2D *******************************************
3312 :   SI (EGA &DIME 2);
3313 :     TITR CHA1;
3314 :     SI (EGA IINTE 99);
3315 :       REPETER BB1 (DIME IND1);
3316 :          MOT1 = IND1.&BB1; PT = PROG; PG = PROG;
3317 :          IND2 = INDE (SUPTAB.'RESULTATS'.MOT1);
3318 :         REPETER BB2 (DIME IND2);
3319 :            P1 = &BB2 - 1;
3320 :            PT = PT ET (PROG SUPTAB.'SOLUTION_PASAPAS'.'TEMPS'.P1);
3321 :            PG = PG ET (PROG SUPTAB.'RESULTATS'.MOT1.P1);
3322 :         FIN BB2;
3323 :          E1 = EVOL 'MANU' 'TEMPS' PT MOTTI PG;
3324 :          SUPTAB.'EVOLUTION_RESULTATS'.MOT1 = E1;
3325 :       FIN BB1;
3326 :     SINON;
3327 :        PT = PROG; PG = PROG;
3328 :       REPETER BB1 (DIME IND1);
3329 :          P1 = &BB1 - 1;
3330 :          PT = PT ET (PROG SUPTAB.'SOLUTION_PASAPAS'.'TEMPS'.P1);
3331 :          PG = PG ET (PROG SUPTAB.'RESULTATS'.P1);
3332 :       FIN BB1;
3333 :        E1 = EVOL 'MANU' 'TEMPS' PT MOTTI PG;
3334 :        SUPTAB.'EVOLUTION_RESULTATS' = E1;
3335 :     FINSI;
3336 :   FINSI;
3337 : 
3338 : * Cas 3D *******************************************
3339 : * Rem BP :
3340 : * !!! Attention passer ici lorsque la fissure propage doit conduire
3341 : *     a des erreurs, puisque les points du front changent ...!!!
3342 : *     Il faudrait a terme supprimer cette mise en forme des resultats
3343 : *     ainsi que celle utilisant la table RESULTATS (cf g_calcul).
3344 : *     Seul CHPO_RESULTATS semblent perenne pour la porpagation.
3345 :   SI ((EGA &DIME 3) ET (NON ICOQU));
3346 :     SI (EGA IINTE 99);
3347 :       REPETER BB1 (DIME IND1);
3348 :          MOT1 = IND1.&BB1;
3349 :          IND2 = INDE (SUPTAB.'RESULTATS'.MOT1);
3350 :          IND3 = INDE (SUPTAB.'RESULTATS'.MOT1.(IND2.1));
3351 :         REPETER BB2 (DIME IND3);
3352 :            PM = IND3.&BB2; PT = PROG; PG = PROG;
3353 :           SI (EGA &BB2 (DIME IND3));
3354 :              CHA2 = CHAI ' (Global)';
3355 :           SINON;
3356 : *              CHA2 = CHAI ' (Pt ' ('NOEUD' PM) ')';
3357 :              CHA2 = CHAI ' (Pt ' (&BB2) ')';
3358 :           FINSI;
3359 :           'TITR' (CHAI CHA1 CHA2);
3360 :           REPETER BB3 (DIME IND2);
3361 :              P1 = &BB3 - 1;
3362 :              PT = PT ET (PROG SUPTAB.'SOLUTION_PASAPAS'.'TEMPS'.P1);
3363 :              PG = PG ET (PROG SUPTAB.'RESULTATS'.MOT1.P1.PM);
3364 :           FIN BB3;
3365 :            E1 = EVOL 'MANU' 'TEMPS' PT MOTTI PG;
3366 :            SUPTAB.'EVOLUTION_RESULTATS'.MOT1.PM = E1;
3367 :         FIN BB2;
3368 :       FIN BB1;
3369 :     SINON;
3370 :        IND2 = INDE (SUPTAB.'RESULTATS'.(IND1.1));
3371 :       REPETER BB1 (DIME IND2);
3372 :          PM = IND2.&BB1; PT = PROG; PG = PROG;
3373 :         SI (EGA &BB1 (DIME IND2));
3374 :            CHA2 = CHAI ' (Global)';
3375 :         SINON;
3376 : *            CHA2 = CHAI ' (Pt ' ('NOEUD' PM) ')';
3377 :            CHA2 = CHAI ' (Pt ' (&BB1) ')';
3378 :         FINSI;
3379 :         'TITR' (CHAI CHA1 CHA2);
3380 :         REPETER BB2 (DIME IND1);
3381 :            P1 = &BB2 - 1;
3382 :            PT = PT ET (PROG SUPTAB.'SOLUTION_PASAPAS'.'TEMPS'.P1);
3383 :            PG = PG ET (PROG SUPTAB.'RESULTATS'.P1.PM);
3384 :         FIN BB2;
3385 :          E1 = EVOL 'MANU' 'TEMPS' PT MOTTI PG;
3386 :          SUPTAB.'EVOLUTION_RESULTATS'.PM = E1;
3387 :       FIN BB1;
3388 :     FINSI;
3389 :   FINSI;
3390 : 
3391 : * Cas coque *******************************************
3392 :   SI ICOQU;
3393 :      IND2 = INDE (SUPTAB.'RESULTATS'.(IND1.1));
3394 :     REPETER BB1 (DIME IND2);
3395 :        PM = MOT IND2.&BB1; PT = PROG; PG = PROG;
3396 : *        CHA2 = CHAI ' (' PM ')';
3397 :        CHA2 = CHAI ' (Pt ' (&BB1) ')';
3398 :       'TITR' (CHAI CHA1 CHA2);
3399 :       REPETER BB2 (DIME IND1);
3400 :          P1 = &BB2 - 1;
3401 :          PT = PT ET (PROG SUPTAB.'SOLUTION_PASAPAS'.'TEMPS'.P1);
3402 :          PG = PG ET (PROG SUPTAB.'RESULTATS'.P1.PM);
3403 :       FIN BB2;
3404 :        E1 = EVOL 'MANU' 'TEMPS' PT MOTTI PG;
3405 :        SUPTAB.'EVOLUTION_RESULTATS'.PM = E1;
3406 :     FIN BB1;
3407 :   FINSI;
3408 : 
3409 : FINSI;
3410 : 
3411 : 
3412 : ********************************************************
3413 : **** STOCKAGE POUR UNE EVENTUELLE REPRISE DE CALCUL ****
3414 : ********************************************************
3415 : SI IPAP;
3416 :    SUPTAB.'OBJ1' = MOT SUPTAB.'OBJECTIF' ;
3417 :    SUPTAB.'IABC' = IABC ;
3418 :    SUPTAB.'MAT1' = MAT1 ;
3419 :    SUPTAB.'END1' = ENERM;
3420 :   SI (((DIME MODPLA) '>' 0) ET ITHER);
3421 :      SUPTAB.'ENV1' = WVMIS;
3422 :   FINSI;
3423 :   SI (EGA ITYPEF 2);
3424 :      SUPTAB.'VDI1' = VDI1 ;
3425 :   FINSI;
3426 : FINSI;
3427 : 
3428 : ********************************************************
3429 : ************ ENLEVEMENT DES OBJETS INUTILES ************
3430 : ********************************************************
3431 : si(exis SUPTAB 'COUCHE');
3432 :       SUPTAB.'COU1' = SUPTAB.'COUCHE' ;
3433 : fins;
3434 : *    SUPTAB.'CHAMP_THETA' = ELTETA ;
3435 : *bp: soyons coherent... CHAMP_THETA est deja renseigné
3436 : SUPTAB.'CHAMP_THET1' = SUPTAB.'CHAMP_THETA' ;
3437 : OTER SUPTAB 'CHAMP_THETA';
3438 : * UTILTETA sera pratique pour la propagation xfem
3439 : SUPTAB.'UTILTET1' = UTILTETA;
3440 : * ELTETA sera pratique pour les tests de reprise
3441 : SUPTAB.'ELTET1' = ELTETA ;
3442 : * SUPTAB = 'ENLE' SUPTAB 'MAILLAGE';
3443 : * si (exis SUPTAB 'FISSURE');
3444 : *     SUPTAB = 'ENLE' SUPTAB 'FISSURE';
3445 : * fins;
3446 : * SUPTAB = 'ENLE' SUPTAB 'CHAMP_THETA';
3447 : * SI (EXIS SUPTAB 'CHAMP_PI');
3448 : *    SUPTAB = 'ENLE' SUPTAB 'CHAMP_PI';
3449 : * FINSI;
3450 : * SI (EXIS SUPTAB 'EPAISSEUR');
3451 : *    SUPTAB = 'ENLE' SUPTAB 'EPAISSEUR';
3452 : * FINSI;
3453 : * FINPROC SUPTAB;
3454 : *bp: ENLE crée une 2nde table jamais récupérée:
3455 : *    mieux vaut 1 seule SUPTAB bien définie!
3456 : FINPROC;
3457 : 
3458 : 
3459 : 
3460 : 
3461 : 
3462 : 
3463 : 
3464 :  
3465 :  
3466 :  
3467 :  

© Cast3M 2003 - All rights reserved.
Disclaimer