Accueil | Cast3M

* fichier :  mma_01.dgibi
************************************************************************
************************************************************************
 
************************************************************************
*      Test de l'opérateur MMA : Méthode des Asymptotes Mobiles        *
*    Application a une fonction a 2 variables et sans contraintes      *
*                                                                      *
* Le problème s'écrit :                                                *
*    Minimiser :  f0(x) = (x1-50)² + (x2-25)² + 25                     *
*      sur x                                                           *
*         avec :  1 <= xi <= 100   pour i=1,2                          *
*                                                                      *
* La solution (evidente) est : x1 = 50   x2 = 25                       *
************************************************************************
 
* Options
OPTI 'ECHO' 0 ;
itrac = FAUX ;
 
* Procedures pour calculer f0 ainsi que les derivees partielles df0/dx
* Juste les valeurs de la fonction
DEBP F00 x1 x2 ;
  f0 = ((x1 - 50.) ** 2) + ((x2 - 25.) ** 2) + 25. ;
FINP f0 ;
 
* La fonction et ses derivees
DEBP FONC lx*'LISTREEL' ;
  x1 = EXTR lx 1 ;
  x2 = EXTR lx 2 ;
  f0 = F00 x1 x2 ;
  f  = PROG 0. ;
  df0dx = 2. * (lx - (PROG 50. 25.)) ;
  dfdx = TABL ;
  dfdx . 1 = PROG 0. 0. ;
FINP f0 df0dx f dfdx ;
 
* Solution de reference
x1ref = 50. ;
x2ref = 25. ;
xref  = PROG x1ref x2ref ;
f0ref df0dxref fval0 dfdxref = FONC xref ;
MESS ;
MESS 'Solution de reference' ;
MESS '     x1          x2          f0' ;
MESS (CHAI 'FORMAT' '(F10.5)' x1ref /4 x2ref /16 f0ref /28) ;
 
* Choix d'un point initial pour les inconnues
x0 = PROG 1. 1. ;
 
* Calcul de f0(x0), df0dx(x0), fval(x0), dfdx(x0)
f0 df0dx fval dfdx = FONC x0 ;
 
* Initialisation de la table pour la MMA
t = TABL ;
t . 'X'     = x0 ;
t . 'XMIN'  = 1. ;
t . 'XMAX'  = 100. ;
t . 'F0VAL' = f0 ;
t . 'DF0DX' = df0dx ;
t . 'FVAL'  = fval ;
t . 'DFDX'  = dfdx ;
t . 'A0'    = 1. ;
t . 'A'     = PROG 0. ;
t . 'C'     = PROG 1.E7 ;
t . 'D'     = PROG 1. ;
t . 'MOVE'  = 0.1 ;
 
* Iterations de la MMA
x01 = EXTR x0 1 ;
x02 = EXTR x0 2 ;
lx1 = PROG x01 ;
lx2 = PROG x02 ;
lit = PROG 0. ;
lf0 = PROG f0 ;
f1  = EXTR fval 1 ;
li  = PROG (MAXI 0. f1) ;
MESS ;
MESS 'Optimisation par MMA' ;
MESS 'It   x1          x2          f0         kktnorm' ;
MESS (CHAI 'FORMAT' '(F10.5)' 0 x01 /4 x02 /16 f0 /28) ;
* Boucle d'optimisation
REPE loop 25 ;
* Appel a MMA
  lit = lit ET &loop ;
  MMA t ;
  xmma = t . 'X' ;
  x1 = EXTR xmma 1 ;
  x2 = EXTR xmma 2 ;
  lx1 = lx1 ET x1 ;
  lx2 = lx2 ET x2 ;
* Mise a jour des valeurs des fonctions
  f0 df0dx fval dfdx = FONC xmma ;
  lf0  = lf0 ET f0 ;
  f1   = EXTR fval 1 ;
  li   = li ET (MAXI 0. f1) ;
  t . 'F0VAL'  = f0 ;
  t . 'DF0DX'  = df0dx ;
  t . 'FVAL'   = fval ;
  t . 'DFDX'   = dfdx ;
* Calcul du residu pour les conditions KKT
  res kkt2 kktinf = KKT_MMA t ;
