Download psury.dgibi

Back to the list

* fichier :  psury.dgibi
*****************************************************
*                                                   * 
*        TEST DE VALIDATION DE LA LOI               * 
*        D'ENDOMMAGEMENT TRIAXIAL P/Y               *  
*                                                   *  
*  TEST:                                            * 
*    UN BARREAU EST CHARGE EN TRACTION              * 
*    LE CHARGEMENT EST DES DEPLACEMENTS IMPOSES     *  
*    CALCUL EN MASSIF 3D                            *  
*                                                   *  
*  MAILLAGE:                                        * 
*    UN BARREAU DE SECTION CARREE DE .5 M SUR .01 M * 
*    50 ELEMENTS                                    *  
*                                                   * 
*****************************************************
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 ;  
*  
*  Deplacement maximal 
* 
depmax = .1 ;
*   
*  Courbe d'ecrouissage 
*
beta0 = 125. ; 
sig0 = 410.E6 ; 
xm = .1 ;       
youn0 = 2.1E11 ; 
epsp0 = 'PROG' 0. ; 
deb0 = 'PROG' 0. ;  
epsmax0 = depmax / l1 ;  
pas0 = epsmax0 / 50 ; 
epsp1 = 'PROG' 0. 'PAS' pas0 epsmax0 ;  
dim00 = 'DIME' epsp1 ;  
un = 'PROG' dim00 '*' 1. ; 
Y0 = sig0 ;  
epsp2 = beta0 * epsp1 ; 
epsp3 = epsp2 + un ; 
epsp31 = log epsp3 ; 
epsp32 = epsp31 * xm ;   
epsp4 = exp epsp32 ; 
Y1 = Y0 * epsp4 ; 
*  
epsp1 = epsp0 'ET' epsp1 ;  
Y1 = deb0 'ET' Y1 ;   
epsp2 = epsp1 + ( Y1 / youn0 ) ;   
*  
*    Courbe d'evolution de l'endommagement
*
n10 = 5. ;
alp0 = 'PROG' 0. 'PAS' .01 1. ;
*dd10 = 'SIN' ( 90. * alp0 ) ;
dd10 = ( alp0 ) ** ( 1. / n10 ) ;
evol0 = 'EVOL' 'MANU' 'defo' alp0 'sigm' dd10 ;
* 
*    Modele de calcul 
* 
trac0 = 'EVOL' 'MANU' 'defo' epsp2 'sigm' Y1 ; 
endo0 = 'EVOL' 'MANU' 'defo' ( 'PROG' .05 .1 .25 .5 ) 
        'sigm' ( 'PROG' .6 .3 .15 .1 ) ;  
mod0 = 'MODE' vol1 'MECANIQUE' 'ELASTIQUE' 'ISOTROPE' 
       'PLASTIQUE_ENDOM' 'PSURY' ;  
mat0 = 'MATE' mod0 'YOUN' youn0 'NU' 0.3 'RHO' 7800.     
       'TRAC' trac0 'EVOL' endo0 'COMP' evol0 ; 
*  
*    Conditions aux limites 
*  
cl1 = 'BLOQ' sur1 'UX' ; 
cl2 = 'BLOQ' sur2 'UX' ; 
cl3 = 'BLOQ' p1 'UX' 'UY' 'UZ' ;  
cl4 = 'BLOQ' p5 'UX' 'UY' ; 
cl5 = 'BLOQ' p4 'UX' 'UZ' ; 
cl0 = cl1 'ET' cl2 'ET' cl3 'ET' cl4 'ET' cl5 ;  
*  
*    Chargement 
*  
dep1 = 'DEPI' cl2 depmax ; 
ev0 = 'EVOL' 'MANU' temps ( 'PROG' 0. 10. ) y ( 'PROG' 0. 1. ) ; 
cha0 = 'CHAR' 'DIMP' dep1 ev0 ; 
*  
*  Temps du calcul  
*  
dt0 = 1. ; 
tfin0 = 10. ;  
*   
*   
*  Resolution par PASAPAS  
*---------------------------------------
*    
ta1 = 'TABLE' ; 
ta1 .'MODELE' = mod0 ; 
ta1 .'CARACTERISTIQUES' = mat0 ; 
ta1 .'BLOCAGES_MECANIQUES' = cl0 ; 
ta1 .'CHARGEMENT' = cha0 ; 
ta1 .'TEMPS_CALCULES' = 'PROG' 0. 'PAS' dt0 tfin0 ;  
*
TMASAU=table;
ta1 . 'MES_SAUVEGARDES'=TMASAU;
TMASAU .'DEFTO'=VRAI;
TMASAU .'DEFIN'=VRAI;  
PASAPAS ta1 ; 
*  
*  
*  Post-traitement 
*-----------------------------------------
* 
err = faux ; 
coeff = 1. ; 
rho_thf = 0. ; 
epsp1 = 'PROG' 0. ; 
epst_x1 = 'PROG' 0. ; 
def_1 = 'PROG' 0. ; 
sig1 = 'PROG' 0. ; 
sig_x1 = 'PROG' 0. ; 
eps_th1 = 'PROG' 0. ; 
con_th1 = 'PROG' 0. ; 
err_e1 = 'PROG' 0. ; 
err_s1 = 'PROG' 0. ; 
err_g1 = 'PROG' 0. ; 
g1 = 'PROG' 1. ;
g_th1 = 'PROG' 1. ;   
eps_th0 = 0. ; 
dim0 = 'DIME' ta1 .'TEMPS' - 1 ; 
mod0 = ta1 .'MODELE' ;   
'REPETER' bloc1 dim0 ; 
  i1 = &bloc1 ; 
