Télécharger la plaquette Cast3M


fr en

Télécharger cl_E_6.dgibi

Retour à la liste

Numérotation des lignes : 

  1. * fichier :  cl_E_6.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. *C***********************************************************************
  5. *C*                                                                     *
  6. *C* PROJET      : Opérateur CLMI                                        *
  7. *C* NOM         : cl_E_6  .dgibi                                        *
  8. *C* DESCRIPTION : Jeu de données pour le calcul de la couche limite     *
  9. *C*               turbulente en utilisant les relations de fermeture de *
  10. *C*               Cousteix - Cas de la plaque plane avec gradient de    *
  11. *C*               pression positif exponentiel -Méthode de Head         *
  12. *C* LANGAGE     : Gibiane                                               *
  13. *C* AUTEUR      : Guillaume VENCO - DRN/DMT/SEMT/LTMF                   *
  14. *C*                                                                     *
  15. *C***********************************************************************
  16. *C*                                                                     *
  17. *C* APPELES     : Opérateur CLMI                                        *
  18. *C*                                                                     *
  19. *C***********************************************************************
  20. *C*                                                                     *
  21. *C*                                                                     *
  22. *C***********************************************************************
  23. *C*                                                                     *
  24. *C* VERSION     : 19/05/2000                                            *
  25. *C* CREATION    : 07/04/2000                                            *
  26. *C*                                                                     *
  27. *C***********************************************************************
  28. 'OPTION' dime 2 elem seg2 ;
  29. GRAPH=FAUX;
  30. COMPLET=FAUX;
  31. *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
  32. *>>>>>>>>>>>>>Maillage et domaine>>>>>>>>>>>>>>>>>>>>>>>>>>>
  33. *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
  34. P1 = 0.15 0.;
  35. p2 =  10 0. ;
  36. nb1 = 10;
  37. p1p2 = 'DROIT' p1 p2 'DINI' 0.0001 'DFIN' 1. ;
  38. *p1p2 = 'DROIT' p1 p2 100;
  39. mt = p1p2 ;
  40. * 'TRACER'  mt ;
  41.  
  42. *Définition du modèle
  43. Mmt = 'CHANGER'  mt 'QUAF' ;
  44. Mail = 'CHANGER' Mmt 'POI1';
  45. Pt0 = Mail 'ELEM' 1;
  46. $mt = 'MODELISER'  Mmt 'NAVIER_STOKES' 'LINE' ;
  47.  
  48. *Données 
  49. Ue0 = 1.;
  50. NU = 1.d-6 ;
  51. *Valeurs initiales données par la loi de puissance
  52.  
  53. D20=4.335d-4;
  54. Q10=6.37*D20;
  55. *Pas de temps et nombre de pas
  56. DT = 1.d-1;
  57.  
  58. NBIT =10;
  59.  
  60. Si COMPLET ;
  61. NBIT =100;
  62. FINSI;
  63.  
  64. DT2= 1.;
  65.  
  66. NBI2=10;
  67. Si COMPLET ;
  68. NBI2=100;
  69. FINSI;
  70.  
  71. *Champ de vitesse à l'extérieur
  72.  
  73. * vitesse exponentielle
  74. m1= 'DOMA' $mt 'SOMMET';
  75. X = 'COORDONNEE' 1 m1;
  76. X2 = (-2.)*(X-0.15);
  77. X1 = 'EXP' (X2);
  78. UEX = 1+X1;
  79. UE = 'KCHT' $mt scal 'SOMMET' UEX;
  80. DUEX = (-2.)*X1;
  81. DUE = 'KCHT' $mt scal 'SOMMET' DUEX;
  82.  
  83. *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
  84. *>>>>>>>>>>>>Création de la table RV et de la table RK>>>>>>
  85. *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
  86.  
  87. *Dans CLMI le 1er coef correspond au type de fermeture
  88. *le 2ème coef correspond au n°de l'équation à traiter
  89. *1=QDM 2=Energie cinétique
  90.  
  91. rv = 'EQEX'  'OMEGA' 1.  'NITER' 1  'ITMA' 0   'FIDT' 1
  92.      'OPTI' 'EF' 'IMPL' 'CENTREE' 'BDF2'
  93.      'ZONE' $mt
  94.      'OPER' 'DFDT' 1. 'D2NM' 'D2N2' 'DT'  'INCO' 'D2' ;
  95. *     'OPER' 'DFDT' 1. 'D2NM' DT  'INCO' 'D2' ;
  96.  
  97. rv = 'EQEX' rv
  98.      'OPTI' 'EF' 'IMPL' 'SUPG'
  99.       'ZONE'  $mt
  100.      'OPER'  'CLMI' 5 1 UE DUE 'D2NM' 'Q1NM'
  101.      'INCO' 'D2' ;
  102.  
  103. rv = 'EQEX' rv
  104.      'OPTI' 'EF' 'IMPL' 'CENTREE' 'BDF2'
  105.      'ZONE' $mt
  106.      'OPER' 'DFDT' 1. 'Q1NM' 'Q1N2' 'DT'  'INCO' 'Q1' ;
  107. *     'OPER' 'DFDT' 1. 'D1NM' DT  'INCO' 'D1' ;
  108.  
  109. rv = 'EQEX' rv
  110.      'OPTI' 'EF' 'IMPL' 'SUPG'
  111.      'ZONE'  $mt
  112.      'OPER'  'CLMI' 5 3 UE DUE 'Q1NM' 'D2NM'
  113.      'INCO' 'Q1' ;
  114.  
  115. rv = 'EQEX' rv
  116.       'CLIM' 'D2' 'TIMP' pt0 D20
  117.       'CLIM' 'Q1' 'TIMP' pt0 Q10;
  118.  
  119.  
  120. rv.inco = 'TABLE'  inco ;
  121. rv.inco.'D2'=  'KCHT'  $mt scal sommet D20;
  122. rv.inco.'D2NM'=  'KCHT'  $mt scal sommet D20;
  123. rv.inco.'D2N2'=  'KCHT'  $mt scal sommet D20;
  124. rv.inco.'Q1'=  'KCHT'  $mt scal sommet Q10;
  125. rv.inco.'Q1NM'=  'KCHT'  $mt scal sommet Q10;
  126. rv.inco.'Q1N2'=  'KCHT'  $mt scal sommet Q10;
  127. rv.inco.'DT' = dt;
  128. ltps = 'PROG' 0. ;
  129. lteta = 'PROG' 0. ;
  130. TPS = 0. ;
  131. nupt = 0 ;
  132. mts1 = 'DOMA'  $mt 'SOMMET' ;
  133. pt1 = 'POIN'  mts1 2 ;
  134. 'LISTE'  pt1 ;
  135.  
  136. *'OPTION'  donn 5 ;
  137.  
  138. 'REPETER'  BCLT nbit ;
  139. nupt=nupt+1;
  140. tps=tps+dt ;
  141. mess ' Pas de temps DT=' dt ' Temps = ' tps ;
  142. ltps=ltps et (prog tps) ;
  143. EXEC  rv ;
  144. teta='EXTRAIRE'  (rv.inco.'D2') 'SCAL' pt1 ;
  145. mess ' Teta=' teta ;
  146. lteta= lteta et ('PROG' teta);
  147.  
  148.  
  149. YNM = rv.'INCO'.'D2NM';
  150. ZNM = rv.'INCO'.'Q1NM';
  151.  
  152. rv.'INCO'.'D2N2' = rv.'INCO'.'D2NM';
  153. rv.'INCO'.'D2NM' = rv.'INCO'.'D2';
  154. rv.'INCO'.'Q1N2' = rv.'INCO'.'Q1NM';
  155. rv.'INCO'.'Q1NM' = rv.'INCO'.'Q1';
  156.  
  157. 'FIN'  bclt ;
  158.  
  159.  
  160. *******************************************************
  161.  
  162. 'REPETER'  BCLT2 nbi2 ;
  163. rv.inco.'DT'=1.;
  164. DT=1.; 
  165. nupt=nupt+1;
  166. tps=tps+dt2 ;
  167. mess ' Pas de temps DT=' dt ' Temps = ' tps ;
  168. ltps=ltps et (prog tps) ;
  169. EXEC  rv ;
  170. teta='EXTRAIRE'  (rv.inco.'D2') 'SCAL' pt1 ;
  171. mess ' Teta=' teta ;
  172. lteta= lteta et ('PROG' teta);
  173.  
  174.  
  175. YNM = rv.'INCO'.'D2NM';
  176. ZNM = rv.'INCO'.'Q1NM';
  177.  
  178. rv.'INCO'.'D2N2' = rv.'INCO'.'D2NM';
  179. rv.'INCO'.'D2NM' = rv.'INCO'.'D2';
  180. rv.'INCO'.'Q1N2' = rv.'INCO'.'Q1NM';
  181. rv.'INCO'.'Q1NM' = rv.'INCO'.'Q1';
  182.  
  183. 'FIN'  bclt2 ;
  184.  
  185. li4 = 'PROG' 'EXPO' 'A' -1. 'B' 0. ltps;
  186. evt2 = 'EVOL' 'ROUGE' manu ltps li4 ;
  187. evt= evol manu ltps lteta;
  188.  
  189.  
  190.  
  191. *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
  192. *>>>>>>>>>>>>>>>>>>>>>>> résultats>>>>>>>>>>>>>>>>>>>>>>>>>>
  193. *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
  194. *évolution de la vitesse extérieure et du gradient
  195. vit = 'EVOL' 'CHPO' UE mt;
  196. gvit = 'EVOL' 'ROUG' 'CHPO' DUE mt;
  197. *quantité de mouvement
  198. D2S = (rv.inco.'D2');
  199. RD2S = UE*D2S/NU;
  200. ev1s = 'EVOL' 'ROUGE' 'CHPO' RD2S mt;
  201. *différence entre (d-d1)
  202. Q1S = (rv.inco.'Q1');
  203. RQ1S = UE*Q1S/NU;
  204. ev2s = 'EVOL' 'ROUGE' 'CHPO' RQ1S mt;
  205. *Epaisseur de déplacement
  206. D1S = (rv.'INCO'.'D1');
  207. RD1S = UE*D1S/NU;
  208. ev3s = 'EVOL' 'ROUGE' 'CHPO' RD1S mt;
  209. *évolution du frottement à la paroi
  210. CfS = rv.'INCO'.'CF';
  211. ev4s = 'EVOL' 'ROUGE' 'CHPO' CfS mt;
  212. *évolution du facteur de forme H final
  213. HS = rv.'INCO'.'H';
  214. ev5s = 'EVOL' 'ROUGE' 'CHPO' HS mt;
  215. *évolution du facteur H*
  216. HHS = Q1S*(D2S**(-1));
  217. ev6s = 'EVOL' 'ROUGE' 'CHPO' HHS mt;
  218.  
  219. *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
  220. *solution théorique cf:SCHLICHTING >>>>>>>>>>>>>>>>>>>>>>>>
  221. *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
  222. *RX = X*UE*(1./NU);
  223. RX = 1.*X*(1./NU);
  224. *Quantité de mvt
  225. *----Rteta donné par Nikuradze
  226. RD2N = (0.01738/1.304)*(RX**(0.861)); 
  227. ev1n =  'EVOL' 'BLEU' 'CHPO' RD2N mt;
  228. *----Rteta par la loi de puissance
  229. RD2P = (1./1.e-6)*(0.036*X*(RX**(-1./5)))*UE;
  230. ev1p = 'EVOL' 'VERT' 'CHPO' RD2P mt;
  231. *différence (d-d1) théorique
  232. Q1P=(0.37-(1.2857*0.036))*(X*(RX**(-1./5)));
  233. RQ1P=UE*Q1P/NU;
  234. ev2p='EVOL' 'VERT' 'CHPO' RQ1P mt;
  235. *Epaisseur de déplacement
  236. *----Rdelta1 par la Nikuradze
  237. RD1N = 1.304*RD2N;
  238. ev3n = 'EVOL' 'BLEU' 'CHPO' RD1N mt;
  239. *----Rdelta1 par la loi de puissance
  240. RD1P = (9./7.)*RD2P;
  241. ev3p = 'EVOL' 'VERT' 'CHPO' RD1P mt;
  242. *Coefficient de frottement
  243. *----Formule par Nikuradze
  244. CfN = 0.02296*(RX**(-0.139));
  245. ev4n = 'EVOL' 'BLEU' 'CHPO' CfN mt;
  246. *----Formule par loi de puissance
  247. CfP = 0.0576*(RX**(-1./5.));
  248. ev4p = 'EVOL' 'VERT' 'CHPO' CfP mt;
  249. *facteur de forme
  250. *----H par Nikuradze
  251. HN = 'KCHT' $mt scal sommet 1.304;
  252. ev5n = 'EVOL' 'BLEU' 'CHPO' HN mt;
  253. *----H par loi de puissance
  254. HP = 'KCHT' $mt scal sommet 1.2857;
  255. ev5p = 'EVOL' 'VERT' 'CHPO' HP mt;
  256. *coef H* théorique
  257. HHP='KCHT' $mt scal sommet 7.89;
  258. ev6p= 'EVOL' 'BLEU' 'CHPO' HHP mt;
  259.  
  260. *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
  261. *Courbes>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
  262. *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
  263. Si GRAPH ;
  264.  
  265. 'DESSIN'  (evt) 'TITR' 'Convergence de la valeur au dernier noeud'
  266.            'TITY' 'Valeurs de Y' 'TITX' 'Temps';
  267.  
  268. 'DESSIN' (vit 'ET' gvit)
  269.      'TITRE' 'Vitesse et gradient de vitesse extérieure';
  270.  
  271. TAB1 = 'TABLE' ;
  272. TAB1.'TITRE' = 'TABLE' ;
  273. TAB1.'TITRE' . 1 = MOT ' Simulation';
  274. TAB1.'TITRE' . 2 = MOT ' Nikuradze';
  275. TAB1.'TITRE' . 3 = MOT ' Loi de puissance';
  276.  
  277. TAB2 = 'TABLE' ;
  278. TAB2.'TITRE' = 'TABLE' ;
  279. TAB2.'TITRE' . 1 = MOT ' Simulation';
  280. TAB2.'TITRE' . 2 = MOT ' Loi de puissance';
  281.  
  282. 'DESSIN' (ev1s 'ET' ev1n 'ET' ev1p)
  283.        'TITRE' 'Epaisseur de quantité de mouvement'
  284.        'TITY' 'Rteta'
  285.        'TITX' 'Abscisses'
  286.        'GRIL'
  287.        'LEGE' TAB1;       
  288.  
  289. 'DESSIN' (ev2s 'ET' ev2p) 'TITRE' 'différence entre delta et delta1' 
  290.        'TITY' 'Rdelta-Rdelta1'
  291.        'TITX' 'Abscisses'
  292.        'GRIL'
  293.        'LEGE' TAB2;              
  294.  
  295. 'DESSIN' (ev3s 'ET' ev3n 'ET' ev3p)
  296.        'TITRE' 'Epaisseur de déplacement'
  297.        'TITY' 'Rdelta1'
  298.        'TITX' 'Abscisses'
  299.        'GRIL'
  300.        'LEGE' TAB1;
  301.  
  302. 'DESSIN' (ev4s 'ET' ev4n 'ET' ev4p)
  303.            'TITRE' 'Coefficient de frottement à la paroi'
  304.            'TITX' 'Abscisses'
  305.            'TITY' 'Cf'
  306.            'GRIL'
  307.            'LEGE' TAB1;
  308.  
  309. 'DESSIN' (ev5s 'ET' ev5n 'ET' ev5p)
  310.         'TITRE' 'Evolution du facteur de forme'
  311.         'TITX' 'Abscisses'
  312.         'TITY' 'H=delta1/theta'
  313.         'LEGE' TAB1;
  314.  
  315. 'DESSIN' (ev6s 'ET' ev6p) 'TITRE' 'Evolution du facteur H*'
  316.         'TITX' 'Abscisses'
  317.         'TITY' 'H'
  318.         'LEGE' TAB1;
  319. Finsi ;
  320.  
  321. *calcul de l'erreur au sens de la norme L2
  322. diag = 'DOMA' $mt 'XXDIAGSI';
  323. r = (RD2S-RD2N)**2;
  324. r2= (RD2S-RD2P)**2;
  325. EL2N = 'SOMT'(diag*r);
  326. EL2P = 'SOMT'(diag*r2);
  327. 'MESSAGE' 'Erreur sur teta au sens de la norme L2=';
  328. 'LISTE'  EL2N;
  329. 'LISTE'  EL2P;
  330. r3=(RQ1S-RQ1P)**2;
  331. EL3P= 'SOMT'(diag*r3);
  332. 'MESSAGE' 'Erreur sur (delta-delta1) au sens de la norme L2=';
  333. 'LISTE'  EL3P;
  334.  
  335. *SAUVEGARDE
  336. *'OPTION' 'SAUVER' '/u2/ttmf2/guillaume/c20/turb6H_sauv';
  337. *SAUVER;
  338. **************************************************************************
  339. **************************************************************************
  340. *Facteur de forme
  341. x = 'EXTRAIRE' ev5s 'ABSC';
  342. hs = 'EXTRAIRE' ev5s 'ORDO';
  343. hn = 'EXTRAIRE' ev5n 'ORDO';
  344. hp = 'EXTRAIRE' ev5p 'ORDO';
  345. nlist = 'DIME' x;
  346. *'OPTION' echo 0;
  347. *'OPTION' impr 'turb6H_H.data';
  348. 'REPETER' BLIST nlist;
  349.    xi = 'EXTRAIRE' x &BLIST;
  350.    hsi = 'EXTRAIRE' hs &BLIST;
  351.    hni = 'EXTRAIRE' hn &BLIST;
  352.    hpi = 'EXTRAIRE' hp &BLIST;
  353. *  'MESSAGE' ('CHAINE' xi ' ' hsi ' ' hni ' ' hpi);
  354. 'FIN' BLIST;
  355. *'OPTION' impr 'Bidon';
  356. *Epaisseur de quantité de mouvement
  357. x = 'EXTRAIRE' ev1s 'ABSC';
  358. d2s = 'EXTRAIRE' ev1s 'ORDO';
  359. d2n = 'EXTRAIRE' ev1n 'ORDO';
  360. d2p = 'EXTRAIRE' ev1p 'ORDO';
  361. nlist = 'DIME' x;
  362. *'OPTION' echo 0;
  363. *'OPTION' impr 'turb6H_D2.data';
  364. 'REPETER' BLIST nlist;
  365.    xi = 'EXTRAIRE' x &BLIST;
  366.    d2si = 'EXTRAIRE' d2s &BLIST;
  367.    d2ni = 'EXTRAIRE' d2n &BLIST;
  368.    d2pi = 'EXTRAIRE' d2p &BLIST;
  369. *  'MESSAGE' ('CHAINE' xi ' ' d2si ' ' d2ni ' ' d2pi);
  370. 'FIN' BLIST;
  371. *'OPTION' impr 'Bidon';
  372. *Epaisseur de déplacement
  373. x = 'EXTRAIRE' ev3s 'ABSC';
  374. d1s = 'EXTRAIRE' ev3s 'ORDO';
  375. d1n = 'EXTRAIRE' ev3n 'ORDO';
  376. d1p = 'EXTRAIRE' ev3p 'ORDO';
  377. nlist = 'DIME' x;
  378. *'OPTION' echo 0;
  379. *'OPTION' impr 'turb6H_D1.data';
  380. 'REPETER' BLIST nlist;
  381.    xi = 'EXTRAIRE' x &BLIST;
  382.    d1si = 'EXTRAIRE' d1s &BLIST;
  383.    d1ni = 'EXTRAIRE' d1n &BLIST;
  384.    d1pi = 'EXTRAIRE' d1p &BLIST;
  385. *  'MESSAGE' ('CHAINE' xi ' ' d1si ' ' d1ni ' ' d1pi);
  386. 'FIN' BLIST;
  387. *'OPTION' impr 'Bidon';
  388.  
  389. *Coeff de frottement
  390. x = 'EXTRAIRE' ev4s 'ABSC';
  391. cfs = 'EXTRAIRE' ev4s 'ORDO';
  392. cfn = 'EXTRAIRE' ev4n 'ORDO';
  393. cfp = 'EXTRAIRE' ev4p 'ORDO';
  394. nlist = 'DIME' x;
  395. *'OPTION' echo 0;
  396. *'OPTION' impr 'turb6H_CF.data';
  397. 'REPETER' BLIST nlist;
  398.    xi = 'EXTRAIRE' x &BLIST;
  399.    cfsi = 'EXTRAIRE' cfs &BLIST;
  400.    cfni = 'EXTRAIRE' cfn &BLIST;
  401.    cfpi = 'EXTRAIRE' cfp &BLIST;
  402. *  'MESSAGE' ('CHAINE' xi ' ' cfsi ' ' cfni ' ' cfpi);
  403. 'FIN' BLIST;
  404. *'OPTION' impr 'Bidon';
  405.  
  406. 'FIN' ;
  407.  
  408.  
  409.  
  410.  
  411.  
  412.  
  413.  
  414.  
  415.  
  416.  
  417.  
  418.  
  419.  
  420.  
  421.  
  422.  
  423.  
  424.  

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