Author's explanation 01.04.09
----------------------------------------
This case is a development of case 290 which concerns turbulent
flow over and downstream of a backward-facing step but uses a grid
which is far too coarse for accuracy.
It has therefore been modified so as to enable the effects of the
grid-defining parameters NX, NY, XPOWER and YPOWER to be easily varied.
This is done by changing the name of the STEP Patch to .STEP .
The '.' which starts its name is a signal to EARTH that locating
arguments are geometric rather than indicial, albeit non-
dimensionalised by reference to the domain dimensions.
[The following statement was true when the .patch was first introduced,
but it the multiply-by-1000 practice is no longer needed:
A multiplier of 1000 is used (and matched by division on arrival
in EARTH) because, although the arguments are real numbers, the
SATELLITE reduces them to whole numbers.]
Although it is not necessary, because they lie at domain
boundaries, the inlet, outlet and wall boundary conditions are
also expressed via dot-patches so as to test and illustrate the
technique.
DBS
-----------------------------------------------------------------
Setting the grid fineness and distribution interactively
NX=20; ny=15 ! defaults
real(xpower,ypower)
xpower=1.0; ypower=1.0
if(iqalib.ne.0) then
+ ans=y
endif
label message
mesgm(nx=:nx: ny=:ny: xpower=:xpower: ypower=:ypower: If ok, RETURN
readvdu(ans,char,y)
if(:ans:.eq.y) then
goto begin
else
mesg(insert required NX, or RETURN
readvdu(nx,int,nx)
mesg(insert required NY, or RETURN
readvdu(ny,int,ny)
mesg(insert required XPOWER, or RETURN
readvdu(xpower,real,xpower)
mesg(insert required YPOWER, or RETURN
readvdu(ypower,real,ypower)
goto message
endif
label begin
PHOTON USE
p;;;;
use patgeo
gr ou z 1
msg contours of velocity, and velocity vectors
con u1 z 1 fi;0.01;vec z 1
msg press RETURN for length scale
pause;con off;vec off;red
msg contours of turbulence length scale
con len1 z 1 fi;0.01;pause;con off;red
msg contours of effective viscosity
con enut z 1 fi;0.01
enduse
GROUP 1. Run title and other preliminaries
TEXT(Backward-Facing Step; K-E or LVEL
TITLE
DISPLAY
This case involves steady, turbulent flow over a backward-facing
step, of height h, in a two-dimensional channel of width 3*h.
The calculations are started at a distance 4h upstream of the step
and terminated at a distance 16h downstream.
//////////////////////// wall ///////////////////////
-----------------------------------------------------
Pressure
Inlet -------> ---------> fixed at
zero
/________________
| ////////////////| Exit
| wall /| Recirculation
| /| <---- ----->
\ constant /|____________________________________
prescribed mass //////////////////////////////////////
inflow rate wall
ENDDIS
REAL(HEIGHT,WIDTH,CLEN,SLEN,UIN,TKEIN,EPSIN,GENUT)
UIN=1.0
** Calculation of domain specifications
HEIGHT=0.0381;
WIDTH=3.*HEIGHT
SLEN=4.*HEIGHT;
CLEN=20.*HEIGHT
XULAST=CLEN
YVLAST=WIDTH
GRDPWR(X,NX,XULAST,XPOWER)
GRDPWR(Y,NY,YVLAST,YPOWER)
GROUP 7. Variables stored, solved & named
SOLVE(P1,U1,V1);SOLUTN(P1,Y,Y,Y,N,N,N);STORE(ENUT,LEN1)
mesg(Turbulence model is KEMODL. Use LVEL instead? (y/n)
CHAR(char1)
READVDU(ANS,CHAR,N)
IF(:ANS:.EQ.Y) THEN
TURMOD(LVEL)
char1=LVEL
store(wgap)
ELSE
TURMOD(KEMODL)
char1=KE-EPS
KELIN=3
mesgb(kelin=3; no under-relaxation is used for ke or ep0
ENDIF
GROUP 9. Properties of the medium (or media)
GROUP 11. Initialization of variable or porosity fields
** Calculation of KE (where fric=0.018)...
TKEIN=0.25*UIN*UIN*0.018
** Calculation of EP (where lmix=0.09 x h)...
EPSIN=TKEIN**1.5*0.1643/(0.09*height)
** Initial values
FIINIT(U1)=UIN;FIINIT(P1)=1.3E-4;FIINIT(KE)=TKEIN;FIINIT(EP)=EPSIN
** Initialization of variables in blocked region
EGWF=T;STORE(PRPS)
declarations and settings of dot-patch arguments
REAL(RXF,RXL,RYF,RYL,RZF,RZL)
The next three de-activated statements illustrate the now discarded
multiplication-by-1000 practice.
RXF=0.0; RXL=SLEN/XULAST*1000
RYF=0.0; RYL=HEIGHT/YVLAST*1000
RZF=0.0; RZL=ZWLAST/ZWLAST*1000
RXF=0.0; RXL=SLEN/XULAST
RYF=0.0; RYL=HEIGHT/YVLAST
RZF=0.0; RZL=ZWLAST/ZWLAST
PATCH(.STEP,INIVAL,RXF,RXL,RYF,RYL,RZF,RZL,1,1) ! Important warning
! Although it might seem reasonable to do so, DO NOT set the
! last but one arguments above to 0 because doing so confuses the
! indexing system of the satellite, with unforeseeable results.
! Note also that 1,1000 is satisfactory for the last two
! arguments of the group 13 dot-patches; but the equally
! reasonable 0,1000, although OK for the x, y and z arguments,
! causes the aforesaid damage.
COVAL(.STEP,PRPS,0.0,198.0)
TEXT( :char1: model; EGWF= :EGWF:
TITLE
GROUP 13. Boundary conditions and special sources
** Inlet
PATCH(INLET,WEST,1,1,1,NY,1,1,1,1) ! the no-dot patch
PATCH(.INLET,WEST,0,0,0,1000,0,1000,1,1000) ! the discarded
! multiply-by-1000 practice
PATCH(.INLET,WEST,0,0,0,1,0,0,0,1)
COVAL(.INLET,P1,FIXFLU,RHO1*UIN)
COVAL(.INLET,U1,ONLYMS,UIN)
COVAL(.INLET,KE,ONLYMS,TKEIN)
COVAL(.INLET,EP,ONLYMS,EPSIN)
** Exit
PATCH(OUTLET,EAST,NX,NX,1,NY,1,1,1,1) ! the no-dot patch
PATCH(.OUTLET,EAST,1000,1000,0,1000,0,1000,1,1000) ! the discarded
! multiply-by-1000 practice
PATCH(.OUTLET,EAST,1,1,0,1,0,1,0,1)
COVAL(.OUTLET,P1,1.0E5,0.0)
** N-wall
WALL (WFNN,NORTH,1,NX,NY,NY,1,1,1,1) ! the no-dot patch
WALL (.WFNN,NORTH,0,1000,1000,1000,0,1000,1,1000) ! the discarded
! multiply-by-1000 practice
WALL (.WFNN,NORTH,0,1,1,1,0,1,1,1)
** S-wall
WALL (.WFNS,SOUTH,0,1000,0,0,0,1000,1,1000) ! the discarded
! multiply-by-1000 practice
WALL (.WFNS,SOUTH,0,1,0,0,0,1,0,1)
GROUP 15. Termination of sweeps
LSWEEP=1000;SELREF=T;RESFAC=0.01
GROUP 16. Termination of iterations
LITER(U1)=20;LITER(V1)=20
GROUP 17. Under-relaxation devices
RELAX(U1,FALSDT,SLEN/UIN);RELAX(V1,FALSDT,SLEN/UIN)
relax(ke,linrlx,0.5);relax(ep,linrlx,0.5)
GROUP 22. Monitor print-out
IYMON=NY/2;IXMON=NX/2;NPRMON=100;TSTSWP=12;UWATCH=T
GROUP 23. Field print-out and plot control
NPLT=1;IPLTL=LSWEEP;ITABL=1
PATCH(MAP,CONTUR,1,NX,1,NY,1,1,1,1)
PLOT(MAP,P1,0.0,10.0);PLOT(MAP,U1,0.0,10.0);PLOT(MAP,V1,0.0,10.0)
tstswp=-1
lsweep=1
debug=t
dbindx=t
search=t
idbf0=414
flag=t
stop