* fichier :  steinb.dgibi
**************************************************
*                                                *  
*  TEST DE VALIDATION D'UNE LOI DE COMPORTEMENT  *  
*  --------------------------------------------  * 
*              DE MATERIAU                       * 
*              -----------                       *   
*                                                *  
*  LOI DE COMPORTEMENT DE:                       * 
*   STEINBERG-COCHRAN-GUINAN                     * 
*                                                *  
*  MAILLAGE:                                     *  
*   UNE BARRE DE SECTION CARREE                  * 
*   LONGUEUR L=.5 M                              *  
*   LARGEUR  l=.01 M                             * 
*                                                *  
*  CHARGEMENT:                                   * 
*   ESSAI DE TRACTION                            * 
*   DEPLACEMENTS IMPOSES                         * 
*                                                * 
**************************************************  
* 
* 
*  
opti echo 0 dime 3 elem cub8 ;
* ajout de option epsilon lineaire pour la precision des test!
OPTION epsilon lineaire; 
*  
*  
*  Maillage 
*  
l1 = .5 ; 
l2 = .01 ; 
l3 = .01 ;  
n1 = 50 ; 
n2 = 1 ; 
n3 = 1 ; 
p1 = 0. 0. 0. ; 
p2 = l1 0. 0. ; 
p3 = l1 l2 0. ; 
p4 = 0. l2 0. ; 
p5 = 0. 0. l3 ; 
p6 = l1 0. l3 ; 
p7 = l1 l2 l3 ; 
p8 = 0. l2 l3 ; 
d1 = d p1 n3 p5 ; 
d2 = d p5 n2 p8 ; 
d3 = d p8 n3 p4 ; 
d4 = d p4 n2 p1 ; 
sur1 = 'DALL' d1 d2 d3 d4 ; 
d5 = d p2 n3 p6 ; 
d6 = d p6 n2 p7 ;  
d7 = d p7 n3 p3 ; 
d8 = d p3 n2 p2 ; 
sur2 = 'DALL' d5 d6 d7 d8 ; 
vol1 = sur1 'VOLU' n1 sur2 ;    
*  
*    Modele de calcul 
* 
INCO_GT = 'EVOL' 'MANU' 'T' ( 'PROG' -273. 27. 1027. ) 
          'GT' ( 'PROG' 5133.E6 0. -17110.E6 ) ;   
mod1 = 'MODE' vol1 'MECANIQUE' 'ELASTIQUE' 'ISOTROPE' 
       'PLASTIQUE' 'STEINBERG' ; 
mat1 = 'MATE' mod1 'YOUN' 2.1E11 'NU' 0.3 
       'RHO' 7.8E3 'ALPHA' 0. 'SIGY' 410.E6  
       'BETA' 125 'N' .1 'EPSI' 0.  
       'GP' 1.79 'GT' INCO_GT 'YMAX' 680.E6 'TMO' 1220. 
       'MU' 1.95 ;
*  
*  Modele thermique  
*  
mod2 = 'MODE' vol1 'THERMIQUE' 'ISOTROPE' ; 
mat2 = 'MATE' mod2 'K' 50. 'C' 140. 'RHO' 7.8E3 ;  
*  
mod3 = 'MODE' ( 'ENVE' vol1 ) thermique 'CONVECTION' 
        ; 
