Procédures | Cast3M

* NELDMEAD  PROCEDUR  FD218221  26/09/22    21:15:06     12654          
DEBP NELDMEAD proc1*'PROCEDUR' ;
 
 
 
*---------------- V A L E U R S   P A R   D E F A U T -----------------*
*                     D E S   P A R A M E T R E S                      *
 
lnompar  = MOTS 'ALPH' 'BETA' 'GAMM' 'DELT' 'CRIT' 'PAS' ;
valpar   = PROG     1.     2.    0.5    0.5  1.E-4    1. ;
maxit    = 200 ;
 
 
 
*------------- L E C T U R E   D E S   A R G U M E N T S --------------*
 
* Simplexe initial (obligatoire, mais 2 manieres possibles)
icas     = 1 ;
* cas 1 : donnee d'un LISTREEL
ARGU x0/'LISTREEL' ;
SI (EGA (TYPE x0) 'LISTREEL') ;
* dimension
  n        = DIME x0 ;
  np1      = n + 1 ;
SINON ;
* cas 2 : donnee d'un LISTOBJE de LISTREEL
  ARGU s*'LISTOBJE' ;
  icas     = 2 ;
* test du type d'objets (doit etre LISTREEL)
  SI (NEG (EXTR s 'TYPE') 'LISTREEL') ;
    ERRE 1156 'AVEC' 'LISTREEL' ;
  FINSI ;
* dimension
  np1      = DIME s ;
  n        = np1 - 1 ;
* test de la dimension des LISTREELs
  REPE b1 np1 ;
    lr1      = EXTR s &b1 ;
    SI (NEG (DIME lr1) n) ;
      ERRE 1158 'AVEC' n ;
    FINSI ;
  FIN b1 ;
FINSI ;
 
* Nombre d'iterations (facultatifs)
ARGU i1/'ENTIER' ;
SI (EGA (TYPE i1) 'ENTIER') ;
  maxit    = i1 ;
FINSI ;
 
* Parametres utilisateur (tous facultatifs)
REPE b1 (DIME lnompar) ;
  ARGU mot2/'MOT' ;
  SI (EGA (TYPE mot2) 'MOT') ;
    pos      = POSI mot2 'DANS' lnompar ;
    SI (pos > 0) ;
      ARGU val1*'FLOTTANT' ;
      REMP valpar pos val1 ;
    FINSI ;
  SINON ;
    QUIT b1 ;
  FINSI ;
FIN b1 ;
 
 
 
*-------------------- I N I T I A L I S A T I O N ---------------------*
 
* Si cas 1 (donnee d'un LISTREEL + pas), calcul automatique du simplexe
* initial selon le point initial et le pas d'initialisation
* (Spendley 1962)
SI (EGA icas 1) ;
  tau      = EXTR valpar 6 ;
  p        = tau * ((np1 ** 0.5) + n - 1) / (n * (2. ** 0.5)) ;
  q        = tau * ((np1 ** 0.5)     - 1) / (n * (2. ** 0.5)) ;
  xq       = PROG n*q ;
  s        = ENUM x0 ;
  REPE b1 n ;
    x1       = x0 + xq ;
    vali     = EXTR x0 &b1 ;
    REMP x1 &b1 (vali + p) ;
    s        = s  ET x1 ;
  FIN b1 ;
FINSI ;
 
** Evaluation de la fonction aux points du simplexe
fs       = PROG ;
REPE b1 np1 ;
  xi       = EXTR s &b1 ;
  fi       = PROC1 xi ;
  fs       = fs ET fi ;
FIN b1 ;
 
 
*-------------- B O U C L E   D' O P T I M I S A T I O N --------------*
 
* Parametres d'evolution du simlexe
alpha    = EXTR valpar 1 ;
beta     = EXTR valpar 2 ;
gamma    = EXTR valpar 3 ;
delta    = EXTR valpar 4 ;
crit0    = EXTR valpar 5 ;
 
* Boucle
REPE b1 maxit ;
* Tri des points du simplexe selon les valeurs de fs
  fs s     = ORDO fs s ;
* Points utiles pour decider de l'evolution
  x1       = EXTR s 1 ;
  xnp1     = EXTR s np1 ;
  f1       = EXTR fs 1 ;
  fn       = EXTR fs n ;
  fnp1     = EXTR fs np1 ;
* Centre de gravite de {x1 x2 ... xn}
  x0       = PROG n*0. ;
  REPE b2 n ;
    x0       = x0 + (EXTR s &b2) ;
  FIN b2 ;
  x0       = x0 / n ;
* Reflexion de xnp1 par rapport a x0
  xr       = x0 + (alpha * (x0 - xnp1)) ;
  fr       = PROC1 xr ;
* Si f(x1) <= f(xr) < f(xn), remplacement de xnp1 par xr
  SI ((f1 &lt;EG fr) ET (fr < fn)) ;
    REMP fs np1 fr ;
    REMP s  np1 xr ;
  SINON ;
*   Si f(xr) < f(x1), expansion du simplexe
    SI (fr < f1) ;
      xe       = x0 + (beta * (xr - x0)) ;
      fe       = PROC1 xe ;
*     Si f(xe) <= f(xr), remplacement de xnp1 par xe
      SI (fe &lt;EG fr) ;
        REMP fs np1 fe ;
        REMP s  np1 xe ;
*     Sinon, remplacement de xnp1 par xr
      SINON ;
        REMP fs np1 fr ;
        REMP s  np1 xr ;
      FINSI ;
*   Si f(xn) < f(xr) : contraction du simplexe
    SINON ;
      xc       = x0 + (gamma * (xnp1 - x0)) ;
      fc       = PROC1 xc ;
*     Si f(xc) <= f(xnp1), remplacement de xnp1 par xc
      SI (fc < fnp1) ;
        REMP fs np1 fc ;
        REMP s  np1 xc ;
*     Sinon, homothetie du simplexe autour de x1
      SINON ;
        REPE b2 np1 ;
          xi       = EXTR s &b2 ;
          xi       = x1 + (delta * (xi - x1)) ;
          fi       = PROC1 xi ;
          REMP fs &b2 fi ;
          REMP s  &b2 xi ;
        FIN b2 ;
      FINSI ;
    FINSI ;
  FINSI ;
* Critere d'arret : dispersion des valeurs de la fonction objectif
*                   aux sommets du simplexe
  crit1    = AGRE 'ECTY' fs ;
  SI (crit1 < crit0) ;
    QUIT b1 ;
  FINSI ;
FIN b1 ;
 
 
** On rend le dernier simplexe trie et les valeurs de la fonction associees
fs s     = ORDO fs s ;
FINP s fs ;