*  
  mat0 = ta1 .'CARACTERISTIQUES' ; 
  youn0 = 'MAXI' ( 'EXCO' mat0 'YOUN' ) ;  
  rho0 = 'MAXI' ( 'EXCO' mat0 'RHO' ) ; 
  depl0 = ta1 .'DEPLACEMENTS'.i1 ; 
  sigm0 = ta1 .'CONTRAINTES'.i1 ; 
  depeq0 = ta1 .'VARIABLES_INTERNES'.i1 ; 
  def0 = ta1 .'DEFORMATIONS_INELASTIQUES'.i1 ;  
  eps0 = 'EPSI' mod0 depl0 ; 
* 
  def_0 = 'MAXI' ( 'EXCO' def0 'EIXX' 'SCAL' ) ; 
  g0 = 'MAXI' ( 'EXCO' depeq0 'VHWG' 'SCAL' ) ;
  epsp0 = 'MAXI' ( 'EXCO' depeq0 'EPSE' 'SCAL' ) ; 
  rho1 = 'MAXI' ( 'EXCO' depeq0 'VHHO' 'SCAL' ) ; 
  rho_f = 'MAXI' ( 'EXCO' depeq0 'VHOF' 'SCAL' ) ; 
  sig_0 = 'MAXI' ( 'VMIS' mod0 sigm0 ) ; 
  sig_x0 = 'MAXI' ( 'EXCO' sigm0 'SMXX' 'SCAL' ) ; 
  epst_x0 = 'MAXI' ( 'EXCO' eps0 'EPXX' 'SCAL' ) ; 
  epst_y0 = 'MAXI' ( 'EXCO' eps0 'EPYY' 'SCAL' ) ; 
  epst_z0 = 'MAXI' ( 'EXCO' eps0 'EPZZ' 'SCAL' ) ; 
  tr_eps0 = epst_x0 + epst_y0 + epst_z0 ; 
*  
*  
*  Endommagement theorique 
*--------------------------------------
*
*  Calcul de la densite rho
*
  rho = rho0 / ( 1. + tr_eps0 ) ;
*
*  deformation plastique de debut d'endommagement  
*  
absc0 = 'EXTR' endo0 'ABSC' ; 
ordo0 = 'EXTR' endo0 'ORDO' ; 
'REPETER' bloc00 ( 'DIME' absc0 ) ; 
  i00 = &bloc00 ; 
  ab1 = 'EXTR' absc0 i00 ; 
  or1 = 'EXTR' ordo0 i00 ; 
  'SI' ( i00 'EGA' ( 'DIME' absc0 ) ) ; 
     eps1 = ab1 ; 
  'SINON' ; 
     ab2 = 'EXTR' absc0 ( i00 + 1 ) ; 
     or2 = 'EXTR' ordo0 ( i00 + 1 ) ; 
     'SI' ( ( or1 '>EG' ( 1. / 3. ) ) 'ET' 
          ( or2 '<' ( 1. / 3. ) ) ) ;  
         eps1 = ( ab2 - ab1 ) / ( or2 - or1 ) ; 
         eps1 = eps1 * ( ( 1. / 3. ) - or1 ) ; 
         eps1 = eps1 + ab1 ; 
         'QUITTER' bloc00 ; 
     'FINSI' ; 
   'FINSI' ; 
'FIN' bloc00 ; 
* 
*        
*  Courbe theorique 
*--------------------------------------  
* 
m0 = xm ;  
'SI' ( ( youn0 * epst_x0 ) '>EG' sig0 ) ; 
x0 = 0. ; 
x2 = x0 ;   
y2 = ( sig0 * ( ( 1 + ( beta0 * x2 ) ) ** m0 ) ) ;
y2 = y2 + ( youn0 * x2 ) ; 
y2 = y2 - ( youn0 * epst_x0 ) ; 
*
x1 = epst_x0 ; 
x2 = x1 ; 
y2 = ( sig0 * ( ( 1 + ( beta0 * x2 ) ) ** m0 ) ) ;
y2 = y2 + ( youn0 * x2 ) ; 
y2 = y2 - ( youn0 * epst_x0 ) ; 
*  
'REPETER' bloc2 ; 
*  
x2 = ( x0 + x1 ) / 2. ; 
y2 = ( sig0 * ( ( 1 + ( beta0 * x2 ) ) ** m0 ) ) ;
y2 = y2 + ( youn0 * x2 ) ; 
y2 = y2 - ( youn0 * epst_x0 ) ; 
*  
'SI' ( y2 > 1.E-4 ) ; 
  x1 = x2 ; 
