TALK=F;RUN( 1, 1)
TEXT(USP. Test 11

title
  DISPLAY 
    This case solves a two-dimensional steady hydrodynamics
    problem about laminar flow in plane channel.
    Length of channel is 1 m, height of channel is 0,1 m
    The viscosity varies with Y-coordinate.
    The minimum of viscosity is on walls. The maximum viscosity
    is in centre of channel.
  ENDDIS
  
REAL(X1,X2,Y1,Y2,DTF,VISC,UAVE,RHO)
REAL(REYNO)

       Initial values
NX = 100
NY = 20
NZ = 1
X1=0.;X2=1.0
Y1=0.;Y2=0.1
UAVE = 0.01
RHO = 1.
REYNO = 100

LSWEEP=800

RSET(M,NX,NY,NZ)
RSET(D,DOM,X2,Y2,1.E-3)
 


SOLVE(P1,U1,V1)
STORE(VISL)
PARSOL = F

REAL(VISMAX,VISMIN,VIS,DVIS)
VISMAX = 1.e-4
VISMIN = 1.e-5
mesg(Default value of minimum viscosity is :VISMIN:
mesg(Default value of maximum viscosity is :VISMAX:
mesg(Do you want to changed them? (y/n)
readvdu(ans,char,n)
if(:ans:.eq.y)then
mesg(Enter minimum viscosity
readvdu(VISMIN,real,VISMIN)
mesg(New value of minimum viscosity is :VISMIN:
mesg(Enter maximum viscosity
readvdu(VISMAX,real,VISMAX)
mesg(New value of maximum viscosity is :VISMAX:
endif
      Average velocity
VISC = (VISMAX+VISMIN)/2
UAVE = REYNO*VISC/Y2
ENUL=VISC
RHO1=RHO



DVIS = 2*(VISMAX-VISMIN)/NY

TERMS (U1 ,Y,Y,Y,N,N,N)
TERMS (V1 ,Y,Y,Y,N,N,N)

REAL(XCOR,DXC,YCOR,DYC,U1C,P1C,HY,P1AVE)
DYC=2./NY
DXC=X2/NX
HY=Y2*0.5


    Following commands set the initial fields and inlet boundary source
VIS = VISMIN + 0.5*DVIS
DO JJ=1,NY
+ YCOR=0.5*DYC+(JJ-1)*DYC-1
+ U1C = 1.5*UAVE*(1-YCOR*YCOR)
+ PATCH(IV:JJ:,INIVAL,1,NX,:JJ:,:JJ:,1,NZ,1,1)
+ INIT(IV:JJ:,VISL,0.0,VIS)
+ PATCH(IN:JJ:,WEST,1,1,:JJ:,:JJ:,1,NZ,1,1)
+ COVAL(IN:JJ:,P1,FIXFLU,U1C*RHO)
+ COVAL(IN:JJ:,U1,ONLYMS,U1C)
+ COVAL(IN:JJ:,V1,ONLYMS,0.)
+IF(JJ.LT.NY/2)THEN
+ VIS = VIS + DVIS
+ELSE
+ IF(JJ.GT.NY/2)THEN
+ VIS = VIS - DVIS
+ ENDIF
+ENDIF
ENDDO

PATCH(OUTLET,EAST,NX,NX,1,NY,1,1,1,1)
COVAL(OUTLET,P1, 1.0E+3 ,0.0)
NONORT=T
     ========= Walls =====================
PATCH(WALL1,SWALL,1,NX,1,1,1,1,1,1)
COVAL(WALL1,U1,1.,0.0)
COVAL(WALL1,V1,1.,0.0)

PATCH(WALL2,NWALL,1,NX,NY,NY,1,1,1,1)
COVAL(WALL2,U1,1.,0.0)
COVAL(WALL2,V1,1.,0.0)

DTF=1.E+01
RELAX(P1,LINRLX,1.)
RELAX(U1,FALSDT,DTF)
RELAX(V1,FALSDT,DTF)


USP    = T
UAUTO  = F
USPDBG = F
UTCPLT = T
USPVTK = T
USPIMB = F
MXLEV  = 0
MYLEV  = 0
MZLEV  = 0
DOMAT  = -1
MINPRP = -1
MAXPRP = 250
CELLST = 10
FACEST = 10

mesg(Do you want to use collocated arrangement (y) or staggered one (n)? (y/n)
readvdu(ans,char,n)
if(:ans:.eq.y)then
 SPEDAT(SET,USP,METHOD,I,1)
 RELAX(P1  ,LINRLX, 6.000000E-01)
 RELAX(U1  ,FALSDT, 1.000000E-01)
 RELAX(V1  ,FALSDT, 1.000000E-01)
mesg(Do you want to use SIMPLEST (y) or SIMPLE (n)? (y/n)
readvdu(ans,char,n)
if(:ans:.eq.y)then
 RELAX(P1  ,LINRLX, 1.000000E-00)
 SPEDAT(SET,USP,SIMPLEST,L,T)
endif

endif

mesg(Do you want to view results in the centres of cells? (y/n)
readvdu(ans,char,n)
if(:ans:.eq.y)then
SPEDAT(SET,USPIO,VERTCENT,L,F)
endif




SELREF = T; RESFAC =1.0E-7
ECHO=F;IXMON=NX-1;IYMON=NY/2
TSTSWP=-1
NXPRIN= 1; NYPRIN=1
STOP