* Bilan de l'iteration
  MESS (CHAI 'FORMAT' '(F10.5)' &loop x1 /4 x2 /16 f0 /28 kkt2 /40) ;
FIN loop ;
nit = (DIME lf0) - 1 ;
MESS ;
 
* Verification du resultat
errmax = MAXI 'ABS' ((xmma - xref) / xref) ;
MESS 'Erreur relative max (sur x)' ;
MESS errmax ;
MESS ;
 
* Évolutions temporelles de f, de l'infaisabilité et des variables
* en fonction des iterations d'optimisation
SI itrac ;
  OPTI 'DIME' 2 'ELEM' 'QUA8' ;
  msh   = (DROI 20 ((t . 'XMIN') (t . 'XMIN')) ((t . 'XMAX') (t . 'XMIN'))) TRAN 20 (0. ((t . 'XMAX') - (t . 'XMIN'))) ;
  cmsh  = CONT msh ;
  x y   = COOR msh ;
  f0msh = F00 x y ;
  f0msh = ((x - 50.) ** 2) + ((y - 25.) ** 2) + 25. ;
  path  = QUEL 'SEG2' lx1 lx2 ;
  pdep  = x01 x02 ;
  pfin  = x1 x2 ;
  cdep  = (CERC 10 'ROTA' 360. (pdep PLUS (1. 0.)) pdep) COUL 'VIOL' ;
  cfin  = (CERC 10 'ROTA' 360. (pfin PLUS (1. 0.)) pfin) COUL 'ROUG' ;
  annd  = ANNO 'ETIQ' pdep 'VIOL' 'NO' 10. VRAI 'Depart' ; 
  annf  = ANNO 'ETIQ' pfin 'ROUG' 'NE' 10. VRAI 'Arrivee' ; 
  TRAC f0msh msh (cmsh ET path ET cdep ET cfin) 25 (annd ET annf) 'TITR' 'Isovaleurs de la fonction objectif F0 et chemin au cours de l''optimisation' ;
  lit   = PROG 0. 'PAS' 1. nit ;
  evf0  = EVOL 'ROUG' 'MANU' 'Iteration' lit 'F0'  lf0 ;
  evfr  = EVOL 'ROUG' 'MANU' 'Iteration' (PROG 0. nit) 'F0' (PROG f0ref f0ref) ;
  tl    = TABL ;
  tl . 2 = 'TIRR' ;
  tl . 'TITRE' = TABL ;
  tl . 'TITRE' . 1  = 'F0 calculee' ;
  tl . 'TITRE' . 2  = 'F0 ref.' ;
  DESS (evf0 ET evfr) 'TITR' 'Fonction objectif F0(x) VS Iterations' 'LEGE' tl ;
  evi   = EVOL 'VERT' 'MANU' 'Iteration' lit 'Inf' li ;
  tl . 'TITRE' . 1  = 'Infaisabilite' ;
  DESS evi 'TITR' 'Infaisabilite VS Iterations' 'LEGE' tl ;
  evx1  = EVOL 'ROUG' 'MANU' 'Iteration' lit 'x' lx1 ;
  evx1r = EVOL 'ROUG' 'MANU' 'Iteration' (PROG 0. nit) 'x' (PROG x1ref x1ref) ;
  evx2  = EVOL 'ORAN' 'MANU' 'Iteration' lit 'x' lx2 ;
  evx2r = EVOL 'ORAN' 'MANU' 'Iteration' (PROG 0. nit) 'x' (PROG x2ref x2ref) ;
  tl . 4  = 'TIRR' ;
  tl . 'TITRE' . 1  = 'x1 calcule' ;
  tl . 'TITRE' . 2  = 'x1 ref.' ;
  tl . 'TITRE' . 3  = 'x2 calcule' ;
  tl . 'TITRE' . 4  = 'x2 ref.' ;
  DESS (evx1 ET evx1r ET evx2 ET evx2r) 'TITR' 'Valeurs x VS Iterations' 'LEGE' tl ;
FINSI ;
 
* Sortie en erreur si l'écart aux valeurs de références est trop important
SI (errmax > 1.E-8) ;
  ERRE 'Erreur dans le calcul d''optimisation' ;
SINON ;
  MESS 'Cas test passe avec succes !' ;
FINSI ;
 
 
FIN ;