'SINON' ; 
  'SI' ( y2 <  -1.E-4 ) ; 
     x0 = x2 ; 
  'SINON' ; 
     eps_th0 = x2 ; 
     con_th0 = youn0 * ( epst_x0 - x2 ) ; 
    'QUITTER' bloc2 ;  
  'FINSI' ; 
'FINSI' ; 
*  
'FIN' bloc2 ;
*  
'SINON' ; 
  eps_th0 = 0. ; 
  con_th0 = youn0 * epst_x0 ; 
'FINSI' ; 
*    
'SI' ( eps_th0 '>EG' eps1 ) ; 
*
*  Calcul de rho_thf ( debut d'endommagement)   
*  
  'SI' ( rho_thf '<EG' 0. ) ;  
      tr_f =  tr_eps0 ; 
      sig_f = ( 1. + ( beta0 * eps1 ) ) ** ( m0 ) ; 
      sig_f = sig0 * sig_f ; 
      tr_f = sig_f / ( 3. * 1.75E11 ) ;  
      rho_thf = rho0 / ( 1. + tr_f ) ; 
      err_rho = ( 'ABS' ( rho_f - rho_thf ) ) / rho_thf ; 
      err_rho = err_rho * 100. ; 
      'SI' ( err_rho '>' 1.E-3 ) ; 
         err = vrai ; 
      'FINSI' ; 
  'FINSI' ; 
*
alpha = 1. - ( rho / rho_f ) ;
absc0 = 'EXTR' evol0 'ABSC' ;
ordo0 = 'EXTR' evol0 'ORDO' ;
'REPETER' bloc01 ( 'DIME' absc0 ) ;
  i01 = &bloc01 ;
  ab1 = 'EXTR' absc0 i01 ;
  or1 = 'EXTR' ordo0 i01 ;
  'SI' ( i01 'EGA' ( 'DIME' absc0 ) ) ;
     eps1 = ab1 ;
  'SINON' ;
     ab2 = 'EXTR' absc0 ( i01 + 1 ) ;
     or2 = 'EXTR' ordo0 ( i01 + 1 ) ;
     'SI' ( ( alpha '>EG' or1 ) 'ET'
          ( alpha '<' or2 ) ) ;
         D_end0 = ( or2 - or1 ) / ( ab2 - ab1 ) ;
         D_end0 = D_end0 * ( alpha - ab1 ) ;
         D_end0 = D_end0 + or1 ;
         'QUITTER' bloc01 ;
     'FINSI' ;
   'FINSI' ;
'FIN' bloc01 ;
*
  g_th0 = 1. - D_end0 ; 
'SINON' ; 
  g_th0 = 1. ; 
'FINSI' ;
* 
*  Calcul des contraintes endommagées 
* 
con_th0 = con_th0 * g_th0 ; 
*  
*  'MESS' ' Fin du calcul de restitution pour le pas' ta1 .'TEMPS'.i1 ;  
  max0 = ( 'ABS' ( con_th0 - sig_x0 ) + con_th0 + sig_x0 ) ; 
  max0 = 0.5 * max0 ; 
  max1 = ( 'ABS' ( eps_th0 - def_0 ) + eps_th0 + def_0 ) ; 
  max1 = 0.5 * max1 ;  
  max2 = ( 'ABS' ( g_th0 - g0 ) + g_th0 + g0 ) ;  
  max2 = 0.5 * max2 ; 
  'SI' ( ( 'ABS' max0 ) > 0. ) ; 
     err_s0 = ( 'ABS' ( con_th0 - sig_x0 ) ) / max0 ;  
     err_s0 = err_s0 * 100. ; 
  'SINON' ;  
     err_s0 = 0. ; 
  'FINSI' ; 
  'SI' ( ( 'ABS' max1 ) > 0. ) ; 
     err_e0 = ( 'ABS' ( eps_th0 - def_0 ) ) / max1 ;  
     err_e0 = err_e0 * 100. ; 
  'SINON' ;  
     err_e0 = 0. ; 
  'FINSI' ; 
  'SI' ( ( 'ABS' max2 ) > 0. ) ; 
     err_g0 = ( 'ABS' ( g_th0 - g0 ) ) / max2 ; 
     err_g0 = err_g0 * 100. ; 
  'SINON' ; 
     err_g0 = 0. ; 
  'FINSI' ; 
  'SI' ( ( err_s0 '>' .5 ) 'OU' ( err_e0 '>' .5 ) 'OU' 
         ( err_g0 '>' 1.E-2 ) ) ;  
       err = vrai ; 
  'FINSI' ; 
*  
*  
'FIN' bloc1 ; 
* 
'SI' ( err ) ; 
   'ERRE' 5 ; 
'SINON' ; 
   'ERRE' 0 ; 
'FINSI' ; 
*  
*
'FIN' ;   
 

 

 

 

 

 

© Cast3M 2003 - All rights reserved.
Disclaimer