Télécharger pret1.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier : pret1.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. **********************************************************************
  5. **** APPROCHE VF "Cell-Centred Formulation" pour la solution des ****
  6. **** Equations d'Euler pour un gaz parfait. ****
  7. **** ****
  8. **** OPERATEUR PRET ****
  9. **** Operateur qui 'recontruit les variables primitives aux faces ****
  10. **** Cas gaz monoespece "thermally perfect" ****
  11. **** Differents cas tests ****
  12. **** ****
  13. **** A. BECCANTINI DRN/DMT/SEMT/TTMF AOUT 1998 ****
  14. **********************************************************************
  15.  
  16. 'OPTION' 'DIME' 2 ;
  17. 'OPTION' 'ELEM' QUA4 ;
  18. 'OPTION' 'ECHO' 1 ;
  19. 'OPTION' 'TRAC' 'PS' ;
  20.  
  21. *
  22. *** GRAPH
  23. *
  24.  
  25. * GRAPH = VRAI ;
  26. GRAPH = FAUX ;
  27.  
  28. ***************************
  29. *** PROCEDURE RUO1 ****
  30. ***************************
  31. *
  32. *** Precedure pour la rotation d'un tenseur de premier
  33. * ordre
  34. *
  35.  
  36. 'DEBPROC' RUO1 ;
  37. 'ARGUMENT' ALPHA*'FLOTTANT' UN*'FLOTTANT' UT*'FLOTTANT';
  38.  
  39. SINA = 'SIN' ALPHA ;
  40. COSA = 'COS' ALPHA ;
  41.  
  42. UX = (UN * COSA ) '-' (UT * SINA) ;
  43. UY = (UN * SINA ) '+' (UT * COSA) ;
  44.  
  45.  
  46. 'FINPROC' UX UY;
  47.  
  48. ***************************
  49. *** PROCEDURE RUO2 ****
  50. ***************************
  51. *
  52. *** Precedure pour la rotation d'un tenseur de deuxieme
  53. * ordre
  54. *
  55.  
  56. 'DEBPROC' RUO2 ;
  57. 'ARGUMENT' ALPHA*'FLOTTANT' UNN*'FLOTTANT' UNT*'FLOTTANT'
  58. UTN*'FLOTTANT' UTT*'FLOTTANT';
  59. *
  60. **** (n,t) -> (x,y)
  61. *
  62. * n = CA x '+' SA y ;
  63. * t = -SA x '+' CA y ;
  64. *
  65. * UNT = DUN/DT
  66. *
  67.  
  68. SA = 'SIN' ALPHA ;
  69. CA = 'COS' ALPHA ;
  70. CA2 = CA * CA ;
  71. CASA = CA * SA ;
  72. SA2 = SA * SA;
  73.  
  74. UXX = (CA2 * UNN) '-' (CASA * (UNT '+' UTN)) '+' (SA2 *UTT) ;
  75. UYX = (CASA * (UNN '-' UTT)) '-' (SA2 * UNT) '+' (CA2 * UTN ) ;
  76. UXY = (CASA * (UNN '-' UTT)) '+' (CA2 * UNT) '-' (SA2 * UTN ) ;
  77. UYY = (SA2 * UNN) '+' (CASA * (UNT '+' UTN)) '+' (CA2 *UTT) ;
  78.  
  79. 'FINPROC' UXX UXY UYX UYY;
  80.  
  81.  
  82. ***************************
  83. ***** DOMAINE SPATIAL ****
  84. ***************************
  85.  
  86. *
  87. **** Deux carre
  88. *
  89.  
  90. A1 = 0.0D0 0.0D0;
  91. A2 = 1.0D0 0.0D0;
  92. A3 = 2.0D0 0.0D0;
  93. A4 = 2.0D0 1.0D0;
  94. A5 = 1.0D0 1.0D0;
  95. A6 = 0.0D0 1.0D0;
  96.  
  97. L12 = A1 'DROIT' 1 A2;
  98. L23 = A2 'DROIT' 1 A3;
  99. L34 = A3 'DROIT' 1 A4;
  100. L45 = A4 'DROIT' 1 A5;
  101. L56 = A5 'DROIT' 1 A6;
  102. L61 = A6 'DROIT' 1 A1;
  103. L25 = A2 'DROIT' 1 A5;
  104.  
  105.  
  106. DOM10 = 'DALL' L12 L25 L56 L61
  107. 'PLANE';
  108. DOM20 = 'DALL' L23 L34 L45 ('INVERSE' L25)
  109. 'PLANE';
  110.  
  111. *
  112. *** Point face entre le deux carre, ou on fait les controlles
  113. *
  114.  
  115. P10 = 1.0 0.5;
  116.  
  117. DOM1 = DOM10 ;
  118. DOM2 = DOM20 ;
  119.  
  120.  
  121. 'ELIMINATION' (DOM1 ET DOM2) 1D-6;
  122. DOMTOT = DOM1 ET DOM2;
  123.  
  124. $DOMTOT = 'MODELISER' DOMTOT 'EULER';
  125.  
  126. $DOM1 = 'MODELISER' DOM1 'EULER';
  127. $DOM2 = 'MODELISER' DOM2 'EULER';
  128.  
  129. TDOMTOT = 'DOMA' $DOMTOT 'VF';
  130.  
  131. TDOM1 = 'DOMA' $DOM1 'VF';
  132. TDOM2 = 'DOMA' $DOM2 'VF';
  133.  
  134. MDOM1 = TDOM1 . 'QUAF' ;
  135. MDOM2 = TDOM2 . 'QUAF' ;
  136. MDOMTOT = TDOMTOT . 'QUAF' ;
  137.  
  138. 'ELIMINATION' (MDOMTOT ET MDOM1) 0.0001 ;
  139. 'ELIMINATION' (MDOMTOT ET MDOM2) 0.0001 ;
  140.  
  141. FACTOT = 'DOMA' $DOMTOT 'FACEL' ;
  142.  
  143.  
  144. 'SI' GRAPH;
  145. 'TRACER' (DOMTOT 'ET' FACTOT
  146. 'ET' P10) 'TITRE' 'Domaine et FACEL';
  147. 'FINSI' ;
  148.  
  149.  
  150. *
  151. **** TEST1: si les limiteurs sont nuls ou les gradients son nuls
  152. * le deuxieme ordre en espace 'et premier ordre
  153. * en temps degenere en premier ordre en espace
  154. *
  155.  
  156.  
  157. RN = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' 1.0D0;
  158.  
  159. PN = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' 2.0D0;
  160.  
  161. VIT = 'KCHT' $DOMTOT 'VECT' 'CENTRE' (4.0 5.0);
  162.  
  163. GAMMA = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' 1.5D0;
  164.  
  165.  
  166. *
  167. **** Limiteurs nuls, gradients non-nuls
  168. *
  169.  
  170.  
  171.  
  172. GRADR = 'KCHT' $DOMTOT 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  173. (1.0D0 1.0D0);
  174.  
  175.  
  176.  
  177. GRADP = 'KCHT' $DOMTOT 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  178. (1.0D0 1.0D0);
  179.  
  180.  
  181.  
  182. GRADVX = 'KCHT' $DOMTOT 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  183. (1.0D0 1.0D0);
  184.  
  185.  
  186. GRADVY = 'KCHT' $DOMTOT 'VECT' 'CENTRE' 'COMP' 'P2DX' 'P2DY'
  187. (1.0D0 1.0D0);
  188.  
  189. GRADV = GRADVX 'ET' GRADVY ;
  190.  
  191.  
  192.  
  193. ALR = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' 'COMP' 'P1' 0.0D0 ;
  194.  
  195.  
  196. ALP = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' 'COMP' 'P1' 0.0D0 ;
  197.  
  198.  
  199. ALV = 'KCHT' $DOMTOT 'VECT' 'CENTRE' 'COMP' 'P1' 'P2'
  200. (0.0 0.0D0) ;
  201.  
  202. *
  203.  
  204. ORDESP = 1;
  205. ORDTEM = 1;
  206.  
  207.  
  208. ROF1 VITF1 PF1 GAMF1 = 'PRET' 'PERFMONO' ORDESP ORDTEM
  209. $DOMTOT RN
  210. VIT
  211. PN
  212. GAMMA ;
  213.  
  214.  
  215.  
  216. ORDESP = 2;
  217. ORDTEM = 1;
  218.  
  219.  
  220. ROF VITF PF GAMF = 'PRET' 'PERFMONO' ORDESP ORDTEM
  221. $DOMTOT RN GRADR ALR
  222. VIT GRADV ALV
  223. PN GRADP ALP
  224. GAMMA ;
  225.  
  226.  
  227.  
  228. ERRO = 'MAXIMUM' ('PROG'
  229. ('MAXIMUM' (ROF '-' ROF1) 'ABS')
  230. ('MAXIMUM' (PF '-' PF1) 'ABS')
  231. ('MAXIMUM' (VITF '-' VITF1) 'ABS')
  232. ) 'ABS' ;
  233.  
  234. 'SI' (ERRO > 1.0D-15);
  235. 'ERREUR' 5;
  236. 'FINSI' ;
  237.  
  238.  
  239.  
  240. *
  241. **** Limiteurs non-nuls, gradients nuls
  242. *
  243.  
  244.  
  245. GRADR = 'KCHT' $DOMTOT 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  246. (0.0D0 0.0D0);
  247.  
  248. GRADP = 'KCHT' $DOMTOT 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  249. (0.0D0 0.0D0);
  250.  
  251. GRADVX = 'KCHT' $DOMTOT 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  252. (0.0D0 0.0D0);
  253. GRADVY = 'KCHT' $DOMTOT 'VECT' 'CENTRE' 'COMP' 'P2DX' 'P2DY'
  254. (0.0D0 0.0D0);
  255. GRADV = GRADVX 'ET' GRADVY ;
  256.  
  257. ALR = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' 'COMP' 'P1' 0.4D0 ;
  258.  
  259.  
  260. ALP = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' 'COMP' 'P1' 0.8D0 ;
  261.  
  262.  
  263. ALV = 'KCHT' $DOMTOT 'VECT' 'CENTRE' 'COMP' 'P1' 'P2'
  264. (0.9 0.9D0) ;
  265.  
  266.  
  267. ORDESP = 1;
  268. ORDTEM = 1;
  269.  
  270. ROF1 VITF1 PF1 GAMF1 = 'PRET' 'PERFMONO' ORDESP ORDTEM
  271. $DOMTOT RN
  272. VIT
  273. PN
  274. GAMMA ;
  275.  
  276. ORDESP = 2;
  277. ORDTEM = 1;
  278.  
  279. ROF VITF PF GAMF = 'PRET' 'PERFMONO' ORDESP ORDTEM
  280. $DOMTOT RN GRADR ALR
  281. VIT GRADV ALV
  282. PN GRADP ALP
  283. GAMMA ;
  284.  
  285.  
  286. ERRO = 'MAXIMUM' ('PROG'
  287. ('MAXIMUM' (ROF '-' ROF1) 'ABS')
  288. ('MAXIMUM' (PF '-' PF1) 'ABS')
  289. ('MAXIMUM' (VITF '-' VITF1) 'ABS')
  290. ) 'ABS' ;
  291.  
  292. 'SI' (ERRO > 1.0D-15);
  293. 'ERREUR' 5;
  294. 'FINSI' ;
  295.  
  296.  
  297. *
  298. **** TEST2 : Deuxieme ordre en espace, premiere ordre en
  299. * temps: invariance par rotation
  300. *
  301. * Il faut remarquer que:
  302. * densite, pression sont des tenseurs d'ordre 0,
  303. * donc leurs gradients sont des tenseurs d'ordre 1
  304. * la vitesse est un tenseur d'ordre 1,
  305. * donc sont gradient est un tenseur d'ordre 2
  306. *
  307. **** On considere un repere solidal avec la surface, (n,t),
  308. * et un repere (x,y)
  309. *
  310. *** Etats gauche et droite
  311. *
  312.  
  313. gam = 1.4D0;
  314.  
  315. rog = 1.00 ;
  316. grong = 1.0D0 ;
  317. grotg = 2.0D0;
  318. pg = 2.0D0 ;
  319. gpng = 3.0D0 ;
  320. gptg = 4.0D0;
  321. ung = 3.0D0 ;
  322. gunng = 5.0D0 ;
  323. guntg = 6.0D0 ;
  324. utg = 4.0D0 ;
  325. gutng = 7.0D0 ;
  326. guttg = 8.0D0;
  327.  
  328. rod = 1.50 ;
  329. grond = 1.5D0 ;
  330. grotd = 2.5D0 ;
  331. pd = 2.5D0 ;
  332. gpnd = 3.5D0 ;
  333. gptd = 4.5D0 ;
  334. und = 3.5D0 ;
  335. gunnd = 5.5D0 ;
  336. guntd = 6.5D0 ;
  337. utd = 4.5D0 ;
  338. gutnd = 7.5D0 ;
  339. guttd = 8.5D0;
  340.  
  341. *
  342. *** Etats a gauche et a droite du point face P10
  343. *
  344. * ro, un, ut, p
  345.  
  346. ETATG = 'PROG' 1.25000 3.625 4.875 2.1875;
  347. ETATD = 'PROG' 1.12500 2.8125 3.5625 2.28125;
  348.  
  349.  
  350. *
  351. *** Rotation
  352. *
  353.  
  354. ANGLE = 7.0D0;
  355. DANGLE = 85;
  356. ORIG = 0.0D0 0.0D0;
  357.  
  358.  
  359. 'REPETER' BL1 5;
  360.  
  361. ANGLE = ANGLE '+' DANGLE ;
  362.  
  363. uxg uyg = RUO1 ANGLE UNG UTG ;
  364. uxd uyd = RUO1 ANGLE UND UTD ;
  365.  
  366. groxg groyg = RUO1 ANGLE grong grotg;
  367. groxd groyd = RUO1 ANGLE grond grotd;
  368.  
  369. gpxg gpyg = RUO1 ANGLE gpng gptg;
  370. gpxd gpyd = RUO1 ANGLE gpnd gptd;
  371.  
  372. guxxg guxyg guyxg guyyg = RUO2 ANGLE gunng guntg gutng guttg;
  373. guxxd guxyd guyxd guyyd = RUO2 ANGLE gunnd guntd gutnd guttd;
  374.  
  375. 'MESSAGE' ;
  376. 'MESSAGE' (CHAIN 'Angle de rotation= ' ANGLE);
  377. 'MESSAGE' ;
  378.  
  379. DOM1 = DOM10 'TOURNER' ANGLE ORIG;
  380. DOM2 = DOM20 'TOURNER' ANGLE ORIG;
  381. P1 = P10 'TOURNER' ANGLE ORIG;
  382.  
  383. 'ELIMINATION' (DOM1 ET DOM2) 1D-6;
  384. DOMTOT = DOM1 ET DOM2;
  385.  
  386. $DOMTOT = 'MODELISER' DOMTOT 'EULER';
  387.  
  388. $DOM1 = 'MODELISER' DOM1 'EULER';
  389. $DOM2 = 'MODELISER' DOM2 'EULER';
  390.  
  391. TDOMTOT = 'DOMA' $DOMTOT 'VF';
  392.  
  393. TDOM1 = 'DOMA' $DOM1 'VF';
  394. TDOM2 = 'DOMA' $DOM2 'VF';
  395.  
  396. MDOM1 = TDOM1 . 'QUAF' ;
  397. MDOM2 = TDOM2 . 'QUAF' ;
  398. MDOMTOT = TDOMTOT . 'QUAF' ;
  399.  
  400. 'ELIMINATION' (MDOMTOT ET MDOM1) 0.0001 ;
  401. 'ELIMINATION' (MDOMTOT ET MDOM2) 0.0001 ;
  402.  
  403. FACTOT = 'DOMA' $DOMTOT 'FACEL' ;
  404.  
  405. 'SI' GRAPH;
  406. 'TRACER' (DOMTOT 'ET' FACTOT
  407. 'ET' P1) 'TITRE' 'Domaine et FACEL';
  408. 'FINSI' ;
  409.  
  410. *
  411. **** Redefinition de P1 dans $DOMTOT . 'FACE'
  412. *
  413.  
  414. P1 = (TDOMTOT . 'FACE') 'POIN' 'PROC' P1;
  415.  
  416.  
  417. ***********************
  418. **** Les CHPOINTs ****
  419. ***********************
  420.  
  421. GAMMA = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' gam;
  422.  
  423. RN1 = 'KCHT' $DOM1 'SCAL' 'CENTRE' rog;
  424. RN2 = 'KCHT' $DOM2 'SCAL' 'CENTRE' rod;;
  425. RN = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' (RN1 'ET' RN2);
  426.  
  427.  
  428. VIT1 = 'KCHT' $DOM1 'VECT' 'CENTRE' (uxg uyg);
  429. VIT2 = 'KCHT' $DOM2 'VECT' 'CENTRE' (uxd uyd);
  430. VIT = 'KCHT' $DOMTOT 'VECT' 'CENTRE' (VIT1 'ET' VIT2);
  431.  
  432. PN1 = 'KCHT' $DOM1 'SCAL' 'CENTRE' pg;
  433. PN2 = 'KCHT' $DOM2 'SCAL' 'CENTRE' pd;
  434. PN = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' (PN1 'ET' PN2);
  435.  
  436.  
  437.  
  438. *
  439. **** On impose les gradients et le limiteurs
  440. *
  441.  
  442.  
  443. ALR = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' 'COMP' 'P1' 0.5D0 ;
  444.  
  445.  
  446. ALP = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' 'COMP' 'P1' 0.125D0 ;
  447.  
  448.  
  449. ALV = 'KCHT' $DOMTOT 'VECT' 'CENTRE' 'COMP' 'P1' 'P2'
  450. (0.25 0.25D0) ;
  451.  
  452.  
  453. GRADR1 = 'KCHT' $DOM1 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  454. (groxg groyg) ;
  455.  
  456. GRADR2 = 'KCHT' $DOM2 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  457. (groxd groyd) ;
  458.  
  459. GRADR = 'KCHT' $DOMTOT 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  460. (GRADR1 'ET' GRADR2);
  461.  
  462.  
  463.  
  464. GRADP1 = 'KCHT' $DOM1 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  465. (gpxg gpyg) ;
  466.  
  467. GRADP2 = 'KCHT' $DOM2 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  468. (gpxd gpyd) ;
  469.  
  470. GRADP = 'KCHT' $DOMTOT 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  471. (GRADP1 'ET' GRADP2);
  472.  
  473.  
  474.  
  475. GRADVX1 = 'KCHT' $DOM1 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  476. (guxxg guxyg);
  477.  
  478. GRADVX2 = 'KCHT' $DOM2 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  479. (guxxd guxyd);
  480.  
  481. GRADVY1 = 'KCHT' $DOM1 'VECT' 'CENTRE' 'COMP' 'P2DX' 'P2DY'
  482. (guyxg guyyg);
  483.  
  484. GRADVY2 = 'KCHT' $DOM2 'VECT' 'CENTRE' 'COMP' 'P2DX' 'P2DY'
  485. (guyxd guyyd);
  486.  
  487. GRADVX = GRADVX1 'ET' GRADVX2 ;
  488. GRADVY = GRADVY1 'ET' GRADVY2 ;
  489. GRADV = GRADVX 'ET' GRADVY ;
  490.  
  491. ORDESP = 2;
  492. ORDTEM = 1;
  493.  
  494. ROF VITF PF GAMF = 'PRET' 'PERFMONO' ORDESP ORDTEM
  495. $DOMTOT RN GRADR ALR
  496. VIT GRADV ALV
  497. PN GRADP ALP
  498. GAMMA ;
  499.  
  500.  
  501. *********************************************************
  502. *** Control des etats sur la surface qui contient P1 ****
  503. *********************************************************
  504.  
  505.  
  506. GEOP1 = (TDOMTOT. 'FACEL') 'ELEM' 'APPUYE' 'LARGEMENT' P1;
  507. GEOP2 = (TDOMTOT. 'FACE') 'ELEM' 'APPUYE' 'LARGEMENT' P1;
  508.  
  509.  
  510.  
  511. ROGEOP1 = 'REDU' ROF GEOP1;
  512. VGEOP1 = 'REDU' VITF GEOP1;
  513. PGEOP1 = 'REDU' PF GEOP1;
  514. GAMGEOP1 = 'REDU' GAMF GEOP1;
  515. REFGEOP1 = 'REDU' VITF GEOP2;
  516.  
  517.  
  518. rofag = 'EXTRAIRE' ROGEOP1 'SCAL' 1 1 1;
  519. rofad = 'EXTRAIRE' ROGEOP1 'SCAL' 1 1 3;
  520. unfag = 'EXTRAIRE' VGEOP1 'UN ' 1 1 1;
  521. unfad = 'EXTRAIRE' VGEOP1 'UN ' 1 1 3;
  522. utfag = 'EXTRAIRE' VGEOP1 'UT ' 1 1 1;
  523. utfad = 'EXTRAIRE' VGEOP1 'UT ' 1 1 3;
  524. pfag = 'EXTRAIRE' PGEOP1 'SCAL' 1 1 1;
  525. pfad = 'EXTRAIRE' PGEOP1 'SCAL' 1 1 3;
  526.  
  527.  
  528. *
  529. **** Orientation de la normal n de castem par raport a la
  530. * notre; t est par consequence
  531. *
  532.  
  533. NCOS = 'EXTRAIRE' REFGEOP1 'NX' 1 1 1;
  534. NSIN = 'EXTRAIRE' REFGEOP1 'NY' 1 1 1;
  535.  
  536. 'SI' (('ABS' NCOS) > ('ABS' NSIN));
  537. ORIENT = ('COS' ANGLE) '/' NCOS;
  538. 'SINON';
  539. ORIENT = ('SIN' ANGLE) '/' NSIN;
  540. 'FINSI' ;
  541.  
  542. ORIENT = 'SIGN' ORIENT;
  543.  
  544.  
  545. *
  546. **** ORIENT = -1 -> Mon etat gauche est son etat droite
  547. *
  548.  
  549. 'SI' (ORIENT > 0);
  550.  
  551. ERRLIG = 'PROG' rofag (unfag '*' ORIENT)
  552. (utfag '*' ORIENT) pfag ;
  553.  
  554. ERRLID = 'PROG' rofad (unfad '*' ORIENT)
  555. (utfad '*' ORIENT) pfad ;
  556.  
  557. 'SINON' ;
  558.  
  559. ERRLID = 'PROG' rofag (unfag '*' ORIENT)
  560. (utfag '*' ORIENT) pfag ;
  561.  
  562. ERRLIG = 'PROG' rofad (unfad '*' ORIENT)
  563. (utfad '*' ORIENT) pfad ;
  564.  
  565. 'FINSI' ;
  566.  
  567. ERRO = 'MAXIMUM' ('PROG'
  568. ('MAXIMUM' (ETATG '-' ERRLIG) 'ABS')
  569. ('MAXIMUM' (ETATD '-' ERRLID) 'ABS')
  570. );
  571.  
  572.  
  573. 'SI' (ERRO > 1.0D-14)
  574. 'MESSAGE' 'Ordre en espace = 2';
  575. 'MESSAGE' 'Ordre en temps = 1';
  576. 'ERREUR' 5 ;
  577. 'FINSI' ;
  578.  
  579. 'FIN' BL1;
  580.  
  581.  
  582.  
  583.  
  584. *
  585. **** TEST3: deuxieme ordre en espace et en temps
  586. * Invariance par rotation
  587. *
  588. *
  589. **** On considere un repere solidal avec la surface, (n,t),
  590. * et un repere (x,y)
  591. *
  592. *** Etats gauche et droite
  593. *
  594.  
  595. gam = 1.4D0;
  596.  
  597. rog = 1.00 ;
  598. grong = 1.0D0 ;
  599. grotg = 2.0D0;
  600. pg = 2.0D0 ;
  601. gpng = 3.0D0 ;
  602. gptg = 4.0D0;
  603. ung = 3.0D0 ;
  604. gunng = 5.0D0 ;
  605. guntg = 6.0D0 ;
  606. utg = 4.0D0 ;
  607. gutng = 7.0D0 ;
  608. guttg = 8.0D0;
  609.  
  610. rod = 1.50 ;
  611. grond = 1.5D0 ;
  612. grotd = 2.5D0 ;
  613. pd = 2.5D0 ;
  614. gpnd = 3.5D0 ;
  615. gptd = 4.5D0 ;
  616. und = 3.5D0 ;
  617. gunnd = 5.5D0 ;
  618. guntd = 6.5D0 ;
  619. utd = 4.5D0 ;
  620. gutnd = 7.5D0 ;
  621. guttd = 8.5D0;
  622.  
  623.  
  624. *
  625. *** Rotation
  626. *
  627.  
  628. ANGLE = 7.0D0;
  629. DANGLE = 85;
  630. ORIG = 0.0D0 0.0D0;
  631.  
  632.  
  633. 'REPETER' BL1 5;
  634.  
  635. ANGLE = ANGLE '+' DANGLE ;
  636.  
  637. uxg uyg = RUO1 ANGLE UNG UTG ;
  638. uxd uyd = RUO1 ANGLE UND UTD ;
  639.  
  640. groxg groyg = RUO1 ANGLE grong grotg;
  641. groxd groyd = RUO1 ANGLE grond grotd;
  642.  
  643. gpxg gpyg = RUO1 ANGLE gpng gptg;
  644. gpxd gpyd = RUO1 ANGLE gpnd gptd;
  645.  
  646. guxxg guxyg guyxg guyyg = RUO2 ANGLE gunng guntg gutng guttg;
  647. guxxd guxyd guyxd guyyd = RUO2 ANGLE gunnd guntd gutnd guttd;
  648.  
  649. 'MESSAGE' ;
  650. 'MESSAGE' (CHAIN 'Angle de rotation= ' ANGLE);
  651. 'MESSAGE' ;
  652.  
  653. DOM1 = DOM10 'TOURNER' ANGLE ORIG;
  654. DOM2 = DOM20 'TOURNER' ANGLE ORIG;
  655. P1 = P10 'TOURNER' ANGLE ORIG;
  656.  
  657. 'ELIMINATION' (DOM1 ET DOM2) 1D-6;
  658. DOMTOT = DOM1 ET DOM2;
  659.  
  660. $DOMTOT = 'MODELISER' DOMTOT 'EULER';
  661.  
  662. $DOM1 = 'MODELISER' DOM1 'EULER';
  663. $DOM2 = 'MODELISER' DOM2 'EULER';
  664.  
  665. TDOMTOT = 'DOMA' $DOMTOT 'VF';
  666.  
  667. TDOM1 = 'DOMA' $DOM1 'VF';
  668. TDOM2 = 'DOMA' $DOM2 'VF';
  669.  
  670. MDOM1 = TDOM1 . 'QUAF' ;
  671. MDOM2 = TDOM2 . 'QUAF' ;
  672. MDOMTOT = TDOMTOT . 'QUAF' ;
  673.  
  674. 'ELIMINATION' (MDOMTOT ET MDOM1) 0.0001 ;
  675. 'ELIMINATION' (MDOMTOT ET MDOM2) 0.0001 ;
  676.  
  677. FACTOT = 'DOMA' $DOMTOT 'FACEL' ;
  678.  
  679.  
  680. 'SI' GRAPH;
  681. 'TRACER' (DOMTOT 'ET' FACTOT
  682. 'ET' P1) 'TITRE' 'Domaine et FACEL';
  683. 'FINSI' ;
  684.  
  685. *
  686. **** Redefinition de P1 dans $DOMTOT . 'FACE'
  687. *
  688.  
  689. P1 = (TDOMTOT . 'FACE') 'POIN' 'PROC' P1;
  690.  
  691.  
  692. ***********************
  693. **** Les CHPOINTs ****
  694. ***********************
  695.  
  696. GAMMA = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' gam;
  697.  
  698. RN1 = 'KCHT' $DOM1 'SCAL' 'CENTRE' rog;
  699. RN2 = 'KCHT' $DOM2 'SCAL' 'CENTRE' rod;;
  700. RN = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' (RN1 'ET' RN2);
  701.  
  702.  
  703. VIT1 = 'KCHT' $DOM1 'VECT' 'CENTRE' (uxg uyg);
  704. VIT2 = 'KCHT' $DOM2 'VECT' 'CENTRE' (uxd uyd);
  705. VIT = 'KCHT' $DOMTOT 'VECT' 'CENTRE' (VIT1 'ET' VIT2);
  706.  
  707. PN1 = 'KCHT' $DOM1 'SCAL' 'CENTRE' pg;
  708. PN2 = 'KCHT' $DOM2 'SCAL' 'CENTRE' pd;
  709. PN = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' (PN1 'ET' PN2);
  710.  
  711.  
  712.  
  713. *
  714. **** On impose les gradients et le limiteurs
  715. *
  716.  
  717.  
  718. ALR = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' 'COMP' 'P1' 0.5D0 ;
  719.  
  720.  
  721. ALP = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' 'COMP' 'P1' 0.125D0 ;
  722.  
  723.  
  724. ALV = 'KCHT' $DOMTOT 'VECT' 'CENTRE' 'COMP' 'P1' 'P2'
  725. (0.25 0.25D0) ;
  726.  
  727.  
  728. GRADR1 = 'KCHT' $DOM1 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  729. (groxg groyg) ;
  730.  
  731. GRADR2 = 'KCHT' $DOM2 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  732. (groxd groyd) ;
  733.  
  734. GRADR = 'KCHT' $DOMTOT 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  735. (GRADR1 'ET' GRADR2);
  736.  
  737.  
  738.  
  739. GRADP1 = 'KCHT' $DOM1 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  740. (gpxg gpyg) ;
  741.  
  742. GRADP2 = 'KCHT' $DOM2 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  743. (gpxd gpyd) ;
  744.  
  745. GRADP = 'KCHT' $DOMTOT 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  746. (GRADP1 'ET' GRADP2);
  747.  
  748.  
  749.  
  750. GRADVX1 = 'KCHT' $DOM1 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  751. (guxxg guxyg);
  752.  
  753. GRADVX2 = 'KCHT' $DOM2 'VECT' 'CENTRE' 'COMP' 'P1DX' 'P1DY'
  754. (guxxd guxyd);
  755.  
  756. GRADVY1 = 'KCHT' $DOM1 'VECT' 'CENTRE' 'COMP' 'P2DX' 'P2DY'
  757. (guyxg guyyg);
  758.  
  759. GRADVY2 = 'KCHT' $DOM2 'VECT' 'CENTRE' 'COMP' 'P2DX' 'P2DY'
  760. (guyxd guyyd);
  761.  
  762. GRADVX = GRADVX1 'ET' GRADVX2 ;
  763. GRADVY = GRADVY1 'ET' GRADVY2 ;
  764. GRADV = GRADVX 'ET' GRADVY ;
  765.  
  766. ORDESP = 2;
  767. ORDTEM = 2;
  768.  
  769. *
  770. *** L = 1
  771. *
  772.  
  773. DTCFLG = (gam '*' pg) '/' rog ;
  774. DTCFLG = DTCFLG '**' 0.5 ;
  775. DTCFLG = DTCFLG '+' ung ;
  776. DTCFLG = 1 '/' DTCFLG ;
  777.  
  778. DTCFLD = (gam '*' pd) '/' rod ;
  779. DTCFLD = DTCFLD '**' 0.5 ;
  780. DTCFLD = DTCFLD '+' und ;
  781. DTCFLD = 1 '/' DTCFLD ;
  782.  
  783. DTCFLD = 'MINIMUM' ('PROG' DTCFLD DTCFLG);
  784.  
  785.  
  786. ROF VITF PF GAMF = 'PRET' 'PERFMONO' ORDESP ORDTEM
  787. $DOMTOT RN GRADR ALR
  788. VIT GRADV ALV
  789. PN GRADP ALP
  790. GAMMA (1D-3 '*' DTCFLD) ;
  791.  
  792.  
  793. *********************************************************
  794. *** Control des etats sur la surface qui contient P1 ****
  795. *********************************************************
  796.  
  797.  
  798. GEOP1 = (TDOMTOT. 'FACEL') 'ELEM' 'APPUYE' 'LARGEMENT' P1;
  799. GEOP2 = (TDOMTOT. 'FACE') 'ELEM' 'APPUYE' 'LARGEMENT' P1;
  800.  
  801.  
  802.  
  803. ROGEOP1 = 'REDU' ROF GEOP1;
  804. VGEOP1 = 'REDU' VITF GEOP1;
  805. PGEOP1 = 'REDU' PF GEOP1;
  806. GAMGEOP1 = 'REDU' GAMF GEOP1;
  807. REFGEOP1 = 'REDU' VITF GEOP2;
  808.  
  809.  
  810. rofag = 'EXTRAIRE' ROGEOP1 'SCAL' 1 1 1;
  811. rofad = 'EXTRAIRE' ROGEOP1 'SCAL' 1 1 3;
  812. unfag = 'EXTRAIRE' VGEOP1 'UN ' 1 1 1;
  813. unfad = 'EXTRAIRE' VGEOP1 'UN ' 1 1 3;
  814. utfag = 'EXTRAIRE' VGEOP1 'UT ' 1 1 1;
  815. utfad = 'EXTRAIRE' VGEOP1 'UT ' 1 1 3;
  816. pfag = 'EXTRAIRE' PGEOP1 'SCAL' 1 1 1;
  817. pfad = 'EXTRAIRE' PGEOP1 'SCAL' 1 1 3;
  818.  
  819.  
  820. *
  821. **** Orientation de la normal n de castem par raport a la
  822. * notre; t est par consequence
  823. *
  824.  
  825. NCOS = 'EXTRAIRE' REFGEOP1 'NX' 1 1 1;
  826. NSIN = 'EXTRAIRE' REFGEOP1 'NY' 1 1 1;
  827.  
  828. 'SI' (('ABS' NCOS) > ('ABS' NSIN));
  829. ORIENT = ('COS' ANGLE) '/' NCOS;
  830. 'SINON';
  831. ORIENT = ('SIN' ANGLE) '/' NSIN;
  832. 'FINSI' ;
  833.  
  834. ORIENT = 'SIGN' ORIENT;
  835.  
  836.  
  837. *
  838. **** ORIENT = -1 -> Mon etat gauche est son etat droite
  839. *
  840.  
  841. 'SI' (ORIENT > 0);
  842.  
  843. ERRLIG = 'PROG' rofag (unfag '*' ORIENT)
  844. (utfag '*' ORIENT) pfag ;
  845.  
  846. ERRLID = 'PROG' rofad (unfad '*' ORIENT)
  847. (utfad '*' ORIENT) pfad ;
  848.  
  849. 'SINON' ;
  850.  
  851. ERRLID = 'PROG' rofag (unfag '*' ORIENT)
  852. (utfag '*' ORIENT) pfag ;
  853.  
  854. ERRLIG = 'PROG' rofad (unfad '*' ORIENT)
  855. (utfad '*' ORIENT) pfad ;
  856.  
  857. 'FINSI' ;
  858.  
  859. *
  860. **** 'SI' &BL1 > 1, on controle qui rien n'a change!
  861. *
  862.  
  863. 'SI' (&BL1 > 1);
  864. ERRO = 'MAXIMUM' ('PROG'
  865. ('MAXIMUM' (ETATG '-' ERRLIG) 'ABS')
  866. ('MAXIMUM' (ETATD '-' ERRLID) 'ABS')
  867. );
  868.  
  869.  
  870. 'SI' (ERRO > 1.0D-14)
  871. 'MESSAGE' 'Ordre en espace = 2';
  872. 'MESSAGE' 'Ordre en temps = 2';
  873. 'ERREUR' 5 ;
  874. 'FINSI' ;
  875. 'FINSI' ;
  876.  
  877. ETATG = ERRLIG ;
  878. ETATD = ERRLID ;
  879.  
  880. 'FIN' BL1;
  881.  
  882. 'FIN' ;
  883.  
  884.  
  885.  
  886.  
  887.  
  888.  
  889.  
  890.  
  891.  
  892.  
  893.  
  894.  
  895.  
  896.  
  897.  

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