mat3 = 'MATE' mod3 'H' 750. ; 
*  
*  Modele total 
*  
mod0 = mod1 'ET' mod2 'ET' mod3 ; 
mat0 = mat1 'ET' mat2 'ET' mat3 ;  
*
*    Conditions aux limites 
*  
cl1 = 'BLOQ' sur1 'UX' ; 
cl2 = 'BLOQ' sur2 'UX' ; 
cl3 = 'BLOQ' p1 'UX' 'UY' 'UZ' ;  
cl4 = 'BLOQ' p5 'UX' 'UY' ; 
cl0 = cl1 'ET' cl2 'ET' cl3 'ET' cl4 ;  
* 
*  Chargement  
*   
dep1 = 'DEPI' cl2 .1 ; 
tex0 = 'MANU' 'CHPO' ( 'ENVE' vol1 ) 1 'T' 500. ;  
ev0 = 'EVOL' 'MANU' temps ( 'PROG' 0. 1000. ) y ( 'PROG' 0. 1. ) ; 
cha1 = 'CHAR' 'DIMP' dep1 ev0 ;  
cha2 = 'CHAR' 'TECO' tex0 ev0 ; 
cha0 = cha1 'ET' cha2 ;  
* 
*  Temps du calcul  
*  
dt0 = 10. ; 
tfin0 = 100. ;  
*   
*   
*  Resolution par PASAPAS  
*---------------------------------------
*    
ta1 = 'TABLE' ; 
ta1 .'MODELE' = mod0 ; 
ta1 .'CARACTERISTIQUES' = mat0 ; 
ta1 .'BLOCAGES_MECANIQUES' = cl0 ; 
ta1 .'CHARGEMENT' = cha0 ; 
ta1 .'PROCEDURE_THERMIQUE' = LINEAIRE ; 
ta1 .'TEMPS_CALCULES' = 'PROG' 0. 'PAS' dt0 tfin0 ; 
*  
TMASAU=table;
ta1 . 'MES_SAUVEGARDES'=TMASAU;
TMASAU .'DEFTO'=VRAI;
TMASAU .'DEFIN'=VRAI;
PASAPAS ta1 ; 
*  
*  
*  Post-traitement 
*-----------------------------------------
* 
*  Quelques donnees du calcul 
*  
youn1 = 'EXTR' ( 'EXCO' mat0 'YOUN' 'SCAL' ) 'SCAL' 1 1 1 ; 
xnu0 = 'EXTR' ( 'EXCO' mat0 'NU' 'SCAL' ) 'SCAL' 1 1 1 ;
sig0 = 'EXTR' ( 'EXCO' mat0 'SIGY' 'SCAL' ) 'SCAL' 1 1 1 ; 
m0 = 'EXTR' ( 'EXCO' mat0 'N' 'SCAL' ) 'SCAL' 1 1 1 ; 
beta0 = 'EXTR' ( 'EXCO' mat0 'BETA' 'SCAL' ) 'SCAL' 1 1 1 ; 
G0 = youn1 / ( 2. * ( 1. + xnu0 ) ) ; 
GP0 = 'EXTR' ( 'EXCO' mat0 'GP' 'SCAL' ) 'SCAL' 1 1 1 ;  
eps_el0 = sig0 / youn1 ;
*  
*  Post-traitement 
*-----------------------------------------
* 
sigx1 = 'PROG' 0. ; 
sigx2 = 'PROG' 0. ; 
depx1 = 'PROG' 0. ;  
depx2 = 'PROG' 0. ; 
sig_1 = 'PROG' 0. ; 
def_1 = 'PROG' 0. ;   
eps_th1 = 'PROG' 0. ; 
con_th1 = 'PROG' 0. ; 
err_e1 = 'PROG' 0. ; 
err_s1 = 'PROG' 0. ; 
tp_max1 = 'PROG' 0. ; 
tp_min1 = 'PROG' 0. ;  
A1 = 2. * ( 1. + xnu0 ) ; 
B1 = A1 / ( 3. * ( 1. - ( 2. * xnu0 ) ) ) ; 
dim0 = 'DIME' ta1 .'TEMPS' - 1 ; 
mod0 = ta1 .'MODELE' ; 
'REPETER' bloc1 dim0 ; 
  i1 = &bloc1 ; 
  mat0 = ta1 .'CARACTERISTIQUES' ;
  depl0 = ta1 .'DEPLACEMENTS'.i1 ; 
  sigm0 = ta1 .'CONTRAINTES'.i1 ; 
  depeq0 = ta1 .'VARIABLES_INTERNES'.i1 ;  
  def0 = ta1 .'DEFORMATIONS_INELASTIQUES'.i1 ;  
  cht0 = ta1 .'TEMPERATURES'.i1 ; 
  GT0 = 'EXCO' ( 'VARI' 'NUAG' mod1 mat1 cht0 ) 'GT' 'SCAL' ;  
  GT0 = 'EXTR' GT0 'SCAL' 1 1 1 ;  
  t1 = 'EXTR' cht0 'T' ( 'NOEUD' 1 ) ; 
*  'MESS' ' Temperature au temps ' i1 ' = ' t1 ; 
  tp_max0 = 'MAXI' cht0 ; 
  tp_min0 = 'MINI' cht0 ;  
  tp_max1 = tp_max1 'ET' ( 'PROG' tp_max0 ) ; 
  tp_min1 = tp_min1 'ET' ( 'PROG' tp_min0 ) ;   
  sigx0 = 'MAXI' ( 'VMIS' mod1 sigm0 ) ;     
  sigx00 = 'MINI' ( 'VMIS' mod1 sigm0 ) ; 
  depx0 = 'MAXI' ( 'EXCO' depeq0 'EPSE' 'SCAL' ) ;   
  depx00 = 'MINI'( 'EXCO' depeq0 'EPSE' 'SCAL' ) ;  
  sig_00 = 'EXCO' sigm0 'SMXX' 'SCAL' ; 
  sig_0 = 'EXTR' sig_00 'SCAL' 1 1 1 ;    
  def_00 = 'EXCO' def0 'EIXX' 'SCAL' ; 
  def_0 = 'EXTR' def_00 'SCAL' 1 1 1 ; 
  sigx1 = sigx1 'ET' ( 'PROG' sigx0 ) ; 
  sigx2 = sigx2 'ET' ( 'PROG' sigx00 ) ;  
  depx1 = depx1 'ET' ( 'PROG' depx0 ) ;  
  depx2 = depx2 'ET' ( 'PROG' depx00 ) ; 
  sig_1 = sig_1 'ET' ( 'PROG' sig_0 ) ; 
  def_1 = def_1 'ET' ( 'PROG' def_0 ) ;  
