Télécharger 15wedge.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier : 15wedge.dgibi
  2. ******************************************************************
  3. * CALCUL DE L'ECOULEMENT SUPERSONIQUE STATIONNAIRE DANS UN CANAL *
  4. * AVEC RAMPE INCLINEE A 15deg. *
  5. * POWELL AND DEZEEUW, AIAA-91-1542-CP *
  6. * FORMULATION VF COMPRESSIBLE EXPLICITE/IMPLICITE *
  7. * *
  8. * H. PAILLERE/P. GALON TTMF AOUT 1997 *
  9. * *
  10. * MODIF BECCANTINI MARS 97 (reconstruction lineaure exacte) *
  11. * MODIF BECCANTINI MARS 01 (implicite) *
  12. * MODIF, 10/07/01, syntaxe de PENT changée *
  13. ******************************************************************
  14.  
  15.  
  16. 'MESSAGE' 'A mettre a jours' ;
  17. 'FIN' ;
  18.  
  19. TYEL = 'QUA4' ;
  20.  
  21. 'OPTION' 'DIME' 2 'ELEM' TYEL 'ISOV' 'SULI'
  22. 'ECHO' 1 'TRAC' 'X' ;
  23.  
  24. *****************************************************
  25. *****************************************************
  26. ***** PROCEDURE PROLIM *****
  27. *****************************************************
  28. *****************************************************
  29. *
  30. *
  31. **** Cas Euler mono-espèce
  32. *
  33.  
  34. 'DEBPROC' PROLIM ;
  35. 'ARGUMENT' RV * TABLE ;
  36.  
  37. UN = RV . 'UN' ;
  38. UNLIM = RV . 'UNLIM' ;
  39. MAILLIM = 'EXTRAIRE' UNLIM 'MAILLAGE' ;
  40. UNLEV = 'REDU' UN MAILLIM '*' (-1.0) ;
  41.  
  42. UN = UN '+' UNLEV '+' UNLIM ;
  43.  
  44. 'FINPROC' UN ;
  45.  
  46. *********************************************************************
  47. **** Procedure CALC *************************************************
  48. *********************************************************************
  49. *
  50. * Cette procédure, utilisée comme un opérateur, calcule la différence
  51. * entre deux pas de temps toutes les (RVX . 'FCALC') itérations.
  52. * L'evolution de cette différence (erreur absolue) au cours du temps
  53. * est conservée en vue d'un post-traitement.
  54. *
  55. *
  56. 'DEBP' CALC ;
  57. 'ARGU' RVX*'TABLE' ;
  58. *
  59. RV = RVX . 'EQEX' ;
  60. UN = RV . 'UN' ;
  61. NOMDEN = 'EXTRAIRE' (RV . 'LISTCONS') 1 ;
  62. MOMMOM = 'EXTRAIRE' (RV . 'LISTCONS') ('LECT' 2 3) ;
  63. RN = 'EXCO' UN NOMDEN ;
  64. GN = 'EXCO' UN NOMMOM ('MOTS' 'UX' 'UY') ;
  65.  
  66. *
  67. 'SI' ( 'EXIS' RVX 'COMPT') ;
  68. RVX . 'COMPT' = RVX . 'COMPT' + 1 ;
  69. 'SINON' ;
  70. RVX . 'COMPT' = 1 ;
  71. RVX . 'IT' = 'PROG' ;
  72. RVX . 'TPS' = 'PROG' ;
  73. RVX . 'ER1' = 'PROG' ;
  74. RVX . 'ERINF' = 'PROG' ;
  75. RVX . 'GN0' = 'COPI' GN ;
  76. RVX . 'CFLX' = 'PROG' ;
  77. RVX . 'CFLY' = 'PROG' ;
  78. RVX . 'RESLINY' = 'PROG' ;
  79. 'FINSI' ;
  80. *
  81. DD = RVX . 'COMPT' ;
  82. COMBIEN = RVX . 'FCALC' ;
  83. NN = DD '/' COMBIEN ;
  84. LO = (DD '-' (COMBIEN '*' NN)) 'EGA' 0 ;
  85. 'SI' LO ;
  86. ERR = GN '-' (RVX . 'GN0') ;
  87. ERRINF = 'MAXIMUM' ERR 'ABS' ;
  88. ERR = ('XTX' ERR) '**' 0.5 ;
  89. ERR = ERR '/' ('NBNO' ('EXTRAIRE' GN 'MAILLAGE')) ;
  90. RVX . 'GN0' = 'COPI' GN ;
  91. EL1 = ('LOG' (ERR + 1.0E-50))/('LOG' 10.) ;
  92. ELINF = ('LOG' (ERRINF + 1.0E-50))/('LOG' 10.) ;
  93. 'MESSAGE' ;
  94. 'MESSAGE' 'Erreur de convergence' ;
  95. 'MESS' 'ITERATION ' (RVX . 'COMPT') ' LOG10 ERR_1(GN)' EL1 ;
  96. 'MESS' 'ITERATION ' (RVX . 'COMPT') ' LOG10 ERR_I(GN)' ELINF ;
  97. 'MESSAGE' ;
  98. IT = 'PROG' RVX . 'COMPT' ;
  99. ER1 = 'PROG' EL1 ;
  100. ERINF = 'PROG' ELINF ;
  101. TPS = 'PROG' (RV . 'PASDETPS' . 'TPS') ;
  102. RVX . 'IT' = (RVX . 'IT') 'ET' IT ;
  103. RVX . 'ER1' = (RVX . 'ER1') 'ET' ER1 ;
  104. RVX . 'ERINF' = (RVX . 'ERINF') 'ET' ERINF ;
  105. RVX . 'TPS' = (RVX . 'TPS') 'ET' TPS ;
  106. 'FINS' ;
  107. *
  108. * On peut augmunter le CFL (et MATTAB . 'RESID')
  109. *
  110. 'SI' VRAI ;
  111. COMBIEN = RVX . 'FCFL' ;
  112. NN = DD '/' COMBIEN ;
  113. LO = (DD '-' (COMBIEN '*' NN)) 'EGA' 0 ;
  114. 'SI' LO ;
  115. CFL = PKONTAB . 'ALPHA' ;
  116. CFL = CFL '*' 2. ;
  117. CFL = 'MINIMUM' ('PROG' CFL (RVX . 'ALPMAX')) ;
  118. PKONTAB . 'ALPHA' = CFL ;
  119. 'MESSAGE' ;
  120. 'MESSAGE' ('CHAINE' 'CFL =' ((PKONTAB . 'ALPHA') '/' 2)) ;
  121. * We decrease the linear iterations error
  122. MATTAB . 'RESID' = 'MAXIMUM' ('PROG'
  123. (MATTAB . 'RESID' '/' 3.1623) 1.0D-14) ;
  124. 'MESSAGE' ('CHAINE' 'MATTAB . RESID = ' MATTAB . 'RESID') ;
  125. 'FINSI' ;
  126. RVX . 'CFLX' = RVX . 'CFLX' 'ET' ('PROG' DD) ;
  127. RVX . 'CFLY' = RVX . 'CFLY' 'ET' ('PROG'
  128. ((PKONTAB . 'ALPHA') '/' 2)) ;
  129. RVX . 'RESLINY' = (RVX . 'RESLINY') 'ET'
  130. ('PROG' MATTAB . 'RESID') ;
  131. 'FINSI' ;
  132.  
  133. IRESU IJACO ='KOPS' 'MATRIK' ;
  134. IALPDT = 'PROG' ;
  135.  
  136. 'MENAGE' ;
  137.  
  138. 'RESPRO' IRESU IALPDT ;
  139.  
  140. 'FINP' ;
  141.  
  142. ************************************************************************
  143. ************************************************************************
  144. ***************** FIN PARTIE PROCEDURES ********************************
  145. ************************************************************************
  146. ************************************************************************
  147.  
  148. GRAPH = VRAI ;
  149. * GRAPH = FAUX ;
  150.  
  151. ***********************************************************
  152. *********************** MAILLAGE **************************
  153. ***********************************************************
  154.  
  155. NRAFF = 8 ;
  156.  
  157. NX1 = 1 '*' NRAFF ;
  158. NX2 = 1 '*' NRAFF ;
  159. NX3 = 4 '*' NRAFF ;
  160. NY = 2 '*' NRAFF ;
  161. DY = 2.0 '/' NY ;
  162.  
  163. H1 = 0.268 ;
  164.  
  165. A1 = (1. '-' DY) 0.0 ;
  166. A2 = 1.0 0.0 ;
  167. A3 = 2.0 0.0 ;
  168. A4 = 3.0 H1 ;
  169. A5 = 7.0 H1 ;
  170. A6 = (7.0 '+' DY) H1 ;
  171.  
  172. B1 = (1.0 '-' DY) 2.0 ;
  173. B2 = 1.0 2.0 ;
  174. B3 = 2.0 2.0 ;
  175. B4 = 3.0 2.0 ;
  176. B5 = 7.0 2.0 ;
  177. B6 = (7.0 '+' DY) 2.0 ;
  178.  
  179. BAS1 = A2 'DROI' NX1 A3 ;
  180. BAS2 = A3 'DROI' NX2 A4 ;
  181. BAS3 = A4 'DROI' NX3 A5 ;
  182. HAU3 = B5 'DROI' NX3 B4 ;
  183. HAU2 = B4 'DROI' NX2 B3 ;
  184. HAU1 = B3 'DROI' NX1 B2 ;
  185.  
  186. COTEG = B2 'DROI' NY A2 ;
  187. COTED = A5 'DROI' NY B5 ;
  188.  
  189. FRONTG = 'DALLER' (A1 'DROI' 1 A2) ('INVE' COTEG)
  190. (B2 'DROI' 1 B1) (B1 'DROI' NY A1) ;
  191. FRONTG = FRONTG 'COUL' ROUG ;
  192.  
  193. FRONTD = 'DALLER' (A5 'DROI' 1 A6) (A6 'DROI' NY B6)
  194. (B6 'DROI' 1 B5) ('INVE' COTED) ;
  195. FRONTD = FRONTD 'COUL' VERT ;
  196.  
  197. FRONT = FRONTG 'ET' FRONTD ;
  198. 'ELIMINATION' FRONT 1D-6;
  199.  
  200. DOMINT = 'DALLER' (BAS1 ET BAS2 ET BAS3) COTED
  201. (HAU3 ET HAU2 ET HAU1) COTEG 'PLAN' ;
  202.  
  203. DOMTOT = DOMINT ET FRONT ;
  204. 'ELIMINATION' DOMTOT (1D-6 '/' NRAFF) ;
  205.  
  206. $DOMTOT = 'DOMA' DOMTOT ;
  207. $DOMINT = 'DOMA' DOMINT 'INCL' $DOMTOT (1D-6 '/' NRAFF) ;
  208. $FRONT = 'DOMA' FRONT 'INCL' $DOMTOT (1D-6 '/' NRAFF) ;
  209.  
  210. 'DOMA' $DOMINT 'XXNORMAF' ;
  211. CHPONOR = $DOMINT . 'XXNORMAF' ;
  212.  
  213. *****************************************************************
  214. * CONDITIONS INITIALES ET LIMITES : ECOULEMENT SUPERSONIQUE M=2 *
  215. * Adimensionalisation cf. These HP. p. 152 *
  216. *****************************************************************
  217.  
  218. gamscal = 1.4 ;
  219. GAMN = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' gamscal;
  220.  
  221. Minf = 2.0D0 ;
  222. alpha = 0.D0 ;
  223.  
  224. 'SI' ('NEG' ALPHA 0.0) ;
  225. 'MESSAGE' 'Attention aux C.L.' ;
  226. 'ERREUR' 21 ;
  227. 'FINSI' ;
  228.  
  229. rhoinf = 1.0D0 ;
  230. rhouinf = 'COS' alpha ;
  231. rhovinf = 'SIN' alpha ;
  232. rhoEtinf = 2.0D0 '+' (gamscal '*' Minf '*' Minf
  233. '*' (gamscal '-' 1.0D0)) ;
  234. rhoEtinf = rhoEtinf '/' (2.D0 '*' gamscal '*' Minf
  235. '*' Minf '*' (gamscal '-' 1.0D0)) ;
  236.  
  237. RN0 = 'MANUEL' 'CHPO' ($DOMTOT . 'CENTRE') 1 'SCAL' rhoinf ;
  238. GN0 = 'MANUEL' 'CHPO' ($DOMTOT . 'CENTRE') 2 'UX' rhouinf
  239. 'UY' rhovinf ;
  240. RETN0 = 'MANUEL' 'CHPO' ($DOMTOT . 'CENTRE') 1 'SCAL' rhoEtinf ;
  241.  
  242. VN0 PN0 = 'PRIM' 'PERFMONO' RN0 GN0 RETN0 GAMN ;
  243.  
  244. VN2 = 'PSCAL' VN0 VN0 ('MOTS' 'UX' 'UY') ('MOTS' 'UX' 'UY') ;
  245. C2 = GAMN '*' ( PN0 '/' RN0) ;
  246.  
  247. MACH2 = VN2 '/' C2;
  248. MACHN0 = MACH2 '**' 0.5;
  249.  
  250. MOD1 = 'MODELISER' ($DOMTOT . 'MAILLAGE' ) 'THERMIQUE';
  251.  
  252. 'SI' GRAPH ;
  253.  
  254. CHM_RN = 'KCHA' $DOMTOT 'CHAM' RN0 ;
  255. CHM_PN = 'KCHA' $DOMTOT 'CHAM' PN0 ;
  256. CHM_VN = 'KCHA' $DOMTOT 'CHAM' VN0 ;
  257. CHM_MACH = 'KCHA' $DOMTOT 'CHAM' MACHN0 ;
  258.  
  259. 'TRACER' CHM_RN MOD1
  260. 'TITR' ('CHAINE' 'RN at t=' 0.0);
  261. 'TRACER' CHM_PN MOD1
  262. 'TITR' ('CHAINE' 'PN at t=' 0.0);
  263. 'TRACER' CHM_MACH MOD1
  264. 'TITR' ('CHAINE' 'MACHN at t=' 0.0);
  265. 'TRACER' CHM_VN MOD1
  266. 'TITR' ('CHAINE' 'VN at t=' 0.0);
  267.  
  268. 'FINSI' ;
  269.  
  270. NOMDEN = 'RN' ;
  271. NOMMOM = 'MOTS' 'RUX' 'RUY' ;
  272. NOMENE = 'RETN' ;
  273. LISTCONS = ('MOTS' NOMDEN) 'ET' NOMMOM 'ET' ('MOTS' NOMENE) ;
  274.  
  275. UN =
  276. ('NOMC' NOMDEN RN0 'NATURE' 'DISCRET') 'ET'
  277. ('EXCO' ('MOTS' 'UX' 'UY') GN0 NOMMOM
  278. 'NATURE' 'DISCRET') 'ET'
  279. ('NOMC' NOMENE RETN0 'NATURE' 'DISCRET') ;
  280.  
  281. *****************************************************
  282. *****************************************************
  283. *****************************************************
  284. *****************************************************
  285. **************** La table RV **********************
  286. *****************************************************
  287. *****************************************************
  288. *****************************************************
  289. *****************************************************
  290. RV = 'TABLE' ;
  291. * Constant time step
  292. * RV . 'DTIMP' = 1.191 ;
  293. *
  294. RV . 'UN' = UN ;
  295. * Primitive formulation
  296. RV . 'FORMPRIM' = FAUX ;
  297. RV . 'LISTCONS' = LISTCONS ;
  298. RV . 'LISTPRIM' = LISTCONS ;
  299. **** Le domaine interne
  300. RV . 'DOMINT' = $DOMINT . 'CENTRE' ;
  301. *
  302. RV . 'CLIM' = VRAI ;
  303. RV . 'MAIFAN' = $FRONT . 'CENTRE' ;
  304. *
  305. RV . 'INST' = VRAI ;
  306. RV . 'CONS' = FAUX ;
  307. *
  308. 'SI' (RV . 'INST') ;
  309. RV . 'TFINAL' = 10.0D5 ;
  310. RV . 'NITMAEX' = 100 ;
  311. 'FINSI' ;
  312. *
  313. RV . 'NITMAIN' = 1 ;
  314. RV . 'NITMIIN' = 1 ;
  315. RV . 'EPSINT' = 3.0D-10 ;
  316. *
  317. RV . 'FEXT' = 20 ;
  318. RV . 'FINT' = 1 ;
  319. *
  320. RV . 'LISTOPER' = 'MOTS' 'CALC' 'PKON' ;
  321. RV . '1CALC' = 'TABLE' ;
  322. CALCTAB = RV . '1CALC' ;
  323. RV . '2PKON' = 'TABLE' ;
  324. PKONTAB = RV . '2PKON' ;
  325. *
  326. RV . 'MATIDE' = 'KOPS' 'MATIDE' LISTCONS ($DOMTOT . 'CENTRE') ;
  327. *
  328. **** Parametres de la procedure PROLIM
  329. *
  330. RV . 'UNLIM' = 'REDU' (RV . 'UN') (RV . 'MAIFAN') ;
  331. * Gas model (Procedure EXEXPRIM, EXEXCONV)
  332. RV . 'GAMMA' = GAMN ;
  333. * Bas Mach (Procedures EXEXDT, PKON)
  334. RV . 'ORDTPS' = 1 ;
  335. RV . 'PRECTEMP' = FAUX ;
  336. * RV . DUALTEXP = VRAI ;
  337. * RV . DUALTIMP = VRAI ;
  338. RV . CFLDTAU = 1.0 ;
  339. RV . 'DIAM' = ('DOMA' $DOMTOT 'DIAMIN') '*' 1.0 ;
  340. * RV . 'CO1' = 'MANUEL' 'CHPO' ($DOMTOT . 'CENTRE')
  341. * 1 'SCAL' (0.1 * U_G) '+'
  342. * ('MANUEL' 'CHPO' ($FRONT . 'CENTRE')
  343. * 1 'SCAL' (0. * U_G)) ;
  344. * RV . 'CO2' = 'MANUEL' 'CHPO' ($DOMTOT . 'CENTRE')
  345. * 1 'SCAL' (0. * U_G) ;
  346. *****************
  347. ** CALCUL *******
  348. *****************
  349. CALCTAB . 'EQEX' = RV ;
  350. CALCTAB . 'FCALC' = 5 ;
  351. CALCTAB . 'IMPL' = FAUX ;
  352. * Doubling frequence for CFL
  353. CALCTAB . 'FCFL' = 5 ;
  354. * Maximum of CFL
  355. CALCTAB . 'ALPMAX' = 10000.0 ;
  356.  
  357. *****************
  358. * PKON **********
  359. *****************
  360. PKONTAB . 'EQEX' = RV ;
  361. PKONTAB . 'GAZ' = 'PERFMONO' ;
  362. PKONTAB . 'GAMN' = GAMN ;
  363. PKONTAB . 'DOMA' = $DOMTOT ;
  364. PKONTAB . 'LISTCONS' = LISTCONS ;
  365. * PKONTAB . 'LISTPRIM' = LISTPRIM ;
  366. PKONTAB . 'METHODE' = 'VLH' ;
  367. *
  368. *
  369. PKONTAB . 'IMPL' = VRAI ;
  370. 'SI' (PKONTAB . 'IMPL') ;
  371. PKONTAB . 'TYPEJACO' = 'VLH' ;
  372. 'FINSI' ;
  373. *
  374. PKONTAB . 'ORDREESP' = 2 ;
  375. 'SI' ((PKONTAB . 'ORDREESP') 'EGA' 2) ;
  376. PKONTAB . 'LIMITEUR' = 'LIMITEUR' ;
  377. CHPVID MATVID = 'KOPS' 'MATRIK' ;
  378. PKONTAB . 'VLIM' = CHPVID ;
  379. CACCA PIPI GRADSCAL = 'PENTE' $DOMTOT 'CENTRE' 'EULESCAL'
  380. 'NOLIMITE' ('EXCO' NOMDEN (RV . 'UN')) ;
  381. CACCA PIPI GRADVEC = 'PENTE' $DOMTOT 'CENTRE' 'EULEVECT'
  382. 'NOLIMITE' ('EXCO' NOMMOM (RV . 'UN')
  383. ('MOTS' 'UX' 'UY')) 'CLIM' (PKONTAB . 'VLIM') ;
  384. PKONTAB . 'GRADRN' = GRADSCAL ;
  385. PKONTAB . 'GRADPN' = GRADSCAL ;
  386. PKONTAB . 'GRADVN' = GRADVEC ;
  387. 'FINSI' ;
  388. *
  389. ELP1 = 'MANUEL' 'POI1' (123.56 124.87) ;
  390. PKONTAB . 'FACELIM' = 'DIFF' ELP1 ELP1 ;
  391. PKONTAB . 'ALPHA' = 10. ;
  392. ***********************************
  393. * Inversion de la matrice *********
  394. ***********************************
  395. * RV . 'LINEHIST' = historique de convergence des iterations internes
  396. * RV . 'MATUPDAT' updating
  397. * RV . 'MUPINI' we update de matrix if the external iterations
  398. * are less than (RV . 'MUPINI')
  399. * RV . 'MUPEXT' external iterations updating
  400. * RV . 'MUPINT' internal iterations updating
  401. * RV . 'MATACC' Broyden method?
  402. * RV . 'MUPLIN' We update the matrix if in the previous
  403. * internal iteration the number of linear
  404. * iterations to solve the system were bigger
  405. * than (RV . 'MUPLIN')
  406. RV . 'LINEHIST' = 'PROG' ;
  407. RV . 'MATUPDAT' = VRAI ;
  408. 'SI' ('NON' ( RV . 'MATUPDAT')) ;
  409. RV . 'MUPEXT' = 4000 ;
  410. RV . 'MUPINT' = 200 ;
  411. RV . 'MATACC' = FAUX ;
  412. RV . 'MUPINI'= 0 ;
  413. RV . 'MUPLIN' = 80 ;
  414. 'FINSI' ;
  415. RV . 'MATINV' = 'TABLE' 'METHINV' ;
  416. MATTAB = RV . 'MATINV' ;
  417. MATTAB . 'TYPINV' = 5 ;
  418. MATTAB . 'IMPINV' = 0 ;
  419. *
  420. * Matrice pour assurer que la matrice à inverser est correctement
  421. * assemblé
  422. *
  423. * MATTAB . 'MATASS' definie en EXEXIM
  424. * MATTAB . 'MAPREC' "
  425. * MATTAB . 'XINIT' "
  426. *
  427.  
  428. * Methode de numerotation de DDL
  429. MATTAB . 'TYRENU' = 'RIEN' ;
  430. MATTAB . 'PCMLAG' = 'APR2' ;
  431. MATTAB . 'GMRESTRT' = 500 ;
  432. * ILU 3
  433. * ILUT (dual truncation) 5
  434. * Dans le cas ILUT il faut definir
  435. * ILUTLFIL
  436. * ILUTDTOL
  437. MATTAB . 'PRECOND' = 3 ;
  438. * MATTAB . 'ILUTLFIL' = 32. ;
  439. * MATTAB . 'ILUTDTOL' = 0. ;
  440. MATTAB . 'NITMAX' = 1000 ;
  441. MATTAB . 'RESID' = 1.D-3 ;
  442. * bcgstab
  443. MATTAB . 'BCGSBTOL' = 1.0d-40 ;
  444.  
  445. 'OPTION' DONN 5 ;
  446.  
  447. 'TEMPS' 'ZERO' ;
  448. EXEXIM RV ;
  449. 'TEMPS' ;
  450.  
  451. *
  452. **** Les variables conservatives
  453. *
  454.  
  455. 'SI' (RV . 'INST') ;
  456. TFIN = RV . 'PASDETPS' . 'TPS' ;
  457. 'SINO' ;
  458. TFIN = 'INF' ;
  459. 'FINSI' ;
  460.  
  461. RN = 'EXCO' 'RN' (RV . 'UN') ;
  462. GN = 'EXCO' ('MOTS' 'RUX' 'RUY') (RV . 'UN') ('MOTS' 'UX' 'UY') ;
  463. RETN = 'EXCO' 'RETN' (RV . 'UN') ;
  464.  
  465. *
  466. **** Les variables primitives
  467. *
  468.  
  469. VN PN = 'PRIM' 'PERFMONO'
  470. RN GN RETN GAMN ;
  471.  
  472. CSON2 = (GAMN '*' PN) '/' RN ;
  473. VN2 = 'PSCAL' VN VN ('MOTS' 'UX' 'UY') ('MOTS' 'UX' 'UY') ;
  474. MACHN2 = VN2 '/' CSON2 ;
  475. MACHN = MACHN2 '**' 0.5 ;
  476.  
  477. *
  478. *** GRAPHIQUE DES SOLUTIONS
  479. *
  480.  
  481. MOD1 = 'MODELISER' $DOMINT 'THERMIQUE';
  482.  
  483. 'SI' GRAPH ;
  484.  
  485. CHM_RN ='KCHA' $DOMINT 'CHAM' ('REDU' RN ($DOMINT . 'CENTRE')) ;
  486. CHM_PN ='KCHA' $DOMINT 'CHAM' ('REDU' PN ($DOMINT . 'CENTRE')) ;
  487. CHM_VN = 'KCHA' $DOMINT 'CHAM' ('REDU' VN ($DOMINT . 'CENTRE')) ;
  488. CHM_MACH= 'KCHA' $DOMINT 'CHAM' ('REDU' MACHN ($DOMINT . 'CENTRE'));
  489.  
  490.  
  491. 'TRACER' CHM_RN MOD1 ('CONTOUR' DOMINT)
  492. 'TITR' ('CHAINE' 'METO =' METO ' NBEL =' NTOTEL
  493. ' RN at t=' TFIN ' OE =' 2);
  494. 'TRACER' CHM_PN MOD1 ('CONTOUR' DOMINT)
  495. 'TITR' ('CHAINE' 'METO =' METO ' NBEL =' NTOTEL
  496. ' PN at t=' TFIN ' OE =' 2);
  497. 'TRACER' CHM_MACH MOD1 ('CONTOUR' DOMINT)
  498. 'TITR' ('CHAINE' 'METO =' METO ' NBEL =' NTOTEL
  499. ' MACHN at t=' TFIN ' OE =' 2);
  500. 'TRACER' CHM_VN MOD1 ('CONTOUR' DOMINT)
  501. 'TITR' ('CHAINE' 'METO =' METO ' NBEL =' NTOTEL
  502. ' VN at t=' TFIN ' OE =' 2);
  503.  
  504.  
  505. RNC = 'ELNO' $DOMINT ('REDU' RN ($DOMINT . 'CENTRE'));
  506. 'TRACER' DOMINT RNC ('CONTOUR' DOMINT)
  507. 'TITR' ('CHAINE' 'METO =' METO ' NBEL =' NTOTEL
  508. ' RN at t=' TFIN ' OE =' 2) 15 ;
  509.  
  510. 'FINSI' ;
  511.  
  512. *
  513. *** Test
  514. *
  515.  
  516. PCEN1 = ($DOMTOT . 'CENTRE') 'POIN' 'PROC' (2.5 (H1 '/' 2.0)) ;
  517. 'SI' GRAPH ;
  518. 'TRACER' (DOMTOT 'ET' PCEN1) 'TITRE'
  519. 'Test 1: solution derrier le premier choc' ;
  520. 'FINSI' ;
  521. ro1 = 'EXTRAIRE' RN 'SCAL' PCEN1 ;
  522.  
  523. erro = ((ro1 '-' 1.7289) 'ABS') '/' 1.7289 ;
  524. 'SI' (erro > 1.0D-1) ;
  525. 'ERREUR' 5 ;
  526. 'FINSI' ;
  527.  
  528. p1 = 'EXTRAIRE' PN 'SCAL' PCEN1 ;
  529.  
  530. erro = ((p1 '-' 0.3919) 'ABS') '/' 0.3919 ;
  531. 'SI' (erro > 1.0D-1) ;
  532. 'ERREUR' 5 ;
  533. 'FINSI' ;
  534.  
  535. m1 = 'EXTRAIRE' MACHN 'SCAL' PCEN1 ;
  536.  
  537. erro = ((m1 '-' 1.44572) 'ABS') '/' 1.44572 ;
  538. 'SI' (erro > 1.0D-1) ;
  539. 'ERREUR' 5 ;
  540. 'FINSI' ;
  541.  
  542. * Derrier le mach stem on a le mach min
  543. * MMIN = 0.577350
  544. * Il faut avoir un maillage assez fin pour l'avoir
  545.  
  546. MMIN = 'MINIMUM' MACHN ;
  547. 'MESSAGE' ('CHAINE' 'MMIN_ex =' 0.577350) ;
  548. 'MESSAGE' ('CHAINE' 'MMIN =' MMIN) ;
  549.  
  550. 'SI' (MMIN > 1.5) ;
  551. 'ERREUR' 5 ;
  552. 'FINSI' ;
  553.  
  554. 'FIN' ;
  555.  
  556.  
  557.  
  558.  
  559.  
  560.  
  561.  
  562.  
  563.  
  564.  
  565.  

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