Download hasofer.procedur

Back to the list

   1 : * HASOFER   PROCEDUR  AM        09/12/07    21:15:37     6578           
   2 : DEBPROC HASOFER TAB1*TABLE;
   3 : 
   4 : *-----------------------------------------------------------------
   5 : * Calcul de l'indice de fiabilite d'Hasofer-Lind
   6 : * Entree :
   7 : *       - nombre de variables                : tab1.nbre_variables
   8 : *       - parametres des variables NATAF    : tab1.param_va
   9 : *       - matrice de covariance              : tab1.matcov;
  10 : *       - coefficient de majoration 1        : tab1.major1
  11 : *       - coefficient de majoration 2        : tab1.major2
  12 : *       - precision relative du BETA         : tab1.prec
  13 : *       - nombre max d'iterations            : tab1.itmax
  14 : *       - point de depart P0                 : tab1.depart
  15 : *       - parametres modele MECA             : tab1.parametres
  16 : * Sortie :
  17 : *       - point d'arrivee P*                 : tab1.arrivee
  18 : *       - indices de fiabilite               : tab1.beta
  19 : *       - nombres d'appels à la FEL          : tab1.appels
  20 : *       - points P*(iter)                    : tab1.beta_point
  21 : *       - valeurs de G(iter)                 : tab1.glim
  22 : *       - ecarts |UP*i(iter+1)-UP*i(iter)|   : tab1.ecu
  23 : *       - increments DELTA_X(iter)           : tab1.delx
  24 : *       - indice d'erreur                    : err1
  25 : *-----------------------------------------------------------------
  26 : 
  27 : 
  28 : 
  29 : * Recuperation des entrees
  30 : *-------------------------
  31 : 
  32 : nva1=tab1.nbre_variables;
  33 : para1=tab1.param_va;
  34 : matc1=tab1.matcov;
  35 : vecu1=tab1.depart;
  36 : SI (EXISTE tab1.parametres);
  37 :   para2=tab1.parametres;
  38 : FINSI;
  39 : 
  40 : * Controle des entrees
  41 : *---------------------
  42 : 
  43 : err1=0;
  44 : SI ((DIME vecu1) NEG nva1);
  45 :  MESS 'HASOFER  : incoherence DIME depart / nva';
  46 :  err1=1;
  47 :  QUIT HASOFER ;
  48 : FINSI;
  49 : SI ((DIME para1) NEG nva1);
  50 :  MESS 'HASOFER  : incoherence DIME param_va / nva';
  51 :  err1=1;
  52 :  QUIT HASOFER ;
  53 : FINSI;
  54 : naux1 = 0;
  55 : REPE bouc1 nva1;
  56 :  naux1=naux1 + &bouc1;
  57 : FIN bouc1;
  58 : SI ((DIME matc1) NEG naux1);
  59 :  MESS 'HASOFER  : incoherence DIME matcov / nva';
  60 :  err1=1;
  61 :  QUIT HASOFER ;
  62 : FINSI;
  63 : 
  64 : * Declarations et initialisations
  65 : *--------------------------------
  66 : 
  67 : * Coefficient d'amplification c1
  68 : maj1=10.;
  69 : 
  70 : * Coefficient d'amplification c2
  71 : maj2=10.;
  72 : 
  73 : * Precision sur beta (inutile en pratique) epsilon_beta
  74 : prec1=0.005;
  75 : 
  76 : * Precision sur la distance entre deux P* consecutifs epsilon_p
  77 : prec2=0.05;
  78 : 
  79 : * Nombre iteration maximal
  80 : itmax1=100;
  81 : 
  82 : tab1.beta=TABLE;
  83 : tab1.beta_pointu=TABLE;
  84 : tab1.beta_pointx=TABLE;
  85 : tab1.glim=TABLE;
  86 : tab1.ecu=TABLE;
  87 : tab1.delx=TABLE;
  88 : tab1.distance=TABLE;
  89 : tab1.gradu=TABLE;
  90 : tab1.gradx=TABLE;
  91 : tab1.cal_grad=TABLE;
  92 : tabx1=TABLE;
  93 : tabx1.nbre_variables=nva1;
  94 : SI (EXISTE para2);
  95 :   tabx1.parametres=para2;
  96 : FINSI;
  97 : tabg1=TABLE;
  98 : tabg1.nbre_variables=nva1;
  99 : SI (EXISTE para2);
 100 :   tabg1.parametres=para2;
 101 : FINSI;
 102 : 
 103 : tab1.appels=0;
 104 : ch1=MOTS;
 105 : id_cal1=LECT;
 106 : REPE bouc1 nva1;
 107 :  id_cal1=id_cal1 ET (LECT 1);
 108 :  ch1=ch1 ET (MOTS (CHAIN 'V' &bouc1));
 109 : FIN bouc1;
 110 : 
 111 : * Calcul de BETA
 112 : beta1 = NORV1 vecu1;
 113 : 
 114 : *------------------------------------------------------------
 115 : * Boucle principale
 116 : *------------------
 117 : 
 118 : REPE boucp1;
 119 : 
 120 :  SI (&boucp1 > itmax1);
 121 :    mess 'HASOFER  : depassement de itmax';
 122 :    err1=1;
 123 :    QUIT boucp1;
 124 :  FINSI;
 125 : 
 126 : * Point courant en espace physique
 127 : *---------------------------------
 128 : 
 129 :  tabn1=TABLE;
 130 :  tabn1.transformation_directe=FAUX;
 131 :  tabn1.points_espace_reference=vecu1;
 132 :  tabn1.noms_des_variables=ch1;
 133 :  tabn1.matcov=matc1;
 134 :  tabn1.param_va=para1;
 135 :  NATAF tabn1;
 136 :  vecx1=tabn1.points_espace_physique;
 137 : 
 138 : * Valeur de la fonction d'EL au point courant
 139 : *--------------------------------------------
 140 : 
 141 :  tabx1.variables=vecx1;
 142 :  err1=FLIM1 tabx1;
 143 :  SI (EGA err1 0);
 144 :    g1=tabx1.valeur_G;
 145 :    precr1=tabx1.precision_R;
 146 :    precs1=tabx1.precision_S;
 147 :    r1=tabx1.valeur_R;
 148 :    s1=tabx1.valeur_S;
 149 :  SINON;
 150 :    QUIT boucp1;
 151 :  FINSI;
 152 :  tab1.appels=tab1.appels+1;
 153 : 
 154 : * Valeurs de la fonction d'EL aux bornes de la perturbation delta_X
 155 : *------------------------------------------------------------------
 156 : 
 157 : * Increments sur les axes en espace physique
 158 : *-------------------------------------------
 159 : 
 160 :  SI (EGA &boucp1 1);
 161 :    delx1=maj1*(precr1+precs1)*vecx1;
 162 :  SINON;
 163 :    delx1=delx2;
 164 :  FINSI;
 165 :  tabg1.vecx=vecx1;
 166 :  tabg1.delx=delx1;
 167 :  tabg1.calcul=id_cal1;
 168 :  err1=GDFLIM1 tabg1;
 169 :  SI (EGA err1 0);
 170 :    delg1=tabg1.valeurs_G;
 171 :    appel1=tabg1.appels;
 172 :  SINON;
 173 :    QUIT boucp1;
 174 :  FINSI;
 175 :  tab1.appels=tab1.appels + appel1;
 176 : 
 177 : * Gradient en espace physique X
 178 : *------------------------------
 179 : 
 180 :  gradx1=PROG;
 181 :  REPE bouc2 nva1;
 182 :    dg1=EXTR &bouc2 delg1;
 183 :    dx1=EXTR &bouc2 delx1;
 184 :    SI (EGA (EXTR &bouc2 id_cal1) 1);
 185 :      gradx1=gradx1 ET (PROG ((dg1-g1)/dx1));
 186 :    SINON;
 187 :      gradx1=gradx1 ET (PROG (EXTR &bouc2 gradx2));
 188 :    FINSI;
 189 :  FIN bouc2;
 190 : 
 191 : * Gradient en espace standard U
 192 : *------------------------------
 193 : 
 194 : * Increment delta_Ui pour calculer le vecteur gradient de X par rapport a Ui
 195 :  delu1=1.E-7;
 196 :  gradu1=PROG;
 197 : 
 198 :  REPE bouc3 nva1;
 199 : 
 200 : * Gradient de X par rapport a Ui
 201 : *-------------------------------
 202 :    ui1=EXTR &bouc3 vecu1;
 203 :    v_delu1=INSERER vecu1 &bouc3 (ui1+delu1);
 204 :    v_delu1=ENLEVER v_delu1 (&bouc3+1);
 205 :    tabn1=TABLE;
 206 :    tabn1.transformation_directe=FAUX;
 207 :    tabn1.points_espace_reference=v_delu1;
 208 :    tabn1.noms_des_variables=ch1;
 209 :    tabn1.matcov=matc1;
 210 :    tabn1.param_va=para1;
 211 :    NATAF tabn1;
 212 :    v_delx1=tabn1.points_espace_physique;
 213 :    gradxu1=(v_delx1 - vecx1)/delu1;
 214 : 
 215 : * Gradient de G par rapport a Ui
 216 : *-------------------------------
 217 :    gradui1=LTL gradx1 gradxu1;
 218 : 
 219 :    gradu1=gradu1 ET (PROG gradui1);
 220 :  FIN bouc3;
 221 : 
 222 : * Nouveau point en espace standard (Rackwitz-Fiessler)
 223 : *-----------------------------------------------------
 224 : 
 225 :  alpu1=(-1.)*gradu1/(NORV1 gradu1);
 226 :  vecu2=((LTL vecu1 alpu1)+(g1/(NORV1 gradu1)))*alpu1;
 227 : 
 228 :  beta2=NORV1 vecu2;
 229 : 
 230 : * Nouveau point en espace physique
 231 : *---------------------------------
 232 : 
 233 :  tabn1=TABLE;
 234 :  tabn1.transformation_directe=FAUX;
 235 :  tabn1.points_espace_reference=vecu2;
 236 :  tabn1.noms_des_variables=ch1;
 237 :  tabn1.matcov=matc1;
 238 :  tabn1.param_va=para1;
 239 :  NATAF  tabn1;
 240 :  vecx2=tabn1.points_espace_physique;
 241 : 
 242 : * Stockage des resultats
 243 : *-----------------------
 244 : 
 245 :  tab1.beta.&boucp1=beta2;
 246 :  tab1.beta_pointu.&boucp1=vecu2;
 247 :  tab1.beta_pointx.&boucp1=vecx2;
 248 :  tab1.glim.&boucp1=g1;
 249 :  tab1.ecu.&boucp1=ABS (vecu2-vecu1);
 250 :  tab1.distance.&boucp1=NORV1 (vecu2-vecu1);
 251 :  tab1.gradu.&boucp1=gradu1;
 252 :  tab1.gradx.&boucp1=gradx1;
 253 :  tab1.iterations=&boucp1;
 254 :  tab1.delx.&boucp1=delx1;
 255 :  tab1.cal_grad.&boucp1 = id_cal1;
 256 : 
 257 : * Verification du critere de convergence entre 2 P* consecutifs
 258 : * (crtitere plus severe que sur 2 BETA consecutifs)
 259 : *--------------------------------------------------------------
 260 : 
 261 :  indc=1;
 262 :  REPETER boucc1 nva1;
 263 :    SI ((EXTR &boucc1 tab1.ecu.&boucp1) < prec2);
 264 :      indc=indc*1;
 265 :    SINON;
 266 :      indc=indc*0;
 267 :    FINSI;
 268 :  FIN boucc1;
 269 : 
 270 :  SI (EGA indc 1);
 271 : 
 272 : *  Convergence : fin
 273 : ********************
 274 : 
 275 :    tab1.arrivee = vecu2;
 276 :    QUIT boucp1;
 277 :  SINON;
 278 : 
 279 : * Non convergence
 280 : *****************
 281 : 
 282 : * Increments optimaux sur les axes en espace physique
 283 : *----------------------------------------------------
 284 : 
 285 : * Distance entre 2 points P* consecutifs
 286 : *---------------------------------------   
 287 :    d1_x12=PROG;
 288 :    REPE bouc4 nva1;
 289 :      x1=(EXTR &bouc4 vecx1);
 290 :      x2=(EXTR &bouc4 vecx2);
 291 :      SI (EGA (EXTR &bouc4 id_cal1) 1);
 292 :        d1_x12=d1_x12 ET (PROG (x2 - x1));
 293 :      SINON;
 294 :        d1_x12=d1_x12 ET (PROG ((x2 - x1) + (EXTR &bouc4 d2_x12)));
 295 :      FINSI;
 296 :    FIN bouc4;
 297 : 
 298 : * Indice de calcul des derivees partielles: 1 calcul, 0 pas de calcul
 299 : *-------------------------------------------------------------------------
 300 :    id_cal1=LECT;
 301 :    REPE bouc5 nva1;
 302 :      dx1=(EXTR &bouc5 delx1);
 303 :      SI (OU ((ABS (EXTR &bouc5 d1_x12)) > dx1) (EGA &boucp1 1));
 304 :        id_cal1=id_cal1 ET (LECT 1);
 305 :      SINON;
 306 :        id_cal1=id_cal1 ET (LECT 0);
 307 :      FINSI;
 308 :    FIN bouc5;
 309 : 
 310 : * Calcul des gradients et increments
 311 : *-----------------------------------
 312 :    gradx1=PROG;
 313 :    delx2=PROG;
 314 :    REPE bouc7 nva1;
 315 :      dx1=(EXTR &bouc7 delx1);
 316 :      SI (EGA (EXTR &bouc7 id_cal1) 1);
 317 :        dg1=(EXTR &bouc7 delg1);
 318 :        gradx1=gradx1 ET (PROG ((dg1-g1)/dx1));
 319 :      SINON;
 320 :        gradx1=gradx1 ET (PROG (EXTR &bouc7 gradx2));
 321 :      FINSI;
 322 :      dx1=maj2*((precr1*r1)+(precs1*s1))/(EXTR &bouc7 gradx1);
 323 :      delx2=delx2 ET (ABS (PROG dx1));
 324 :    FIN bouc7;
 325 : 
 326 : * Mises à jour pour l'iteration suivante
 327 : *---------------------------------------
 328 : 
 329 :    gradx2=gradx1;
 330 :    beta1=beta2;
 331 :    vecu1=vecu2;
 332 :    d2_x12=d1_x12;
 333 : 
 334 :  FINSI;
 335 : 
 336 : FIN boucp1;
 337 : 
 338 : *------------------
 339 : * Fin Boucle principale
 340 : *------------------------------------------------------------
 341 : 
 342 : FINPROC err1;
 343 : 
 344 :  

© Cast3M 2003 - All rights reserved.
Disclaimer