* fichier : condmixtesVF.dgibi
*
********************* CAS TEST : transport1.dgibi ********************
*
GRAPH = 'N' ;
'OPTI' 'ECHO' 0 ;
CRIT1 = 2.D-2 ;
'SAUT' 'PAGE' ;
*
*---------------------------------------------------------------------
* TEST TRANSPORT - conditions mixtes
*
* CALCUL DARCY ISOTROPE TRANSPORT.
* Transport d'un front.
*
* dT
* -- + div (uT - Kgrad(T)) = 0
* dt
*
* Ce test permet de vérifier le bon fonctionnement des opérateurs
* utilisés afin de résoudre l'équation de transport par diffusion
* convection d'un champ scalaire passif par la méthode d'éléments
* finis mixtes hybrides implanté dans CASTEM2000.
*
* Le maillage ne comporte qu'un élément dans la direction transverse
* à l'écoulement puisque le phénomène étudié est 1D. Dans le sens de
* l'écoulement, le pas de discrétisation est constant.
*
* Les conditions initiales correspondent à l'entrée du front dans le
* domaine. Comme le champ de vitesse choisi est positif, la trace de
* concentration est égale à 1 à gauche du domaine.
*
* Les conditions aux limites imposent la trace de concentration à
* gauche du domaine égale à un via une condition mixte au lieu
* de conditions de Dirichlet; sur les autres frontières ce sont
* les conditions aux limites naturelles.
*
* Le schéma en temps utilisé est le schéma d'Euler implicite
* (teta-méthode avec teta=1.0).
*
* Données physiques :
* VK : Coefficient de diffusion ................. 1.D-2
* VX1 : Vitesse suivant x......................... 1.
* VY1 : Vitesse suivant y......................... 0.
* T0 : Concentration initiale.................... 0.
* T1 : Concentration en entrée................... 1.
*
* Données temporelles :
* TMIN : Temps initial............................. 0.
* TMAX : Temps final............................... 30.
* DELTAT: Pas de temps.............................. 1.
* T0 : Concentration initiale.................... 0.
* T1 : Concentration en entrée................... 1.
*
* Parametres de maillage :
* L : longueur du domaine ....................... 100.
* H : hauteur du domaine ........................ 1.
* INUMX : Nombre d'éléments en x .................... 100
* INUMY : Nombre d'éléments en y .................... 1 car 1D
*
*---------------------------------------------------------------------
*
*------------------
* Options generales
*------------------
*
'OPTI' 'DIME' 2 'ELEM' 'QUA4' ;
'OPTI' 'ISOV' 'SURF' ;
*
*
*=========
* MAILLAGE
*=========
*
*
*- Création des points supports du contour du domaine, et des droites
*- passant par les centres et les faces pour le post-traitement.
*
L = 100.D0 ;
LS2 = L / 2.D0 ;
H = 10.D0 ;
HS2 = H / 2.D0 ;
X0 = 0.D0 ;
X1 = X0 + L ;
Y0 = 0.D0 ;
Y1 = Y0 + H ;
INUMX = 100 ;
INUMY = 1 ;
INUM1 = INUMX - 1 ;
Y01 = Y0 + Y1 * 0.5D0 ;
DX = X1 - X0 / INUMX ;
DX1 = DX / 2.D0 ;
XG = X0 + DX1 ;
XD = X1 - DX1 ;
*
A1 = X0 Y0 ;
A3 = X1 Y0 ;
D1 = X0 Y1 ;
D3 = X1 Y1 ;
B1 = X0 Y01 ;
B3 = X1 Y01 ;
C1 = XG Y01 ;
C3 = XD Y01 ;
P6 = LS2 Y01 ;
*
*- Création des DROITES frontieres
*
DRBAS = A3 'DROI' INUMX A1 ;
DRGAU = A1 'DROI' INUMY D1 ;
DRHAU = D1 'DROI' INUMX D3 ;
DRDRO = D3 'DROI' INUMY A3 ;
PELIM = DX1 / (5. * INUMX) ;
*
*- Creation maillage GEOMETRIQUE
*
PTOT1 = 'DALL' DRBAS DRGAU DRHAU DRDRO ;
PTOT2 = 'ORIE' PTOT1 ;
*
*- Creation maillage HYBRIDE y compris sous-objets (cond. limites)
*
DRMID = B1 'DROI' INUMX B3 ;
DRMIC = C1 'DROI' INUM1 C3 ;
EXT1 = 'MANU' 'POI1' B1 ;
*
QFTOT= 'CHAN' PTOT2 QUAF ;
QFGAU= 'CHAN' DRGAU QUAF ;
QFDRO= 'CHAN' DRDRO QUAF ;
ELIM PELIM (QFTOT ET QFGAU ET QFDRO ET DRMID ET DRMIC ET EXT1) ;
*
*================
* INITIALISATIONS
*================
*
* ----------------
* = MODELISATION =
* ----------------
MODHYB = 'MODL' QFTOT 'DARCY' 'ANISOTROPE' ;
MODDRO = 'MODL' QFDRO 'DARCY' 'ANISOTROPE' ;
MODGAU = 'MODL' QFGAU 'DARCY' 'ANISOTROPE' ;
CHYB1 = 'DOMA' MODHYB 'SURFACE' ;
CHYB2 = 'DOMA' MODHYB 'NORMALE' ;
*
* ---------------------
* = Donnees physiques =
* ---------------------
*
T0 = 0.D0 ;
* concentration à gauche (pour le posttraitement)
CT1 = 1.D0;
* on prépare les conditions mixts équivalentes à C = 1 à gauche
T1 = 1.D0 ;
T1 = T1 * H / INUMY;
VX1 = 1.D0 ; VY1 = 0.D0 ;
VK = 1.D-2 ;
MATI2 = 'MATE' MODHYB DIRECTION (1. 0.)
'K11' VK 'K21' 0.D0 'K22' VK;
MATI2 = KCHA MODHYB MATI2 'CHPO';
*
SPEED = 'MANU' 'CHPO' ('DOMA' MODHYB 'FACE') 2
'VX' VX1 'VY' VY1 'NATURE' 'DISCRET' ;
*
* -----------------------
* = Donnees transitoire =
* -----------------------
* TETA : Parametre de le theta-méthode
* TMAX : Temps final
* TSUP : Temps pour conditions aux limites
* DELTAT : Pas de temps
*
TETA = 0.00D0 ;
TMIN = 0.D0 ;
TMAX2 = 15.00D0 ;
TMAX = 30.00D0 ;
TSUP = 1.2D0 * TMAX ;
DELTAT = 0.5D0 ;
*
LICALC = 'PROG' TMIN 'PAS' (DELTAT/2.D0) TMAX2 ;
LICALC = LICALC ET ('PROG' (TMAX2 + DELTAT) 'PAS' DELTAT TMAX) ;
LISAUV = 'PROG' TMAX ;
*
* ------------------------
* = Conditions initiales =
* ------------------------
TFHY0 = 'MANU' 'CHPO' ('DOMA' MODHYB 'FACE') 1 'TH' T0
'NATURE' 'DISCRET' ;
TFHY1 = 'MANU' 'CHPO' EXT1 1 'TH' T1 'NATURE' 'DISCRET' ;
TFHYB = TFHY0 'ET' TFHY1 ;
TCHYB = 'MANU' 'CHPO' ('DOMA' MODHYB 'CENTRE') 1 'H' T0
'NATURE' 'DISCRET' ;
*
* --------------
* = T imposée =
* --------------
BBGAU = 'BLOQ' ('DOMA' MODGAU 'CENTRE') 'TH' ;
EEGAU = 'DEPI' BBGAU T1 ;
TT1 = 'KCHT' MODGAU SCAL CENTRE (-1.* T1);
TTT1 = 'CHAR' TT1 ('EVOL' 'MANU' ('PROG' 0. TSUP) ('PROG' 1. 1.)) ;
CHAGAU = 'CHAR' EEGAU ('EVOL' 'MANU' ('PROG' 0. TSUP) ('PROG' 1. 1.)) ;
*
* ---------------------------
* = Table DARCY_TRANSITOIRE =
* ---------------------------
*
*
*
*-- Table de transport :
Transp = 'TABLE';
Transp . 'MODELE' = MODHYB ;
Transp.'TEMPS' = 'TABLE';
Transp.'CONCENTRATION' = 'TABLE';
Transp.'FLUXDIFF' = 'TABLE';
Transp.'FLUXCONV' = 'TABLE';
Transp.'CARACTERISTIQUES' = MATI2;
*Transp.'POROSITE' = MaPor;
Transp.'CONVECTION' = SPEED ;
Transp.'VITELEM' = kcht MODHYB VECT CENTRE (1.D0 0.);
Transp.'VITELEM' = nomc (mots UX UY) (mots VX VY)
Transp . 'VITELEM';
* Conditions initiales :
Transp.'TEMPS'. 0 = TMIN ;
Transp.'CONCENTRATION'. 0 = TCHYB;
Transp.'FLUXDIFF'. 0 = 0.D0 * TFHY0;
Transp.'FLUXCONV'. 0 = 0.D0 * TFHY0;
* Conditions aux limites :
*Transp . 'TRACE_IMPOSE' = TTT1 ;
*Transp . 'MIXTES' = TABLE;
Transp . 'MIXTES' = TTT1 ;
Transp . 'MIXCOFA' = 'KCHT' MODGAU SCAL CENTRE 1.D0;
Transp . 'MIXCOFB' = 'KCHT' MODGAU SCAL CENTRE -1.D0;
*Transp . 'FLUXTOT_IMP' = TTT1 ;
* Paramètres numériques :
Transp.'THETA_DIFF' = 1.D0;
Transp.'THETA_CONVECTION' = 1.0D0;
Transp.'LUMP' = FAUX;
Transp.'TYPDISCRETISATION' = 'VF';
Transp.'DECENTR' = FAUX;
TABRES = table METHINV;
TABRES . 'TYPINV' = 1;
TABRES . 'PRECOND' = 3;
Transp . 'METHINV' = TABRES;
Transp.'TEMPS_CALCULES' = LiCalc;
Transp.'TEMPS_SAUVES' = LiSauv;
Transp . INTCONC = TABLE;
Transp . INTCONC . 0 = 0.D0 * TCHYB;
* ==========
*Transp . 'DECROISSANCE' = 0.01D0;
* | CALCUL |
* ==========
*=======================
* Resolution transitoire
*=======================
*
TRANSGEN TRANSP ;
*
*
*=================
* POST-TRAITEMENT
*=================
TRANS2 = TABLE TRANSP;
*
*--------------------------------------------------------------------
* Dans chaque cas on trace
* La trace de concentration le long de DRMID
* La concentration le long de DRMIC
*--------------------------------------------------------------------
* Tests de NON-REGRESSION :
* Principe du maximum
* Position du centre de gravité du champ de concentration
*--------------------------------------------------------------------
*
* Critères numériques
CFL = VX1 * DELTAT / DX ;
PEK = VX1 * DX / (2. * VK) ;
FOU = 2. * VK * DELTAT / (DX * DX) ;
'SAUT' 1 'LIGNE' ;
'MESS' ' Critères numériques ' ;
'MESS' ' CFL ' CFL ;
'MESS' ' PECLET ' PEK ;
'MESS' ' FOURIER ' FOU ;
'MESS' ' ' ;
*
XC YC = 'COOR' (DOMA MODHYB 'CENTRE') ;
ISOR1 = INDEX ( TRANS2 . 'TEMPS') ;
NTSOR = DIME ISOR1 ;
NTSO1 = NTSOR - 1 ;
IOK = FAUX ;
IRESU = 1 ;
*
*-----------------------
REPETER VISURESU NTSO1 ;
*-----------------------
*
IRESU = IRESU + 1 ;
INDI1 = ISOR1.IRESU ;
TTRA = TRANS2 . 'TEMPS' . INDI1 ;
* Principe du maximum
EPR1 = TRANS2 . 'CONCENTRATION' . INDI1 ;
PMA = 'MAXI' EPR1 ; PMI = 'MINI' EPR1 ;
VTEST = PMA + PMI / 1.D0 ;
'SAUT' 1 'LIGNE' ;
'MESS' ' Principe du maximum ' ;
'MESS' ' Max et min théorique ' CT1 T0 ;
'MESS' ' Max et min par maille ' PMA PMI ;
'MESS' ' Précision demandée ' CRIT1 ;
'MESS' ' ' ;
* Centre de gravité
M0 = 'KCHA' MODHYB 'CHAM' EPR1 ;
M1 = 'INTG' MODHYB M0 ;
*
XV1 = XC * EPR1 ;
MXV1 = 'KCHA' MODHYB 'CHAM' XV1 ;
XV2 = 'INTG' MODHYB MXV1 ;
XV3 = XV2 / M1 ;
YV1 = YC * EPR1 ;
MYV1 = 'KCHA' MODHYB 'CHAM' YV1 ;
YV2 = 'INTG' MODHYB MYV1 ;
YV3 = YV2 / M1 ;
XTEST = TTRA / 2.D0 ;
*
'SAUT' 1 'LIGNE' ;
'MESS' ' Centre de gravité du nuage ' ;
'MESS' ' THEORIQUE ' XTEST HS2 ;
'MESS' ' OBTENU ' XV3 YV3 ;
'MESS' ' TAILLE D UNE MAILLE ' DX ;
'MESS' ' ' ;
'SAUT' 1 'LIGNE' ;
* Tests de non regression
'SI' ( VTEST 'NEG' CT1 CRIT1 ) ;
IOK = VRAI ;
mess 'VTEST';
'FINS' ;
'SI' ( XV3 'NEG' XTEST DX ) ;
IOK = VRAI ;
mess 'XV3';
'FINS' ;
'SI' ( YV3 'NEG' HS2 CRIT1 ) ;
IOK = VRAI ;
mess 'YV3';
'FINS' ;
*
'SI' ('NEG' GRAPH 'N' ) ;
LTI2 = 'CHAINE' 'Front 1D-h temps ' TTRA ;
'TITR' LTI2 ;
AV2 = 'EVOL' 'ROUG' 'CHPO' EPR1 'H' DRMIC ;
'DESS' AV2 'MIMA' ;
'FINS' ;
*
*-------------
FIN VISURESU ;
*-------------
*
*
'SI' ( IOK ) ;
mess 'erreur' ((VTEST - CT1) / CT1);
mess 'erreur' ((XV3 - XTEST) / XTEST);
mess 'erreur' ((YV3 - HS2) / HS2);
'ERRE' 5 ;
'SINO' ;
'ERRE' 0 ;
'FINSI' ;
*
'FIN' ;
© Cast3M 2003 - All rights reserved.
Disclaimer