Télécharger benchmark_imst.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier : benchmark_imst.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. ***************************************************************
  5. * *
  6. * NOM : benchmark_imst.dgibi *
  7. * DESCRIPTION : fiche de validation CASTEM2000 *
  8. * Mecanique des Fluides *
  9. * Convection naturelle laminaire *
  10. * en cavite rectangulaire (2D) *
  11. * FONCTIONS *
  12. * TESTEES : 'NS' algorithme semi-implicite *
  13. * option CENTREE *
  14. * approximation de Boussinesq *
  15. * elements QUA8 'MACRO' *
  16. * *
  17. * AUTEUR : Stephane GOUNAND (CEA/DRN/DMT/SEMT/TTMF) *
  18. * e-mail : gounand@semt2.smts.cea.fr *
  19. * *
  20. ***************************************************************
  21. * *
  22. * VERSION : v1, 19/9/97, version initiale *
  23. * HISTORIQUE : v1, 19/9/97, creation *
  24. * HISTORIQUE : 15/06/2014 passage EQPR -> EQEX *
  25. * HISTORIQUE : 11/07/2023 fixer la pression en 1 point *
  26. * HISTORIQUE : *
  27. * HISTORIQUE : *
  28. * *
  29. ***************************************************************
  30. * *
  31. * Priere de completer la partie HISTORIQUE apres modification *
  32. * du jeu de donnees afin de faciliter la maintenance *
  33. * *
  34. ***************************************************************
  35. 'OPTI' echo 0 ;
  36. *
  37. * DESCRIPTION DETAILLEE :
  38. *
  39. * -- IMST -- Convection naturelle laminaire à Prandtl (Pr) = 0
  40. * Algorithme semi implicite
  41. * operateur NS (option CENTREE)
  42. *
  43. * phi = 0
  44. * ------------------------
  45. * | |
  46. * T1 = 0 | H | T2 = A
  47. * | L |
  48. * ------------------------
  49. * phi = 0
  50. *
  51. * A = L/H = 4 l'allongement
  52. * U = 0 sur les parois
  53. * Le profil de temperature est impose T = T1 + (T2-T1)*(x/L)
  54. *
  55. * La structure de l'ecoulement ne depend que
  56. * du nombre de Grashof Gr :
  57. * Gr = g beta (T2 - T1) H**3 A / nu**2
  58. *
  59. * ^
  60. * | Psi max (valeur max de la fonction de courant)
  61. * |
  62. * |
  63. * -| o o
  64. * | unsteady -> o v
  65. * | (3 cells) . v
  66. * | ----------- o v . o
  67. * | steady -> .^ v o
  68. * | (3 cells) . ^ v .
  69. * | ----------- o ^ o
  70. * -| ^ .
  71. * | steady -> . ^ . <-- steady (2 cells)
  72. * | o
  73. * | (1 cell) .
  74. * |
  75. * | .
  76. * | o
  77. * -|-------------------------------------------> Gr
  78. * I I I
  79. * 1.E4 2.E4 3.E4
  80. **************************************************************
  81. * References : - Gounand, S. 1997
  82. * Simulation numerique de cellules
  83. * thermoconvectives.
  84. * Rapport CEA/DMT 97/357.
  85. * - Le Garrec, S. 1988
  86. * Convection naturelle en cavite pour des
  87. * fluides à faible nombre de Prandtl.
  88. * Rapport CEA/DMT 88/216.
  89. * - Le Garrec, S. & Magnaud, J.-P.
  90. * Numerical Simulation of oscillatory
  91. * convection in low Prandtl Fluids with
  92. * TRIO-EF.
  93. * A GAMM Workshop Edited by Bernard Roux
  94. * Notes on Numerical Fluid Mechanics,Volume 27
  95. * pp. 189-199. (Vieweg, Braunschweig 1990)
  96. *
  97. ***************************************************************
  98. *
  99. * JEU DE DONNEE :
  100. *
  101. prg = ('CHAINE' 'Benchmark Imst') ;
  102. **** Options
  103. * court = VRAI : le cas-test s'effectue en une dizaine
  104. * de secondes.
  105. * Gr = 5000 => une cellule stationnaire
  106. * la comparaison s'effectue sur le profil
  107. * Uy(x) le long de la mediane horizontale.
  108. * On verifie egalement que l'ecoulement
  109. * est centrosymetrique.
  110. * court = FAUX : maillage raffine ;
  111. * Gr = 40000 => bifurcation 3->2 cellules
  112. * environ 1h30 CPU avec les graphiques
  113. *
  114. * graph = VRAI : on trace les graphiques (historiques, champs)
  115. * On a plus de calculs avec cette option
  116. * (fonction de courant, convergence...)
  117. * On mettra donc graph = faux pour la fiche
  118. * de validation
  119. *
  120. * inta = VRAI si on est en interactif ; 0 sinon
  121. *
  122. * discconv : option de discretisation des termes de
  123. * convection pour les equations de Navier-Stokes
  124. * choix possibles : 'CENTREE', 'SUPG'...
  125. * voir notice de 'NS'
  126. court = VRAI ;
  127. graph = FAUX ;
  128. inta = FAUX ;
  129. typelem = 'QUA8' ;
  130. discconv = 'CENTREE' ;
  131.  
  132. DISCR = 'MACRO' ;
  133. KPRESS = 'CENTRE' ;
  134. * betap : parametre de stabilisation de la pression
  135. BETAP = 1. ;
  136.  
  137. 'OPTION' 'DIME' 2 'ELEM' typelem ;
  138.  
  139. 'SI' ('NON' inta) ;
  140. 'OPTION' 'TRAC' 'PSC';
  141. 'SINON' ;
  142. **od = 'TEXTE' 'OPTION DONN 3' ;
  143. 'FINSI' ;
  144.  
  145. *************************************************************
  146. **** Definition de quelques procedures :
  147. ** -> SUMMARY : affiche un resume des parametres
  148. ** -> CALCPSI : calcul de la fonction de courant
  149. ** (appelee par 'EXEC' via la table 'EQEX')
  150. ** -> MAJHIST : mise à jour d'historiques
  151. ** (appelee par 'EXEC' via la table 'EQEX')
  152. ** -> TRHIST : trace des historiques
  153. ** -> TRCHAMP : trace des differents champs
  154. ** (appelee par 'EXEC' via la table 'EQEX')
  155. ** NB : pour simplifier, on ne transmet pas les parametres à
  156. ** SUMMARY. C'est de la mauvaise programmation...
  157. ** Les autres procedures sont correctes
  158. ** sous toutes reserves
  159. *
  160. ** SUMMARY
  161. 'DEBPROC' SUMMARY ;
  162. 'MESSAGE' ('CHAINE' 'Programme :' ' ' prg) ;
  163. 'SAUTER' 'LIGNE' ;
  164. 'MESSAGE' ('CHAINE' '*** GEOMETRIE :') ;
  165. 'MESSAGE' ('CHAINE' 'Boite carree :' ' ' l ' x' h) ;
  166. 'MESSAGE' ('CHAINE' ' soient :' ' ' nl 'x' nh '=' (nl '*' nh) ' elements' ' ' typelem) ;
  167. 'MESSAGE' ('CHAINE' ' soient :' ' ' (nl '*' nh '*' 4) ' elements lineaires') ;
  168. 'SAUTER' 'LIGNE' ;
  169. 'MESSAGE' ('CHAINE' '*** PHYSIQUE :') ;
  170. 'MESSAGE' ('CHAINE' 'Rapport d_aspect : ' ' ' ('ENTIER' aspect)) ;
  171. 'MESSAGE' ('CHAINE' 'Nombre de Grashof :' ' '('ENTIER' Gr) ) ;
  172. 'MESSAGE' ('CHAINE' 'Temperature de reference :' ' ' tref ' comprise entre 0 et ' ('ENTIER' aspect)) ;
  173. 'SAUTER' 'LIGNE' ;
  174. 'MESSAGE' ('CHAINE' '*** NUMERIQUE :') ;
  175. 'MESSAGE' ('CHAINE' 'Discretisation' ' ' discconv ' des termes de convection.');
  176. 'MESSAGE' ('CHAINE' ' Alfa (mult. pasdt.) =' ' ' alpha) ;
  177. 'MESSAGE' ('CHAINE' ' Beta (stab. pression) =' ' ' betap) ;
  178. 'SAUTER' 'LIGNE' ;
  179. rvpdt = rv . 'PASDETPS' ;
  180. 'MESSAGE' ('CHAINE' '*** TEMPS :') ;
  181. 'MESSAGE' ('CHAINE' 'Nb. iterations :' ' ' ((rvpdt . 'NUPASDT') '-' 1)) ;
  182. 'MESSAGE' ('CHAINE' 'Temps physique :' ' ' (rvpdt . 'TPS')) ;
  183. 'MESSAGE' ('CHAINE' 'Temps CPU :' ' ' ('ENTIER' (thist . 'tpscpu')) ' ms') ;
  184. 'FINPROC' ;
  185.  
  186.  
  187. ** CALCPSI
  188. 'DEBPROC' CALCPSI ;
  189. 'ARGUMENT' trx*'TABLE ' ;
  190. as2 ama1 = 'KOPS' 'MATRIK' ;
  191. trv = (trx . 'EQEX') ;
  192. tinco = (trv . 'INCO') ;
  193. $dom = (trx . 'DOMZ') ;
  194. cdom = 'CONTOUR' ('DOMA' $dom 'MAILLAGE') ;
  195. vitesse = 'KCHT' $dom 'VECT' 'SOMMET' (tinco . 'UN') ;
  196. sw = 'KOPS' vitesse 'ROT' $dom ;
  197. rk = 'EQEX' $dom 'OPTI' 'EF' 'IMPL' 'ZONE' $dom 'OPER' 'LAPN' 1.0D0
  198. 'INCO' 'PSI' 'ZONE' $dom 'OPER' 'FIMP' sw
  199. 'INCO' 'PSI' 'CLIM' 'PSI' 'TIMP' cdom 0.0D0 ;
  200. rk . 'INCO' = 'TABLE' 'INCO' ;
  201. rk . 'INCO' . 'PSI' = 'KCHT' $dom 'SCAL' 'SOMMET' 0.0D0 ;
  202. *-*-*-*- A remplacer par EXEC qd ca marchera...
  203. EXEC rk ;
  204. *-*-*-*-
  205. tinco . 'PSI' = (rk . 'INCO' . 'PSI') ;
  206.  
  207. 'FINPROC' as2 ama1 ;
  208.  
  209. ** MAJHIST
  210. 'DEBPROC' MAJHIST ;
  211. 'ARGUMENT' trx*'TABLE ' ;
  212. as2 ama1 = 'KOPS' 'MATRIK' ;
  213. trv = (trx . 'EQEX') ;
  214. tinco = (trv . 'INCO') ;
  215. th = (trv . 'HISTOIRE') ;
  216. tp = (trv . 'PASDETPS') ;
  217. vit = (tinco . 'UN') ; vitm1 = (tinco . 'UN-1') ;
  218. vitref = (tinco . 'UREF') ;
  219. psi = (tinco . 'PSI') ;
  220. npdt = (tp . 'NUPASDT') ;
  221. * Calcul de l'erreur relative (norme Linfini) entre deux
  222. * champs de vitesse successifs
  223. errrel = ('MAXIMUM' (vit '-' vitm1) 'ABS') '/' vitref ;
  224. 'SI' ('MULT' npdt (trv . 'FIDT')) ;
  225. 'MESSAGE' ('CHAINE' 'Pdt : ' npdt ' |Un - Un_1| Linf = ' errrel) ;
  226. 'FINSI' ;
  227. tinco . 'UN-1' = 'COPIER' vit ;
  228. * Mise à jour des differentes listes
  229. mpsi = 'MAXIMUM' psi 'ABS' ;
  230. th . 'ltps' = ((th . 'ltps') 'ET' ('PROG' (tp . 'TPS'))) ;
  231. th . 'lpsi' = ((th . 'lpsi') 'ET' ('PROG' mpsi)) ;
  232. th . 'lerr' = ((th . 'lerr') 'ET' ('PROG' errrel)) ;
  233. thl = (th . 'lpdt') ;
  234. thl . 'conv' = ((thl . 'conv') 'ET' ('PROG' (tp . 'DTCONV'))) ;
  235. thl . 'diff' = ((thl . 'diff') 'ET' ('PROG' (tp . 'DTDIFU'))) ;
  236. thl . 'min' = ((thl . 'min') 'ET' ('PROG' (tp . 'DELTAT-1'))) ;
  237. 'FINPROC' as2 ama1 ;
  238.  
  239.  
  240. ** TRHIST
  241. 'DEBPROC' TRHIST ;
  242. 'ARGUMENT' th*'TABLE ' ;
  243. * On enleve le premier indice des listes
  244. hltps = 'ENLEVER' (th . 'ltps') 1 ;
  245. hlpsi = 'ENLEVER' (th . 'lpsi') 1 ;
  246. hlerr = 'ENLEVER' (th . 'lerr') 1 ;
  247. thl = (th . 'lpdt') ;
  248. hlpdt = 'ENLEVER' (thl . 'min') 1 ;
  249. hlconv = 'ENLEVER' (thl . 'conv') 1 ;
  250. hldiff = 'ENLEVER' (thl . 'diff') 1 ;
  251. * Pas de temps = fn(t)
  252. evpdt = 'EVOL' 'MANU' 'TEMPS' hltps 'PDT' hlpdt ;
  253. 'DESSIN' evpdt 'MIMA' 'TITR' 'Pas de temps' ;
  254. * LOG10(pas de temps) = fn(t)
  255. ltco = ('LOG' (hlconv '*' alpha)) '/' ('LOG' 10.0D0) ;
  256. ltdi = ('LOG' (hldiff '*' alpha)) '/' ('LOG' 10.0D0) ;
  257. ltmi = ('LOG' hlpdt) '/' ('LOG' 10.0D0) ;
  258. evco = 'COULEUR' ('EVOL' 'MANU' 'TEMPS' hltps 'dtconv' ltco) 'BLEU' ;
  259. evdi = 'COULEUR' ('EVOL' 'MANU' 'TEMPS' hltps 'dtdiff' ltdi) 'ROUG' ;
  260. evmi = 'COULEUR' ('EVOL' 'MANU' 'TEMPS' hltps 'dt' ltmi) 'BLAN' ;
  261. tabt = 'TABLE' ; tabt . 'TITRE' = 'TABLE' ;
  262. tabt . 1 = 'TIRC' ;
  263. tabt . 2 = 'TIRL' ;
  264. tabt . 3 = 'NOLI' ;
  265. tabt . 'TITRE' . 1 = 'CHAINE' 'Pdt convection' ;
  266. tabt . 'TITRE' . 2 = 'CHAINE' 'Pdt diffusion' ;
  267. tabt . 'TITRE' . 3 = 'CHAINE' 'Pdt min(conv, diff)' ;
  268. 'DESSIN' (evco 'ET' evdi 'ET' evmi) 'MIMA' 'TITR' 'Pas de temps' 'LEGE' tabt ;
  269. * Max |Fn de courant| = fn(t)
  270. evpsi = 'EVOL' 'MANU' 'TEMPS' hltps 'PSI' hlpsi ;
  271. 'DESSIN' evpsi 'MIMA' 'TITR' 'Max |Fn de courant|' ;
  272. * Max |Erreur Linf Vit.| = fn(t) (convergence)
  273. everr = 'EVOL' 'MANU' 'TEMPS' hltps 'ERRLINF' (('LOG' hlerr) '/' ('LOG' 10.0D0)) ;
  274. 'DESSIN' everr 'MIMA' 'TITR' 'Log Erreur Linf Vit.' ;
  275. 'FINPROC' ;
  276.  
  277. ** TRCHAMP
  278. 'DEBPROC' TRCHAMP ;
  279. 'ARGUMENT' trx*'TABLE ' ;
  280. as2 ama1 = 'KOPS' 'MATRIK' ;
  281. trv = (trx . 'EQEX') ;
  282. 'SI' ('NON' ('MULT' (trv . 'PASDETPS' . 'NUPASDT') (trv . 'NPTRC'))) ;
  283. 'QUITTER' TRCHAMP ;
  284. 'FINSI' ;
  285. tinco = (trv . 'INCO') ;
  286. $dom = (trx . 'DOMZ') ;
  287. dom = ('DOMA' $dom 'MAILLAGE') ;
  288. cdom = 'CONTOUR' dom ;
  289. tphy = (trv . 'PASDETPS' . 'TPS') ;
  290. * Vitesses
  291. vitesse = (tinco . 'UN') ; vref = (tinco . 'UREF') ;
  292. titv = ('CHAINE' 'Champ de vitesse (VECT SOMMET) à t = ' tphy) ;
  293. vv = 'VECTEUR' vitesse (0.2D0 '/' vref) 'UX' 'UY' 'JAUNE' ;
  294. 'TRACER' vv cdom 'TITR' titv ;
  295. * Fonction de courant
  296. 'OPTION' 'ISOV' 'SULI' ;
  297. psi = (tinco . 'PSI') ;
  298. titpsi = ('CHAINE' 'Fonction de courant (SCAL SOMMET) à t = ' tphy) ;
  299. 'TRACER' psi dom cdom 14 'TITR' titpsi ;
  300. * Pression
  301. * NB : On n'utilise pas 'ELNO' : on dessine une pression
  302. * constante sur un element.
  303. pression = (trv . 'INCO' . 'PRESSION') ;
  304. modelp = 'MODELISER' dom 'THERMIQUE' ;
  305. elemp = 'KCHA' $dom 'CHAM' pression ;
  306. titpres = ('CHAINE' 'Pression (SCAL CENTRE) à t = ' tphy) ;
  307. 'TRACER' elemp modelp cdom 'TITR' titpres ;
  308. 'FINPROC' as2 ama1 ;
  309. *
  310. ** Fin de la definition des procedures
  311. **************************************************************
  312.  
  313. **************************************************************
  314. ** Debut du 'main'
  315. *
  316. **** Definition des constantes
  317. * l, h, nl, nh : donnees pour la boite.
  318. * egeom : erreur geometrique pour 'ELIMINATION'
  319. * aspect : rapport d'aspect de la boite
  320. * nu, g, rho : parametres du fluide
  321. * tg, td : temperatures sur les murs gauches et droits
  322. h = 1.0D0 ;
  323. 'SI' court ; nh = 3 ; 'SINON' ; nh = 5 ; 'FINSI' ;
  324. aspect = 4.0D0 ;
  325. l = h '*' aspect ;
  326. 'SI' court ;
  327. nl = 'ENTIER' (nh '*' 3) ;
  328. 'SINON' ;
  329. nl = 'ENTIER' (nh '*' aspect) ;
  330. 'FINSI' ;
  331. egeom = 1.0D-4 ;
  332. nu = 1.0D0 ; g = (0.D0 -1.D0) ;
  333. tg = 0.D0 ; td = aspect ;
  334.  
  335. **** Creation des maillages
  336. lp = l '/' 2.0D0 ; hp = h '/' 2.0D0 ;
  337. mlp = (-1.0D0) '*' lp ; mhp = (-1.0D0 '*' hp) ;
  338. pmil = lp hp ;
  339. pA = 0.0D0 0.0D0 ; pB = l 0.0D0 ;
  340. pC = l h ; pD = 0.0D0 h ;
  341. pDA = 0.0D0 hp ; pBC = l hp ;
  342. bas = pA 'DROIT' nl pB ; haut = pD 'DROIT' nl pC ;
  343. rwall = pB 'DROIT' nh pC ; lwall = pA 'DROIT' nh pD ;
  344. median = pDA 'DROIT' nl pBC ;
  345. mt = 'DALLER' haut rwall bas lwall ;
  346. mt = chan mt QUAF ;
  347. $mt= MODE mt 'NAVIER_STOKES' DISCR ; DOMA $MT 'IMPR' ;
  348. mt = (DOMA $mt 'MAILLAGE') ; 'TASSER' mt ;
  349. haut = 'CHANGER' haut SEG2 ; bas = 'CHANGER' bas SEG2 ;
  350. lwall = 'CHANGER' lwall SEG2 ; rwall = 'CHANGER' rwall SEG2 ;
  351. median = 'CHANGER' median SEG2 ;
  352. cmt = (lwall 'ET' bas 'ET' rwall 'ET' haut) ;
  353. 'ELIMINATION' (mt 'ET' bas 'ET' rwall 'ET' haut 'ET' lwall 'ET' cmt 'ET' median) egeom ;
  354.  
  355. **** Initialisation de la table des historiques
  356. *
  357. thist = 'TABLE' 'HIST' ;
  358. TABTPS= 'TEMPS' 'NOEC';
  359. thist . 'tpscpu' = TABTPS.'TEMPS_CPU'.'INITIAL' ;
  360. * Listes
  361. thist . 'ltps' = 'PROG' ;
  362. thist . 'lpsi' = 'PROG' ;
  363. thist . 'lerr' = 'PROG' ;
  364. thist . 'lpdt' = 'TABLE' 'LPDT' ; thl = (thist . 'lpdt') ;
  365. thl . 'conv' = 'PROG' ; thl . 'diff' = 'PROG' ;
  366. thl . 'min' = 'PROG' ;
  367.  
  368. **** Initialisation des tables rv et rvp
  369. rv = 'EQEX' $mt 'OPTI' discconv 'ZONE' $mt 'OPER' 'NS' 'NU' 'GB' 'TN' 'TREF' 'INCO' 'UN' 'OPTI' 'CENTREE' 'ZONE' $mt 'OPER' 'DFDT' 1. 'UN' 'DELTAT' 'INCO' 'UN' ;
  370. * Conditions aux limites
  371. rv = 'EQEX' rv 'CLIM' 'UN' 'UIMP' cmt 0.0D0 'UN' 'VIMP' cmt 0.0D0 ;
  372. 'SI' graph ;
  373. rv = 'EQEX' rv 'ZONE' $mt 'OPER' 'CALCPSI' 'ZONE' $mt 'OPER' 'MAJHIST' 'ZONE' $mt 'OPER' 'TRCHAMP' ;
  374. 'FINSI' ;
  375.  
  376. RVP = EQEX 'OPTI' 'EF' KPRESS
  377. 'ZONE' $MT 'OPER' 'KBBT' -1. betap 'INCO' 'UN' 'PRES' ;
  378. pfpres = 'POIN' ('DOMA' $MT KPRESS) 'PROC' pmil ;
  379. ppfpres = 'MANU' 'POI1' pfpres ;
  380. rvp = 'EQEX' rvp 'CLIM' 'PRES' 'TIMP' ppfpres 0.0D0 ;
  381.  
  382.  
  383. rvp.'METHINV'.TYPINV=1 ;
  384. rvp.'METHINV'.IMPINV=0 ;
  385. rvp.'METHINV'.NITMAX=300;
  386. rvp.'METHINV'.PRECOND=3 ;
  387. rvp.'METHINV'.RESID =1.e-8 ;
  388. rvp.'METHINV' . 'FCPRECT'=100 ;
  389. rvp.'METHINV' . 'FCPRECI'=100 ;
  390.  
  391. RV.'PROJ' =RVP ;
  392.  
  393. rv . 'INCO' = 'TABLE' 'INCO' ;
  394. rv . 'HISTOIRE' = thist ;
  395.  
  396. **** Initialisation des inconnues
  397. * V -> 0 ; T gradient constant avec tg sur lwall, td sur rwall
  398. rv . 'INCO' . 'UN-1' = 'KCHT' $mt 'VECT' 'SOMMET' (0.D0 0.D0) ;
  399. rv . 'INCO' . 'UN' = 'KCHT' $mt 'VECT' 'SOMMET' (0.D0 0.D0) ;
  400. rv . 'INCO' . 'PRES' = 'KCHT' $mt 'SCAL' KPRESS 0.D0 ;
  401. rv . 'INCO' . 'TN' = 'KCHT' $mt 'SCAL' 'SOMMET' (tg '+' ((('COORDONNEE' 1 mt) '/' l) '*' (td '-' tg)) ) ;
  402.  
  403. **** Valeurs des parametres :
  404. *
  405. ** Physiques
  406. *
  407. * Gr : Nombre de Grashof
  408. * tref : temperature de reference (approx. de Boussinesq)
  409. * Il faut la dissymetriser pour observer la transition
  410. * 3 -> 2 cellules
  411. * uref : echelle de vitesse (cf. adimensionnement)
  412. * ampvit : facteur d'amplification pour le trace des vitesses
  413. *
  414. 'SI' court ;
  415. Gr = 5.0D3 ; tref = (tg '+' td) '/' 2.0D0 ;
  416. 'SINON' ;
  417. Gr = 40.0D3 ; tref = tg ;
  418. 'FINSI' ;
  419. uref = (Gr '**' 0.5D0) ;
  420. *
  421. ** Numeriques
  422. *
  423. * nbiter : nombre d'iterations
  424. * nptrc : frequence pour le trace des champs (.ps ou à l'ecran)
  425. * fidt : frequence d'affichage des messages à l'ecran
  426. * alpha : multiplicateur du pas de temps
  427. * pour les elements MACRO
  428. 'SI' court ;
  429. nbiter = 250 ; nptrc = 125 ; fidt = 50 ;
  430. 'SINON' ;
  431. nbiter = 15000 ; nptrc = 2500 ; fidt = 100 ;
  432. 'FINSI' ;
  433. alpha = 0.9D0 ;
  434. *
  435. ** Remplissage des tables
  436. *
  437. rv . 'INCO' . 'NU' = nu ;
  438. rv . 'INCO' . 'GB' = g '*' Gr ;
  439. rv . 'INCO' . 'TREF' = tref ;
  440. rv . 'INCO' . 'UREF' = uref ;
  441. rv . 'ALFA' = alpha ;
  442. rv . 'FIDT' = fidt ;
  443. rv . 'NPTRC' = nptrc ;
  444. rv . 'ITMA' = nbiter ;
  445. *
  446. **
  447. ***
  448. **** Execution
  449. SUMMARY ;
  450. 'SI' inta ; 'OPTION' 'DONN' 5 ; 'FINSI' ;
  451. EXEC rv ;
  452. TABTP1= 'TEMPS' 'NOEC';
  453. thist . 'tpscpu' = (thist . 'tpscpu') '+' TABTP1.'TEMPS_CPU'.'INITIAL' ;
  454. 'SI' graph ; TRHIST thist ; 'FINSI' ;
  455. ****
  456. ***
  457. **
  458. *
  459. ** Le cas-test donne-t-il un resultat correct ??
  460. 'SI' court ;
  461. *
  462. ** Verification de la centrosymetrie sur les vitesse
  463. *
  464. vit = (rv . 'INCO' . 'UN') ;
  465. tvit = vit 'TOURNER' 180.0D0 pmil ;
  466. vit = 'KCHT' $mt 'VECT' 'SOMMET' vit ;
  467. 'ELIMINATION' (('EXTRAIRE' tvit 'MAIL') 'ET' mt) egeom ;
  468. tvit = 'KCHT' $mt 'VECT' 'SOMMET' tvit ;
  469. dv = 'KOPS' vit '-' tvit ;
  470. critere1 = ('MAXIMUM' dv 'ABS') '/' uref ;
  471. * Valeur de critere1 au 2/10/97 : 1.718764099177966D-15
  472. limite1 = 1.0D-12 ;
  473. *
  474. ** Evolution vitesse verticale le long de la mediane
  475. *
  476. evvm = 'EVOL' 'CHPO' vit 'UY' median ;
  477. evym = 'EXTRAIRE' evvm 'ORDO' ;
  478. list (evym *1.e-6);
  479. *
  480. * Resultat obtenu le 2/10/97
  481. * correspondant à critere2 = 1.507288760336422D-16
  482. evvref = 'PROG' 0.0D0 -.2205372973655706D+02 -.1537657415971809D+02 -.9886959157032784D+01 -.4994589224716669D+01 -.1325630316341117D+01 .5065207872228873D+00 .9297558065189220D+00 .5908058903982639D+00 .2945842774751989D-13 -.5908058903982122D+00 -.9297558065189337D+00 -.5065207872229142D+00 .1325630316341034D+01 .4994589224716784D+01 .9886959157032821D+01 .1537657415971808D+02 .2205372973655705D+02 0.0D0 ;
  483. * Resultat obtenu le 15/06/2014
  484. evvref = 'PROG' -7.24189E-39 -22.049 -15.375 -9.8821 -4.9918 -1.3199 0.50926 0.93264 0.59143 -1.04084E-12 -0.59143 -0.93264 -0.50926 1.3199 4.9918 9.8821 15.375 22.049 7.24188E-39 ;
  485. critere2 = ('MAXIMUM' (evym '-' evvref) 'ABS') '/' uref ;
  486. limite2 = 1.0D-4 ;
  487. * pour tenir compte de l'affichage de 5 chiffres significatifs
  488. * Tps CPU (le 2/10/97) ; 1009 centiemes de secondes
  489. * (court = VRAI ; graph = FAUX ; inta = FAUX ;)
  490. TABTP1= 'TEMPS' 'NOEC';
  491. thist . 'tpscpu'= (thist . 'tpscpu') '+' TABTP1.'TEMPS_CPU'.'INITIAL';
  492. *
  493. ** Affichage des messages et des erreurs eventuelles
  494. *
  495. 'SAUTER' 2 'LIGNE' ;
  496. 'MESSAGE' 'Debriefing :' ;
  497. 'MESSAGE' '------------' ;
  498. 'SAUTER' 'LIGNE' ;
  499. 'MESSAGE' ('CHAINE' 'Erreur/Limite de centrosymetrie du champ de vitesses :') ;
  500. 'MESSAGE' ('CHAINE' 'FORMAT' '(D24.16)' critere1 ' / ' limite1) ;
  501. 'MESSAGE' ('CHAINE' 'Erreur/Limite Vy sur la mediane horizontale') ;
  502. 'MESSAGE' ('CHAINE' 'FORMAT' '(D24.16)' critere2 ' / ' limite2) ;
  503. 'MESSAGE' ('CHAINE' 'Temps CPU / Temps CPU (30/9/97) : ') ;
  504. 'MESSAGE' ('CHAINE' ('ENTIER' (thist . 'tpscpu')) ' / ~1000' ' centieme(s) de secondes') ;
  505. mess ' ->' critere1 ' 'limite1' 'critere2' 'limite2 ;
  506. 'SI' ('OU' (critere1 '>' limite1) (critere2 '>' limite2)) ;
  507. 'ERREUR' 5 ;
  508. 'FINSI' ;
  509. 'MESSAGE' 'Tout s_est bien passe' ;
  510. 'OPTION' 'ECHO' 1 ;
  511. 'FINSI' ;
  512. SUMMARY ;
  513. 'SI' inta ; 'OPTION' 'DONN' 5 ; 'FINSI' ;
  514. 'FIN' ;
  515. *
  516. ** Fin du 'main'
  517. **************************************************************
  518.  
  519.  
  520.  

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