* fichier : pente.dgibi
***********************************************************
**** Finite Volume, "Cell-Centred Formulation". ****
**** PENT, operator to compute gradients and limiters ****
**** ****
**** A. BECCANTINI, TTMF JANUARY 2001 ****
**** Modif, 10/07/01, syntaxe de PENT changée ****
**** Modif, 25/02/04, ajoute 'MODE' 'AXIS' ****
***********************************************************
'OPTION' 'ECHO' 0 'DIME' 2
'ELEM' 'TRI3'
'TRAC' 'X' ;
GRAPH = FAUX ;
* GRAPH = VRAI ;
COMPLET = VRAI ;
*
*** MESH
*
P0 = 0.0D0 0.0D0 ;
P1 = 3.0D0 0.0D0 ;
P2 = 3.0D0 3.0D0 ;
P3 = 0.0D0 3.0D0 ;
P4 = 6.0D0 0.0D0 ;
P5 = 6.0D0 3.0D0 ;
'SI' COMPLET ;
N1 = 80 ;
N2 = 60 ;
N3 = 40 ;
N4 = 90 ;
N5 = 90 ;
N6 = 10 ;
N7 = 50 ;
'SINON' ;
N1 = 8 ;
N2 = 7 ;
N3 = 6 ;
N4 = 9 ;
N5 = 9 ;
N6 = 7 ;
N7 = 8 ;
'FINSI' ;
LINEXT1 = ((P0 'DROIT' N1 P1) 'ET' (P1 'DROIT' N2 P2) 'ET'
(P2 'DROIT' N3 P3) 'ET' (P3 'DROIT' N4 P0)) ;
LINEXT2 = ((P1 'DROIT' N2 P2) 'ET' (P2 'DROIT' N5 P5) 'ET'
(P5 'DROIT' N6 P4) 'ET' (P4 'DROIT' N7 P1)) ;
'OPTION' 'ELEM' QUA4 ;
DOM1 = 'SURFACE' LINEXT1 'PLAN' ;
'OPTION' 'ELEM' TRI3 ;
DOM2 = 'SURFACE' LINEXT2 'PLAN' ;
DOMTOT = DOM1 'ET' DOM2;
'ELIMINATION' 1D-6 DOMTOT;
LTYPE = 'MOTS' 'TRI3' 'QUA4' ;
LELEM = 'NBEL' DOMTOT LTYPE ;
'MESSAGE' ;
'MESSAGE' 'Le maillage :';
'MESSAGE' ('CHAINE' 'Nombre TRI3: ' ('EXTRAIRE' LELEM 1));
'MESSAGE' ('CHAINE' 'Nombre QUA4: ' ('EXTRAIRE' LELEM 2));
*
**** Internal and external cells
*
DOMLIM = 'CONTOUR' DOMTOT 'COULEUR' 'VERTE' ;
MAICON = DOMTOT 'ELEM' 'APPUYE' 'LARGEMENT' DOMLIM ;
DOMINT = DOMTOT 'DIFF' MAICON 'COULEUR' 'ROUGE' ;
$DOMTOT = 'MODE' DOMTOT 'EULER' ;
$DOMINT = 'MODE' DOMINT 'EULER' ;
$DOM1 = 'MODE' DOM1 'EULER' ;
TDOMTOT = 'DOMA' $DOMTOT 'VF' ;
TDOMINT = 'DOMA' $DOMINT 'VF' ;
TDOM1 = 'DOMA' $DOM1 'VF' ;
MDOMTOT = TDOMTOT . 'QUAF' ;
MDOMINT = TDOMINT . 'QUAF' ;
MDOM1 = TDOM1 . 'QUAF' ;
'ELIM' (MDOMTOT 'ET' MDOMINT 'ET' MDOM1) 1D-6;
'SI' GRAPH ;
'TRACER' MAICON 'TITRE' 'Elements sur le contour';
'TRACER' DOMINT 'TITRE' 'Elements internes' ;
'FINSI' ;
*
**** The interface points (for limit condition)
*
POIN0 = DOMLIM 'POIN' 1;
X1 Y1 = 'COORDONNEE' POIN0;
'REPETER' BLLIM (('NBNO' DOMLIM)) ;
X0 = X1 ;
Y0 = Y1 ;
POIN0 = DOMLIM 'POIN' (&BLLIM '+' 1) ;
X1 Y1 = 'COORDONNEE' POIN0 ;
XFAC = (X0 '+' X1) '/' 2 ;
YFAC = (Y0 '+' Y1) '/' 2 ;
PFAC = ('DOMA' $DOMTOT 'FACE') 'POIN' 'PROC' (XFAC YFAC);
'SI' (&BLLIM 'EGA' 1) ;
GEOLIM = 'MANUEL' 'POI1' PFAC 'BLEU' ;
'SINON' ;
GEOLIM = GEOLIM 'ET' ('MANUEL' 'POI1' PFAC 'BLEU') ;
'FINSI' ;
'FIN' BLLIM ;
*
* 'ELIMINATION' GEOLIM 1.0D-4 ($DOMTOT . 'FACE') ;
*
'SI' GRAPH ;
'TRACER' (DOMTOT 'ET' GEOLIM)
'TITRE' 'DOMTOT et points faces aux bords';
'FINSI' ;
***********************************************
********* CHP: scalar linear field *********
********* EULESCAL *********
***********************************************
*
* (A11 '*' XX) '+' (A12 '*' YY) '+' A0;
*
A11 = 12.01517 ;
A12 = 13.1421 ;
A0 = -3.21 ;
B11 = 22.01517 ;
B12 = 43.1421 ;
B0 = -5.21 ;
XX YY = 'COORDONNEE' ('DOMA' $DOMTOT 'CENTRE') ;
CHP1 = (A11 '*' XX) '+' (A12 '*' YY) '+' A0;
CHP2 = (B11 '*' XX) '+' (B12 '*' YY) '+' B0;
CHP1 = 'NOMC' 'C1' CHP1 'NATU' 'DISCRET' ;
CHP2 = 'NOMC' 'C2' CHP2 'NATU' 'DISCRET' ;
CHP = CHP1 'ET' CHP2 ;
MOTC = ('MOTS' 'C1' 'C2') ;
GRCHP LIMCH GRG = 'PENT' $DOMTOT 'CENTRE'
'EULESCAL' 'NOLIMITE' MOTC CHP ;
'LISTE' MOTC ;
LIMCH1 = 'MANU' 'CHPO' ('DOMA' $DOMTOT 'CENTRE')
2 'P1' 1.0 'P2' 1.0 ;
ERRO = 'MAXIMUM' (LIMCH1 '-' LIMCH) 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
*
* TEST1 : we check that we have the same values if we use
* GRG to compute gradients
*
MOTC = ('MOTS' 'C1' 'C2') ;
GRCHP1 LIMCH = 'PENT' $DOMTOT 'CENTRE'
'EULESCAL' 'LIMITEUR' MOTC
CHP 'GRADGEO' GRG ;
'LISTE' MOTC ;
ERRO = 'MAXIMUM' (GRCHP '-' GRCHP1) 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
*
* TEST2 : we check that the computation is exact on the internal domain
*
GRCHPI = 'REDU' GRCHP ('DOMA' $DOMINT 'CENTRE') ;
GR1DX = 'MANUEL' 'CHPO' ('DOMA' $DOMINT 'CENTRE') 1 'P1DX' A11
'NATU' 'DISCRET' ;
GR1DY = 'MANUEL' 'CHPO' ('DOMA' $DOMINT 'CENTRE') 1 'P1DY' A12
'NATU' 'DISCRET';
GR2DX = 'MANUEL' 'CHPO' ('DOMA' $DOMINT 'CENTRE') 1 'P2DX' B11
'NATU' 'DISCRET';
GR2DY = 'MANUEL' 'CHPO' ('DOMA' $DOMINT 'CENTRE') 1 'P2DY' B12
'NATU' 'DISCRET';
GREXA = GR1DX 'ET' GR1DY 'ET' GR2DX 'ET' GR2DY ;
ERRCH = GREXA '-' GRCHPI ;
LMOT = 'MOTS' 'P1DX' 'P1DY' 'P2DX' 'P2DY' ;
MERRCH = ('PSCAL' ERRCH ERRCH LMOT LMOT) '**' 0.5 ;
*
*** Error graphics (we change the error CHAMPOINT into a CHAMELEM
*
'SI' GRAPH ;
MOD1 = 'MODELISER' ('DOMA' $DOMINT 'MAILLAGE') 'THERMIQUE' ;
CHM_MER = 'KCHA' $DOMINT 'CHAM' MERRCH ;
'TRACER' CHM_MER MOD1 'TITR' ('CHAINE' 'EULESCAL: erreur');
'FINSI' ;
ERRO = 'MAXIMUM' ERRCH 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
*
* TEST3 : we impose the linear field values on the border.
* In that case error must be "zero" everywhere.
*
XLIM YLIM = 'COORDONNEE' GEOLIM ;
CHPLI1 = (A11 '*' XLIM) '+' (A12 '*' YLIM) '+' A0;
CHPLI2 = (B11 '*' XLIM) '+' (B12 '*' YLIM) '+' B0;
CHPLI1 = 'NOMC' 'C1' CHPLI1 'NATU' 'DISCRET' ;
CHPLI2 = 'NOMC' 'C2' CHPLI2 'NATU' 'DISCRET' ;
CHPLI = CHPLI1 'ET' CHPLI2 ;
MOTC = ('MOTS' 'C1' 'C2') ;
GRCHPL LIMCH GRGL = 'PENT' $DOMTOT 'CENTRE'
'EULESCAL' 'NOLIMITE' MOTC CHP 'CLIM' CHPLI ;
MOTC = ('MOTS' 'C1' 'C2') ;
GRCHP1 LIMCH = 'PENT' $DOMTOT 'CENTRE'
'EULESCAL' 'NOLIMITE' MOTC
CHP 'CLIM' CHPLI 'GRADGEO' GRGL ;
ERRO = 'MAXIMUM' (GRCHPL '-' GRCHP1) 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
GR1DX = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'P1DX' A11
'NATU' 'DISCRET' ;
GR1DY = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'P1DY' A12
'NATU' 'DISCRET';
GR2DX = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'P2DX' B11
'NATU' 'DISCRET';
GR2DY = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'P2DY' B12
'NATU' 'DISCRET';
GREXA = GR1DX 'ET' GR1DY 'ET' GR2DX 'ET' GR2DY ;
ERRCH = GREXA '-' GRCHPL ;
LMOT = 'MOTS' 'P1DX' 'P1DY' 'P2DX' 'P2DY' ;
MERRCH = ('PSCAL' ERRCH ERRCH LMOT LMOT) '**' 0.5 ;
'SI' GRAPH ;
MOD1 = 'MODELISER' ('DOMA' $DOMTOT 'MAILLAGE') 'THERMIQUE' ;
CHM_MER = 'KCHA' $DOMTOT 'CHAM' MERRCH ;
'TRACER' CHM_MER MOD1 'TITR'
('CHAINE' 'EULESCAL avec C.L. : Erreur');
* CHM_GRA = 'KCHA' $DOMTOT 'CHAM' GRCHPL ;
* 'TRACER' CHM_GRA MOD1 'TITR' ('CHAINE' 'Gradient');
'FINSI' ;
ERRO = 'MAXIMUM' ERRCH 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
*
**** TEST4: boudary condition for a scalar field (EULESCAL)
*
'OPTION' 'ELEM' 'QUA4' ;
L0 = (P0 'DROIT' P3 'DINI' 0.3 'DFIN' 1 ) ;
DOMC0 = L0 'TRANSLATION' 1 (-1.0 0) ;
DOMC1 = L0 'TRANSLATION' 5 (5.0 0) 'COULEUR' 'ROUGE' ;
DOMC2 = DOMC0 'ET' DOMC1 ;
'ELIMINATION' DOMC2 0.0002 ;
'SI' GRAPH ;
'TRACER' (DOMC1 'ET' DOMC2) 'TITRE' 'Traitement de C.L.' ;
'FINSI' ;
$DOMC1 = 'MODE' DOMC1 'EULER' ;
$DOMC2 = 'MODE' DOMC2 'EULER' ;
TDOMC1 = 'DOMA' $domc1 'VF' ;
TDOMC2 = 'DOMA' $domc2 'VF' ;
MDOMC1 = TDOMC1 . 'QUAF';
MDOMC2 = TDOMC2 . 'QUAF';
'ELIM' (MDOMC1 'ET' MDOMC2) 1D-6;
*
* CHP is symmetric with respect to the line L0
*
XX = 'COORDONNEE' 1 ('DOMA' $DOMC2 'CENTRE') ;
CHP = ('NOMC' 'C1' (XX '*' XX) 'NATU' 'DISCRET') 'ET'
('NOMC' 'C2' (4 '*' XX '*' XX) 'NATU' 'DISCRET') ;
GRCHP LIMCHP GRG = 'PENT' $DOMC2 'CENTRE'
'EULESCAL' 'NOLIMITE' ('MOTS' 'C1' 'C2') CHP ;
GRCHP1 LIMCHP GRG1 = 'PENT' $DOMC1 'CENTRE'
'EULESCAL' 'NOLIMITE' ('MOTS' 'C1' 'C2')
('REDU' CHP ('DOMA' $DOMC1 'CENTRE')) ;
ERRCH = GRCHP '-' ('REDU' GRCHP1 ('DOMA' $DOMC1 'CENTRE'));
LMOT = 'MOTS' 'P1DX' 'P1DY' 'P2DX' 'P2DY' ;
MERRCH = ('PSCAL' ERRCH ERRCH LMOT LMOT) '**' 0.5 ;
'SI' GRAPH ;
MOD1 = 'MODELISER' ('DOMA' $DOMC1 'MAILLAGE') 'THERMIQUE';
CHM_MER = 'KCHA' $DOMC1 'CHAM' MERRCH ;
'TRACER' CHM_MER MOD1
'TITR' ('CHAINE' 'EULESCAL sur les bords: erreur');
'FINSI' ;
ERRO = 'MAXIMUM' ERRCH 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
***************************************************
********* CHP: vectorial linear field ************
********* EULEVECT ****************
***************************************************
*
* (A11 '*' XX) '+' (A12 '*' YY) '+' A0;
*
A11 = 11.3 ;
A12 = 0.23 ;
A0 = 1.0 ;
B11 = 77.0 ;
B12 = 0.6 ;
B0 = 9.11 ;
XX YY = 'COORDONNEE' ('DOMA' $DOMTOT 'CENTRE') ;
CHP1 = (A11 '*' XX) '+' (A12 '*' YY) '+' A0;
CHP2 = (B11 '*' XX) '+' (B12 '*' YY) '+' B0;
CHP1 = 'NOMC' 'UX ' CHP1 'NATU' 'DISCRET' ;
CHP2 = 'NOMC' 'UY ' CHP2 'NATU' 'DISCRET' ;
CHP = CHP1 'ET' CHP2 ;
MOTC = ('MOTS' 'UX' 'UY') ;
GRCHP LIMCHP GRG = 'PENT' $DOMTOT 'CENTRE'
'EULEVECT' 'NOLIMITE' MOTC
CHP ;
'LISTE' MOTC ;
*
* TEST1 : we check that we have the same values if we use
* GRG to compute gradients
*
GRCHP1 LIMCHP = 'PENT' $DOMTOT 'CENTRE'
'EULEVECT' 'NOLIMITE' ('MOTS' 'UX' 'UY')
CHP 'GRADGEO' GRG ;
ERRO = 'MAXIMUM' (GRCHP '-' GRCHP1) 'ABS' ;
'SI' (ERRO > 5.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
*
* TEST2 : we check that the computation is exact on the internal domain
*
GRCHPI = 'REDU' GRCHP ('DOMA' $DOMINT 'CENTRE') ;
GR1DX = 'MANUEL' 'CHPO' ('DOMA' $DOMINT 'CENTRE') 1 'P1DX' A11
'NATU' 'DISCRET' ;
GR1DY = 'MANUEL' 'CHPO' ('DOMA' $DOMINT 'CENTRE') 1 'P1DY' A12
'NATU' 'DISCRET';
GR2DX = 'MANUEL' 'CHPO' ('DOMA' $DOMINT 'CENTRE') 1 'P2DX' B11
'NATU' 'DISCRET';
GR2DY = 'MANUEL' 'CHPO' ('DOMA' $DOMINT 'CENTRE') 1 'P2DY' B12
'NATU' 'DISCRET';
GREXA = GR1DX 'ET' GR1DY 'ET' GR2DX 'ET' GR2DY ;
ERRCH = GREXA '-' GRCHPI ;
LMOT = 'MOTS' 'P1DX' 'P1DY' 'P2DX' 'P2DY' ;
MERRCH = ('PSCAL' ERRCH ERRCH LMOT LMOT) '**' 0.5 ;
*
*** CREATION DE 'MODEL' POUR GRAPHIQUER LE CHAMELEM ***
*
'SI' GRAPH ;
MOD1 = 'MODELISER' ('DOMA' $DOMINT 'MAILLAGE') 'THERMIQUE';
CHM_MER = 'KCHA' $DOMINT 'CHAM' MERRCH ;
'TRACER' CHM_MER MOD1 'TITR' ('CHAINE' 'EULEVECT: Erreur');
'FINSI' ;
ERRO = 'MAXIMUM' ERRCH 'ABS' ;
'SI' (ERRO > 5.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
*
* TEST3 : we impose the linear field values on the border.
* In that case error must be "zero" everywhere.
*
XLIM YLIM = 'COORDONNEE' GEOLIM ;
CHPLI1 = (A11 '*' XLIM) '+' (A12 '*' YLIM) '+' A0;
CHPLI2 = (B11 '*' XLIM) '+' (B12 '*' YLIM) '+' B0;
CHPLI1 = 'NOMC' 'UX' CHPLI1 'NATU' 'DISCRET' ;
CHPLI2 = 'NOMC' 'UY' CHPLI2 'NATU' 'DISCRET' ;
CHPLI = CHPLI1 'ET' CHPLI2 ;
GRCHPL LIMCHP GRGL = 'PENT' $DOMTOT 'CENTRE'
'EULEVECT' 'NOLIMITE' ('MOTS' 'UX' 'UY')
CHP 'CLIM' CHPLI ;
GRCHP1 LIMCHP = 'PENT' $DOMTOT 'CENTRE'
'EULEVECT' 'NOLIMITE' ('MOTS' 'UX' 'UY')
CHP 'CLIM' CHPLI 'GRADGEO' GRGL ;
ERRO = 'MAXIMUM' (GRCHPL '-' GRCHP1) 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
GR1DX = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'P1DX' A11
'NATU' 'DISCRET' ;
GR1DY = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'P1DY' A12
'NATU' 'DISCRET';
GR2DX = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'P2DX' B11
'NATU' 'DISCRET';
GR2DY = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'P2DY' B12
'NATU' 'DISCRET';
GREXA = GR1DX 'ET' GR1DY 'ET' GR2DX 'ET' GR2DY ;
ERRCH = GREXA '-' GRCHPL ;
LMOT = 'MOTS' 'P1DX' 'P1DY' 'P2DX' 'P2DY' ;
MERRCH = ('PSCAL' ERRCH ERRCH LMOT LMOT) '**' 0.5 ;
'SI' GRAPH ;
MOD1 = 'MODELISER' ('DOMA' $DOMTOT 'MAILLAGE') 'THERMIQUE';
CHM_MER = 'KCHA' $DOMTOT 'CHAM' MERRCH ;
'TRACER' CHM_MER MOD1 'TITR'
('CHAINE' 'EULEVECT avec C.L.: erreur');
* CHM_GRA = 'KCHA' $DOMTOT 'CHAM' GRCHPL ;
* 'TRACER' CHM_GRA MOD1 'TITR' ('CHAINE' 'Gradient');
'FINSI' ;
ERRO = 'MAXIMUM' ERRCH 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
*
**** TEST4: boudary condition for a vectorial field (EULESCAL)
*
'OPTION' 'ELEM' 'QUA4' ;
L0 = (P0 'DROIT' P3 'DINI' 0.3 'DFIN' 1 ) ;
DOMC0 = L0 'TRANSLATION' 1 (-1.0 0) ;
DOMC1 = L0 'TRANSLATION' 5 (5.0 0) 'COULEUR' 'ROUGE' ;
DOMC2 = DOMC0 'ET' DOMC1 ;
'ELIMINATION' DOMC2 0.0002 ;
'SI' GRAPH ;
'TRACER' DOMC2 'TITRE' 'Traitement de C.L.' ;
'FINSI' ;
$DOMC1 = 'MODE' DOMC1 'EULER' ;
$DOMC2 = 'MODE' DOMC2 'EULER' ;
TDOMC1 = 'DOMA' $domc1 'VF' ;
TDOMC2 = 'DOMA' $domc2 'VF' ;
MDOMC1 = TDOMC1 . 'QUAF';
MDOMC2 = TDOMC2 . 'QUAF';
'ELIM' (MDOMC1 'ET' MDOMC2) 1D-6;
*
* L0: symmetry line
*
XX = 'COORDONNEE' 1 ('DOMA' $DOMC2 'CENTRE') ;
CHP = ('NOMC' 'UX' (1.15 '*' XX '*' XX '*' XX) 'NATU' 'DISCRET') 'ET'
('NOMC' 'UY' (0 '*' XX '*' XX) 'NATU' 'DISCRET') ;
GRCHP LIMCHP GRG = 'PENT' $DOMC2 'CENTRE'
'EULEVECT' 'NOLIMITE' ('MOTS' 'UX' 'UY') CHP ;
GRCHP1 LIMCHP GRG1 = 'PENT' $DOMC1 'CENTRE'
'EULEVECT' 'NOLIMITE' ('MOTS' 'UX' 'UY')
('REDU' CHP ('DOMA' $DOMC1 'CENTRE')) ;
ERRCH = GRCHP1 '-' ('REDU' GRCHP ('DOMA' $DOMC1 'CENTRE'));
LMOT = 'MOTS' 'P1DX' 'P1DY' 'P2DX' 'P2DY' ;
MERRCH = ('PSCAL' ERRCH ERRCH LMOT LMOT) '**' 0.5 ;
'SI' GRAPH ;
MOD1 = 'MODELISER' ('DOMA' $DOMC1 'MAILLAGE') 'THERMIQUE';
CHM_MER = 'KCHA' $DOMC1 'CHAM' MERRCH ;
'TRACER' CHM_MER MOD1
'TITR' ('CHAINE' 'EULEVECT sur les bords: erreur');
* CHM_GRA = 'KCHA' $DOMC1 'CHAM' GRCHP1 ;
* 'TRACER' CHM_GRA MOD1 'TITR' ('CHAINE' 'Gradient');
'FINSI' ;
ERRO = 'MAXIMUM' ERRCH 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
****************************
**** Diamond method ********
****************************
*
* 2D
*
*** MESH
*
P0 = 0.0D0 0.0D0 ;
P1 = 3.0D0 0.0D0 ;
P2 = 3.0D0 3.0D0 ;
P3 = 0.0D0 3.0D0 ;
P4 = 6.0D0 0.0D0 ;
P5 = 6.0D0 3.0D0 ;
'SI' COMPLET ;
N1 = 80 ;
N2 = 60 ;
N3 = 40 ;
N4 = 90 ;
N5 = 90 ;
N6 = 10 ;
N7 = 50 ;
'SINON' ;
N1 = 8 ;
N2 = 6 ;
N3 = 4 ;
N4 = 9 ;
N5 = 9 ;
N6 = 3 ;
N7 = 5 ;
'FINSI' ;
LINEXT1 = ((P0 'DROIT' N1 P1) 'ET' (P1 'DROIT' N2 P2) 'ET'
(P2 'DROIT' N3 P3) 'ET' (P3 'DROIT' N4 P0)) ;
LINEXT2 = ((P1 'DROIT' N2 P2) 'ET' (P2 'DROIT' N5 P5) 'ET'
(P5 'DROIT' N6 P4) 'ET' (P4 'DROIT' N7 P1)) ;
'OPTION' 'ELEM' QUA4 ;
DOM1 = 'SURFACE' LINEXT1 'PLAN' ;
'OPTION' 'ELEM' TRI3 ;
DOM2 = 'SURFACE' LINEXT2 'PLAN' ;
DOMTOT = DOM1 'ET' DOM2;
'ELIMINATION' 1D-6 DOMTOT;
LTYPE = 'MOTS' 'TRI3' 'QUA4' ;
LELEM = 'NBEL' DOMTOT LTYPE ;
'MESSAGE' ;
'MESSAGE' 'Le maillage :';
'MESSAGE' ('CHAINE' 'Nombre TRI3: ' ('EXTRAIRE' LELEM 1));
'MESSAGE' ('CHAINE' 'Nombre QUA4: ' ('EXTRAIRE' LELEM 2));
*
**** Internal and external cells
*
DOMLIM = 'CONTOUR' DOMTOT 'COULEUR' 'VERTE' ;
$DOMTOT = 'MODE' DOMTOT 'EULER' ;
TDOMTOT = 'DOMA' $DOMTOT 'VF' ;
MDOMTOT = TDOMTOT . 'QUAF';
'ELIM' MDOMTOT 1D-6;
'SI' GRAPH ;
'TRACER' DOMTOT 'TITRE' 'Domaine' ;
'FINSI' ;
A11 = 12.01517 ;
A12 = 13.1421 ;
A0 = -3.21 ;
B11 = 22.01517 ;
B12 = 43.1421 ;
B0 = -5.21 ;
XX YY = 'COORDONNEE' ('DOMA' $DOMTOT 'CENTRE') ;
CHP1 = (A11 '*' XX) '+' (A12 '*' YY) '+' A0;
CHP2 = (B11 '*' XX) '+' (B12 '*' YY) '+' B0;
CHP1 = 'NOMC' 'C1' CHP1 'NATU' 'DISCRET' ;
CHP2 = 'NOMC' 'C2' CHP2 'NATU' 'DISCRET' ;
CHP = CHP1 'ET' CHP2 ;
GRCHP GRG = 'PENT' $DOMTOT 'FACE'
'DIAMANT' CHP ;
*
* TEST1 : we check that we have the same values if we use
* GRG to compute gradients
*
GRCHP1 = 'PENT' $DOMTOT 'FACE'
'DIAMANT' CHP 'GRADGEO' GRG ;
ERRO = 'MAXIMUM' (GRCHP '-' GRCHP1) 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
*
* TEST2 : we check that the computation is exact on the internal domain
*
GR1DX = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'FACE') 1 'P1DX' A11
'NATU' 'DISCRET' ;
GR1DY = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'FACE') 1 'P1DY' A12
'NATU' 'DISCRET';
GR2DX = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'FACE') 1 'P2DX' B11
'NATU' 'DISCRET';
GR2DY = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'FACE') 1 'P2DY' B12
'NATU' 'DISCRET';
GREXA = GR1DX 'ET' GR1DY 'ET' GR2DX 'ET' GR2DY ;
ERRCH = GREXA '-' GRCHP1 ;
LMOT = 'MOTS' 'P1DX' 'P1DY' 'P2DX' 'P2DY' ;
MERRCH = ('PSCAL' ERRCH ERRCH LMOT LMOT) '**' 0.5 ;
ERRO = 'MAXIMUM' ERRCH 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
*
* TEST3 : we impose the linear field values on the border.
* In that case error must be "zero" everywhere.
*
XLIM YLIM = 'COORDONNEE' DOMLIM ;
CHPLI1 = (A11 '*' XLIM) '+' (A12 '*' YLIM) '+' A0;
CHPLI2 = (B11 '*' XLIM) '+' (B12 '*' YLIM) '+' B0;
CHPLI1 = 'NOMC' 'C1' CHPLI1 'NATU' 'DISCRET' ;
CHPLI2 = 'NOMC' 'C2' CHPLI2 'NATU' 'DISCRET' ;
CHPLI = CHPLI1 'ET' CHPLI2 ;
GRCHPL GRGL = 'PENT' $DOMTOT 'FACE'
'DIAMANT' CHP 'CLIM' CHPLI ;
GRCHP1 = 'PENT' $DOMTOT 'FACE'
'DIAMANT' CHP 'CLIM' CHPLI 'GRADGEO' GRGL ;
ERRO = 'MAXIMUM' (GRCHPL '-' GRCHP1) 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
GR1DX = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'FACE') 1 'P1DX' A11
'NATU' 'DISCRET' ;
GR1DY = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'FACE') 1 'P1DY' A12
'NATU' 'DISCRET';
GR2DX = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'FACE') 1 'P2DX' B11
'NATU' 'DISCRET';
GR2DY = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'FACE') 1 'P2DY' B12
'NATU' 'DISCRET';
GREXA = GR1DX 'ET' GR1DY 'ET' GR2DX 'ET' GR2DY ;
ERRCH = GREXA '-' GRCHPL ;
LMOT = 'MOTS' 'P1DX' 'P1DY' 'P2DX' 'P2DY' ;
MERRCH = ('PSCAL' ERRCH ERRCH LMOT LMOT) '**' 0.5 ;
ERRO = 'MAXIMUM' ERRCH 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
****************************
**** Diamond method ********
****************************
* The BC on flux are taken into account
*
* 2D
*
*** MESH
*
'OPTION' 'ELEM' QUA4 ;
P0 = 0.0D0 0.0D0 ;
P1 = 3.0D0 0.0D0 ;
P2 = 3.0D0 3.0D0 ;
P3 = 0.0D0 3.0D0 ;
P4 = 6.0D0 0.0D0 ;
P5 = 6.0D0 3.0D0 ;
'SI' COMPLET ;
N1 = 80 ;
N2 = 60 ;
N3 = 40 ;
N4 = 90 ;
N5 = 90 ;
N6 = 10 ;
N7 = 50 ;
'SINON' ;
N1 = 8 ;
N2 = 6 ;
N3 = 4 ;
N4 = 9 ;
N5 = 9 ;
N6 = 3 ;
N7 = 5 ;
'FINSI' ;
LINEXT1 = ((P0 'DROIT' N1 P1) 'ET' (P1 'DROIT' N2 P2) 'ET'
(P2 'DROIT' N3 P3) 'ET' (P3 'DROIT' N4 P0)) ;
LINEXT2 = ((P1 'DROIT' N2 P2) 'ET' (P2 'DROIT' N5 P5) 'ET'
(P5 'DROIT' N6 P4) 'ET' (P4 'DROIT' N7 P1)) ;
'OPTION' 'ELEM' QUA4 ;
DOM1 = 'SURFACE' LINEXT1 'PLAN' ;
'OPTION' 'ELEM' TRI3 ;
DOM2 = 'SURFACE' LINEXT2 'PLAN' ;
DOMTOT = DOM1 'ET' DOM2;
'ELIMINATION' 1D-6 DOMTOT;
* BC
DOMLIM = 'CONTOUR' DOMTOT 'COULEUR' 'VERTE' ;
DOMLIM1 = (DOMLIM 'INTERSECTION' LINEXT1) ;
DOMLIM2 = (DOMLIM 'INTERSECTION' LINEXT2) ;
* MODELS and MESHES
*
MDOMTOT = 'MODE' DOMTOT 'EULER' ;
MDOMLIM = 'MODE' DOMLIM 'EULER' ;
MDOMLIM1 = 'MODE' DOMLIM1 'EULER' ;
MDOMLIM2 = 'MODE' DOMLIM2 'EULER' ;
*
TDOMTOT = 'DOMA' MDOMTOT 'VF' ;
TDOMLIM = 'DOMA' MDOMLIM 'VF' ;
TDOMLIM1 = 'DOMA' MDOMLIM1 'VF' ;
TDOMLIM2 = 'DOMA' MDOMLIM2 'VF' ;
QDOMTOT = TDOMTOT . 'QUAF' ;
QDOMLIM = TDOMLIM . 'QUAF' ;
QDOMLIM1 = TDOMLIM1 . 'QUAF' ;
QDOMLIM2 = TDOMLIM2 . 'QUAF' ;
'ELIMINATION' (QDOMTOT 'ET' QDOMLIM 'ET' QDOMLIM1 'ET' QDOMLIM2)
1.0D-6 ;
'SI' GRAPH ;
'TRACER' (DOMTOT 'ET' ('DOMA' MDOMTOT 'CENTRE')
'ET' ('DOMA' MDOMTOT 'FACE')) 'TITRE' 'Centres et faces' ;
'FINSI' ;
* CHAMPOINT
A11 = 12.01517 ;
A12 = 13.1421 ;
A0 = -3.21 ;
B11 = 22.01517 ;
B12 = 43.1421 ;
B0 = -5.21 ;
XX YY = 'COORDONNEE' ('DOMA' MDOMTOT 'CENTRE') ;
CHP1 = (A11 '*' XX) '+' (A12 '*' YY) '+' A0;
CHP2 = (B11 '*' XX) '+' (B12 '*' YY) '+' B0;
CHP1 = 'NOMC' 'C1' CHP1 'NATU' 'DISCRET' ;
CHP2 = 'NOMC' 'C2' CHP2 'NATU' 'DISCRET' ;
CHP = CHP1 'ET' CHP2 ;
*
* TEST1 : Dirichlet BC everywhere
*
XXF YYF = 'COORDONNEE' ('DOMA' MDOMLIM 'CENTRE') ;
CHP1 = (A11 '*' XXF) '+' (A12 '*' YYF) '+' A0;
CHP2 = (B11 '*' XXF) '+' (B12 '*' YYF) '+' B0;
CHP1 = 'NOMC' 'C1' CHP1 'NATU' 'DISCRET' ;
CHP2 = 'NOMC' 'C2' CHP2 'NATU' 'DISCRET' ;
CHPL1 = CHP1 'ET' CHP2 ;
CACCA PIPI = 'KOPS' 'MATRIK' ;
CHPL2 = 'COPIER' CACCA ;
*
MOTC = ('MOTS' 'C1' 'C2') ;
MOTC1 = ('MOTS' 'C1DX' 'C1DY' 'C2DX' 'C2DY') ;
GRCHP GRG = 'PENT' MDOMTOT 'FACE'
'DIAMAN2' MOTC MOTC1
CHP CHPL1 CHPL2 ;
'LISTE' MOTC ;
'LISTE' MOTC1 ;
GRCHP1 = 'PENT' MDOMTOT 'FACE'
'DIAMAN2' ('MOTS' 'C1' 'C2') ('MOTS' 'C1DX' 'C1DY' 'C2DX' 'C2DY')
CHP CHPL1 CHPL2 'GRADGEO' GRG ;
ERRO = 'MAXIMUM' (GRCHP '-' GRCHP1) 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
GREXA = 'MANUEL' 'CHPO' ('DOMA' MDOMTOT 'FACE') 4 'C1DX' A11
'C1DY' A12 'C2DX' B11 'C2DY' B12 'NATU' 'DISCRET' ;
ERRCH = GREXA '-' GRCHP ;
ERRO = 'MAXIMUM' ERRCH 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
*
* TEST2 : von Neumann BC everywhere
*
CHPL2 = 'MANUEL' 'CHPO' ('DOMA' MDOMLIM 'CENTRE') 4
'C1DX' A11 'C1DY' A12 'C2DX' B11 'C2DY' B12 ;
CACCA PIPI = 'KOPS' 'MATRIK' ;
CHPL1 = CACCA ;
GRCHP GRG = 'PENT' MDOMTOT 'FACE'
'DIAMAN2' ('MOTS' 'C1' 'C2') ('MOTS' 'C1DX' 'C1DY' 'C2DX' 'C2DY')
CHP CHPL1 CHPL2 ;
ERRO = 'MAXIMUM' (GRCHP '-' GREXA) 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
*
* TEST3 : Dirichlet BC and von Neumann BC
*
XXF YYF = 'COORDONNEE' ('DOMA' MDOMLIM1 'CENTRE') ;
CHP1 = (A11 '*' XXF) '+' (A12 '*' YYF) '+' A0;
CHP2 = (B11 '*' XXF) '+' (B12 '*' YYF) '+' B0;
CHP1 = 'NOMC' 'C1' CHP1 'NATU' 'DISCRET' ;
CHP2 = 'NOMC' 'C2' CHP2 'NATU' 'DISCRET' ;
CHPL1 = CHP1 'ET' CHP2 ;
CHPL2 = 'MANUEL' 'CHPO' ('DOMA' MDOMLIM2 'CENTRE') 4
'C1DX' A11 'C1DY' A12 'C2DX' B11 'C2DY' B12 ;
GRCHP GRG = 'PENT' MDOMTOT 'FACE'
'DIAMAN2' ('MOTS' 'C1' 'C2') ('MOTS' 'C1DX' 'C1DY' 'C2DX' 'C2DY')
CHP CHPL1 CHPL2 ;
ERRCH = GREXA '-' GRCHP ;
ERRO = 'MAXIMUM' ERRCH 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
*
***** Elementary mesh
*
P0 = 0.0D0 0.0D0 ;
P1 = 1.0D0 0.0D0 ;
P2 = 1.0D0 1.0D0 ;
P3 = 0.0D0 1.0D0 ;
LINEXT1 = ((P0 'DROIT' 1 P1) 'ET' (P1 'DROIT' 1 P2) 'ET'
(P2 'DROIT' 1 P3) 'ET' (P3 'DROIT' 1 P0)) ;
'OPTION' 'ELEM' 'QUA4' ;
DOMTOT = 'SURFACE' LINEXT1 'PLAN' ;
* BC
DOMLIM = 'CONTOUR' DOMTOT 'COULEUR' 'VERTE' ;
* MODELS and MESHES
MDOMTOT = 'MODE' DOMTOT 'EULER' ;
MDOMLIM = 'MODE' DOMLIM 'EULER' ;
TDOMTOT = 'DOMA' MDOMTOT 'VF' ;
TDOMLIM = 'DOMA' MDOMLIM 'VF' ;
* QUAF
QDOMTOT = TDOMTOT . 'QUAF' ;
QDOMLIM = TDOMLIM . 'QUAF' ;
'ELIMINATION' (QDOMTOT 'ET' QDOMLIM)
1.0D-6 ;
*
'SI' GRAPH ;
'TRACER' (DOMTOT 'ET' ('DOMA' MDOMTOT 'CENTRE')
'ET' ('DOMA' MDOMTOT 'FACE')) 'TITRE' 'Centres et faces' ;
'FINSI' ;
* CHAMPOINT
A11 = 12.01517 ;
A12 = 13.1421 ;
A0 = -3.21 ;
B11 = 22.01517 ;
B12 = 43.1421 ;
B0 = -5.21 ;
XX YY = 'COORDONNEE' ('DOMA' MDOMTOT 'CENTRE') ;
CHP1 = (A11 '*' XX) '+' (A12 '*' YY) '+' A0;
CHP2 = (B11 '*' XX) '+' (B12 '*' YY) '+' B0;
CHP1 = 'NOMC' 'C1' CHP1 'NATU' 'DISCRET' ;
CHP2 = 'NOMC' 'C2' CHP2 'NATU' 'DISCRET' ;
CHP = CHP1 'ET' CHP2 ;
*
* TEST1 : Dirichlet BC everywhere
*
XXF YYF = 'COORDONNEE' ('DOMA' MDOMLIM 'CENTRE') ;
CHP1 = (A11 '*' XXF) '+' (A12 '*' YYF) '+' A0;
CHP2 = (B11 '*' XXF) '+' (B12 '*' YYF) '+' B0;
CHP1 = 'NOMC' 'C1' CHP1 'NATU' 'DISCRET' ;
CHP2 = 'NOMC' 'C2' CHP2 'NATU' 'DISCRET' ;
CHPL1 = CHP1 'ET' CHP2 ;
CACCA PIPI = 'KOPS' 'MATRIK' ;
CHPL2 = 'COPIER' CACCA ;
*
GRCHP GRG = 'PENT' MDOMTOT 'FACE'
'DIAMAN2' ('MOTS' 'C1' 'C2') ('MOTS' 'C1DX' 'C1DY' 'C2DX' 'C2DY')
CHP CHPL1 CHPL2 ;
GRCHP1 = 'PENT' MDOMTOT 'FACE'
'DIAMAN2' ('MOTS' 'C1' 'C2') ('MOTS' 'C1DX' 'C1DY' 'C2DX' 'C2DY')
CHP CHPL1 CHPL2 'GRADGEO' GRG ;
ERRO = 'MAXIMUM' (GRCHP '-' GRCHP1) 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
GREXA = 'MANUEL' 'CHPO' ('DOMA' MDOMTOT 'FACE') 4 'C1DX' A11
'C1DY' A12 'C2DX' B11 'C2DY' B12 'NATU' 'DISCRET' ;
ERRCH = GREXA '-' GRCHP ;
ERRO = 'MAXIMUM' ERRCH 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
*
* TEST2 : von Neumann BC everywhere
*
CHPL2 = 'MANUEL' 'CHPO' ('DOMA' MDOMLIM 'CENTRE') 4
'C1DX' A11 'C1DY' A12 'C2DX' B11 'C2DY' B12 ;
CACCA PIPI = 'KOPS' 'MATRIK' ;
CHPL1 = CACCA ;
GRCHP GRG = 'PENT' MDOMTOT 'FACE'
'DIAMAN2' ('MOTS' 'C1' 'C2') ('MOTS' 'C1DX' 'C1DY' 'C2DX' 'C2DY')
CHP CHPL1 CHPL2 ;
ERRO = 'MAXIMUM' (GRCHP '-' GREXA) 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
*
* TEST3 : Dirichlet BC and von Neumann BC
*
DOMLIM1 = 'MANUEL' 'POI1' ('POIN' 1 ('DOMA' MDOMLIM 'CENTRE')) ;
XXF YYF = 'COORDONNEE' DOMLIM1 ;
CHP1 = (A11 '*' XXF) '+' (A12 '*' YYF) '+' A0;
CHP2 = (B11 '*' XXF) '+' (B12 '*' YYF) '+' B0;
CHP1 = 'NOMC' 'C1' CHP1 'NATU' 'DISCRET' ;
CHP2 = 'NOMC' 'C2' CHP2 'NATU' 'DISCRET' ;
CHPL1 = CHP1 'ET' CHP2 ;
DOMLIM2 = 'DIFF' ('DOMA' MDOMLIM 'CENTRE') DOMLIM1 ;
CHPL2 = 'MANUEL' 'CHPO' DOMLIM2 4
'C1DX' A11 'C1DY' A12 'C2DX' B11 'C2DY' B12 ;
GRCHP GRG = 'PENT' MDOMTOT 'FACE'
'DIAMAN2' ('MOTS' 'C1' 'C2') ('MOTS' 'C1DX' 'C1DY' 'C2DX' 'C2DY')
CHP CHPL1 CHPL2 ;
ERRCH = GREXA '-' GRCHP ;
ERRO = 'MAXIMUM' ERRCH 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
************************************
************************************
************************************
******** EULESCAL AXI **************
************************************
************************************
************************************
************************************
'OPTION' 'ECHO' 0 'DIME' 2
'ELEM' 'TRI3'
'TRAC' 'X'
'MODE' 'AXIS' ;
GRAPH = FAUX ;
* GRAPH = VRAI ;
COMPLET = VRAI ;
*
*** MESH
*
P0 = 1.0D0 0.0D0 ;
P1 = 3.0D0 0.0D0 ;
P2 = 3.0D0 3.0D0 ;
P3 = 1.0D0 3.0D0 ;
P4 = 6.0D0 0.0D0 ;
P5 = 6.0D0 3.0D0 ;
'SI' COMPLET ;
N1 = 80 ;
N2 = 60 ;
N3 = 40 ;
N4 = 90 ;
N5 = 90 ;
N6 = 10 ;
N7 = 50 ;
'SINON' ;
N1 = 8 ;
N2 = 7 ;
N3 = 6 ;
N4 = 9 ;
N5 = 9 ;
N6 = 7 ;
N7 = 8 ;
'FINSI' ;
LINEXT1 = ((P0 'DROIT' N1 P1) 'ET' (P1 'DROIT' N2 P2) 'ET'
(P2 'DROIT' N3 P3) 'ET' (P3 'DROIT' N4 P0)) ;
LINEXT2 = ((P1 'DROIT' N2 P2) 'ET' (P2 'DROIT' N5 P5) 'ET'
(P5 'DROIT' N6 P4) 'ET' (P4 'DROIT' N7 P1)) ;
'OPTION' 'ELEM' QUA4 ;
DOM1 = 'SURFACE' LINEXT1 'PLAN' ;
'OPTION' 'ELEM' TRI3 ;
DOM2 = 'SURFACE' LINEXT2 'PLAN' ;
DOMTOT = DOM1 'ET' DOM2;
'ELIMINATION' 1D-6 DOMTOT;
LTYPE = 'MOTS' 'TRI3' 'QUA4' ;
LELEM = 'NBEL' DOMTOT LTYPE ;
'MESSAGE' ;
'MESSAGE' 'Le maillage :';
'MESSAGE' ('CHAINE' 'Nombre TRI3: ' ('EXTRAIRE' LELEM 1));
'MESSAGE' ('CHAINE' 'Nombre QUA4: ' ('EXTRAIRE' LELEM 2));
*
**** Internal and external cells
*
DOMLIM = 'CONTOUR' DOMTOT 'COULEUR' 'VERTE' ;
MAICON = DOMTOT 'ELEM' 'APPUYE' 'LARGEMENT' DOMLIM ;
DOMINT = DOMTOT 'DIFF' MAICON 'COULEUR' 'ROUGE' ;
$DOMTOT = 'MODE' DOMTOT 'EULER' ;
$DOMINT = 'MODE' DOMINT 'EULER' ;
$DOM1 = 'MODE' DOM1 'EULER' ;
TDOMTOT = 'DOMA' $DOMTOT 'VF' ;
TDOMINT = 'DOMA' $DOMINT 'VF' ;
TDOM1 = 'DOMA' $DOM1 'VF' ;
MDOMTOT = TDOMTOT . 'QUAF' ;
MDOMINT = TDOMINT . 'QUAF' ;
MDOM1 = TDOM1 . 'QUAF' ;
'ELIM' (MDOMTOT 'ET' MDOMINT 'ET' MDOM1) 1D-6;
'SI' GRAPH ;
'TRACER' MAICON 'TITRE' 'Elements sur le contour';
'TRACER' DOMINT 'TITRE' 'Elements internes' ;
'FINSI' ;
*
**** The interface points (for limit condition)
*
POIN0 = DOMLIM 'POIN' 1;
X1 Y1 = 'COORDONNEE' POIN0;
'REPETER' BLLIM (('NBNO' DOMLIM)) ;
X0 = X1 ;
Y0 = Y1 ;
POIN0 = DOMLIM 'POIN' (&BLLIM '+' 1) ;
X1 Y1 = 'COORDONNEE' POIN0 ;
XFAC = (X0 '+' X1) '/' 2 ;
YFAC = (Y0 '+' Y1) '/' 2 ;
PFAC = ('DOMA' $DOMTOT 'FACE') 'POIN' 'PROC' (XFAC YFAC);
'SI' (&BLLIM 'EGA' 1) ;
GEOLIM = 'MANUEL' 'POI1' PFAC 'BLEU' ;
'SINON' ;
GEOLIM = GEOLIM 'ET' ('MANUEL' 'POI1' PFAC 'BLEU') ;
'FINSI' ;
'FIN' BLLIM ;
*
* 'ELIMINATION' GEOLIM 1.0D-4 ($DOMTOT . 'FACE') ;
*
'SI' GRAPH ;
'TRACER' (DOMTOT 'ET' GEOLIM)
'TITRE' 'DOMTOT et points faces aux bords';
'FINSI' ;
***********************************************
********* CHP: scalar linear field *********
********* EULESCAL *********
***********************************************
*
* (A11 '*' XX) '+' (A12 '*' YY) '+' A0;
*
A11 = 12.01517 ;
A12 = 13.1421 ;
A0 = -3.21 ;
B11 = 22.01517 ;
B12 = 43.1421 ;
B0 = -5.21 ;
XX YY = 'COORDONNEE' ('DOMA' $DOMTOT 'CENTRE') ;
CHP1 = (A11 '*' XX) '+' (A12 '*' YY) '+' A0;
CHP2 = (B11 '*' XX) '+' (B12 '*' YY) '+' B0;
CHP1 = 'NOMC' 'C1' CHP1 'NATU' 'DISCRET' ;
CHP2 = 'NOMC' 'C2' CHP2 'NATU' 'DISCRET' ;
CHP = CHP1 'ET' CHP2 ;
*
* TEST3 : we impose the linear field values on the border.
* In that case error must be "zero" everywhere.
*
XLIM YLIM = 'COORDONNEE' GEOLIM ;
CHPLI1 = (A11 '*' XLIM) '+' (A12 '*' YLIM) '+' A0;
CHPLI2 = (B11 '*' XLIM) '+' (B12 '*' YLIM) '+' B0;
CHPLI1 = 'NOMC' 'C1' CHPLI1 'NATU' 'DISCRET' ;
CHPLI2 = 'NOMC' 'C2' CHPLI2 'NATU' 'DISCRET' ;
CHPLI = CHPLI1 'ET' CHPLI2 ;
GRCHPL LIMCH GRGL = 'PENT' $DOMTOT 'CENTRE'
'EULESCAL' 'NOLIMITE' ('MOTS' 'C1' 'C2') CHP 'CLIM' CHPLI ;
GRCHP1 LIMCH = 'PENT' $DOMTOT 'CENTRE'
'EULESCAL' 'NOLIMITE' ('MOTS' 'C1' 'C2')
CHP 'CLIM' CHPLI 'GRADGEO' GRGL ;
ERRO = 'MAXIMUM' (GRCHPL '-' GRCHP1) 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
GR1DX = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'P1DX' A11
'NATU' 'DISCRET' ;
GR1DY = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'P1DY' A12
'NATU' 'DISCRET';
GR2DX = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'P2DX' B11
'NATU' 'DISCRET';
GR2DY = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'P2DY' B12
'NATU' 'DISCRET';
GREXA = GR1DX 'ET' GR1DY 'ET' GR2DX 'ET' GR2DY ;
ERRCH = GREXA '-' GRCHPL ;
LMOT = 'MOTS' 'P1DX' 'P1DY' 'P2DX' 'P2DY' ;
MERRCH = ('PSCAL' ERRCH ERRCH LMOT LMOT) '**' 0.5 ;
'SI' GRAPH ;
MOD1 = 'MODELISER' ('DOMA' $DOMTOT 'MAILLAGE') 'THERMIQUE' ;
CHM_MER = 'KCHA' $DOMTOT 'CHAM' MERRCH ;
'TRACER' CHM_MER MOD1 'TITR'
('CHAINE' 'EULESCAL avec C.L. : Erreur');
* CHM_GRA = 'KCHA' $DOMTOT 'CHAM' GRCHPL ;
* 'TRACER' CHM_GRA MOD1 'TITR' ('CHAINE' 'Gradient');
'FINSI' ;
ERRO = 'MAXIMUM' ERRCH 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
************************************
************************************
**** Diamond method with BC ********
************************************
************************************
*
* 3D
*
*** MESH
*
*
'OPTION' 'DIME' 3 'ELEM' 'CUB8' 'ECHO' 0 ;
'MESSAGE' ;
'MESSAGE' '3D' ;
'MESSAGE' 'PENTE DIAMAN2' ;
'MESSAGE' ;
P1 = 0.0 0.0 0.0 ;
P2 = 3.0 0.0 0.0 ;
P3 = 3.0 3.0 0.0 ;
P4 = 0.0 3.0 0.0 ;
BAS = 'SURFACE' (P1 'DROIT' 4 P2 'DROIT' 3 P3 'DROIT' 5 P4 'DROIT' 6
P1) 'PLAN' ;
DOMTOT = BAS 'VOLUME' 'TRANSLATION' 2 (0.0 0.0 2.0) ;
CONDOM = 'ENVELOPPE' DOMTOT ;
RESTEN = CONDOM 'DIFF' BAS ;
MDOMTOT = 'MODELISER' DOMTOT 'EULER' ;
MCONDOM = 'MODELISER' CONDOM 'EULER' ;
MBAS = 'MODELISER' BAS 'EULER' ;
MRESTEN = 'MODELISER' RESTEN 'EULER' ;
TDOMTOT = 'DOMA' MDOMTOT 'VF' ;
TCONDOM = 'DOMA' MCONDOM 'VF' ;
TBAS = 'DOMA' MBAS 'VF' ;
TRESTEN = 'DOMA' MRESTEN 'VF' ;
QDOMTOT = TDOMTOT . 'QUAF' ;
QCONDOM = TCONDOM . 'QUAF' ;
QBAS = TBAS . 'QUAF' ;
QRESTEN = TRESTEN . 'QUAF' ;
'ELIMINATION' (QDOMTOT 'ET' QCONDOM 'ET' QBAS 'ET' QRESTEN) 1.0D-6 ;
'SI' GRAPH ;
'TRACER' (DOMTOT 'ET' ('DOMA' MDOMTOT 'CENTRE')
'ET' ('DOMA' MDOMTOT 'FACE')) 'TITRE' 'Centres et faces' ;
'FINSI' ;
* CHAMPOINT
A11 = 12.01517 ;
A12 = 13.1421 ;
A13 = 11.7 ;
A0 = -3.21 ;
B11 = 22.01517 ;
B12 = 43.1421 ;
B13 = 7.01 ;
B0 = -5.21 ;
XX YY ZZ = 'COORDONNEE' ('DOMA' MDOMTOT 'CENTRE') ;
CHP1 = (A11 '*' XX) '+' (A12 '*' YY) '+' (A13 '*' ZZ) '+' A0;
CHP2 = (B11 '*' XX) '+' (B12 '*' YY) '+' (B13 '*' ZZ) '+' B0;
CHP1 = 'NOMC' 'C1' CHP1 'NATU' 'DISCRET' ;
CHP2 = 'NOMC' 'C2' CHP2 'NATU' 'DISCRET' ;
CHP = CHP1 'ET' CHP2 ;
*
* TEST1 : Dirichlet BC everywhere
*
XXF YYF ZZF = 'COORDONNEE' ('DOMA' MCONDOM 'CENTRE') ;
CHP1 = (A11 '*' XXF) '+' (A12 '*' YYF) '+' (A13 '*' ZZF) '+' A0;
CHP2 = (B11 '*' XXF) '+' (B12 '*' YYF) '+' (B13 '*' ZZF) '+' B0;
CHP1 = 'NOMC' 'C1' CHP1 'NATU' 'DISCRET' ;
CHP2 = 'NOMC' 'C2' CHP2 'NATU' 'DISCRET' ;
CHPL1 = CHP1 'ET' CHP2 ;
CACCA PIPI = 'KOPS' 'MATRIK' ;
CHPL2 = 'COPIER' CACCA ;
*
GRCHP GRG = 'PENT' MDOMTOT 'FACE'
'DIAMAN2' ('MOTS' 'C1' 'C2') ('MOTS' 'C1DX' 'C1DY' 'C1DZ'
'C2DX' 'C2DY' 'C2DZ') CHP CHPL1 CHPL2 ;
GRCHP1 = 'PENT' MDOMTOT 'FACE'
'DIAMAN2' ('MOTS' 'C1' 'C2') ('MOTS' 'C1DX' 'C1DY' 'C1DZ'
'C2DX' 'C2DY' 'C2DZ') CHP CHPL1 CHPL2 'GRADGEO' GRG ;
ERRO = 'MAXIMUM' (GRCHP '-' GRCHP1) 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
GREXA = 'MANUEL' 'CHPO' ('DOMA' MDOMTOT 'FACE') 6 'C1DX' A11
'C1DY' A12 'C1DZ' A13 'C2DX' B11 'C2DY' B12 'C2DZ' B13
'NATU' 'DISCRET' ;
ERRCH = GREXA '-' GRCHP ;
ERRO = 'MAXIMUM' ERRCH 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
*
* TEST2 : von Neumann BC everywhere
*
CHPL2 = 'MANUEL' 'CHPO' ('DOMA' MCONDOM 'CENTRE') 6
'C1DX' A11 'C1DY' A12 'C1DZ' A13 'C2DX' B11 'C2DY'
B12 'C2DZ' B13 'NATU' 'DISCRET' ;
CACCA PIPI = 'KOPS' 'MATRIK' ;
CHPL1 = CACCA ;
GRCHP GRG = 'PENT' MDOMTOT 'FACE'
'DIAMAN2' ('MOTS' 'C1' 'C2') ('MOTS' 'C1DX' 'C1DY' 'C1DZ'
'C2DX' 'C2DY' 'C2DZ') CHP CHPL1 CHPL2 ;
ERRO = 'MAXIMUM' (GRCHP '-' GREXA) 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
*
* TEST2 : von Neumann and Dirichlet BC
*
XXF YYF ZZF = 'COORDONNEE' ('DOMA' MBAS 'CENTRE') ;
CHP1 = (A11 '*' XXF) '+' (A12 '*' YYF) '+' (A13 '*' ZZF) '+' A0;
CHP2 = (B11 '*' XXF) '+' (B12 '*' YYF) '+' (B13 '*' ZZF) '+' B0;
CHP1 = 'NOMC' 'C1' CHP1 'NATU' 'DISCRET' ;
CHP2 = 'NOMC' 'C2' CHP2 'NATU' 'DISCRET' ;
CHPL1 = CHP1 'ET' CHP2 ;
CHPL2 = 'MANUEL' 'CHPO' ('DOMA' MRESTEN 'CENTRE') 6
'C1DX' A11 'C1DY' A12 'C1DZ' A13 'C2DX' B11 'C2DY'
B12 'C2DZ' B13 'NATU' 'DISCRET' ;
GRCHP GRG = 'PENT' MDOMTOT 'FACE'
'DIAMAN2' ('MOTS' 'C1' 'C2') ('MOTS' 'C1DX' 'C1DY' 'C1DZ'
'C2DX' 'C2DY' 'C2DZ') CHP CHPL1 CHPL2 ;
ERRO = 'MAXIMUM' (GRCHP '-' GREXA) 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
*
*
*****************************
*****************************
***** Maillage regulier *****
*****************************
*****************************
*
* Géométrie
*
pA = 0. 0. 0. ;
pB = 1. 0. 0. ;
pC = 1. 1. 0. ;
pD = 0. 1. 0. ;
pE = 0. 0. 1. ;
pF = 1. 0. 1. ;
pG = 1. 1. 1. ;
pH = 0. 1. 1. ;
*
* Paramètres de la discrétisation de base
*
nAB = 1 ;
nBC = 1 ;
nCD = 1 ;
nDA = 1 ;
nH = 1 ;
*
* Géométrie discrétisée
*
bas = 'DROIT' nAB pA pB ;
droite = 'DROIT' nBC pB pC ;
haut = 'DROIT' nCD pC pD ;
gauche = 'DROIT' nDA pD pA ;
spourt = bas 'ET' droite 'ET' haut 'ET' gauche ;
sbas = 'SURFACE' spourt 'PLAN' ;
mt = 'VOLUME' sbas 'TRAN' nH pE ;
*
cnt = 'ENVE' mt 'COULEUR' 'VERT' ;
'SI' (graph ET FAUX) ;
'TRACER' (mt 'ET' cnt) 'TITRE' 'Enveloppe' ;
'FINSI' ;
*
* Creation des modèles
*
* MODELS
$mt = 'MODE' mt 'EULER' ;
$cnt = 'MODE' cnt 'EULER' ;
tmt = 'DOMA' $mt 'VF' ;
tcnt = 'DOMA' $cnt 'VF' ;
* QUAF
Mmt = tmt . 'QUAF' ;
Mcnt = tcnt . 'QUAF' ;
'ELIM' (Mmt 'ET' Mcnt) 1.E-8 ;
* MAILLAGES
'SI' graph ;
'TRACER' (mt 'ET' (('DOMA' $mt 'FACE') 'COULEUR' 'ROUG') 'ET'
(('DOMA' $MT 'CENTRE') 'COUL' 'VERT'))
'TITRE' 'Centres et faces' ;
'FINSI' ;
*
* solex = (alpha '*' x) '+' (beta '*' y) '+' (gamma '*' z) '+'
* (eta '*' x '*' y '*' z) '+' delta
*
alpha = (** 2. 0.5) ;
beta = (** 3. 0.5) ;
gamma = (** 5. 0.5) ;
* eta = (** 7. 0.5) ;
eta = 0. ;
delta = (** 9. 0.5) ;
*
* Solution exacte bilinéaire (sur le centre du contour)
*
xx yy zz = 'COORDONNEE' ('DOMA' $cnt 'CENTRE') ;
solex = (xx '*' alpha) '+' (yy '*' beta) '+'
(zz '*' gamma) '+'
(xx '*' yy '*' zz '*' eta) '+' delta ;
*
* Solution exacte aux centres
*
xxc yyc zzc = 'COORDONNEE' ('DOMA' $mt 'CENTRE') ;
solexc = (xxc '*' alpha) '+' (yyc '*' beta) '+'
(zzc '*' gamma) '+'
(xxc '*' yyc '*' zzc '*' eta) '+' delta ;
*
* Gradient exact (aux faces)
*
xxf yyf zzf = 'COORDONNEE' ('DOMA' $mt 'FACE') ;
gradtx = 'NOMC' 'UX' (alpha '+' (yyf '*' zzf '*' eta)) 'NATU'
'DISCRET' ;
gradty = 'NOMC' 'UY' (beta '+' (xxf '*' zzf '*' eta)) 'NATU'
'DISCRET' ;
gradtz = 'NOMC' 'UZ' (gamma '+' (xxf '*' yyf '*' eta)) 'NATU'
'DISCRET' ;
gradex = gradtx 'ET' gradty 'ET' gradtz ;
*
* BC
*
*
* TEST1 : Dirichlet BC everywhere
*
cdir = ('DOMA' $cnt 'CENTRE') ;
* cvn = set of FACE centers where we impose von Neumann boundary
* conditions
cvn = ('DIFF' ('DOMA' $cnt 'CENTRE') cdir) 'COUL' 'TURQ' ;
'SI' graph ;
'TRACER' (cnt et cvn) 'TITRE' 'Von Neumann B.C.' ;
'TRACER' (cnt et cdir) 'TITRE' 'Dirichlet B.C.' ;
'FINSI' ;
tlim1 = 'REDU' solex cdir ;
gradtl = 'REDU' (gradtx 'ET' gradty 'ET' gradtz) cvn ;
tlim2 = 'NOMC' ('MOTS' 'UX' 'UY' 'UZ') ('MOTS' 'UX' 'UY' 'UZ')
gradtl ;
gradt mchamt = 'PENT' $mt 'FACE' 'DIAMAN2' ('MOTS' 'SCAL')
('MOTS' 'UX' 'UY' 'UZ') solexc tlim1 tlim2 ;
ERRO = 'MAXIMUM' (gradt '-' gradex) 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
*
* TEST2 : von Neumann BC everywhere
*
cdir = 'DIFF' ('DOMA' $cnt 'CENTRE') ('DOMA' $cnt 'CENTRE') ;
* cvn = set of FACE centers where we impose von Neumann boundary
* conditions
cvn = ('DIFF' ('DOMA' $cnt 'CENTRE') cdir) 'COUL' 'TURQ' ;
'SI' graph ;
'TRACER' (cnt et cvn) 'TITRE' 'Von Neumann B.C.' ;
'TRACER' (cnt et cdir) 'TITRE' 'Dirichlet B.C.' ;
'FINSI' ;
tlim1 = 'REDU' solex cdir ;
gradtl = 'REDU' (gradtx 'ET' gradty 'ET' gradtz) cvn ;
tlim2 = 'NOMC' ('MOTS' 'UX' 'UY' 'UZ') ('MOTS' 'UX' 'UY' 'UZ')
gradtl ;
gradt mchamt = 'PENT' $mt 'FACE' 'DIAMAN2' ('MOTS' 'SCAL')
('MOTS' 'UX' 'UY' 'UZ') solexc tlim1 tlim2 ;
ERRO = 'MAXIMUM' (gradt '-' gradex) 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
*
* TEST3 : von Neumann and Dirichlet bc
*
cdir = 'MANUEL' 'POI1' ('POIN' 1 ('DOMA' $cnt 'CENTRE')) ;
* cvn = set of FACE centers where we impose von Neumann boundary
* conditions
cvn = ('DIFF' ('DOMA' $cnt 'CENTRE') cdir) 'COUL' 'TURQ' ;
'SI' graph ;
'TRACER' (cnt et cvn) 'TITRE' 'Von Neumann B.C.' ;
'TRACER' (cnt et cdir) 'TITRE' 'Dirichlet B.C.' ;
'FINSI' ;
tlim1 = 'REDU' solex cdir ;
gradtl = 'REDU' (gradtx 'ET' gradty 'ET' gradtz) cvn ;
tlim2 = 'NOMC' ('MOTS' 'UX' 'UY' 'UZ') ('MOTS' 'UX' 'UY' 'UZ')
gradtl ;
gradt mchamt = 'PENT' $mt 'FACE' 'DIAMAN2' ('MOTS' 'SCAL')
('MOTS' 'UX' 'UY' 'UZ') solexc tlim1 tlim2 ;
ERRO = 'MAXIMUM' (gradt '-' gradex) 'ABS' ;
'SI' (ERRO > 1.0D-8) ;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'Probleme dans PENT');
'MESSAGE' ;
'ERREUR' 5;
'FINSI' ;
'FIN' ;
© Cast3M 2003 - All rights reserved.
Disclaimer