Procédures | Cast3M

* @SOLVMEC  PROCEDUR  FD218221  21/05/07    21:15:01     11001          
************************************************************************
*       Procedure minimaliste de solveur implicite en mecanique        *
*        - comportements lineaires et non lineaires (appel a COMP)     *
*        - chargements en forces ou deplacements imposes               *
*        - hypothese petits deplacements                               *
*                 et grands deplacements en lagrangien reactualise     *
*        - acceleration de convergence (appel a ACT3)                  *
************************************************************************
 
DEBP @SOLVMEC t*'TABLE' ;
 
************************************************************************
*                 Partie 1 : initialisations diverses                  *
************************************************************************
* nombre max d'iteration de minimisation du residu (sur un pas de temps)
  nmaxit = 1000 ;
 
* objets fournis par l'utilisateur
  mod = t . 'MODELE' ;
  mat = t . 'CARACTERISTIQUES' ;
  blo = t . 'BLOCAGES_MECANIQUES' ;
  cha = t . 'CHARGEMENT' ;
  ltc = t . 'TEMPS_CALCULES' ;
 
* lecture des options
  igrdepl = FAUX ;
  SI (EXIS t 'GRANDS_DEPLACEMENTS') ;
    igrdepl = t . 'GRANDS_DEPLACEMENTS' ;
  FINSI ;
  iaccconv = FAUX ;
  SI (EXIS t 'ACCELERATION_CONVERGENCE') ;
    iaccconv = t . 'ACCELERATION_CONVERGENCE' ;
  FINSI ;
 
* initialisations des champs solution a tps = 0
  tps0 = 0. ;
  u0   = CHAN 'CHPO' (ZERO mod 'DEPLACEM') mod ;
  fr0  = REAC u0 blo ;
  eps0 = ZERO mod 'DEFORMAT' ;
  epi0 = ZERO mod 'DEFINELA' ;
  sig0 = ZERO mod 'CONTRAIN' ;
  fi0  = BSIG sig0 mod ;
  vi0  = ZERO mod 'VARINTER' ;
  SI (EXIS cha 'MECA') ;
    fex0 = TIRE cha tps0 'MECA' ;
  SINON ;
    fex0 = VIDE 'CHPOINT'/'DISCRET' ;
  FINSI ;
  SI (EXIS cha 'DIMP') ;
    flx0 = TIRE cha tps0 'DIMP' ;
  SINON ;
    flx0 = VIDE 'CHPOINT'/'DISCRET' ;
  FINSI ;
 
* initialisation des indices de sortie de la table
  t . 'TEMPS'                  = TABL ;
  t . 'TEMPS' . 0              = tps0 ;
  t . 'DEPLACEMENTS'           = TABL ;
  t . 'DEPLACEMENTS' . 0       = u0 ;
  t . 'REACTIONS'              = TABL ;
  t . 'REACTIONS' . 0          = fr0 ;
  t . 'CONTRAINTES'            = TABL ;
  t . 'CONTRAINTES' . 0        = sig0 ;
  t . 'VARIABLES_INTERNES'     = TABL ;
  t . 'VARIABLES_INTERNES' . 0 = vi0 ;
 
* [grands deplacements] configuration initiale
  SI igrdepl ;
    conf_ini = FORM ;
  FINSI ;
 
* [petits deplacements] calcul de la rigidite elastique initiale
  SI (NON igrdepl) ;
    k = RIGI mod mat ;
    k = k ET blo ;
  FINSI ;
 
 
************************************************************************
*          Partie 2 : debut de la boucle sur les pas de temps          *
************************************************************************
  ntc = DIME ltc ;
  i = 0 ;
  REPE btemps ntc ;
    iconv = FAUX ;
*   instant de calcul, fin du pas de temps
    tps1 = EXTR ltc &btemps ;
    SI ((ABS tps1) < 1.E-15) ;
      ITER btemps ;
    FINSI ;
 
*   interpolation du second membre a la fin du pas de temps
*   forces imposees (composantes FX, FY, FZ, ...)
    SI (EXIS cha 'MECA') ;
      fex1 = TIRE cha tps1 'MECA' ;
    SINON ;
      fex1 = VIDE 'CHPOINT'/'DISCRET' ;
    FINSI ;
*   deplacements imposes (composante FLX)
    SI (EXIS cha 'DIMP') ;
      flx1 = TIRE cha tps1 'DIMP' ;
    SINON ;
      flx1 = VIDE 'CHPOINT'/'DISCRET' ;
    FINSI ;
 
*   initialisation du residu (avec les FLX des deplacements imposes)
*   et norme de convergence
    resi0 = fex1 + fr0 - fi0 + flx1 - flx0 ;
    fref = MAXI 'ABS' (fex1 + fr0) ;
    SI (fref < 1.E-12) ;
      fref = 1. ;
    FINSI ;
 
*   [acceleration] initialisation de la correction
    SI iaccconv ;
      cor1 = 0. ;
    FINSI ;
 
*   affichages
    MESS '----- Resolution a l''instant :' tps1 ' (' &btemps ' /' ' ' ntc ') -----' ;
    MESS ' Effort de reference =' fref ;
    MESS ' Iter.    Residu' ;
 
*   [grands deplacements] passage sur la configuration debut de pas et
*                         calcul de la rigidite elastique
    SI igrdepl ;
      conf_0 = FORM u0 ;
      k = (RIGI mod mat) ET (KSIG mod sig0) ;
      k = k ET blo ;
    FINSI ;
 
 
************************************************************************
*       Partie 3 : debut de la boucle de minimisation du residu        *
************************************************************************
    u00 = u0 * 1. ;
    REPE bminim nmaxit ;
