TEXT(Inclined Supersonic Flow In A Duct
TITLE
  DISPLAY
  The problem considered is inclined supersonic flow in a plane
  duct. The flow enters at 6o to the horizontal at M=3, and
  Pin/Ptot=0.02722 with RHOin/RHOtot=0.07623. A compression wave is
  formed from the top corner of the duct inlet, and an expansion
  from bottom corner of the duct inlet. The calculations are made
  using the parabolic solver option with IPARAB=4. The flow is
  modelled as isentropic. The comparison between PHOENICS and the
  analytical results is as follows:
 
                           Expansion wave
                      Mach      P/Po       RHO/RHOo
        PHOENICS       3.3     0.0166      0.1047
        Analytical     3.3     0.0167      0.1045
 
                           Compression wave
                      Mach      P/Po       RHO/RHOo
        PHOENICS       2.7     0.0424       0.0537
        Analytical     2.7     0.0425       0.0542
 
  ENDDIS
   So as to allow a direct computation of dimensionless flow
   variables, the flow equations are normalised such that the flow
   variables can be interpreted as: P/Po; RHO/RHOo; T/To; U/Uref;
   H/Href;
 
   Here: Po, RHOo and To are the inlet total pressure, density
   and temperature; Uref=Ao/SQRT(GAMMA); and Href=(gam-1.)*Ho/gam,
   where Ho is the inlet total enthalpy and gam is the specific
   heat ratio. Ao is the acoustic velocity at To (see Palacio et al,
   Int.J.Heat Mass Transfer, Vol.33, No.6, p1193, [1990] ).
 
  PHOTON USE
  p
  parphi
 
 
 
  msg The grid used
  gr x 1
  msg            -
  msg Press  to continue
  pause
  gr off;red;gr ou x 1
  msg pressure distribution
  con p1 x 1 fi;0.01
  int 20
  gr ou x 1
  msg            -
  msg Press  to continue
  pause
  con off;r
  ed;vec x 1 sh
  msg velocity vectors
  msg            -
  msg Press  to continue
  pause;cl
  con on;con mach x 1 fi;.01
  msg Mach number contours
  msg
  msg Press e to END
  enduse
 
REAL(GASCON,GAMMA,PTOTAL,TTOTAL,RHOTOT,MACHI,AGAM1,RGAM)
REAL(PIN,TIN,VIN,QTOT,POWER,WIN,RHOIN,PI,COSA,SINA)
REAL(ANG1,ZLEN,WIDTH,AIN,FLOWIN)
PI=3.1415927
GASCON=1.0;GAMMA=1.4
 
PTOTAL=1.0;TTOTAL=1.0;RHOTOT=1.0;MACHI=3.0
   ** Inlet flow angle
ANG1=6.0*PI/180.
ZLEN=1.0;WIDTH=1.0
 
   ** Calculation of inlet velocity
AGAM1=GAMMA-1.;RGAM=1./GAMMA;POWER=GAMMA/AGAM1
PIN=PTOTAL/(1.+AGAM1*MACHI*MACHI/2.)**POWER
RHOIN=RHOTOT/(PTOTAL/PIN)**RGAM
WIN=MACHI*(GAMMA*PIN/RHOIN)**0.5
QTOT=WIN
   ** Calculation of Inlet Temperature
TIN=PIN/(GASCON*RHOIN)
 
COSA=COS(ANG1);SINA=SIN(ANG1)
WIN=QTOT*COSA;VIN=QTOT*SINA
    GROUP 1. Run title and other preliminaries
    GROUP 2. Transience; time-step specification
PARAB=T;IPARAB=4
    GROUP 3. X-direction grid specification
    GROUP 4. Y-direction grid specification
GRDPWR(Y,40,WIDTH,1.0)
    GROUP 5. Z-direction grid specification
GRDPWR(Z,40,ZLEN,1.0)
    GROUP 6. Body-fitted coordinates or grid distortion
    GROUP 7. Variables stored, solved & named
SOLVE(P1,W1,V1)
  ** Provide storage for the density.
STORE(RHO1,MACH)
    GROUP 8. Terms (in differential equations) & devices
TERMS(W1,Y,Y,N,Y,Y,Y);TERMS(V1,Y,Y,N,Y,Y,Y)
    GROUP 9. Properties of the medium (or media)
   ** Use Isentropic Density Law
RHO1=COMPRESS;RHO1A=RHOTOT/PTOTAL**RGAM;RHO1B=RGAM
RHO1C=0.;PRESS0=0.;DRH1DP=COMPRESS
    GROUP 10. Inter-phase-transfer processes and properties
    GROUP 11. Initialization of variable or porosity fields
FIINIT(P1)=PIN;FIINIT(W1)=WIN;FIINIT(RHO1)=RHOIN
FIINIT(V1)=VIN
    GROUP 13. Boundary conditions and special sources
INLET(IN,LOW,1,NX,1,NY,1,1,1,1)
VALUE(IN,P1,RHOIN*WIN);VALUE(IN,W1,WIN);VALUE(IN,V1,VIN)
    GROUP 16. Termination of iterations
SELREF=F
AIN=YVLAST;FLOWIN=WIN*RHOIN*AIN
RESREF(P1)=1.E-12*FLOWIN;RESREF(V1)=RESREF(P1)*WIN
RESREF(W1)=RESREF(V1);LITER(P1)=20
    GROUP 17. Under-relaxation devices
DENPCO=T
REAL(DTF);DTF=ZWLAST/WIN
RELAX(P1,LINRLX,1.0)
RELAX(V1,FALSDT,DTF);RELAX(W1,FALSDT,DTF)
    GROUP 18. Limits on variables or increments to them
    GROUP 19. Data communicated by satellite to GROUND
    GROUP 20. Preliminary print-out
    GROUP 21. Print-out of variables
    GROUP 22. Spot-value print-out
IYMON=1
    GROUP 23. Field print-out and plot control
ITABL=3;NPLT=2;NYPRIN=1;NZPRIN=NZ/4;TSTSWP=-1
IF(NZ.GT.1) THEN
+ IDISPA=1;IDISPB=1;IDISPC=NZ
ENDIF
LITHYD=10