* fichier : sine_bumpBM.dgibi ****************************************************************** * CALCUL DE L'ECOULEMENT SUBSONIQUE STATIONNAIRE DANS UN CANAL * * AVEC SINE-SHAPED BUMP * * FORMULATION VF COMPRESSIBLE EXPLICITE/IMPLICIT * * * * H. PAILLERE/P. GALON TTMF AOUT 1997 * * * * MODIF BECCANTINI MARS 97 * * * * BECCANTINI NOVEMBRE 98 * * MODIF BECCANTINI AUGUST 2021 * * * * Compressible AUSM+ scheme is used for computation * ****************************************************************** GRAPH = VRAI ; COMPLET = VRAI ; GRAPH = FAUX ; COMPLET = FAUX ; 'SI' complet ; RAF = 8 ; NITER = 20 ; 'SINON' ; RAF = 5 ; NITER = 20 ; 'FINSI' ; TYEL = 'QUA4' ; 'OPTION' 'DIME' 2 'ELEM' TYEL 'ISOV' 'SULI' 'ECHO' 1 'TRAC' 'X'; * *** C.L. et initiales * RO_INF = 1.4D0 ; P_INF = 1.0D0 ; U_INF = 0.50 ; GAMSCAL = 1.4D0; * * ************************************************************************ ************************************************************************ ***************** PARTIE PROCEDURES ************************************ ************************************************************************ ************************************************************************ ************************************************************************ * * * Put it in utilproc * Procedures will finish in FIN PARTIE PROCEDURES * ***** $$$$ EXEXCH * ***************************************************** ***************************************************** ** PROCEDURE EXEX POUR FORMULATION VF COMPRESSIBLE ** ***************************************************** ***************************************************** * * We check wether the table RV is complete * * 'DEBPROC' EXEXCH ; 'ARGUMENT' RV*TABLE ; 'MESSAGE' ; 'MESSAGE' 'PROCEDURE EXEXCH' ; * **** Initialisation d'un CHPOINT 'ET' d'une MATRIK vide * CHPVID MATVID = 'KOPS' 'MATRIK' ; RESPRO CHPVID MATVID ; * **** Les inconnues * 'SI' ('NON' ('EXISTE' RV 'UN')) ; 'MESSAGE' ; 'MESSAGE' 'UN = ???' ; 'ERREUR' 21 ; 'FINSI' ; * **** We check the compatibility between the variables and their names * 'SI' ('NON' ('EXISTE' RV 'LISTCONS')) ; 'MESSAGE' ; 'MESSAGE' 'LISTCONS = ???' ; 'ERREUR' 21 ; 'FINSI' ; 'MESSAGE' ; 'MESSAGE' 'Conservative variables check' ; UNCELL = 'EXCO' (RV . 'LISTCONS') (RV . 'UN') (RV . 'LISTCONS') ; ERRO = 'MAXIMUM' (UNCELL '-' (RV . 'UN')) 'ABS' ; ERRO = ERRO '/' (('MAXIMUM' UNCELL 'ABS') '+' 1.0D-15) ; 'SI' (ERRO > 1.0D-6) ; 'MESSAGE' ; 'MESSAGE' 'UN = ???' ; 'MESSAGE' 'LISTCONS = ???' ; 'ERREUR' 21 ; 'SINON' ; 'MESSAGE' 'OK' ; 'FINSI' ; * **** On which variables the error has to be computed? * 'SI' ('NON' ('EXISTE' RV 'ERROMOTS')) ; 'MESSAGE' ; 'MESSAGE' 'ERROMOTS = ???' ; 'ERREUR' 21 ; 'FINSI' ; * **** Definition of the table containing some results * 'SI' ('NON' ('EXISTE' RV 'RESULTATS')) ; RV . 'RESULTATS' = 'TABLE' ; 'FINSI' ; * **** Existence d'une procédure pour imposer le conditions aux limites * à chaque iteration (interne ou externe) * 'SI' ('NON' ('EXISTE' RV 'CLIM')) ; 'MESSAGE' ; 'MESSAGE' 'CLIM = ???' ; 'ERREUR' 21 ; 'FINSI' ; * Par securité, on les impose tout de suite 'SI' (RV . 'CLIM') ; 'SI' ('NON' ('EXISTE' RV 'MAIFAN')) ; 'MESSAGE' ; 'MESSAGE' 'MAIFAN = ???' ; 'ERREUR' 21 ; 'FINSI' ; RV . 'UN' = PROLIM RV ; * DUCLIM à imposer dans l'inversion matricielle: * pas d'increment sur les cellules fantomes DUCLIM = 0.0 '*' ('REDU' (RV . 'UN') (RV . 'MAIFAN')) ; 'SINON' ; DUCLIM = 'COPIER' CHPVID ; 'FINSI' ; RESPRO DUCLIM ; * * ****** Mise a jour de la matrice à inverser pendant les iterations internes? * * RV . 'MATUPDAT' * RV . 'MUPINT' * 'SI' ('NON' ('EXISTE' RV 'MATUPDAT')) ; 'MESSAGE' ; 'MESSAGE' 'MATUPDAT = ???' ; 'ERREUR' 21 ; 'FINSI' ; 'SI' (RV . 'MATUPDAT') ; 'MESSAGE' ; 'MESSAGE' 'We always update the matrix to inverse' ; LOGMAOLD = FAUX ; 'SINON' ; 'MESSAGE' ; 'MESSAGE' 'We don t always update the matrix to inverse' ; LOGMAOLD = VRAI ; 'FINSI' ; * RESPRO LOGMAOLD ; * **** Iterations externes et/ou temps final * La table 'PASDETPS' * 'SI' ('NON' (('EXISTE' RV 'NITMAEX') 'OU' ('EXISTE' RV 'TFINAL'))) ; 'MESSAGE' ; 'MESSAGE' 'NITMAEX = ???' ; 'MESSAGE' 'ou' ; 'MESSAGE' 'TFINAL = ???' ; 'ERREUR' 21 ; 'FINSI' ; * **** Gas model * 'SI' ('NON' ('EXISTE' RV 'PGAZ')) ; 'MESSAGE' ; 'MESSAGE' 'PGAZ = ???' ; 'ERREUR' 21 ; 'FINSI' ; * **** BAS MACH * 'SI' ('NON' ('EXISTE' RV 'BASMACH')) ; 'MESSAGE' ; 'MESSAGE' 'BASMACH = ???' ; 'ERREUR' 21 ; 'SINON' ; 'SI' (RV . 'BASMACH') ; 'SI' (('NON' ('EXISTE' RV 'CO1')) 'OU' ('NON' ('EXISTE' RV 'CO2'))) ; 'MESSAGE' ; 'MESSAGE' 'CO1 = ???' ; 'MESSAGE' 'CO2 = ???' ; 'ERREUR' 21 ; 'FINSI' ; 'FINSI' ; 'FINSI' ; * *** La table PASDETPS * 'SI' ('NON' ('EXISTE' RV 'PASDETPS')) ; 'MESSAGE' ; 'MESSAGE' 'We create the table PASDETPS' ; RV . 'PASDETPS' = 'TABLE' ; 'SINON' ; 'MESSAGE' ; 'MESSAGE' 'We use the table PASDETPS that already exists' ; 'FINSI' ; KTPS = RV . 'PASDETPS' ; * RESPRO KTPS ; * Initialisation éventuelle de la table 'SI' (('NON' ('EXISTE' KTPS 'NUPASDT')) 'OU' ('NON' ('EXISTE' KTPS 'TPS')) 'OU' ('NON' ('EXISTE' KTPS 'TPSM'))) ; 'MESSAGE' ; 'MESSAGE' 'Table PASDETPS' ; 'MESSAGE' 'NUPASDT = 0' ; 'MESSAGE' 'TPS = 0.0' ; 'MESSAGE' ; KTPS . 'NUPASDT' = 0 ; KTPS . 'TPSM' = 0.0D0 ; KTPS . 'TPS' = 0.0D0 ; 'FINSI' ; * * KTPS contient * * KTPS . 'NUPASDT' = numero de pas de TPS actuel dans les itérations * internes * * KTPS . 'TPSM' = le TPS aprés (KTPS . 'NUPASDT' '-' 2) itérations * KTPS . 'TPS' = le TPS aprés (KTPS . 'NUPASDT' '-' 1) itérations * KTPS . 'TPSP' = le TPS aprés (KTPS . 'NUPASDT') itérations, i.e. à * la fin de l'itération actuelle * * * * N.B. 'XINIT' is used in the iterative methods * We have decided to keep it 0 RV . 'MATINV' . 'XINIT' = 0.0 '*' (RV . 'UN') ; 'FINPROC' ; **** $$$$ EXEXIM * ***************************************************** ***************************************************** ** PROCEDURE EXEX POUR FORMULATION VF COMPRESSIBLE ** ***************************************************** ***************************************************** * * RV . 'UN' = les variables conservatives * * RV . 'ANOM' = LOGIQUE (anomalie detectée ?) * * RV . 'DTIMP' = pas de tmps imposé (facultatif) * * RV . 'LISTCONS'= noms des variables conservatives * * RV . 'ERROMOTS'= noms des variables sur lesquelles on calcule Linf * * RV . 'CLIM' = logique qui me dit si existe une procédure pour le * calcul de conditions limites (procédure PROLIM) * * RV . 'MAIFAN' = (a définir si RV . 'CLIM') * les cellules fantômes * * RV . 'TFINAL' = le temps final * * RV . 'NITMAEX' = le nombre d'itération externes * * N.B. Si RV . 'TFINAL' et RV . 'NITMAEX' sont les deux spécifiés, on * s'arrête quand un des deux critères est satisfait * * RV . 'MATUPDAT' = si VRAI, on ne met pas toujours à jour la matrice * à inverser * * RV . 'MUPEXT' = si (RV . 'MATUPDAT'), indice de mise à jours de la * matrice jacobienne pendant les itérations externes * * RV . 'MUPINT' = si (RV . 'MATUPDAT'), indice de mise à jours de la * matrice jacobienne pendant les itérations internes * * RV . 'MUPINI' = si (RV . 'MATUPDAT'), nombre initial de pas de temps * pendant lesquels on met à jour la matrice à chaque * itération externe * * RV . 'MUPLIN' = si (RV . 'MATUPDAT') et on utilise un solveur * itératif pour la résolution du système linéaire, * si le nombre d'itérations linéaires est > que * RV . 'MUPLIN', alors on met à jour le * preconditionneur * * RV . 'MATACC' = si (RV . 'MATUPDAT'), logique qui nous dit si on * veut utiliser la méthode de Broyden. * * RV . 'LISTOPER' = liste des opérateurs (ou des procédures) qui * interviennent dans le calcul (chaque opérateur a une * table associée), qui s'appelle &NOMPER ou * & = position de l'opérateur dedans cette liste * NOMPER = noms de l'opérateur ou de la procédure * * RV . 'NITMAIN' = (à définir dans le cas d'un calcul implicite) * le nombre max d'itération internes * RV . 'NITMIIN' = (à définir dans le cas d'un calcul implicite) * le nombre min d'itération internes * * RV . 'EPSINT' = (à définir dans le cas d'un calcul implicite) * l'erreur pour le critère de convergence sur les * itérations internes * * RV . 'MATHINV' = (a définir sans la cas d'un calcul implicite) * table de SOUSTYPE 'TYPINV' pour l'inversion de * MATRIK (pour l'opérateur 'KRES') * * RV . 'RESULTATS' = table qui contient des resultats * 'NITERLIN' = nombre de iterations lineaires * dans le solveur iteratif * 'ERROLIN' = residu de l'erreur du solveur * iteratif * 'ITERIN' = iteration interne (0, 1, ...) * ************************************************************************* * MODIF ************************************************************************* * Message en anglais! * 'DEBPROC' EXEXIM ; 'ARGUMENT' RV*TABLE ; 'MESSAGE' ; 'MESSAGE' 'PROCEDURE EXEXIM' ; * **** Initialisation and controls into exexch * CHPVID MATVID DUCLIM LOGMAOLD KTPS = EXEXCH RV ; ************************************************************************* ************************************************************************* ************************************************************************* ************************************************************************* ************** ITÉRATIONS EXTERNES ************************************ ************************************************************************* ************************************************************************* ************************************************************************* ************************************************************************* LOGEXP = VRAI ; * LOGEXP = variable logique qui me dit si on est en explicite ; LOGQIE = FAUX ; * **** Boucle qui s'arrête quand LOGQIE = VRAI ; i.e. * * (KTPS . NUPASDT) = (RV . 'NITMAEX') * ou * (KTPS . 'TPS') = (RV . 'TFINAL'); * * 'SI' ('EXISTE' RV 'DTIMP') ; ALPDT = RV . 'DTIMP' ; 'SINON' ; ALPDT = 0.0 ; 'FINSI' ; NITEX = 0 ; 'REPETER' BLEX -1 ; NITEX = NITEX '+' 1 ; * ****** Iteration interne * RV . 'RESULTATS' . 'ITERIN' = 0 ; * ****** Mise à jour de la matrice à inverser * 'SI' ((NITEX 'EGA' 1) 'OU' (RV . 'ANOM')) ; * A la premiere iteration externe on doit calculer la matrice LOGUPEX = VRAI ; 'SINON' ; 'SI' LOGMAOLD ; MUPEXT = RV . 'MUPEXT' ; * On met a jour la matrice toutes les MUPEXT iterations externes LOGUPEX = ((NITEX '/' MUPEXT) '*' MUPEXT) 'EGA' NITEX ; 'FINSI' ; 'FINSI' ; * ****** Trop d'iterations externes? * KTPS . 'NUPASDT' = KTPS . 'NUPASDT' '+' 1 ; 'SI' ('EXISTE' RV 'NITMAEX') ; 'SI' ( (KTPS . 'NUPASDT') '>EG' (RV . 'NITMAEX')) ; LOGQIE = VRAI ; 'FINSI' ; 'FINSI' ; * ****** Impression * 'MESSAGE' ; 'MESSAGE' ('CHAINE' 'PASDETPS = ' (KTPS . 'NUPASDT') ' TPS = ' (KTPS . 'TPS') ' DTM1 = ' ALPDT) ; * *** ALGORITHM A RESOUDRE * * F(U^{n+1}) = 0 * * if (RV . INST) then * * F(U^{n+1}) = -1. * (U^{n+1} - U^{n}) '/' (\delta t) * '+' \sum_k Res_k(U^{n+1}) '+' \sum_kexpl Res_kexpl(U^{n}) * * else * * F(U^{n+1}) = -\sum_k Res_k(U^{n+1}) * * endif * * * Avec une methode de type Newton-Raphson * * a) U^{n+1,0} = U^{n} * * b) for l=0,1,2,... * * MAT(U^{n+1,l}) \delta U^{n+1,l} = -F(U^{n+1,l}) * * U^{n+1,l+1} = U^{n+1,l} '+' \delta U^{n+1,l} * * if || \delta U^{n+1,l} ||_{\inf} < \epsilon goto c) * * c) * * U^{n+1} = U^{n+1,l+1} * * endif * *** LES VARIABLES * * RESIMP = \sum_k Res_k(U^{n+1,l}) * RESEXP = \sum_kexpl Res_kexpl(U^{n}) * MATASS = - \sum_k JACR_k(U^{n+1,l}) * UNM = conserved variables at t^{n} * UNM = conserved variables at t^{n-1} * LREEDT = LISTREEL de (\delta t)_k t.q. * \delta t = min (\delta t)_k * RESIMP = 'COPIER' CHPVID ; RESEXP = 'COPIER' CHPVID ; 'SI' ('EXISTE' RV 'UNM2') ; RV . 'UNM3' = 'COPIER' (RV . 'UNM2') ; 'FINSI' ; 'SI' ('EXISTE' RV 'UNM') ; RV . 'UNM2' = 'COPIER' (RV . 'UNM') ; 'FINSI' ; RV . 'UNM' = 'COPIER' (RV . 'UN') ; MATASS = 'KOPS' 'MULT' 0.0 MATVID ; LREEDT = 'PROG' ; * *********************************************** ********* Boucle sur les operateurs *********** *********************************************** *** On calcule: LREEDT * RESIMP * RESEXP * (MATASS) * 'REPETER' BLOP ('DIME' (RV . 'LISTOPER')) ; NOMPER = 'EXTRAIRE' &BLOP (RV . 'LISTOPER') ; NOTABLE = 'MOT' ('TEXTE' ('CHAINE' &BLOP NOMPER) ) ; 'SI' (RV . NOTABLE . 'IMPL') ; 'SI' ('EGA' NITEX 1) ; 'MESSAGE' ; 'MESSAGE' ('CHAINE' ('MOT' NOMPER) ': implicit') ; 'FINSI' ; * JACO = objet de type MATRIK (éventuellement vide) * RESIDU = " RESIDU " * ALPHADT = " LISTREEL " LOGEXP = FAUX ; JACO RESIDU ALPHADT = ('TEXTE' NOMPER) (RV . NOTABLE) ; 'SI' (RV . NOTABLE . 'ANOM') ; RV . 'ANOM' = VRAI ; RV . NOTABLE . 'ANOM' = FAUX ; 'QUITTER' BLOP ; 'SINON' ; MATASS = MATASS 'ET' ('KOPS' 'MULT' -1.0D0 JACO ) ; RESIMP = RESIMP 'ET' RESIDU ; 'FINSI' ; 'SINON' ; 'SI' ('EGA' NITEX 1) ; 'MESSAGE' ; 'MESSAGE' ('CHAINE' ('MOT' NOMPER) ': explicit') ; 'FINSI' ; RESIDU ALPHADT = ('TEXTE' NOMPER) (RV . NOTABLE) ; 'SI' (RV . NOTABLE . 'ANOM') ; RV . 'ANOM' = VRAI ; 'QUITTER' BLOP ; 'SINON' ; RESEXP = RESEXP 'ET' RESIDU ; 'FINSI' ; 'FINSI' ; LREEDT = LREEDT 'ET' ALPHADT ; 'FIN' BLOP ; *********************************************** ***** Fin de la boucle sur les operateurs **** *********************************************** * ******* Anomalie detectée (CTRL+S 9999) * On arrete ici * 'SI' ('NON' (RV . 'ANOM')) ; * ******* On controlle la compatibilité E/S (si NITEX = 1) * 'SI' ('EGA' NITEX 1) ; * **** Dans le cas implicite, on verifie l'existence * des parametres pour les itérations internes * 'SI' ('NON' LOGEXP) ; 'SI' ('NON' (('EXISTE' RV 'NITMAIN') 'ET' ('EXISTE' RV 'EPSINT') 'ET' ('EXISTE' RV 'NITMIIN'))); 'MESSAGE' ; 'MESSAGE' 'NITMAIN = ??? ' ; 'MESSAGE' 'NITMIIN = ??? ' ; 'MESSAGE' 'EPSINT = ??? ' ; 'ERREUR' 21 ; 'FINSI' ; 'FINSI' ; 'FINSI' ; * ******* Fin contrôle compatibilité E/S * *** Mise a jour de la table RV . 'PASDETPS' au debu du calcul * 'SI' ('EXISTE' RV 'DTIMP') ; ALPDT = RV . 'DTIMP' ; 'SINON' ; ALPDT = 'MINIMUM' LREEDT ; 'FINSI' ; ALPDT1 = (RV . 'TFINAL') '-' (KTPS . 'TPS') ; 'SI' ((ALPDT 'EGA' ALPDT1 (ALPDT '*' 1.0D-6)) 'OU' (ALPDT > ALPDT1)) ; KTPS . 'TPSP' = (RV . 'TFINAL') ; ALPDT = ALPDT1 ; LOGQIE = VRAI ; * Dans ce cas, il faut metre a jour la matrice à inverser * car elle peut etre tres differnte pas rapport à celle * calculé precedentment LOGMAOLD = FAUX ; 'SINO' ; KTPS . 'TPSP' = (KTPS . 'TPS') '+' ALPDT ; 'FINSI' ; * 'SI' LOGEXP ; ************************************************************************* ********** Cas explicite pour les variables conservatives *************** ************************************************************************* 'SI' ((NITEX 'EGA' 1) 'OU' ((RV . 'ORDTPS') 'EGA' 1)) ; RV . 'UN' = (RV . 'UNM') '+' (ALPDT '*' RESEXP) ; 'SINON' ; RV . 'UN' = ((4. '/' 3.) '*' (RV . 'UNM')) '+' ((-1. '/' 3.) '*' (RV . 'UNM2')) '+' (((2. * ALPDT) '/' 3.) '*' RESEXP) ; 'FINSI' ; 'SI' (RV . 'CLIM') ; * 'MESSAGE' 'PROLIM after an explicit iteration' ; RV . 'UN' = PROLIM RV ; 'FINSI' ; 'SINON' ; ************************************************************************* ********** Cas implicite pour les variables conservatives *************** ************************************************************************* * ****** Initialisation of RV . 'RESULTATS' . 'ERRONLIN' ; * RV . 'RESULTATS' . 'ERRONLIN' = 'PROG' ; * ************************************************************************* ************** Les iterations internes ********************************** ************************************************************************* * 'REPETER' BLINT (RV . 'NITMAIN') ; * RV . 'RESULTATS' . 'ITERIN' = &BLINT ; * ****** Implicite * 'SI' ((NITEX 'EGA' 1) 'OU' ((RV . 'ORDTPS') 'EGA' 1)) ; * *************** Implicit Euler * MAT1 = 'KOPS' 'MULT' (1. '/' ALPDT) (RV . 'MATIDE') ; SMB = ((RV . 'UNM') '-' (RV . 'UN')) '/' ALPDT ; 'SINON' ; * **************** BDF2 * MAT1 = 'KOPS' 'MULT' (1.5 '/' ALPDT) (RV . 'MATIDE') ; SMB = ((4. '*' (RV . 'UNM')) '+' (-3. '*' (RV . 'UN')) '+' (-1. '*' (RV . 'UNM2'))) '/' (2. '*' ALPDT) ; 'FINSI' ; * ********* Bas Mach * 'SI' (RV . 'BASMACH') ; MAT1 = MAT1 'ET' (RV . 'GAMSDTAU') ; 'FINSI' ; * MATTOT = MAT1 'ET' MATASS ; RESTOT = SMB 'ET' RESIMP 'ET' RESEXP ; * ********* LOGOLD = VRAI * Si methode directe, on utilise la meme parametrisation LU * Si methode iterative, on utilise le meme preconditionneur * LOGOLD = VRAI ; 'SI' LOGMAOLD ; 'SI' (LOGUPEX 'ET' (&BLINT 'EGA' 1)) ; * LOGUPEX = VRAI -> At this external iteration we can * update the matrix to inverse 'OUBLIER' MATOLD ; 'MENAGE' ; MATOLD = MATTOT ; LOGOLD = FAUX ; 'MESSAGE' ; 'MESSAGE' 'We update the matrix' ; 'SINON' ; * We don't care about LOGUPEX if (&BLINT > 1) * In this case we update the matrix each MUPINT-th * iteration 'SI' (((&BLINT '/' (RV . 'MUPINT')) '*' (RV . 'MUPINT')) 'EGA' &BLINT) ; 'OUBLIER' MATOLD ; 'MENAGE' ; MATOLD = MATTOT ; LOGOLD = FAUX ; 'MESSAGE' ; 'MESSAGE' 'We update the matrix' ; 'FINSI' ; 'FINSI' ; 'SINON' ; * * We always update the matrix to inverse * 'OUBLIER' MATOLD ; 'MENAGE' ; LOGOLD = FAUX ; MATOLD = MATTOT ; 'FINSI' ; * 'SI' ('EXISTE' (RV . 'MATINV') 'CONVINV') ; 'SI' LOGOLD ; NLIT = ('DIME' (RV . 'MATINV' . 'CONVINV')) ; 'SI' (NLIT > RV . 'MUPLIN') ; 'MESSAGE' ; 'MESSAGE' 'We update the matrix' ; 'OUBLIER' MATOLD ; 'MENAGE' ; MATOLD = MATTOT ; 'FINSI' ; 'FINSI' ; 'FINSI' ; RV . 'MATINV' . 'MATASS' = MATOLD ; RV . 'MATINV' . 'MAPREC' = MATOLD ; * 'SI' LOGOLD ; * * LOGOLD vrai si on utilise la meme * parametrisation LU ou le meme preconditionneur * pour calculer la solution du systeme lineaire * 'SI' ((RV . 'MATINV' . 'TYPINV') 'EGA' 1) ; * Meme parametrisation LU que MATOLD DELTAU = 'KRES' MATOLD 'TYPI' (RV . 'MATINV') 'CLIM' DUCLIM 'SMBR' RESTOT 'IMPR' 0 ; 'SINON' ; * Meme preconditionneur que MATOLD DELTAU = 'KRES' MATTOT 'TYPI' (RV . 'MATINV') 'CLIM' DUCLIM 'SMBR' RESTOT 'IMPR' 0 ; 'FINSI' ; 'SINON' ; DELTAU = 'KRES' MATOLD 'TYPI' (RV . 'MATINV') 'CLIM' DUCLIM 'SMBR' RESTOT 'IMPR' 0 ; 'FINSI' ; RV . 'UN' = (RV . 'UN') '+' DELTAU ; 'SI' (RV . 'CLIM') ; UNCELL = RV . 'UN' ; * 'PROLIM after the matrix inversion' ; RV . 'UN' = PROLIM RV ; * We redefine DELTAU to compute the error DELTAU = ((RV . 'UN') '-' UNCELL) '+' DELTAU ; 'FINSI' ; * ********* Boucle sur les operateurs implicites pour calculer RESIMP * et MATASS de RV . 'UN' RESIMP = 'COPIER' CHPVID ; MATASS = 'KOPS' 'MULT' 0.0 MATVID ; 'REPETER' BLOP ('DIME' (RV . 'LISTOPER')) ; NOMPER = 'EXTRAIRE' &BLOP (RV . 'LISTOPER') ; NOTABLE = 'MOT' ('TEXTE' ('CHAINE' &BLOP NOMPER) ) ; 'SI' (RV . NOTABLE . 'IMPL') ; JACO RESIDU ALPHADT = ('TEXTE' NOMPER) (RV . NOTABLE) ; 'SI' (RV . NOTABLE . 'ANOM') ; RV . 'ANOM' = VRAI ; 'QUITTER' BLINT ; 'SINON' ; MATASS = MATASS 'ET' ('KOPS' 'MULT' -1.0D0 JACO ) ; RESIMP = RESIMP 'ET' RESIDU ; 'FINSI' ; 'FINSI' ; 'FIN' BLOP ; * ********* Test de convergence * ERRO = 'MAXIMUM' DELTAU 'ABS' 'AVEC' (RV . 'ERROMOTS') ; 'SI' (&blint 'EGA' 1) ; ERRO0 = ERRO ; 'FINSI' ; RV . 'RESULTATS' . 'ERRONLIN' = (RV . 'RESULTATS' . 'ERRONLIN') 'ET' ('PROG' ERRO) ; * ********* Impression * 'MESSAGE' ; 'MESSAGE' ('CHAINE' 'PASDETPS = ' (KTPS . 'NUPASDT') ' TPS = ' (KTPS . 'TPS') ' DT = ' ALPDT) ; 'MESSAGE' ('CHAINE' 'ITERIN = ' &BLINT) ; 'REPETER' BLERRO ('DIME' (RV . 'LISTCONS')) ; MOTERR = 'EXTRAIRE' (RV . 'LISTCONS') &BLERRO ; ERRBLE = 'MAXIMUM' DELTAU 'ABS' 'AVEC' ('MOTS' MOTERR) ; 'MESSAGE' (CHAIN 'ERREUR on ' MOTERR ' = ' ERRBLE) ; 'FIN' BLERRO ; 'SI' (('EXISTE' (RV . 'MATINV') 'CONVINV')) ; NLIT = ('DIME' (RV . 'MATINV' . 'CONVINV')) ; 'MESSAGE' ('CHAINE' 'Linear iterations = ' NLIT ); RV . 'RESULTATS' . 'NITERLIN' = NLIT ; REIT = 'EXTR' (RV . 'MATINV' . 'CONVINV') NLIT ; 'MESSAGE' ('CHAINE' 'Linear residuum = ' REIT ); RV . 'RESULTATS' . 'ERROLIN' = REIT ; 'FINSI' ; * ************ Criteres de sortie * 'SI' (&BLINT > (RV . 'NITMIIN')) ; 'SI' ('EXISTE' RV 'EPSREL') ; 'SI' (ERRO < ((RV . 'EPSREL') '*' ERRO0)) ; 'QUITTER' BLINT ; 'FINSI' ; 'FINSI' ; 'SI' (ERRO '<' (RV . 'EPSINT')) ; 'QUITTER' BLINT ; 'FINSI' ; 'FINSI' ; * 'SI' (ERRO > (1.0D3 '*' ERRO0)) ; * RV . 'ANOM' = VRAI ; * 'QUITTER' BLINT ; * 'FINSI' ; * 'SI' (&BLINT '>EG' (RV . 'NITMAIN')) ; * RV . 'ANOM' = VRAI ; * 'QUITTER' BLINT ; * 'FINSI' ; 'FIN' BLINT ; ************************************************************************* **************** Fin boucle iterations internes ************************* ************************************************************************* *********************************************** ************* Fin cas implicite *************** *********************************************** 'FINSI' ; * ('SI' (RV . 'ANOM') ;) 'FINSI' ; * 'SI' (RV . 'ANOM') ; LOGQIE = FAUX ; NITEX = NITEX '-' 1 ; KTPS . 'NUPASDT' = KTPS . 'NUPASDT' '-' 1 ; KTPS . 'TPSP' = KTPS . 'TPS' ; RV . 'UN' = 'COPIER' RV . 'UNM' ; 'SI' ('EXISTE' RV 'UNM2') ; RV . 'UNM' = 'COPIER' (RV . 'UNM2') ; 'FINSI' ; 'SINON' ; * *** Mise a jour de la table RV . 'PASDETPS' à la fin du calcul * KTPS . 'TPSM' = KTPS . 'TPS' ; KTPS . 'TPS' = KTPS . 'TPSP' ; 'FINSI' ; * **** Dernier commande de BLEX * 'SI' LOGQIE ; 'QUITTER' BLEX ; 'FINSI' ; 'FIN' BLEX ; 'FINPROC' ; ***************************************************** ***************************************************** ** FIN PROCEDURE EXEX ** ***************************************************** ***************************************************** *****$$$$ PKON ***************************************************** ***************************************************** ** PROCÉDURE PKON ** ***************************************************** ***************************************************** * Il faut définir: * * *KKONV . 'EQEX' = table générale, qui contient toutes les * informations sur le calcul qu'on va faire. * Dans ce table, on ne prend que: * - la table (KKONV . 'EQEX' . 'PASDETPS') * - les inconnues du problème * (KKONV . 'EQEX' . 'UN') * - le maillage fantome * (KKONV . 'EQEX' . 'MAIFAN') * - les vitesses de cut-off pour le bas Mach * (KKONV . 'EQEX' . 'CO1') * (KKONV . 'EQEX' . 'CO2') * * *KKONV . 'GAZ' = le modelé de gaz qu'on considère * - si KKONV . 'GAZ' = 'PERFMONO', alors * on considère un gaz parfait mono-espèce * polytropique * * *KKONV . 'MODELE' = objet modele * * *KKONV . 'LISTCONS' = les noms des variables conservatives, i.e. * densité, q.d.m., energie totale par unité de * volume * * *KKONV . 'IMPL' = calcul implicite ou non ? * * *KKONV . 'METHODE' = méthode pour le calcul du flux convectif * * *KKONV . 'TYPEJACO' = 'VLH' * 'AUSMPLUS' * 'AUSMPLM' * (à donner si (KKONV . 'IMPL')) * * *KKONV . 'ORDREESP' = ordre en espace (1 ou 2) ; * * *KKONV . 'GRADRN' coeff. pour calculer le gradient de la densité * *KKONV . 'GRADVN' coeff. pour calculer le gradient de la vitesse * *KKONV . 'GRADPN' coeff. pour calculer le gradient de la pression * *KKONV . 'VLIM' vitesses au bord imposées * (à donner si (KKONV . 'ORDREESP') = 2) * * *KKONV . 'LIMITEUR' = limiteur utilisé * (à donner si (KKONV . 'ORDREESP') = 2) * * *KKONV . 'NFROZLIM' = entier * On gele les limiteurs quand * RV . 'PASDETPS' . 'NUPASDT' = KKONV . 'NFROZLIM' * On le met dedans * *KKONV . 'FROZALR' * *KKONV . 'FROZALP' * *KKONV . 'FROZALV' * * *KKONV . 'ALPHA' = coefficient de securité pour le quel on * multiplie le pas de tps determiné par un * condition de type CFL. * Il peut etre: * - un flottant * - une mot que veut 'INF' * * *KKON . 'FACELIM' = maillage de centres de face ou on ne calcule pas * le flux convectif et le jacobien * * *KKON . 'CLIM' = LOGIQUE * si vrai on doit definir une procedure PKOCLI * * *KKON . 'DT' = output. Donne les pas de temps plus petit * (r_elem/(u+c)) * 'DEBPROC' PKON ; 'ARGUMENT' KKONV * TABLE ; * CHPVID MATVID = 'KOPS' 'MATRIK' ; * METO = KKONV . 'METHODE' ; * RV = KKONV . 'EQEX' ; * LOGIMP = KKONV . 'IMPL'; LOGEXP = 'NON' LOGIMP ; * ORDESP = KKONV . 'ORDREESP' ; * 'SI' (('NEG' (KKONV . 'ALPHA') 'INF') 'ET' ('NEG' ('TYPE' (KKONV . 'ALPHA')) FLOTTANT)) ; 'MESSAGE' 'PKON . ALPHA = ???' ; 'ERREUR' 21 ; 'FINSI' ; * UN = 'REDU' (RV . 'UN') ('DOMA' (KKONV . 'MODELE') 'CENTRE') ; * 'SI' ('NEG' (KKONV . 'GAZ') 'PERFMONO') ; * ******** EULER, monoespece, "calorically perfect" (cv = constant) * 'MESSAGE' ; 'MESSAGE' 'PKON . GAZ = ???' ; 'ERREUR' 21 ; 'FINSI' ; * **** Type de jacobien (Bas Mach ?) * 'SI' LOGIMP ; ITJACO = 0 ; 'SI' (('EGA' (KKONV . 'TYPEJACO') 'AUSMPLM') 'OU' ('EGA' (KKONV . 'TYPEJACO') 'RUSANOLM') 'OU' ('EGA' (KKONV . 'TYPEJACO') 'ROELM') 'OU' ('EGA' (KKONV . 'TYPEJACO') 'HLLCLM')) ; ITJACO = -1 ; 'FINSI' ; 'FINSI' ; * ITMETO = 0 ; 'SI' (('EGA' METO 'AUSMPLM') 'OU' ('EGA' METO 'RUSANOLM') 'OU' ('EGA' METO 'ROELM') 'OU' ('EGA' METO 'HLLCLM')) ; ITMETO = -1 ; 'FINSI' ; * ***** Les variables conservatives * MOT1 = 'EXTRAIRE' (KKONV . 'LISTCONS') 1 ; 'SI' ('EGA' ('VALE' 'DIME') 2) ; NOMMOM = 'EXTRAIRE' (KKONV . 'LISTCONS') ('LECT' 2 3 ) ; NOMVEL = 'MOTS' 'UX' 'UY' ; MOT2 = 'EXTRAIRE' (KKONV . 'LISTCONS') 4 ; NOMGP = 'MOTS' 'P1DX' 'P1DY' ; 'SINON' ; NOMMOM = 'EXTRAIRE' (KKONV . 'LISTCONS') ('LECT' 2 3 4) ; NOMVEL = 'MOTS' 'UX' 'UY' 'UZ' ; MOT2 = 'EXTRAIRE' (KKONV . 'LISTCONS') 5 ; NOMGP = 'MOTS' 'P1DX' 'P1DY' 'P1DZ' ; 'FINSI' ; RN = 'EXCO' MOT1 UN 'SCAL' ; GN = 'EXCO' NOMMOM UN NOMVEL ; RETN = 'EXCO' MOT2 UN 'SCAL' ; * ***** On calcule les variables primitive * GAMN = (RV . 'PGAZ' . 'GAMN') ; VN PN = 'PRIM' 'PERFMONO' RN GN RETN GAMN ; * VN PN = PPRIM RN GN RETN GAMN ; * 'SI' ((('MINIMUM' PN) < 0) 'OU' (('MINIMUM' RN) < 0)) ; * 'MESSAGE' 'Negative density or pressure' ; * KKONV . 'ANOM' = VRAI ; * RN = ('ABS' RN) '+' (('MAXIMUM' RN 'ABS') '*' 1.001) ; * PN = ('ABS' PN) '+' (('MAXIMUM' PN 'ABS') '*' 1.001) ; * 'FINSI' ; * **** Conditions aux limits ??? * 'SI' (KKONV . 'CLIM') ; * TABLIM = 'TABLE' definie dans le programme principal ; KKONV . 'TABLIM' . 'RN' = RN ; KKONV . 'TABLIM' . 'VN' = VN ; KKONV . 'TABLIM' . 'PN' = PN ; KKONV . 'TABLIM' . 'GAMN' = GAMN ; CHPLIM RESLIM JACLIM = PKOLIM (KKONV . 'TABLIM') ; MAILLIM = ('EXTRAIRE' CHPLIM 'MAILLAGE') 'ET' (KKONV . 'FACELIM') ; 'SINON' ; MAILLIM = KKONV . 'FACELIM' ; RESLIM JACLIM = 'KOPS' 'MATRIK' ; CHPLIM = 'COPIER' RESLIM ; 'FINSI' ; MAILLIM = 'CHANGER' 'POI1' MAILLIM ; * **** La gravite * * RESGRA = 'FIMP' 'VF' 'GRAVMONO' 'RESI' (KKONV . 'LISTCONS') RN GN (KKONV . 'GRAVITE') ; * JACGRA = 'FIMP' 'VF' 'GRAVMONO' 'JACOCONS' (KKONV . 'LISTCONS') RN GN (KKONV . 'GRAVITE') ; * * Fin contribution gravité * * ***** On calcule les variables aux faces * ORDTPS = 1 ; * 'SI' (ORDESP 'EGA' 1) ; * ********* Ordre 1 en espace * ROF VITF PF GAMF = 'PRET' 'PERFMONO' 1 ORDTPS (KKONV . 'MODELE') RN VN PN GAMN ; 'SINON' ; * ********* Ordre 2 en espace * 'SI' (KKONV . 'CLIM') ; GRADV ALV = 'PENT' (KKONV . 'MODELE') 'CENTRE' 'EULEVECT' (KKONV . 'LIMITEUR') NOMVEL VN 'CLIM' ((KKONV . 'VLIM') '+' ('EXCO' NOMVEL CHPLIM)) 'GRADGEO' (KKONV . 'GRADVN') ; GRADR ALR = 'PENT' (KKONV . 'MODELE') 'CENTRE' 'EULESCAL' (KKONV . 'LIMITEUR') ('MOTS' 'SCAL') RN 'CLIM' ('EXCO' 'RN' CHPLIM) 'GRADGEO' (KKONV . 'GRADRN') ; * GRADP ALP = 'PENT' (KKONV . 'MODELE') 'CENTRE' 'EULESCAL' (KKONV . 'LIMITEUR') ('MOTS' 'SCAL') PN 'CLIM' ('EXCO' 'PN' CHPLIM) 'GRADGEO' (KKONV . 'GRADPN') ; * 'SINON' ; GRADV ALV = 'PENT' (KKONV . 'MODELE') 'CENTRE' 'EULEVECT' (KKONV . 'LIMITEUR') NOMVEL VN 'CLIM' (KKONV . 'VLIM') 'GRADGEO' (KKONV . 'GRADVN') ; GRADR ALR = 'PENT' (KKONV . 'MODELE') 'CENTRE' 'EULESCAL' (KKONV . 'LIMITEUR') ('MOTS' 'SCAL') RN 'GRADGEO' (KKONV . 'GRADRN') ; * GRADP ALP = 'PENT' (KKONV . 'MODELE') 'CENTRE' 'EULESCAL' (KKONV . 'LIMITEUR') ('MOTS' 'SCAL') PN 'GRADGEO' (KKONV . 'GRADPN') ; 'FINSI' ; * ****** Limiters = 0 on ghostcells * 'SI' (('EXISTE' RV 'MAIFAN')) ; ALRLEV = ('REDU' ALR (RV . 'MAIFAN')) ; ALPLEV = ('REDU' ALP (RV . 'MAIFAN')) ; ALVLEV = ('REDU' ALV (RV . 'MAIFAN')) ; CELL = ALR ; ALR = ALR '-' ALRLEV ; 'DETR' CELL ; CELL = ALP ; ALP = ALP '-' ALPLEV ; 'DETR' CELL ; CELL = ALV ; ALV = ALV '-' ALVLEV ; 'DETR' CELL ; 'FINSI' ; * ****** Frozen Limiters * * External iterations * 'SI' ('EXISTE' KKONV 'NFROZLIM') ; ITEREX = RV . 'PASDETPS' . 'NUPASDT' ; 'SI' (ITEREX 'EGA' (KKONV . 'NFROZLIM')) ; 'MESSAGE' ; 'MESSAGE' 'External iterations' ; 'MESSAGE' 'We froze the limiters' ; KKONV . 'FROZALR' = 'COPIER' ALR ; KKONV . 'FROZALP' = 'COPIER' ALP ; KKONV . 'FROZALV' = 'COPIER' ALV ; 'SINON' ; 'SI' (ITEREX > (KKONV . 'NFROZLIM')) ; * * Min * DAR = (ALR '-' (KKONV . 'FROZALR')) 'ABS' ; DAP = (ALP '-' (KKONV . 'FROZALP')) 'ABS' ; DAV = (ALV '-' (KKONV . 'FROZALV')) 'ABS' ; ALR = (ALR '+' (KKONV . 'FROZALR')) '-' DAR ; ALP = (ALP '+' (KKONV . 'FROZALP')) '-' DAP ; ALV = (ALV '+' (KKONV . 'FROZALV')) '-' DAV ; ALR = ALR '/' 2.0 ; ALV = ALV '/' 2.0 ; ALP = ALP '/' 2.0 ; KKONV . 'FROZALR' = 'COPIER' ALR ; KKONV . 'FROZALP' = 'COPIER' ALP ; KKONV . 'FROZALV' = 'COPIER' ALV ; 'FINSI' ; 'FINSI' ; 'FINSI' ; * * Internal iterations * 'SI' ('EXISTE' KKONV 'NFROZLII') ; ITERIN = RV . 'RESULTATS' . 'ITERIN' ; 'SI' (ITERIN 'EGA' (KKONV . 'NFROZLII')) ; 'MESSAGE' ; 'MESSAGE' 'Internal iterations' ; 'MESSAGE' 'We froze the limiters' ; KKONV . 'FROZALR' = 'COPIER' ALR ; KKONV . 'FROZALP' = 'COPIER' ALP ; KKONV . 'FROZALV' = 'COPIER' ALV ; 'SINON' ; 'SI' (ITERIN > (KKONV . 'NFROZLII')) ; DAR = (ALR '-' (KKONV . 'FROZALR')) 'ABS' ; DAP = (ALP '-' (KKONV . 'FROZALP')) 'ABS' ; DAV = (ALV '-' (KKONV . 'FROZALV')) 'ABS' ; ALR = (ALR '+' (KKONV . 'FROZALR')) '-' DAR ; ALP = (ALP '+' (KKONV . 'FROZALP')) '-' DAP ; ALV = (ALV '+' (KKONV . 'FROZALV')) '-' DAV ; ALR = ALR '/' 2.0 ; ALV = ALV '/' 2.0 ; ALP = ALP '/' 2.0 ; KKONV . 'FROZALR' = 'COPIER' ALR ; KKONV . 'FROZALP' = 'COPIER' ALP ; KKONV . 'FROZALV' = 'COPIER' ALV ; 'FINSI' ; 'FINSI' ; 'FINSI' ; * * ********* Ordre 1 en temps * ROF VITF PF GAMF = 'PRET' 'PERFMONO' ORDESP ORDTPS (KKONV . 'MODELE') RN GRADR ALR VN GRADV ALV PN GRADP ALP GAMN ; 'FINSI' ; * 'SI' (ITMETO 'EGA' (-1)) ; RESIDU DELTAT = 'KONV' 'VF' 'PERFMONO' 'RESI' METO (KKONV . 'LISTCONS') (KKONV . 'MODELE') ROF VITF PF GAMF MAILLIM (RV . 'CO1') (RV . 'CO2') ; 'SINON' ; RESIDU DELTAT = 'KONV' 'VF' 'PERFMONO' 'RESI' METO (KKONV . 'LISTCONS') (KKONV . 'MODELE') ROF VITF PF GAMF MAILLIM ; 'FINSI' ; * 'SI' ('EGA' ('TYPE' (KKONV . 'ALPHA')) FLOTTANT) ; ALPDT = 'PROG' ((KKONV . 'ALPHA') '*' DELTAT) ; 'SINON' ; ALPDT = 'PROG' ; 'FINSI' ; KKONV . 'DT' = DELTAT ; * **** Low Mach * 'SI' (RV . 'BASMACH') ; MATBM = 'KONV' 'VF' 'PERFMONO' 'GAMMCONS' (KKONV . 'LISTCONS') ('DOMA' (KKONV . 'MODELE') 'CENTRE') ('DOMA' (KKONV . 'MODELE') 'DIAMIN') RN VN PN GAMN (RV . 'CO1') (RV . 'CO2') ; RV . 'GAMSDTAU' = 'KOPS' 'MULT' (1. '/' (KKONV . 'CFLDTAU')) MATBM ; 'FINSI' ; * 'SI' LOGIMP ; 'SI' (ITJACO 'EGA' 0) ; JACO = 'KONV' 'VF' 'PERFMONO' 'JACOCONS' (KKONV . 'MODELE') (KKONV . 'LISTCONS') MAILLIM (KKONV . 'TYPEJACO') RN VN PN GAMN ; 'SINON' ; JACO = 'KONV' 'VF' 'PERFMONO' 'JACOCONS' (KKONV . 'MODELE') (KKONV . 'LISTCONS') MAILLIM (KKONV . 'TYPEJACO') RN VN PN GAMN (RV . 'CO1') (RV . 'CO2') ; 'FINSI' ; 'RESPRO' (JACO 'ET' JACLIM 'ET' JACGRA) (RESIDU '+' RESLIM '+' RESGRA) ALPDT ; 'SINON' ; 'RESPRO' (RESIDU '+' RESLIM '+' RESGRA) ALPDT ; 'FINSI' ; * **** On detrui les choses qui ne servent plus * 'DETR' UN ; 'OUBL' UN ; 'DETR' RN ; 'DETR' GN ; 'DETR' RETN ; 'DETR' VN ; 'DETR' PN ; 'OUBL' RN ; 'OUBL' GN ; 'OUBL' RETN ; 'OUBL' VN ; 'OUBL' PN ; * **** Les MCHAML faces * 'DETR' ROF ; 'DETR' VITF ; 'DETR' PF ; 'DETR' GAMF ; 'OUBL' ROF ; 'OUBL' VITF ; 'OUBL' PF ; 'OUBL' GAMF ; * **** Les pentes et les limiteurs * 'SI' (ORDESP 'EGA' 2); * 'DETR' GRADR ; 'DETR' GRADP ; 'DETR' GRADV ; 'DETR' ALR ; 'DETR' ALP ; 'DETR' ALV; * 'OUBL' GRADR ; 'OUBL' GRADP ; 'OUBL' GRADV ; 'OUBL' ALR ; 'OUBL' ALP ; 'OUBL' ALV; * 'FINSI' ; 'FINPROC' ; ***************************************************** ***************************************************** ** FIN PROCEDURE PKON ** ***************************************************** ***************************************************** ********************************************************************* **** Procedure PKOLIM *********************************************** ********************************************************************* 'DEBPROC' PKOLIM ; 'ARGUMENT' TABLIM*'TABLE' ; LISTP = 'MOTS' 'RN' 'UX' 'UY' 'PN' ; RCHLIM1 RCHRES1 = 'KONV' 'VF' 'PERFMONO' 'CLIM' 'RESI' MDOMINT MLIGG LISTCONS LISTP (TABLIM . 'RN') (TABLIM . 'VN') (TABLIM . 'PN') (TABLIM . 'GAMN') (TABLIM . 'CHINRI') 'INRI' ; RJACO1 = 'KONV' 'VF' 'PERFMONO' 'CLIM' 'JACOCONS' MDOMINT MLIGG LISTCONS LISTP (TABLIM . 'RN') (TABLIM . 'VN') (TABLIM . 'PN') (TABLIM . 'GAMN') (TABLIM . 'CHINRI') 'INRI' ; RCHLIM2 RCHRES2 = 'KONV' 'VF' 'PERFMONO' 'CLIM' 'RESI' MDOMINT MLIGD LISTCONS LISTP (TABLIM . 'RN') (TABLIM . 'VN') (TABLIM . 'PN') (TABLIM . 'GAMN') (TABLIM . 'CHOUTR') 'OUTRI' ; RJACO2 = 'KONV' 'VF' 'PERFMONO' 'CLIM' 'JACOCONS' MDOMINT MLIGD LISTCONS LISTP (TABLIM . 'RN') (TABLIM . 'VN') (TABLIM . 'PN') (TABLIM . 'GAMN') (TABLIM . 'CHOUTR') 'OUTRI' ; CHPLIM = RCHLIM1 '+' RCHLIM2 ; RESLIM = RCHRES1 '+' RCHRES2 ; JACLIM = RJACO1 'ET' RJACO2 ; 'RESPRO' CHPLIM RESLIM JACLIM ; 'FINPROC' ; ********************************************************************* **** Procedure CALC ************************************************* ********************************************************************* * * Personal procedure * 'DEBP' CALC ; 'ARGU' RVX*'TABLE' ; * RV = RVX . 'EQEX' ; UN = RV . 'UN' ; MCALC = RVX . 'MODELE' ; * * During the external iteration (except the first; in which * UNM = UN of the previous step ) * 'SI' ('NEG' 0 (RV . 'RESULTATS' . 'ITERIN')) ; RN = 'REDU' ('DOMA' MCALC 'CENTRE') ('EXCO' UN NOMDEN ('MOTS' 'SCAL')) ; GN = 'REDU' ('DOMA' MCALC 'CENTRE') ('EXCO' NOMMOM UN NOMVEL) ; RETN = 'REDU' ('DOMA' MCALC 'CENTRE') ('EXCO' NOMRET UN ('MOTS' 'SCAL')) ; VN PN = 'PRIM' 'PERFMONO' RN GN RETN (RV . 'PGAZ' . 'GAMN') ; * 'SI' ('EGA' (RVX . 'COMPT') 0) ; RVX . 'IT' = 'PROG' ; RVX . 'ER' = 'PROG' ; 'SINON' ; CHPVOL = 'DOMA' MCALC 'VOLUME' ; VOLTOT = 'XTY' CHPVOL ('MANUEL' CHPO ('DOMA' MCALC 'CENTRE') 1 'SCAL' 1.) ('MOTS' 'SCAL') ('MOTS' 'SCAL'); ERRO = 'ABS' (PN '-' (RVX . 'PN0')) ; ERRO = 'XTY' ERRO CHPVOL ('MOTS' 'SCAL') ('MOTS' 'SCAL') ; ERRO = ERRO '/' VOLTOT ; RVX . 'IT' = (RVX . 'IT') 'ET' ('PROG' (RVX . 'COMPT')) ; RVX . 'ER' = (RVX . 'ER') 'ET' ('PROG' ERRO) ; 'FINSI' ; RVX . 'PN0' = PN ; RVX . 'COMPT' = (RVX . 'COMPT') '+' 1 ; 'FINSI' ; * * We change the cut-off during the internal iterations * 'SI' FAUX ; 'SI' ((RV . 'RESULTATS' . 'ITERIN') > 0) ; CUMAX = RV . 'CO1MAX' ; CUMIN = RV . 'CO1MIN' ; DCU = (CUMIN '-' CUMAX) '/' 5 ; CU = CUMAX '+' (DCU * (RV . 'RESULTATS' . 'ITERIN')) ; * Take the max between CU and CUMIN A = CUMIN '+' CU ; B = 'ABS' (CUMIN '-' CU) ; CU = 0.5 '*' (A '+' B) ; * RV . 'CO1' = CU ; RV . 'CO2' = CU ; 'FINSI' ; 'FINSI' ; * IALPDT = 'PROG' ; IRESU IJACO ='KOPS' 'MATRIK' ; * 'MENAGE' ; * 'SI' (RVX . 'IMPL') ; 'RESPRO' IJACO IRESU IALPDT ; 'SINON' ; 'RESPRO' IRESU IALPDT ; 'FINSI' ; * 'FINP' ; * ************************************************************************ ************************************************************************ ***************** FIN PARTIE PROCEDURES ******************************** ************************************************************************ ************************************************************************ ************************************************************************ ************************************************************************ ************************************************************************ **************************** MESH ************************************** ************************************************************************ ************************************************************************ NY = 5 '*' RAF ; NX1 = 4 '*' RAF ; NX2 = 2 '*' NX1 ; NX3 = NX1 ; NX = (NX1 '+' NX2 '+' NX3) ; DX = (4.0 '/' NX) ; A0 = -2.0 0.0 ; A1 = -1.0 0.0 ; A2 = 1.0 0.0 ; A3 = 2.0 0.0 ; A4 = 2.0 1.0 ; A5 = -2.0 1.0 ; * **** LIGB * LIGB1 = A0 'DROIT' NX1 A1 ; * LIGB2 (On utilise un propriete de 'ET' ; si 'ET' change ?) xcel = ('COORDONNEE' 1 A1) '+' DX ; ycel = 0.1 '*' ( 1.0 '+' ('COS' (180 '*' xcel))); ACEL = xcel ycel ; LIGB2 = A1 'DROIT' 1 ACEL ; 'REPETER' BL1 (NX2 '-' 2) ; ACEL0 = ACEL ; xcel = xcel '+' DX ; ycel = 0.1 '*' ( 1.0 '+' ('COS' (180 '*' xcel))); ACEL = xcel ycel ; LIGB2 = LIGB2 'ET' (ACEL0 'DROIT' 1 ACEL) ; 'FIN' BL1; LIGB2 = LIGB2 'ET' (ACEL 'DROIT' 1 A2) ; LIGB3 = A2 'DROIT' NX3 A3 ; LIGB = LIGB1 'ET' LIGB2 'ET' LIGB3 ; * **** LIGH * LIGH = A4 'DROIT' NX A5 ; * **** DOMINT * DOMINT = LIGH 'REGLER' NY ('INVERSE' LIGB) ; LIGCON = 'CONTOUR' DOMINT ; * *** LIGG * LIGG = LIGCON 'ELEM' 'COMP' A5 A0 ; * **** LIGD * LIGD = LIGCON 'ELEM' 'COMP' A3 A4 ; * **** MODEL OBJECTS * MDOMINT = 'MODELISER' DOMINT 'EULER' ; MLIGG = 'MODELISER' LIGG 'EULER' ; MLIGD = 'MODELISER' LIGD 'EULER' ; * **** Creation of DOMAINE tables via the MODEL object * TDOMINT = 'DOMA' MDOMINT 'VF' ; TLIGG = 'DOMA' MLIGG 'VF' ; TLIGD = 'DOMA' MLIGD 'VF' ; QDOMINT = TDOMINT . 'QUAF' ; QLIGD = TLIGD . 'QUAF' ; QLIGG = TLIGG . 'QUAF' ; 'ELIMINATION' QDOMINT (DX '/' 100.) QLIGG ; 'ELIMINATION' QDOMINT (DX '/' 100.) QLIGD ; 'SI' GRAPH ; 'TRACER' (DOMINT 'ET' (LIGG 'COULEUR' 'ROUG') 'ET' (LIGD 'COULEUR' 'BLEU')) 'TITRE' 'Domaine total' ; 'FINSI' ; ******************************************* ****** LIGNE de Post-traitement *********** ******************************************* LIGB = LIGCON 'ELEM' 'COMPRIS' A0 A3 ; POIN0 = LIGB 'POIN' 1; X1 Y1 = 'COORDONNEE' POIN0 ; 'REPETER' BLLIM (('NBNO' LIGB) '-' 1) ; X0 = X1 ; Y0 = Y1 ; POIN0 = LIGB 'POIN' (&BLLIM '+' 1) ; X1 Y1 = 'COORDONNEE' POIN0 ; XFAC = (X0 '+' X1) '/' 2 ; YFAC = (Y0 '+' Y1) '/' 2 ; PFAC = (TDOMINT . 'FACE') 'POIN' 'PROC' (XFAC YFAC); GEOFAC1 = (TDOMINT . 'FACEL') 'ELEM' 'APPUYE' 'LARGEMENT' PFAC ; * **** Tranformation en POI1 * GEO1POI1 = 'CHANGER' 'POI1' GEOFAC1 ; PCEL11 = GEO1POI1 'POIN' 1 ; PCEL12 = GEO1POI1 'POIN' 2 ; * **** Il faur verifier que PFAC = PCEL12 = PCEL22 * ('NBEL' GEO1POI1) = ('NBEL' GEO2POI1) = 2 * 'SI' ( ('NBEL' GEO1POI1) 'NEG' 2); 'MESSAGE' ; 'MESSAGE' 'Probleme dans la creation de la ligne pour le post.'; 'MESSAGE' ; 'ERREUR' 21 ; 'FINSI' ; 'SI' ( PCEL12 'NEG' PFAC); 'MESSAGE' ; 'MESSAGE' 'Probleme dans la creation de la ligne pour le post.'; 'MESSAGE' ; 'ERREUR' 21 ; 'FINSI' ; * *** Creation d'un maillage SEG2 * 'SI' (&BLLIM 'EGA' 1); PCEN0 = PCEL11 ; 'SINON' ; 'SI' (&BLLIM 'EGA' 2); LIGPOST = 'MANUEL' 'SEG2' PCEN0 PCEL11 ; PCEN0 = PCEL11 ; 'SINO' ; LIGPOST = LIGPOST 'ET' ('MANUEL' 'SEG2' PCEN0 PCEL11) ; PCEN0 = PCEL11 ; 'FINSI' ; 'FINSI' ; 'FIN' BLLIM ; 'SI' GRAPH ; 'TRACER' (DOMINT 'ET' (LIGPOST 'COULEUR' 'VERT')) 'TITRE' 'LIGPOST' ; 'FINSI' ; *************************************************************** *************************************************************** *************************************************************** ************** Initial conditions ***************************** *************************************************************** *************************************************************** *************************************************************** * *** C.L. et initiales * C_INF = GAMSCAL * P_INF '/' RO_INF ; C_INF = C_INF '**' 0.5 ; M_INF = U_INF '/' C_INF ; * * Names of conserved variables and others * NOMDEN = 'MOTS' 'RN' ; NOMMOM = 'MOTS' 'RUX' 'RUY' ; NOMRET = 'MOTS' 'RETN' ; NOMVEL = 'MOTS' 'UX' 'UY' ; * NOMPRE = 'MOTS' 'NOMPRE' ; LISTCONS = NOMDEN 'ET' NOMMOM 'ET' NOMRET ; * GAMN = 'MANUEL' 'CHPO' (TDOMINT . 'CENTRE') 1 'SCAL' gamscal; rek_inf = 0.5D0 * u_inf * u_inf * ro_inf ; rei_inf = P_INF '/' (gamscal '-' 1.0D0) ; RN0 = 'MANUEL' 'CHPO' (TDOMINT . 'CENTRE') 1 'SCAL' ro_inf ; GN0 = 'MANUEL' 'CHPO' (TDOMINT . 'CENTRE') 2 'UX' (ro_inf * u_inf) 'UY' 0.0 ; RETN0 = 'MANUEL' 'CHPO' (TDOMINT . 'CENTRE') 1 'SCAL' (rei_inf '+' rek_inf) ; UNCONS = ('NOMC' 'RN' RN0 'NATU' 'DISCRET') 'ET' ('NOMC' ('MOTS' 'UX' 'UY') GN0 ('MOTS' 'RUX' 'RUY') 'NATU' 'DISCRET') 'ET' ('NOMC' 'RETN' RETN0 'NATU' 'DISCRET') ; VN0 PN0 = 'PRIM' 'PERFMONO' RN0 GN0 RETN0 GAMN ; ERRO = ('MAXIMUM' (PN0 '-' p_inf) 'ABS') / p_inf ; 'SI' (ERRO > 1.0D-16) ; 'ERREUR' 21 ; 'FINSI' ; aa = 'EXCO' 'UX' VN0 ; ERRO = ('MAXIMUM' (aa '-' u_inf) 'ABS') / u_inf ; 'SI' (ERRO > 1.0D-16) ; 'ERREUR' 21 ; 'FINSI' ; VN2 = 'PSCAL' VN0 VN0 ('MOTS' 'UX' 'UY') ('MOTS' 'UX' 'UY') ; C2 = GAMN '*' ( PN0 '/' RN0) ; CN0 = C2 '**' 0.5 ; MACH2 = VN2 '/' C2; MACHN0 = MACH2 '**' 0.5; * **** Plot of IC * 'SI' GRAPH ; CHM_RN = 'KCHA' MDOMINT 'CHAM' RN0 ; CHM_PN = 'KCHA' MDOMINT 'CHAM' PN0 ; CHM_VN = 'KCHA' MDOMINT 'CHAM' VN0 ; CHM_MACH = 'KCHA' MDOMINT 'CHAM' MACHN0 ; 'TRACER' CHM_RN MDOMINT 'TITR' ('CHAINE' 'RN at t=' 0.0); 'TRACER' CHM_PN MDOMINT 'TITR' ('CHAINE' 'PN at t=' 0.0); 'TRACER' CHM_VN MDOMINT 'TITR' ('CHAINE' 'VN at t=' 0.0); 'TRACER' CHM_MACH MDOMINT 'TITR' ('CHAINE' 'MACH at t=' 0.0); 'FINSI' ; ***************************************************** ***************************************************** ***************************************************** **************** La table RV ********************** ***************************************************** ***************************************************** ***************************************************** * RV = 'TABLE' ; RV . 'ANOM' = FAUX ; * Constant time step CFLREF = 100.0 ; UREF = u_inf ; UREF = c_inf ; RV . 'DTIMP' = CFLREF * DX '/' UREF ; * RV . 'ORDTPS' = 1 ; * RV . 'UN' = UNCONS ; RV . 'LISTCONS' = LISTCONS ; RV . 'ERROMOTS' = NOMRET ; * RV . 'CLIM' = FAUX ; * RV . 'CONS' = FAUX ; * RV . 'TFINAL' = 1.0E6 ; RV . 'NITMAEX' = NITER ; * RV . 'NITMAIN' = 6 ; RV . 'NITMIIN' = 6 ; RV . 'EPSINT' = 1.0D-16 ; * * RV . 'FEXT' = 20 ; * RV . 'FINT' = 1 ; * RV . 'LISTOPER' = 'MOTS' 'CALC' 'PKON' ; * CALCTAB = 'TABLE' ; RV . '1CALC' = CALCTAB ; PKONTAB = 'TABLE' ; RV . '2PKON' = PKONTAB ; * RV . 'MATIDE' = 'KOPS' 'MATIDE' LISTCONS ('DOMA' MDOMINT 'CENTRE') 'MATRIK' ; * **** Parametres de la procedure PROLIM * * Gas model * Gas properties RV . 'PGAZ' = 'TABLE' ; RV . 'PGAZ' . 'R' = Rair ; RV . 'PGAZ' . 'GAMN' = gamn ; * RV . 'PGAZ' . 'MU' = 0.0 ; * RV . 'PGAZ' . 'LAMBDA' = 0.0 ; * Bas Mach (PKON) RV . 'BASMACH' = FAUX ; RV . 'CO1' = 'MANUEL' 'CHPO' ('DOMA' MDOMINT 'CENTRE') 1 'SCAL' (10 * U_INF) ; RV . 'CO1MAX' = 'MANUEL' 'CHPO' ('DOMA' MDOMINT 'CENTRE') 1 'SCAL' (10 * U_INF) ; RV . 'CO1MIN' = 'MANUEL' 'CHPO' ('DOMA' MDOMINT 'CENTRE') 1 'SCAL' (10 * U_INF) ; RV . 'CO2' = 'COPIER' ( RV . 'CO1') ; * RV . 'CO2' can be Modified into PDIF ***************** ** CALCUL ******* ***************** * Personal procedure * CALCTAB . 'ANOM' = FAUX ; CALCTAB . 'EQEX' = RV ; CALCTAB . 'MODELE' = MDOMINT ; * We call this procedure during the external iterations -> CALCTAB . 'IMPL' = VRAI ; CALCTAB . 'COMPT' = 0 ; * ***************** * PKON ********** ***************** CHPVID MATVID = 'KOPS' 'MATRIK' ; PKONTAB . 'CFLDTAU' = 10000. ; PKONTAB . 'ANOM' = FAUX ; PKONTAB . 'EQEX' = RV ; PKONTAB . 'GAZ' = 'PERFMONO' ; PKONTAB . 'MODELE' = MDOMINT ; PKONTAB . 'LISTCONS' = LISTCONS ; * PKONTAB . 'METHODE' = 'CENTERED' ; * PKONTAB . 'METHODE' = 'RUSANOLM' ; * PKONTAB . 'METHODE' = 'AUSMPLM' ; PKONTAB . 'METHODE' = 'AUSMPLUS' ; * PKONTAB . 'METHODE' = 'VLH' ; * PKONTAB . 'METHODE' = 'ROELM' ; * PKONTAB . 'METHODE' = 'ROE' ; * PKONTAB . 'METHODE' = 'HLLCLM' ; * PKONTAB . 'IMPL' = VRAI ; 'SI' (PKONTAB . 'IMPL') ; PKONTAB . 'TYPEJACO' = 'AUSMPLUS' ; * PKONTAB . 'TYPEJACO' = 'VLH' ; * PKONTAB . 'TYPEJACO' = 'RUSANOLM' ; * PKONTAB . 'TYPEJACO' = 'AUSMPLM' ; 'FINSI' ; * PKONTAB . 'ORDREESP' = 2 ; 'SI' ((PKONTAB . 'ORDREESP') 'EGA' 2) ; PKONTAB . 'VLIM' = CHPVID ; CHPENLIM = 'MANUEL' 'CHPO' (('DOMA' MLIGG 'CENTRE') 'ET' ('DOMA' MLIGD 'CENTRE')) 1 'SCAL' 0.0 ; CHPE2LIM = 'MANUEL' 'CHPO' (('DOMA' MLIGG 'CENTRE') 'ET' ('DOMA' MLIGD 'CENTRE')) 2 'UX' 1.0 'UY' 0.0 ; TOTO TITI MCHRCON = 'PENTE' MDOMINT 'CENTRE' 'EULESCAL' 'LIMITEUR' ('MOTS' 'SCAL') RN0 'CLIM' CHPENLIM ; TOTO TITI MCHVCON = 'PENTE' MDOMINT 'CENTRE' 'EULEVECT' 'LIMITEUR' ('MOTS' 'UX' 'UY') GN0 'CLIM' CHPE2LIM ; PKONTAB . 'LIMITEUR' = 'NOLIMITE' ; PKONTAB . 'GRADRN' = MCHRCON ; PKONTAB . 'GRADPN' = MCHRCON ; PKONTAB . 'GRADVN' = MCHVCON ; 'FINSI' ; PKONTAB . 'NFROZLIM' = 50 ; * Gravité * PKONTAB . 'CORRGRA' = FAUX ; PKONTAB . 'GRAVITE' = 'MANUEL' 'CHPO' ('DOMA' MDOMINT 'CENTRE') 2 'UX' 0.0 'UY' 0.0 ; * MAIVID = 'DIFF' ('DOMA' MLIGG 'CENTRE') ('DOMA' MLIGG 'CENTRE') ; PKONTAB . 'FACELIM' = MAIVID ; PKONTAB . 'ALPHA' = 'INF' ; PKONTAB . 'CLIM' = VRAI ; PKONTAB . 'TABLIM' = TABLE ; * * Normalisation du profile parabolique, pour obtenir que * le flux rentrante vaut MINJ * PKONTAB . 'TABLIM' . 'CHINRI' = 'MANUEL' 'CHPO' ('DOMA' MLIGG 'CENTRE') 4 'RN' ro_inf 'UX' U_INF 'UY' 0.0 'PN' p_inf ; PKONTAB . 'TABLIM' . 'CHOUTR' = 'MANUEL' 'CHPO' ('DOMA' MLIGD 'CENTRE') 4 'RN' ro_inf 'UX' U_INF 'UY' 0.0 'PN' p_inf ; * *********************************** * Inversion de la matrice ********* *********************************** * * RV . 'MATUPDAT' updating * RV . 'MUPEXT' external iterations updating * RV . 'MUPINT' internal iterations updating * RV . 'MUPLIN' We update the matrix if in the previous * internal iteration the number of linear * iterations to solve the system were bigger * than (RV . 'MUPLIN') RV . 'MATUPDAT' = FAUX ; 'SI' ('NON' ( RV . 'MATUPDAT')) ; RV . 'MUPEXT' = 1 ; RV . 'MUPINT' = 50 ; RV . 'MUPLIN' = 4000 ; 'FINSI' ; RV . 'MATINV' = 'TABLE' 'METHINV' ; MATTAB = RV . 'MATINV' ; * MATTAB . 'TYPINV' = 1 : methode exact * MATTAB . 'TYPINV' = 5 ; GMRES MATTAB . 'TYPINV' = 1 ; MATTAB . 'IMPINV' = 0 ; * * Matrice pour assurer que la matrice à inverser est correctement * assemblé * * MATTAB . 'MATASS' definie en EXEXIM * MATTAB . 'MAPREC' " * MATTAB . 'XINIT' " * * Methode de numerotation de DDL MATTAB . 'TYRENU' = 'RIEN' ; MATTAB . 'PCMLAG' = 'APR2' ; MATTAB . 'GMRESTRT' = 500 ; MATTAB . 'SCALING' = 1 ; * ILU 3 * ILUT (dual truncation) 5 * ILUTP 9 * Dans le cas ILUT il faut definir * ILUTLFIL * ILUTDTOL MATTAB . 'PRECOND' = 7 ; MATTAB . 'OUBMAT' = 0 ; MATTAB . 'ILUTPPIV' = 0.1 ; MATTAB . 'ILUTPMBL' = 0 ; MATTAB . 'ILUTLFIL' = 4. ; MATTAB . 'ILUTDTOL' = 0. ; MATTAB . 'NITMAX' = 3000 ; MATTAB . 'RESID' = 1.D-6 ; * **** Save results into a file * SI FAUX ; FICHER1 = ('CHAINE' FICHER 'main' RAF (PKONTAB . 'METHODE') 'OE' (PKONTAB . 'ORDREESP') 'OT' (RV . 'ORDTPS') '.sauv') ; 'OPTION' 'SAUVER' FICHER1 ; 'MESSAGE' 'We save results into' ; 'MESSAGE' FICHER ; 'MESSAGE' ; 'FINSI' ; * **** Exexcution EXEXIM * 'SI' FAUX ; LISTTPS = ('PROG' 100.0) ; 'REPETER' BL1 ('DIME' LISTTPS) ; RV . 'TFINAL' = 'EXTRAIRE' LISTTPS &BL1 ; 'TEMPS' 'ZERO' ; EXEXIM RV ; 'TEMPS' 'IMPR' ; * **** To save memory storage * 'OUBLIER' MATTAB MATASS ; 'OUBLIER' MATTAB MAPREC ; 'MENAGE' ; 'SI' FAUX ; 'SAUVER' 'LABEL' ('CHAINE' &BL1) ; 'FINSI' ; 'FIN' BL1 ; SINON ; * 'TEMPS' 'ZERO' ; EXEXIM RV ; * 'TEMPS' 'IMPR' ; 'FINSI' ; * *************************************************************** *************************************************************** *************************************************************** ************** Plot of the solutions ************************** *************************************************************** *************************************************************** *************************************************************** * * GRAPH = VRAI ; * 'SI' FAUX ; FICHER = './D_sauv/main16VLHOE2OT1.sauv' ; 'OPTION' 'RESTITUER' FICHER ; 'RESTITUER' ; 'FINSI' ; * * LOGGNP = logical. If VRAI, files for gnuplot are created * LOGGNP = VRAI ; VALECH = 1 ; 'OPTION' 'ECHO' VALECH 'TRAC' 'X' ; TPS = RV . 'PASDETPS' . 'TPS' ; METO = PKONTAB . 'METHODE' ; ORD_ESP = PKONTAB . 'ORDREESP' ; ORD_TPS = RV . 'ORDTPS' ; NTOTEL = 'NBNO' ('DOMA' MDOMINT 'CENTRE') ; * **** Les variables conservatives * GAMN = 'REDU' GAMN ('DOMA' MDOMINT 'CENTRE') ; RN = 'EXCO' NOMDEN (RV . 'UN') ('MOTS' 'SCAL') ; GN = 'EXCO' NOMMOM (RV . 'UN') ('MOTS' 'UX' 'UY') ; RETN = 'EXCO' NOMRET (RV . 'UN') ('MOTS' 'SCAL') ; * **** Les variables primitives * VN PN = 'PRIM' 'PERFMONO' RN GN RETN GAMN ; CSON2 = (GAMN '*' PN) '/' RN ; VN2 = 'PSCAL' VN VN ('MOTS' 'UX' 'UY') ('MOTS' 'UX' 'UY') ; MACHN2 = VN2 '/' CSON2 ; MACHN = MACHN2 '**' 0.5 ; TITOLO = ('CHAINE' METO ', OE = ' ORD_ESP ', OT =' ORD_TPS ', NBEL = ' NTOTEL ' ,tps =' TPS) ; EVO1 = 'EVOL' 'MANU' 'it' (CALCTAB . 'IT') 'err' ((CALCTAB . 'ER') '+' 1.0D-16) ; * *** GRAPHIQUE DES SOLUTIONS * 'SI' GRAPH ; 'TRACER' DOMINT 'TITRE' ('CHAINE' 'Maillage') ; CHM_RN = 'KCHA' MDOMINT 'CHAM' RN ; CHM_PN = 'KCHA' MDOMINT 'CHAM' PN ; CHM_VN = 'KCHA' MDOMINT 'CHAM' VN ; CHM_MACH = 'KCHA' MDOMINT 'CHAM' MACHN ; 'TRACER' CHM_RN MDOMINT ('CONTOUR' DOMINT) 'TITR' ('CHAINE' ' RN : ' TITOLO) ; RNV = 'ELNO' TDOMINT ('REDU' RN TDOMINT . 'CENTRE') ; 'OPTION' 'ISOV' 'LIGN' ; 'TRACER' DOMINT RNV ('CONTOUR' DOMINT) 'TITRE' ('CHAINE' 'RN : ' TITOLO) 15 ; 'OPTION' 'ISOV' 'SULI' ; 'TRACER' CHM_PN MDOMINT ('CONTOUR' DOMINT) 'TITR' ('CHAINE' 'PN : ' TITOLO) ; PNV = 'ELNO' TDOMINT ('REDU' PN TDOMINT . 'CENTRE') ; 'OPTION' 'ISOV' 'LIGN' ; 'TRACER' DOMINT PNV ('CONTOUR' DOMINT) 'TITRE' ('CHAINE' 'PN : ' TITOLO) 15 ; 'OPTION' 'ISOV' 'SULI' ; 'TRACER' CHM_MACH MDOMINT ('CONTOUR' DOMINT) 'TITR' ('CHAINE' 'MACHN : ' TITOLO); MACHNV = 'ELNO' TDOMINT ('REDU' MACHN TDOMINT . 'CENTRE') ; 'OPTION' 'ISOV' 'LIGN' ; 'TRACER' DOMINT MACHNV ('CONTOUR' DOMINT) 'TITRE' ('CHAINE' 'Mach : ' TITOLO) 15 ; 'OPTION' 'ISOV' 'SULI' ; 'TRACER' CHM_VN MDOMINT ('CONTOUR' DOMINT) 'TITR' ('CHAINE' 'VN : ' TITOLO) ; VECN = 'VECTEUR' VN (8. '/' RAF) 'UX' 'UY' 'JAUNE' ; 'TRACER' DOMINT VECN ('CONTOUR' DOMINT) 'TITRE' ('CHAINE' 'VN : ' TITOLO) ; * 'OPTION' 'ISOV' 'SULI' ; 'DESSIN' EVO1 LOGY 'TITRE' 'Error' ; * 'FINSI' ; * **** Test de convergence * AA = 'EXTRAIRE' (CALCTAB . 'ER') ('DIME' (CALCTAB . 'ER')) ; 'SI' (AA > 1.0D-8) ; 'MESSAGE' 'Probleme' ; 'ERREUR' 5 ; 'FINSI' ; * XLIGPOST = 'COORDONNEE' 1 LIGPOST ; EVCOOR = 'EVOL' 'CHPO' XLIGPOST 'SCAL' LIGPOST ; LISTX = 'EXTRAIRE' EVCOOR 'ORDO' ; EVMACH = 'EVOL' 'CHPO' MACHN 'SCAL' LIGPOST ; LISTMACH = 'EXTRAIRE' EVMACH 'ORDO' ; EVMACH = 'EVOL' 'MANU' 'x' LISTX 'Mach' LISTMACH ; 'SI' GRAPH ; 'DESSIN' EVMACH 'MIMA' 'TITRE' TITOLO ; 'FINSI' ; MMAX = 'MAXIMUM' MACHN ; MMIN = 'MINIMUM' MACHN ; MMAX_REF = 0.93624 ; MMIN_REF = 0.42153 ; ERRO = 'ABS' (MMAX_REF '-' MMAX) '/' M_INF ; 'SI' (ERRO > 5.0D-2) ; 'MESSAGE' 'Probleme' ; 'ERREUR' 5 ; 'FINSI' ; * ERRO = 'ABS' (MMIN_REF '-' MMIN) '/' M_INF ; 'SI' (ERRO > 1.0D-2) ; 'MESSAGE' 'Probleme' ; 'ERREUR' 5 ; 'FINSI' ; 'FIN' ;