*     correction de l'increment de deplacements correspondant au residu
*     et a la raideur elastique
      ddu = RESO k resi0 ;
 
*     deplacement a la fin du pas de temps et increment
      u1 = u0 + ddu ;
      du = u1 - u00 ;
 
*     reactions a la fin du pas de temps
      fr1 = REAC blo u1 ;
 
*     increment et deformation totale a la fin du pas de temps
*     [petits deplacements] on utilise la deformation lineaire
*     [grands deplacements] on utilise la deformation quadratique
      SI igrdepl ;
        deps = EPSI mod du 'QUAD' ;
      SINON ;
        deps = EPSI mod du 'LINE' ;
      FINSI ;
      eps1 = eps0 + deps ;
 
*     application de la loi de comportment pour l'increment de deformation sur le pas de temps
      etat_0 = (ADET 'NOUV' mod 'TEMP' tps0) ET eps0 ET mat ET sig0 ET vi0 ET epi0 ;
      etat_1 = (ADET 'NOUV' mod 'TEMP' tps1) ET eps1 ET mat ;
      ch1 = COMP mod etat_0 etat_1 ;
 
*     contraintes, deformations, variables internes a la fin du pas de temps
      sig1 = CHAN (EXCO ch1 (EXTR mod 'CONT') 'NOID') 'TYPE' 'CONTRAINTES' ;
      eps1 = CHAN (EXCO ch1 (EXTR mod 'DEFO') 'NOID') 'TYPE' 'DEFORMATIONS' ;
      epi1 = CHAN (EXCO ch1 (EXTR mod 'DEIN') 'NOID') 'TYPE' 'DEFORMATIONS INELASTIQUES' ;
      vi1  = CHAN (EXCO ch1 (EXTR mod 'VARI') 'NOID') 'TYPE' 'VARIABLES INTERNES' ;
 
*     [grands deplacements] transport des contraintes sur la configuration fin de pas
      SI igrdepl ;
        sig1 = PICA mod sig1 du ;
      FINSI ;
 
*     forces interieures a la fin du pas de temps
*     [grands deplacements] l'integration est faite sur la configuration fin de pas
      SI igrdepl ;
        conf_1 = FORM du ;
      FINSI ;
      fi1 = BSIG mod sig1 ;
*     [grands deplacements] retour a la configuration debut de pas
      SI igrdepl ;
        FORM conf_0 ;
      FINSI ;
 
*     nouveau residu et sa norme
      resi1 = fex1 + fr1 - fi1 ;
      crit1 = (MAXI 'ABS' resi1) / fref ;
 
*     affichages
      logmess = FAUX ;
      SI ((&bminim < 11) OU ((&bminim < 101)  ET (EGA (@MOD &bminim 10)  0))
                         OU ((&bminim < 1001) ET (EGA (@MOD &bminim 100) 0))) ;
        logmess = VRAI ;
        MESS &bminim'/'2 crit1'/'10 ;
      FINSI ;
 
*     test de convergence
      SI (crit1 < 1.E-4) ;
        iconv = VRAI ;
        SI (NON logmess) ;
          MESS &bminim'/'2 crit1'/'10 ;
        FINSI ;
        MESS ' Convergence atteinte !' ;
        MESS ' ' ;
        QUIT bminim ;
      FINSI ;
 
*     [acceleration de convergence] a partir des 4 derniers residus et
*                                   seulement toutes les 4 iterations
      SI iaccconv ;
        cor1p  = cor1 ;
        cor1   = 0. ;
        acfp1  = resi1 ;
        acfep1 = acfp1 - cor1p ;
        SI ((&bminim MULT 4) ET (&bminim > 3)) ;
          MESS ' Acceleration' ;
          cor1  = ACT3 acfep3 acfep2 acfep1 acfp4 acfp3 acfp2 acfp1 ;
          resi1 = resi1 - cor1 ;
        FINSI ;
        acfp4  = acfp3 ;
        acfp3  = acfp2 ;
        acfp2  = acfp1 ;
        acfep3 = acfep2 ;
        acfep2 = acfep1 ;
      FINSI ;
 
*     preparation de la prochaine iteration de minimisation
      resi0 = resi1 ;
      u0 = u1 ;
    FIN bminim ;
************************************************************************
 
*   [grands deplacements] remise en configuration initiale
    SI igrdepl ;
      FORM conf_ini ;
    FINSI ;
 
*   si convergence : preparation du prochain pas de temps
    SI iconv ;
      tps0 = tps1 ;
      u1 = EXCO u1 (EXTR mod 'DEPL') ;
      u0 = u1 ;
      fr0 = fr1 ;
      fi0 = fi1 ;
      sig0 = sig1 ;
      flx0 = flx1 ;
      eps0 = eps1 ;
      epi0 = epi1 ;
      vi0  = vi1 ;
*     remplissage de la table avec les resultats
      i = i + 1 ;
      t . 'TEMPS'              . i = tps1 ;
      t . 'DEPLACEMENTS'       . i = u1 ;
      t . 'REACTIONS'          . i = fr1 ;
      t . 'CONTRAINTES'        . i = sig1 ;
      t . 'VARIABLES_INTERNES' . i = vi1 ;
*   si non convergence : sortie en erreur
    SINON ;
      MESS '***** NOMBRE MAXI D''ITERATIONS ATTEINT *****' ;
      MESS '*****          NON CONVERGENCE          ***** !' ;
      ERRE 'Echec du calcul, on quitte la procedure' ;
    FINSI ;
  FIN btemps ;
************************************************************************
 
FINP ;