Télécharger la plaquette Cast3M


fr en

Télécharger unpas.procedur

Retour à la liste

Numérotation des lignes : 

  1. * UNPAS     PROCEDUR  MB234859  24/07/08    21:15:03     11958          
  2. 'DEBPROC' UNPAS PRECED*'TABLE';
  3. *----------------------------------------------------------------------*
  4. *                          PROCEDURE UNPAS                             *
  5. *                                                                      *
  6. * Calcul d'un increment de solution en grand deplacement plastique par *
  7. * la methode des residus.                                              *
  8. *                                                                      *
  9. * Les differentes configurations qui interviennent sont :              *
  10. * WTAB.'FOR0'    : configuration de debut de calcul                    *
  11. * WTAB.'GE0_DEB' :                                                     *
  12. *      -> en GRANDS_DEPLACEMENTS ou FEFP_FORMULATION (1) :             *
  13. *         configuration de debut de pas                                *
  14. *      -> dans les autres cas (2) :                                    *
  15. *         configuration de debut de calcul (= WTAB.'FOR0)              *
  16. *  Dans le cas (1), la configuration WTAB.'GE0_DEB' est actualisee     * 
  17. *  dans PASAPAS avec GEOM2.                                            *
  18. *  Dans le cas (2), il n'y a pas d'actualisation de la configuration   *
  19. *  au cours du calcul.                                                 * 
  20. *                                                                      *
  21. *  Pour les calculs DYNAMIQUE, un schema de Newmark (implicite)        *
  22. *  avec gamma = 1/2 et beta = 1/4 est utilise.                         *
  23. *                                                                      *
  24. *----------------------------------------------------------------------*
  25. *CB215821 : Recuperation de XPETIT (07/12/2016)
  26. XPETIT   = 'VALE' 'PETI'    ;
  27. XGRAND   = 'VALE' 'GRAND'   ;
  28. *
  29. WTAB  = PRECED.'WTABLE'      ;
  30. conti = PRECED.'CONTINUATION';
  31. estim = PRECED.'ESTIMATION'  ;
  32. *
  33. * Liste de composantes utiles
  34. MXMYMZ = 'MOTS' 'MX'  'MY'  'MZ'  'MT'  'FP'  'FPQ'  'FTP'
  35.                 'IMX' 'IMY' 'IMZ' 'IMT' 'IFP' 'IFPQ' 'IFTP' ;
  36. MXMFLX = 'MOTS' 'MX'  'MY'  'MZ'  'MT'  'FLX' 'FP'  'FPQ'  'FTP'
  37.                 'IMX' 'IMY' 'IMZ' 'IMT' 'FLX' 'IFP' 'IFPQ' 'IFTP' ;
  38. MNPRIM = 'MOTS' 'UX'  'UY'  'UZ'  'UR'  'UT'   'RX'   'RY'  'RZ'  'RT'
  39.                       'P'   'PQ'  'TP'   'ALFA' 'BETA'
  40.                 'IUX' 'IUY' 'IUZ' 'IUR'  'IUT'  'IRX'  'IRY' 'IRZ' 'IRT'
  41.                 'IP' 'IPQ'  'ITP' 'IALF' 'IBET' ;
  42. MLPRIM = MNPRIM 'ET' ('MOTS' 'LX');
  43. MNDUAL = 'MOTS' 'FX'  'FY'  'FZ'  'FR'   'FT'   'MX'   'MY'  'MZ'  'MT'
  44.                       'FP'  'FPQ' 'FTP'  'FALF' 'FBET'
  45.                 'IFX' 'IFY' 'IFZ' 'IFR'  'IFT'  'IMX'  'IMY' 'IMZ' 'IMT'
  46.                 'IFP' 'IFPQ' 'IFTP' 'IFAL' 'IFBE' ;
  47. MLDUAL = MNDUAL 'ET' ('MOTS' 'FLX');
  48. MLDEPL = 'MOTS' 'UX' 'UY' 'UZ' 'UR' 'UT'  'IUX' 'IUY' 'IUZ' 'IUR' 'IUT'
  49.                 'ALFA' 'BETA' 'IALF' 'IBET' ;
  50. MLROTA = 'MOTS' 'RX' 'RY' 'RZ' 'RT' 'P' 'PQ' 'TP'
  51.                 'IRX' 'IRY' 'IRZ' 'IRT' 'IP' 'IPQ' 'ITP';
  52. MLDEFOR = 'MOTS'  'EPXX' 'EPYY' 'EPZZ' 'EPSS' 'EPTT' 'EPRR'
  53.                   'GAXY' 'GAXZ' 'GAYZ' 'GAST' 'GASN' 'GATN'
  54.                   'GARZ' 'GART' 'GAZT' 'GXY '
  55. *                 'CX  ' 'CY  ' 'CZ  '
  56.                   'EPSE' 'EPS' ;
  57. * Pour MOD_LIAISON :
  58. MLPRIM_LIA = MLPRIM ;
  59. MVPRIM_LIA = 'MOTS' 'VTX' 'VTY' 'VTZ' 'VTR' 'VTT' 'VWX' 'VWY' 'VWZ' 'VWT'
  60.                     'VVP' 'VVPQ' 'VVTP' 'VALF' 'VBET'
  61.                     'IVTX' 'IVTY' 'IVTZ' 'IVTR' 'IVTT' 'IVWX' 'IVWY' 'IVWZ' 'IVWT'
  62.                     'IVVP' 'IVPQ' 'IVTP' 'IVAL' 'IVBE' 'VLX';
  63. MLDUAL_LIA = MLDUAL ;
  64.  
  65. *=DEB===== Formulation HHO = Definition des DDLs ===============================
  66. * Definition des DDLs primaux : inconnues de deplacements des cellules et des Faces
  67. MNPRIM_HHO = 'MOTS'
  68. *                'UXC0' 'UXC1' 'UXC2' 'UXC3' 'UXC4' 'UXC5' 'UXC6' 'UXC7' 'UXC8' 'UXC9'
  69. *                'UYC0' 'UYC1' 'UYC2' 'UYC3' 'UYC4' 'UYC5' 'UYC6' 'UYC7' 'UYC8' 'UYC9'
  70. *                'UZC0' 'UZC1' 'UZC2' 'UZC3' 'UZC4' 'UZC5' 'UZC6' 'UZC7' 'UZC8' 'UZC9'
  71.                 'UXF0' 'UXF1' 'UXF2' 'UXF3' 'UXF4' 'UXF5' 'UXF6' 'UXF7' 'UXF8' 'UXF9'
  72.                 'UYF0' 'UYF1' 'UYF2' 'UYF3' 'UYF4' 'UYF5' 'UYF6' 'UYF7' 'UYF8' 'UYF9'
  73.                 'UZF0' 'UZF1' 'UZF2' 'UZF3' 'UZF4' 'UZF5' 'UZF6' 'UZF7' 'UZF8' 'UZF9' ;
  74. MNDUAL_HHO = 'MOTS'
  75. *                'FXC0' 'FXC1' 'FXC2' 'FXC3' 'FXC4' 'FXC5' 'FXC6' 'FXC7' 'FXC8' 'FXC9'
  76. *                'FYC0' 'FYC1' 'FYC2' 'FYC3' 'FYC4' 'FYC5' 'FYC6' 'FYC7' 'FYC8' 'FYC9'
  77. *                'FZC0' 'FZC1' 'FZC2' 'FZC3' 'FZC4' 'FZC5' 'FZC6' 'FZC7' 'FZC8' 'FZC9'
  78.                 'FXF0' 'FXF1' 'FXF2' 'FXF3' 'FXF4' 'FXF5' 'FXF6' 'FXF7' 'FXF8' 'FXF9'
  79.                 'FYF0' 'FYF1' 'FYF2' 'FYF3' 'FYF4' 'FYF5' 'FYF6' 'FYF7' 'FYF8' 'FYF9'
  80.                 'FZF0' 'FZF1' 'FZF2' 'FZF3' 'FZF4' 'FZF5' 'FZF6' 'FZF7' 'FZF8' 'FZF9' ;
  81. * Pour l'istant on ne considere que les inconnues (DEPL et FORC) de face pour
  82. * le calcul de l'equilibre et l'acceleration de convergence
  83. * D'ou les lignes commentees ci-dessus dans la definition de MNPRIM_HHO et MNDUAL_HHO
  84. MNPRIM = MNPRIM 'ET' MNPRIM_HHO ;
  85. MLPRIM = MLPRIM 'ET' MNPRIM_HHO ;
  86. MNDUAL = MNDUAL 'ET' MNDUAL_HHO ;
  87. MLDUAL = MLDUAL 'ET' MNDUAL_HHO ;
  88. MLDEPL = MLDEPL 'ET' MNPRIM_HHO ;
  89. *A Ameliorer : ne faire que si presence HHO et n'ajouter que les DDLs utilises par HHO !
  90. *=FIN===== Formulation HHO =====================================================
  91. *
  92. * Pour stocker des informations necessaires aux calculs d'usure
  93. WTAB.'POST_COFR' = 'TABL' ;
  94. *----------------------------------------------------------------------*
  95. *             Options de pasapas - parametres du calcul                *
  96. *----------------------------------------------------------------------*
  97. * 1- rigidite
  98. IKSIA  = WTAB.'K_SIGMA';
  99. IKTAN  = WTAB.'K_TANGENT' ;
  100. ZKTASYM = 'TEXTE' '    ' ;
  101. 'SI' WTAB.'K_TANGENT_SYME' ;
  102.   ZKTASYM = 'MOT' 'SYME' ;
  103. 'FINS';
  104. IMPLP   = WTAB.'LIAISON_PERSISTANTE';
  105. AUTAUG  = WTAB.'AUTOAUGM';
  106. IRAUGLU = WTAB.'RAIDAUGM';
  107. IRAUG   = FAUX;
  108. 'SI' IRAUGLU ;
  109.   RIG_AUG = WTAB.'RIGIDITE_AUGMENTEE';
  110.   'SI' ('NON' AUTAUG) ; 
  111.     IRAUG = VRAI ;
  112.   'FINS';
  113. 'FINS';
  114. IRCON = WTAB.'RAIDCONST';
  115. 'SI' IRCON;
  116.   RIG_CONS = WTAB.'RIGIDITE_CONSTANTE';
  117.   MAI_CONS ='EXTR' RIG_CONS 'MAIL'    ;
  118. 'FINS';
  119. *  declenchement du recalcul de la matrice
  120. DLTAIT = WTAB.'DELTAITER' ;
  121. ITRCLC =  1  '*' DLTAIT ;
  122. 'SI' (ITRCLC <  20) ; ITRCLC =  20; 'FINSI';
  123. *
  124. * 2- Type de formulation
  125. IFEFP  = WTAB.'FEFP_FORMULATION' ;
  126. LAG_TOT= WTAB.'LAG_TOT';
  127. ISSTE  = WTAB.'SUBSTEPPING';
  128. IFEFPUL= WTAB.'UPDATE_LAGRANGIAN';
  129. LNLOC  = WTAB.'NLOC';
  130. 'SI' (LNLOC 'ET' ('EGA' WTAB.'NON_LOCAL' 'HELM'));
  131.   TAHELM = WTAB.'HELMHOLTZ' ;
  132.   NHELM =  TAHELM . 'N_VARI_NL' ;
  133. 'FINSI' ;
  134. *
  135. * 3- Type de materiau
  136. ICERAM = WTAB.'CERAMIQUE' ;
  137. IENDOM = WTAB.'ENDOMMAGEMENT';
  138. IPLAVI = WTAB.'IPLAVI';
  139. POR1   = WTAB.'POR1' ;
  140. IVIEXT = WTAB.'VISCO_EXTERNE';
  141. IVIDOM = WTAB.'VISCODOMMAGE';
  142. IVISCO = WTAB.'VISCOPLASTIQUE';
  143. *
  144. * 4- Critere de convergence/non-convergence
  145. ZMAXIT = WTAB.'MAXITERATION' ;
  146. NITMA  = WTAB.'NITERINTER_MAX';
  147. ZPRECM = WTAB.'PRECFLEX' ;
  148. ZPREK  = WTAB.'PRECISINTER' ;
  149. ZPREC  = WTAB.'PRECISION' ;
  150. ZPRECHPP = ZPREC ;
  151. IFTOL  = 'NEG' WTAB.'FTOL' 'INCONNU' ;
  152. 'SI' IFTOL ;
  153.   ZFTOL = 'ABS' WTAB.'FTOL' ;
  154. 'FINS';
  155. IMTOL  = 'NEG' WTAB.'MTOL' 'INCONNU';
  156. 'SI' IMTOL ;
  157.   ZMTOL = 'ABS' WTAB.'MTOL' ;
  158. 'FINS';
  159. 'SI' ('OU' IVISCO IVIDOM IVIEXT); ZPREK = 5.E-7 ; 'FINS';
  160. 'SI' IENDOM;                      ZPREK = ZPREC ; 'FINS';
  161. NSOINCR = WTAB.'SOUS_INCREMENT' ;
  162. *
  163. * 5- Type de calcul
  164. IPILOT  = WTAB.'AUTOMATIQUE';
  165. ISOL    = WTAB.'CONSOLIDATION';
  166. IDYN    = WTAB.'DYNAMIQUE';
  167. IGRD    = WTAB.'GRANDS_DEPLACEMENTS';
  168. HYPDEF  = 'MOT' WTAB.'HYPOTHESE_DEFORMATIONS';
  169. IJAUMA  = 'EGA' HYPDEF 'JAUMANN' ;
  170. HPP_EPS = 'EGA' WTAB.'PREDICTEUR' 'HPP';
  171. *
  172. * 6- Chargement particulier
  173. LOGDEF = WTAB.'CHAR_DEFI';
  174. LOGPRE = WTAB.'CHAR_PRES' ;
  175. LOGPIL = WTAB.'CHAR_PILO' ;
  176. ITHER  = WTAB.'ITHER' ;
  177. *
  178. * 7- Acceleration de convergence
  179. ZNACCE = VRAI;
  180. *
  181. * 8- Instant de debut et de fin de pas
  182. TEMPS0 = WTAB.'TEMPS0';
  183. TI     = WTAB.'T_FINAL';
  184. *
  185. * 9- Blocages mecaniques (BLOM)
  186. 'SI' WTAB.'CHAR_BLOM' ;
  187.   BLOM1 = 'TIRE' PRECED.'BLOCAGES_MECANIQUES' 'BLOM' TI ;
  188.   WTAB.'BLOCAGES_MECANIQUES' = BLOM1 ;
  189. 'FINSI' ;
  190. ZCLIM0 = WTAB.'BLOCAGES_MECANIQUES';
  191. ZCLIM  = ZCLIM0 ;
  192. *----------------------------------------------------------------------*
  193. *           Verification de la compatibilite entre options             *
  194. *----------------------------------------------------------------------*
  195. 'SI' (IKTAN 'ET' ('NON' IPLAVI)) ;
  196.   'MESS' 'IPLAVI faux : pas de matrice tangente ->'
  197.          ' on utilise la rigidite elastique' ;
  198.   IKTAN = FAUX ;
  199. 'FINS' ;
  200. IPERT = WTAB.'K_TANGENT_PERT' 'ET' ('NON' LNLOC) 'ET' IPLAVI ;
  201. * Matrice tangente : pas d'acceleration en cas de modele FEFP ou SSTE
  202. 'SI' IKTAN ;
  203.   'SI' (IFEFP 'OU' ISSTE) ; ZNACCE = FAUX; 'FINS';
  204. 'FINS' ;
  205. *
  206. * Configurations de reference et de debut de pas
  207. GEOREF0 = WTAB.'FOR0' ;
  208. GEOM1   = WTAB.'GE0_DEB' ;
  209. *----------------------------------------------------------------------*
  210. *           Modele complet                                             *
  211. *                                                                      *
  212. * Notations utilisees :                                                *
  213. *   ZMODLI : modeles mecanique + poreux (non //)                       *
  214. *   -> aucune parallelisation       alors ZMODL = ZMODLI               *
  215. *   -> parallelisation comportement alors ZMODL = ZMODLI               *
  216. *                                         MODRELOC est //              *
  217. *   -> parallelisation automatique  alors ZMODL est //                 *
  218. *----------------------------------------------------------------------*
  219. * PAS_MODL : mise a jour des indices de PRECED.WTABLE relatifs aux 
  220. *            modeles si PRECED.WTABLE.MODELE a ete modifie.
  221. PAS_MODL PRECED ;
  222. 'SI' WTAB.'FOR_MECA' ;
  223.   PARALLEL = FAUX ; PARTLOCA = FAUX ;
  224.   ZMODLI = WTAB.'MO_TOT' ;
  225.   ZMODL  = ZMODLI ;
  226.   'SI' ('EGA' WTAB.'PROCESSEURS' 'COMPORTEMENT') ;
  227.     PARALLEL = VRAI ; PARTLOCA = VRAI ;
  228.     MODRELOC = 'PART' 'ARLE' ZMODLI WTAB.'NBPART' ;
  229.   'FINS' ;
  230.   'SI' ('EGA' WTAB.'PROCESSEURS' 'AUTOMATIQUE') ;
  231.     PARALLEL = VRAI ; PARTLOCA = FAUX ;
  232.     ZMODL = 'PART' 'ARLE' ZMODLI WTAB.'NBPART' ;
  233.     'OPTI' 'PARA' VRAI ;
  234.   'FINS' ;
  235. *
  236. * Si le modele a change, adapter egalement les MCHAMLS
  237.   'SI' ('NEG' ZMODLI WTAB.'MO_TOT_PREC') ;
  238. *
  239.     CONTI.'CONTRAINTES' =('REDU' CONTI. 'CONTRAINTES' ZMODLI) '+' ('ZERO' ZMODLI 'CONTRAINTES') ;
  240.     CONTI.'DEFORMATIONS'=('REDU' CONTI.'DEFORMATIONS' ZMODLI) '+' ('ZERO' ZMODLI 'DEFORMAT'   ) ;
  241.     'SI' IPLAVI ;
  242.       CONTI.'DEFORMATIONS_INELASTIQUES'=('REDU' CONTI.'DEFORMATIONS_INELASTIQUES' ZMODLI) '+' ('ZERO' ZMODLI 'DEFINELA') ;
  243.       CONTI.'VARIABLES_INTERNES'       =('REDU' CONTI.'VARIABLES_INTERNES'        ZMODLI) '+' ('ZERO' ZMODLI 'VARINTER') ;
  244.     'FINS';
  245. *   
  246.     WTAB.'ETAT1' = PAS_ETAT PRECED TEMPS0      ;
  247.     WTAB.'MAT1'  = PAS_MATE PRECED WTAB.'ETAT1';
  248. *  
  249.     'SI' ITHER ;
  250.       WTAB.'ETHER1'    =('REDU' WTAB.'ETHER1'   ZMODLI) '+' ('ZERO' ZMODLI 'DEFORMAT') ;
  251.       'SI' POR1;
  252.         WTAB.'MSRTHER1'=('REDU' WTAB.'MSRTHER1' ZMODLI) '+' ('ZERO' ZMODLI 'DEFORMAT') ;
  253.       'FINS';
  254.     'FINS';
  255. *  
  256.     'SI' LOGDEF ; 
  257.       WTAB.'DEFOR1' = ('REDU' WTAB.'DEFOR1' ZMODLI) '+' ('ZERO' ZMODLI 'DEFORMAT') ;
  258.       WTAB.'DEFOR2' = ('REDU' WTAB.'DEFOR2' ZMODLI) '+' ('ZERO' ZMODLI 'DEFORMAT') ;
  259.     'FINS';
  260. *
  261.   'FINS';
  262. 'FINSI';
  263. *----------------------------------------------------------------------*
  264. *           Champ materiau et caracteristiques a l'instant TEMPS0      *
  265. *                                                                      *
  266. * Notations utilisees :                                                *
  267. *   ZMAT11 : materiau a l'instant TEMPS0 (non //)                      *
  268. *   -> aucune parallelisation       alors ZMAT1 = ZMAT11               *
  269. *   -> parallelisation comportement alors ZMAT1 = ZMAT11               *
  270. *   -> parallelisation automatique  alors ZMAT1 est //                 *
  271. *----------------------------------------------------------------------*
  272. ZMAT11 = WTAB.'MAT1' ;
  273. ZMAT1  = ZMAT11;
  274. 'SI' WTAB.'FOR_MECA' ;
  275.   'SI' ('EGA' WTAB.'PROCESSEURS' 'AUTOMATIQUE') ;
  276.     ZMAT1 = 'REDU' ZMAT11 ZMODL ;
  277.   'FINS';
  278. *
  279. * Distinguer les modeles poreux des modeles mecaniques
  280.   'SI' POR1;
  281.     MO_PORI = 'EXTR' ZMODLI 'FORM' 'POREUX';
  282.     MO_POR  = 'EXTR' ZMODL  'FORM' 'POREUX';
  283. *
  284. **  kich ma_por0 intervient si ISOL
  285. **  initialement     MA_POR0= 'REDU' MO_POR MAT1  ;
  286.     MA_POR0 = 'REDU' ZMAT11 MO_POR ;
  287.   'FINS';
  288. 'FINSI';
  289. *----------------------------------------------------------------------*
  290. *           Champ materiau et caracteristiques a l'instant TI          *
  291. *                                                                      *
  292. * Notations utilisees :                                                *
  293. *   ZMAT22 : materiau a l'instant TI (non //)                          *
  294. *   -> aucune parallelisation       alors ZMAT2 = ZMAT22               *
  295. *   -> parallelisation comportement alors ZMAT2 = ZMAT22               *
  296. *   -> parallelisation automatique  alors ZMAT2 est //                 *
  297. *----------------------------------------------------------------------*
  298. ZETAT2 = PAS_ETAT PRECED TI    ;
  299. ZMAT2I = PAS_MATE PRECED ZETAT2;
  300. *
  301. 'SI' WTAB.'FOR_MECA' ;
  302.   ZMAT22 = 'REDU' ZMAT2I ZMODLI ;
  303.   'SI' WTAB.'ITCAR';
  304.     'FORM' GEOREF0 ;
  305.     GEOBID ZMAT22 = 'FORM' CONTI.'DEPLACEMENTS' ZMODLI ZMAT22;
  306.     'FORM' GEOM1 ;
  307.   'FINS';
  308.   ZMAT2 = 'REDU' ZMAT22 ZMODL ;
  309.   'SI' POR1;
  310.     MA_POR = 'REDU' ZMAT2I MO_POR ;
  311.   'FINS';
  312. 'FINSI' ;
  313. *----------------------------------------------------------------------*
  314. *           Chargement a imposer a l'instant TI      (F^ext_n+1)       *
  315. *----------------------------------------------------------------------*
  316. ZFEXT2= 'VIDE' 'CHPOINT' / 'DISCRET';
  317. 'SI' ('EXIS' WTAB.'CHARGEMENT' 'MECA');
  318.   ZFEXT2 = 'TIRE' WTAB.'CHARGEMENT' 'MECA' TI;
  319.   TYP_2  = 'TYPE' ZFEXT2;
  320.   'SI' ('NEG' TYP_2 'CHPOINT ');
  321.     'MESS' '*** Erreur dans la definition du chargement ( MECA ) ***';
  322.     'ERRE' '*** Le type du champ n est pas CHPOINT ***';
  323.   'FINS';
  324. 'FINS';
  325. *
  326. 'SI' ('EXIS' WTAB.'CHARGEMENT' 'FORC');
  327.   F2_FOR = 'TIRE' WTAB.'CHARGEMENT' 'FORC' TI;
  328.   TYP_2 = 'TYPE' F2_FOR;
  329.   'SI' ('NEG' TYP_2 'CHPOINT ');
  330.     'MESS' '*** Erreur dans la definition du chargement ( FORC ) ***';
  331.     'ERRE' '*** Le type du champ n est pas CHPOINT ***';
  332.   'FINS';
  333.   ZFEXT2 = ZFEXT2 '+' F2_FOR ;
  334. 'FINS';
  335. *
  336. 'SI' (LOGPRE 'ET' ('NON' IGRD)) ;
  337.   MOP   = WTAB.'MOD_PRE' ;
  338.   ZPEXT ='TIRE' WTAB.'CHARGEMENT' 'PRES' TI ;
  339.   'SI' ('EXIS' WTAB 'MAT_PRE') ;
  340.     ZFPEXT = 'BSIG' MOP ZPEXT WTAB.'MAT_PRE' ;
  341.   'SINON' ;
  342.     ZFPEXT = 'BSIG' MOP ZPEXT ;
  343.   'FINS' ;
  344.   ZFEXT2 = ZFEXT2 '+' ZFPEXT ;
  345. 'FINS' ;
  346. *
  347. 'SI' ('EXIS' WTAB.'CHARGEMENT' 'DIMP');
  348.   F2_DEP = 'TIRE' WTAB.'CHARGEMENT' 'DIMP' TI;
  349.   ZFEXT2 = ZFEXT2 + F2_DEP;
  350. 'FINS';
  351. *
  352. 'SI' ('EXIS' WTAB.'CHARGEMENT' 'DINC');
  353.   F2_INC = 'TIRE' WTAB.'CHARGEMENT' 'DINC' TI;
  354.   F2_mail = extr f2_inc maillage;
  355.   zclim_mail = zclim extrai 'MAILLAGE';
  356.   F2_Mcmp = zclim_mail elem appu larg f2_mail;
  357.   zclim_inc = zclim redu f2_mcmp;
  358.   f2_base = (zclim_inc * conti.'DEPLACEMENTS') exco 'FLX' 'NOID' 'FLX';
  359.   ZFEXT2 = ZFEXT2 + F2_INC + F2_base; ;
  360. 'FINS';
  361. *
  362. 'SI' LOGDEF ;
  363.   WTAB.'DEFOR2' = 'TIRE' WTAB.'CHARGEMENT' 'DEFI' TI;
  364. 'FINS';
  365. *
  366. 'SI' WTAB.'FOR_MECA' ;
  367.   'SI' ITHER ;
  368.     'SI' WTAB.'CHAR_THE';
  369.       WTAB.'TET2' = 'TIRE' WTAB.'CHARGEMENT' 'T   ' TI;
  370.     'FINSI';
  371.     'SI' WTAB.'FOR_THER';
  372.       WTAB.'TET2' = estim.'TEMPERATURES' ;
  373.     'FINSI';
  374.  
  375. *   Champs de temperature en debut (ZTEMP1) et fin de pas (ZTEMP2)
  376. *   Si les maillages mecanique/thermique sont differents,
  377. *   on projette le champ de temperature sur le modele mecanique
  378.     'SI' WTAB.'PROJECTION';
  379.       'FORM' GEOREF0 ;
  380.       'SI' ('EGA' ('TYPE' WTAB.'TET1') 'CHPOINT') ;
  381.         ZTEMP1 = 'PROI' ('EXTR' ZMODL 'MAILLAGE') ('CHAN' 'CHAM' WTAB.'TET1' WTAB.'MOD_THE' 'NOEUD') ;
  382.         ZTEMP2 = 'PROI' ('EXTR' ZMODL 'MAILLAGE') ('CHAN' 'CHAM' WTAB.'TET2' WTAB.'MOD_THE' 'NOEUD') ;
  383.       'SINON' ;
  384.         ZTEMP1 = 'PROI' ('EXTR' ZMODL 'MAILLAGE') WTAB.'TET1' ;
  385.         ZTEMP2 = 'PROI' ('EXTR' ZMODL 'MAILLAGE') WTAB.'TET2' ;
  386.       'FINSI' ;
  387.       'FORM' GEOM1 ;
  388. *   Sinon, on vient seulement recuperer les champs de temperature
  389.     'SINON' ;
  390.       ZTEMP1 = WTAB.'TET1';
  391.       ZTEMP2 = WTAB.'TET2';
  392.     'FINSI' ;
  393.     ETT0 ='REDU' WTAB.'ETHER1' ZMODL;     'COMM' 'REDU au cas ou le modele soit parallele';
  394.     'SI' POR1;
  395.       MSRTT0='REDU' WTAB.'MSRTHER1' ZMODL; 'COMM' 'REDU au cas ou le modele soit parallele';
  396.     'FINS';
  397.   'SINON' ;
  398.     che3 = 'MANU' 'CHML' ZMODL 'T' 20. 'STRESSES' ;
  399.     che4 = 'MANU' 'CHML' ZMODL 'T' 20. 'STRESSES' ;
  400.   'FINSI';
  401. 'FINSI';
  402. *----------------------------------------------------------------------*
  403. *         Dynamique : ajout au second membre de                        *
  404. *                     F0 + 4/DT*M*V0 - B0*SIG0                         *
  405. *----------------------------------------------------------------------*
  406. * Matrice de masse
  407. 'SI' (IDYN 'OU' WTAB.'FREQUENTIEL') ;
  408.   ICALMAS = FAUX ;
  409.   'SI' ('NON' ('EXIS' WTAB 'MASSE'));
  410.     ICALMAS = VRAI ;
  411.   'FINSI' ;
  412.   'SI' (IGRD 'ET' IDYN) ;
  413.     ICALMAS = VRAI ;
  414.   'FINSI' ;
  415.   'SI' ICALMAS ;
  416.     WTAB.'MASSE' = 'MASS' WTAB.'MO_TOT' ('CHAN' 'MASSE'
  417.                          ZMAT22 WTAB.'MO_TOT' ) ;
  418.     'SI' (IDYN 'ET' WTAB.'MASSCONST') ;
  419.       WTAB.'MASSE' = WTAB.'MASSE' 'ET' WTAB.'MASSE_CONSTANTE';
  420.     'FINS';
  421.   'FINS';
  422. 'FINS';
  423. *
  424. 'SI' IDYN ;
  425.   UNSURH = 1.D0 '/' WTAB.'DT' ;
  426.   'SI' ('EGA' WTAB.'FREA1' 'INCONNU');
  427.     F1 = 'VIDE' 'CHPOINT' / 'DISCRET';
  428. *
  429. *   Forces exterieures a TEMPS0
  430.     'SI' ('EXIS' WTAB.'CHARGEMENT' 'MECA');
  431.       F1M = 'TIRE' WTAB.'CHARGEMENT' 'MECA' WTAB.'TEMPS0' ;
  432.       F1  = F1 '+' F1M ;
  433.     'FINS';
  434. *
  435.     'SI' ('EXIS' WTAB.'CHARGEMENT' 'FORC');
  436.       F1F = 'TIRE' WTAB.'CHARGEMENT' 'FORC' WTAB.'TEMPS0' ;
  437.       F1  = F1 '+' F1F ;
  438.     'FINS';
  439. *
  440.     'SI' WTAB.'PROCEDURE_CHARMECA';
  441.       TFF1 = CHARMECA PRECED TEMPS0 ;
  442.       'SI' ('EXIS' TFF1 'ADDI_SECOND');
  443.         FF1=TFF1 .'ADDI_SECOND';
  444.       'FINS';
  445.       F1 = F1 '+' FF1;
  446.     'FINS';
  447. *
  448.     'SI' LOGPRE ;
  449.       MOP = WTAB.'MOD_PRE' ;
  450.       ZPEXT0 = 'TIRE' WTAB.'CHARGEMENT' 'PRES' TEMPS0 ;
  451.       'SI' ('EXIS' WTAB 'MAT_PRE') ;
  452.         FF1 = 'BSIG' MOP ZPEXT0 ('REDU' WTAB.'MAT1' MOP) ;
  453.       'SINON' ;
  454.         FF1 = 'BSIG' MOP ZPEXT0 ;
  455.       'FINS' ;
  456.       F1 = F1 '+' FF1;
  457.     'FINS';
  458. *
  459. *   Forces interieures a TEMPS0
  460. *HHO : Modifications appel a BSIGMA (ajout du champ de deplacements necessaire en HHO)
  461.     LAF0 = 'BSIG' WTAB.'MO_TOT' conti.'CONTRAINTES' ZMAT22 conti.'DEPLACEMENTS' ;
  462.     'SI'  IRCON  ;
  463.       LAF0 = LAF0 'ET'
  464.                   ('REDU' MAI_CONS (RIG_CONS '*' conti.'DEPLACEMENTS'));
  465.     'FINS';
  466. *
  467. *   FREA1 : (masse*acceleration initiale)+(amortissement*vitesse initiale)
  468.     WTAB.'FREA1' = F1 '-' LAF0 ;
  469.   'FINS';
  470. *
  471.   'SI' IMPLP ;
  472. *  forces d'acceleration au debut du pas
  473.     'SI' ('NEG' WTAB.'AMORTISSEMENT' 'INCONNU');
  474.       FF4 = WTAB.'AMORTISSEMENT' '*' conti.'VITESSES';
  475.       WTAB.'FMAN' = WTAB.'FREA1' - FF4 ;
  476.     'SINON' ;
  477.       WTAB.'FMAN' = WTAB.'FREA1' ;
  478.     'FINS';
  479.   'FINS';
  480. *
  481.   FF = WTAB.'MASSE' '*' conti.'VITESSES';
  482.   FF4 = 4. '/' WTAB.'DT' '*' FF; 'DETR' FF;
  483. *
  484. * partie du second membre qui ne depend que des informations du pas prec
  485.   WTAB.'FREA1' = FF4 '+' WTAB.'FREA1';
  486. 'FINS';
  487. *
  488. * Advection
  489. 'SI' WTAB.'NVSTNL'                                        ;
  490.   'SI' ('EXISTE' WTAB.'CHARGEMENT' 'ADVE')                ;
  491.      UADV = 'TIRE' WTAB.'CHARGEMENT' 'ADVE' TI            ;
  492.      LOGADV = VRAI                                        ;
  493.   'SINON'                                                 ;
  494.      LOGADV = FAUX                                        ; 
  495.   'FINSI'                                                 ; 
  496. 'FINSI'                                                   ;
  497. *----------------------------------------------------------------------*
  498. *     Consolidation : ajout au second membre de                        *
  499. *                     -B0*SIG0 + DT*(1-TETA)*FI0 + DT*H*P              *
  500. *----------------------------------------------------------------------*
  501. 'SI' ISOL ;
  502.   FF  = 'BSIG' WTAB.'MOD_POR' conti.'CONTRAINTES' ;
  503.   FF4 = 'EXCO' WTAB.'MOT_POR' FF
  504.                WTAB.'MOT_POR' 'NOID' 'NATURE' 'DISCRET' ;
  505.   ZFEXT2 = ZFEXT2 '+' FF4  ; 'DETR' FF ;
  506. *
  507. *           ---- traitement des flux si besoin ----
  508.   'SI' ( 'EXIS' WTAB.'CHARGEMENT' 'FLUX' ) ;
  509.     FLUXT0= 'TIRE' WTAB.'CHARGEMENT' 'FLUX' TEMPS0;
  510.     FLUXTI= 'TIRE' WTAB.'CHARGEMENT' 'FLUX' TI  ;
  511.     FACFLU =  -1. * WTAB.'DT';
  512.     FLUXT =   ( FACFLU * (1 - WTAB.'TETA') * FLUXT0 )
  513.             + ( FACFLU * WTAB.'TETA' *FLUXTI ) ;
  514.     ZFEXT2 = ZFEXT2 '+' FLUXT ;
  515.     'DETR' FLUXT ; 'DETR' FLUXT0; 'DETR' FLUXTI;
  516.   'FINS' ;
  517. 'FINS';
  518. *----------------------------------------------------------------------*
  519. *                  Pilotage indirect : preparation                     *
  520. *----------------------------------------------------------------------*
  521. 'SI' LOGPIL ;
  522.   ZFPILIN = WTAB.'FORCES_PILOTEES' 'ET' WTAB.'DEPLACEMENTS_PILOTES' ;
  523.   DTAU = ('IPOL' WTAB.'PARAMETRE_DE_PILOTAGE' TI) -
  524.          ('IPOL' WTAB.'PARAMETRE_DE_PILOTAGE' TEMPS0) ;
  525.   ETA0 = 'EXTR' PRECED.'COEFFICIENT_DE_PILOTAGE' 
  526.                              ('DIME' PRECED.'COEFFICIENT_DE_PILOTAGE') ;
  527.   D_ETA=0.;                                                   
  528. 'FINS' ;
  529. *----------------------------------------------------------------------*
  530. *         Partie constante du second membre sur le pas de temps        *
  531. *----------------------------------------------------------------------*
  532. ZFCONST1 = ZFEXT2 ;
  533. 'SI' IDYN ;
  534.   ZDYFEXT = 'ENLE' ZFEXT2 'FLX';
  535.   ZFEXT2  = ZFEXT2 '+' WTAB.'FREA1';
  536. 'FINSI' ;
  537. ZFLXB    = 'EXCO' ZFEXT2 'FLX' 'NOID' 'FLX' ;
  538. *----------------------------------------------------------------------*
  539. *   Recuperation de valeurs de WTAB pour initialisations eventuelles   *
  540. *----------------------------------------------------------------------*
  541. ZRIBLO_M = 'MOT' 'INCONNU' ;
  542. 'SI' WTAB.'FOR_MECA' ;
  543. *
  544.   ZXDENO  = WTAB.'XDENO'  ;
  545.   ZXDENOM = WTAB.'XDENOM' ;
  546. * XDENOo  = WTAB.'XDENO';
  547. * XDENOMo = WTAB.'XDENOM';
  548.   ZETAT1 = WTAB.'ETAT1';
  549.   'SI' LOGDEF ;
  550.     ZDEFOR1 = WTAB.'DEFOR1' ;
  551.     ZDEFOR2 = WTAB.'DEFOR2' ;
  552.   'FINS';
  553.   'SI' ITHER ;
  554.     ZTET1 = WTAB.'TET1' ;
  555.     ZTET2 = WTAB.'TET2' ;
  556.   'FINS';
  557.   'SI' IPILOT ;
  558.     ZAUTOCOEF = WTAB.'AUTOCOEF' ;
  559.     ZAUTOREDU = WTAB.'AUTOREDU' ;
  560.     ZAUTORED1 = WTAB.'AUTORED1' ;
  561.   'FINS';
  562.   ZLASTKTAN = WTAB.'LASTKTAN' ;
  563. *
  564.   'SI' ('NEG' WTAB.'RIBLO_M' 'INCONNU') ;
  565.     ZRIBLO_M = WTAB.'RIBLO_M' ;
  566.     ZLISEA_M = WTAB.'LISEA_M' ;
  567.   'FINS' ;
  568. *
  569.   ZFNONL = WTAB.'FNONL' ;
  570.   'SI' ('NEG' WTAB.'INCREMENT' 'INCONNU');
  571.     INCRPREC = WTAB.'INCREMENT';
  572.   'FINS';
  573. *
  574. * teste t'on les moments  ?
  575.   TSTMOM   = 'EGA' ('VALE' 'MODE') 'PLANGENE' ;
  576.   TSTMOM = TSTMOM 'OU' ('EXIS' ZMAT22 'EPAI') 'OU' ('EXIS' ZMAT22 'INRY') 
  577.                   'OU' ('EXIS' ZMAT22 'MODS') 'OU'
  578.            (('EXIS' ZMODLI 'ELEM' 'JOI1') 'ET' ('EXIS' ZMAT22 'V1X ')
  579.                                           'ET' ('EXIS' ZMAT22 'V1Y '));
  580. * teste t'on les POREUX   ?
  581. 'SI' POR1 ; TSTMOM=VRAI ; 'FINS';
  582. *
  583.   IKLFFF=VRAI;
  584.   'SI'TSTMOM;    
  585.     'SI' IFTOL; 'SI' IMTOL;
  586.        IKLFFF=FAUX;
  587.     'FINS'; 'FINS';
  588.   'SINON' ;
  589.     'SI' IFTOL; IKLFFF=FAUX; 'FINS';
  590.   'FINS';
  591. *
  592. * Jacobien du modele pour ponderer les champs de act3
  593.   CH1  = 'MANU' 'CHML' zmodli 'SCAL' 1 'GRAVITE';
  594.   zjac = 'INTG' 'ELEM' CH1  zmodli;
  595.   zjac = zjac zmodli chan 'CHPO' 'SOMM' + xpetit;
  596.   zjac = zjac **-0.5;
  597. *
  598. * Initialisation CHArgement SANS T :
  599.   CHASANST= WTAB.'CHARGEMENT' 'ENLE' 'T' ;
  600. *
  601. * isoucomp : Indicateur SOUci dans COMPortement, FAUX par defaut 
  602.   isoucomp = faux ;
  603.   PASUNIL = FAUX;
  604.   residc  = 'VIDE' 'CHPOINT' / 'DISCRET';
  605.   RFNS    = 'VIDE' 'RIGIDITE' / 'RIGIDITE' ;
  606.   RTRSF   = 'VIDE' 'RIGIDITE' / 'RIGIDITE' ;
  607.   ZNITE   = 0 ;
  608.   ISOUSPPP = WTAB.'ISOUSPAS' ;
  609.   WTAB.'ISOUSPAS' = 0;
  610. 'FINSI' ;
  611. *
  612. RED_URG = 0;
  613. resmul  = 1;
  614. augmult = 0.30000000;
  615. XCONV   = 0.;
  616. DEPSTDM = 0.;
  617.  dpsmax = 0;
  618. WTAB.'CONV' = VRAI;
  619. KNOCONV     = 0 ;
  620. *----------------------------------------------------------------------*
  621. *              Valeurs des champs en debut de pas TEMPS0               *
  622. *----------------------------------------------------------------------*
  623. 'SI' WTAB.'FOR_MECA'  ;
  624.   ZDEP0   = CONTI.'DEPLACEMENTS' ;
  625.   ZDEF0   = 'REDU' CONTI.'DEFORMATIONS' ZMODL ;
  626.   ZSIG0   = 'REDU' CONTI.'CONTRAINTES'  ZMODL ;
  627.   com_sig = 'EXTR' ZMODLI 'CONT';
  628.   'SI' IPLAVI ;
  629.     ZDEI0 = 'REDU' CONTI.'DEFORMATIONS_INELASTIQUES' ZMODL ;
  630.     ZVAR0 = 'REDU' CONTI.'VARIABLES_INTERNES'        ZMODL ;
  631.     'SI' ('NEG' WTAB.'MOVA' 'RIEN') ;
  632.       ACC0 = 'EXCO' ZVAR0 WTAB.'MOVA' ;
  633.     'FINS' ;
  634.     com_var = 'EXTR' ZMODLI 'VARI' ;
  635.     com_dei = 'EXTR' ZMODLI 'DEIN' ;
  636.     'SI' ISOL  ;
  637.       com_maa = 'EXTR' ZMODLI 'MATE' ;
  638.     'FINS' ;
  639.   'FINS' ;
  640.   ZGRDU0 = 'MOT' 'INCONNU' ;
  641.   'SI' ('NEG' WTAB.'MOD_MEC' 'INCONNU') ;
  642.     'SI' (IGRD 'ET' WTAB.'UTILISATEUR');
  643.       'FORM' GEOREF0 ;
  644.       ZGRDU0 ='GRAD' ZMODL ZMAT2 ZDEP0 ;
  645.       'FORM' GEOM1 ;
  646.     'FINS' ;
  647.   'FINS' ;
  648. 'FINSI' ;
  649. *######################################################################*
  650. *----------------------------------------------------------------------*
  651. *       Boucle de non convergence BONOCONV (Configuration GEOM1)       *
  652. *----------------------------------------------------------------------*
  653. 'REPETER' BONOCONV WTAB.'MAXSOUSPAS' ;
  654. *
  655. * Initialisation objectif de fin de pas. Forces et deplacements
  656.   ZFCONSTA = ZFEXT2 ;
  657. *
  658.   'SI' IGRD ;
  659.     'FORM' GEOM1;
  660.   'FINSI' ;
  661. *
  662.   excfconv = 0;
  663.   'SI' ('EGA' RED_URG 0);
  664.     'SI' ((XCONV < ZPREC) 'ET' (DEPSTDM < ZPREC) 'ET' (resmul > 0.99) 'ET' (AUGMULT < 100.));
  665.       augmult = 0.30000000;
  666. *pv resmul = 1;
  667.       augauto = augmult;
  668.       augk    = augmult;
  669.  
  670.       'SI' autaug ; 
  671.         RE_CNTRL = 'VIDE' 'RIGIDITE' / 'RIGIDITE';
  672.         DE_CNTRL  = FAUX;
  673.         'SI' (dpsmax  < zprec);
  674.           iraug = faux;
  675.         'FINSI';
  676.       'FINSI';
  677.       znacce=VRAI;
  678.       'SI' (IFEFP 'OU' ISSTE) ; ZNACCE = FAUX; 'FINS';
  679.     'SINON';
  680.       augmult = augmult / 1.8;
  681.       si (augmult < 1d-3); augmult = 1d-3; finsi;
  682.     'FINSI';
  683.   'SINON';
  684.     augmult = augmult * 1.5;
  685.     'SI' autaug; IRAUG = VRAI; 'FINSI';
  686.   'FINSI';
  687. *
  688.   'SI' WTAB.'FOR_MECA' ;
  689.     KNOCONV = KNOCONV+1 ;
  690.     'SI' (knoconv '>' 1) ;
  691.       HPP_EPS = FAUX ;
  692.     'FINSI' ;
  693. *
  694.     ZDEP0N = -1. '*' ZDEP0 ;
  695.     ZD0SLX = 'ENLE' ZDEP0 'LX';
  696.     ZDEPT  = ZDEP0 '-' ZD0SLX ;
  697.     DEPST  = 'ZERO' ZMODL 'DEFORMATIONS' ;
  698.     ZSIGF  = ZSIG0 ;
  699. *  
  700.     'SI' ITHER  ;
  701.       DTETD = ZTEMP2 '-' ZTEMP1 ;
  702.     'FINS' ;
  703. *  
  704. *   --------------------------------------------------------------------
  705. *   Reevaluer les materiaux en convergence forcee (TEMPS0 =TI)
  706.     'SI' (knoconv '>' 1) ;
  707.       ZMAT11 = PAS_MATE PRECED ZETAT1 ;
  708.       'SI' WTAB.'ITCAR';
  709.         'FORM' GEOREF0 ;
  710.         GEOBID ZMAT11 = 'FORM' ZDEP0 ZMODLI ZMAT11;
  711.         'FORM' GEOM1 ;
  712.       'FINS';
  713.       ZMAT1 = 'REDU' ZMAT11 ZMODL ;
  714. *  
  715.       MCHC   = PAS_ETAT PRECED TI  ;
  716.       ZMAT22 = PAS_MATE PRECED MCHC;
  717.       'SI' (IGRD 'ET' ('EGA' LAG_TOT 1)) ;
  718.         ZMAT2I = ZMAT22 ;
  719.       'FINSI' ;
  720.       'SI' WTAB.'ITCAR';
  721.         'FORM' GEOREF0 ;
  722.         GEOBID ZMAT22 = 'FORM' ZDEP0 ZMODLI ZMAT22;
  723.         'FORM' GEOM1 ;
  724.       'FINS';
  725.       ZMAT2 = 'REDU' ZMAT22 ZMODL ;
  726.     'FINS' ;
  727. *  
  728. *   --------------------------------------------------------------------
  729. *   Rigidite a la fin du pas 
  730.     'SI' ('OU' WTAB.'RECARI' IRAUG ('NON' ('EXIS' WTAB 'RRRR')));
  731. *
  732.       'SI' ('OU' IENDOM IVIDOM ICERAM);
  733.         HOOKENDO = 'HOOK' ZMODL ZMAT2    ZVAR0 ;
  734.         RH       = 'RIGI' ZMODL HOOKENDO ZMAT2 ;
  735.         'DETR' HOOKENDO;
  736.       'SINON';
  737.         'SI' (IGRD 'ET' ('EGA' LAG_TOT 1)) ;
  738.           'FORM' GEOREF0;
  739.           HOOKRH  = 'HOOK' ZMODL ('REDU' ZMAT2I ZMODL) ;
  740.           'FORM' GEOM1;
  741.           HOOKRH2 = 'CAPI' HYPDEF HOOKRH ZMODL ZDEP0N;
  742.           RH = 'RIGI' HOOKRH2 ZMODL ZMAT2;
  743.           'DETR' HOOKRH; 'DETR' HOOKRH2;
  744.         'SINON';
  745.           RH  = 'RIGI' ZMODL ZMAT2 ;
  746.         'FINSI';
  747.       'FINS';
  748. *
  749. *     RH peut contenir des CL
  750.       ZCL = 'EXTR' RH 'RIGI' 'MULT' ;
  751.       'SI' ('NEG' ('DIME' ZCL) 0) ;
  752.         ZCLIM = ZCLIM0 'ET' ZCL ;
  753.       'FINSI' ;
  754.       RRRR = RH 'ET' ZCLIM0 ;
  755. *
  756. *     Option 'RIGIDITE_CONSTANTE'
  757.       'SI' IRCON;
  758.         RRRR = RRRR 'ET' RIG_CONS;
  759.       'FINS';
  760. *
  761. *     Option 'RIGIDITE_AUGMENTEE'
  762.       'SI' IRAUGLU ;
  763.         'SI' ('NON' (AUTAUG 'ET' IRAUG)) ;
  764.           RRRR = RRRR 'ET' RIG_AUG  ;
  765.         'FINSI';
  766.       'FINSI';
  767. *  
  768. *     Option 'AUGMENTATION_AUTOMATIQUE'
  769.       'SI' (AUTAUG 'ET' ('NON' IRAUGLU));
  770.         RIG_AUG = 'MASSE' ZMODL ZMAT2 ;
  771. ***     mess 'actualisation rig_aug';
  772.       'FINSI';
  773.       'SI' (AUTAUG 'ET' IRAUG);
  774.         RRRR = RRRR 'ET' (RIG_AUG * augauto) 'ET' (RH   * augk) ;
  775.       'FINS';
  776. *  
  777. *     Stockage de la rigidite pour eviter de la recalculer
  778.       'SI' ('EGA' ('DIME' ZCL) 0) ;
  779.         WTAB.'RRRR'=RRRR;
  780.       'FINSI';
  781. *
  782.     'FINS';
  783.     ZRAID = WTAB.'RRRR' ;
  784.  
  785.   'FINSI' ;
  786. *FOR_MECA 
  787. *
  788.   'SI' WTAB.'NVSTNL'                                        ;
  789.      'SI' ('NON' ('EXISTE' WTAB 'RNSL'))                    ;
  790.    'SI' (WTAB.'DT' '&lt;EG' 0.) ; 'MESS' 'PAS DE TEMPS !! ' WTAB.'DT' ;
  791.                                         'ERRE' 21 ; 'FINSI' ;
  792.         WTAB.'MAT_NSL' = WTAB.'MAT_NSL' 'ET' 
  793.            ('MATE' WTAB.'MOD_NSL' 'FREQ' (WTAB.'DT'**(-1))) ;
  794.         RINSL = 'RIGI' WTAB.'MOD_NSL' WTAB.'MAT_NSL'        ;
  795.         ndiv = 'RIGI' 'DIVE' WTAB.'MOD_NSL' WTAB.'MAT_NSL'  ;
  796.         WTAB.'NSMA' = 'MASSE' WTAB.'MOD_NSL' WTAB.'MAT_NSL' ;
  797.         WTAB.'RNSL' = RINSL 'ET' ZCLIM0 'ET' ndiv           ;
  798.      'FINSI'                                                ;
  799.   'FINSI'                                                   ;
  800. *
  801. *------------ consolidation ou dynamique faut-il recalculer l'operateur?
  802.   'SI' (ISOL 'OU' IDYN);
  803.     'SI' ( '>' (WTAB.'DTPREC' '*' 0.9999) WTAB.'DT') ;
  804.       WTAB.'RECAOP' = VRAI;
  805.     'FINS';
  806.     'SI' ( '<' (WTAB.'DTPREC' '*' 1.0001) WTAB.'DT') ;
  807.       WTAB.'RECAOP' = VRAI;
  808.     'FINS';
  809.     'SI' WTAB.'MATVAR';
  810.       WTAB.'RECAOP' = VRAI;
  811.     'FINS';
  812.     'SI' ('NON' WTAB.'RECAOP') ;
  813.       'SI' ('NEG' WTAB.'OPERATEUR' 'INCONNU');
  814.         ZRAID = WTAB.'OPERATEUR';
  815.       'FINS';
  816.     'FINS';
  817.   'FINS';
  818. *
  819. *------------ operateur amortissement en frequentiel
  820.   'SI' WTAB.'FREQUENTIEL' ;
  821.     RRR2 = 'AMOR' ZMODL ZMAT2 ;
  822.     RR2  = 'CHAN' 'INCO' RRR2
  823.             ('MOTS' 'ALFA' 'BETA') ('MOTS' 'IALF' 'IBET')
  824.             ('MOTS' 'FALF' 'FBET') ('MOTS' 'FALF' 'FBET') 'QUEL' ;
  825.     RR3  = 'CHAN' 'INCO' RRR2
  826.               ('MOTS' 'ALFA' 'BETA') ('MOTS' 'ALFA' 'BETA')
  827.               ('MOTS' 'FALF' 'FBET') ('MOTS' 'IFAL' 'IFBE') 'QUEL' ;
  828.     RRR2 = RR2 'ET' RR3 ;
  829.  
  830.     RR1  = ZRAID ;
  831.     OMEGI= 2.* PI * TI ;
  832.     RRR1 = OMEGI * OMEGI * (-1.) * WTAB.'MASSE' ;
  833.     RR1  = ZRAID 'ET' RRR1 ;
  834.     RR4  = 'CHAN' 'INCO' (RR1 '*' (-1.D0))
  835.             ('MOTS' 'ALFA' 'BETA') ('MOTS' 'IALF' 'IBET')
  836.             ('MOTS' 'FALF' 'FBET') ('MOTS' 'IFAL' 'IFBE') 'QUEL';
  837.     ZRAID= RR1   'ET' RR4  ;
  838.     RR5  = OMEGI '*'  RRR2 ;
  839.     ZRAID= ZRAID 'ET' RR5  ;
  840.   'FINS' ;
  841.  
  842. *--------------- et la perméabilité  ----------------------------------
  843.   'SI' ISOL ;
  844.     'SI' (IGRD 'OU' WTAB.'MATVAR') ;
  845.       WTAB.'PERMEABILITE' = 'PERM' WTAB.'MOD_POR' ZMAT2;
  846.       WTAB.'RECAOP'       =  VRAI                      ;
  847.     'FINS';
  848.   'FINS';
  849.  
  850. *------------- Cas de la consolidation ou de la dynamique -------------
  851. *------------- il faut recalculer l'operateur d'iteration -------------
  852.   'SI' WTAB.'RECAOP' ;
  853. *
  854.     'SI' IDYN;
  855.       ZRMAS = 4.D0 '/' (WTAB.'DT' '**' 2) '*' WTAB.'MASSE' ;
  856.       ZRAID = ZRMAS 'ET' ZRAID;
  857.       'SI' ('NEG' WTAB.'AMORTISSEMENT' 'INCONNU');
  858.         ZRAID = WTAB.'AMORTISSEMENT' '*' (2.D0 '/' WTAB.'DT') 'ET' ZRAID;
  859.       'FINS';
  860.     'FINS' ;
  861. *
  862.     'SI' ISOL ;
  863.       ZRAID =-1.* WTAB.'DT'* WTAB.'TETA'* WTAB.'PERMEABILITE' 'ET' ZRAID ;
  864.     'FINS' ;
  865. *
  866.     WTAB.'OPERATEUR'= ZRAID ;
  867.   'FINS';
  868. *
  869.   'SI' WTAB.'FOR_MECA' ;
  870. *
  871. *   --------------------------------------------------------------------
  872. *   Conditions de contact-frottement
  873.     'SI' WTAB.'CONTACT';
  874.       MODCON = WTAB.'MODCONTA';
  875.       'SI' ('EXIS' WTAB 'MATCONTA') ;
  876.         WTAB.'MATCONTA' = 'REDU' MODCON WTAB.'MATCONTA';
  877.         CJEU CRR = 'RFCO' MODCON  WTAB.'CONV' WTAB.'MATCONTA' ;
  878.       'SINON' ;
  879.         CJEU CRR = 'RFCO' MODCON  WTAB.'CONV'                 ;
  880.       'FINS';
  881. *
  882. *     Separer ce qui concerne l'adherence (FADH) et les jeux (FLX)
  883.       'SI' ('NEG' CJEU 0) ;
  884.         FADHE = 'EXCO' CJEU 'FADH' 'NOID' 'FADH' ;
  885.         CDAP  = 'EXCO' CJEU 'FLX'  'NOID' 'FLX' ;
  886.       'FINSI' ;
  887. *
  888.       'SI' WTAB.'MODAL' ;
  889.         CRR  = 'PJBA' CRR   ZMODL ZMAT2 ;
  890.         MCRR = 'EXTR' CRR  'MAIL' 'MULT' ;
  891.         MCDAP= 'EXTR' CDAP 'MAIL'        ;
  892.         'REPETER' BCDA ('NBNO' MCDAP) ;
  893.           PBCDA ='POIN' MCDAP &BCDA ;
  894.           PCRR  ='POIN' MCRR  &BCDA ;
  895.           CHCR  ='MANU' 'CHPO' PCRR 1 'FLX' ('EXTR' CDAP 'FLX' PBCDA)
  896.                  'NATURE' 'DISCRETE' ;
  897.           'SI' ('EGA' 1 &BCDA) ;
  898.             CCDA = CHCR ;
  899.           'SINON' ;
  900.             CCDA = CHCR 'ET' CCDA ;
  901.           'FINS' ;
  902.         'FIN' BCDA ;
  903.         CDAP = CCDA ;
  904.       'FINS' ;
  905. *
  906.       ZFCONT = 'VIDE' 'CHPOINT' / 'DISCRET' ; 
  907.       'SI' ('NEG' CRR 0) ;
  908.         ZFCONT = CRR '*' ZD0SLX 'EXCO' 'FLX' 'NOID' 'FLX';
  909.         ZCLIM  = CRR 'ET' ZCLIM ;
  910.         ZRAID  = CRR 'ET' ZRAID ;
  911.       'FINS';
  912. *
  913.       'SI' ('NEG' CJEU 0) ;
  914.         ZFCONT = ZFCONT '+' CDAP ;
  915.       'FINS';
  916. *
  917.       'SI' (('NEG' CRR 0) 'OU' ('NEG' CJEU 0)) ;
  918.         ZFCONSTA = ZFCONSTA '+' ZFCONT ;
  919.       'FINS';
  920. *
  921.     'FINS';
  922. *
  923. *   --------------------------------------------------------------------
  924. *   Matrice tangente et fefp
  925.     'SI' IKTAN ;
  926.       'SI' IFEFP ;
  927.         IKT_SAUV = VRAI ;
  928.         'SI' ('NEG' WTAB.'LASTKTAN' 'INCONNU') ;
  929.           'MESS' 'FEFP: Start with LASTKTAN' ;
  930.           ZRIKTA = ZLASTKTAN ;
  931.           ZRAID = ZCLIM 'ET' ZRIKTA ;
  932.         'SINON' ;
  933.           'MESS' 'FEFP: Previous KTAN not available' ;
  934.           ZRAID = ZRAID 'ET' ('KSIGMA' ZMODL ZSIG0 ZMAT2) ;
  935.         'FINS' ;
  936.       'SINON' ;
  937.         'SI' ('NEG' WTAB.'LASTKTAN' 'INCONNU') ;
  938.           IKT_SAUV = VRAI ;
  939.           'SI' IPERT ;
  940.             'MESS' 'Matrice tangente par perturbation - '
  941.             'Demarrage avec KTAN = LASTKTAN' ;
  942.           'SINON' ;
  943.             'MESS' 'Matrice tangente "coherente" - '
  944.             'Demarrage avec KTAN = LASTKTAN' ;
  945.           'FINS' ;
  946.           ZRIKTA = ZLASTKTAN ;
  947.           ZRAID = ZCLIM 'ET' ZRIKTA ;
  948.         'SINON' ;
  949.           'SI' IPERT ;
  950.             IKT_SAUV='NEG' WTAB.'K_TANGENT_ITER0'
  951.                                      'MAT_ELASTIQUE';
  952.             'MESS' 'Matrice tangente par perturbation - '
  953.             'Demarrage avec KTAN = rigidite elastique' ;
  954.           'SINON' ;
  955.             IKT_SAUV = ('NEG' WTAB.'K_TANGENT_ITER0' 'MAT_ELASTIQUE')
  956.                  'ET' ('NEG' WTAB.'K_TANGENT_ITER0' 'MAT_TANGENTE') ;
  957.             'SI' ('EGA' WTAB.'K_TANGENT_ITER0' 'MAT_ELASTIQUE') ;
  958.               'MESS' 'Matrice tangente "coherente" - '
  959.                               'Demarrage avec KTAN = rigidite elastique' ;
  960.             'SINON' ;
  961.               'MESS' 'Matrice tangente "coherente" - '
  962.                                  'Demarrage avec KTAN (DTTAN = 0.)' ;
  963.               DTTAN = 0. ;
  964.               ZRIKTA = 'KTAN' ZMODL ZSIG0 ZVAR0 ZMAT2
  965.                               'PREC' ZPREK 'DT  ' DTTAN ZKTASYM ;
  966.               ZRAID = ZCLIM 'ET' ZRIKTA ;
  967.             'FINS' ;
  968.           'FINS' ;
  969.         'FINS' ;
  970.       'FINS' ;
  971.     'FINS' ;
  972. *
  973. *   --------------------------------------------------------------------
  974. *   Raideur geometrique (option K_SIGMA)
  975.     'SI' (IGRD 'ET' ('NON' HPP_EPS)) ;
  976.       'SI' (IKSIA 'ET' ('NON' IFEFP)) ;
  977. *       Passage de GEOREF0 -> GEOM1
  978.         'FORM' GEOM1;
  979.         ZSIGKS   = 'CAPI' HYPDEF ZSIG0 ZMODL ZDEP0N ;
  980.         KSIGTC   = 'KSIGMA' ZMODL ZSIGKS ZMAT2;
  981.         ZRAIDINI = ZRAID           ;
  982.         ZRAID    = ZRAID 'ET' KSIGTC;
  983.       'FINS' ;
  984.     'FINS' ;
  985. *----------------------------------------------------------------------*
  986. *                  Pilotage automatique 
  987. *----------------------------------------------------------------------*
  988.     COEPI = 1.D0 ; CORPREC = 1. ;
  989. *
  990.     'SI' IPILOT ;
  991.       AL1 = 1. ; COEPI0 = 1.d0;
  992.       'SI' WTAB.'AUTODEUX' ;
  993.         COEPI = 'ABS' ZAUTOCOEF ;
  994.         COEPI = COEPI / (1.-COEPI);
  995.         'SI' (COEPI > 1.D0) ; COEPI=1.D0;'FINS';
  996.         COEPI0=COEPI;ZAUTOCOEF =COEPI;
  997.       'SINON';
  998.         ZAUTOCOEF = 1.D0;
  999.       'FINS' ;
  1000.       RED1 = 1. ;
  1001.       RED2 = 0 ;
  1002.       'SI' ('NEG' ZAUTORED1 'INCONNU') ;
  1003.         'SI' (ZAUTORED1 > 0);
  1004.           ZAUTORED1 = ZAUTORED1 - 1;
  1005.           'SI' (ZAUTORED1 'EGA' 0) ;
  1006.             COEPI = 3 *  COEPI ;
  1007.             ZAUTOREDU = ZAUTOREDU / 3.;
  1008.             'SI' (ZAUTOREDU > 1.1 );
  1009. *             on travaille encore avec un critere reduit
  1010.               ZAUTORED1 = 4 ;
  1011.               RED1 = 3. ;
  1012.             'FINS';
  1013.             'MESS' 'On multiplie le critere de pilotage par 3';
  1014.           'FINS';
  1015.         'FINS';
  1016.       'FINS';
  1017.       'SI' (COEPI > 1d0);
  1018.         RED1 = RED1 / COEPI ;COEPI =1d0;
  1019.       'FINS';
  1020.       'SI' ('NEG' WTAB.'NBPLAS' 'INCONNU') ;
  1021.         'SI' (WTAB.'NBPLAS' 'EGA' 0) ;
  1022.           'SI' (COEPI < 0.) ;COEPI = COEPI * -2.;'FINS';
  1023.         'FINS';
  1024.       'FINS';
  1025.       COEPI = 'ABS' COEPI ;
  1026.       COEPI0 = COEPI;
  1027.     'SINON';
  1028.       ZAUTOCOEF = 1.D0 ;
  1029.     'FINS';
  1030. *
  1031.   'FINSI' ;
  1032. *FOR_MECA 
  1033. *----------------------------------------------------------------------*
  1034. *             Partie variable du chargement a imposer a TI             *
  1035. * (geometrie, materiaux, chargement en deformations)                   *
  1036. *----------------------------------------------------------------------*
  1037.   XDENDEF = 0. ; FDENDEF = 0. ; 
  1038.   DEPST0  = 0. ; DMSRT0  = 0. ; 
  1039.   DFTHE   = 0. ; DFDEF   = 0. ;
  1040.   ZFSUIV  = 0. ;
  1041. *
  1042. * - Corrections en cas de materiaux variables
  1043.   'SI' WTAB.'MATVAR' ;
  1044.     'SI' ('OU' IENDOM IVIDOM ICERAM);
  1045.       XXX3 ='ELAS' ZMODL ZSIG0 ZMAT2 ZVAR0;
  1046.       XXX4 ='ELAS' ZMODL ZSIG0 ZMAT1 ZVAR0;
  1047.     'SINON';
  1048.       XXX3 ='ELAS' ZMODL ZSIG0 ZMAT2      ;
  1049.       XXX4 ='ELAS' ZMODL ZSIG0 ZMAT1      ;
  1050.     'FINS';
  1051.     DDEF0  = XXX4 - XXX3;
  1052.     'DETR' XXX3; 'DETR' XXX4;
  1053.     DEPST0 = -1.* DDEF0;
  1054.   'FINS';
  1055. *
  1056. * - En cas de chargement thermiques
  1057.   'SI' ITHER ;
  1058.     'SI' ('EGA' ('TYPE' ZTEMP2) 'CHPOINT');
  1059.       MCHTETA2 ='CHAN' 'CHAM' ZTEMP2 ZMODL 'STRESSES' 'TEMPERATURES' ;
  1060.     'SINON';
  1061.       MCHTETA2 = ZTEMP2;
  1062.     'FINS';
  1063.     ETREF = 'REDU' WTAB.'ETREF' ZMODL ; 'COMM' 'REDU au cas ou le modele soit parallele';
  1064.     ETT   =('EPTH' ZMODL ZMAT2 MCHTETA2 ) - ETREF ;
  1065.     DTT   =  ETT '-' ETT0 ;
  1066.     XDENDEF = XDENDEF '+' ETT ;
  1067.  
  1068.     'SI' WTAB.'POR1';
  1069. *------   Cas du milieu poreux avec chargement thermique ----------
  1070. *         cas isotrope seulement pour le moment
  1071. *         et on  ne s'occupe pas du alpha-reference !!
  1072.           MSRTT = PAS_EPTH PRECED ZMODL ZMAT2 MCHTETA2 ;
  1073.           DMSRT0= MSRTT '-' MSRTT0 ;
  1074.     'FINS' ;
  1075.  
  1076.     DEPST0 = DTT '+' DEPST0 ;
  1077.   'FINS';
  1078. *
  1079. * - Contraintes et forces equivalentes associees a ces deformations
  1080.   'SI' (DEPST0 'NEG' 0);
  1081.     'SI' ('OU'  IENDOM IVIDOM ICERAM);
  1082.       DSIGT0 = 'ELAS' ZMODL DEPST0 ZMAT2 ZVAR0;
  1083.       'SI' (XDENDEF 'NEG' 0) ;
  1084.         XDENDEF = 'ELAS' ZMODL XDENDEF ZMAT2 ZVAR0;
  1085.       'FINSI' ;
  1086.     'SINON';
  1087.       DSIGT0 = 'ELAS' ZMODL DEPST0 ZMAT2 ;
  1088.       'SI' (XDENDEF 'NEG' 0) ;
  1089.         XDENDEF = 'ELAS' ZMODL XDENDEF ZMAT2 ;
  1090.       'FINSI' ;
  1091.     'FINS' ;
  1092.  
  1093.     'SI' (DMSRT0 'NEG' 0);
  1094.       DSIGT0 = DSIGT0 '+' DMSRT0 ;
  1095.       XDENDEF = XDENDEF '+' MSRTT ;
  1096.     'FINS';
  1097.     DFTHE ='BSIG' ZMODL DSIGT0 ZMAT2 ;
  1098.     'SI' (XDENDEF 'NEG' 0) ;
  1099.       FDENDEF = 'BSIG' ZMODL XDENDEF ZMAT2 ;
  1100.     'FINSI' ;
  1101.   'FINS';
  1102. *
  1103. * - Chargement en deformations imposees
  1104.   'SI' LOGDEF;
  1105.     DDEFOR0 = 'REDU' (ZDEFOR2 - ZDEFOR1) ZMODL;
  1106.     XDENDEF = 'REDU' ZDEFOR2 ZMODL ;
  1107.     'SI' ('OU' IENDOM IVIDOM ICERAM);
  1108.       DSI1 = 'ELAS' ZMODL DDEFOR0 ZMAT2 ZVAR0;
  1109.       XDENDEF = 'ELAS' ZMODL XDENDEF ZMAT2 ZVAR0;
  1110.     'SINON';
  1111.       DSI1 = 'ELAS' ZMODL DDEFOR0 ZMAT2    ;
  1112.       XDENDEF = 'ELAS' ZMODL XDENDEF ZMAT2 ;
  1113.     'FINS';
  1114.     DFDEF = 'BSIG' ZMODL DSI1 ZMAT2 ;
  1115.     FDENDEF = FDENDEF '+' ('BSIG' ZMODL XDENDEF ZMAT2) ;
  1116.   'FINS';
  1117. *
  1118. * - Chargement de pression suiveuse en grands deplacements
  1119.   'SI' (LOGPRE 'ET' IGRD) ;
  1120.     MOP = WTAB.'MOD_PRE' ;
  1121.     ZPEXTF = 'TIRE' WTAB.'CHARGEMENT' 'PRES' TI     ;
  1122.     'SI' ('EXIS' WTAB 'MAT_PRE') ;
  1123.       ZFPEXTF = 'BSIG' MOP ZPEXTF WTAB.'MAT_PRE' ;
  1124.     'SINON' ;
  1125.       ZFPEXTF = 'BSIG' MOP ZPEXTF ;
  1126.     'FINS' ;
  1127.     ZFSUIV = ZFSUIV '+' ZFPEXTF   ;
  1128.   'FINS';
  1129. *
  1130. * - Procedure UTILISATEUR : y a-t-il des forces non conservatives
  1131.   'SI' WTAB.'PROCEDURE_CHARMECA';
  1132. *   on ajoute l indice ADDI_MATRICE pour signaler a charmeca qu on
  1133. *   souhaite aussi l operateur linearisé des Forces NL de charmeca
  1134.     PRECED.'ADDI_MATRICE' = vrai;
  1135.     TFP22  = CHARMECA PRECED TI ;
  1136.     PRECED.'ADDI_MATRICE' = faux;
  1137. *
  1138. *   FP22 = F^suiv_n+1
  1139.     ADDISEC2 = 'EXIS' TFP22 'ADDI_SECOND' ;
  1140.     'SI' ADDISEC2 ;
  1141.       FP22   = TFP22.'ADDI_SECOND' ;
  1142.       FP022  = 'COPIER' FP22       ;
  1143.       ZFSUIV = ZFSUIV '+' FP22     ;
  1144.     'FINS';
  1145.     'SI' ('EXIS' TFP22 'ADDI_MATRICE');
  1146.       ZRAID = ZRAID 'ET' TFP22.'ADDI_MATRICE';
  1147.     'FINS';
  1148.   'FINS';
  1149. *
  1150. * A-t-on des C.L. unilaterales ?
  1151.   XXX1   = 'EXTR' ZCLIM 'MAIL' 'UNIL';
  1152.   IMPO12 = ('NBEL' XXX1) '>' 0 ;
  1153. *----------------------------------------------------------------------*
  1154. *                      Second membre    RESIDU                         *
  1155. *                                                                      *
  1156. * Calcul du premier residu : desequilibre entre les forces externes et *
  1157. * le calcul B*SIGMA. Le sigma qui sert est celui qui existerait si     *
  1158. *   le champ de deplacement ne changeait pas (ZDEP0).                  *
  1159. *                                                                      *
  1160. *          | forces exterieures sans reactions - forces interieures |  *
  1161. *          |   F^ext_n+1 + DF^suiv + DF^ther + DF^defi - F^int_n    |  *
  1162. *          |                         |-------------------------|    |  *
  1163. *          |                                  -1.*ZFPLO             |  *
  1164. *          |          (composantes de forces FX FY FZ ...)          |  *  
  1165. * RESIDU = |                                                        |  *
  1166. *          |      increment des relations imposees Du^imp           |  *
  1167. *          |          (composantes de depl. FLX)                    |  *
  1168. *                                                                      *
  1169. *  A F^ext peuvent s'ajouter des termes supplementaires p.ex. en       *
  1170. *   dynamique ou en poreux.                                            *
  1171. *----------------------------------------------------------------------*
  1172. * ZFCONSTA = [ F^ext_n+1 ; u^imp_n+1 ]
  1173.   ZFEXT = 'ENLE' ZFCONSTA 'FLX';
  1174.   ZFLX1 = 'EXCO' ZFCONSTA 'FLX' 'NOID' 'FLX' 'NATURE' 'DISCRET';
  1175. *
  1176.   'SI' WTAB.'FOR_MECA' ;
  1177. *
  1178. *   Forces externes deja equilibrees au debut du pas par B*sigma
  1179. *   ZF1 = F^int_n = B*sigma_n + K^cst*u_n
  1180. *HHO : Modifications appel a BSIGMA (ajout du champ de deplacements necessaire en HHO)
  1181.  
  1182.  
  1183. *  transport de zsig0 sur config debut de pas
  1184.     ZSIGBS = ZSIG0; 
  1185. *   Passage de GEOREF0 -> GEOM1
  1186.     'SI' IGRD ;
  1187.       'FORM' GEOM1;
  1188.       ZSIGBS = 'CAPI' HYPDEF ZSIG0 ZMODL ZDEP0N ;
  1189.     'FINSI';
  1190.     ZF1 = 'BSIG' ZMODL ZSIGBS ZMAT1 ZDEP0 ;
  1191.     'SI' IRCON;
  1192.       ZF1 = ZF1 '+' ('REDU' MAI_CONS (RIG_CONS '*' ZDEP0));
  1193.     'FINS';
  1194. *
  1195.     'SI' IDYN ; 
  1196.       FFDYN = 'COPIER' ZF1; 
  1197.     'FINS';
  1198. *
  1199.     'SI' ISOL ;
  1200.       GRAP0= 'GRAD' MO_POR ZDEP0 MA_POR0 'CONS' ;
  1201.       XXX1 = 'GRAD' MO_POR ZDEP0 MA_POR  'CONS' ;
  1202.       XXXS =((1.- WTAB. 'TETA' )*GRAP0)+ (WTAB. 'TETA' * XXX1);
  1203.       XXX2 = WTAB.'DT' '*' ('GNFL' MO_POR XXXS) ;
  1204.       XXX3 = ZF1  ;
  1205.       ZF1  = XXX3 - XXX2;
  1206.       'DETR' XXX3; 'DETR' XXX2 ;
  1207.     'FINS';
  1208. *
  1209. *   [ FLXINI ; FREAP ] = [ u^imp_n ; -F^reac_n ]
  1210.     FZDEP0 = ZDEP0 '*' ZCLIM;
  1211.     FLXINI = 'EXCO' FZDEP0 'FLX' 'NOID' 'FLX' 'NATURE' 'DISCRET';
  1212.     FREAP  = 'ENLE' FZDEP0 'FLX';
  1213. *
  1214. *   XXX1 = [ F^ext_n+1 ; u^imp_n+1 - u^imp_n ]
  1215.     XXX1 = ZFCONSTA '-' FLXINI ;
  1216. *
  1217. *   FEXT0 : chargement externe (sans reactions) au pas precedent
  1218. *   FEXT0 = F^int_n - F^reac_n = F^ext_n
  1219.     FEXT0 = ZF1 '+' FREAP;
  1220. *
  1221. *   RESIDU : forces exterieures sans reactions (avec des termes 
  1222. *   supplementaires le cas echeant p.ex. en dynamique ou en poreux)
  1223. *   - forces interieures et increment des relations imposees
  1224. *
  1225. *   RESIDU = [DF^tot ; Du^imp]
  1226. *   RESIDU = [F^ext_n+1 + F^suiv_n+1 + DF^ther + DF^defi - F^int_n ; Du^imp]
  1227.     RESIDU = XXX1 '+' DFTHE '+' DFDEF '-' ZF1;
  1228.     RESIDU = RESIDU '+' ZFSUIV ;
  1229. *
  1230. *   ZFPLO  = -1.* [- F^int_n + DF^ther + DF^defi ]
  1231.     ZFPLO  = ZF1 '-' DFTHE '-' DFDEF ;
  1232. *
  1233. *   ZDFINI = [F^ext_n+1 + DF^ther + DF^defi - (F^int_n - F^reac_n) ; Du^imp]
  1234.     ZDFINI = XXX1 '-' FEXT0 ;
  1235.     'SI' IPILOT ;
  1236.       ZLX = 'EXCO' ZDEP0 'LX' 'NOID' 'LX' 'NATURE' 'DIFFUS';
  1237.       DFEXT0F = 'ENLE' ZDFINI 'FLX';
  1238.       DFEXT0L = 'EXCO' ZDFINI 'FLX' 'NOID' 'FLX';
  1239.       'SI' IMPO12 ;
  1240.         DUUNIL = 'REDU' DFEXT0L XXX1 ;
  1241.         XXX2   = 'DIFF' ('EXTR' ZCLIM 'MAIL') XXX1 ;
  1242.         DUIMPO = 'REDU' DFEXT0L XXX2 ;
  1243.       'FINS';
  1244.     'FINSI' ;
  1245.     'SI' (ITHER 'OU' WTAB.'MATVAR');
  1246.       ZDFINI = ZDFINI '+' DFTHE;
  1247.     'FINS';
  1248.     'SI' LOGDEF;
  1249.       ZDFINI = ZDFINI '+' DFDEF;
  1250.     'FINS';
  1251. *
  1252.     ZSDMBR = RESIDU '*' 1.D0;
  1253. *
  1254.   'FINSI' ;
  1255. * FOR_MECA
  1256. *
  1257.   'SI' WTAB.'NVSTNL'                          ;
  1258.     ZVIFL = conti.'VITESSES_FLUIDE'           ;
  1259.     movifl = 'EXTR' WTAB.'MOD_NSL' 'DEPL'     ;
  1260.     UN = 'EXCO' ZVIFL movifl                  ;
  1261.     ch_un = 'MANU' 'CHML' WTAB.'MOD_NSL' 'ADVE' un 'TYPE' 'CARACTERISTIQUES';
  1262.     nugrad = 'ADVE' WTAB.'MOD_NSL' ch_un      ;
  1263.     'SI' ('EXIS' PRECED 'OPTION_BDF2') ;
  1264.       ZRNSL = WTAB.'RNSL' 'ET' nugrad 'ET' (WTAB.'NSMA' '*' 1.5)  ;
  1265.       unm = 'EXCO' conti.'VITESSES_FLUIDE_0' movifl ;
  1266.       fmass = WTAB.'NSMA' '*' ((2. '*' un) '-' (0.5 '*' unm))  ;
  1267.     'SINON' ;
  1268.       ZRNSL = WTAB.'RNSL' 'ET' nugrad 'ET' WTAB.'NSMA' ;
  1269.       fmass = WTAB.'NSMA' * un                ;
  1270.     'FINSI' ;
  1271.     ftnsl  = ZFCONSTA 'ET' fmass              ;
  1272.  
  1273.     'SI' LOGADV                               ;
  1274.       ch_ad = 'MANU' 'CHML' WTAB.'MOD_NSL' 'ADVE' UADV 'TYPE' 'CARACTERISTIQUES';
  1275.       nuadv = 'ADVE' WTAB.'MOD_NSL' ch_ad     ;
  1276.       ZRNSL = ZRNSL 'ET' nuadv                ;
  1277.       fadv = nugrad * UADV                    ;
  1278.       mess 'advec / inertie' (maxi (abs fadv)) (maxi (abs ftnsl)) (maxi (abs UADV)) ;
  1279.       ftnsl = ftnsl '-' fadv                  ;
  1280.     'FINSI'                                   ;
  1281.  
  1282.     resnsl = ftnsl '-' (ZRNSL * ZVIFL)        ;
  1283.  
  1284. * resolution 
  1285.     ZVIF1  = 'KRES' ZRNSL resnsl              ; 
  1286.     ZVIFL = ZVIFL '+' ZVIF1                   ;
  1287.   'FINSI'                                     ;
  1288. *  
  1289. *----------------------------------------------------------------------*
  1290. *           1ere Resolution  - ZDEP1 = [ Du^0 ; LX_n+1 ]               *
  1291. *                                                                      *
  1292. *           -> estimation a TI a partir de l'instant TEMPS0            *
  1293. *                                                                      *
  1294. *             K*Du^0 + At*LX_n+1 = F^ext_n+1 - Bsigma_n                *
  1295. *             A*Du^0             = d1 - A*u_n                          *
  1296. *             |----------------|  |-------------------|                *
  1297. *               ZRAID * ZDEP1    =       RESIDU                        *
  1298. *----------------------------------------------------------------------*
  1299.   'SI' WTAB.'FOR_MECA' ;
  1300.     FEXCI = 'VIDE' 'CHPOINT' / 'DISCRET' ;
  1301.     'SI' WTAB.'ADHERENCE' ;
  1302.       'SI' ('EXIS' WTAB 'MAIL_BLOM') ;
  1303.         FADRE = 'REDU' FADHE WTAB.'MAIL_BLOM' ;
  1304.       'SINON' ;
  1305.         FADRE = FADHE ;
  1306.       'FINSI' ;
  1307.       FEXCI = FEXCI 'ET' FADRE ;
  1308.     'FINSI' ;
  1309. *  
  1310. *   Modification du residu si automatique
  1311.     'SI' LOGPIL ;      
  1312.       RESIDU = RESIDU 'ET' (ETA0 '*' ZFPILIN);        
  1313.     'FINSI' ;
  1314. *  
  1315. *   Force limite de frottement
  1316.     'SI' (WTAB.'CAFROTTE' 'ET' IMPO12);
  1317.       'SI' WTAB.'FROCABL' ;
  1318.         zsigfT = 'REDU' zsigf ZMODLI;
  1319.         FFROT = 'EXCC' ZCLIM ZDEPT MODCON WTAB.'MATCONTA'  ZSIGFT;
  1320.         RTRSF = 'VIDE' 'RIGIDITE' / 'RIGIDITE' ;
  1321.         RFNS  = 'VIDE' 'RIGIDITE' / 'RIGIDITE' ;
  1322.         FEXCI = FEXCI '+' FFROT ;
  1323.       'FINSI';
  1324.       'SI' WTAB.'FROCOUL' ;
  1325.         excfconv = excfconv + 1 ;  
  1326.         RTRSF*'RIGIDITE' RFNS*'RIGIDITE'  FFROT = 'EXCF' ZCLIM ZDEPT MODCON WTAB.'MATCONTA'  excfconv;
  1327.         FEXCI = FEXCI '+' FFROT ;
  1328.         MFROT = 'EXTR' FFROT 'MAIL' ;
  1329.       'FINSI';
  1330.     'FINS';
  1331. *  
  1332.     'SI' ('NEG' ZRIBLO_M 'INCONNU') ;
  1333.       'SI' ((WTAB.'CAFROTTE' 'OU' WTAB.'ADHERENCE') 'ET' IMPO12);
  1334.         ZDEP1 BID BID BID BID = 'RESO' ZRAID 'SOUC' RESIDU 'INIB'
  1335.                                        ZRIBLO_M ZLISEA_M FEXCI RFNS RTRSF ;
  1336.       'SINON';
  1337.         ZDEP1 BID BID BID     = 'RESO' ZRAID 'SOUC' RESIDU 'INIB'
  1338.                                               ZRIBLO_M ZLISEA_M ;
  1339.       'FINS';
  1340.       ZRIBLO_M = 'MOT' 'INCONNU' ;
  1341.     'SINON';
  1342.       'SI' ((WTAB.'CAFROTTE' 'OU' WTAB.'ADHERENCE') 'ET' IMPO12);
  1343.         ZDEP1 = 'RESO' ZRAID 'SOUC' RESIDU FEXCI RFNS RTRSF          ;
  1344.       'SINON';
  1345.         ZDEP1 = 'RESO' ZRAID 'SOUC' RESIDU;
  1346.       'FINS';
  1347.     'FINS';
  1348. *  
  1349.     'SI' LOGPIL ; 
  1350.       ZDEPII  = 'RESO' ZRAID ZFPILIN ; 
  1351.       ZDEPILO = 0. * ZDEPII ;     
  1352. *     Calcul de D_eta par appel a la procedure PILOINDI
  1353.       D_ETA = PILOINDI PRECED ZDEP0 ZDEPILO ZDEP1 ZDEPII DTAU; 
  1354. *     Mise à jour  
  1355.       ZDEP1 = ZDEP1 '+' (D_ETA '*' ZDEPII) ;       
  1356.       ETA = ETA0 ;  
  1357.     'FINSI' ;
  1358. *
  1359. *   Pour initialiser le statut des CL unil lors du prochain appel a RESO
  1360.     ZRAID_T = 'EXTR' ZRAID 'CONT';
  1361.     'SI' ('EXIS' ZRAID_T 'NITER');
  1362. *     Conserver les conditions unilaterales pour VITEUNIL
  1363.       'SI' IDYN;
  1364.         WTAB.'ZRAIDV' = ZRAID;
  1365.       'FINS';
  1366.       ZRIBLO_M = ZRAID_T. 7 ;
  1367.       ZLISEA_M = ZRAID_T. 6 ;
  1368.     'FINS';
  1369. *
  1370. *   Pour eventuel calcul de augauto
  1371.     zdep1d = zdep1 'ENLE' 'LX';
  1372. *   et initialisation zdepl (sert pour xnum)
  1373.     ZDEPL = 'EXCO' ZDEP1 'LX' 'NOID' 'LX' ;     
  1374. *----------------------------------------------------------------------*
  1375. *              Calcul d'une norme pour la convergence                  *
  1376. *----------------------------------------------------------------------*
  1377.     XXX1 = ZFEXT '+' ZFSUIV ;
  1378.     'SI' LOGPIL ;      
  1379.       XXX1 = XXX1 '+' ((ETA '+' D_ETA) '*' ZFPILIN);  
  1380.     'FINSI' ;
  1381. *
  1382.     ZDEP1P50 = ZDEP1       ;
  1383.     XDENO='XTY' ZDEP1P50  ( XXX1 -( RESIDU  'EXCO'
  1384.         'FLX' 'NOID' 'FLX' 'NATURE' 'DISCRET')) MLPRIM MLDUAL;
  1385.     'SI' ('VERI' xdeno) ; 'SINON'; xdeno = 1; 'FINSI';
  1386.     MZDEP1M = 'MAXI' ZDEP1P50 'ABS' 'AVEC' MLDEPL;
  1387.     MZFM = 'MAXI' (FDENDEF + DFTHE + DFDEF + ZF1) 'ABS' 'AVEC' MLDUAL;
  1388.     XDENO1 = 'ABS' XDENO  + (MZFM * MZDEP1M);
  1389.     MZDEP1M = MZDEP1M + XPETIT;
  1390.     XDENO=XDENO1/MZDEP1M;
  1391.     XDENO = XDENO + MZFM;
  1392.     XDENO = 'MAXI' XDENO ('MAXI' XXX1 'ABS' 'AVEC' MLDUAL);
  1393.     'SI' (XDENO < XPETIT); XDENO = 1.; 'FINSI';
  1394.     XDENOM=XDENO;
  1395.     'SI' TSTMOM ;
  1396.       ZDEP1P50 = ZDEP1 + XPETIT;          ;
  1397.       XDENOM = XDENO1/('MAXI' ZDEP1P50 'ABS' 'AVEC' MLROTA);
  1398.       XDENOM = XDENOM + XPETIT       ;
  1399.     'FINS' ;
  1400. *
  1401. **  mess 'xdeno propose' ' ' xdeno ' xdeno precedent' ' ' xdenoo;
  1402. **     XDENO = MAXI XDENOo XDENO;
  1403. **     XDENOo = XDENO;
  1404. **     XDENOM = MAXI XDENOMo XDENOM;
  1405. **     XDENOMo = XDENOM;
  1406. *
  1407.     'SI' WTAB.'CONV';
  1408.       ZINCREMENT = ZDFINI '+' ZFSUIV ; 
  1409.     'SINON' ;
  1410. *     'MESS' 'non convergence on garde xdeno xdenom ' xdeno xdenom;
  1411.       'SI' (XDENO < ZXDENO); XDENO = ZXDENO ; 'FINSI';
  1412.       'SI' (XDENOM < ZXDENOM); XDENOM = ZXDENOM ; 'FINSI';
  1413.     'FINS';
  1414. *
  1415.     RESIDNOR = 'COPIER' RESIDU ;
  1416. *
  1417.     'SI' IPILOT;
  1418.       'SI' WTAB.'AUTODEUX' ;
  1419.         XDENO  = ZXDENO ;
  1420.         XDENOM = ZXDENOM ;
  1421.       'FINS';
  1422.       DFEXT = DFEXT0F '*' COEPI ;
  1423.       ZFEXT = DFEXT   '+' FEXT0 ;
  1424.       XXX1  = RESIDU  '*' COEPI; 
  1425.       XXX1  = XXX1 '+' ((1. '-' COEPI) '*' FREAP) ; 
  1426.       'SI' IMPO12;
  1427.         RESIDU = 'ENLE' XXX1 'FLX' ;
  1428.         XJEUX  = DUIMPO '*' COEPI ;
  1429.         RESIDU = RESIDU '+' XJEUX '+' DUUNIL ;
  1430.       'SINON' ;
  1431.         RESIDU = XXX1 ;
  1432.       'FINS';
  1433.     'FINS';
  1434. *----------------------------------------------------------------------*
  1435. *            Corriger le residu a partir du pas precedent              *
  1436. *----------------------------------------------------------------------*
  1437. *   petite correction du residu pour esperer gagner du temps ***********
  1438.     INIT = FAUX ;
  1439.     'SI'(('NEG' ZFNONL 'INCONNU') 'ET' WTAB.'INITIALISATION');
  1440.       'SI' IPILOT;
  1441.         'SI' ( WTAB.'AUTODEUX' 'ET' (COEPI 'NEG' 1.D0)) ;
  1442.           'MESS' 'Initialisation a partir du pas precedent 'COEPI;
  1443.           INIT = VRAI;
  1444.           XXX1 = RED1 * ZFNONL ;
  1445.           XXX2 = XXX1 + RESIDU ;
  1446.           'DETR' RESIDU ;'DETR' XXX1 ;
  1447.           RESIDU = XXX2;
  1448.         'FINS' ;
  1449.       'SINON';
  1450. * on fait la correction si le pas precedent a converge sans non convergence
  1451.         'SI' (WTAB.'CONV' 'ET' (ISOUSPPP 'EGA' 0));
  1452. * on fait la correction si le pas precedent etait non lineaire.
  1453.           'SI' (('MAXI' 'ABS' ZFNONL) > (ZPREC * XDENO)) ;
  1454. * on enleve le residu du pas precedent pour recuperer l'increment
  1455. * nominal du second membre e imposer f1.
  1456. * f2 et l'increment du second membre du pas precedent
  1457.             ZINCREMENT = ZINCREMENT - WTAB.'RESIDU' ;
  1458.             zdeps = WTAB.'ZDEP1' + zdep1;
  1459.             FFNO= 'XTY' ZFNONL zdeps  MNDUAL MNPRIM;   
  1460.             F1 = ZINCREMENT ;
  1461.             F2 = INCRPREC ;
  1462.             f12 = 'XTY' f1 zdeps MNDUAL MNPRIM;                
  1463.             f22 = 'XTY' f2 zdeps MNDUAL MNPRIM;                
  1464.             AMPL = f12/(f22 + XPETIT); AMPLT = 0;
  1465.             'SI' ('NEG' WTAB.'DTPREC' 0);
  1466.               DTPREC = WTAB.'DTPREC' ;
  1467.               'SI' (DTPREC > XPETIT);
  1468.                 AMPLT = WTAB.'DT' '/' DTPREC;
  1469.               'FINS';
  1470.             'FINS';
  1471. *  Le chargement n'est il pas de fluage ou de thermique ?';
  1472.             XDCOMP = ('XTY' ZDEP1 F1 MNPRIM MNDUAL) ;  
  1473.             'SI' (('ABS' XDCOMP) < (ZPREC * XDENO * mzdep1m)
  1474.                       'OU' (('ABS' FFNO) > (('ABS' F22) * 2.e2   )));
  1475.               DTPREC = WTAB.'DTPREC' ;
  1476.               AMPL=AMPLT;
  1477. * la decharge est-elle significative
  1478.               'SI'((F12/(f22 + XPETIT)) < -0.05);ampl=0.;'FINS';
  1479. **            'MESS' 'F12 F22_' f12 f22;
  1480.               'MESS' 'Pas d increment de charge, initialisation calculee avec le temps'     ;
  1481.             'SINON';
  1482.               AMPL = F12 / (F22 + XPETIT);
  1483.             'FINS';
  1484. * changement de modele on n'initialise pas
  1485.             'SI' ('NEG' ZMODLI WTAB.'MO_TOT_PREC');  AMPL = 0; 'FINSI';
  1486.               AMPL = MINI (prog AMPL AMPLT);
  1487.               'SI' ((AMPL > 0) 'ET' (AMPL < 2e1)) ;
  1488.       'MESS' 'Initialisation a partir de la solution precedente Coeff'AMPL;
  1489.                 XXX1 = AMPL * ZFNONL ;
  1490.                 XXX2 =RESIDU+ XXX1;
  1491.                 'DETR' XXX1;'DETR' RESIDU;
  1492.                 RESIDU = XXX2;
  1493.                 INIT = VRAI;
  1494.             'FINS';
  1495.           'FINS';
  1496.         'FINS';
  1497.       'FINS';
  1498.     'FINS';
  1499.     WTAB.'ZDEP1'=zdep1;
  1500. *
  1501. *   Initialisation pour l'acceleration de convergence (on peut mettre 
  1502. *   n'importe quoi c'est pour ne pas faire de tests dans la boucle)
  1503.     ACFP1  = 'COPIER' ZFEXT2 *0. ;
  1504.     ACFP2  = ACFP1 ;
  1505.     ACFP3  = ACFP1 ;
  1506.     ACFEP1 = ACFP1 ;
  1507.     ACFEP2 = ACFP1 ;
  1508.     FCORF  = 'COPIER' FREAP ;
  1509.     CORREC = 0;
  1510. *
  1511.     zdept = (zdept 'EXCO' 'LX' 'NOID' 'LX') + (ZDEPT * 0.);
  1512. *----------------------------------------------------------------------*
  1513. *       Initialisation des messages pour le pas de temps courant       *
  1514. *----------------------------------------------------------------------*
  1515.     'SI' IKLFFF ;
  1516.       'MESSAGE'
  1517. 'Iter'*13   'Nplas'*26    'Critere'*39        'Deps.max'*52       'Eps.max'*65      'Crit.flex'*78
  1518. ;
  1519.     'SINON';
  1520.       'MESSAGE'
  1521. 'Iter'*13   'Nplas'*26     'Fresidu'*39        'Deps.max'*52       'Eps.max'*65      'Mresidu'*78
  1522. ;
  1523.     'FINS';
  1524.   'FINSI' ;
  1525. *FOR_MECA
  1526. *
  1527.   NONCONV = FAUX; ZICONV = VRAI; PASTEST = FAUX; PASREINI = VRAI;
  1528.   ITACC = 0 ; RECA_N = 0 ; IT_RECA = 0 ;
  1529.   RED_URG = 0 ; DPSMAXP = 0; DEPSTDM = 0. ; MMCMAX = 0 ; 
  1530.   IT = 0 ; MMC = 0 ; XCONV = 0. ; DPSMAX = xpetit ; EPSM = 0. ;
  1531.   IPREM = VRAI ; IRATE = FAUX; URG = FAUX ;
  1532.   IKT = FAUX ; zprecnc=1e-5 ; TABCONV = 'TABL';
  1533.   RE_CNTRL = 'VIDE' 'RIGIDITE' / 'RIGIDITE'; DE_CNTRL = FAUX;
  1534. *######################################################################*
  1535. *----------------------------------------------------------------------*
  1536. *           Boucle de convergence ETIQ (Configuration GEOM1)           *
  1537. *----------------------------------------------------------------------*
  1538. * preparation pour la gestion du pas avant etiq
  1539.   depstp zdeptp zsigfp fcorfp correcp = depst zdept zsigf fcorf 0;
  1540. * initialisation sous decoupage pour epsi
  1541.   n_sous_dec = 1 ;
  1542. *
  1543.   'REPETER' ETIQ ;
  1544. *
  1545. *   IT est le compteur de ETIQ
  1546.     IT = IT + 1 ;
  1547. *
  1548.     'SI' WTAB.'FOR_MECA' ;
  1549. *
  1550. *     A tester mais ne fonctionne pas en l'etat
  1551. *           nsoincr = n_sous_dec;
  1552. *
  1553.       depst zdept zsigf fcorf correc = depstp zdeptp zsigfp fcorfp correcp;
  1554. *     Cas ou l'iteration precedente qui ne s'est pas faite entierement
  1555. *     -> partir du bon etat (mis a jour si le sous increment est acceptable) 
  1556.       'SI' IRATE;
  1557.         HPP_EPS = FAUX ;
  1558.         URG     = VRAI ;
  1559.         ITACC   = 4 ;
  1560.         'SI' AUTAUG ; IRAUG = VRAI ; 'FINSI' ;
  1561.       'FINSI' ;
  1562.       IRATE = VRAI;
  1563. *
  1564. *     Champ materiau final
  1565.       ZMAT2F = ZMAT2 ;
  1566. *              
  1567.       'SI' autaug;
  1568.         resmul = 2. * resmul;
  1569.         'SI' ((dpsmax <   1d-3) 'ET' (resmul > 1D-20) 'ET' ('NON' PASTEST)); 
  1570.           resmul = resmul / (dpsmax + xpetit) * 1d-3; 
  1571.         'FINSI';
  1572.         'SI' (resmul > 1.0); resmul = 1.0; 'SINON'; itacc = 4; 'FINSI';
  1573.       'FINSI';            
  1574. *
  1575.       zdep1 = (zdept 'EXCO' 'LX' 'NOID' 'LX') + (ZDEPT * 0.);
  1576. *
  1577. *     Pas de temps courant (WTAB.'DT') eventuellement sous-decoupe via COEPI
  1578.       'SI' (knoconv '>' 1) ;
  1579.         ZDT = 0. ;
  1580.       'SINON' ;
  1581.         ZDT = WTAB.'DT' '*' COEPI ;
  1582.       'FINS' ;
  1583. *
  1584. *     ITACC doit etre =< 0 pour qu'on accelere
  1585.       ITACC = ITACC - 1;
  1586. *
  1587.       PASTEST = FAUX;
  1588.       RECA_K  = FAUX;
  1589. *
  1590. *     Pour conserver les criteres de convergence
  1591.       tabconv. it = 1;
  1592. *----------------------------------------------------------------------*
  1593. *                Recalcul de la rigidite si necessaire                 *
  1594. *----------------------------------------------------------------------*
  1595.       'SI' IGRD ;
  1596.         'FORM' GEOM1;
  1597. *
  1598.         'SI' ((IT > 1) 'ET' (URG 'OU' ('MULT' (IT - IT_RECA) ITRCLC))) ;
  1599. *
  1600.           URG     = FAUX ;
  1601.           IT_RECA = IT;
  1602.           RECA_K  = VRAI;
  1603.           RECA_N  = RECA_N + 1;
  1604.           HPP_EPS = FAUX;
  1605.           'SI' (RECA_N > 20) ;
  1606.             nonconv = vrai;
  1607.           'FINS';
  1608.           PASREINI = FAUX;
  1609. *
  1610.           GEOR ZMATTEMP = 'FORM' ZDEPT ZMODLI ZMAT22 ;
  1611. *         A cause de FORM, ZMATTEMP n'est plus parallele...
  1612.           ZMATTEMP = 'REDU' ZMATTEMP ZMODL;
  1613. *
  1614.           ZDEPTR = ZDEP0 '+' ZDEPT ;
  1615. *
  1616. *         ----------------------------------------------------------------
  1617. *         Recalcul de la rigidite a la fin du pas 
  1618.           txt_k ='CHAI' ' Recalcul de K (= K^el';
  1619.           'SOUC' 0;
  1620.           'SI' ('EGA' LAG_TOT 1) ;
  1621.             'FORM' GEOREF0;
  1622.             HOOKRH  = 'HOOK' ZMODL ('REDU' ZMAT2I ZMODL) ;
  1623.             'FORM' GEOR;
  1624.             HOOKRH2 = 'CAPI' HYPDEF HOOKRH ZMODL (-1. * ZDEPTR);
  1625.             RITC = 'RIGI' HOOKRH2 ZMODL zmattemp 'NOER';
  1626.             'DETR' HOOKRH; 'DETR' HOOKRH2;
  1627.           'SINON';
  1628.             RITC = 'RIGI' ZMODL ZMATTEMP 'NOER' ;
  1629.           'FINSI';
  1630.           'FORM' GEOM1;
  1631.           ZMATTEMP = 0 ;
  1632.           'SI' (SOUCI) ;
  1633.             'MESS' 'rigi rate. on reessaye avec la configuration geom1';
  1634.             RITC = 'RIGI' ZMODL ZMAT2F ;
  1635.           'FINSI';
  1636. *
  1637. *         RITC peut contenir des CL
  1638.           ZRI = 'EXTR' RITC 'RIGI' 'NOMU' ;
  1639.           ZCL = 'EXTR' RITC 'RIGI' 'MULT' ;
  1640.           'SI' ('NEG' ('DIME' ZCL) 0) ;
  1641.             ZCLIM = ZCLIM0 'ET' ZCL ;
  1642.           'SINON' ;
  1643.             ZCLIM = ZCLIM0 ;
  1644.           'FINSI' ;
  1645. *
  1646. *         ----------------------------------------------------------------
  1647. *         Option 'AUGMENTATION_AUTOMATIQUE' : actualiser la rig augm
  1648.           'SI' (AUTAUG 'ET' ('NON' IRAUGLU));
  1649.             'FORM' GEOR ;
  1650.             RIG_AUG = 'MASSE' ZMODL ZMAT2F ;
  1651.             RH      = RITC;
  1652.            'FORM' GEOM1;
  1653. ***         mess 'actualisation rig_aug';
  1654.           'FINSI';
  1655. *
  1656. *         ----------------------------------------------------------------
  1657. *         Ajout de relations issues de la procedure DEFO_IMP
  1658.           ZFCONT   = 0 ;
  1659.           IACTURES = FAUX ;
  1660.           CDEPSLX  = 'ENLE' (ZDEP0 '+' ZDEPT) 'LX' ;
  1661.  
  1662. **        DE_CNTRL = FAUX;
  1663.           'SI' (AUTAUG 'ET' DE_CNTRL);
  1664.             IACTURES = VRAI ;
  1665.             ZFCONT = RE_CNTRL '*' CDEPSLX 'EXCO' 'FLX' 'NOID' 'FLX'  ;
  1666.             ZCLIM  = RE_CNTRL 'ET' ZCLIM ;
  1667.             txt_k  = 'CHAI' txt_k ' + K^ctr';
  1668.           'FINSI';
  1669. *
  1670. *         ----------------------------------------------------------------
  1671. *         Recalcul des conditions de contact-frottement ds la conf finale
  1672.           'SI' WTAB.'CONTACT';
  1673.             'FORM' GEOR ;
  1674. *
  1675.             ZFLX = 'EXCO' (CRR '*' ZDEPT) 'FLX' 'NOID' 'FLX' ;
  1676. *     
  1677.             MODCON  = wtab.'MODCONTA';
  1678.             'SI' ('EXIS' WTAB 'MATCONTA');
  1679.               CJEU CRR = 'RFCO' MODCON  WTAB.'CONV'  WTAB.'MATCONTA';
  1680.             'SINON' ;
  1681.               CJEU CRR = 'RFCO' MODCON  WTAB.'CONV'                ;
  1682.             'FINS' ; 
  1683. *
  1684. *           Separer ce qui concerne l'adherence (FADH) et les jeux (FLX)
  1685.             'SI' ('NEG' CJEU 0) ;
  1686.               FADHE = 'EXCO' CJEU 'FADH' 'NOID' 'FADH' ;
  1687.               CDAP  = 'EXCO' CJEU 'FLX'  'NOID' 'FLX' ;
  1688. *
  1689. *             Ajout du glissement deja parcouru sur les conditions de frot.
  1690.               'SI' WTAB.'CAFROTTE';
  1691.               'SI' ('EGA' ('TYPE' MFROT) 'MAILLAGE');
  1692.                 CGLI = 'REDU' ZFLX MFROT ;
  1693.                 CDAP = CDAP ET (-1*CGLI) ;
  1694.               'FINSI'; 
  1695.               'FINSI'; 
  1696.             'FINSI'; 
  1697. *
  1698.             'SI' WTAB.'MODAL' ;
  1699.               CRR  = 'PJBA' CRR  ZMODL ZMAT2F ;
  1700.               MCRR = 'EXTR' CRR 'MAIL' 'MULT';
  1701.               MCDAP = 'EXTR' CDAP 'MAIL' ;
  1702.  
  1703.               'REPETER' BCDA ('NBNO' MCDAP) ;
  1704.                  PBCDA = MCDAP 'POINT' &BCDA ;
  1705.                  PCRR  ='POINT' MCRR &BCDA   ;
  1706.                  CHCR  ='MANU' 'CHPO' PCRR 1 'FLX' ('EXTR' CDAP 'FLX' PBCDA)
  1707.                         'NATURE' 'DISCRETE' ;
  1708.                 'SI' ('EGA' 1 &BCDA) ;
  1709.                    CCDA = CHCR ;
  1710.                 'SINON' ;
  1711.                    CCDA = CCDA 'ET' CHCR ;
  1712.                 'FINS' ;
  1713.               'FIN' BCDA ;
  1714.  
  1715.               CDAP = CCDA ;
  1716.             'FINS' ;
  1717. *
  1718.             'SI' ('NEG' CRR 0) ;
  1719.               IACTURES = VRAI ;
  1720.               CHPZ     = CRR '*' CDEPSLX 'EXCO' 'FLX' 'NOID' 'FLX'  ;
  1721.               'SI' ('NEG' ZFCONT 0) ;
  1722.               ZFCONT = ZFCONT '+' CHPZ ;
  1723.               'SINON' ;
  1724.               ZFCONT = CHPZ ;
  1725.               'FINSI' ;
  1726.               ZCLIM  = CRR 'ET' ZCLIM ;
  1727.             'FINS';
  1728. *
  1729.             'SI' ('NEG' CJEU 0) ;
  1730.               IACTURES = VRAI ;
  1731.               ZFCONT = ZFCONT '+' CDAP ;
  1732.             'FINS';
  1733.           'FINS';
  1734. *
  1735. *         Mise a jour des jeux dans le residu (partie force inchangee)
  1736.           ZFCONSTA = ZFEXT2 ;
  1737.           ZFLX1    = ZFLXB ;
  1738.           'SI' IACTURES ;
  1739.             ZFCONSTA = ZFCONSTA '+' ZFCONT ;
  1740.             RFORCE = 'ENLE' RESIDU 'FLX' ;
  1741.             FLXPRE = 'EXCO' (ZCLIM '*' (ZDEP0 '+' ZDEPT)) 'FLX' 'NOID' 'FLX';
  1742.             RJEUX  = ZFCONT '+' ZFLXB '-' FLXPRE ;
  1743.             RESIDU = RFORCE 'ET' RJEUX ;
  1744.             ZFLX1  = ZFLX1 'ET' ZFCONT ;
  1745.           'FINS';
  1746. *
  1747.           ZRAID = ZCLIM 'ET' ZRI;
  1748. *         Mise a jour de FREAP (faut-il faire cette actualisation?)
  1749.           FREAP = (ZDEP0 'EXCO' 'LX' 'LX')  * ZCLIM;
  1750. *
  1751. *         ----------------------------------------------------------------
  1752. *         Ajout de la partie dynamique
  1753.           'SI' IDYN;
  1754. *           bp : en toute rigueur, il faudrait aussi recalculer la MASSE ...
  1755. *                et ajouter l'amortissement le cas échéant ...
  1756.             ZRAID = ZRAID 'ET' ZRMAS ;
  1757.           'FINS';
  1758. *
  1759. *         ----------------------------------------------------------------
  1760. *         Recalcul de la raideur geometrique (option K_SIGMA)
  1761.           'SI' IKSIA ;
  1762.             'FORM' GEOR;
  1763.             ZSIGKS = 'CAPI' HYPDEF ZSIGF ZMODL (-1 * ZDEPTR) ;
  1764.             KSIGTC = 'KSIGMA' ZMODL ZSIGKS ZMAT2;
  1765.             ZRAID  = ZRAID 'ET' ksigtc ;
  1766.             txt_k  = 'CHAI' txt_k ' + K^sig';
  1767.             'FORM' GEOM1;
  1768.           'FINS';
  1769. *
  1770. *         ----------------------------------------------------------------
  1771. *         Prise en compte d'eventuelles 'RIGIDITE_CONSTANTE'
  1772.           'SI' IRCON;
  1773.             ZRAID = ZRAID 'ET' RIG_CONS;
  1774.             txt_k = 'CHAI' txt_k ' + K^cst';
  1775.           'FINS';
  1776. *
  1777. *         ----------------------------------------------------------------
  1778. *         Traitement en de la rigidite augmentee
  1779.           'SI' IRAUG;
  1780.             'SI' AUTAUG ;
  1781.               'SI'  iprem; zdep1d = ZDEP0 'ENLE' 'LX'; 'FINSI';
  1782.               xkx = (xtmx zdep1d rh)  + xpetit;
  1783.               xmx = (xtmx zdep1d rig_aug) + xpetit;
  1784.               'SI' (('VERI' xkx) 'ET' ('VERI' xmx) 'NON');
  1785.                 zu1l = ZDEP0 'ENLE' 'LX';
  1786.                 xkx = (xtmx zu1l rh)  + xpetit;
  1787.                 xmx = (xtmx zu1l rig_aug) + xpetit;
  1788.               'FINSI';
  1789.               augauto = xkx / xmx; augk = 1.d-3;
  1790.               augauto = abs augauto;
  1791.               augk    = abs augk;
  1792.               augmult = augmult * 1.02;
  1793.               'SI' (augmult > 1D4 ); augmult=1D4 ; 'FINSI';
  1794.               augauto = augauto * augmult ;
  1795.               augk = augk * augmult;
  1796.  
  1797.               ZRAID = ZRAID 'ET' (RIG_AUG * augauto) 'ET' (RH * augk) ;
  1798.  
  1799.               'MESS' 'multiplicateur d augmentation masse' augauto ' raideur' augk ; 
  1800.             'SINON';
  1801.               ZRAID = ZRAID 'ET' RIG_AUG ;
  1802.             'FINSI';
  1803.             txt_k = 'CHAI' txt_k ' + K^aug';
  1804. ***         ZNACCE = FAUX;
  1805.           'SINON';
  1806.             augmult = augmult * 0.55 ;
  1807.             si (augmult < 1d-3); augmult = 1d-3; finsi;
  1808.           'FINSI';
  1809.           'FORM' GEOM1;
  1810. *
  1811.           XXX1   = 'EXTR' ZCLIM 'MAIL' 'UNIL';
  1812.           IMPO12 = ('NBEL' XXX1) '>' 0 ;
  1813.           'SI' (IPILOT 'ET' IMPO12) ;
  1814.             DUUNIL = 'REDU' DFEXT0L XXX1 ;
  1815.             XXX2   = 'DIFF' ('EXTR' ZCLIM 'MAIL') XXX1 ;
  1816.             DUIMPO = 'REDU' DFEXT0L XXX1 ;
  1817.           'FINS';
  1818. *  
  1819. *         ----------------------------------------------------------------
  1820. *         Procedure UTILISATEUR
  1821.           'SI' WTAB.'PROCEDURE_CHARMECA';
  1822.             PRECED.'ADDI_MATRICE' = vrai;
  1823.             TFP22= CHARMECA PRECED TI ;
  1824.             PRECED.'ADDI_MATRICE' = faux;
  1825.             'SI' (EXIS TFP22 'ADDI_MATRICE');
  1826.               zraid = zraid 'ET' TFP22.'ADDI_MATRICE';
  1827.               txt_k = 'CHAI' txt_k ' + K^cent';
  1828.             'FINS';
  1829.           'FINS';
  1830. *
  1831.           'MESS' ('CHAI' txt_k ' ) dans config deformee ') ;
  1832. *
  1833. *         on impose le recalcul de K a la prochaine iteration si it=2 (pq ?)
  1834.           'SI' (IT 'EGA' 2) ;
  1835.             ITACC = 4 ;
  1836.           'FINS' ;
  1837.         'FINS';
  1838.       'FINS';
  1839. *----------------------------------------------------------------------*
  1840. *            Evaluation de la matrice tangente si demandee             *
  1841. *----------------------------------------------------------------------*
  1842.       'SI' IKTAN ;
  1843.         'SI' ('NON' IFEFP) ;
  1844.           'SI' (('MULT' IT WTAB.'NITER_KTANGENT') 'ET' (IT > 1)) ;
  1845. *
  1846.             IKT = VRAI ;
  1847.             'SI' (IGRD 'ET' ('NON' HPP_EPS));
  1848.               'FORM' GEOR ;
  1849.             'FINS' ;
  1850.  
  1851.             'SI' ('NON' ISSTE) ;
  1852.               'SI' IPERT ;
  1853.                 'SI' PARTLOCA ;
  1854.                   cha1  = 'ASSIS' 'TOUS' 'REDU' MODRELOC Z1COMP ;
  1855.                   cha2  = 'ASSIS' 'TOUS' 'REDU' MODRELOC Z2COMP ;
  1856.                   zktap = 'ASSIS' 'TOUS' 'KTAN' 'PERT' MODRELOC cha1 cha2
  1857.                            'C1' WTAB.'K_TANG_PERT_C1' 
  1858.                            'C2' WTAB.'K_TANG_PERT_C2'
  1859.                            ZKTASYM ;
  1860.                   ZRIKTA = 'ET' zktap ;
  1861.                 'SINON';
  1862.                   ZRIKTA = 'KTAN' 'PERT' ZMODL Z1COMP Z2COMP
  1863.                            'C1' WTAB.'K_TANG_PERT_C1'
  1864.                            'C2' WTAB.'K_TANG_PERT_C2'
  1865.                            ZKTASYM ;
  1866.                 'FINS';
  1867.               'SINON' ;
  1868.                 'SI' ('NON' (IVISCO 'OU' IVIDOM)) ;
  1869.                   DTTAN = 0. ;
  1870.                 'SINON' ;
  1871.                   DTTAN = ZDT ;
  1872.                 'FINS' ;
  1873.                 ZRIKTA = 'KTAN' ZMODL ZSIGF ZVARF ZMAT2F
  1874.                                 'PREC' ZPREK 'DT  ' DTTAN ZKTASYM ;
  1875.               'FINS' ;
  1876.             'FINS' ;
  1877.  
  1878.             'SI' IKSIA ;
  1879.               ZKSIG  = 'KSIGMA' ZMODL ZSIGF ZMAT2F ;
  1880.               ZRIKTA = ZRIKTA 'ET' ZKSIG ;
  1881.             'FINS' ;
  1882.  
  1883.             ZRAID = ZRIKTA 'ET' ZCLIM ;
  1884.  
  1885.             'SI' IGRD ; 'FORM' GEOM1 ; 'FINS' ;
  1886.  
  1887.           'FINS' ;
  1888.         'FINS' ;
  1889.       'FINS' ;
  1890.     'FINSI' ;
  1891. *FOR_MECA
  1892. *
  1893.     'SI' WTAB.'NVSTNL'                          ;
  1894.       UN = 'EXCO' ZVIFL movifl                  ;
  1895.       ch_un = 'MANU' 'CHML' WTAB.'MOD_NSL' 'ADVE' un 'TYPE' 'CARACTERISTIQUES';
  1896.       nugrad = 'ADVE' WTAB.'MOD_NSL' ch_un      ;
  1897.       'SI' ('EXIS' PRECED 'OPTION_BDF2') ;
  1898.         ZRNSL = WTAB.'RNSL' 'ET' nugrad 'ET' (WTAB.'NSMA' '*' 1.5)  ;
  1899.       'SINON'                                   ;
  1900.         ZRNSL = WTAB.'RNSL' 'ET' nugrad 'ET' WTAB.'NSMA' ;
  1901.       'FINSI' ;
  1902.       ftnsl  = ZFCONSTA 'ET' fmass              ;
  1903.  
  1904.       'SI' LOGADV                               ;
  1905.         ch_ad = 'MANU' 'CHML' WTAB.'MOD_NSL' 'ADVE' UADV 'TYPE' 'CARACTERISTIQUES';
  1906.         NUADV = 'ADVE' WTAB.'MOD_NSL' ch_ad     ;
  1907.         ZRNSL = ZRNSL 'ET' NUADV                ;
  1908.         ftnsl = ftnsl  '-' fadv                 ;
  1909.       'FINSI'                                   ;
  1910.  
  1911.       resnsl = ftnsl '-' (ZRNSL * ZVIFL)        ;
  1912.  
  1913.     'FINSI'                                       ;
  1914. *----------------------------------------------------------------------*
  1915. *                Acceleration de convergence effective                 *
  1916. *----------------------------------------------------------------------*
  1917.     'SI' WTAB.'FOR_MECA' ;
  1918.       CORRECA = CORREC ;
  1919.       CORREC = 0;
  1920.       ACFEP0 = RESIDU - FCORF;
  1921.       ACFEP0 = ACFEP0 'ENLE' 'FLX';
  1922.       ACFEPT = ACFEP0;
  1923.       ACFP0  = ACFEP0 * ZJAC ;
  1924.       ACFEP0 = ACFEP0 - CORRECA;
  1925.       ACCEL = FAUX;
  1926.       'SI' ((ITACC '&lt;EG' 0) 'ET' (IT '>' 3) 'ET' ZNACCE);
  1927.         ITACC = 2;
  1928.         ACCEL = VRAI;
  1929.         CORREC = 'ACT3' ACFEP2 ACFEP1 ACFEP0
  1930.                         ACFP3  ACFP2  ACFP1  ACFP0 ;
  1931. *
  1932. *       verif que l'acceleration ne renvoie pas en arriere
  1933.         'SI' WTAB.'STABILITE' ;
  1934.           acc_ref = xty acfept zdep1d MNDUAL MNPRIM;
  1935.           acc_ref = acc_ref + xpetit;
  1936.  
  1937.           acc_dir = xty (acfept - correc) zdep1d MNDUAL MNPRIM;
  1938.           acc_rap = acc_dir/acc_ref;
  1939.           acc_lim = 256.; 'SI' (acc_ref < 0.); acc_lim = 1.; 'FINSI';
  1940.           'SI' (acc_rap  '>' acc_lim) ;
  1941. *       en cas d'acceleration trop grande, on la limite
  1942.             'MESS' 'Limitation acceleration ' ' 'acc_rap;  CORREC = correc * (acc_lim/acc_rap);
  1943.           'FINSI';
  1944.           'SI' (acc_rap '&lt;EG' 0.) ;
  1945. *       en cas d'acceleration retrograde, on n'accelere pas
  1946.             'MESS' 'Annulation acceleration: retrograde' ' ' acc_rap;
  1947.             correc =  0.;
  1948.           'FINSI';
  1949.         'FINSI';
  1950. *
  1951.         RESIDU = RESIDU '-' CORREC;
  1952.       'FINS';
  1953. *
  1954.       ACFP3  = ACFP2  ;
  1955.       ACFP2  = ACFP1  ;
  1956.       ACFP1  = ACFP0  ;
  1957.       ACFEP2 = ACFEP1 ;
  1958.       ACFEP1 = ACFEP0 ;
  1959. *
  1960.       'SI' (resmul < 0.99);
  1961.         'MESS' 'Reduction du chargement. Coefficient: ' resmul;
  1962.         'SI' (resmul < 1.d-40); 'ERREUR' 996; 'FINSI';
  1963.         residu = ((residu -fcorf) * resmul) +  fcorf;
  1964.         ACFP0 = (RESIDU  - FCORF) ;
  1965.       'FINSI';
  1966. *----------------------------------------------------------------------*
  1967. *                 Systeme a resoudre a l'iteration (i)                 *
  1968. *                                                                      *
  1969. *           | K  At | |  du   |   | Fsur + Fvol - Bsigma^(i) |         *
  1970. *           |       | |       | = |                          |         *
  1971. *           | A  0  | | LX(i) |   | d1 - A*u^(i)             |         *
  1972. *----------------------------------------------------------------------*
  1973.       'SOUCI' 0;
  1974.       FEXCI = 'VIDE' 'CHPOINT' / 'DISCRET' ;
  1975.       'SI' WTAB.'ADHERENCE' ;
  1976.         'SI' ('EXIS' WTAB 'MAIL_BLOM') ;
  1977.           FADRE = 'REDU' FADHE WTAB.'MAIL_BLOM' ;
  1978.         'SINON' ;
  1979.           FADRE = FADHE ;
  1980.         'FINSI' ;
  1981.         FEXCI = FEXCI 'ET' FADRE ;
  1982.       'FINSI' ;
  1983. *
  1984. *     Force limite de frottement
  1985.       'SI' (WTAB.'CAFROTTE' 'ET' IMPO12);
  1986.         'SI' WTAB.'FROCABL' ;
  1987.           zsigfT = 'REDU' zsigf ZMODLI;
  1988.           FFROT = 'EXCC' ZCLIM ZDEPT MODCON WTAB.'MATCONTA' ZSIGFT;
  1989.           RTRSF = 'VIDE' 'RIGIDITE' / 'RIGIDITE';
  1990.           RFNS  = 'VIDE' 'RIGIDITE' / 'RIGIDITE';
  1991.           FEXCI = FEXCI '+' FFROT ;
  1992.         'FINSI';
  1993.         'SI' WTAB.'FROCOUL' ;
  1994.           excfconv = excfconv + 1 ;  
  1995.           RTRSF*'RIGIDITE' RFNS*'RIGIDITE'  FFROT = 'EXCF' ZCLIM ZDEPT MODCON WTAB.'MATCONTA' excfconv;
  1996.           FEXCI = FEXCI '+' FFROT ;
  1997.           MFROT = 'EXTR' FFROT 'MAIL' ;
  1998.         'FINSI';
  1999.       'FINS';
  2000. *
  2001.       'SI' ( ('NEG' ZRIBLO_M 'INCONNU') 'ET' (IPREM 'OU' RECA_K) ) ;
  2002.         'SI' ((WTAB.'CAFROTTE' 'OU' WTAB.'ADHERENCE') 'ET' IMPO12);
  2003.           ZDEP1 BID BID BID BID = 'RESO' ZRAID 'SOUC' RESIDU 'INIB'
  2004.                                    ZRIBLO_M ZLISEA_M
  2005.                                    FEXCI RFNS RTRSF;
  2006.         'SINON';
  2007.           ZDEP1 BID BID BID     = 'RESO' ZRAID 'SOUC' RESIDU 'INIB'
  2008.                                    ZRIBLO_M ZLISEA_M ;
  2009.         'FINS';
  2010.       'SINON';
  2011.         'SI' ((WTAB.'CAFROTTE' 'OU' WTAB.'ADHERENCE') 'ET' IMPO12);
  2012.           ZDEP1 = 'RESO' ZRAID 'SOUC' RESIDU FEXCI       RFNS RTRSF;
  2013.         'SINON';
  2014.           ZDEP1 = 'RESO' ZRAID 'SOUC' RESIDU 'NOID';
  2015.           'SI' (WTAB.'MAN' 'ET' IPREM);
  2016.             ORDRE      = WTAB.'ORDRE' ;
  2017.             ZDEP2 IOUT = CORMAN ZRAIDINI ZMODL ZMAT2F ORDRE
  2018.                           ZDEP0 ZSIG0      RESIDNOR  WTAB ;
  2019.             'SI' (IOUT 'EGA' 1) ;
  2020.               ZDEP1=ZDEP2;
  2021.             'FINS';
  2022.           'FINSI' ;
  2023.         'FINS';
  2024.       'FINS';
  2025. *
  2026. *     verif stabilite matrice elastique
  2027.       'SI' (wtab.'STABILITE' 'ET' ('NON' HPP_EPS) 'ET' AUTAUG);
  2028.         RESIDD = RESIDU - (ZRAID * (ZDEP1 'EXCO' 'LX' 'LX'));
  2029.         acc_ref= xty (RESIDD       ) zdep1 MNDUAL MNPRIM;
  2030. ** diag refactorise la matrise. On evite donc provisoirement de l'appeler
  2031. **      nbng = 'DIAG' ZRAID;
  2032.         nbng = 0;
  2033.         'SI' ((acc_ref <   0.) 'OU' (nbng > 0));
  2034.           'MESS' 'Raideur negative' ' ' acc_ref ;
  2035. **        resmul = resmul * 0.25;
  2036. **          augmult = augmult * 1.4;
  2037.           'SI' FAUX;
  2038. *  controler la deformation max
  2039.             DEPST = 'EPSI' 'LINE' ZMODL (ZDEPT + ZDEP1) ZMAT2F;
  2040.             DPSMAX = 'MAXI' 'ABS' DEPST 'AVEC' MLDEFOR;
  2041.             RE_CNTRL2*'RIGIDITE'   = DEFO_IMP DEPST  (zdept + zdep1)  ZRAID ;
  2042.             DE_CNTRL = FAUX;
  2043.             ZFCONT = RE_CNTRL2 '*' (        ZDEP1) 'EXCO'  'FLX' 'NOID' 'FLX'  ;
  2044.             cmul = MINI  1. ( 1d-3 / DPSMAX);
  2045.             ZFEXT2   = ZFEXT2   + (ZFCONT * cmul);
  2046.             ZFLXB    = 'EXCO' ZFEXT2 'FLX' 'NOID' 'FLX' ;
  2047.             ZFCONSTA = ZFEXT2 ;
  2048.             RE_CNTRL = RE_CNTRL 'ET' RE_CNTRL2;
  2049.           'FINSI';
  2050. **        'ITERER' ETIQ;
  2051. ** essai en inversant zdep1 pour avoir la positivite du travail
  2052.           zdep1 = ((zdep1 'ENLE' 'LX') * -1.) 'ET' (zdep1 'EXCO' 'LX' 'LX');
  2053.         'FINSI';
  2054.       'FINSI';
  2055. *
  2056.       'SI' ('SOUCI') ; 
  2057.         'SI' AUTAUG ; 
  2058.           'MESS' 'Annulation du sous-increment: souci dans la resolution';
  2059.           resmul = resmul * 0.25;
  2060.           augmult = augmult * 1.4;
  2061.           'ITERER' ETIQ;
  2062.         'SINON' ;
  2063.           'MESS' 'Souci':  (vale 'SOUCI');
  2064.           pastest = VRAI;
  2065. **        IKSIA = FAUX;
  2066.           HPP_EPS = FAUX;
  2067.           URG     = VRAI;
  2068.           ITACC   = 4;
  2069. ***       'ITERER' ETIQ;
  2070.         'FINSI';
  2071.       'FINSI';
  2072. *
  2073.       'SI' LOGPIL ; 
  2074. *       Calcul de D_eta par appel a la procedure PILOINDI
  2075.         D_ETA = PILOINDI PRECED ZDEP0 ZDEPT ZDEP1 ZDEPII DTAU ;       
  2076. *       Mise à jour  
  2077.         ZDEP1 = ZDEP1 '+' (D_ETA '*' ZDEPII) ;
  2078.         ETA = ETA '+' D_ETA ;
  2079.       'FINSI' ;
  2080. *
  2081. *  verif sens ----------------------------------------------------------
  2082. *
  2083.       'SI' (resmul < 1d-2) ; souci 0; 'FINSI';
  2084.       zdep1s = zdep1 'ENLE' 'LX';
  2085.       pasunil = 'SOUCI' ;
  2086. *
  2087.       ZRAID_T = 'EXTR' ZRAID 'CONT';
  2088.       'SI' ('EXIS' ZRAID_T 'NITER');
  2089. *       Conserver les conditions unilaterales pour VITEUNIL
  2090.         'SI' IDYN ; WTAB.'ZRAIDV' = ZRAID; 'FINS';
  2091.         ZRIBLO_M = ZRAID_T.7;
  2092.         ZLISEA_M = ZRAID_T.6;
  2093.         WTAB.'MAIL_BLOM' = 'EXTR' ZRIBLO_M 'MAIL' 'MULT' ;
  2094.       'FINS';
  2095.     'FINSI' ;
  2096. *FOR_MECA
  2097.  
  2098.     'SI' WTAB.'NVSTNL'             ;
  2099.       ZVIF1  = 'KRES' ZRNSL Resnsl ;  
  2100.       ZVIFL = ZVIFL '+' ZVIF1      ;
  2101.     'FINSI'                        ;
  2102. *----------------------------------------------------------------------*
  2103. *            Increment de deplacements et mult. de Lagrange            *
  2104. *                                                                      *
  2105. * ZDEP1 est un increment d'increment de deplacements, note du, et les  *
  2106. * multiplicateurs de Lagrange sont resolus "en total".                 *
  2107. * La solution de la boucle ETIQ a l'iteration i a fourni :             *
  2108. *                ZDEP1  = [ du  ; LX(i) ]                              *
  2109. * Le cumul est realise de maniere a avoir :                            *
  2110. *                ZDEPT  = [ Du^(i) ; LX(i) ]                           *
  2111. * avec Du^(i) = Du^(i-1) + du                                          *
  2112. * Le champ de deplacements total s'ecrit :                             *
  2113. *                ZDETOT = ZDEP0 '+' ZDEPT                              *
  2114. *----------------------------------------------------------------------*
  2115.     'SI' WTAB.'FOR_MECA' ;
  2116.       'SI' IPREM ;
  2117.         ZDEPT = 'COPIER' ZDEP1 ;
  2118.         ZDELA = 'COPIER' ZDEPT ;
  2119.       'SINON';
  2120.         XXX1  = ZDEPT 'ENLEVER' 'LX' ;
  2121.         ZDEPT = XXX1 '+' ZDEP1 ;
  2122.         'DETR' XXX1 ;
  2123.       'FINS' ;
  2124. *
  2125.       'SI' WTAB.'MODAL' ;
  2126.         'SI' ('EXIS' WTAB 'MODCONTA') ;
  2127.           'SI' IPREM ;
  2128.             mamoco1 = 'EXTR' (ZDEPT 'ENLEVER' 'LX') 'MAIL' ;
  2129.           'SINON' ;
  2130. *           mettre les point materiels dans zdept
  2131.             zdeptu1 =  'REDU' zdept   mamoco1      ;
  2132.             ch_dco  =  'RECO' zdeptu1 ZMODL ZMAT22 ;
  2133.             ZDEPT   = ('EXCO' zdept 'LX' 'LX') 'ET' zdeptu1 'ET' ch_dco ;
  2134.           'FINS';
  2135.         'FINS';
  2136.       'FINS' ;
  2137. *
  2138. *     Pilotage automatique
  2139.       'SI' IPILOT;
  2140.         'SI' ACCEL ;
  2141.           RED1 = ZAUTOREDU ;
  2142.  
  2143. *         reduction du critere de pilotage red2 sert d'incateur si
  2144. *         pas accelere tous les 2 pas
  2145.  
  2146. *         MODIFICATION CB215821 : 18/06/2015
  2147.           'SI' (RED2 < (IT '/' 20));
  2148.             ZAUTOREDU = ZAUTOREDU '*' 3.D0 ;
  2149.             RED2 = RED2 '+' 1 ;
  2150.             ITACC= 3;
  2151.           'FINS';
  2152.  
  2153.           RED1 = RED1 '/' ZAUTOREDU ;
  2154.           'SI' (RED1 '<' 0.9) ;  ZAUTORED1 = 10 ; 'FINS' ;
  2155.  
  2156.           'SI' (ZAUTOREDU '>' 1.D0 ) ;
  2157.             'MESS' 'On divise le critere de pilotage par 'ZAUTOREDU;
  2158.           'FINS' ;
  2159.  
  2160.           OO   = WTAB.'AUTOCRIT' '/' ZAUTOREDU ;
  2161.           U1MA = AUTOPILO ZDEPT (COEPI'*'ZDELA) ZMODLI ZMAT22 WTAB;
  2162.           AL1  = OO '/' U1MA ;
  2163. *         al1 est le coefficient de normalisation
  2164.           'SI' (( al1 '>EG' 1.d0) 'ET' (coepi '>EG' 1.d0)) ;
  2165. *         pour eviter d'aller au dela de alpha=1 on ignore le critere
  2166.             AL1 = 1.d0 ;
  2167.           'FINS' ;
  2168.           'SI'((AL1 '>' 1.D0) 'ET' (COEPI '>' 0.D0));
  2169.             'SI' ( AL1 '>' (1.D0 '/' COEPI));
  2170.               AL1 = 1.D0 '/' COEPI;
  2171.             'FINS';
  2172.           'FINS';
  2173. *         normalisation
  2174.           PASTEST = VRAI;
  2175.           XXX1    = ZLX   '*' ( 1.d0 '-' AL1);
  2176.           XXX3    = ZDEPT '*' AL1 ;
  2177.           ZDEPT   = XXX3  '+' XXX1;
  2178.           'DETR' XXX3;
  2179.         'SINON';
  2180. *         le calcul d'al1 est fait pour eviter la convergence
  2181. *         s'il est  externe a l'interval 0.5 --- 2.
  2182.           OO   = WTAB.'AUTOCRIT' '/' ZAUTOREDU ;
  2183.           U1MA = AUTOPILO ZDEPT (COEPI'*'ZDELA) ZMODLI ZMAT22 WTAB;
  2184.           AL1  = OO '/' U1MA ;
  2185.           'SI' (( AL1 '>EG' 1.d0) 'ET' (COEPI '>EG' 1.d0)) ;
  2186.             AL1 = 1.d0 ;
  2187.           'FINS' ;
  2188.         'FINS' ;
  2189. *
  2190.       'SINON' ;
  2191. *       On part dans le decor : redemarrer a 0
  2192.         'SI' ((XCONV > 1E8) 'ET' PASREINI) ;
  2193.           'MESS' 'Reinitialisation du schema';
  2194.           PASREINI=FAUX;
  2195.           'ITERER' ETIQ;
  2196.         'FINS';
  2197.       'FINS';
  2198. *
  2199. *     garder les reactions pour le test de convergence
  2200.       ZDEPLP = ZDEPL                         ;
  2201.       ZDEPL  = ZDEPT 'EXCO' 'LX' 'NOID' 'LX' ;
  2202.       ZGRDUF = 'MOT' 'INCONNU';
  2203. *
  2204. *----------------------------------------------------------------------
  2205. *     le nouveau champ est fixe on va tester l'equilibre(convergence)
  2206. *     et calculer la force motrice pour l'iteration suivante
  2207. *----------------------------------------------------------------------
  2208.       'DETRUIRE' FCORF;
  2209. *
  2210. *     Champ de deplacement total
  2211.       ZDETOT = ZD0SLX '+' ZDEPT ;
  2212. *
  2213. *     XXX1 = [ FCORF ; FCORU ] = [ -Freac_it ; u^imp_it ]
  2214.       XXX1  = ZCLIM '*' ZDETOT;
  2215.       FCORF = 'ENLE' XXX1 'FLX';
  2216.       FCORU = 'EXCO' XXX1 'FLX' 'NOID' 'FLX';
  2217.       'DETR' XXX1 ;
  2218. *----------------------------------------------------------------------*
  2219. *              Increment de deformations totales DEPST                 *
  2220. *----------------------------------------------------------------------*
  2221. *     Update or total lagrangian ---------------------------------------
  2222.       'SI' IFEFP;
  2223.         'SI' IFEFPUL ;
  2224. *         mess ' update lagrangian ZRIKTA';
  2225.           GEOM2 = 'FORM' ZDEPT ;
  2226.           ZRIKTA ZSIGF ZVARF ZDEIF = 'ECFEFP'
  2227.                 ZMODL ZDEI0 ZVAR0 ZDEPT ZMAT2F ZPREK NITMA 1 ;
  2228.         'SINON';
  2229. *         mess ' total lagrangian ZRIKTA';
  2230.           chp_z = ZDEPT + ZDEP0 ;
  2231.           GEOM2 = 'FORM' chp_z ;
  2232.           ZRIKTA ZSIGF ZVARF ZDEIF = 'ECFEFP'
  2233.                 ZMODL ZDEI0 ZVAR0 chp_z ZMAT2F ZPREK NITMA ;
  2234.           chp_z = 1 ;
  2235.         'FINS';
  2236.         FEQU2 = 'BSIG' ZMODL ZSIGF ZMAT2F ;
  2237.         ZRAID = ZRIKTA 'ET' ZCLIM ;
  2238.         'FORM' GEOM1 ;
  2239. *
  2240.         XXX1 = 'EXCO' (WTAB.'MOVA') ZVARF ;
  2241.         EPSM = 'MAXI' XXX1 'AVEC' MLDEFOR ;
  2242.         ACC  = 'ABS' ( XXX1 - ACC0 ) ;
  2243.         MMC = 'MASQUE' ACC 'SUPERIEUR' 1.D-10 'SOMME' ;
  2244.         'SI' (MMC '>' MMCMAX) ; MMCMAX = MMC ; 'FINS' ;
  2245.         DPSMAX  = 'MAXI' ACC 'AVEC'  MLDEFOR ;
  2246.         DEPST   = 'CHAN' 'TYPE' ZDEIF 'DEFORMATIONS' ;
  2247.         ZDEFF   = ZDEF0 '+' DEPST ;
  2248. *
  2249. *     cas standard -----------------------------------------------------
  2250.       'SINON';
  2251.         'SI' IGRD;
  2252.           zdepm  = zdept '*' 0.5;
  2253.           'FORM' GEOM1 ;
  2254.           'SI' WTAB.'ITCAR';
  2255.             GEOM2M ZMATTEMP = 'FORM' ZDEPM ZMODLI ZMAT22 ;
  2256. *           A cause de FORM, ZMATTEMP n'est plus parallele...
  2257.             ZMATTEMP        = 'REDU' ZMATTEMP ZMODL ;
  2258.           'SINON';
  2259.             GEOM2M   = 'FORM' GEOM1 ZDEPM ;
  2260.             ZMATTEMP = ZMAT2F;
  2261.           'FINS';
  2262.           GEOM2 = 'FORM' GEOM1 ZDEPT ;
  2263.           'FORM' GEOM1 ;
  2264.         'FINS';
  2265. *
  2266. *       Pour les hypotheses de deformations non lineaires,
  2267. *       on calcule l'increment de deformation sur le pas avec EPSI 'LINE'
  2268. *       mais en se placant sur une configuration intermediaire (mi pas)
  2269.         iDEFOK = VRAI ; PASOK = vrai;
  2270.         'SI' (('NEG' HYPDEF 'LINEAIRE') 'ET' IGRD);
  2271.          'SOUC' 0;
  2272.          zdepti = zdept * (1.0 / n_sous_dec);
  2273.          depst = 0;
  2274.          'REPE' b_eps n_sous_dec;
  2275.            teps = (&b_eps * (1. / n_sous_dec)) - (1. / n_sous_dec /2.);
  2276.            zdepa  = zdept '*' teps;
  2277.  
  2278.            geom2a = form zdepa;
  2279.           'SI' ('OU' IENDOM IVIDOM ICERAM);
  2280.             HOOKENDO = 'HOOK'        ZMODL ZMAT2F ZVAR0 ;
  2281.              DEPSTat  = 'EPSI' 'LINE' 'NOER' ZMODL ZDEPTi  HOOKENDO ZMAT2F;
  2282.             'DETR' HOOKENDO;
  2283.           'SINON';
  2284.             DEPSTat    = 'EPSI' 'LINE' 'NOER' ZMODL ZDEPti ZMATTEMP ;
  2285.           'FINSI';
  2286.            'SI' ('SOUCI');  'FORM' GEOM1;  IDEFOK = faux; 'QUITTER' b_eps; 'FINSI';
  2287.            FORM geom1;
  2288.             DEPSTa = 'CAPI' HYPDEF DEPSTat ZMODL ZDEPa ;
  2289.            depst = depst + depsta;
  2290.            'DETR' geom2a;               'DETR' depsta; 'DETR' depstat;
  2291.          'FIN' b_eps;
  2292.  
  2293.           PASOK = 'SOUCI' 'OU' ('NON' IDEFOK);
  2294.  
  2295.           DPSMAX = 1001;
  2296.           'SI' ('NON' PASOK);
  2297.             DPSMAX = 'MAXI' 'ABS' DEPST 'AVEC' MLDEFOR;
  2298.           'FINSI';
  2299.           'SI' (PASOK 'OU' (DPSMAX > 1000))  ;
  2300.             MSGZ = 'CHAI' 'Deformation non lineaire non calculable' ;
  2301.             'SI' autaug;  
  2302.               mess ('CHAI' MSGZ ' : annulation sous-increment'); 
  2303.               augmult = augmult * 1.4;
  2304.               resmul = resmul * 0.25;
  2305.               'ITERER' ETIQ;
  2306.             'FINSI';
  2307.             MESS ('CHAI' MSGZ ' : deformation lineaire'); 
  2308.             pastest = vrai;
  2309.             iDEFOK = FAUX ;
  2310. *
  2311.           'SINON' ;
  2312.             'FORM' GEOM1 ;
  2313.             DPSMAX = 'MAXI' 'ABS' DEPST 'AVEC' MLDEFOR;
  2314.             'SI' (DPSMAX '>' 1000);
  2315.               si faux;
  2316.                 RE_CNTRL2*'RIGIDITE'  = DEFO_IMP DEPST (zdept + zdep1)  ZRAID;
  2317.                 DE_CNTRL = FAUX;
  2318.                 ZFCONT = RE_CNTRL2 '*' (      ZDEP1) 'EXCO'  'FLX' 'NOID' 'FLX'  ;
  2319.                 cmul = MINI  1. ( 1d-3 / DPSMAX);
  2320.                 ZFEXT2   = ZFEXT2   + (ZFCONT * cmul);
  2321.                 ZFLXB    = 'EXCO' ZFEXT2 'FLX' 'NOID' 'FLX' ;
  2322.                 ZFCONSTA = ZFEXT2 ;
  2323.                 RE_CNTRL = RE_CNTRL 'ET' RE_CNTRL2;
  2324.               finsi;
  2325.               augmult = augmult * 1.4;
  2326.               'ITERER' ETIQ;
  2327.             'FINS';
  2328.           'FINS';
  2329.  
  2330.         'FINS';
  2331. *
  2332. ** on n'a pas pu calcule DEPST : on utilise les deformations lineaires
  2333.         'SI' PASOK;
  2334.           'SI' ('OU' IENDOM IVIDOM  ICERAM);
  2335.             HOOKENDO = 'HOOK'        ZMODL ZMAT2F ZVAR0         ;
  2336.             DEPST    = 'EPSI' 'LINE' ZMODL ZDEPT HOOKENDO ZMAT2F;
  2337.             'DETR' HOOKENDO;
  2338.           'SINON';
  2339.             DEPST = 'EPSI' 'LINE' ZMODL ZDEPT ZMAT2F;
  2340.           'FINSI';
  2341.           'SI' (('NEG' HYPDEF 'LINEAIRE') 'ET' IGRD);
  2342.             'MESS' 'Attention utilisation des deformations lineaires';
  2343.           'FINS';
  2344.         'FINS';
  2345. *
  2346. *  ramener l'increment de deformation sur la config initiale pour diminuer les tests ensuite
  2347.          'SI' IGRD ;
  2348.            'FORM' GEOREF0;
  2349.            DEPST = 'CAPI' HYPDEF DEPST ZMODL ZDEP0 ;
  2350.          'FINSI';
  2351.          DPSMAX = 'MAXI' 'ABS' DEPST 'AVEC' MLDEFOR;
  2352. *
  2353. *       DEPST ZSIG0 et ZDEF0 sont sur la configuration initiale
  2354. *       Les exprimer sur la configuration adequate
  2355.         DEPSTa = DEPST;
  2356.         ZSIG0a = ZSIG0; 
  2357.         ZDEF0a = ZDEF0;
  2358.         ZDETOTN = -1. '*' ZDETOT ;
  2359.         'SI' (IGRD 'ET' iDEFOK);
  2360. *         Lagrangien total : sur GEOREF0
  2361.           'SI' ('EGA' LAG_TOT 1);
  2362.             'FORM' GEOREF0;
  2363.           'FINSI';
  2364. *         Lagrangien fin pas : GEOREF0 -> GEOM2
  2365.           'SI' ('EGA' LAG_TOT 2);
  2366.             'FORM' GEOM2;
  2367.             DEPSTa = 'CAPI' HYPDEF DEPST ZMODL ZDETOTN ;
  2368.             ZSIG0a = 'CAPI' HYPDEF ZSIG0 ZMODL ZDETOTN ;
  2369.             ZDEF0a = 'CAPI' HYPDEF ZDEF0 ZMODL ZDETOTN ;
  2370.           'FINS';
  2371. *         Lagrangien mi pas : GEOREF0 -> GEOM2M
  2372.           'SI' ('EGA' LAG_TOT 3);
  2373.             'FORM' GEOM2M;
  2374.             ZDEPMN = -1. '*' (ZDEP0 '+' ZDEPM) ;
  2375.             DEPSTa = 'CAPI' HYPDEF DEPST ZMODL ZDEPMN ;
  2376.             ZSIG0a = 'CAPI' HYPDEF ZSIG0 ZMODL ZDEPMN ;
  2377.             ZDEF0a = 'CAPI' HYPDEF ZDEF0 ZMODL ZDEPMN ;
  2378.           'FINS';
  2379.         'FINS';
  2380. *
  2381. *       Verifier si les deformations obtenues ne sont pas trop
  2382. *       importantes
  2383.         'SI' (IGRD 'ET' (('NON' HPP_EPS) 'OU' (IT > 0 )));
  2384.           ref = 20;
  2385.           n_sous_dec = MAXI 1 (n_sous_dec - 3) (((dpsmax/5d-3) + 1) ENTIER);
  2386.           n_sous_dec = MINI 10 n_sous_dec                          ;
  2387. **        mess 'sous decoupage pour epsi : ' n_sous_dec;
  2388.           coefmul = 1.d0 / (DPSMAX + XPETIT);
  2389.           'SI' (((coefmul < ref) 'OU'  PASUNIL) 'ET' AUTAUG);
  2390.             'MESS' 'annulation du sous-increment: trop grand ' dpsmax;
  2391.               si faux;
  2392.               RE_CNTRL2*'RIGIDITE'   = DEFO_IMP DEPSTa  (zdept + zdep1)  ZRAID ;
  2393. **pv          DE_CNTRL = VRAI ;
  2394.               ZFCONT = RE_CNTRL2 '*' ( ZDEP1) 'EXCO'  'FLX' 'NOID' 'FLX'  ;
  2395.               cmul = MINI  1. ( 1d-3 / DPSMAX);
  2396.               ZFEXT2   = ZFEXT2   + (ZFCONT * cmul);
  2397.               ZFLXB    = 'EXCO' ZFEXT2 'FLX' 'NOID' 'FLX' ;
  2398.               ZFCONSTA = ZFEXT2 ;
  2399.               RE_CNTRL = RE_CNTRL 'ET' RE_CNTRL2;
  2400.               finsi;
  2401. *          il faut reduire d'urgence les deplacements
  2402.             coefmul = 0.0;
  2403.             'SI' ((red_urg > 3) 'ET' ('NON' URG)); 
  2404.               resmul = resmul * ('MINI' 0.25 (1d-2 / dpsmax));
  2405.             'FINSI';
  2406.             augmult = augmult * 1.5;
  2407. * on garde zdep1d pour le calcul des coeff d'acceleration
  2408.             zdep1d = zdep1 'ENLE' 'LX';
  2409.             'SI' ((dpsmaxp < 5e-4) 'OU' (dpsmaxp > 2e-2 )) ;
  2410.               PASTEST=VRAI;
  2411.             'FINSI';
  2412.             'SI' (RED_URG > 2 )  ;
  2413.               'SI' ('NON' nonconv);  
  2414.                 'MESS' ' non convergence detectee 1' ;
  2415.                 ZMAXIT = IT+5;
  2416.               'FINS';
  2417.               nonconv = vrai;
  2418.             'FINS';
  2419.             'SI' (resmul < 0.99); 
  2420.               resmul = resmul * ('MINI' 0.25 (1d-2 / dpsmax));
  2421.             'FINSI';
  2422.             RED_URG = RED_URG + 1;
  2423.             'SI' ((dpsmaxp > 1.8d-2) 'ET' (augmult > 100));
  2424.               zprecnc= 2.1d-2;
  2425.               'SI' (zmaxit > it) ;
  2426.                 zmaxit = it - 1;
  2427.               'FINSI';
  2428.             'FINSI';
  2429.             'SI' ((zmaxit > it) 'ET' (RED_URG > 4 )); ZMAXIT = IT-1; 'FINS';
  2430.             'ITERER' etiq;
  2431.           'SINON';
  2432. **          RE_CNTRL = 'VIDE' 'RIGIDITE' / 'RIGIDITE';
  2433. **          DE_CNTRL = FAUX;
  2434.           'FINS';
  2435.         'FINS';
  2436. *----------------------------------------------------------------------*
  2437. *              Loi de comportement (operateur ELAS ou COMP)            *
  2438. *               sur la configuration definie par LAG_TOT               *
  2439. *----------------------------------------------------------------------*
  2440. *       Comportement elastique lineaire
  2441.         'SI' ('NON' IPLAVI) ;
  2442. *
  2443. *         Increment de contraintes
  2444.           DSIGTa = 'ELAS' ZMODL DEPSTa ZMAT2F ;
  2445. *
  2446. *         Soustraire les parties thermique et deformations imposees
  2447.           'SI' (ITHER 'OU' WTAB.'MATVAR') ;
  2448.             XXX1 = DSIGT0 '*' COEPI;
  2449.             XXX2 = DSIGTa  '-' XXX1;
  2450.             DSIGTa =  XXX2  ;
  2451.           'FINS';
  2452. *
  2453.           'SI' LOGDEF  ;
  2454.             XXX1 = DSI1  '*' COEPI;
  2455.             XXX2 = DSIGTa '-' XXX1;
  2456.             DSIGTa =  XXX2 ;
  2457.           'FINS';
  2458. **        mess 'depsta en elastique: ' (maxi abs depsta);
  2459. **        mess 'dsigta en elastique: ' (maxi abs dsigta);
  2460. *
  2461. *       Comportement non lineaire : integration avec l'operateur COMP
  2462.         'SINON';
  2463. *
  2464. *  dans le cas des modeles endommageables de Lemaitre, on ecoule
  2465. *  en tenant compte, dans les iterations internes, de la variation du
  2466. *  materiau avec la temperature
  2467. *
  2468. *         Soustraire les parties thermique et deformations imposees
  2469.           'SI' (ITHER 'OU' WTAB.'MATVAR') ;
  2470.             XXX1   = DEPST0 '*' COEPI;
  2471.             DEPST  = DEPST  '-' XXX1 ;
  2472.             DEPSTa = DEPSTa '-' XXX1 ;
  2473.           'FINS' ;
  2474. *
  2475.           'SI' LOGDEF  ;
  2476.             XXX1   = DDEFOR0 '*' COEPI;
  2477.             DEPST  = DEPST  '-' XXX1 ;
  2478.             DEPSTa = DEPSTa '-' XXX1 ;
  2479.           'FINS' ;
  2480. *
  2481.           ZMATT = ZMAT2F ;
  2482.           'SI' ('ET' ('ET' ITHER ('OU' IENDOM IVIDOM )) WTAB.'MATVAR') ;
  2483. * on recupere certain materiau avec les parametres fct de la temperature
  2484. * voir PAS_mate (il ne faut que la dependance thermique)
  2485.             ZMATT = 'REDU' WTAB.'MA_COMP' ZMODL;
  2486.           'FINS';
  2487. *
  2488. *         ...cas SSTE...
  2489.           'SI' ISSTE;
  2490.             ZRIKTA ZSIGFa ZVARF ZDEIF =
  2491.               'SSTE' ZMODL ZSIG0a ZVAR0 DEPSTa ZMATT
  2492.                      ZPREK WTAB.'NMAXSUBSTEPS' NITMA;
  2493.             zvar0 = ('EXCO' zvar0 ('ENLE' com_var ('DIME' com_var))) 
  2494.                       'ET' ('EXCO' zvarf ssii);
  2495.             zvar0 = 'CHANGER' 'TYPE' zvar0 'VARIABLES INTERNES';
  2496. *
  2497. *         ...cas non SSTE...
  2498.           'SINON';
  2499. *
  2500. *           Gradient des deplacements (total) a la fin du pas de temps 
  2501. *           par rapport a la configuration de reference GEOREF0
  2502.             'SI' IGRD;
  2503.               'SI' ('NEG' ZGRDU0 'INCONNU');
  2504.                 'SI' ('NEG' LAG_TOT 1); 
  2505.                   geot = 'FORM';
  2506.                   'FORM' GEOREF0 ; 
  2507.                 'FINS' ;
  2508.                 ZGRDUF = 'GRAD' ZMODL ZMAT2F ZDETOT ;
  2509.                 D_GR_U = ZGRDUF '-' ZGRDU0 ;
  2510.                 'SI' ('NEG' LAG_TOT 1); 'FORM' GEOT; 'FINS';
  2511.               'FINS';
  2512.             'FINS';
  2513. *
  2514.             ZSIG01 = ZSIG0a;
  2515.             ZVAR01 = ZVAR0;
  2516.             ZDEI01 = ZDEI0;
  2517.  
  2518.             'SI' (nsoincr 'NEG' 1);
  2519.               ZSIGM = ZSIG0a '/' 2;
  2520.             'FINS';
  2521. *
  2522.             'REPETER' sousinc nsoincr;
  2523. *
  2524.               codeb  = 'FLOT' (&sousinc '-' 1) '/' nsoincr ;
  2525.               cofin  = 'FLOT'  &sousinc        '/' nsoincr ;
  2526. *
  2527.               tsodeb = WTAB.'TEMPS0' '+' (ZDT '*' codeb) ;
  2528.               tsofin = WTAB.'TEMPS0' '+' (ZDT '*' cofin) ;
  2529.               'SI' (knoconv '>' 1) ;
  2530.                 tsodeb = tsofin ;
  2531.               'FINS' ;
  2532.               che1 = 'ADET' 'NOUV' ZMODL CHASANST tsodeb 'TEMP' tsodeb;
  2533.               che2 = 'ADET' 'NOUV' ZMODL CHASANST tsofin 'TEMP' tsofin;
  2534. *
  2535.               'SI' ('EXIS' WTAB 'MOD_LIA') ;
  2536.                  vite1 =  'CHAN' 'COMP' conti.'VITESSES' MLPRIM_LIA MVPRIM_LIA ;
  2537.                  che1  = che1
  2538.                    'ET' ('CHAN' 'CHAM' ZMODL vite1                'STRESSES')
  2539.                    'ET' ('CHAN' 'CHAM' ZMODL conti.'DEPLACEMENTS' 'STRESSES');
  2540.                  che2  = che2
  2541.                    'ET' ('CHAN' 'CHAM' ZMODL  ZDETOT              'STRESSES')
  2542.                    'ET' ('CHAN' 'CHAM' ZMODL  ZFCONST1            'STRESSES');
  2543.               'FINS';
  2544.  
  2545.               'SI' ither ;
  2546.                 TETDE = ZTET1 + (DTETD '*' codeb) ;
  2547.                 TETDF = ZTET1 + (DTETD '*' cofin) ;
  2548.                 che3 = 'CHAN' 'CHAM' ZMODL TETDE 'STRESSES' ;
  2549.                 che4 = 'CHAN' 'CHAM' ZMODL TETDF 'STRESSES' ;
  2550.               'FINS' ;
  2551.  
  2552.               che1 = che1 'ET' che3 ;
  2553.               che2 = che2 'ET' che4 ;
  2554.               aa5  = DEPSTa '*' codeb;
  2555.               che5 = ZDEF0a '+'   aa5;
  2556.               cc5  = DEPSTa '*' cofin;
  2557.               che6 = ZDEF0a '+'  CC5;
  2558. *On met les deformations en tete des champs pour COMP ('KTAN' 'PERT')
  2559.               che1 = che5 'ET' che1 ;
  2560.               che2 = che6 'ET' che2 ;
  2561. * pour les materiaux on garde toujours la valeur a la fin du pas
  2562. * et au debut du pas
  2563.               che1 = che1 'ET' ZMAT1 ;
  2564. *
  2565. *  Modele NON LINEAIRE UTILISATEUR + GRANDES DEFORMATIONS
  2566. *  On ajoute les gradients de deplacements qui seront transformes en
  2567. *  gradients de transformation avant appel a UMAT (dans WKUMA1.ESO)
  2568.               'SI' (( 'NEG' ZGRDU0 'INCONNU') 'ET' igrd) ;
  2569.                 gru1 = ZGRDU0 + (D_GR_U '*' codeb) ;
  2570.                 che1 = che1 'ET' gru1 ;
  2571.                 gru2 = ZGRDU0 + (D_GR_U '*' cofin) ;
  2572.                 che2 = che2 'ET' gru2 ;
  2573.               'FINS';
  2574. *
  2575.               'SI' LNLOC ;
  2576.                  ISTEP = 1 ;
  2577.                  chm_z = 'ADET' 'NOUV' ZMODL CHASANST tsodeb 'STEP' ISTEP ;
  2578.                  Che1a = che1 ;
  2579.                  che1  = che1 'ET' chm_z ;
  2580.                  chm_z = 'ADET' 'NOUV' ZMODL CHASANST tsofin 'STEP' ISTEP ;
  2581.                  Che2a = che2 ;
  2582.                  che2  = che2 'ET' chm_z ;
  2583.               'FINS';
  2584. *
  2585.               che11 = che1 'ET' ZSIG01 'ET' ZVAR01 'ET' ZDEI01 ;
  2586.               che22 = che2 'ET' ZMATT ;
  2587.               'SI' PARTLOCA;
  2588.                 cha11 = 'ASSIS' 'TOUS' 'REDU' MODRELOC che11 ;
  2589.                 cha22 = 'ASSIS' 'TOUS' 'REDU' MODRELOC che22 ;
  2590.                 souci 0 ;
  2591.                 cho22 = 'ASSIS' 'TOUS' 'COMP' MODRELOC cha11 cha22 ;
  2592.                 isoucomp = souci ;
  2593.                 cho2  = 'REDU' ZMODL ('ETG' cho22) ;
  2594.                 cha11 = 0 ; cha22 = 0 ; cho22 = 0 ;
  2595.               'SINON';
  2596.                 souci 0 ;
  2597.                 cho2  = 'COMP' ZMODL che11 che22 ;
  2598.                 isoucomp = souci ;
  2599.               'FINS';
  2600. *
  2601. *++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  2602. *             DEBUT MODELE NON LOCAL
  2603.               'SI' LNLOC ;
  2604.                 ZVARF = 'CHANGER' ('EXCO' cho2 com_var 'NOID') 'TYPE' 'VARIABLES INTERNES' ;
  2605. * NLOC ne traitant pas des champs //, on reduit tout au modele initial.
  2606.                 'SI' (PARALLEL 'ET' ('NON' PARTLOCA)) ;
  2607.                   ZVARF = 'REDU' ZVARF ZMODLI ;
  2608.                 'FINS' ;
  2609.  
  2610. ********** boucle sous-iterations Helmholtz **************
  2611.  
  2612. *               donnees pour acceleration Helmholtz 
  2613. *               periode de recours a  l acceleration ACT3 pour Helmholtz     
  2614.                 ZNACCEHL = 3;
  2615. *               NBR min de sous-iterations d Helmholtz  a faire avant d appeler ACT3 pour la 1ere fois
  2616.                 ITDEPHL  = 4;
  2617.  
  2618.                 'REPE' BIHL WTAB.'MAXSOUSITERATION' ;
  2619.  
  2620.                   'SI' ('NEG'  WTAB.'NON_LOCAL' 'HELM');
  2621. *                   cas non-local MOYE ou SB
  2622.                     'SI' ('EGA' WTAB.'NON_LOCAL' 'SB') ;
  2623.                       MOD_SB =  WTAB.'NLOC_MODL' ;
  2624.                       CONTP = 'PRIN' ('REDU' ZSIG01 MOD_SB) MOD_SB ;
  2625.                       ZVARF = ZVARF '+' CONTP '+' WTAB.'NLOC_SB_REGU' ;
  2626.                     'FINS' ;
  2627.                     ZVARN = 'NLOC'  ZVARF  WTAB.'CONN'  ;
  2628.                     ZVARN = 'CHANGER' ('EXCO' ZVARN com_var) 'TYPE' 'VARIABLES INTERNES' ;
  2629.                     LIERR1='PROG' 0 ; 
  2630.                     LOGHLIN = VRAI ;
  2631.                   'SINON' ;
  2632.  
  2633. *                   cas non-local HELM Sellier Millard
  2634.                     MOD_HELM =  WTAB.'NLOC_MODL' ;
  2635.                     ZVARN = 'COPI' ZVARF ;
  2636.                     ZVARF1 = 'REDU' ZVARN MOD_HELM ;
  2637.                     LICOHELM = 'EXTR' ZVARF1 'COMP' ;
  2638. *                   (re) lecture du materiau si init ou evolutif
  2639.                     ZMATHL = PAS_MATE PRECED cho2 ; 
  2640. *                   (re) calcul des matrices de rigidite de Helmholtz Sellier        
  2641.                     PAS_HELM PRECED ZMATHL;
  2642. *                   variable pour la liste des residus de convergence sur Helmholtz
  2643.                     LIERR1='PROG'; 
  2644.  
  2645. *                   ----- debut boucle sur les formulations de Helmholtz ---
  2646. *                   on reste en lineaire si toutes les vari sont lineaires (sinon on aura des istep a 3)
  2647.                     LOGHLIN = VRAI ;
  2648.                     'REPE' BH NHELM ;
  2649.                       LOGHLIN=LOGHLIN 'ET' TAHELM . &BH. 'LINEAIRE' ;
  2650. *                     nom de la variable a delocaliser
  2651.                       LEMO = TAHELM . &BH. 'NOM' ;
  2652. *                     recuperation de la variable sur le modele mecanique
  2653.                       ZVAUX = 'EXCO' ZVARF1 LEMO 'SCAL';
  2654.                       ZVAUX2 = 'CHAN' 'NOEUD' MOD_HELM ZVAUX;
  2655. *                     projection sur le modele de Helmholtz associe
  2656.                       ZVAUX2 = 'PROI' TAHELM . &BH . 'H_MODELE'  ZVAUX2 'MINI';
  2657. *                     TRAC zvaux2 TAHELM . &BH . 'H_MODELE' titre ('CHAI' 'Vari a diffuser:'LEMO );
  2658. *                     la variable devient la source
  2659.                       FSOUR = 'SOUR' (TAHELM . &BH . 'H_MODELE') ZVAUX2 ;
  2660. *                     trac fsour (extr TAHELM . &BH . 'H_MODELE' 'MAIL') titre 'Source';
  2661. *                     trac (TAHELM . &BH . 'H_OPER')  titre 'H_Oper';
  2662. *                     ***  pour istep=1 *** 1ere etape non locale ***             
  2663.                       'SI' ('EGA' &BIHL 1)  ;
  2664. *                       &BIHL est le compteur de sous iteration de Helmholtz
  2665. *                       il vaut forcement 1 si istep vaut 1, on fait la resolution classique
  2666.                         ZVNEW = 'RESO' (TAHELM . &BH . 'H_OPER') FSOUR ;
  2667. *                       sauvegarde pour possible acceleration ulterieure 
  2668.                         TAHELM . &BH . 'XN' = ZVNEW ;
  2669.                       'SINON';
  2670. *                       recuperation de la solution de la sous iteration de Helmholtz precedente
  2671.                         ZVNEW = TAHELM . &BH . 'XN' ;
  2672.                       'FINSI';
  2673. *                     *** traitement non lineaire de la formulation ***
  2674.                       'SI' ( 'NON'  TAHELM . &BH . 'LINEAIRE' ) ;
  2675. *                       calcul du residu sur l evolution de la source         
  2676.                         RESN = FSOUR - (TAHELM . &BH . 'H_OPER' '*' TAHELM . &BH . 'XN' );              
  2677. *                       stockage du residu pour la mesure d erreur                   
  2678.                         RESN0='NOMC' RESN 'SCAL' ;
  2679. *                       normalisation des residus par l inverse de la capacite 
  2680. *                       (calculee dans pas_helm.procedur)
  2681.                         ZERR0 = RESN0 '*'  (TAHELM . &BH . 'INV_CAPA' );
  2682. *                       erreur residu
  2683.                         ZERR2 = 'MAXI' ( 'ABS' (ZERR0) );
  2684. *                       'SI' ('EGA' &BIHL 1);
  2685. *                         TAHELM . &BH . 'RES1' = ZERR2;
  2686. *                       'FINSI';   
  2687. *                       carre de norme du residu 
  2688. *                       ZERR1 = 'XTX' ZERR0; 
  2689. *                       Residu 
  2690. *                       ZERR2 = (ZERR1)**0.5;
  2691. *                       normalisation de la source par la capacite
  2692.                         ZERR3 = FSOUR '*' (TAHELM . &BH . 'INV_CAPA' )('MOTS' 'Q') ('MOTS' 'SCAL') ('MOTS' 'SCAL');
  2693. *                       calcul erreur
  2694.                         ZERR6 = 'MAXI' ( 'ABS' (ZERR3) );
  2695. *                       carre de la norme de la source
  2696. *                       ZERR5 = 'XTX' ZERR3; 
  2697. *                       norme de la source
  2698. *                       ZERR6 = ZERR5 **0.5;
  2699. *                       valeur maximale du champs
  2700. *                       ZERR3='MAXI' (TAHELM . &BH . 'XN') 'ABS';
  2701. *                       mess 'variable Helmholtz', LEMO;
  2702. *                       mess 'residu', ZERR2;
  2703. *                       mess 'norme de la source ', ZERR5;                     
  2704. *                       mess 'norme Source', ZERR3;
  2705. *                       on teste si la valeur maximale est capable de normaliser
  2706.                         'SI' (ZERR6 '>' 0.);
  2707.                           ZERR4 = ZERR2 '/' ZERR6;
  2708.                         'SINON';
  2709.                           ZERR4 = ZERR2;
  2710.                         'FINS';
  2711. *                       preparation acceleration                              
  2712.                         'SI' ( 'EGA' &BIHL 1 ) ;         
  2713.                           TAHELM . &BH . 'H_RESO' = TAHELM . &BH . 'H_OPER' ;            
  2714.                           FWOR = FSOUR * 0. ;
  2715.                           TAHELM . &BH . 'ACFP1'  = FWOR ;
  2716.                           TAHELM . &BH . 'ACFP2'  = FWOR ;
  2717.                           TAHELM . &BH . 'RESNP'  = FWOR ;
  2718.                           TAHELM . &BH . 'ACFP3'  = FWOR ;
  2719.                           TAHELM . &BH . 'ACFEP1' = FWOR ;
  2720.                           TAHELM . &BH . 'ACFEP2' = FWOR ;
  2721.                           TAHELM . &BH . 'CORREC' = 0.;
  2722.                         'SINON' ;
  2723. *                         acceleration de convergence act3
  2724.                           CORRECPHL = TAHELM . &BH . 'CORREC';
  2725.                           CORRECHL = 0;
  2726.                           ACFP0HL = RESN - TAHELM . &BH . 'RESNP' ;
  2727.                           ACFEP0HL = ACFP0HL - CORRECPHL ;
  2728.                           'SI' ('MULT' &BIHL ZNACCEHL)   ;
  2729.                             'SI' (&BIHL > ITDEPHL);
  2730.                               CORRECHL   = 'ACT3' 
  2731.                               TAHELM . &BH . 'ACFEP2' 
  2732.                               TAHELM . &BH . 'ACFEP1'  ACFEP0HL 
  2733.                               TAHELM . &BH . 'ACFP3' 
  2734.                               TAHELM . &BH . 'ACFP2' 
  2735.                               TAHELM . &BH . 'ACFP1'   ACFP0HL ;
  2736.                               RESN = RESN - CORRECHL ;
  2737.                             'FINS';
  2738.                           'FINS';
  2739.                           TAHELM . &BH . 'RESNP'  = RESN ; 
  2740.                           TAHELM . &BH . 'ACFP3'  = TAHELM . &BH . 'ACFP2'  ; 
  2741.                           TAHELM . &BH . 'ACFP2'  = TAHELM . &BH . 'ACFP1'  ; 
  2742.                           TAHELM . &BH . 'ACFP1'  = ACFP0HL ;
  2743.                           TAHELM . &BH . 'ACFEP2' = TAHELM . &BH . 'ACFEP1' ; 
  2744.                           TAHELM . &BH . 'ACFEP1' = ACFEP0HL ;    
  2745. *                         calcul correction delta-sigma                               
  2746.                           DZVN = 'RESO' TAHELM . &BH . 'H_RESO' RESN;                 
  2747.                           ZVNEW = ZVNEW + DZVN ;
  2748.                           TAHELM . &BH . 'XN'   = ZVNEW ;
  2749.                         'FINSI' ;                          
  2750.                       'SINON';
  2751. *                       cette variable est traitee de facon lineaire                       
  2752.                         ZERR4=0.;
  2753.                       'FINSI';
  2754.                       ZVNEW1 = 'NOMC' ZVNEW 'SCAL' ;              
  2755.                       ZVNEW = 'CHAN' 'CHAM' ZVNEW1 MOD_HELM 'STRESSES' 'VARIABLES INTERNES';
  2756. *                     trac zvnew mod_helm titre 'ZVNEW';
  2757.                       DZVN = ZVNEW '-' ZVAUX ;
  2758.                       DZVN = 'NOMC' DZVN LEMO ;
  2759.                       ZVARN = ZVARN '+' DZVN ;
  2760.  
  2761. *                     mess 'Erreur retenue', ZERR4;
  2762. *                     stockage de l erreur normalisee dans la liste des erreurs
  2763.                       LIERR1 = LIERR1 'ET' ('PROG' ZERR4);
  2764.  
  2765. *                     cas particulier TWL2 Sellier theorie WL2
  2766.                       'SI' ('EGA' LEMO 'TWL2') ;
  2767.                         ZSMAX = 'EXCO' ZVARF1 'SMAX' 'SMAX' ;
  2768.                         ZSMAN = 'MANU' 'CHML' MOD_HELM 'SMAX' ('MAXI' ZSMAX) 'TYPE' 'VARIABLES INTERNES' 'STRESSES' ;
  2769.                         ZVARN = ZVARN '+' ZSMAN '-' ZSMAX ;
  2770.                       'FINS' ;
  2771. *                     ---- cas particulier SER Sellier renfort avec gradient
  2772. *                     calcul GRADIENT des contraintes sur les renforts (SERi)
  2773. *                     stockage dans la VARI HELM1(&BH,2) du gradient selon l absc curviligne 
  2774.                       MOPRE = 'EXTR' LEMO 1 3 ;
  2775.                       'SI'  ('EGA' MOPRE 'SHR'); 
  2776.                         NUME='EXTR' LEMO 4; 
  2777.                         NDIM1='VALE' 'DIME';
  2778.                         'SI' ( 'EGA' NDIM1 3 );
  2779. *                         valeur SERi du pas precedent dans Helm(&bh,1)                             
  2780.                           LEMO1=CHAI 'H' &BH '1'; 
  2781. *                         valeur SERi du pas actuel 
  2782.                           ZAUX1='EXCO' ZVARN LEMO 'SCAL';                 
  2783.                           ZVNEW1='NOMC' ZAUX1 'T' ;          
  2784.                           ZVAUX1 = 'CHAN' 'NOEUD' MOD_HELM ZVNEW1;              
  2785.                           ZVNEW2='PROI' TAHELM . &BH . 'H_MODELE' ZVAUX1 'MINI';
  2786.                           GRA1=GRAD TAHELM . &BH . 'H_MODELE' ZVNEW2;                 
  2787.                           GRA2=CHAN 'CHAM' 'NOEUD' GRA1 TAHELM . &BH . 'H_MODELE';
  2788.                           GRA3=PROI GRA2 MOD_HELM; 
  2789.                           LMO3=EXTR GRA3 COMP; 
  2790.                           NOM1=EXTR LMO3 1; 
  2791.                           NOM2=EXTR LMO3 2;
  2792.                           NOM3=EXTR LMO3 3;
  2793.                           GRA11=EXCO NOM1 GRA3 'SCAL';
  2794.                           GRA22=EXCO NOM2 GRA3 'SCAL';
  2795.                           GRA33=EXCO NOM3 GRA3 'SCAL';                  
  2796. *                         recuperation vecteur direction de projection
  2797.                           LEMON1=CHAI 'VR' NUME '1';
  2798.                           NX1=EXCO LEMON1 ZMATHL 'SCAL';
  2799.                           NX1=CHAN 'NOEUD' NX1 MOD_HELM;
  2800.                           NX11=CHAN NX1 TYPE 'SCALAIRE';
  2801.                           LEMON2=CHAI 'VR' NUME '2';
  2802.                           NX2=EXCO LEMON2 ZMATHL 'SCAL';
  2803.                           NX2=CHAN 'NOEUD' NX2 MOD_HELM;                  
  2804.                           NX22=CHAN NX2 TYPE 'SCALAIRE';
  2805.                           LEMON3=CHAI 'VR' NUME '3';
  2806.                           NX3=EXCO LEMON3 ZMATHL 'SCAL'; 
  2807.                           NX3=CHAN 'NOEUD' NX3 MOD_HELM;
  2808.                           NX33=CHAN NX3 TYPE 'SCALAIRE';
  2809. *                         produit scalaire                  
  2810.                           ZVAUX3=(NX11*GRA11)+(NX22*GRA22)+(NX33*GRA33);
  2811. *                         stockage du gradient de epse dans YU2j
  2812.                           mot4=CHAI 'H' NUME '2'; 
  2813.                           zvaux5= 'NOMC' zvaux3 mot4;
  2814.                           zvaux6= 'CHAN' 'CHAM' zvaux5 MOD_HELM 'STRESSES' 'VARIABLES INTERNES';
  2815. *                         trac zvaux6 MOD_HELM titre '3D H2';                                             
  2816.                           ZVAUX7 = 'EXCO' ZVARF1 MOT4 ;
  2817.                           DZVN = ZVAUX6 '-' ZVAUX7;
  2818.                           ZVARN = ZVARN '+' DZVN ;        
  2819.                         'SINON';
  2820. *                         valeur SERi du pas precedent dans Helm(&bh,1)                             
  2821.                           LEMO1=CHAI 'H' &BH '1';  
  2822.                           ZAUX1='EXCO' ZVARF1 LEMO1 'SCAL';                 
  2823.                           ZVNEW1='NOMC' ZAUX1 'T' ;          
  2824.                           ZVAUX1 = 'CHAN' 'NOEUD' MOD_HELM ZVNEW1;              
  2825.                           ZVNEW2='PROI' TAHELM . &BH . 'H_MODELE' ZVAUX1 'MINI';
  2826.                           GRA1=GRAD TAHELM . &BH . 'H_MODELE' ZVNEW2;                 
  2827.                           GRA2=CHAN 'CHAM' 'NOEUD' GRA1 TAHELM . &BH . 'H_MODELE';
  2828.                           GRA3=PROI GRA2 MOD_HELM; 
  2829.                           LMO3=EXTR GRA3 COMP; 
  2830.                           NOM1=EXTR LMO3 1; 
  2831.                           NOM2=EXTR LMO3 2;
  2832.                           GRA11=EXCO NOM1 GRA3 'SCAL';
  2833.                           GRA22=EXCO NOM2 GRA3 'SCAL';                 
  2834. *                         recuperation vecteur direction de projection
  2835.                           LEMON1=CHAI 'VR' NUME '1';
  2836.                           NX1=EXCO LEMON1 ZMATHL 'SCAL';
  2837.                           NX1=CHAN 'NOEUD' NX1 MOD_HELM;
  2838.                           NX11=CHAN NX1 TYPE 'SCALAIRE';
  2839.                           LEMON2=CHAI 'VR' NUME '2';
  2840.                           NX2=EXCO LEMON2 ZMATHL 'SCAL';
  2841.                           NX2=CHAN 'NOEUD' NX2 MOD_HELM;                  
  2842.                           NX22=CHAN NX2 TYPE 'SCALAIRE';
  2843. *                         produit scalaire                  
  2844.                           ZVAUX3=(NX11*GRA11)+(NX22*GRA22);
  2845. *                         stockage du gradient de epse dans YU2j
  2846.                           mot4=CHAI 'H' NUME '2'; 
  2847.                           zvaux5= 'NOMC' zvaux3 mot4;
  2848.                           zvaux6= 'CHAN' 'CHAM' zvaux5 MOD_HELM 'STRESSES' 'VARIABLES INTERNES';
  2849. *                         trac zvaux6 MOD_HELM titr '2DH2';                                            
  2850.                           ZVAUX7 = 'EXCO' ZVARF1 MOT4 ;
  2851.                           DZVN = ZVAUX6 '-' ZVAUX7;
  2852.                           ZVARN = ZVARN '+' DZVN ;        
  2853.                         'FINSI';     
  2854.                       'FINSI';
  2855. *                        ------fin du calcul des gradients  et Hessiens--------
  2856.                     'FIN' BH ;
  2857. *                 ------fin boucle sur les formulations de Helmholtz ---
  2858.                   'FINS' ;
  2859.  
  2860. *                 On rend ZVARN // si besoin :
  2861.                   'SI' (PARALLEL 'ET' ('NON' PARTLOCA)) ;
  2862.                     ZVARN = 'REDU' ZVARN ZMODL ;
  2863.                   'FINS' ;
  2864.  
  2865. *                 Recherche du residu maxi parmi les formulations non locales 
  2866.                   ISTEP0=ISTEP;
  2867. *                 residu maxi de tous les Helmholtz actifs
  2868.                   ZERR6='MAXI' LIERR1;
  2869.  
  2870.                   'SI' ('EGA' &BIHL WTAB.'MAXSOUSITERATION') ;
  2871.                     'MESSAGE' ' SOUS ITERATION HELMHOLTZ CONVERGENCE FORCEE';
  2872. *                   on indique au modele de mecanique que la sous 
  2873. *                   iteration de Helmholtz est achevee                 
  2874.                     ISTEP = 2 ;
  2875.                   'SINON' ;
  2876.                     'SI' ((( ZERR6 '<' WTAB.'PRECSOUSITERATION') 'ET' ('NEG' ISTEP 1)) 'OU' LOGHLIN ) ;
  2877. *                     on indique au modele de mecanique que la sous 
  2878. *                     iteration de Helmholtz est achevee 
  2879.                       ISTEP = 2 ;
  2880.                     'SINON';
  2881. *                     on indique au modele de mecanique que les sous
  2882. *                     iterations de Helmholtz continuent
  2883.                       ISTEP = 3 ;
  2884.                     'FINS';
  2885.                   'FINS' ;
  2886. *
  2887.                   'SI' ('EGA' ISTEP0 1);
  2888.                     mess1=chai ' SOUS ITERATION HELMHOLTZ:' &BIHL  ' ISTEP:' ISTEP0 ;
  2889.                   'SINON';
  2890.                     mess1=chai ' SOUS ITERATION HELMHOLTZ:' &BIHL  ' ISTEP:' ISTEP0 ' CRITERE:' ZERR6;
  2891.                   'FINSI';
  2892.                   'MESSAGE' mess1 ;           
  2893. *
  2894.                   chm_z = 'ADET' 'NOUV' ZMODL CHASANST tsodeb 'STEP' ISTEP ;
  2895.                   che11 = CHE1A 'ET' chm_z ;
  2896.                   che11 = che11 'ET' ZSIG01 'ET' ZVARN 'ET' ZDEI01 ;
  2897.  
  2898.                   chm_z = 'ADET' 'NOUV' ZMODL CHASANST tsofin 'STEP' ISTEP ;
  2899.                   che22 = CHE2A 'ET' chm_z ;
  2900.                   che22 = che22 'ET' ZMATT ;
  2901.  
  2902.                   'SI' PARTLOCA ;
  2903.                     cha11 = 'ASSIS' 'TOUS' 'REDU' MODRELOC che11 ;
  2904.                     cha22 = 'ASSIS' 'TOUS' 'REDU' MODRELOC che22 ;
  2905.                     souci 0 ;
  2906.                     cho22 = 'ASSIS' 'TOUS' 'COMP' MODRELOC cha11 cha22 ;
  2907.                     isoucomp = souci ;
  2908.                     cho2  = 'REDU' ZMODL ('ETG' cho22) ;
  2909.                     cha11 = 0 ; cha22 = 0 ; cho22 = 0 ;
  2910.                   'SINON';
  2911.                     souci 0 ;
  2912.                     cho2  = 'COMP' ZMODL che11 che22 ;
  2913.                     isoucomp = souci ;
  2914.                   'FINS' ;
  2915. *                         
  2916.                   'SI' ('EGA'  ISTEP 2 );
  2917.                     'QUITTER' BIHL;
  2918.                   'FINS';
  2919. *                 fin de boucle sous-iterations 
  2920.                   ZVARF = 'CHANGER' ('EXCO' cho2 com_var 'NOID') 'TYPE' 'VARIABLES INTERNES' ;
  2921. *                 NLOC ne traitant pas des champs //, on reduit tout au modele initial.
  2922.                   'SI' (PARALLEL 'ET' ('NON' PARTLOCA)) ;
  2923.                     ZVARF = 'REDU'  ZVARF ZMODLI ;
  2924.                   'FINS' ;                
  2925.                 'FIN' BIHL ; 
  2926. *      ------ fin des sous iterations de Helmholtz ------
  2927.  
  2928. *               Un peu de menage
  2929.                 CHE1A = 1 ;
  2930.                 CHE2A = 1 ;
  2931.                 ZVARN = 1 ;
  2932.                 ZVARF = 1 ;
  2933.                 chm_z = 1 ;
  2934.                 LIERR1 = 1 ;      
  2935.                 ZMATHL = 1; 
  2936.               'FINS' ;
  2937. *             FIN MODELE NON LOCAL
  2938. *++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  2939. *
  2940.               si isoucomp ;
  2941.                 mess ' Souci dans le comportement' ;
  2942.               fins ;
  2943.               ZSIGFa = 'CHANGER' ('EXCO' cho2 com_sig 'NOID')'TYPE' 'CONTRAINTES' ;
  2944.               ZVARF  = 'CHANGER' ('EXCO' cho2 com_var 'NOID')'TYPE' 'VARIABLES INTERNES' ;
  2945.               ZDEIF  = 'CHANGER' ('EXCO' cho2 com_dei 'NOID')'TYPE' 'DEFORMATIONS INELASTIQUES' ;
  2946. *
  2947.               'SI' ('EXIS' WTAB 'MOD_LIA') ;
  2948.                 ZFLIA = 'EXCO' cho2 'FLIA' 'NOID';
  2949.               'FINS' ;
  2950. *
  2951.               'SI' (nsoincr 'NEG' 1);
  2952.                 'SI' (&sousinc 'EGA' nsoincr);
  2953.                   aaa1=zsigfa/2;
  2954.                   aaa2=aaa1+ zsigm;
  2955.                   'DETR' aaa1; 'DETR' zsigm;
  2956.                   zsigm=aaa2;
  2957.                 'SINON';
  2958.                   aaa2=zsigfa+ zsigm;
  2959.                   'DETR' zsigm;
  2960.                   zsigm=aaa2;
  2961.                 'FINS';
  2962.               'FINS';
  2963. *
  2964.               zsig01 = zsigfa ;
  2965.               ZVAR01 = ZVARF ;
  2966.               ZDEI01 = ZDEIF ;
  2967.             'FIN' sousinc;
  2968. *
  2969. * Matrice tangente par perturbation evaluee pour la derniere iteration calculee
  2970. * A voir : cas grand deplacement ZSIGF et ZMAT2, cas poreux et thermique ZSIGF
  2971.             'SI' (IKTAN 'ET' IPERT) ;
  2972.               Z1COMP = che11 ;
  2973.               Z2COMP = che22 'ET' ZSIGFa ;
  2974.             'FINS' ;
  2975. *
  2976.             'SI' ((nsoincr 'NEG' 1) 'ET' IDYN);
  2977. * pour tenir compte de ce que le travail de la correction est 1/2 FU,
  2978. * on la multiplie par 2
  2979.               CONT  =(ZSIGFa + ZSIG0a)          ;
  2980.               ZSIGM = ZSIGM*(2. /nsoincr)     ;
  2981.               ZSIGM = ZSIGM - CONT            ;
  2982.               BZSIGM='BSIG' ZMODL ZSIGM ZMAT2F ;
  2983.             'FINS';
  2984.             zsig01=1;ZVAR01=1;ZDEI01=1;
  2985.  
  2986. *     le fin d'en dessous est le fin de "si isste ... sinon ... finsi"
  2987. *     on est encore dans "si iplavi ... finsi"
  2988.           'FINS';
  2989. *...fin du si ISSTE sinon ...
  2990. *
  2991. *         cas particulier poreux et thermique
  2992.           'SI' ITHER;
  2993.             'SI' POR1 ;
  2994.                ZSIGFa = ZSIGFa - ( COEPI * DMSRT0 ) ;
  2995.             'FINS';
  2996.           'FINS';
  2997. *
  2998. *  max de epse pendant l'iteration
  2999. *  nbre de points qui ont une evolution non lineaire
  3000.           'SI' ('NEG' WTAB.'MOVA' 'RIEN') ;
  3001.             XXX1 = 'EXCO' (WTAB.'MOVA') ZVARF ;
  3002.             ACC = 'ABS' ( XXX1 - ACC0 ) ;
  3003.             MMC = 'MASQUE' ACC 'SUPERIEUR' 1.D-10 'SOMME' ;
  3004.             'SI' (MMC > MMCMAX) ; MMCMAX = MMC ; 'FINS' ;
  3005.           'FINS';
  3006.  
  3007.           'SI' ('NON' WTAB.'MODAL') ;
  3008.             DSIGTa = ZSIGFa - ZSIG0a;
  3009.           'FINSI' ;
  3010. **        mess 'depsta en plastique: ' (maxi abs depsta);
  3011. **        mess 'dsigta en plastique: ' (maxi abs dsigta);
  3012. *  le FINSI suivant est celui de SI IPLAVI
  3013.         'FINS';
  3014. *
  3015. *       Exprimer les contraintes sur la configuration initiale
  3016. *         zdmax = maxi abs zsigf;
  3017. *         mess 'zdmax en 3014:' zdmax  ;
  3018.         'SI' ('NON' WTAB.'MODAL') ;
  3019.           DSIGT = DSIGTa ;
  3020.           'SI' (IGRD 'ET' iDEFOK) ;
  3021.             'FORM' GEOREF0 ;
  3022. *           Lagrangien total : deja sur GEOREF0
  3023.             'SI' ('EGA' LAG_TOT 1);
  3024.             'FINS';
  3025. *           Lagrangien fin pas : GEOM2 -> GEOREF0
  3026.             'SI' ('EGA' LAG_TOT 2);
  3027.               DSIGT = 'CAPI' HYPDEF DSIGTa ZMODL ZDETOT ;
  3028.             'FINS';
  3029. *           Lagrangien mi pas : GEOM2M -> GEOREF0
  3030.             'SI' ('EGA' LAG_TOT 3);
  3031.               DSIGT = 'CAPI' HYPDEF DSIGTa ZMODL (ZDEPM +ZDEP0) ;
  3032.             'FINS';
  3033.           'FINS';
  3034.           ZSIGF = ZSIG0 + DSIGT;
  3035.         'SINON' ;
  3036.           ZSIGF = ZSIGFa ;
  3037.         'FINS';
  3038. *
  3039.         DPSMAX = 'MAXI' 'ABS' DEPST 'AVEC' MLDEFOR ;
  3040.         ZDEFF  = ZDEF0 '+' DEPST ;
  3041.         EPSM   = 'MAXI' 'ABS' ZDEFF 'AVEC' MLDEFOR ;
  3042. *
  3043. * temporaire verification zsigf
  3044. ***     zsigtmp = 'ELAS' ZMODL zdeff  ZMAT2F ;
  3045. ***     zdiff = (zsigf - zsigtmp);
  3046. ***     zdmax = maxi abs zdiff;
  3047. ***     zdref= maxi abs zsigf;
  3048. ***     mess 'zdmax:' zdmax ' zdref:' zdref;
  3049. **        mess 'iraug rate' ' ' iraug ' ' irate;
  3050. **        zsigf = zsigtmp;
  3051. *
  3052. *----------------------------------------------------------------------*
  3053. *          Forces nodales equivalentes au champ de contraintes         *
  3054. *----------------------------------------------------------------------*
  3055.         'SI' IGRD ;
  3056.           'SI' ('NON' HPP_EPS) ;
  3057.             'SOUC' 0;
  3058.             'SI' WTAB.'ITCAR';
  3059.               'FORM' GEOM1 ;
  3060.               GEOM2 ZMAT2F = 'FORM' ZDEPT  ZMODLI ZMAT22;
  3061.               ZMAT2F       = 'REDU' ZMAT2F ZMODL  ;
  3062.             'SINON';
  3063.               GEOM2 = 'FORM' GEOM1 ZDEPT ;
  3064.             'FINS';
  3065. *           Changer de GEOREF0 -> GEOM2
  3066.             'SI' IJAUMA; 
  3067.               'FORM' GEOM1;
  3068.               ZSIGG = 'CAPI' HYPDEF ZSIGF ZMODL ZDEP0N ;
  3069.               'FORM' GEOM2;
  3070.               ZSIGG = 'CAPI' HYPDEF ZSIGG ZMODL (zdept * -1) ;
  3071.             'SINON';
  3072.               ZSIGG = 'CAPI' HYPDEF ZSIGF ZMODL ZDETOTN ;
  3073.              'FINSI';
  3074. *HHO : Modifications appel a BSIGMA (ajout du champ de deplacements necessaire en HHO)
  3075.             FEQU2 = 'BSIG' 'NOER' ZMODL ZSIGG ZMAT2F ZDETOT ;
  3076. *
  3077.             'SI' IJAUMA;
  3078.               'FORM' GEOREF0 ;
  3079.               ZSIGF = 'CAPI' HYPDEF ZSIGG ZMODL ZDETOT ;
  3080.             'FINSI' ;
  3081. *
  3082.             'FORM' GEOM1;
  3083.             'SI' (SOUCI) ;
  3084.               ZSIGG = 'CAPI' HYPDEF ZSIGG ZMODL ZDEPT ;
  3085.               MESS 'BSIG impossible : tentative sur geom1' ;
  3086. *HHO : Modifications appel a BSIGMA (ajout du champ de deplacements necessaire en HHO)
  3087.               FEQU2 = 'BSIG' 'NOER' ZMODL ZSIGG ZMAT2F ZDETOT ;
  3088.               pastest = vrai;
  3089.             'FINSI';
  3090.           'SINON';
  3091. *HHO : Modifications appel a BSIGMA (ajout du champ de deplacements necessaire en HHO)
  3092. *    config debut de pas pour predicteur hpp
  3093.             'FORM' GEOM1 ;
  3094.             ZSIGG = 'CAPI' HYPDEF ZSIGF ZMODL ZDEP0N;
  3095.             FEQU2 = 'BSIG' ZMODL ZSIGG ZMAT2F ZDETOT;
  3096.             GEOM2 = 'FORM' GEOM1 ZDEPT;
  3097.           'FINSI';
  3098.         'SINON';
  3099. *         config unique car HPP (igrd faux)                     
  3100.           FEQU2 = 'BSIG' ZMODL ZSIGF ZMAT2F ZDETOT ;
  3101.         'FINS';
  3102. *
  3103.         'SI' IRCON;
  3104.           FEQU2 = FEQU2 '+' ('REDU' MAI_CONS (RIG_CONS '*' ZDETOT));
  3105.         'FINS';
  3106. *----------------------------------------------------------------------*
  3107. *          Mise a jour des chargements dependants de la geometrie      *
  3108. *----------------------------------------------------------------------*
  3109.         iACTU = FAUX ;
  3110.         'SI' (LOGPRE 'ET' IGRD 'ET' ('NON' HPP_EPS)) ;
  3111.           'FORM' GEOM2 ;
  3112.           'SI' ('EXIS' WTAB 'MAT_PRE') ;
  3113.             ZFPEXTF = 'BSIG' MOP ZPEXTF WTAB.'MAT_PRE' ;
  3114.           'SINON' ;
  3115.             ZFPEXTF = 'BSIG' MOP ZPEXTF ;
  3116.           'FINS' ;
  3117.           ZFSUIV = ZFPEXTF ;
  3118.           iACTU  = VRAI ;
  3119.         'FINS' ;
  3120.         'FORM' GEOM1 ;
  3121. *
  3122.         'SI' WTAB.'PROCEDURE_CHARMECA';
  3123.           'SI' (IGRD 'ET' ('NON' HPP_EPS)) ;
  3124.             'FORM' GEOM2 ;
  3125.           'FINSI' ;
  3126.           TFP22    = CHARMECA PRECED TI ;
  3127.           ADDISEC2 = 'EXIS' TFP22 'ADDI_SECOND' ;
  3128.           'SI' ADDISEC2 ;
  3129.             FP22   = TFP22.'ADDI_SECOND' ;
  3130.             'SI' iACTU ;
  3131.               ZFSUIV = ZFSUIV '+' FP22 ;
  3132.             'SINON';
  3133.               ZFSUIV = FP22 ;
  3134.             'FINS';
  3135.           'FINS';
  3136.         'FINS' ;
  3137.         'FORM' GEOM1 ;
  3138. *
  3139.       'FINS';
  3140. * le finsi au dessus est le finsi de si ifefp
  3141. *
  3142.       'SI' ('EXIS' WTAB 'MOD_LIA') ;
  3143.         FLIAI = 0.d0 ;
  3144.         'SI' ('EGA' ('TYPE' ZFLIA) 'MCHAML  ') ;
  3145.           NZLIA = 'EXTR' ZFLIA 'NBZO' ;
  3146.           'SI' (NZLIA > 0) ;
  3147.             'REPETER' BZLIA NZLIA ;
  3148. * un point support par zone - 2010 kich
  3149.               FLIAI= FLIAI +
  3150.   ('EXCO' ('EXTR' ZFLIA 'FLIA' &BZLIA 1 1) MLDUAL_LIA 'NOID') ;
  3151.             'FIN' BZLIA ;
  3152.           'FINS' ;
  3153.           FEQU2 = FEQU2 + FLIAI ;
  3154.         'FINS' ;
  3155.       'FINS' ;
  3156. *
  3157.       'SI' ISOL  ;
  3158.         'SI' IPLAVI ;
  3159.           ZMATPF = 'CHANGER' ('EXCO' cho2 com_maa 'NOID')
  3160.           'TYPE' 'CARACTERISTIQUES';
  3161.           MA_POR = 'REDU' ZMATPF MO_POR ;
  3162.         'FINS';
  3163.         XXX1 = 'GRAD' MO_POR ZDETOT MA_POR 'CONS'  ;
  3164.         XXXS = (1. - WTAB.'TETA' )*GRAP0 + (WTAB.'TETA'*XXX1) ;
  3165.         XXX2 = ZDT * ('GNFL' MO_POR XXXS) ;
  3166.         XXX3 = FEQU2 ;
  3167.         FEQU2 = XXX3 - XXX2;'DETR' XXX3; 'DETR' XXX2 ;
  3168.       'FINS' ;
  3169. *
  3170. * si on a fait de la sous incrementation ( a cause du dynamique
  3171.       'SI' ((nsoincr > 1) 'ET' IDYN);
  3172.         FEQU2 = FEQU2 + BZSIGM ;
  3173.       'FINS';
  3174.     'FINSI' ;
  3175. *FOR_MECA
  3176. *----------------------------------------------------------------------*
  3177. *          Equilibre (RESIDC) et nouveau second membre (RESIDU)        *
  3178. *----------------------------------------------------------------------*
  3179. *   Ajout dans le residu de contributions propres a la dynamique :
  3180. *   [ 4/h2*m + 2/h *c ] * du
  3181.     'SI' IDYN;
  3182.       'DETR' FFDYN;
  3183.       FFDYN = 'COPIER' FEQU2 ;
  3184.       XXX1 = WTAB.'MASSE' * ZDEPT;
  3185.       XXX3 = 4. * UNSURH * UNSURH  * XXX1 ;
  3186.       'DETR' XXX1;
  3187.       'SI' ('EXIS' WTAB 'MOD_LIA');
  3188.         'SI' (NZLIA > 0) ;
  3189.           XVIT2 = 0.d0;
  3190.           'REPETER' BZLIA NZLIA ;
  3191.             XVIT2 = XVIT2 +
  3192.        ('EXCO' ('EXTR' ZFLIA 'FLIA' &BZLIA 1 1) MVPRIM_LIA MLPRIM_LIA 'NOID') ;
  3193.           'FIN' BZLIA ;
  3194.         'FINS' ;
  3195.       'FINS' ;
  3196.       'SI' ( 'NEG' WTAB.'AMORTISSEMENT' 'INCONNU') ;
  3197.         XXX1 = WTAB.'AMORTISSEMENT' * ZDEPT;
  3198.         XXX2 = 2. * UNSURH * XXX1 ;
  3199.         XXX4 = XXX3 + XXX2 ;
  3200.         XXX3 = XXX4 ;
  3201.       'FINS' ;
  3202.       XXX4 = FEQU2 + XXX3;
  3203.       'DETR' XXX3 ; 'DETR' FEQU2 ;
  3204.       FEQU2 = XXX4 ;
  3205. *
  3206. *         forces correctrices en cas de liaison persistante :
  3207. *  on veut avoir (forces inertielles + forces visqueuses) compatibles
  3208. *  avec accelerations et vitesses relatives nulles aux points de contact
  3209. *  (pendant le contact). On modifie fequ2 ---> residu et l'iteration sui
  3210. *  fournira les bonnes reactions
  3211. *
  3212.       'SI' IMPLP;
  3213.         VADD XXX3 = VITETFOR ZRAID_T WTAB FCORU ;
  3214.         XXX4      = FEQU2 - XXX3;
  3215.         'DETR' XXX3 ; 'DETR' FEQU2 ;
  3216.         FEQU2     = XXX4;
  3217.       'FINS' ;
  3218.     'FINS';
  3219. *
  3220.     'SI' WTAB.'FOR_MECA' ;
  3221.       XXX1 = ZFEXT ;
  3222.       XXX1 = XXX1 '+' (ZFSUIV '*' COEPI) ;
  3223. *
  3224. *     Modification du chargement
  3225.       'SI' LOGPIL ;      
  3226.         XXX1 = XXX1 '+' (ETA '*' ZFPILIN);      
  3227.       'FINSI' ;
  3228. *
  3229. *     Forces exterieures (+ autres termes p.ex. en dynamique ou poreux)
  3230. *     sans reactions - forces interieures
  3231.       RESIDU = XXX1 '-' FEQU2;
  3232. *
  3233. *     Forces exterieures (+ autres termes p.ex. en dynamique ou poreux)
  3234. *     + reactions
  3235.       FEQUI = XXX1 '-' FCORF;
  3236. *
  3237. *     Forces exterieures (+ autres termes p.ex. en dynamique ou poreux)
  3238. *     + reactions - forces interieures
  3239.       RESIDC = RESIDU '-' FCORF;
  3240. *
  3241. *     Mise a jour des jeux (FLX) qui travaillent en incremental :
  3242. *     - pour les depl imposes cela revient a imposer un increment nul
  3243. *     - pour les rela unilaterales cela revient a mettre a jour le jeu
  3244.       'SI' (IPILOT 'ET' IMPO12);
  3245.         RESIDU = 'ENLE' RESIDU 'FLX' ;
  3246.         XJEUX  = (DUIMPO '*' COEPI) '+' DUUNIL ;
  3247.         XXX3   = XJEUX   '-' FCORU ;
  3248.         RESIDU = RESIDU '+' XXX3  ;
  3249.       'SINON' ;
  3250.         XXX1   = RESIDU '+' ZFLX1; 'DETR' RESIDU;
  3251.         RESIDU = XXX1   '-' FCORU; 'DETR' XXX1;
  3252.       'FINS';
  3253. *----------------------------------------------------------------------*
  3254. * Calcul du nouveau alpha pour le pilotage
  3255. *----------------------------------------------------------------------*
  3256.       'SI' (IPILOT 'ET' PASTEST) ;
  3257. *
  3258.         'SI' ('OU' IVISCO IVIDOM ITHER LOGDEF);
  3259.           'MESS' 'ALPHA calcule avec la norme de l increment' ;
  3260.           COEPI0 = COEPI;
  3261.           COEPI  = 'MINI' ('PROG' 1. (ZAUTOCOEF * al1) );
  3262. *         les acceleration auront lieu toutes les 4 pas
  3263.           COEINC = COEPI0 - COEPI ;
  3264.         'SINON' ;
  3265. *         Les multiplicateurs sont affectes d'un signe -
  3266.           XXX4 = ZDEPT '-' (2.D0 '*' ('EXCO' ZDEPT 'LX' 'NOID' 'LX')); 
  3267.           XXX5 = RESIDNOR* 1.D0;          
  3268.           XXX6 = DFEXT0L ;         
  3269.           'SI' IMPO12 ;
  3270.             XXX5 = XXX5 - DUUNIL ;
  3271.             XXX6 = DUIMPO ;
  3272.           'FINS';
  3273.           XXX3 = COEPI * ( 1 - AL1 ) * XXX6 ;          
  3274.           XXX2 = XXX3 + RESIDC ; 'DETR' XXX3;
  3275.  
  3276.           XX1 ='XTY' XXX2 XXX4 MLDUAL MLPRIM;
  3277.           XX2 ='XTY' XXX5 XXX4 MLDUAL MLPRIM;
  3278.           COEINC =  XX1 '/' XX2;
  3279.           'DETR' XXX4;'DETR' XXX2;'DETR' XXX5;
  3280. *
  3281.           COEPI0 = COEPI;
  3282.           COEPI = 'MINI' ('PROG' (COEPI0 - COEINC) 1.d0 ) ;
  3283.           'SI' ((COEPI0 'EGA' 1d0) 'ET' (AL1 'EGA' 1D0));
  3284.             COEPI = 1d0;
  3285.           'FINS';
  3286.           COEINC = COEPI0 - COEPI;
  3287.           'SI' ((COEPI 'EGA' 1d0) 'ET' (COEPI0 'NEG' 1d0));ITACC=3; 'FINS';
  3288.         'FINS' ;
  3289.  
  3290. *       Mise a jour des termes force de l'acceleration de convergence
  3291.         XXXZ = RESIDNOR ;
  3292.         'SI' WTAB.'PROCEDURE_CHARMECA';
  3293.           XXXZ = RESIDNOR '+' (FP22 '-' FP022) ;
  3294.         'FINS';
  3295.         XXX1 = 'ENLE' (COEINC '*' XXXZ) 'FLX' ;
  3296.         XXX3 = ACFEP2 - XXX1 ;
  3297.         XXX4 = ACFEP1 - XXX1 ;
  3298.         'DETR' ACFEP2; ACFEP2= XXX3;
  3299.         'DETR' ACFEP1; ACFEP1= XXX4;
  3300.         'DETR' XXX1;
  3301. *
  3302. *   mise a jour de ZFEXT en fonction de coepi on ne met a jour
  3303. *   de maniere definitive que les forces, les modifs sur les FLX
  3304. *   pour reimposer les conditions de deplacements ne
  3305. *   sont faites que pour le premier residu attention au cas des jeux
  3306. *   sortis de IMPO qu'il ne faut pas toucher
  3307. *
  3308.         DFEXT = DFEXT0F '*' COEPI ;
  3309.         ZFEXT = DFEXT   '+' FEXT0 ;
  3310.         XXX2  = ZFEXT   '+' (COEPI '*' ZFSUIV) ;
  3311.         XXX1  = XXX2 - FEQU2;
  3312. *
  3313.         XXX3  = COEPI * DFEXT0L;
  3314.         XXX3  = XXX3 + FLXINI - FCORU;
  3315.         RESIDU = XXX1 + XXX3;
  3316.         'SI' IMPO12;
  3317.           RESIDU = 'ENLE' RESIDU 'FLX' ;
  3318.           XXX3 = COEPI '*' DUIMPO ;
  3319.           RESIDU = RESIDU + XXX3 '+' DUUNIL ;
  3320.         'FINS';
  3321. *
  3322.         ZAUTOCOEF = COEPI ;
  3323. *
  3324. *        test si on est sur une branche descendante
  3325. *              corprec est la correction de precision
  3326.         CORPREC = 1. ;
  3327.         'SI' (COEPI < 0) ;
  3328.           CORPREC = WTAB.'PRECDECHARGE' / ZPREC;
  3329.           'MESS' 'Chargement decroissant. Critere divise par ' CORPREC;
  3330.         'FINS';
  3331.         'MESS' 'nouvelle valeur de ALPHA' COEPI;
  3332.       'FINS' ;
  3333.       'DETR' FCORU;
  3334. *----------------------------------------------------------------------*
  3335. *             Critere pour statuer de la convergence ou non            *
  3336. *----------------------------------------------------------------------*
  3337.       XNUMF = 'MAXI' RESIDC 'ABS' 'SANS' MXMFLX;
  3338.       XNUMG = 'MAXI' ((ZDEPL '-' ZDEPLP) '*' ZCLIM) 'ABS';
  3339.       XNUMF = 'MAXI' ('PROG' XNUMF XNUMG);
  3340.       'SI' TSTMOM ;
  3341.         XNUMM = 'MAXI' RESIDC 'ABS' 'AVEC' MXMYMZ;
  3342.       'FINS';
  3343. *
  3344.       ZFAU1 = FEQUI '-' ZFPLO ;
  3345.       XAUXF = 'MAXI' ZFAU1 'ABS' 'SANS' MXMFLX;
  3346.       'SI' TSTMOM ;
  3347.         XAUXM = 'MAXI' ZFAU1 'ABS' 'AVEC' MXMYMZ;
  3348.       'FINS';
  3349. *
  3350.       'SI' IFTOL;
  3351.         ZPREC = ZFTOL ;
  3352.         XDENO = 1.;
  3353.       'FINS';
  3354.       'SI' IMTOL;
  3355.         ZPRECM = ZMTOL ;
  3356.         XDENOM = 1.;
  3357.       'FINS';
  3358. *
  3359.       XCONV  = XNUMF '/' XDENO ;
  3360.       XCONVM = XCONV ;
  3361. *
  3362.       'SI' TSTMOM ;
  3363.         'SI' ('<' XAUXM XPETIT) ;
  3364.           XCONVM = 0.;
  3365.         'SINON';
  3366.           'SI' ('<' (XAUXF / XAUXM) 1.D-12) ;
  3367.             XCONV = 0.;
  3368.           'FINS';
  3369.           XCONVM = XNUMM / XDENOM ;
  3370.         'FINS';
  3371.       'FINS' ;
  3372.  
  3373. *  residu etant calcule on met a jour zdep1d qui doit etre coherent avec residu
  3374.       zdep1d = zdep1 'ENLE' 'LX';
  3375. **
  3376.     'FINSI' ;
  3377. *FOR_MECA
  3378. *
  3379.     'SI' WTAB.'NVSTNL'                               ;
  3380.         ZFLCONV = FAUX                               ; 
  3381.         mofofl = 'EXTR' WTAB.'MOD_NSL' 'FORC   '     ;
  3382.           dvit2 = 'XTX'  ('EXCO' ZVIF1 movifl)       ;
  3383.           err  = dvit2 '**' 0.5                      ;
  3384.           res2 = 'XTX' ('EXCO' Resnsl mofofl 'NOID') ; 
  3385.           errres2  = res2 **0.5                      ;
  3386.       'MESSAGE' ' ite:' IT*26 ' errdu:' err*52 ' errres:' errres2*78  ; 
  3387.        'SI' (err < 1.d-9); ZFLCONV = VRAI ;   'FINSI'; 
  3388.        'SI' (ZFLCONV 'OU' (it '>EG' WTAB.'MAXSOUSPAS')) ; 
  3389.                              'QUITTER' ETIQ ; 'FINSI' ;
  3390.     'FINSI' ;
  3391. *----------------------------------------------------------------------*
  3392. *                resume de l'iteration                                 *
  3393. *----------------------------------------------------------------------*
  3394.     'SI' WTAB.'FOR_MECA' ;
  3395.       'MESS'  IT*13 MMC*26 XCONV*39  DPSMAX*52 EPSM*65 XCONVM*78  ;
  3396. *       sauver depstp pour le test increment de deformation
  3397.       depstps = depstp;
  3398. *       acceptation de l'etat courant pour le sous pas suivant
  3399.       depstp zdeptp zsigfp fcorfp correcp = depst zdept zsigf fcorf correc;
  3400.       IRATE = FAUX;
  3401. *       premiere iteration prise en compte
  3402.       iprem = faux;
  3403. *----------------------------------------------------------------------*
  3404. *          test de convergence : equilibre de la structure             *
  3405. *----------------------------------------------------------------------*
  3406.       TABCONV.IT=XCONV;
  3407. *
  3408.       'SI' IPILOT;'SI' ((MMC 'EGA' 0) 'ET' (MMCMAX > 0.) 'ET' (COEPI < 0.));
  3409. *       on refuse de converger si on est elastique et en decharge
  3410.         COEPI = 'ABS' COEPI;PASTEST = VRAI;
  3411.       'FINS'; 'FINS';
  3412. *
  3413. *     Si frottement, faire au moins 2 iterations
  3414.       'SI' (WTAB.'CAFROTTE' 'ET' IMPO12 'ET' (IT < 2 ));
  3415.         PASTEST = VRAI;
  3416.       'FINS';
  3417. *
  3418.       PASUNIL = FAUX;
  3419.       'SI' ('EXIS' ZRAID_T 'OK');
  3420.         PASUNIL = 'NON' ZRAID_T.'OK' ;
  3421.       'FINS';
  3422. *
  3423.       'SI' ('NON' PASTEST) ;
  3424. *
  3425. *       Variation de Despi entre 2 itérés
  3426.         DEPSTD  = DEPST - DEPSTPS ;
  3427.         DEPSTDM = 'MAXI' 'ABS' DEPSTD 'AVEC' MLDEFOR;
  3428. *
  3429. *       Pas de test sur la convergence a la premiere iteration
  3430. *       ou apres une initialisation a partir de la solution precedente
  3431.         'SI' ((IT > 1) 'OU' ('NON' INIT)) ;
  3432. *
  3433. *         Si les criteres (deplacements + moments) sont < precision
  3434. *         et la variation sur Depsi est < precision souhaitee : on a convergé!
  3435.           ZPRECHPP = ZPREC; 
  3436.           'SI' HPP_EPS; ZPRECHPP = ZPREC * 1d-1; 'FINSI';
  3437.  
  3438.           'SI' ( XCONV < (ZPRECHPP * CORPREC )) ;
  3439.             'SI' ( XCONVM < (ZPRECM * CORPREC )) ;
  3440.               'SI' ( DEPSTDM < (ZPRECHPP * CORPREC )) ;
  3441. *
  3442.                 'SI' PASUNIL ;
  3443.                   'MESS' ' ****** NON CONVERGENCE DE CONTACT A L''ITERATION' ' ' IT ' SOUS-PAS' ' ' WTAB.'ISOUSPAS';
  3444.                   ZICONV = FAUX;
  3445.                   'QUITTER' ETIQ;
  3446.                 'FINS';
  3447. *
  3448.                 'SI' (IRAUG 'ET' AUTAUG 'ET' ((DPSMAX > ZPREC) 'OU' (augmult > 0.6e0)) 'OU' 
  3449.                       (resmul < 0.99) 'OU' DE_CNTRL);
  3450.                   'MESS' ' ****** NON CONVERGENCE DUE A L AUGMENTATION A L''ITERATION' ' ' IT ' SOUS-PAS' ' ' WTAB.'ISOUSPAS';
  3451.                   ZICONV = FAUX;
  3452.                   'QUITTER' ETIQ;
  3453.                 'FINS';
  3454. *
  3455.                 iRECA = FAUX ;
  3456.                 'SI' (WTAB.'CONTACT' 'ET' ('NON' HPP_EPS) 'ET' ('NON' NONCONV));
  3457.                   'SI' ('NON' (WTAB.'MODAL' 'OU' WTAB.'FROCABL')) ;
  3458.                     MODCON = WTAB.'MODCONTA';
  3459.                     'FORM' ZDEPT ;
  3460.                     'SI' ('EXIS' WTAB 'MATCONTA') ;
  3461.                       WTAB.'MATCONTA' = 'REDU' MODCON WTAB.'MATCONTA';
  3462.                       CJEU CRRN = 'RFCO' MODCON  WTAB.'CONV' WTAB.'MATCONTA' ;
  3463.                     'SINON' ;
  3464.                       CJEU CRRN = 'RFCO' MODCON  WTAB.'CONV'                ;
  3465.                     'FINS';
  3466.                     'FORM' GEOM1 ;
  3467.                     ZDLX = 'EXCO' ZDEPT 'LX' 'NOID' 'LX' ;
  3468. * 
  3469.                     REA1 = 'REAC' CRR  ZDLX ;
  3470.                     CMPF = 'EXTR' REA1 'COMP' ;
  3471.                     PSC1 = 'PSCA' REA1 CMPF REA1 CMPF ;
  3472.                     NOR1 = 'MAXI' 'ABS' PSC1 ;
  3473. *
  3474.                     'SI' ('>' NOR1 XDENO) ;
  3475.                       REA2 = 'REAC' CRRN ZDLX ;
  3476.                       DREA = REA2 '-' REA1 ;
  3477.                       PSC2 = 'PSCA' DREA CMPF DREA CMPF ;
  3478. *
  3479.                       RAP4 = ('MAXI' 'ABS' PSC2) '/' NOR1 ;
  3480. *
  3481.                       'SI' ('>EG' RAP4 ZPREC) ;
  3482.                         'MESS' 'Actualiser les conditions de contact-frottement' rap4 nor1;
  3483. *                       En gardant les statuts courants
  3484. *                       ZRIBLO_M = 'MOT' 'INCONNU' ;
  3485.                         iRECA = VRAI ;
  3486.                       'FINSI' ;
  3487.                     'FINSI' ;
  3488.                   'FINSI' ;
  3489.                 'FINSI' ;
  3490. *
  3491. *               Tout semble ok mais il reste encore a verifier que
  3492. *               - le resultat est obtenu sans predicteur HPP
  3493. *               - les relations de cont-fro sont "coherentes" avec la conf finale
  3494.                 'SI' (HPP_EPS 'OU' isoucomp 'OU' iRECA) ;
  3495.                   'SI' (IGRD 'ET' HPP_EPS) ; 
  3496.                     'MESS' ' passage en grands deplacements'; 
  3497.                   'FINSI' ;
  3498. *
  3499.                   HPP_EPS = FAUX;
  3500.                   URG     = VRAI;
  3501.                   ITACC   = 4;
  3502. *
  3503.                   DLTAIT  = WTAB.'DELTAITER' '+' IT ;
  3504.                   PASTEST = VRAI;
  3505.                   'SI' autaug; iraug = faux; de_cntrl = faux ; 'FINSI';
  3506.                 'SINON';
  3507.                   'MESS' ' ****** CONVERGENCE A L ITERATION' ' ' IT '  SOUS-PAS' ' ' WTAB.'ISOUSPAS';
  3508.                   'MESS' ' ';
  3509. *                 on sauve le nombre de points plastique
  3510.                   ZNBPLAS = MMC ;
  3511.                   ZNITE   = ZNITE '+' IT;
  3512.                   'QUITTER' ETIQ ;
  3513.                 'FINS';
  3514.               'FINS';
  3515.             'FINS';
  3516.           'FINS';
  3517.         'FINS';
  3518.       'FINS' ;
  3519. *----------------------------------------------------------------------*
  3520. * test de non convergence
  3521. *----------------------------------------------------------------------*
  3522.       XCONVREF=1E50;
  3523.       'SI' (IT > DLTAIT);
  3524.         XCONVREF=TABCONV.(IT-DLTAIT) * 0.99;
  3525.       'FINS';
  3526. *   si on a depasse le nombre max d'iterations ou si le residu augmente
  3527. *   ou si on aurait du converger et que cela n'est pas le cas :
  3528. *   => non convergence detectee !
  3529.       'SI' (('NON' PASTEST) 'ET' ('NON' NONCONV) 'ET' (IT > 1) 'ET' 
  3530.             ('NON' IPILOT)  'ET' ('NON' PASUNIL)) ;
  3531.         'SI' ((IT '>EG' ZMAXIT) 'OU' (XCONV > XCONVREF) 'OU'
  3532.               ((XCONV < ZPRECHPP) 'ET' (XCONVM < ZPRECHPP) 'ET' (DEPSTDM <
  3533.                ZPRECHPP))) ;
  3534.          'MESS' ' non convergence detectee 2' ;ZMAXIT = 3 * IT;
  3535.          NONCONV = VRAI;
  3536.          PASTEST = VRAI;
  3537. *  test pv
  3538. *
  3539.           HPP_EPS = FAUX;
  3540.           URG     = VRAI;
  3541.           ITACC   = 4;
  3542. *
  3543.           COEPI = WTAB.'RELAXATION_NONCONV';
  3544.         'FINS';
  3545.       'FINS';
  3546. *
  3547. *     Changement de la precision en non convergence
  3548.       'SI' nonconv ;
  3549. *       'SI' (it > (zmaxit *2 /3));
  3550.         'SI' (it > 15); 
  3551.           zprecnc=zprecnc*2;
  3552.         'FINS';
  3553.       'FINS';
  3554. *
  3555.       'SI' (('NON' PASTEST) 'ET' (dpsmax < zprecnc));
  3556.         'SI' ((IT '>EG' ZMAXIT) 'OU' 
  3557.               (NONCONV 
  3558.                'ET' (XCONV '>' XCONVREF)
  3559.                'ET' (DPSMAX '>EG' (DPSMAXP /4 ))) );
  3560.           'MESS' ' ****** NON CONVERGENCE A L''ITERATION' ' ' IT ' SOUS-PAS' ' ' WTAB.'ISOUSPAS';
  3561.           ZICONV = FAUX;
  3562.           'QUITTER' ETIQ;
  3563.         'FINS' ;
  3564.       'FINS' ;
  3565. *----------------------------------------------------------------------*
  3566. *     Menage d'objets temporaires
  3567. *----------------------------------------------------------------------*
  3568.       'DETR' FEQU2; 'DETR' ZDEFF;
  3569. *
  3570. *     Supprimer les configurations intermediaires
  3571.       'SI' (IGRD 'ET' ('NON' HPP_EPS));
  3572.         'DETR' GEOM2 ;
  3573.       'FINS' ;
  3574. *
  3575.       'SI' IFEFP ; 
  3576.         'DETR' ZRIKTA ; 'DETR' GEOM2; 
  3577.       'FINS' ;
  3578. *
  3579.     'FINSI' ;
  3580. *FOR_MECA
  3581. *
  3582.     dpsmaxp = dpsmax;
  3583. *
  3584.   'FIN' ETIQ ;
  3585. *----------------------------------------------------------------------*
  3586. *           Fin de la boucle de convergence ETIQ                       *
  3587. *----------------------------------------------------------------------*
  3588. *######################################################################*
  3589.   'SI' IDYN  ;
  3590.     VITI = UNSURH * 2. * ZDEPT - conti. 'VITESSES' ;
  3591.     VITI = 'ENLE' VITI 'LX' ;
  3592.     ZFP  = ZDYFEXT '+' ZFSUIV '-' FFDYN;
  3593. *
  3594.     'SI' ('EXIS' WTAB 'MOD_LIA');
  3595.       'SI' ('EXIS' WTAB.'MOD_LIA' 'MATE' 'NEWMARK') ;
  3596.         'SI' ('EGA' 'CHPOINT  ' ('TYPE' FLIAI) ) ;
  3597.           MODNEW = 'EXTR' WTAB.'MOD_LIA' 'MATE' 'NEWMARK' ;
  3598.           MAUNEW = 'EXTR' MODNEW 'MAIL' ;
  3599.           'SI' ('EXIS' ZFLIA MODNEW) ;
  3600.             ZFLIANEW = 'REDU' ZFLIA MODNEW ;
  3601.             NZNEW = 'EXTR' ZFLIANEW 'NBZO' ;
  3602.             ZFPNW = 0. ; XVITW = 0. ;
  3603.             'REPETER' BZLIAW NZNEW ;
  3604. *             un point support par zone - 2010 kich
  3605.               ZFPNW = ZFPNW +
  3606.                   ('EXCO' ('EXTR' ZFLIANEW 'FLIA' &BZLIAW 1 1) MLDUAL_LIA 'NOID') ;
  3607.               XVITW = XVITW +
  3608.                   ('EXCO' ('EXTR' ZFLIANEW 'FLIA' &BZLIAW 1 1) MVPRIM_LIA MLPRIM_LIA 'NOID') ;
  3609.             'FIN' BZLIAW ;
  3610.             ZFP = ZFP +  ZFPNW ;
  3611.           'FINS' ;
  3612.         'FINS' ;
  3613.       'FINS' ;
  3614.  
  3615.       'SI' (NZLIA > 0) ;
  3616.         XVIT2 = 0.d0;
  3617.         'REPETER' BZLIA NZLIA ;
  3618.           XVIT2 = XVIT2 +
  3619.                   ('EXCO' ('EXTR' ZFLIA 'FLIA' &BZLIA 1 1) MVPRIM_LIA MLPRIM_LIA 'NOID') ;
  3620.         'FIN' BZLIA ;
  3621.         'SI' ('EGA' ('TYPE' XVIT2) 'CHPOINT  ') ;
  3622.           MAAUT2 = 'DIFF' ('EXTR' VITI 'MAIL') ('EXTR' XVIT2 'MAIL') ;
  3623.           VITI = ('REDU' VITI MAAUT2) + XVIT2 ; 'DETR' XVIT2 ;
  3624.           'SI' ('EGA' ('TYPE' XVITW) 'CHPOINT  ') ;
  3625.             MAAUT3 = 'DIFF' ('EXTR' VITI 'MAIL') ('EXTR' XVITW 'MAIL');
  3626.             VITI = ('REDU' VITI MAAUT3) + XVITW ; 'DETR' XVITW ;
  3627.           'FINS' ;
  3628.         'FINS' ;
  3629.         'DETR' FLIAI ;'DETR' ZFLIA ;
  3630.       'FINS' ;
  3631.     'FINS' ;
  3632. *DMODI_NB
  3633.     ZVITET = VITI - conti. 'VITESSES';
  3634.     ACCEI  = UNSURH * 2. * ZVITET - conti.'ACCELERATIONS';
  3635.     ACCEI  = 'ENLE' ACCEI 'LX' ;
  3636. *FMODI_NB
  3637.   'FINS' ;
  3638. *
  3639.   'SI' WTAB.'NVSTNL' ; 
  3640.     estim.'VITESSES_FLUIDE_0' = estim.'VITESSES_FLUIDE' ;
  3641.     estim.'VITESSES_FLUIDE' = ZVIFL  ;
  3642.     'SI' ('NON' ('EXISTE' WTAB 'QNSL')) ;
  3643.       WTAB.'QNSL' = 'MASSE' 'PRES' WTAB.'MOD_NSL' WTAB.'MAT_NSL' ;
  3644.     'FINSI' ;
  3645.     QNSL = WTAB.'QNSL' ;
  3646.     FPREFL = 'PROI' WTAB.'MOD_NSL' ZVIFL ;
  3647.     ZPREFL = 'RESOUT' QNSL FPREFL ;
  3648.     estim.'PRESSION_FLUIDE' = ZPREFL ;
  3649.     WTAB.'CONV' = ZFLCONV            ; 
  3650.     'SI' ('NON' ZFLCONV) ; 'MESS' 'Pas de convergence NVST NLIN'; 'FINSI' ;
  3651.     'QUITTER' BONOCONV ;
  3652.   'FINSI' ;
  3653. *
  3654. * Redefinir les variables associees a l'instant TEMPS0, pour soit :
  3655. * -> ZICONV = VRAI : initialiser et preparer le pas suivant
  3656. * -> ZICONV = FAUX : redemarrer le pas courant en convergence forcee
  3657.   'SI' WTAB.'FOR_MECA' ;
  3658. *
  3659.     augmult = augmult * 0.55 ;
  3660.     si (augmult < 1d-3); augmult = 1d-3; finsi;
  3661. *
  3662. *   Conserver la partie Forces nonlineaires trouvee pendant ce pas pour 
  3663. *   s'en servir pour l'initialisation du pas suivant
  3664. *   DFNL = K*DU  - DF  - residu
  3665. *   Ktot*dutot ne contient pas les forces de reactions mais les forces
  3666. *   internes dues a un champ de deformation initiales ( thermique)
  3667. *   pour les forces suiveuses on fait delta FP
  3668. *   pour etre plus precis on fait  aussi intervenir le residu
  3669. *
  3670.     XXX1 = ZRAID * ZDEPT;
  3671.     XXX3 = (ZDFINI '+' ZFSUIV) * COEPI;
  3672.     XXX4 = XXX1 - XXX3;
  3673.     XXX3 = XXX4 + RESIDC;
  3674.     XXX1 = XXX3 - FREAP ;
  3675.     XXX5 = XXX1 'ENLEVER' 'FLX ' ;
  3676.     'DETR' XXX1 ; 'DETR' XXX3 ; 'DETR' XXX4 ;
  3677. * en cas de non convergence on cumule les forces non lineaires
  3678.     'SI' WTAB.'CONV';
  3679.       ZFNONL = XXX5 ;
  3680.     'SINON';
  3681.       ZFNONL = XXX5 + ZFNONL ;
  3682.       'DETR' XXX5 ;
  3683.     'FINS';
  3684. *
  3685. *   Garder la derniere matrice KTAN calculee dans ETIQ si necessaire
  3686.     'SI' IKTAN ;
  3687.       'SI' IFEFP ;
  3688. *       'MESS' 'FEFP: Last KTAN is kept for next increment' ;
  3689.         ZLASTKTAN = ZRIKTA ;
  3690.       'SINON' ;
  3691.         'SI' (IKT_SAUV 'ET' IKT) ;
  3692. *         'MESS' 'KTAN : La matrice est conservee pour le pas suivant' ;
  3693.           ZLASTKTAN = ZRIKTA ;
  3694.         'FINS' ;
  3695.       'FINS' ;
  3696.     'FINS' ;
  3697. *
  3698.     'DETR' ZFEXT ; 'DETR' FREAP ;
  3699. *
  3700.     'SI' IGRD ;
  3701.       'FORM' GEOM2;
  3702.     'FINSI';
  3703. *
  3704.     'SI' (IGRD 'OU' IFEFPUL) ;
  3705.       'SI' (knoconv '>' 1) ;
  3706.         'DETR' GEOM1 ;
  3707.       'FINS';
  3708.     'FINS';
  3709. *
  3710. *   Utile pour les calculs d'usure
  3711.     'SI' WTAB.'CAFROTTE' ;
  3712.       ZFLX = 'EXCO' (ZDEPT '*' ZCLIM) 'FLX' 'NOID' 'FLX' ;
  3713.       WTAB.'POST_COFR'.'RIGI_UNILA' = CRR ;
  3714.       WTAB.'POST_COFR'.'MULT_SOLUT' = 'EXCO' ZDEPT 'LX' 'NOID' 'LX' ;
  3715.       WTAB.'POST_COFR'.'GLISSEMENT' = ZFLX ;
  3716.     'FINSI' ;
  3717. *
  3718.     ZXDENO  = XDENO;
  3719.     ZXDENOM = XDENOM;
  3720. *
  3721. *   Faut-il partir en convergence forcee ou ok?
  3722. *   -> si convergence forcee, recalcul du pas avec TEMPS0 = TI
  3723.     WTAB.'CONV' = ZICONV ;
  3724.     'SI' ('NON' ZICONV) ;
  3725. *
  3726. *     -----------------------------------------------------------------
  3727. *     Preparation des donnees pour le nouveau sous-pas
  3728.       TEMPS0 = TI;
  3729. *
  3730.       'SI' (IGRD 'OU' IFEFPUL) ;
  3731.         GEOM1 = GEOM2;
  3732.       'FINS';
  3733. *
  3734.       ZDEP0 = ZDETOT ;
  3735.       ZDEF0 = ZDEFF ;
  3736.       ZSIG0 = ZSIGF ;
  3737.       'SI' IPLAVI ;
  3738.         ZDEI0 = ZDEIF ;
  3739.         ZVAR0 = ZVARF ;
  3740.         'SI' ('NEG' WTAB.'MOVA' 'RIEN') ;
  3741.           ACC0 = 'EXCO' ZVAR0 WTAB.'MOVA' ;
  3742.         'FINS' ;
  3743.       'FINS';
  3744.       ZGRDU0 = ZGRDUF;
  3745. *  
  3746.       ZETAT1 = 'REDU' ZETAT2 ZMODLI ;
  3747.       'SI' ('EXIS' WTAB 'MOD_CHA') ;
  3748.         ETAZ   = 'REDU' ZETAT2 WTAB.'MOD_CHA' ;
  3749.         ZETAT1 = ZETAT1 'ET' ETAZ ;
  3750.         ETAZ   = 1. ;
  3751.       'FINS' ;
  3752. *
  3753.       'SI' ITHER ;
  3754.         ZTET1  = ZTET2 ;
  3755.         ZTEMP1 = ZTEMP2 ;
  3756.         ETT0 = 'REDU' ETT ZMODL ;
  3757.         'SI' WTAB.'POR1' ;
  3758.           MSRTT0 = 'REDU' MSRTT ZMODL ;
  3759.         'FINSI' ;
  3760.       'FINSI' ;
  3761. *
  3762.       'SI' LOGDEF;
  3763.         ZDEFOR1 = ZDEFOR2 ;
  3764.       'FINS';
  3765. *     -----------------------------------------------------------------
  3766.     'FINS' ;
  3767. *
  3768.     'SI' (ZICONV 'OU' ('NON' WTAB.'CONVERGENCE_FORCEE')) ;
  3769.       'QUITTER' BONOCONV;
  3770.     'FINS';
  3771. *
  3772. *   Convergence forcee -> on reinitialise ZCLIM
  3773.     ZCLIM = ZCLIM0 ;
  3774.   'FINSI' ;
  3775. *FOR_MECA
  3776. *
  3777. * Nombre maximum de sous-pas atteint? si oui, arret de pasapas
  3778.   WTAB.'ISOUSPAS' = WTAB.'ISOUSPAS' + 1;
  3779.   'SI' (WTAB.'ISOUSPAS' >EG WTAB.'MAXSOUSPAS');
  3780.     'ERREUR' 996 ;
  3781.   'FINS';
  3782. *
  3783. 'FIN' BONOCONV;
  3784. *----------------------------------------------------------------------*
  3785. *           Fin de la boucle de non convergence BONOCONV               *
  3786. *----------------------------------------------------------------------*
  3787. *######################################################################*
  3788. *----------------------------------------------------------------------*
  3789. *           Preparation du pas suivant TI -> TEMPS0                    *
  3790. *----------------------------------------------------------------------*
  3791. WTAB.'KNOCONV' = KNOCONV;
  3792. *
  3793. 'SI' WTAB.'FOR_MECA' ;
  3794. *
  3795.   estim.'DEPLACEMENTS' = ZDETOT ;
  3796.   estim.'REACTIONS'    = 'REAC' estim.'DEPLACEMENTS' ZCLIM ;
  3797.   estim.'CONTRAINTES'  = 'REDU' ZSIGF ZMODLI;
  3798.   estim.'DEFORMATIONS' = 'REDU' ZDEFF ZMODLI ;
  3799. *
  3800.   'SI'  IPLAVI ;
  3801.     estim.'VARIABLES_INTERNES'        = 'REDU' ZVARF ZMODLI ;
  3802.     estim.'DEFORMATIONS_INELASTIQUES' = 'REDU' ZDEIF ZMODLI ;
  3803.   'FINS';
  3804. *
  3805.   'SI'  IDYN  ;
  3806.     'SI' ('NEG' WTAB.'REAPREC' 'INCONNU');
  3807.       maireac = 'EXTR' estim.'REACTIONS' 'MAIL' ;
  3808.       'SI' ('EGA' ('TYPE' maireac) 'MAILLAGE');
  3809.         reacdif = 'REDU' WTAB.'REAPREC' maireac ;
  3810.         estim.'REACTIONS' = estim.'REACTIONS' - REACDIF;
  3811.       'FINSI' ;
  3812.     'FINSI';
  3813.     estim.'VITESSES' = VITI;
  3814. *
  3815.     'SI' IMPLP; 
  3816. *     Correction des vitesses pour avoir des vitesses relatives nulles
  3817. *     aux points qui sont en contact
  3818.       estim.'VITESSES' = estim.'VITESSES' + VADD ;
  3819.     'SINON' ;
  3820. *     Appuis unilateraux + choc elastique : essai de corriger les 
  3821. *     vitesses fournies par le schema
  3822.       'SI' ('NEG' WTAB.'ZRAIDV' 'INCONNU') ;
  3823.         estim.'VITESSES' = VITEUNIL WTAB.'ZRAIDV' WTAB.'MASSE' 
  3824.                            estim.'VITESSES' ZDEPT conti.'DEPLACEMENTS' 
  3825.                            ZDT ZSDMBR WTAB ;
  3826.         'SI' ('EXIS' WTAB 'RATE_VITEUNIL');
  3827.           'ERREUR' 'Probleme dans viteunil' ; 
  3828.         'FINSI';
  3829.       'FINSI';
  3830.     'FINS';
  3831. *
  3832.     estim.'ACCELERATIONS' = ACCEI ;
  3833.     WTAB.'FOPL'           = ZFP ;
  3834.   'FINS' ;
  3835. *
  3836.   'SI' LOGPIL ; 
  3837.     PRECED.'COEFFICIENT_DE_PILOTAGE' = PRECED.'COEFFICIENT_DE_PILOTAGE' 
  3838.           'ET' ('PROG'ETA) ;                                                                         
  3839.   'FINSI' ;
  3840. *
  3841. *----------------------------------------------------------------------
  3842. *
  3843.   WTAB.'CLIM' = ZCLIM ;
  3844.   'SI' (IGRD 'OU' IFEFPUL) ;
  3845.     WTAB.'FOR' = GEOM2 ;
  3846.   'FINS';
  3847.   WTAB.'MAT1' = 'REDU' ZMAT2F ZMODLI ;
  3848.   'SI' ('EXIS' WTAB 'MOD_CHA') ;
  3849.     MODZ = 'REDU' ZMAT22 WTAB.'MOD_CHA' ;
  3850.     WTAB.'MAT1' = WTAB.'MAT1' 'ET' MODZ ;
  3851.     MODZ = 1.D0 ;
  3852.   'FINS' ;
  3853. *
  3854. * Pour initialisation a partir du pas precedent
  3855.   WTAB.'DTPREC'    = WTAB.'DT';
  3856.   WTAB.'FNONL'     = ZFNONL ;
  3857.   WTAB.'INCREMENT' = ZINCREMENT ;
  3858.   WTAB.'RESIDU'    = RESIDC ;
  3859. *
  3860. * Contact-frottement
  3861.   'SI' ('NEG' ZRIBLO_M 'INCONNU') ;
  3862.     WTAB.'RIBLO_M' = ZRIBLO_M ;
  3863.     WTAB.'LISEA_M' = ZLISEA_M ;
  3864.   'FINSI' ;
  3865. *
  3866. * Pilotage
  3867.   WTAB.'AUTOCOEF' = ZAUTOCOEF ;
  3868.   'SI' IPILOT ;
  3869.     WTAB.'AUTOREDU' = ZAUTOREDU ;
  3870.     WTAB.'AUTORED1' = ZAUTORED1 ;
  3871.     WTAB.'NBPLAS'   = ZNBPLAS ;
  3872.   'FINSI' ;
  3873. *
  3874.   WTAB.'XDENO'             = ZXDENO ;
  3875.   WTAB.'XDENOM'            = ZXDENOM ;
  3876.   WTAB.'LASTKTAN'          = ZLASTKTAN ;
  3877.   WTAB.'NOMBRE_ITERATIONS' = ZNITE ;
  3878. *
  3879.   'SI' ITHER ;
  3880.     WTAB.'ETHER2' = 'REDU' ETT ZMODLI ;
  3881.     'SI' POR1 ;
  3882.       WTAB.'MSRTHER2' = 'REDU' MSRTT ZMODLI;
  3883.     'FINSI' ;
  3884.   'FINSI' ;
  3885. *
  3886. * Mise a jour du MODELE du PAS precedent
  3887.   WTAB.'MO_TOT_PREC' = ZMODLI ;
  3888. *
  3889. 'FINSI' ;
  3890. *FOR_MECA 
  3891. *
  3892. 'OPTI' 'PARA' FAUX;
  3893. *
  3894. IERRMEC = ('NON' WTAB.'CONV') 'ET' ('NON' WTAB.'CONVERGENCE_FORCEE') ;
  3895. *
  3896. 'FINPROC' IERRMEC ;
  3897.  
  3898.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales