Télécharger exec.notice

Retour à la liste thématique

Afficher cette notice en

Numérotation des lignes :
   1 : $$$$ EXEC     NOTICE  GOUNAND   14/06/04    21:15:03     8050           
   2 :                                              DATE     14/06/04
   3 : 
   4 :     Procedure EXEC                           Voir aussi : EQEX MODE
   5 :     --------------                                        DOMA  KCHT  
   6 :                                                           KRES  
   7 :     EXEC  TAB1  ;                                       
   8 : 
   9 : 
  10 : 
  11 : 
  12 : 
  13 :     Objet :
  14 :     _______
  15 : 
  16 :    I/ TRANSPORT/DIFFUSION D'UN SCALAIRE
  17 : 
  18 :     Cette procedure permet d'effectuer un calcul du transport
  19 :    (convection/diffusion), d'un scalaire passif, en transitoire ou en
  20 :    regime permanent, en presence de conditions aux limites variees, valeur
  21 :    imposee, flux, echange avec le milieu exterieur, source etc. Le
  22 :    calcul peut etre lineaire ou non. Les non linearites sont resolues
  23 :    par une methode de point fixe. La description de l'equation a
  24 :    resoudre se fait a l'aide de l'operateur EQEX qui cree la table TAB1.
  25 :    Le scalaire peut aussi bien etre la temperature qu'une concentration
  26 :    ou toute grandeur physique intensive.
  27 :    Les parametres de l'algorithme sont definis dans la table TAB1 a l'aide
  28 :    de EQEX.
  29 : 
  30 : 
  31 :    I.1/ Calcul transitoire explicite.
  32 : 
  33 :    Le pas de temps est soumis a une contrainte de stabilite. Il peut etre
  34 :    impose ou calcule automatiquement.
  35 :    Les non linearites, notamment sur les proprietes physiques,peuvent
  36 :    etre resolues en les actualisant dans une procedure de calcul appelee
  37 :    a chaque pas de temps. Le regime permanent peut etre obtenu comme
  38 :    limite asymptotique du transitoire.
  39 : 
  40 : *.Exemple I.1 :.........................................................
  41 : * Procedure calculant une propriete physique dependant de la temperature
  42 : 
  43 :       'DEBP' CALCUL ;
  44 :       'ARGU' RX*TABLE ;
  45 :       iarg = rx . 'IARG' ;
  46 :       rv   = rx . 'EQEX' ;
  47 : *Lecture des arguments de l'operateur calcul
  48 :       'SI' ( 'NON' ( 'EGA' iarg 1)) ;
  49 :         'MESS' 'Procedure CALCUL : nombre d arguments incorrect ' iarg ;
  50 :         'QUIT' CALCUL ;
  51 :       'FINSI' ;
  52 : 
  53 :       'SI' ( 'EGA' ('TYPE' rx . 'ARG1') 'MOT     ') ;
  54 :         TN = rv . 'INCO' . (rx . 'ARG1') ;
  55 :       'SINON' ;
  56 :         'MESS' 'Procedure CALCUL : type argument invalide ' ;
  57 :         'QUIT' CALCUL ;
  58 :       'FINSI' ;
  59 : 
  60 : * La temperature est exprimee en Kelvin
  61 :       T = TN + 273. ;
  62 : *Viscosite dynamique : loi de Sutherland : Kg/m/s
  63 :       MU = 1.648*(T**1.5) * ('INVE' (T + 0.648))
  64 : *conductivite : loi de Sutherland : W/m/oC
  65 :       LB = 1.368*(T**1.5) * ('INVE' (T + 0.368))
  66 : *J/kg/oC
  67 :       CP = 1015 ;
  68 : * Nombre de Prandtl
  69 :       Pr = MU * CP * ('INVE' LB)
  70 : * le Reynolds est donne
  71 :       Re = 400. ;
  72 :       Pe = Re * Pr ;
  73 :       rv . 'INCO' . 'IPE' = 'INVE' Pe ;
  74 : * La derniere instruction cree des objet vides
  75 : * pour satisfaire la procedure EXEC
  76 :       as2 ama1  = 'KOPS' 'MATRIK' ;
  77 :       'FINPROC' as2 ama1 ;
  78 : 
  79 : * On cree la table RV decrivant le probleme physique
  80 : * On choisit un algorithme explicite (OPTI 'EFM1')
  81 : * Le pas de temps est calcule automatiquement
  82 : * (mot cle 'DELTAT' sur DFDT)
  83 : * On fera 200 pas de temps
  84 : * Le Peclet est calcule dans la procedure CALCUL
  85 : 
  86 :       RV = 'EQEX' 'OMEGA' 1. 'NITER' 1 'ITMA' 200
  87 :       'ZONE' $mt  'OPER' CALCUL 'TN'
  88 :       'OPTI' 'EFM1' 'SUPG'
  89 :       'ZONE' $mt  'OPER' 'TSCA' 'IPE' 'UN' 0. 'INCO' 'TN'
  90 :       'OPTI' 'EFM1' 'CENTREE'
  91 :       'ZONE' $mt  'OPER' 'DFDT' 1. 'TN' 'DELTAT' 'INCO' 'TN'
  92 :       'CLIM' 'TN' 'TIMP' entree 0. 'TN' 'TIMP' paroi  1.
  93 :       ;
  94 :       rv . 'INCO' = 'TABLE' 'INCO' ;
  95 :       rv . 'INCO' .'UN'= 'KCHT' $mt 'VECT' 'SOMMET' (1. 0.) ;
  96 :       rv . 'INCO' .'TN'= 'KCHT' $mt 'SCAL' 'SOMMET' 0. ;
  97 : 
  98 :       EXEC RV ;
  99 : 
 100 : * Les resultats se trouvent dans la table rv . 'INCO'
 101 : 
 102 : *.Fin exemple I.1 ......................................................
 103 : 
 104 :    I.2/ Calcul direct d'un regime permanent.
 105 : 
 106 :    - On peut faire une recherche directe d'un regime permanent avec
 107 :    des iterations internes pour resoudre les non-linearites.
 108 : 
 109 : *.Exemple I.2 :.........................................................
 110 : 
 111 : * On cree la table RV decrivant le probleme physique
 112 : * On choisit un algorithme IMPLICITE (OPTI 'EF' 'IMPL')
 113 : * On fera 10 iterations avec un facteur de relaxation de OMEGA=0.5
 114 : * Le Peclet est calcule dans la procedure CALCUL decrite dans
 115 : * l'exemple 1.
 116 : 
 117 :       RV = 'EQEX' 'OMEGA' 0.5 'NITER' 10 'ITMA' 0
 118 :       'ZONE'  $mt 'OPER' CALCUL 'TN'
 119 :       'OPTI' 'EF' 'SUPG' 'IMPL'
 120 :       'ZONE'  $mt 'OPER' 'TSCA' 'IPE' 'UN' 0. 'INCO' 'TN'
 121 :       'CLIM' 'TN' 'TIMP' entree 0. 'TN' 'TIMP' paroi  1.
 122 :       ;
 123 :       rv . 'INCO' = 'TABLE' 'INCO' ;
 124 :       rv . 'INCO' . 'UN' = 'KCHT' $mt 'VECT' 'SOMMET' (1. 0.) ;
 125 :       rv . 'INCO' . 'TN' = 'KCHT' $mt 'SCAL' 'SOMMET' 0. ;
 126 : 
 127 :       EXEC RV ;
 128 : 
 129 : * Les resultats se trouvent dans la table rv . 'INCO'
 130 : 
 131 : *.Fin exemple I.2 ......................................................
 132 : 
 133 :    I.3/ Calcul transitoire implicite.
 134 : 
 135 :    - On peut faire un calcul transitoire implicite avec ou sans
 136 :    iterations internes a chaque pas de temps.
 137 : 
 138 : *.Exemple I.3 :.........................................................
 139 : 
 140 : * On cree la table RV decrivant le probleme physique
 141 : * On choisit un algorithme IMPLICITE (OPTI 'EF' 'IMPL') ordre 1 en
 142 : * temps ou mieux un schema en temps d'ordre 2 (Crank Nicolson)
 143 : * ('OPTI' 'EF' 'SUPG' 'SEMI' 0.5)
 144 : * On fera 10 pas de temps sans iteration interne
 145 : * Le Peclet est calcule dans la procedure CALCUL decrite dans
 146 : * l'exemple 1.
 147 : 
 148 :       dt = 1. ;
 149 : 
 150 :       RV = 'EQEX' 'OMEGA' 1. 'NITER' 1 'ITMA' 10
 151 :       'ZONE'  $mt 'OPER' CALCUL 'TN'
 152 :       'OPTI' 'EF' 'SUPG' 'SEMI' 0.5
 153 :       'ZONE'  $mt 'OPER' 'TSCA' 'IPE' 'UN' 0. 'INCO' 'TN'
 154 :       'CLIM' 'TN' 'TIMP' entree 0. 'TN' 'TIMP' paroi  1.
 155 :       'OPTI' 'EF' 'CENTREE'
 156 :       'ZONE'  $mt 'OPER' 'DFDT' 1. 'TN' dt 'INCO' 'TN'
 157 :       ;
 158 :       rv . 'INCO' = 'TABLE' 'INCO' ;
 159 :       rv . 'INCO' . 'UN'= 'KCHT' $mt 'VECT' 'SOMMET' (1. 0.) ;
 160 :       rv . 'INCO' . 'TN'= 'KCHT' $mt 'SCAL' 'SOMMET' 0. ;
 161 : 
 162 :       EXEC RV ;
 163 : 
 164 : * Les resultats se trouvent dans la table rv . 'INCO'
 165 : 
 166 : *.Fin exemple I.3 ......................................................
 167 : 
 168 : 
 169 :    II/ NAVIER STOKES INCOMPRESSIBLE
 170 : 
 171 :     La procedure permet de resoudre les equations de Navier_Stokes par
 172 :    une methode d'elements finis (EF) en variables primitives (vitesse et
 173 :    pression), pour un ecoulement incompressible, ou faiblement compre-
 174 :    ssible, en regime permanent ou en transitoire. Les conditions limites
 175 :    peuvent etre variees, vitesse imposee, pression imposee, source de
 176 :    quantite de mouvement, perte de charge ... etc. Le calcul peut etre
 177 :    lineaire ou non. Les non linearites sont resolues par une methode de
 178 :    point fixe. La description de l'equation a resoudre se fait a
 179 :    l'aide de l'operateur EQEX qui cree la table TAB1.
 180 :    Trois algorithmes sont proposes pour resoudre le systeme vitesse-
 181 :    pression.
 182 :    - Un algorithme semi-explicite : implicite pour la pression, explicite
 183 :    pour la vitesse et eventuellement toutes les autres quantites
 184 :    scalaires transportees (a la Gresho).
 185 :    - Un algorithme implicite : resolution directe du systeme en vitesse-
 186 :    pression (a la Taylor-Hood). Cet algorithme permet aussi de faire une
 187 :    recherche directe d'un regime permanent.
 188 :    - Une methode de projection : cet algorithme consiste en deux etapes :
 189 :    une premiere etape de transport diffusion et une seconde de projection
 190 :    sur un espace a divergence nulle (a la Chorin Temam) .
 191 : 
 192 :     Les parametres de l'algorithme sont definis dans la table TAB1
 193 :    Les non linearites, notamment sur les proprietes physiques,peuvent
 194 :    etre resolues, comme precedemment, par l'intermediaire d'une methode de
 195 :    point fixe, en actualisant dans une procedure les non linearites.
 196 :    de calcul appelee a chaque pas de temps ou chaque iteration.
 197 :     Le regime permanent peut etre obtenu comme limite asymptotique du
 198 :    transitoire
 199 : 
 200 : 
 201 :    II.1/ Calcul semi-explicite
 202 : 
 203 :    - On peut faire un calcul transitoire semi explicite. La pression est
 204 :    implicite, la vitesse explicite. De ce fait le pas de temps est soumis
 205 :    a une contrainte de stabilite de type CFL liee a la convection ou a la
 206 :    diffusion. Le pas de temps peut etre impose ou calcule automatiquement
 207 :    (mot cle 'DELTAT' en 3eme argument de DFDT).
 208 :    L'algorithme necessite en general un grand nombre de pas de temps.
 209 :    La construction de la table se fait en deux etapes.
 210 :    - Une premiere etape consiste a creer une table decrivant la partie
 211 :    explicite.
 212 :    - Dans la deuxieme etape on construit la table decrivant les operateurs
 213 :    lies a la pression.
 214 :    Les conditions de vitesse normale ou vitesse tangente en font partie.
 215 :    Elles sont traitees a l'aide de multiplicateurs de Lagrange. Par contre
 216 :    dans nos formulations elements finis (formulation faible) la condition
 217 :    de pression imposee n'en fait pas partie. La pression est imposee comme
 218 :    les contraintes visqueuses sur l'equation de quantite de mouvement.
 219 :    - Enfin on place la deuxieme table a l'indice 'POISSON' de la premiere.
 220 : 
 221 :     L'indice rv.'CALPRE' = VRAI de la premiere table indique que l'on
 222 :    recalcule a chaque pas de temps la matrice de pression. Ceci est
 223 :    necessaire en formulation A.L.E.
 224 :     rv.'CALPRE' = FAUX ou l'absence de cet indice fait que la matrice de
 225 :    pression sera calculee une fois pour toute.
 226 : 
 227 :    Le positionnement de la variable rv.'DETMAT' a VRAI indique que les
 228 :     objets MATRIK seront detruits a la fin de la procedure.
 229 : 
 230 : *.Exemple II.1  :...........................................................
 231 : * Cas de la cavite carree a paroi defilante
 232 : * !!!! ATTENTION : La cavite etant fermee (V.n impose sur toute la frontiere)
 233 : * il est necessaire d'imposer la pression en un point.
 234 : 
 235 : 
 236 :    ro = 400. ;
 237 :    mu = 1. ;
 238 : 
 239 :    RV = 'EQEX' 'OMEGA' 1. 'NITER' 1 'ITMA' 500 'ALFA' 0.5
 240 :    'OPTI' 'EFM1' 'SUPG'
 241 :    'ZONE'  $mt   'OPER' 'NS' (mu/ro) 'INCO' 'UN'
 242 :    'OPTI' 'EFM1' 'CENTREE'
 243 :    'ZONE'  $mt   'OPER' 'DFDT' 1. 'UN' 'DELTAT' 'INCO' 'UN'
 244 :    ;
 245 :    RV = 'EQEX' RV
 246 :    'CLIM'  'UN' 'UIMP' CD  1. 'UN' 'VIMP' CD 0.
 247 :            'UN' 'UIMP' DA  0. 'UN' 'VIMP' DA 0.
 248 :            'UN' 'UIMP' AB  0. 'UN' 'VIMP' AB 0.
 249 :            'UN' 'UIMP' BC  0. 'UN' 'VIMP' BC 0. ;
 250 : * Le choix d'une methode iterative (Bicg stab + preconditionement MILU0)
 251 : * permet de passer des cas un peu plus gros. (Faire INFO KRES ; )
 252 :    rv. 'METHINV' . 'TYPINV'  = 3 ;
 253 :    rv. 'METHINV' . 'IMPINV'  = 0 ;
 254 :    rv. 'METHINV' . 'NITMAX'  = 400 ;
 255 :    rv. 'METHINV' . 'PRECOND' = 3 ;
 256 :    rv. 'METHINV' . 'RESID'   = 1.e-8 ;
 257 :    rv. 'METHINV' . 'FCPRECT' = 1 ;
 258 :    rv. 'METHINV' . 'FCPRECI' = 1 ;
 259 : 
 260 :    betastab=1.e2 ;
 261 : 
 262 :    RVP = 'EQEX'
 263 :    'OPTI' 'EF' 'CENTRE'
 264 :    'ZONE' $mt 'OPER' 'KBBT' -1. betastab 'INCO' 'UN' 'PRES'
 265 :    'CLIM' 'PRES' 'TIMP' bcp 0. ;
 266 : 
 267 :    rvp . 'METHINV' . 'TYPINV'  = 2 ;
 268 :    rvp . 'METHINV' . 'IMPINV'  = 0 ;
 269 :    rvp . 'METHINV' . 'NITMAX'  = 300 ;
 270 :    rvp . 'METHINV' . 'PRECOND' = 3 ;
 271 :    rvp . 'METHINV' . 'RESID'   = 1.e-8 ;
 272 :    rvp . 'METHINV' . 'FCPRECT' = 100 ;
 273 :    rvp . 'METHINV' . 'FCPRECI' = 100 ;
 274 : 
 275 :    rv . 'POISSON' = rvp ;
 276 : 
 277 :    rv . 'INCO' = 'TABLE' 'INCO' ;
 278 :    rv . 'INCO' .'UN'   = 'KCHT' $mt 'VECT' 'SOMMET' (0. 0.) ;
 279 :    rv . 'INCO' .'PRES' = 'KCHT' $mt 'SCAL' 'CENTRE' 0. ;
 280 : 
 281 :     EXEC RV ;
 282 : 
 283 : * Les resultats se trouvent dans la table rv . 'INCO'
 284 : 
 285 : *.Fin exemple II.1 .....................................................
 286 : 
 287 : 
 288 :    II.2/ Calcul implicite.
 289 : 
 290 :    - On peut faire un calcul transitoire implicite avec ou sans
 291 :    iterations internes a chaque pas de temps.
 292 : 
 293 : *.Exemple II.2  :...........................................................
 294 : * Cas de la cavite carree a paroi defilante
 295 : * !!!! ATTENTION : La cavite etant fermee (V.n impose sur toute la frontiere)
 296 : * il est necessaire d'imposer la pression en un point.
 297 : 
 298 : 
 299 :    ro = 400. ;
 300 :    mu = 1. ;
 301 :    dt = 5. ;
 302 : 
 303 :    RV = 'EQEX' 'OMEGA' 1. 'NITER' 1 'ITMA' 20
 304 :    'OPTI' 'EF' 'IMPL' 'SUPG'
 305 :    'ZONE  $mt  'OPER' 'LAPN' mu 'INCO' 'UN'
 306 :    'ZONE  $mt  'OPER' 'KONV' ro 'UN' mu dt 'INCO' 'UN'
 307 :    'OPTI' 'EF' 'CENTREE'
 308 :    'ZONE  $mt  'OPER' 'DFDT' ro 'UN' dt 'INCO' 'UN'
 309 :    'OPTI' 'EF' 'CENTREP1'
 310 :    'ZONE  $mt  'OPER' 'KBBT' 1. 'INCO' 'UN' 'PRES'
 311 :    ;
 312 :    RV = 'EQEX' RV
 313 :    'CLIM' 'PRES' 'TIMP' bcp 0.
 314 :           'UN' 'UIMP' CD  1. 'UN' 'VIMP' CD 0.
 315 :           'UN' 'UIMP' DA  0. 'UN' 'VIMP' DA 0.
 316 :           'UN' 'UIMP' AB  0. 'UN' 'VIMP' AB 0.
 317 :           'UN' 'UIMP' BC  0. 'UN' 'VIMP' BC 0. ;
 318 : * Le choix d'une methode iterative (Bicg stab + preconditionement MILU0)
 319 : * permet de passer des cas un peu plus gros. (Faire INFO KRES ; )
 320 :    rv . 'METHINV' . 'TYPINV'  = 3 ;
 321 :    rv . 'METHINV' . 'IMPINV'  = 0 ;
 322 :    rv . 'METHINV' . 'NITMAX'  = 400;
 323 :    rv . 'METHINV' . 'PRECOND' = 3 ;
 324 :    rv . 'METHINV' . 'RESID'   = 1.e-8 ;
 325 :    rv . 'METHINV' . 'FCPRECT' = 1 ;
 326 :    rv . 'METHINV' . 'FCPRECI' = 1 ;
 327 : 
 328 :    rv . 'INCO' = 'TABLE' 'INCO' ;
 329 :    rv . 'INCO' . 'UN'   = 'KCHT' $mt 'VECT' 'SOMMET' (0. 0.) ;
 330 :    rv . 'INCO' . 'PRES' = 'KCHT' $mt 'SCAL' 'CENTREP1' 0. ;
 331 : 
 332 :    EXEC RV ;
 333 : 
 334 : * Les resultats se trouvent dans la table rv . 'INCO'
 335 : 
 336 : *.Fin exemple II.2 .....................................................
 337 : 
 338 :    II.3/ Calcul direct d'un regime permanent.
 339 : 
 340 :    - On peut faire une recherche directe d'un regime permanent avec
 341 :    des iterations internes pour resoudre les non-linearites.
 342 : 
 343 : *.Exemple II.3  :...........................................................
 344 : * Cas de la cavite carree a paroi defilante
 345 : * !!!! ATTENTION : La cavite etant fermee (V.n impose sur toute la frontiere)
 346 : * il est necessaire d'imposer la pression en un point.
 347 : 
 348 :    ro = 400. ;
 349 :    mu = 1. ;
 350 : 
 351 :    RV = 'EQEX' 'OMEGA' 0.7 'NITER' 10 'ITMA' 0
 352 :    'OPTI' 'EF' 'IMPL' 'SUPG'
 353 :    'ZONE'  $mt  'OPER' 'LAPN' mu 'INCO' 'UN'
 354 :    'ZONE'  $mt  'OPER' 'KONV' ro 'UN' mu 'INCO' 'UN'
 355 :    'OPTI' 'EF'  'CENTREP1'
 356 :    'ZONE'  $mt  'OPER' 'KBBT' 1.         'INCO' 'UN' 'PRES'
 357 :    ;
 358 :    RV = 'EQEX' RV
 359 :    'CLIM' 'PRES' 'TIMP' bcp 0.
 360 :           'UN' 'UIMP' CD  1. 'UN' 'VIMP' CD 0.
 361 :           'UN' 'UIMP' DA  0. 'UN' 'VIMP' DA 0.
 362 :           'UN' 'UIMP' AB  0. 'UN' 'VIMP' AB 0.
 363 :           'UN' 'UIMP' BC  0. 'UN' 'VIMP' BC 0. ;
 364 : * Le choix d'une methode iterative (Bicg stab + preconditionement MILU0)
 365 : * permet de passer des cas un peu plus gros. (Faire INFO KRES ; )
 366 :    rv . 'METHINV' . 'TYPINV'  = 3 ;
 367 :    rv . 'METHINV' . 'IMPINV'  = 0 ;
 368 :    rv . 'METHINV' . 'NITMAX'  = 400;
 369 :    rv . 'METHINV' . 'PRECOND' = 3 ;
 370 :    rv . 'METHINV' . 'RESID'   = 1.e-8 ;
 371 :    rv . 'METHINV' . 'FCPRECT' = 1 ;
 372 :    rv . 'METHINV' . 'FCPRECI' = 1 ;
 373 :    rv . 'INCO' = 'TABLE' 'INCO' ;
 374 :    rv . 'INCO' .'UN'   = 'KCHT' $mt 'VECT' 'SOMMET' (0. 0.) ;
 375 :    rv . 'INCO' .'PRES' = 'KCHT' $mt 'SCAL' 'CENTREP1' 0. ;
 376 : 
 377 :    EXEC RV ;
 378 : 
 379 : * Les resultats se trouvent dans la table rv . 'INCO'
 380 : 
 381 : *.Fin exemple II.3 .....................................................
 382 : 
 383 : 
 384 : 
 385 :    II.4/ Methode de projection
 386 : 
 387 :    - On peut faire un calcul transitoire implicite en decouplant la reso-
 388 :    lution de l'equation de quantite de mouvement de la resolution de l'eq-
 389 :    uation de continuite. Le pas de temps n'est plus soumis a une contrainte
 390 :    de stabilite de type CFL. Cependant en pratique on ne peut pas prendre
 391 :    un pas de temps arbitrairement grand.
 392 :    L'algorithme necessite beaucoup moins de pas de temps que l'algorithme
 393 :    semi explicite. A l'usage il s'avere etre le plus economique pour un
 394 :    calcul transitoire voire meme un calcul permanent.
 395 : 
 396 :    - La mise en oeuvre est assez similaire a la difference que les opera-
 397 :    teurs peuvent etre implicites.
 398 : 
 399 :    La construction de la table se fait en deux etapes.
 400 :    - Une premiere etape consiste a creer une table decrivant l'equation
 401 :    de quantite de mouvement.
 402 :    - Dans la deuxieme etape on construit la table decrivant les operateurs
 403 :    lies a la projection (div U = 0)
 404 :    Les conditions de vitesse normale ou vitesse tangente en font partie,
 405 :    elles sont traitees a l'aide de multiplicateurs de Lagrange.
 406 : 
 407 :    Comme dans l'algorithme semi explicite la condition de pression imposee
 408 :    n'en fait pas partie. La pression est imposee comme les contraintes
 409 :    visqueuses sur l'equation de quantite de mouvement.
 410 : 
 411 :    - Enfin on place la deuxieme table a l'indice 'PROJ' de la premiere.
 412 : 
 413 :     L'indice rv.'CALPRE' = VRAI de la premiere table indique que l'on
 414 :    recalcule a chaque pas de temps la matrice de pression. Ceci est
 415 :    necessaire en formulation A.L.E.
 416 :     rv.'CALPRE' = FAUX ou l'absence de cet indice fait que la matrice de
 417 :    pression sera calculee une fois pour toute.
 418 : 
 419 : *.Exemple II.4  :...........................................................
 420 : * Cas de la cavite carree a paroi defilante
 421 : * !!!! ATTENTION : La cavite etant fermee (V.n impose sur toute la frontiere)
 422 : * il est necessaire d'imposer la pression en un point.
 423 : 
 424 : 
 425 :    ro = 400. ;
 426 :    mu = 1. ;
 427 :    dt = 1. ;
 428 : 
 429 :    RV = 'EQEX' 'OMEGA' 1. 'NITER' 1 'ITMA' 500 'ALFA' 0.5
 430 :    'OPTI' 'EF' 'IMPL' 'SUPG'
 431 :    'ZONE'  $mt 'OPER'  'NS'  (mu/ro) 'INCO' 'UN'
 432 :    'OPTI' 'EFM1' 'CENTREE'
 433 :    'ZONE'  $mt 'OPER' 'DFDT' 1. 'UN' dt 'INCO' 'UN'
 434 :    'CLIM' 'UN' 'UIMP' CD  1. 'UN' 'VIMP' CD 0.
 435 :           'UN' 'UIMP' DA  0. 'UN' 'VIMP' DA 0.
 436 :           'UN' 'UIMP' AB  0. 'UN' 'VIMP' AB 0.
 437 :           'UN' 'UIMP' BC  0. 'UN' 'VIMP' BC 0. ;
 438 : * Le choix d'une methode iterative (Bicg stab + preconditionement MILU0)
 439 : * permet de passer des cas un peu plus gros. (Faire INFO KRES ; )
 440 :     rv . 'METHINV' . 'TYPINV'  = 3 ;
 441 :     rv . 'METHINV' . 'IMPINV'  = 0 ;
 442 :     rv . 'METHINV' . 'NITMAX'  = 400 ;
 443 :     rv . 'METHINV' . 'PRECOND' = 3 ;
 444 :     rv . 'METHINV' . 'RESID'   = 1.e-8 ;
 445 :     rv . 'METHINV' . 'FCPRECT' = 1 ;
 446 :     rv . 'METHINV' . 'FCPRECI' = 1 ;
 447 : 
 448 :     betastab=1.e2 ;
 449 : 
 450 :     RVP = 'EQEX'
 451 :     'OPTI' 'EF' 'CENTRE'
 452 :     'ZONE'  $mt 'OPER' 'KBBT' -1. betastab 'INCO' 'UN' 'PRES'
 453 :     'CLIM' 'PRES' 'TIMP' bcp 0. ;
 454 : 
 455 :     rvp . 'METHINV' . 'TYPINV'  = 2 ;
 456 :     rvp . 'METHINV' . 'IMPINV'  = 0 ;
 457 :     rvp . 'METHINV' . 'NITMAX'  = 300;
 458 :     rvp . 'METHINV' . 'PRECOND' = 3 ;
 459 :     rvp . 'METHINV' . 'RESID'   = 1.e-8 ;
 460 :     rvp . 'METHINV' . 'FCPRECT' = 100 ;
 461 :     rvp . 'METHINV' . 'FCPRECI' = 100 ;
 462 : 
 463 :     rv . 'PROJ' = RVP ;
 464 : 
 465 :     rv . 'INCO' = 'TABLE' 'INCO' ;
 466 :     rv . 'INCO' . 'UN'   = 'KCHT' $mt 'VECT' 'SOMMET' (0. 0.) ;
 467 :     rv . 'INCO' . 'PRES' = 'KCHT' $mt 'SCAL' 'CENTRE' 0. ;
 468 : 
 469 :     EXEC RV ;
 470 : 
 471 : * Les resultats se trouvent dans la table rv . 'INCO'
 472 : 
 473 : *.Fin exemple II.4 .....................................................
 474 : 
 475 :     Remarques :
 476 :     ___________
 477 : 
 478 :     1) Le positionnement de la variable rv.'DETMAT' a VRAI indique que
 479 :        les objets MATRIK seront detruits a la fin de la procedure.
 480 : 
 481 :     2) Deux variables de type LOGIQUE resp. rv.'STOPITER' et 
 482 :        rv.'STOPPDT' initialisees a FAUX peuvent etre mises a VRAI (par
 483 :        exemple a l'aide d'une procedure utilisateur) afin de stopper
 484 :        resp. la boucle de resolution de la non-linearite et la boucle
 485 :        sur les pas de temps. 
 486 :        L'indice rv . 'NUITER' contient le numero de l'iteration non
 487 :        lineaire courante.
 488 : 
 489 :     3) Il est possible d'utiliser une methode de projection algebrique
 490 :        incrementale pour resoudre le systeme NS incompressible de 
 491 :        maniere approchee en gardant la meme syntaxe qu'en calcul
 492 :        transitoire implicite et en rajoutant une table a l'indice 
 493 :        'GPROJ' :
 494 :         rv . 'GPROJ' = 'TABLE' ;
 495 :         rv . 'GPROJ' . 'NOMVIT'  = 'UN' ;
 496 :         rv . 'GPROJ' . 'NOMPRES' = 'CHAINE' 'PRES' ;
 497 :        ou 'UN' est le nom de l'inconnue vitesse et 'PRES' le nom 
 498 :        de l'inconnue pression.
 499 :        Il est egalement possible de donner une table a l'indice :
 500 :         rv . 'GPROJ' . 'METHINV'
 501 :        precisant les options pour l'inversion de la matrice de pression.
 502 :        (cf. notice KRES)
 503 :        On peut également préciser si on veut une méthode de simple ou double  
 504 :        projection :
 505 :        rv . 'GPROJ' . 'dblproj' = FAUX ou VRAI (VRAI par défault)
 506 : 

© Cast3M 2003 - Tous droits réservés.
Mentions légales