Exemples de jeux de données Cast3M

************************************************************************
* Exemple de méthode d'optimisation topologique                        *
* Méthode à densité + pénalisation SIMP                                *
* Algorithme d'optimisation : Method of Moving Asymptotes - Svanberg   *
* Application à une poutre en flexion                                  *
************************************************************************
 
 
** Paramètres globaux
itrac  = FAUX ;
OPTI 'DIME' 2 'MODE' 'PLAN' 'DEFO' 'ELEM' 'QUA4' 'ECHO' 0 ;
l      = 1. ;
h      = 0.2 ;
yo     = 200.E9 ;
nu     = 0.3 ;
 
** Maillage (controlé par sa taille de maille)
den1   = l / 120. ;
OPTI 'DENS' den1 ;
p1     = 0. 0. ;
p2     = l  0. ;
lb     = DROI p1 p2 ;
mail   = lb TRAN (0. h) ;
p3     = mail POIN 'PROC' ((l / 2.) h) ;
con    = CONT mail ;
 
** Modèle
mod    = MODE mail 'MECANIQUE' ;
ma0    = MATE mod 'YOUN' yo 'NU' nu ;
 
** Blocages
blo    = (BLOQ 'UY' p1) ET (BLOQ 'UX' 'UY' p2) ;
 
** Chargement
f      = FORC (0. -1.E5) p3 ;
 
** Paramètres d'optimisation
p      = 3. ;
fv     = 0.4 ;
rmin   = 1.5 * den1 ;
m      = 0.1 ;
xmin   = 0.001 ;
xmax   = 1. ;
 
** Initialisation de la topologie (avec la fraction volumique cible)
x      = MANU 'CHML' mod 'SCAL' fv 'GRAVITE' ;
 
** Matrice de filtrage
un     = MANU 'CHML' mod 'SCAL' 1. 'GRAVITE' ;
vole   = INTG mod un 'ELEM' ;
mcg    = un POIN 'SUPERIEUR' 0. ;
wg     = PROI mcg vole ;
kfil   = MFIL wg rmin 1. 0. ;
 
** Volume plein, initial et cible
v0     = INTG mod un ;
vx     = INTG mod x ;
fvx    = vx / v0 ;
vcib   = v0 * fv ;
chgx   = 0. ;
 
** Initialisation de la table pour la mma
nx     = NBEL mail ;
tmma   = TABL ;
* valeurs initiales de x
tmma . 'VX0' = TABL 'VECTEUR' ;
REPE bx nx ;
  tmma . 'VX0' . &bx = EXTR x 'SCAL' 1 &bx 1 ;
FIN bx ;
* bornes pour les valeurs de x
tmma . 'VXMIN' = TABL 'VECTEUR' ;
tmma . 'VXMAX' = TABL 'VECTEUR' ;
REPE bx nx ;
  tmma . 'VXMIN' . &bx = MAXI (PROG xmin (fv - m)) ;
  tmma . 'VXMAX' . &bx = MINI (PROG xmax (fv + m)) ;
FIN bx ;
* fonction objectif et ses sensibilités
tmma . 'VF' = TABL 'VECTEUR' ;
* fonction contrainte et ses sensibilités
tmma . 'MC' = TABL ;
tmma . 'MC' . 1 = TABL 'VECTEUR' ;
tmma . 'MC' . 1 . 0 = vx ;
REPE bx nx ;
  tmma . 'MC' . 1 . &bx = EXTR vole 'SCAL' 1 &bx 1 ;
FIN bx ;
* bornes pour la fonction contrainte
tmma . 'VCMAX' = TABL 'VECTEUR' ;
tmma . 'VCMAX' . 1 = vcib ;
* paramètres pour la mma
tmma . 'METHODE' = 'STA' ;
tmma . 'T0' = 0.1 ;
 
** Boucle d'optimisation topologique
liso   = PROG 0. 'PAS' 0.05 1. ;
REPE b1 100 ;
* pénalisation de la rigidité
  yop    = (x ** p) * yo ;
  map    = MATE mod 'YOUN' yop 'NU' nu ;
* résolution du problème mécanique
  rip    = RIGI mod map ;
  u      = RESO (rip ET blo) f ;
* fonction objectif : compliance = uT.K.u = Int(sig:eps)dV
  eps    = EPSI 'LINE' mod u ;
  sig    = ELAS mod map eps ;
  psi    = INTG mod (ENER mod sig eps) ;
* sensibilité
  sig0   = ELAS mod ma0 eps ;
  ene0   = ENER mod eps sig0 ;
  ene    = CHAN ene0 mod 'GRAVITE' ;
  dpsi   = (-1. * p * (x ** (p - 1.)) * ene) ;
* filtrage de la sensibilité
  xdpsi  = x * dpsi ;
  xdpsi0 = PROI mcg xdpsi ;
  xdpsi1 = kfil * xdpsi0 ;
  xdpsi  = MANU 'CHML' mod 'REPA' 'SCAL' (EXTR xdpsi1 'VALE') 'TYPE' 'SCALAIRE' 'GRAVITE' ;
  dpsi   = xdpsi / x ;
* infos sur la topologie courante
  info   = CHAI 'It:' (&b1 - 1) / 5 'Obj:' / 10  psi > 1 'Fvol:' > 4 fvx > 1 'Change:' > 4 chgx > 1 ;
  MESS info ;
  SI itrac ;
    def1   = DEFO mail u 1000. ;
    TRAC x mod con def1 liso 'TITR' info 'NCLK' ;
  FINSI ;
* optimisation d'une nouvelle topologie (méthode des asymptotes mobiles)
  tmma . 'VF' . 0 = psi ;
  REPE bx nx ;
    tmma . 'VF' . &bx = EXTR dpsi 'SCAL' 1 &bx 1 ;
  FIN bx ;
  tmma . 'MC' . 1 . 0 = INTG mod x ;
  tmma_new = EXCE tmma ;
  lxnew  = PROG ;
  REPE bx nx ;
    lxnew  = lxnew ET (tmma_new . 'VX0' . &bx) ;
  FIN bx ;
  xnew   = MANU 'CHML' mod 'REPA' 'SCAL' lxnew 'TYPE' 'SCALAIRE' 'GRAVITE' ;
  vxnew  = INTG mod xnew ;
* bilan de l'itération
  fvx    = vxnew / v0 ;
  chgx   = MAXI 'ABS' (x - xnew) ;
  SI (chgx < 0.01) ;
    info   = CHAI 'It:' &b1 / 5 'Obj:' / 10  psi > 1 'Fvol:' > 4 fvx > 1 'Change:' > 4 chgx > 1 ;
    MESS info ;
    QUIT b1 ;
  FINSI ;
* préparation de la nouvelle itération
  x      = xnew ;
  tmma . 'VX0' = tmma_new . 'VX0' ;
  lxmin = BORN (lxnew - m) 'MINIMUM' xmin ;
  lxmax = BORN (lxnew + m) 'MAXIMUM' xmax ;
  REPE bx nx ;
    tmma . 'VXMIN' . &bx = EXTR lxmin &bx ;
    tmma . 'VXMAX' . &bx = EXTR lxmax &bx ;
  FIN bx ;
FIN  b1 ;
 
* Tracé de la dernière itération
SI itrac ;
  def1   = DEFO mail u 1000. ;
  TRAC xnew mod con def1 liso 'TITR' info ;
FINSI ;
 
FIN ;