*  
*  Eta  
*  
  eps00 = 'EPSI' mod1 depl0 ; 
  epsx = 'EXTR' ( 'EXCO' eps00 'EPXX' 'SCAL' ) 'SCAL' 1 1 1 ; 
  epsy = 'EXTR' ( 'EXCO' eps00 'EPYY' 'SCAL' ) 'SCAL' 1 1 1 ; 
  epsz = 'EXTR' ( 'EXCO' eps00 'EPZZ' 'SCAL' ) 'SCAL' 1 1 1 ; 
  treps0 = epsx + epsy + epsz ; 
  eta = 1. / ( 1. + treps0 ) ; 
* 
*  Courbe theorique 
* 
eps0 = epsx ;  
sig1 = sig0 * ( 1. + ( ( GP0 / G0 ) * ( -1. * sig_0 / 3. ) ) + 
( GT0 / G0 ) ) ; 
pente0 = GP0 / ( 3. * ( eta ** ( 1. / 3. ) ) ) ;  
pente1 = G0 + GT0 ; 
*
'SI' ( sigx0 > sig1 ) ; 
x0 = 0. ; 
x2 = x0 ; 
tr1 = treps0 - ( A1 * ( eps0 - x2 ) / ( 3. * B1 ) ) ;
tr1 = tr1 / 3. ;
y2 = ( sig0 * ( ( 1 + ( beta0 * ( x2 - tr1 ) ) ) ** m0 ) ) ; 
y2 = y2 + ( A1 * G0 * x2 ) ;  
y2 = y2 - ( A1 * G0 * eps0 )  ; 
*
x1 = eps0 ; 
x2 = x1 ;   
tr1 = treps0 - ( A1 * ( eps0 - x2 ) / ( 3. * B1 ) ) ;
tr1 = tr1 / 3. ;
y2 = ( sig0 * ( ( 1 + ( beta0 * ( x2 - tr1 ) ) ) ** m0 ) ) ; 
y2 = y2 + ( A1 * G0 * x2 ) ; 
y2 = y2 - ( A1 * G0 * eps0 ) ; 
*  
'REPETER' bloc2 ; 
*  
x2 = ( x0 + x1 ) / 2. ; 
tr1 = treps0 - ( A1 * ( eps0 - x2 ) / ( 3. * B1 ) ) ;
tr1 = tr1 / 3. ;
y2 = ( sig0 * ( ( 1 + ( beta0 * ( x2 - tr1 ) ) ) ** m0 ) ) ; 
y2 = y2 + ( A1 * G0 * x2 ) ; 
y2 = y2 - ( A1 * G0 * eps0 ) ;  
*  
'SI' ( y2 > 1.E-4 ) ; 
  x1 = x2 ; 
'SINON' ; 
  'SI' ( y2 <  -1.E-4 ) ; 
     x0 = x2 ; 
  'SINON' ; 
     eps_th0 = x2 ; 
     con_th0 = ( A1 * pente1 * ( eps0 - x2 ) ) ; 
     con_th0 = con_th0 / ( 1. + ( A1 * pente0 * ( eps0 - x2 ) ) ) ;  
    'QUITTER' bloc2 ;  
  'FINSI' ; 
'FINSI' ; 
*  
'FIN' bloc2 ;
*  
'SINON' ; 
eps_th0 = 0. ; 
con_th0 = ( A1 * pente1 * eps0 ) ;
con_th0 = con_th0 / ( 1. + ( A1 * pente0 * eps0 ) ) ;  
'FINSI' ; 
*    
*  'MESS' ' Fin du calcul de restitution pour le pas' ta1 .'TEMPS'.i1 ;  
  eps_th1 = eps_th1 'ET' ( 'PROG' eps_th0 ) ; 
  con_th1 = con_th1 'ET' ( 'PROG' con_th0 ) ;   
  max0 = ( 'ABS' ( con_th0 - sig_0 ) + con_th0 + sig_0 ) ; 
  max0 = 0.5 * max0 ; 
  max1 = ( 'ABS' ( eps_th0 - def_0 ) + eps_th0 + def_0 ) ; 
  max1 = 0.5 * max1 ;  
  'SI' ( ( 'ABS' max0 ) > 1.E-6 ) ; 
     err_s0 = ( 'ABS' ( con_th0 - sig_0 ) ) / max0 ;  
     err_s0 = err_s0 * 100. ;  
  'SINON' ;  
     err_s0 = 0. ; 
  'FINSI' ; 
  'SI' ( ( 'ABS' max1 ) > ( 1.E-1 * eps_el0 ) ) ; 
     err_e0 = ( 'ABS' ( eps_th0 - def_0 ) ) / max1 ;  
     err_e0 = err_e0 * 100. ;  
  'SINON' ;  
     err_e0 = 0. ; 
  'FINSI' ; 
  'SI' ( ( err_e0 '>' 2.E-2 ) 'OU' ( err_s0 '>' 2.E-2 ) ) ; 
      'ERRE' 5 ; 
  'FINSI' ; 
  err_s1 = err_s1 'ET' ( 'PROG' err_s0 ) ; 
  err_e1 = err_e1 'ET' ( 'PROG' err_e0 ) ;  
