Télécharger la plaquette Cast3M


fr en

Télécharger cl_D_6.dgibi

Retour à la liste

Numérotation des lignes : 

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

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