Accueil | Cast3M

* fichier :  dyna_nl2.dgibi
************************************************************************
************************************************************************
complet = vrai;
*
*************************************************
*  pour calcul complet mettre complet à : vrai;
*************************************************
*************************************************************************       
*
*         dynamique non lineaire geometrique
*
*         reponse forcee d'un oscillateur de type Duffing
*
*         comparaison avec calcul explicite
*               
************************************************************************* 
opti echo 0 ;
graph = faux ;
 
 
OPTI DIME 2 ELEM SEG2 MODE PLAN CONT ; ;
P1 = 0.    0.    ; P2 = 1. 0.  ;
ST = P1 D 1  P2 ;
 
MO = MODE ST  MECANIQUE  BARR     ;                                                                               
MA1 = MATE MO YOUN 1E7  NU 0.3 RHO 0;
CAR1 = CARA MO SECT 1. ;
MACA= MA1 ET CAR1;
MASS1 = MASS UY 250. P2 ;
 
CL1 = BLOQ UX uy  P1 ;
CL2 = BLOQ UX P2 ;
cl3 = bloq uy p2 ;
CL = CL1 ET CL2 ;
 
 
* precontrainte de 500 N 
 
RI1 = RIGI MACA  MO ;
f1 = force fx 500 p2 ;
dep1 = reso (ri1 et cl1 et cl3) f1 ;
sig1 = sigm dep1 mo maca ;
*sig1 = manu 'CHML' mo  'EFFX' 500. 'TYPE' 'CONTRAIN'
*        'STRESSES' ;
 
* calcul dynamique
si complet;
t1 = .4 ;
sinon;
t1 = .1;
finsi;
dt = .005 ;
f1 = force fy 25000. p2 ;
LIX1 = PROG 0. pas dt t1    ;                                                            
LIY1 = (prog (dime lix1)*1. ) - (cos (360*lix1/4/.26479)) ; 
*lix1 = prog 0 8.9 100 ; liy1 = prog 0 1 1 ;                                                           
EV1 = EVOL MANU 'temps' LIX1 'force' LIY1 ;                                                                                                                                                         
CHA1 = CHAR MECA F1   EV1 ;
 
tab2 = table ;                                                                 
 
tab2.grands_deplacements = vrai ;
*tab2.tinit = vrai ;
*tab2.acceleration = 2 ;
tab2.dynamique = vrai ;
tab2.masse_constante = mass1 ;
TAB2.'BLOCAGES_MECANIQUES' = CL  ;
TAB2.'MODELE' = MO;
Tab2.'K_SIGMA' = VRAI ;
TAB2.'CARACTERISTIQUES' = MACA;
TAB2.'CHARGEMENT' = CHA1;
TAB2.'CONTRAINTES' = TABLE;
TAB2.'CONTRAINTES' . 0 = SIG1 ;
LIS11 = PROG 0.  pas dt t1  ;
TAB2.'TEMPS_CALCULES' = LIS11;
PASAPAS TAB2 ;
 
 
NDIM = (DIME ( TAB2 . DEPLACEMENTS )) ;
PGY = PROG ;
REPETER TBOU2 NDIM;
 LEDEP = TAB2 . DEPLACEMENTS.  (&TBOU2 - 1);
 V   = EXTR LEDEP UY P2 ;
 PGY = PGY ET ( PROG V ) ;
FIN TBOU2 ;
ev1 = evol manu 'temps' lis11  'depl' pgy ;
 
*******************************************************************
**
** calcul explicite (differences centrées, 
** on calcule la force de rappel analytiquement)
*
m = 250. ;
x0 = 0. ;
v0 = 0. ;
si complet;
npas = 400 ;
sinon;
npas = 100;
finsi;
dt = .001 ;
f1 =  25000 ;
LIX1 = PROG  pas dt  1.  ;                                                            
lisf1 = f1*((prog (dime lix1)*1. ) - (cos (360*lix1/4/.26479))) ; 
 
n = 0 ;
lisx = prog x0 ;
lisv = prog v0 ;
 repeter toto npas ;
   n = n + 1 ;
   fn = extr lisf1 n ;
    si (n ega 1 ) ;
     xn = x0 ; 
     finT = (500*xn/((1 + (xn**2))**.5))  +
                       (1e7*xn*(1 - (1/((1+(xn**2))**.5)))) ;
     xnm = x0-(dt*v0)- (((dt**2)/2.)*((fint-fn)/m)) ;
    finsi ;
     finT = (500*xn/((1 + (xn**2))**.5))  +
                       (1e7*xn*(1 - (1/((1+(xn**2))**.5)))) ;
     xnp = (-1.)*((dt**2)/m*(finT - fn)) + (2*xn) - xnm ;
       si (n > 1) ;
         vn = (xnp - xnm)/2./dt ;
         lisv = lisv et (prog vn) ;
       finsi ;
     xnm = xn ;
     xn = xnp ;
    lisx = lisx et (prog xnp) ;
 fin toto ;
list1 = (prog 0) et (prog  pas dt (npas*dt)) ;
evx = evol rouge manu  'temps' list1 'depl' lisx ;
*************************************************************************
*
* comparaison entre les deux calculs
*
 
 si graph ;
    dess (evx et ev1) ;
 finsi ;
 
pgy2 = ipol list1 lis11 pgy ;
ddepl = (abs (lisx - pgy2))/ (lisx + (prog (dime lisx)*1e-5)) ;
err1 = maxi ddepl ;
mess 'difference entre pasapas et "explicite_analytique"'(err1*100)'%';
 
 si (err1 > 7);
  erre 5 ;
 sinon ;
  erre 0 ;
 finsi ; 
 
fin;