PHOTON USE
  p
  parphi
  1 1 6
 
 
  up z
  vi -y
  gr ou y 1
  con u1 y 1 fi;0.001 
  msg velocity contours
  msg vertical dimension is time, horizontal is distance 
  pause;cls
  msg pressure contours
  con p1 y 1 fi;0.001
  ENDUSE
    
  GROUP 1. Run title and other preliminaries
TEXT(Shallow-Water Waves After Exit Shut
#cls
TITLE
  DISPLAY
  
  
    The pressure-dependent-porosity method is used to compute the
    transient distributions of velocity and pressure in a one-
    dimensional channel flow.
    
    Initially the water velocity is 0,8 meters per second 
    everywhere; but shutting the exit causes a wave to flow back 
    to the inlet.
    
    Since inflow continues unchecked, a further wave travels  in the
    opposite direction; and so on.
    
    The behaviour can be displayed via PHOTON, because idispa has 
    been set to 1, dictating that information is placed on the parphi 
    file after every second time step.
  ENDDIS
 
    GROUP 2. Transience; time-step specification
STEADY=F;GRDPWR(T,100,500.0,1.0)
 
    GROUP 3. X-direction grid specification
GRDPWR(X,50,1500.0,1.0)
 
    GROUP 4. Y-direction grid specification
YVLAST=500.0
 
    GROUP 5. Z-direction grid specification
ZWLAST=10.0
 
    GROUP 7. Variables stored, solved & named
SOLVE(P1,U1);STORE(EPOR,VPOR)
 
    GROUP 9. Properties of the medium (or media)
RHO1=1000.0;ENUL=1.E-3/RHO1
  **The following statements activate sequences in group 19
    section 3 of GREX3 which make the porosities linearly
    dependent on the pressure at iz=1. For the east cell faces,
      porosity= poria + porib * pressure,
    where porib = 1/(g*rho*z) by hydrostatics.
REAL(BBB);BBB=1./(9.81*RHO1*ZWLAST)
IPORIA=1;PORIA=1.0;PORIB=BBB
  **The dependence of the logarithm of the cell volume on
    pressure is represented by the next statement...
DRH1DP=BBB
    (Note that: d(rho*vol/dt)/dp = (rho*vol/dt)*d(ln vol)/dp.   )
 
    GROUP 11. Initialization of variable or porosity fields
FIINIT(P1)=0.0;FIINIT(U1)=0.08
FIINIT(EPOR)=1.;FIINIT(VPOR)=1.
 
    GROUP 13. Boundary conditions and special sources
INLET(INLET,WEST,1,1,1,1,1,1,1,LSTEP)
VALUE(INLET,P1,RHO1*0.08);VALUE(INLET,U1,0.08)
 
    GROUP 15. Termination of sweeps
LSWEEP=15
SPEDAT(SET,GXMONI,TRANSIENT,L,F) 
    GROUP 21. Print-out of variables
OUTPUT(P1,Y,N,N,N,N,Y);OUTPUT(U1,Y,N,N,N,N,Y)
OUTPUT(EPOR,Y,N,N,N,N,N);OUTPUT(VPOR,Y,N,N,N,N,N)
NTPRIN=5
 
    GROUP 22. Spot-value print-out
ITABL=2;IXMON=NX-1
 
    GROUP 23. Field print-out and plot control
PATCH(WAVE,PROFIL,1,NX,1,1,1,1,1,LSTEP)
PLOT(WAVE,P1,0.0,200.0);PLOT(WAVE,U1,0.,1.1*FIINIT(U1))
PATCH(HALFWAY,PROFIL,NX/2,NX/2,1,1,1,1,1,LSTEP)
PLOT(HALFWAY,P1,0.,200.);PLOT(HALFWAY,U1,0.,1.1*FIINIT(U1))
IDISPA=2