Télécharger la plaquette Cast3M


fr en

Télécharger panachekei.dgibi

Retour à la liste

Numérotation des lignes : 

  1. * fichier : panachekei.dgibi
  2. *********************************************************************
  3. *    JET/PANACHE  2D - SEMI-INFINI                                  *
  4. *    Comparaison K - Epsilon / RNG K - Epsilon                      *
  5. *     il suffit pour cela d'enlever ou remettre                     *
  6. *     l'option RNG                                                  *
  7. *     Mode axisymetrique                                            *
  8. * Le Richardson est parametrable Ri = -5.e-2                        *
  9. *     Gregory Turbelin 29/12/1998                                   *
  10. *                                                                   *
  11. *     jpm, mars 2006, adaptation pour K-epsilon implicite (kei)     *
  12. *     En implicite il est bon de contrôler les réentrées de fluide  *
  13. *     procédure outflow                                             *
  14. *                                                                   *
  15. *********************************************************************
  16.  
  17.  OPTI DIME 2 ELEM QUA8;
  18. *OPTI MODE AXIS; OPTI ISOV SULI;
  19.  OPTI MODE AXIS;
  20. *OPTI TRAC 'PSC' ;
  21.  
  22. COMPLET = VRAI;
  23. COMPLET = FAUX;
  24. graph = faux;
  25. *graph = vrai;
  26. *opti trace 'PSC' ;
  27.  
  28. *--------------------------- Numérique ---------------------------------
  29. Si COMPLET;
  30. NBIT=400;
  31. DT = 1. ;
  32. DISCR = MACRO;
  33. KPRES= CENTREP1;
  34. sinon;
  35. NBIT=40;
  36. DT = 5.   ;
  37. DISCR = LINE;
  38. DISCR = MACRO;
  39. KPRES= CENTRE;
  40. *KPRES= CENTREP1;
  41. finsi;
  42.  
  43.  
  44.  
  45. *ALGOKEPS= mots 'Nut' 'Filtres' ;
  46.  ALGOKEPS= mots 'Nut'  ;
  47. *--------------------------- procédure ---------------------------------
  48. *$$$$ OUTFLOW
  49. 'DEBPROC' OUTFLOW  ;
  50. ARGU RX*TABLE ;
  51. ************************************************************************
  52. * ZONE $Mt OPER 'OUTFLOW' $out Ro Un Muf INCO 'UN'
  53. *
  54. *         UN
  55. *         MUF       viscosité dynamique effective
  56. *                             (tourbillonnaire+moléculaire)
  57. ************************************************************************
  58. rv=rx.'EQEX' ;
  59. iarg=rx.'IARG' ;
  60. *mess 'iarg='iarg;
  61. $mt=rx.'DOMZ' ;
  62.  
  63. * Lecture du 1er Argument objet mmodel
  64. Si(ega ('TYPE' rx.'ARG1') 'MMODEL  ');
  65. $mfront=rx.'ARG1';
  66. Dg=doma $mfront 'XXDIAGSI' ;
  67. Sinon ;
  68. erreur 2;
  69. quitter outflow;
  70. Finsi ;
  71.  
  72. * Lecture du 2ème Argument la densité
  73. Si(ega ('TYPE' rx.'ARG2') 'MOT     ');
  74. Ro = rv.inco.(rx.'ARG2');
  75. Sinon ;
  76. Ro = rx.'ARG2';
  77. Finsi ;
  78.  
  79. * Lecture du 3ème Argument la vitesse
  80. Si(ega ('TYPE' rx.'ARG3') 'MOT     ');
  81. un = rv.inco.(rx.'ARG3');
  82. Sinon ;
  83. un = rx.'ARG3';
  84. Finsi ;
  85.  
  86. * Lecture du 4ème Argument la viscosité
  87. Si(ega ('TYPE' rx.'ARG4') 'MOT     ');
  88. Muf = rv.inco.(rx.'ARG4');
  89. Sinon ;
  90. Muf = rx.'ARG4';
  91. Finsi ;
  92.  
  93. **************** INITIALISATIONs ***********************
  94. Si (Exist rx 'rxtm');
  95. rxtm=rx.'rxtm';
  96. lmx=rxtm.'lmx';
  97. lm0=rxtm.'lm0';
  98. lm=rxtm.'lm';
  99. nc=rxtm.'nc';
  100. MSV = chai rxtm.'msv' ;
  101. MI1 = rxtm.'MI1' ;
  102. KPRES=rxtm.'KPRES';
  103. Sinon;
  104.  
  105. rxtm=table 'KIZX';
  106. rx.'rxtm'=rxtm;
  107. Si (NON (EGA (dime rx.'LISTINCO') 1));
  108.   mess ' Erreur dans la procédure Outflow !';
  109.   mess ' Il doit y avoir une inconnue !';
  110.   quitter outflow ;
  111. Finsi;
  112. MI1=EXTR rx.'LISTINCO' 1 ;
  113. lmx=(extr (rv.inco.MI1) 'COMP');
  114. lm0=(extr un 'COMP');
  115. nc=dime lm0;
  116.  
  117. Si(nc > 1);
  118. *mess 'Outflow Cas vectoriel ' MI1;
  119. MSV=chai 'VECT';
  120. lm=mots (chai 1 MI1) (chai 2 MI1);
  121. Si (EGA (vale dime) 3);
  122. lm=lm et (mots (chai 3 MI1));
  123. Finsi;
  124. kpr=rx.'KOPT'.'KPOIN';
  125. Si(EGA KPR 2);KPRES='CENTRE';Finsi;
  126. Si(EGA KPR 4);KPRES='CENTREP1';Finsi;
  127. Si(EGA KPR 5);KPRES='MSOMMET';Finsi;
  128. Sinon;
  129. *mess 'Outflow Cas scalaire  ' MI1;
  130. MSV=chai 'SCAL';
  131. lm=mots MI1;
  132. Finsi;
  133. rxtm.'lmx'=lmx;
  134. rxtm.'lm0'=lm0;
  135. rxtm.'lm'=lm;
  136. rxtm.'nc'=nc;
  137. rxtm.'msv'=chai MSV;
  138. rxtm.'MI1'=MI1;
  139. rxtm.'KPRES'=KPRES;
  140.  
  141. rxtm.'LISTINCO'=MOTS MI1 ;
  142. rxtm.'DOMZ'=$mfront;
  143. rxtm.'KOPT'=rx.'KOPT';
  144. rxtm.'IARG'=1;
  145. rxtm.'EQEX'=rx.'EQEX';
  146. rxtm.'NOMZONE'=' ';
  147. rxtm.'TDOMZ'=0;
  148. rxtm.'NOMOPER'=mot 'MDIA_T';
  149. Finsi;
  150. **************** Fin Initialisations *******************
  151.  
  152.  mfront=doma $mfront maillage ;
  153.  nj=doma $mt 'NORMALEV';
  154. *unj= vect nj 0.1 ux uy jaune ;
  155. *trace unj mt;
  156.  njf = redu nj mfront ;
  157. *unj= vect njf 0.1 ux uy rouge ;
  158. *trace unj mt TITR ' FRONTIERE';
  159.  xn = rv.inco.MI1;
  160.  ujf = redu un mfront ;
  161.  us=psca ujf lm0 njf lm0;
  162.  ius=masq us 'INFERIEUR' 0.;
  163.  rv.inco.'us'=us;
  164.  rv.inco.'ius'=ius;
  165.  Si(Ega nc 1);
  166.  us=2.*(redu xn mfront);
  167.  Finsi;
  168.  rxtm.'ARG1'=kcht $mfront scal sommet ((-1.)*Ro*us*ius);
  169.  
  170.  st1 mat1 = MDIA rxtm ;
  171. *st1 mat1 = KOPS 'MATRIK' ;
  172.  
  173.  Dgi=Dg*ius;
  174.  Si(nc > 1);
  175.  puj=0.;
  176.  Si(EXIST (rv.inco) 'PRESSION');
  177.  pn=redu (elno $mt rv.inco.'PRESSION' KPRES) mfront;
  178.  puj=Ro*pn*us;
  179.  Finsi;
  180.  gru = redu (Muf*(kops 'GRADS' (exco un 'UX') $mt)) mfront;
  181.  grv = redu (Muf*(kops 'GRADS' (exco un 'UY') $mt)) mfront;
  182.  mgunj = (nomc (Dgi*((psca gru lm0 njf lm0)-puj)) (extr lm 1)) + (nomc (Dgi*((psca grv lm0 njf lm0)-puj)) (extr lm 2)) ;
  183.  Si(Ega (vale dime) 3);
  184.  mgunj = mgunj + (nomc (Dgi*((psca grv lm0 njf lm0)-puj)) (extr lm 3));
  185.  Finsi;
  186.  Sinon;
  187.  grt = redu (Muf*(kops 'GRADS' xn $mt)) mfront ;
  188.  mgunj = (nomc (Dgi*(psca grt lm0 njf lm0)) (extr lm 1)) ;
  189.  Finsi;
  190.  St1 = kcht $mfront MSV sommet comp lm ((-1.)*mgunj);
  191.  
  192. RESPRO St1 mat1 ;
  193. FINPROC ;
  194.  
  195. ** PROCEDURE DE FILTRE DE LA CONCENTRATION
  196.  
  197.  DEBP 4FILTRE ; ARGU RX*TABLE ;
  198.  rv=rx.'EQEX' ;
  199.  cn=rv.inco.'cn';
  200.  cn=KOPS cn '|<' 0.;
  201.  cn=KOPS cn '>|' 1.;
  202.  rv.inco.'cn'=cn ;
  203. as2 ama1 = 'KOPS' 'MATRIK';
  204. 'RESPRO' as2 ama1;
  205.  FINPROC ;
  206. *----------------------- fin procédures --------------------------------
  207.  
  208. ****************************
  209. * CONSTRUCTION DU MAILLAGE *
  210. ****************************
  211.  DJ=1. ;RJ=DJ/2. ;
  212.  
  213. ** POINTS
  214.  P0=0.  0.; P1=RJ  0.;
  215.  P2=8.  0.; P3=8. 15.;
  216.  P4=RJ 15.; P5=0. 15.;
  217.  
  218.  P1b=4. 0. ; P2b=8. 8.;
  219.  P3b=4. 15.; P5b=0. 8.;
  220.  P6=8.   5.; P7=8.  10.;
  221.  P9=0.  10.; P10=0.  5.;
  222.  
  223.  Pe1=20. 0. ; Pe2=20. 8.;
  224.  Pe3=20. 15.; Pe4=20. 50.;
  225.  Pe5=8. 50. ; Pe6=4. 50.;
  226.  Pe7=RJ 50. ; Pe8=0. 50. ;
  227.  
  228. ***** SEGMENTS *****
  229. Si complet ;
  230. N1=5 ; N2=-18 ; N3=-10 ; N4=-12 ; N5=-30 ; N6=-38 ;
  231. N1=3 ; N2=-12; N3=-6 ; N4=-9 ; N5=-21 ; N6=-27;
  232. sinon ;
  233. N1=1 ; N2=-4 ; N3=-2 ; N4=-3 ; N5=-7  ; N6=-9 ;
  234. finsi ;
  235.  
  236.  ENT=P0 DROI N1   P1;
  237.  BA1=P1 DROI N2  P1b DINI (0.008*DJ)  DFIN (0.1*DJ);
  238.  BA2=P1b DROI N3 P2 DINI (1.1*DJ)  DFIN (1.4*DJ);
  239.  BA3=P2 DROI N4 Pe1 DINI (1.6*DJ) DFIN (3.*DJ);
  240.  BAS=BA1 ET BA2 ET BA3;
  241.  BAT=ENT ET BAS;
  242.  
  243.  CO1=Pe1 DROI N5  Pe2 DINI (0.008*DJ)  DFIN (0.1*DJ);
  244.  CO2=Pe2 DROI N2  Pe3 DINI (1.35*DJ)  DFIN (1.5*DJ);
  245.  CO3=Pe3 DROI N6  Pe4 DINI (1.8*DJ) DFIN (3.*DJ);
  246.  COT=CO1 ET CO2 ET CO3;
  247.  
  248.  HT1=Pe4 DROI N4 Pe5 DINI (3.*DJ) DFIN (1.6*DJ);
  249.  HT2=Pe5 DROI N3 Pe6 DINI (1.4*DJ) DFIN (1.1*DJ);
  250.  HT3=Pe6 DROI N2 Pe7 DINI (0.1*DJ) DFIN (0.008*DJ);
  251.  HT4=Pe7 DROI N1 Pe8;
  252.  HTT=HT1 ET HT2 ET HT3 ET HT4;
  253.  
  254.  AX1=Pe8 DROI N6  P5 DINI (3.*DJ) DFIN (1.8*DJ);
  255.  AX2=P5 DROI N2  P5b DINI (1.5*DJ) DFIN (1.35*DJ);
  256.  AX3=P5b DROI N5  P0 DINI (0.1*DJ) DFIN (0.008*DJ);
  257.  AXE=AX1 ET AX2 ET AX3;
  258.  
  259.  ELIM (BAT ET COT ET HTT ET AXE) 1.e-5;
  260.  MT=DALL BAT COT HTT AXE PLAN;
  261.  ORIENTER MT;
  262.  
  263.  Mmt= changer QUAF mt ;
  264.  Mcot = chan QUAF cot ;
  265.  Maxe = chan QUAF axe ;
  266.  Mbat = chan QUAF bat ;
  267.  Mhtt = chan QUAF htt ;
  268.  Ment = chan QUAF ent ;
  269.  Msor = chan QUAF (cot et htt);
  270.  
  271.  Elim (Mmt et Mcot et Maxe et Mbat et Mhtt et Msor) 1.e-5 ;
  272.  
  273. ***** Modèles Navier_Stokes *********
  274.  $MT =mode MMT  'NAVIER_STOKES' DISCR ;
  275.  $COT=mode MCOT 'NAVIER_STOKES' DISCR ;
  276.  $AXE=mode MAXE 'NAVIER_STOKES' DISCR ;
  277.  $BAT=mode MBAT 'NAVIER_STOKES' DISCR ;
  278.  $HTT=mode MHTT 'NAVIER_STOKES' DISCR ;
  279.  $ENT=mode MENT 'NAVIER_STOKES' DISCR ;
  280.  $SOR=mode MSOR 'NAVIER_STOKES' DISCR ;
  281.  
  282.  MT=doma $MT maillage;
  283.  cot = doma $cot maillage ;
  284.  axe = doma $axe maillage ;
  285.  bat = doma $bat maillage ;
  286.  htt = doma $htt maillage ;
  287.  ent = doma $ent maillage ;
  288.  sorp= doma $sor 'MSOMMET';
  289.  
  290. si graph;  TRAC MT; finsi;
  291.  
  292. **********************
  293. * DONNEES PHYSIQUES  *
  294. **********************
  295.  rhoa=2*1.19; rhof=2*0.58;
  296.  mua=1.8e-5 ; nua=mua/rhoa;
  297.  cref=0.    ; Ri=0. -5.e-2;
  298.  cref=0.    ; Ri=0. -1.   ;
  299.  Re=3000.   ; iRe=1./Re;
  300.  db=1.e-5   ; iSc=db/nua;
  301.  RS=iRe*iSc ;  Sct=0.7;
  302.  KA=2.e-5   ; EA=1.E-6;
  303.  
  304. ***************
  305. * EQUATIONS   *
  306. ***************
  307.   rv=EQEX $MT ITMA  NBIT   ALFA 1.  OPTI EF 'SUPG' 'IMPL' 'ZONE' $MT OPER OUTFLOW $sor 1. 'un' 'MUF'  'INCO' 'un' ;
  308.   rv = eqex rv OPTI EF 'SUPG' 'IMPL' 'ZONE' $mt 'OPER' 'KEPSILON' 1. 'un' iRe DT Ri 'cn' 'INCO' 'kn' 'en' ZONE $MT OPER NS 1. 'un' 'MUF' Ri 'cn' cref INCO 'un' ZONE $MT OPER 4FILTRE 'cn' ZONE $MT OPER TSCA 1. 'un' 'MUF' 0.  INCO 'cn' ;
  309.  
  310.   rv = eqex rv OPTI EFM1 'CENTREE' ZONE $MT OPER DFDT 1. 'un' DT 'un' 'MUF' INCO 'un' ZONE $MT OPER DFDT 1. 'cn' DT 'un' 'MUF' INCO 'cn' ;
  311.  
  312.  rv.nomvi='un';
  313.   rv=eqex rv
  314.     CLIM 'un' UIMP ENT 0.
  315.          'un' VIMP ENT 1.
  316.          'un' VIMP BAS 0.
  317.          'un' UIMP AXE 0.
  318.          'un' VIMP COT 0.
  319.          'kn' TIMP ENT KA
  320.          'en' TIMP ENT EA ;
  321.   rv=eqex rv
  322.     CLIM 'cn' TIMP ENT 1.
  323.          'cn' TIMP COT 0.
  324.          'cn' TIMP BAS 0.
  325.          'kn' TIMP COT 1.e-5
  326.          'en' TIMP COT 1.e-6;
  327.  
  328.     RVP= EQEX 'OPTI' 'EF' KPRES 'ZONE' $mt 'OPER' 'KBBT'  -1.  'INCO' 'un' 'PRES'
  329. *  'ZONE' $paroi 'OPER' 'VNIMP' $vc 0.   'INCO' 'UN' 'PRES'
  330.    ;
  331.    Si (EGA KPRES 'MSOMMET');
  332.    RVP= EQEX RVP 'CLIM' 'PRES' TIMP sorp 0.;
  333.    Finsi;
  334.  
  335.     RV.'PROJ'= RVP ;
  336.  
  337.     RV.'ALGO_KEPSILON'=ALGOKEPS;
  338.  
  339. ** INITIALISATIONS
  340.  rv.inco=TABLE INCO;
  341.  rv.inco.'un'=KCHT $MT VECT SOMMET (1.e-7 1.e-7);
  342.  rv.inco.'kn'=KCHT $MT SCAL SOMMET 1.e-3;
  343.  rv.inco.'en'=KCHT $MT SCAL SOMMET 1.e-4;
  344.  rv.inco.'MUF'=KCHT $MT SCAL SOMMET nua;
  345.  rv.inco.'ALFT'=KCHT $MT SCAL CENTRE nua;
  346.  rv.inco.'cn'=KCHT $MT SCAL SOMMET 0.;
  347.  rv.inco.'PRES'=kcht $mt scal KPRES 0. ;
  348.  
  349.  lh=   (POIN MT PROC (0.75 1.))  et (POIN MT PROC (0.75 5.)) et (POIN MT PROC (0.75 10.)) et (POIN MT PROC (0.75 15.))
  350.     et (POIN MT PROC (0.75 20.)) et (POIN MT PROC (0.75 25.)) et (POIN MT PROC (0.75 33.)) et (POIN MT PROC (0.75 40.))
  351.     et (POIN MT PROC (0.75 45.)) et (POIN MT PROC (0.75 50.)); 
  352.   his=khis 'un' 1 lh 'un' 2 lh 'kn'   lh 'cn'   lh 'en'   lh; his.'KFIH'=1;
  353.   rv.'HIST'=his;
  354.  
  355.  EXEC rv;
  356.  
  357.  
  358. ** POST-TRAITEMENT
  359. **----------------
  360. ** HISTORIQUES
  361. **------------
  362. TAB1=TABLE;
  363. TAB1.'TITRE'= table ;
  364. TAB1.1 ='TIRR MARQ CROI';
  365. TAB1.'TITRE' . 1 ='P1';
  366. TAB1.2 ='TIRR MARQ TRIA';
  367. TAB1.'TITRE' . 2 ='P2';
  368. TAB1.3 ='TIRR MARQ CARR';
  369. TAB1.'TITRE' . 3 ='P3';
  370. TAB1.4 ='TIRR MARQ ETOI';
  371. TAB1.'TITRE' . 4 ='P4';
  372. TAB1.5 ='TIRR MARQ PLUS';
  373. TAB1.'TITRE' . 5 ='P5';
  374. TAB1.6 ='TIRR MARQ LOSA';
  375. TAB1.'TITRE' . 6 ='P6';
  376. TAB1.7 ='TIRR MARQ TRIB';
  377. TAB1.'TITRE' . 7 ='P7';
  378.  
  379.  si graph;
  380.  DESS (his.'1un') LEGE  TAB1 ;
  381.  DESS (his.'2un') LEGE  TAB1 ;
  382.  DESS (his.'kn') LEGE TAB1;
  383.  DESS (his.'cn') LEGE TAB1;
  384.  DESS (his.'en') LEGE TAB1;
  385. *dessin his.'TABD' his.'kn' ;
  386. *dessin his.'TABD' his.'en' ;
  387. *dessin his.'TABD' his.'1un' ;
  388. *dessin his.'TABD' his.'2un' ;
  389.  
  390.  finsi ;
  391.  
  392. ** PROFILS SUR LE DOMAINE TOTAL
  393. **------------------------------
  394. vv = vect rv.inco.'un' 5. ux uy jaune ;
  395. kk =  rv.inco.'kn' ;
  396. ee =  rv.inco.'en';
  397. cc =  rv.inco.'cn' ;
  398. nn = rv.inco.'MUF' * (1./iRe);
  399.  
  400.  si graph;
  401. titre 'Vecteur Vitesse';
  402. trac vv MT (CONT MT); ;
  403. titre 'PRODT';
  404. trac (rv.inco.'PRODT') MT (CONT MT);
  405. titre 'Ksi';
  406. trac (rv.inco.'Ksi') MT (CONT MT);
  407. titre 'Energie Cinetique Turbulente';
  408. trac kk MT (CONT MT);
  409. titre 'MUF /mu';
  410. trac nn MT (CONT MT);
  411. titre 'dissipation';
  412. trac ee MT (CONT MT);
  413.  finsi ;
  414.  
  415. ** PROFILS A DIFFERENTES ALTITUDES
  416. **---------------------------------
  417.  LMT=CHAN LIGNE MT;
  418. ** ALTITUDES
  419.  alt1 = 1. ; alt2 = 5.;
  420.  alt3 = 10.; alt4 = 25.;
  421.  
  422. ** LIGNES HORIZONTALES
  423.  lign1=POIN MT DROI (POIN MT PROC (0. alt1)) (POIN MT PROC (20. alt1)) 1.e-3;
  424.  lign2=POIN MT DROI (POIN MT PROC (0. alt2)) (POIN MT PROC (20. alt2)) 1.e-3;
  425.  lign3=POIN MT DROI (POIN MT PROC (0. alt3)) (POIN MT PROC (20. alt3)) 1.e-3;
  426.  lign4=POIN MT DROI (POIN MT PROC (0. alt4)) (POIN MT PROC (20. alt4)) 1.e-3;
  427.  
  428.  L1=ELEM LMT APPUYE lign1;
  429.  L2=ELEM LMT APPUYE lign2;
  430.  L3=ELEM LMT APPUYE lign3;
  431.  L4=ELEM LMT APPUYE lign4;
  432.  
  433.  nut=rv.inco.'MUF' ;
  434.  
  435.  prvy1=EVOL CHPO (rv.inco.'un') uy L1;
  436.  prvx1=EVOL CHPO (rv.inco.'un') ux L1;
  437.  prc1=EVOL CHPO (rv.inco.'cn') scal L1;
  438.  prk1=EVOL CHPO (rv.inco.'kn') scal L1;
  439.  prnut1=EVOL CHPO nut scal L1;
  440.  
  441.  prvy2=EVOL CHPO (rv.inco.'un') uy L2;
  442.  prvx2=EVOL CHPO (rv.inco.'un') ux L2;
  443.  prc2=EVOL CHPO (rv.inco.'cn') scal L2;
  444.  prk2=EVOL CHPO (rv.inco.'kn') scal L2;
  445.  prnut2=EVOL CHPO nut scal L2;
  446.  
  447.  prvy3=EVOL CHPO (rv.inco.'un') uy L3;
  448.  prvx3=EVOL CHPO (rv.inco.'un') ux L3;
  449.  prc3=EVOL CHPO (rv.inco.'cn') scal L3;
  450.  prk3=EVOL CHPO (rv.inco.'kn') scal L3;
  451.  prnut3=EVOL CHPO nut scal L3;
  452.  
  453.  prvy4=EVOL CHPO (rv.inco.'un') uy L4;
  454.  prvx4=EVOL CHPO (rv.inco.'un') ux L4;
  455.  prc4=EVOL CHPO (rv.inco.'cn') scal L4;
  456.  prk4=EVOL CHPO (rv.inco.'kn') scal L4;
  457.  prnut4=EVOL CHPO nut scal L4;
  458.  
  459.  TAB2=TABLE;
  460.  TAB2.'TITRE'= TABLE ;
  461.  TAB2.1='TIRR ';
  462.  TAB2.'TITRE' . 1  = MOT 'Z=1';
  463.  TAB2.2='TIRC';
  464.  TAB2.'TITRE' . 2  = MOT 'Z=5';
  465.  TAB2.3='TIRL';
  466.  TAB2.'TITRE' . 3  = MOT 'Z=10';
  467.  TAB2.4='TIRM';
  468.  TAB2.'TITRE' . 4  = MOT 'Z=25';
  469.  
  470. ** ZOOM SUR L'AXE X
  471.  xbinf = 0.;
  472.  xbsup = 6.;
  473.  xbsup2 =1.2;
  474.  
  475.  si graph;
  476.  dess (prvy1 et prvy2 et prvy3 et prvy4) titr 'PROFIL VITESSE (Uz)' Xbor xbinf xbsup LEGE TAB2 titx 'R/DJ' tity 'Uz';
  477.  dess (prc1 et prc2 et prc3 et prc4) titr 'PROFIL CONCENTRATION' Xbor xbinf xbsup LEGE TAB2 titx 'R/DJ' tity 'Cf';
  478.  dess (prk1 et prk2 et prk3 et prk4) titr 'PROFIL EN CINET TURB' Xbor xbinf xbsup LEGE TAB2 titx 'R/DJ' tity 'K';
  479.  dess (prnut1 et prnut2 et prnut3 et prnut4) titr 'PROFIL VISCO TURB' Xbor xbinf xbsup LEGE TAB2 titx 'R' tity 'Muf';
  480.  finsi ;
  481.  
  482. ** DIAMETRES CARACTERISTIQUES
  483. **----------------------------
  484. DEBP ADIM VISU*EVOLUTIO ;
  485. LV = EXTR VISU 'ORDO' ;
  486. LX = EXTR VISU 'ABSC' ;
  487.  
  488. ** Recherche des grandeurs d'adimensionnalisation
  489. Vaxe = EXTR LV 1 ;
  490. V12 = Vaxe / 2. ;
  491. i = 1 ;
  492. REPETER BOUCLE (dime LV) ;
  493.  Vcour = EXTR LV i ;
  494.  SI (Vcour < V12) ;  quitter boucle   ;
  495.  SINON ; i = i + 1 ;
  496.  FINSI ;
  497. FIN BOUCLE ;
  498.  
  499. ** Determination de X12
  500. Xi = EXTR LX i ;
  501. Xii = EXTR LX (i - 1) ;
  502. Vi = EXTR LV i ;
  503. Vii = EXTR LV (i - 1) ;
  504. X12 = (((Vii - V12) * Xi) + ((V12 -Vi) * Xii)) / (Vii - Vi) ;
  505.  
  506. ** Adimensionnalisation
  507. LVA = LV / (Vaxe + 1e-20) ;
  508. LXA = LX / (X12 + 1e -20) ;
  509. VISUA = EVOL 'MANU' 'X / X1/2' LXA 'V / Vmax' LVA ;
  510. FINP X12 ;
  511.  
  512.  yaxe=COOR 2 axe;
  513.  lyaxe=PROG;
  514.  i=0;
  515.  repe bloc (nbno AXE);
  516.  i=i+1;
  517.  vyaxe=EXTR yaxe SCAL (AXE POIN i);
  518.  lyaxe=lyaxe et (prog vyaxe);
  519.  fin bloc;
  520.  lyaxe = ORDO lyaxe;
  521.  LMT=CHAN LIGNE MT;
  522.  diam2=prog;
  523.  i=0;
  524.  repeter boucl3 (dime lyaxe);
  525.  i=i+1;
  526. *mess 'i vaut' i;
  527.  yi=extr lyaxe i;
  528.  lign6=POIN MT DROI (POIN MT PROC (0. yi)) (POIN MT PROC (10. yi)) 1.e-3;
  529.  L6=elem lmt appuye lign6;
  530.  visu1=evol chpo rv.inco.'un' uy L6;
  531.  XD=adim visu1;
  532.  diam2=diam2 et (prog XD);
  533.  menage;
  534.  fin boucl3;
  535.  dcar1=evol manu diam2 lyaxe;
  536.  
  537.  si graph;
  538.  dess (dcar1)  titx 'R' tity  'Z';
  539.  finsi ;
  540.  
  541. ** GRANDEURS SUR L'AXE ADIMENSIONNEES
  542. **------------------------------------
  543.  
  544. ** Vitesse
  545.  evv4=EVOL CHPO (rv.inco.'un') uy (INVE AXE);
  546.  ev4=EXTR evv4 ORDO ;ev4=ENLEV ev4 1;
  547.  evn4=LOG ev4;
  548.  
  549.  z4=EXTR evv4 ABSC ; z4=ENLEV z4 1;zn4=LOG z4;
  550.  
  551.  nevv4=EVOL MANU  zn4  evn4;
  552.  Si graph;
  553.  dess evv4 titr 'UY SUR L AXE (ADIM)';
  554.  dess nevv4 titr 'UY SUR L AXE (ADIM) ECHELLE LOG';
  555.  Finsi;
  556. ** Concentration
  557. ** Attention, on ne peut calculer la concentration
  558. ** qu'à partir d'un grand nombre de pas de temps
  559. *
  560.   Xa=287.; Xhe=2079.;
  561.   cn1=KOPS (KOPS Xa '*' (rv.inco.'cn')) '/' (KOPS  (KOPS (KOPS 1 '-' (rv.inco.'cn')) '*' Xhe) '+' (KOPS (rv.inco.'cn') '*' Xa));
  562.  
  563.  
  564.   evc1=EVOL CHPO cn1 SCAL (INVE AXE);
  565.   ec1=extr evc1 ORDO; ec1=ENLE ec1 1; ecn1=LOG ec1;
  566.   zc1=extr evc1 ABSC; zc1=ENLE zc1 1; zcn1=LOG zc1;
  567.   nevc1=EVOL MANU  zcn1  ecn1;
  568.  
  569.  Si graph;
  570.  dess evc1 titr 'CONCENTRATION SUR L AXE (ADIM)';
  571.  dess nevc1 titr 'CONCENTRATION SUR L AXE (ADIM) ECHELLE LOG';
  572.  
  573.   TAB3=TABLE;
  574.   TAB3.'TITRE'= TABLE ;
  575.   TAB3.'TITRE' . 1  = MOT 'VITESSE AXIALE';
  576.   TAB3.2='TIRC';
  577.   TAB3.'TITRE' . 2  = MOT 'CONCENTRATION';
  578.  
  579.   dess (evv4 et evc1) LEGE TAB3;
  580.   dess (nevv4 et nevc1) LEGE  TAB3;
  581.  Finsi;
  582.  
  583. ** Energie
  584.  evk1=EVOL CHPO (rv.inco.'kn') SCAL (INVE AXE);
  585.  ek1=EXTR evk1 ORDO;ek1=ENLE ek1 1;ek1=ABS ek1; ekn1=LOG ek1;
  586.  zk1=EXTR evk1 ABSC;zk1=ENLE zk1 1;zk1=ABS zk1; zkn1=LOG zk1;
  587.  nevk1=EVOL MANU 'Z' zkn1 'K' ekn1;
  588.  
  589.  si graph;
  590.  dess evk1 titr 'EN. CINET. TURB. SUR L AXE (ADIM)';
  591.  dess nevk1 titr 'EN. CINET. TURB. SUR L AXE (ADIM) ECHELLE LOG';
  592.  finsi ;
  593.  
  594. **************************
  595. *RESULTATS COMPLEMENTAIRES*
  596. ***************************
  597. ** LIGNES DE COURANT
  598. **-----------------
  599.  corx=COOR 1 MT;
  600.  un=rv.inco.'un';
  601.  un1=EXCO un ux; un1=-1*corx*un1;
  602.  un2=EXCO un uy; un3=corx*un2;
  603.  un1=KCHT $MT SCAL SOMMET un1;
  604.  un2=KCHT $MT SCAL SOMMET un2;
  605.  un3=KCHT $MT SCAL SOMMET un3;
  606.  uncz=NOEL $MT un2;
  607.  uncz=KCHT $MT SCAL CENTRE uncz;
  608.  
  609.  rt2d=KOPS un 'ROT' $MT;
  610.  corx=KCHT $MT SCAL SOMMET corx;
  611.  corm=NOEL $MT corx;
  612.  
  613.  unn1=KCHT $BAT SCAL SOMMET un1;
  614.  unn2=KCHT $COT SCAL SOMMET un3;
  615.  unn3=KCHT $HTT SCAL SOMMET (-1*un1);
  616.  unn4=KCHT $AXE SCAL SOMMET (-1*un3);
  617.  
  618.  unn1=NOEL $BAT unn1;
  619.  unn1=KCHT $BAT SCAL CENTRE unn1;
  620.  unn2=NOEL $COT unn2;
  621.  unn2=KCHT $COT SCAL CENTRE unn2;
  622.  unn3=NOEL $HTT unn3;
  623.  unn3=KCHT $HTT SCAL CENTRE unn3;
  624.  unn4=NOEL $AXE unn4;
  625.  unn4=KCHT $AXE SCAL CENTRE unn4;
  626.  
  627.  sw=2*uncz - (corm*rt2d);
  628.  
  629.  rk=EQEX $MT 'OPTI' 'EF' 'IMPL'
  630.    ZONE $MT OPER LAPN 1. inco 'psi'
  631.    ZONE $MT OPER FIMP sw inco 'psi'
  632.    ZONE $BAT OPER FIMP unn1 inco 'psi'
  633.    ZONE $COT OPER FIMP unn2 inco 'psi'
  634.    ZONE $HTT OPER FIMP unn3 inco 'psi'
  635.    ZONE $AXE OPER FIMP unn4 inco 'psi'
  636.    CLIM 'psi' TIMP BAS 0.;
  637.  rk.inco.'psi'=KCHT $MT SCAL SOMMET 0.;
  638.  exec rk;
  639.  
  640.  si graph;
  641.  psi=rk.inco.'psi';
  642.  OPTI ISOV LIGNE;
  643.  TRAC psi MT  (CONT (MT)) titr 'LIGNES DE COURANT';
  644.  
  645.  psi1=EVOL CHPO psi SCAL COT;
  646.  db1=EXTR psi1 ORDO;db2=EXTR psi1 ABSC;
  647.  npsi1=EVOL MANU 'Z' db2 'DEBIT' db1 ;
  648.  DESS npsi1;
  649.  titr 'DEBIT D ENTRAINEMENT';
  650.  
  651.  FINSI ;
  652.  
  653.  
  654.  Si (NON COMPLET);
  655.  uaxref=prog
  656.    1.0000      0.95603       1.5700       1.9336       2.3534 2.5648       2.8291       2.9057       3.0356       3.0249
  657.    3.0781       3.0320       3.0420       2.9820       2.9616 2.8752       2.8122       2.7166       2.6450       2.5539
  658.    2.4868       2.4086       2.3435       2.2210       2.1198;
  659.  uaxref=uaxref et (prog
  660.    2.0306       1.9548       1.8835       1.8199       1.7601 1.7059       1.6548       1.6088       1.5662       1.5270
  661.    1.4891       1.4534       1.4197       1.3889       1.3608 1.3397 );
  662. zaxe = extr evv4 'ABSC';
  663. list (extr evv4 'ORDO');
  664. evuaxer=evol manu zaxe uaxref ;
  665. m=(INTG evuaxer 'ABSO') ;
  666. delt2= (INTG (evuaxer - evv4) 'ABSO')/m;
  667. mess ' deltaL2 : m=' m ' delt2=' delt2;
  668.   Si (delt2 > 1.e-2);ERREUR 5;Finsi;
  669.  FINSI;
  670.  
  671.  FIN ;
  672.  
  673.  
  674.  
  675.  
  676.  
  677.  
  678.  
  679.  
  680.  
  681.  
  682.  

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