* fichier : mdiavf2.dgibi * * Utilisation de l'operateur MDIA pour creér un jacobien analytique * dans le cas d'injection subsonic d'air dans l'air * * A. Beccantini, SFME/LTMF * GRAPH = VRAI ; GRAPH = FAUX ; ****************** **** MAILLAGE **** ****************** RAF = 8 ; L = 10. ; DX = 1.0 '/' RAF ; NX = 2 ; LINJ = NX '*' DX ; P2 = LINJ 0.0 ; P3 = 0.0 0.0 ; NY = 'ENTIER' ((L '/' DX) '+' 0.5) ; P2P3 = P2 'DROIT' NX P3 ; DOMINT = P2P3 'TRANSLATION' NY (0. L) ; * * Creation of MAILEN (maillage à l'entrée) * LIGAPP = ('CHANGER' 'POI1' P2P3) ; * **** Creation of DOMAINE tables via the MODEL object * MDOMINT = 'MODELISER' DOMINT 'EULER' ; MP2P3 = 'MODELISER' P2P3 'EULER' ; MMAILEN = 'MODELISER' MAILEN 'EULER' ; 'ELIMINATION' (QDOMINT 'ET' QP2P3 'ET' QMAILEN) 0.001 ; 'SI' GRAPH ; 'TRACER' (DOMINT 'ET' (MAILEN 'COULEUR' 'ROUG')) 'TITR' 'Maillage entree'; 'FINSI' ; * * MAIFA1 = internel cells close to ghost cells * MAIFA1 = 'MANUEL' 'POI1' POIN1 ; MAIFA1 = MAIFA1 'ET' ('MANUEL' 'POI1' POIN1) ; ; ; 'FIN' BLLIM ; 'SI' GRAPH ; 'TRACER' (DOMINT 'ET' MAIFA1) 'TITRE' 'Domain and faces -> centres' ; 'FINSI' ; ************** **** I.C. **** ************** ($DOMINT . 'CENTRE')) 'NATU' 'DISCRET') 'ET' ($DOMINT . 'CENTRE')) 'NATU' 'DISCRET') 'ET' ($DOMINT . 'CENTRE')) 'NATU' 'DISCRET') 'ET' ($DOMINT . 'CENTRE')) 'NATU' 'DISCRET') ; ********************************************************************* **** Procedure RESL ************************************************* ********************************************************************* * 'DEBP' RESL ; * **** Euler Monocomponent 2D * Subsonic injection of air * * A definir: * * RVX . 'DOMA' = table domaine des "internal cells close to the * injection interfaces" * RVX . 'DOMA' ne serve que pour MDIA * * RVX . 'VOLUME' = volume des "internal cells close to the injection * interfaces" * * RVX . 'SURFACE'= CHPOINT defini sur les centres des "internal cells * close to the injection interfaces", qui contient * la surface des faces d'injection * |-------| * | | * +---> | * | | * |-------| * * RVX . 'MAILLAGE' = maillage de 'POI1' (centres des "internal cells * close to the injection interfaces") * * RVX . 'GAMLIM' (objet de type FLOTTANT) (gamma du gaz); * RVX . 'ROVINF' (objet de type FLOTTANT) (flux de masse à l'injection) ; * RVX . 'RTINF' (objet de type FLOTTANT) (temperature à l'injection); * * RVX . 'DX' (objet de type FLOTTANT) (longuer pour le calcul du pas * de temps de stabilité pour l'algorithme temporel * explicite) * RVX . 'ALPHA' (objet de type FLOTTANT) (coefficient de securité pour * le calcul du pas de temps de stabilité pour * l'algorithme temporel explicite) * CHVOLUME = RVX . 'VOLUME' ; CHSFACE = RVX . 'SURFACE' ; 'SI' ('NEG' AA 'FLOTTANT') ; 'MESSAGE' 'RV . GAMLIM ???' ; 'ERREUR' 21 ; 'FINSI' ; 'SI' ('NEG' AA 'FLOTTANT') ; 'MESSAGE' 'RV . ROVINF ???' ; 'ERREUR' 21 ; 'FINSI' ; 'SI' ('NEG' AA 'FLOTTANT') ; 'MESSAGE' 'RV . RTINF ???' ; 'ERREUR' 21 ; 'FINSI' ; * Inconnue UN = 'COPIER' (RV . 'UN') ; MAIFA1 = RVX . 'MAILLAGE' ; GAMSCAL = (RVX . 'GAMLIM') ; * LISTREEL pour les conditions limites (flux) . * LISTREEL pour les conditions limites (jacobien) . * Boucle sur maifa1 pour calculer les conditions à imposer RTG = RVX . 'RTINF' ; GG = RVX . 'ROVINF' ; 'REPETER' BL1 NN ; SURCEL = 'EXTRAIRE' CHSFACE 'SCAL' PCEL ; VOLCEL = 'EXTRAIRE' CHVOLUME 'SCAL' PCEL ; ROD = 'EXTRAIRE' RNRED 'SCAL' PCEL ; PD = 'EXTRAIRE' PNRED 'SCAL' PCEL ; UYD = 'EXTRAIRE' VNRED 'UY' PCEL ; UXD = 'EXTRAIRE' VNRED 'UX' PCEL ; UND = UYD ; UTD = UXD ; PG = PD ; ROG = PG '/' RTG ; UNG = GG '/' ROG ; RHTG = ((PG '*' GAMSCAL) '/' (GAMSCAL '-' 1.)) '+' (0.5 '*' ROG '*' UNG '*' UNG) ; * DROUYDP = d(RHOUY)/dP_(primitive) DROUYDP = 1. '-' ((GG '*' GG '*' RTG) '/' (PD '*' PD)) ; * DRETDP = d(RET)/dP_(primitive) DRETDP = -1 '*' ((GG '*' GG '*' GG '*' RTG '*' RTG) '/' (PD '*' PD '*' PD)) ; * * d/dro_(conservative) = DPDRO '*' d()/dP_(primitive) DPDRO = (GAMSCAL '-' 1.) '*' ((UXD ** 2) '+' (UYD ** 2)) '/' 2. ; * * d/droux_(conservative) = DPDUX '*' d()/dP_(primitive) DPDUX = -1 '*' (GAMSCAL '-' 1.) '*' UXD ; * * d/drouy_(conservative) = DPDUY '*' d()/dP_(primitive) DPDUY = -1 '*' (GAMSCAL '-' 1.) '*' UYD ; * * d/dret_(conservative) = DPDRET '*' d()/dP_(primitive) DPDRET = (GAMSCAL '-' 1.) ; * * * Le flux * * ^ * |--> * -------------- * (SURCEL '/' VOLCEL))) ; * * Le jacobien du flux * * Inverse du volume est lié a la syntaxe de MDIA * VOLCEL2 = VOLCEL '*' VOLCEL ; (SURCEL '/' VOLCEL2))) ; (SURCEL '/' VOLCEL2))) ; (SURCEL '/' VOLCEL2))) ; (SURCEL '/' VOLCEL2))) ; (SURCEL '/' VOLCEL2))) ; (SURCEL '/' VOLCEL2))) ; (SURCEL '/' VOLCEL2))) ; (SURCEL '/' VOLCEL2))) ; 'FIN' BL1 ; 'RUY' LRUN 'RETN' LRET 'NATU' 'DISCRET' ; UPC = ((GAMSCAL '*' PG '/' ROG) '**' 0.5) '+' UNG ; THETA = RVX . 'THETA' ; 'SI' (RVX . 'IMPL') ; * **** Jacobian avec MDIA * * La matrice met en relation les couples suivantes: * PRIMALES DUALES *1 densité (q.d.m.)y *2 densité energie totale * *3 (q.d.m.)x (q.d.m.)y *4 (q.d.m.)x energie totale * *5 (q.d.m.)y (q.d.m.)y *6 (q.d.m.)y energie totale * *7 energie totale (q.d.m.)y *8 energie totale energie totale * * **** Test of inlet implicit boundary conditions * * La matrice met en relation les couples suivantes: * PRIMALES DUALES *1 densité (q.d.m.)y *2 densité energie totale * *3 (q.d.m.)x (q.d.m.)y *4 (q.d.m.)x energie totale * *5 (q.d.m.)y (q.d.m.)y *6 (q.d.m.)y energie totale * *7 energie totale (q.d.m.)y *8 energie totale energie totale * * *1 'ZONE' $DOMA 'OPER' 'MDIA' CHRUYR 'INCO' 'RN' 'RUY' *2 'ZONE' $DOMA 'OPER' 'MDIA' CHRETR 'INCO' 'RN' 'RETN' *3 'ZONE' $DOMA 'OPER' 'MDIA' CHRUYRUX 'INCO' 'RUX' 'RUY' *4 'ZONE' $DOMA 'OPER' 'MDIA' CHRETRUX 'INCO' 'RUX' 'RETN' *5 'ZONE' $DOMA 'OPER' 'MDIA' CHRUYRUY 'INCO' 'RUY' 'RUY' ; * *6 'ZONE' $DOMA 'OPER' 'MDIA' CHRETRUY 'INCO' 'RUY' 'RETN' *7 'ZONE' $DOMA 'OPER' 'MDIA' CHRUYRET 'INCO' 'RETN' 'RUY' *8 'ZONE' $DOMA 'OPER' 'MDIA' CHRETRET 'INCO' 'RETN' 'RETN' ; RVDIA . 'INCO' = 'TABLE' 'INCO' ; MAT1 = RVDIA . '1MDIA' . 'MATELM' ; MAT2 = RVDIA . '2MDIA' . 'MATELM' ; MAT3 = RVDIA . '3MDIA' . 'MATELM' ; MAT4 = RVDIA . '4MDIA' . 'MATELM' ; MAT5 = RVDIA . '5MDIA' . 'MATELM' ; MAT6 = RVDIA . '6MDIA' . 'MATELM' ; MAT7 = RVDIA . '7MDIA' . 'MATELM' ; MAT8 = RVDIA . '8MDIA' . 'MATELM' ; IJACO = MAT1 'ET' MAT2 'ET' MAT3 'ET' MAT4 'ET' MAT5 'ET' MAT6 'ET' MAT7 'ET' MAT8 ; 'SINON' ; 'FINSI' ; 'RESPRO' IJACO RESLIM IALPDT THETA ; 'FINP' ; ********************************************************************* **** Fin procedure RESL ********************************************* ********************************************************************* ********************************************************************* **** Procedure JACOVA *********************************************** ********************************************************************* * Derivé partielle du residu en un point par rapport * aux variable dans le meme point * * PPRIM = point ou est localisé la variable primale * PPRIM = point ou est localisé la variable duale * MOTPRI = nom de la composante concernante la variable primale * MOTDUA = nom de la composante concernante la variable duale ELT1 = 'MANUEL' 'POI1' PPRIM ; 'REPETER' BL1 NDIM ; MOTCEL = 'EXTRAIRE' LISTINCO &BL1 ; 0.0 'NATURE' 'DISCRET') ; 'FIN' BL1 ; SCAL = 'EXTRAIRE' D_DMOT PPRIM MOTDUA ; 'FINPROC' SCAL ; ********************************************************************* **** Fin procedure JACOVA ******************************************* ********************************************************************* ********************************************************************* **** Procedure JACNUM *********************************************** ********************************************************************* 'DEBPROC' JACNUM RVX*'TABLE' LISTINCO*'LISTMOTS' PPRIM*'POINT' * PPRIM = point ou est localisé la variable primale * PPRIM = point ou est localisé la variable duale * MOTPRI = nom de la composante concernante la variable primale * MOTDUA = nom de la composante concernante la variable duale * Le valeur dans l'état non-perturbé en PPRIM ; UN0 = 'COPIER' (RV . 'UN') ; THETA0 = RVX . 'THETA' ; RVX . 'THETA' = 0.0 ; IMPL0 = RVX . 'IMPL' ; RVX . 'IMPL' = FAUX ; IJACEL CHPRES0 DT ALPDT = RESL RVX ; VAL0 = 'EXTRAIRE' CHPRES0 PPRIM MOTDUA ; * EPSILON = perturbation * On etabli la variable à perturber 'REPETER' BL1 NDIM ; MOTCEL = 'EXTRAIRE' LISTINCO &BL1 ; 'SI' ('EGA' MOTCEL MOTPRI) ; ICEL = &BL1 ; 'QUITTER' BL1 ; 'FINSI' ; 'FIN' BL1 ; 'SI' (ICEL > NDIM) ; 'MESSAGE' 'Procedure JACNUM' ; 'MESSAGE' 'MOTPRI = ??? '; 'ERREUR' 21 ; 'FINSI' ; ELT1 = 'MANUEL' 'POI1' PPRIM ; 'RN' 1.0 'RUX' 1.0 'RUY' 1.0 'RETN' 1.0 'NATU' 'DISCRET') '+' MOTPRI EPSILON) ; RV . 'UN' = UN0 '*' CHPERT LISTINCO LISTINCO LISTINCO ; DELTATOT = ('EXTRAIRE' (RV . 'UN') PPRIM MOTPRI) '-' ('EXTRAIRE' UN0 PPRIM MOTPRI) ; IJACEL CHPRES0 DT ALPDT = RESL RVX ; VAL1 = 'EXTRAIRE' CHPRES0 PPRIM MOTDUA ; RV . 'UN' = UN0 ; RVX . 'THETA' = THETA0 ; RVX . 'IMPL' = IMPL0 ; 'FINPROC' ((VAL1 '-' VAL0) '/' DELTATOT) ; ********************************************************************* ***** Fin procedure JACNUM ****************************************** ********************************************************************* RV = 'TABLE' ; RV . 'UN' = UN ; * **** Table pour RESL * RV . '3RESL' = 'TABLE' ; RV . '3RESL' . 'THETA' = 1.0 ; RV . '3RESL' . 'IMPL' = VRAI ; RV . '3RESL' . 'VOLUME' = CHVOLU ; RV . '3RESL' . 'SURFACE' = CHSURF ; RV . '3RESL' . 'MAILLAGE' = MAIFA1 ; RV . '3RESL' . 'GAMLIM' = 1.4 ; RV . '3RESL' . 'ROVINF' = 1.0 ; RV . '3RESL' . 'RTINF' = (288. '*' 600.) ; RV . '3RESL' . 'ALPHA' = 1.0 ; ******************** **** Test case ***** ******************** IJAC IRES IALP THE = RESL (RV . '3RESL') ; DELTA = 1.0D-4 ; 'REPETER' BL1 4 ; 'REPETER' BL2 4 ; MOT1 = 'EXTRAIRE' LISTINC1 &BL1 ; MOT2 = 'EXTRAIRE' LISTINC1 &BL2 ; CACCA1 = JACOVA IJAC $DOMINT LISTINC1 POI1 MOT1 MOT2 ; CACCA2 = JACNUM (RV . '3RESL') LISTINC2 POI1 MOT1 MOT2 DELTA ; 'MESSAGE' ('CHAINE' 'Analyt.= ' CACCA1 ' Num.= ' CACCA2 ' Err = ' ('ABS' (CACCA1 '-' CACCA2))) ; 'SI' (('ABS' CACCA1) > 1.0D-4) ; ERRO = 'ABS' ((CACCA1 '-' CACCA2) '/' CACCA1) ; 'SINON' ; ERRO = 'ABS' (CACCA1 '-' CACCA2) ; 'FINSI' ; 'SI' (ERRO > 5.0D-4) ; 'ERREUR' 5 ; 'FINSI' ; 'FIN' BL2 ; 'FIN' BL1 ; 'FIN' ;
© Cast3M 2003 - Tous droits réservés.
Mentions légales