*  
*'MESS' ' Deformation plastique theorique =' eps_th0 ; 
*'MESS' ' Deformation plastique calculee  =' def_0 ;  
* 
*   
'FIN' bloc1 ; 
*  
*  Trace  
*  
trac1 = 'EVOL' 'ROUG' 'MANU' 'epsp' depx1 'sigm' sigx1 ; 
*'DESSIN' trac1 'TITR' 'Contraintes - Deformation plastique' ; 
trac2 = 'EVOL' 'ROUG' 'MANU' 'epsp' depx2 'sigm' sigx2 ;  
TAB = 'TABLE' ; 
TAB.1 = 'MARQ CROI ' ;
TAB.2 = 'MARQ PLUS ' ;
TAB.'TITRE' = 'TABLE' ;
TAB.'TITRE'. 1 = MOT 'MAXI' ;
TAB.'TITRE'. 2 = MOT 'MINI' ;
*'DESSIN' ( trac1 'ET' trac2 )
*'TITR' 'Comparaison mini - maxi de la contrainte equivalente '
*'LEGE' TAB ;
trac3 = 'EVOL' 'ROUG' 'MANU' 'def_pl_xx' def_1 'tp_max' tp_max1 ; 
trac31 = 'EVOL' 'ROUG' 'MANU' 'def_pl_xx' def_1 'tp_max' tp_min1 ;  
TAB = 'TABLE' ;
TAB.1 = 'MARQ CROI ' ;
TAB.2 = 'MARQ PLUS ' ;
TAB.'TITRE' = 'TABLE' ;
TAB.'TITRE'. 1 = MOT 'MAXI' ;
TAB.'TITRE'. 2 = MOT 'MINI' ;
*'DESSIN' ( trac3 'ET' trac31 )
*'TITR' 'Comparaison mini - maxi de la temperature '
*'LEGE' TAB ;
trac4 = 'EVOL' 'ROUG' 'MANU' 'def_pl_xx' def_1 'Sig_xx' sig_1 ;  
trac5 = 'EVOL' 'ROUG' 'MANU' 'eps_pl' eps_th1 'Sig_xx' con_th1 ; 
TAB= 'TABLE' ; 
TAB.1 = 'MARQ CROI ' ;
TAB.2 = 'MARQ PLUS ' ;
TAB.'TITRE' = 'TABLE' ;
TAB.'TITRE'. 1 = MOT 'Calcul' ;
TAB.'TITRE'. 2 = MOT 'Theorie' ;
*'DESSIN' ( trac4 'ET' trac5 )
*'TITR' ' Courbe theorique et calculee'
*'LEGE' TAB ;
trac6 = 'EVOL' 'ROUG' 'MANU' 'epsp' def_1 'erreur %' err_s1 ; 
trac7 = 'EVOL' 'ROUG' 'MANU' 'epsp' def_1 'erreur %' err_e1 ; 
TAB = 'TABLE' ; 
TAB.1 = 'MARQ PLUS ' ;
TAB.2 = 'MARQ CROI ' ;
TAB.'TITRE' = 'TABLE' ;
TAB.'TITRE'. 1 = MOT 'CONTRAINTES' ;
TAB.'TITRE'. 2 = MOT 'DEFORMATIONS' ;
*'DESSIN' ( trac6 'ET' trac7 )
*'TITR' 'Erreur sur les contraintes et la déformation plastiques '
*'LEGE' TAB ;
*  
*'MESS' ' deformation plastique theorique ' eps_th0 ; 
*'MESS' ' deformation plastique calculee  ' def_0 ; 
*'MESS' ' contrainte theorique            ' con_th0 ; 
*'MESS' ' contrainte calculee             ' sig_0 ;  
*'MESS' ' deformation totale finale       ' eps0 ; 
*  
*    
'FIN' ;