Télécharger la plaquette Cast3M


fr en

Télécharger pent3D3.dgibi

Retour à la liste

Numérotation des lignes : 

  1. * fichier :  pent3D3.dgibi
  2. ************************************************************
  3. ************************************************************
  4. **** APPROCHE VF "Cell-Centred Formulation".            ****
  5. **** OPÉRATEUR PENT, pour le calcul des gradients et    ****
  6. **** des limiteurs                                      ****
  7. **** Cas test: calcul du gradient en 3D avec condition  ****
  8. ****           de 'TYPE' mur                            ****
  9. ****                                                    ****
  10. **** A. BECCANTINI, TTMF    MAI     1998                ****
  11. ****                                                    ****
  12. **** A. BECCANTINI, LTMF    MARS    2000                ****
  13. ****                Parametres d'erreurs eccessives pour****
  14. ****                les machines SUN 'ET' HP            ****
  15. ************************************************************
  16. ************************************************************
  17.  
  18.  
  19. 'OPTION'  'ECHO' 1 ;
  20. 'OPTION'  'DIME' 3 ;
  21. 'OPTION'  'ELEM' 'CUB8' ;
  22. 'OPTION'  'TRAC' 'X' ;
  23. GRAPH=FAUX;
  24.  
  25. *
  26. *** MAILLAGE
  27. *
  28.  
  29. A0 = 0.0D0 0.0D0 0.0D0;
  30. A1 = 3.0D0 0.0D0 0.0D0;
  31. A2 = 3.0D0 3.0D0 0.0D0;
  32. A3 = 0.0D0 3.0D0 0.0D0;
  33. B0 = 0.0D0 0.0D0 3.0D0;
  34. B1 = 3.0D0 0.0D0 3.0D0;
  35. B2 = 3.0D0 3.0D0 3.0D0;
  36. B3 = 0.0D0 3.0D0 3.0D0;
  37. N=7;
  38. N1=7;
  39. LIN1 = A0 DROIT N1 A1;
  40. LIN2 = A1 DROIT N1 A2;
  41. LIN3 = A2 DROIT N1 A3;
  42. LIN4 = A3 DROIT N1 A0;
  43. LI1 = B0 DROIT N B1;
  44. LI2 = B1 DROIT N B2;
  45. LI3 = B2 DROIT N B3;
  46. LI4 = B3 DROIT N B0;
  47. L1 = A0 DROIT N B0;
  48. L2 = A1 DROIT N B1;
  49. L3 = A2 DROIT N B2;
  50. L4 = A3 DROIT N B3;
  51. S1='DALL' LIN1 LIN2 LIN3 LIN4 'PLAN';
  52. S2='DALL' LIN1 L2 LI1 L1 'PLAN';
  53. S3='DALL' LIN2 L3 LI2 L2 'PLAN';
  54. S4='DALL' LIN3 L4 LI3 L3 'PLAN';
  55. S5='DALL' LIN4 L1 LI4 L4 'PLAN';
  56. S6='DALL' LI1 LI2 LI3 LI4 'PLAN';
  57. DOMTOT=S1 'VOLU' N S6;
  58. *DOMTOT=(S1 et S2 et S3 et S4 et S5 et S6) 'VOLU';
  59. 'ELIM' 1.0D-12 DOMTOT  (S1 et S2 et S3 et S4 et S5 et S6);
  60. *
  61. **** Le domaine interne
  62. *
  63.  
  64.  GCTOT = DOMTOT 'ELEM'  'APPUYE' 'LARGEMENT'
  65.   (S1 et S2 et S3 et S4 et S5 et S6);
  66.  DOMINT = DOMTOT 'DIFF' GCTOT ;
  67.  
  68.  $DOMTOT = 'MODE' DOMTOT 'EULER'  ;
  69.  $DOMINT = 'MODE' DOMINT 'EULER'  ;
  70.  $GCTOT  = 'MODE' GCTOT  'EULER'  ;
  71.  
  72.  TDOMTOT = 'DOMA' $DOMTOT 'VF'  ;
  73.  TDOMINT = 'DOMA' $DOMINT 'VF'  ;
  74.  TGCTOT  = 'DOMA' $GCTOT  'VF'  ;
  75.  
  76.  MDOMTOT = TDOMTOT . 'QUAF';
  77.  MDOMINT = TDOMINT . 'QUAF';
  78.  MGCTOT  = TGCTOT  . 'QUAF';
  79. 'ELIM' (MDOMTOT 'ET' MDOMINT 'ET' MGCTOT) 1D-6;
  80.  
  81.  $DOMTOT1 = 'MODELISER' MDOMTOT 'NAVIER_STOKES' 'LINE' ;
  82.  $DOMINT1 = 'MODELISER' MDOMINT 'NAVIER_STOKES' 'LINE' ;
  83.  $GCTOT1  = 'MODELISER' MGCTOT  'NAVIER_STOKES' 'LINE' ;
  84.  
  85.  TDOMTOT1 = 'DOMA' $DOMTOT1 'TABLE'  ;
  86.  TDOMINT1 = 'DOMA' $DOMINT1 'TABLE'  ;
  87.  TGCTOT1  = 'DOMA' $GCTOT1  'TABLE'  ;
  88.  
  89. *
  90. *** Calcul du gradient pour un champ lineaire
  91. *
  92.  
  93.  COEF1X = 2.01517  ;
  94.  COEF1Y = 3.1421   ;
  95.  COEF1Z = 1.5;
  96.  COEF2X = -2.7  ;
  97.  COEF2Y = -3.21   ;
  98.  COEF2Z = -0.56;
  99.  XX YY ZZ= 'COORDONNEE' ('DOMA' $DOMTOT 'CENTRE') ;
  100.  CHP1 = (COEF1X '*' XX) '+' (COEF1Y '*' YY) '+' (COEF1Z '*' ZZ) ;
  101.  CHP1 = 'NOMC' 'UX' CHP1 'NATU' 'DISCRET' ;
  102.  CHP2 = (COEF2X '*' XX) '+' (COEF2Y '*' YY) '+' (COEF2Z '*' ZZ);
  103.  CHP2 = 'NOMC' 'UY' CHP2 'NATU' 'DISCRET' ;
  104.  CHP3 = 'BRUIT' 'BLAN' 'UNIF' 0.0D0 1.0D0 ('DOMA' $DOMTOT 'CENTRE');
  105. * CHP3 = (COEF2X '*' XX) '+' (COEF2Y '*' YY) '+' (COEF2Z '*' ZZ);
  106.  CHP3 = 'NOMC' 'UZ' CHP3 'NATU' 'DISCRET' ;
  107.  CHP = CHP1 'ET' CHP2 'ET' CHP3 ;
  108. *
  109. **** Connectivite entre les symetriques
  110. *
  111. 'REPETER' BLLIM ('NBEL' ('ENVE' DOMINT)) ;
  112.     ELE0 = ('ENVE' DOMINT) 'ELEM' &BLLIM;
  113.     ELE0 = ELE0 'CHAN' 'POI1';
  114.     P0 = ELE0 'POINT' 1;
  115.     X0 Y0 Z0= 'COORDONNEE' P0 ;
  116.     P1 = ELE0 'POINT' 2;
  117.     X1 Y1 Z1= 'COORDONNEE' P1 ;
  118.      P2 = ELE0 'POINT' 3;
  119.     X2 Y2 Z2= 'COORDONNEE' P2 ;
  120.     P3 = ELE0 'POINT' 4;
  121.     X3 Y3 Z3= 'COORDONNEE' P3 ;
  122.     XFAC = (X0 '+' X1 '+' X2 '+' X3) '/' 4 ;
  123.     YFAC = (Y0 '+' Y1 '+' Y2 '+' Y3) '/' 4 ;
  124.     ZFAC = (Z0 '+' Z1 '+' Z2 '+' Z3) '/' 4 ;
  125.     PFAC = ('DOMA' $DOMTOT 'FACE') 'POIN' 'PROC' (XFAC YFAC ZFAC);
  126.     GEOFAC1 = ('DOMA' $DOMINT 'FACEL') 'ELEM' 'APPUYE'
  127.              'LARGEMENT' PFAC ;
  128.     GEOFAC2 = ('DOMA' $GCTOT 'FACEL') 'ELEM' 'APPUYE'
  129.              'LARGEMENT' PFAC ;
  130. *
  131. **** Tranformation en POI1
  132. *
  133.     GEO1POI1 = 'CHANGER' 'POI1' GEOFAC1 ;
  134.     PCEL11 = GEO1POI1 'POIN'  1 ;
  135.     PCEL12 = GEO1POI1 'POIN'  2 ;
  136.     GEO2POI1 = 'CHANGER' 'POI1' GEOFAC2 ;
  137.     PCEL21 = GEO2POI1 'POIN' 1 ;
  138.     PCEL22 = GEO2POI1 'POIN' 2 ;
  139. *
  140. **** Il faur verifier que PFAC = PCEL12 = PCEL22
  141. *    ('NBEL' GEO1POI1) = ('NBEL' GEO2POI1) = 2
  142. *
  143.    'SI'  ( ('NBEL' GEO1POI1) 'NEG' 2);
  144.        'MESSAGE' ;
  145.        'MESSAGE'
  146.           'Probleme dans la creation du domaine entree subsonique';
  147.        'MESSAGE' ;
  148.        'ERREUR' 21 ;
  149.    'FINSI' ;
  150.    'SI'  ( ('NBEL' GEO2POI1) 'NEG' 2);
  151.        'MESSAGE' ;
  152.        'MESSAGE' 
  153.           'Probleme dans la creation du domaine entree subsonique';
  154.        'MESSAGE' ;
  155.        'ERREUR' 21 ;
  156.    'FINSI' ;
  157.    'SI'  ( PCEL12 'NEG' PFAC);
  158.        'MESSAGE' ;
  159.        'MESSAGE' 
  160.           'Probleme dans la creation du domaine entree subsonique';
  161.        'MESSAGE' ;
  162.        'ERREUR' 21 ;
  163.    'FINSI' ;
  164.    'SI'  ( PCEL22 'NEG' PFAC);
  165.        'MESSAGE' ;
  166.        'MESSAGE' 
  167.           'Probleme dans la creation du domaine entree subsonique';
  168.        'MESSAGE' ;
  169.        'ERREUR' 21 ;
  170.    'FINSI' ;
  171. *
  172. *** Creation d'un maillage SEG2
  173. *
  174.    'SI' (&BLLIM 'EGA' 1);
  175.        CONBOR = 'MANUEL' 'SEG2' PCEL11 PCEL21 'COULEUR'  'BLEU' ;
  176.        BORD = 'MANU' 'POI1' PFAC;
  177.        CB = PCEL21;
  178.     'SINON' ;
  179.        CONBOR = CONBOR 'ET'
  180.           ( 'MANUEL' 'SEG2' PCEL11 PCEL21 'COULEUR'  'BLEU' );
  181.        BORD = BORD 'ET' ('MANU' 'POI1' PFAC);
  182.        CB = CB 'ET' PCEL21;
  183.     'FINSI' ;   
  184. 'FIN' BLLIM ;
  185. *
  186. **** On doit imposer une condition de bord
  187. *
  188. CHPCONT='KPRO' CHP CONBOR; 
  189. CHP = 'KCHT' $DOMTOT 'VECT' 'CENTRE' CHP CHPCONT;
  190. * 
  191. *
  192. **** N.B. : CHP1 n'as pas comme support geometrique
  193. *           ($DOMTOT 'CENTRE')
  194. *
  195.     GRCHP0 ALCHP0 COEF = 'PENT' $DOMTOT 'CENTRE' 
  196.          'EULESCAL' 'LIMITEUR'
  197.          ('MOTS' 'UX' 'UY' 'UZ') CHP;
  198.  
  199.     GRCHP01 ALCHP01 COEF1 = 'PENT' $DOMTOT1 'CENTRE' 
  200.          'EULESCAL' 'LIMITEUR'
  201.          ('MOTS' 'UX' 'UY' 'UZ') CHP;
  202.  
  203.     GRCHP1 ALCHP1 COEF = 'PENT' $DOMINT 'CENTRE' 'EULESCAL'
  204.     'LIMITEUR' ('MOTS' 'UX' 'UY' 'UZ')
  205.     ('REDU' CHP ('DOMA' $DOMINT 'CENTRE'));
  206.  
  207. GRCHP = 'REDU' GRCHP0 ('DOMA' $DOMINT 'CENTRE');
  208. ALCHP = 'REDU' ALCHP0 ('DOMA' $DOMINT 'CENTRE');
  209. *
  210. **** Le gradient calculé sur le domaine interne est le
  211. **** même que celui calculé avec les symétriques au bord
  212. *
  213.  
  214.  
  215.  ERRO = 'MAXIMUM' (GRCHP1 '-' GRCHP) 'ABS' ;
  216.  'SI' (ERRO  > 5.D-6);
  217.     'MESSAGE' 'Probleme dans Grad';
  218.      'MESSAGE' ERRO;
  219.     'ERREUR' 5 ;
  220.  'FINSI';
  221.  
  222.  
  223.  ERRO = 'MAXIMUM' (ALCHP1 '-' ALCHP) 'ABS';
  224.  
  225. 'SI' (ERRO > 5.D-6);
  226.      'MESSAGE' ;
  227.      'MESSAGE' ('CHAINE' 'Probleme dans Lim');
  228.      'MESSAGE' ERRO;
  229.      'ERREUR' 5;
  230. 'FINSI' ;
  231.  
  232.  
  233. 'SI' GRAPH;
  234. UNCH1=VECT GRCHP1 0.1 'P1DX' 'P1DY' 'P1DZ' ROUGE;
  235. UNCH2=VECT GRCHP1 0.1 'P2DX' 'P2DY' 'P2DZ' VERT;
  236. UNCH3=VECT GRCHP1 0.1 'P3DX' 'P3DY' 'P3DZ' VERT;
  237. UNCH4=VECT GRCHP  0.1 'P1DX' 'P1DY' 'P1DZ' BLEU;
  238. UNCH5=VECT GRCHP  0.1 'P2DX' 'P2DY' 'P2DZ' BLANC;
  239. UNCH6=VECT GRCHP  0.1 'P3DX' 'P3DY' 'P3DZ' BLANC;
  240. TRAC (UNCH1 'ET' UNCH4) DOMINT;
  241. TRAC (UNCH2 'ET' UNCH5) DOMINT;
  242. TRAC (UNCH3 'ET' UNCH6) DOMINT;
  243. 'FINSI';
  244.  
  245.  
  246.  
  247. *
  248. **** On impose le vecteur symetrique au bord
  249. *
  250. NRME = 'REDU' ('DOMA' $DOMTOT 'XXNORMAF') BORD;
  251. GNOLD = CHP;
  252. GNPE = 'KPRO' CHP CONBOR;
  253.  
  254. NL = 'NBEL' BORD;
  255. 'REPETER' BL NL;
  256.   NX = 'EXTR' NRME 'UX' ((BORD 'ELEM' &BL) 'POIN' 1);
  257.   NY = 'EXTR' NRME 'UY' ((BORD 'ELEM' &BL) 'POIN' 1);
  258.   NZ = 'EXTR' NRME 'UZ' ((BORD 'ELEM' &BL) 'POIN' 1);
  259.   GNXOLD = 'EXTR' GNPE 'UX' ((CB 'ELEM' &BL) 'POIN' 1);
  260.   GNYOLD = 'EXTR' GNPE 'UY' ((CB 'ELEM' &BL) 'POIN' 1);
  261.   GNZOLD = 'EXTR' GNPE 'UZ' ((CB 'ELEM' &BL) 'POIN' 1);
  262.   GNN = (NX*GNXOLD) '+' (NY*GNYOLD) '+' (NZ*GNZOLD);
  263.   GNX = GNXOLD '-' (2 '*' NX '*' GNN);
  264.   GNY = GNYOLD '-' (2 '*' NY '*' GNN);
  265.   GNZ = GNZOLD '-' (2 '*' NZ '*' GNN);
  266.   GNNEW= 'MANU' 'CHPO' (CB 'ELEM' &BL) 3 
  267.          'UX' GNX 'UY' GNY 'UZ' GNZ;
  268.   CHP= 'KCHT' $DOMTOT 'VECT' 'CENTRE' GNOLD GNNEW;
  269.   GNOLD = CHP;
  270. 'FIN' BL;
  271.  
  272. *
  273. **** N.B. : CHP1 n'as pas comme support geometrique
  274. *           ($DOMTOT 'CENTRE')
  275. *
  276.     GRCHP0 ALCHP0 COEF = 'PENT' $DOMTOT 'CENTRE' 'EULEVECT' 
  277.                     'LIMITEUR'
  278.                     ('MOTS' 'UX' 'UY' 'UZ') CHP;
  279.  
  280.     GRCHP1 ALCHP1 COEF = 'PENT' $DOMINT 'CENTRE' 'EULEVECT'
  281.                    'LIMITEUR' ('MOTS' 'UX' 'UY' 'UZ')
  282.                    ('REDU' CHP ('DOMA' $DOMINT 'CENTRE'));
  283. GRCHP = 'REDU' GRCHP0 ('DOMA' $DOMINT 'CENTRE');
  284. ALCHP = 'REDU' ALCHP0 ('DOMA' $DOMINT 'CENTRE');
  285.  
  286. *
  287. **** Le bord est bien calcule
  288. *
  289.  
  290.  ERRO = 'MAXIMUM' (GRCHP '-' GRCHP1) 'ABS' ;
  291.  'SI' (ERRO  > 5D-6);
  292.     'MESSAGE' 'Probleme dans Grad';
  293.      'MESSAGE' erro;
  294.     'ERREUR' 5 ;
  295.  'FINSI';
  296.  
  297.  ERRO = 'MAXIMUM' (ALCHP '-' ALCHP1) 'ABS' ;
  298.  
  299. 'SI' (ERRO > 5D-6);
  300.      'MESSAGE' ;
  301.      'MESSAGE' ('CHAINE' 'Probleme dans Lim');
  302.      'MESSAGE' erro;
  303.      'ERREUR' 5;
  304. 'FINSI' ;
  305. 'SI' GRAPH;
  306. UNCH1=VECT GRCHP1 0.1 'P1DX' 'P1DY' 'P1DZ' ROUGE;
  307. UNCH2=VECT GRCHP1 0.1 'P2DX' 'P2DY' 'P2DZ' VERT;
  308. UNCH3=VECT GRCHP1 0.1 'P3DX' 'P3DY' 'P3DZ' VERT;
  309. UNCH4=VECT GRCHP  0.1 'P1DX' 'P1DY' 'P1DZ' BLEU;
  310. UNCH5=VECT GRCHP  0.1 'P2DX' 'P2DY' 'P2DZ' BLANC;
  311. UNCH6=VECT GRCHP  0.1 'P3DX' 'P3DY' 'P3DZ' BLANC;
  312. TRAC (UNCH1 'ET' UNCH4) DOMINT;
  313. TRAC (UNCH2 'ET' UNCH5) DOMINT;
  314. TRAC (UNCH3 'ET' UNCH6) DOMINT;
  315. 'FINSI';
  316.  
  317. 'FIN' ;
  318.  
  319.  
  320.  
  321.  
  322.  
  323.  
  324.  
  325.  
  326.  

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