talk=f;run(1,1)
  PHOTON USE
   p ;;;;;
 
   set vec av off
   msg Velocity vectors:
   vec i 1 sh
   msg Press any key to continue ...
   pause
   cl
   msg Pressure contours:
   con p1 i 1 fi;0.01
   pause
   msg Press any key to continue ...
   cl
   msg Swirl velocity contours:
   con u1 i 1 fi;0.01
   pause
   msg Press any key to continue ...
   cl
   msg W-component contours:
   con w1 i 1 fi;0.01
   msg Press any key to continue ...
  ENDUSE
  **************************************************************
    GROUP 1. Run title and other preliminaries
INTEGER(ICAS)
MESG(
MESG(
MESG(
MESG(            SWIRLING FLOW IN THE EXPANDING PIPE
MESG(
MESG(
MESG( This case illustrates how to use MIGAL with swirling flows. 
MESG( From library case F105.
MESG(
MESG(
MESG(
MESG( Enter required case :  1 - Laminar   Re=100   (default)
MESG(                        2 - Turbulent Re=1.E+5 
MESG(
READVDU(ICAS,INT,1)

BOOLEAN(LTURB)
LTURB= F
IF(ICAS.EQ.2) THEN
+ LTURB=T  
ENDIF

INTEGER(SLV)
MESG( Enter required solver ID :  1 - MIGAL (default)
MESG(                             2 - SIMPLEST 
MESG(
READVDU(SLV,INT,1)
 
IF(SLV.EQ.1) THEN
+ TEXT(Swirling Flow : MIGAL
ENDIF
IF(SLV.EQ.2) THEN
+ TEXT(Swirling Flow : SIMPLEST
ENDIF

  **************************************************************

INTEGER(NZST,NYST,NZ1,NY1,IC)
REAL(REYNO,WIN,UIN,HSTEP,LSTPP,HCHAN,LCHAN,DXX,DTHYD,UCUR)
REAL(TKEIN,EPSIN)
  ** Problem definition:
IF(LTURB) THEN
+ REYNO= 1.0E5
ELSE
+ REYNO= 100.0
ENDIF
HSTEP= 1.0;  LSTPP= 2.0;  HCHAN= 2.0;  LCHAN= 20.0
WIN  = 1.0;  DXX  = 0.1;  UIN  = 3.0;  ENUL = WIN*HSTEP/REYNO
NZST = 6;    NYST = 6;    NZ1  = 27;   NY1  = 12
    GROUP 3. X-direction grid specification
CARTES= F; NX= 1; GRDPWR(X,NX,DXX,1.0)
    GROUP 4. Y-direction grid specification
NREGY = 2; NY = NY1+NYST
IREGY = 1; GRDPWR(Y,-NY1, HCHAN-HSTEP,1.0)
IREGY = 2; GRDPWR(Y,-NYST,HSTEP,      1.0)
    GROUP 5. Z-direction grid specification
NREGZ = 2; NZ = NZST+NZ1
IREGZ = 1; GRDPWR(Z,NZST,LSTPP,      -1.4)
IREGZ = 2; GRDPWR(Z,NZ1, LCHAN-LSTPP, 1.3)

    GROUP 7. Variables stored, solved & named
SOLVE(P1,W1,V1,U1)
IF(LTURB) THEN
+ TURMOD(KEMODL-LOWRE)
+ STORE(ENUT,LEN1)
ENDIF

    GROUP 11. Initialization of variable or porosity fields
INIADD= F
CONPOR(STEP,0.0,CELL,1,NX,NY1+1,NY,1,NZST)
FIINIT(W1)  = 1.E-3
FIINIT(V1)  = 1.E-3
IF(LTURB) THEN
+ TKEIN = (0.05*WIN)**2
+ EPSIN = TKEIN**1.5*0.1643/0.09*HSTEP
+ FIINIT(P1)= 1.3E-4
+ FIINIT(KE)= TKEIN
+ FIINIT(EP)= EPSIN
ENDIF
    GROUP 13. Boundary conditions and special sources
    ** Inlet.
DO II = 1,NY1
+ UCUR = (2*II-1)*UIN/NY1/2
+ INLET(INL:II:,LOW,1,NX,II,II,1,1,1,LSTEP)
+ VALUE(INL:II:,P1,WIN*RHO1)
+ VALUE(INL:II:,W1,WIN)
+ VALUE(INL:II:,V1,0.0)
+ VALUE(INL:II:,U1, UCUR)
+ IF(LTURB) THEN
+  VALUE(INL:II:,KE,TKEIN); VALUE(INL:II:,EP,EPSIN)
+ ENDIF
ENDDO
    ** Outlet.
PATCH(OUT1,HIGH,1,1,1,NY,NZ,NZ,1,1)
COVAL(OUT1,P1,1.E5,0.0)
    ** Walls.
PATCH(WN1,NWALL,1,NX,  NY1,NY1,     1,NZST  ,1,LSTEP)
PATCH(WN2,NWALL,1,NX,   NY, NY,NZST+1,NZ    ,1,LSTEP)
PATCH(WL, LWALL,1,NX,NY1+1, NY,NZST+1,NZST+1,1,LSTEP)
IF(LTURB) THEN
+ COVAL(WN1,KE, 1.0,   0.0)
+ COVAL(WN2,KE, 1.0,   0.0)
+ COVAL(WL, KE, 1.0,   0.0)
+ COVAL(WN1,W1, LOGLAW,0.0)
+ COVAL(WN1,V1, LOGLAW,0.0)
+ COVAL(WN1,U1, LOGLAW,0.0)
+ COVAL(WN2,W1, LOGLAW,0.0)
+ COVAL(WN2,V1, LOGLAW,0.0)
+ COVAL(WN2,U1, LOGLAW,0.0)
+ COVAL(WL ,W1, LOGLAW,0.0)
+ COVAL(WL ,V1, LOGLAW,0.0)
+ COVAL(WL ,U1, LOGLAW,0.0)
ELSE
+ COVAL(WN1,W1, 1.0,0.0)
+ COVAL(WN1,V1, 1.0,0.0)
+ COVAL(WN1,U1, 1.0,0.0)
+ COVAL(WN2,W1, 1.0,0.0)
+ COVAL(WN2,V1, 1.0,0.0)
+ COVAL(WN2,U1, 1.0,0.0)
+ COVAL(WL ,W1, 1.0,0.0)
+ COVAL(WL ,V1, 1.0,0.0)
+ COVAL(WL ,U1, 1.0,0.0)
ENDIF
    GROUP 15. Termination of sweeps

IF(LTURB) THEN
+ LSWEEP = 2000
ELSE
+ LSWEEP = 2000
ENDIF

TSTSWP = -1
    GROUP 16. Termination of iterations
selref=F
resref(P1)=1.E-6
resref(U1)=1.E-6
resref(V1)=1.E-6
resref(W1)=1.E-6
IF(LTURB) THEN
+ resref(KE)=1.E-6
+ resref(EP)=1.E-6
ENDIF
RESFAC = 1.E-3
    GROUP 17. Under-relaxation devices
DTHYD = HCHAN/(NY*WIN)
RELAX(P1,LINRLX,0.25)
RELAX(U1,FALSDT,DTHYD)
RELAX(V1,FALSDT,DTHYD)
RELAX(W1,FALSDT,DTHYD)
IF(LTURB) THEN
+ KELIN= 1
+ RELAX(KE,FALSDT,DTHYD)
+ RELAX(EP,FALSDT,DTHYD)
ENDIF
    GROUP 22. Spot-value print-out
IXMON = NZ/2+1; IYMON = NY/2+1; IXMON=1
    GROUP 23. Field print-out and plot control
NYPRIN= 1; NZPRIN = 1


IF(SLV.EQ.1) THEN
+IF(LTURB) THEN
+ LSWEEP=600
+ SPEDAT(MIGAL, SOLVED1, C, HYDRO)
+ RELAX(W1,FALSDT,DTHYD*10)
+ SPEDAT(MIGAL, SOLVED2, C, U1)
+ SPEDAT(MIGAL, LINRLX2, r, 1.)
+ RELAX(U1,FALSDT,DTHYD*10000)
+ SPEDAT(MIGAL, SOLVED3, C, KE)
+ SPEDAT(MIGAL, LINRLX3, r, 1.)
+ SPEDAT(MIGAL, SOLVED4, C, EP)
+ SPEDAT(MIGAL, LINRLX4, r, 1.)
+ELSE
+ LSWEEP=100
+ SPEDAT(MIGAL, SOLVED1, C, HYDRO)
+ SPEDAT(MIGAL, LINRLX1, r, .9)
+ SPEDAT(MIGAL, SOLVED2, C, U1)
+ SPEDAT(MIGAL, LINRLX2, r, 1.)
+ RELAX(U1,FALSDT,DTHYD*10000)
+ RELAX(V1,FALSDT,DTHYD*10)
+ RELAX(W1,FALSDT,DTHYD*10)
+ENDIF
ENDIF