** MBFGE Test: 3D flow in stirred tank.
  **************************************************************
  DISPLAY
  ----------------------------------------------------------
  The case concerns unsteady 3D laminar or turbulent flow in
  the stirred reactor tank. Geometry of the reactor is  that
  presented in H.Takeda, K.Narasaki, H.Kitajima, S.Sudoh, M.
  Onofusa and S.Iguchi 'Numerical Simulation Of Mixing Flows
  In Agitated Vessels',  Computerts Fluids  Vol.22,  No 2/3,
  pp.223-228, 1993.
 
  The main purpose of the case is to demonstrate the way  to
  set problem using  sliding LINK option  of MBFGE. User  is
  advised to pay attention on the use of the following  PIL-
  commands:
    READCO(MBGR+Y)
    MPATCH(1,MBS1.2,NORTH,...); MPATCH(2,MBS2.1,SOUTH,...)
 
  NOTE!, that the stacking along Y-axis, achieved by the use
  +Y in READCO, provides for XCYCLE=T treatment.
  ----------------------------------------------------------
  ENDDIS
L(PAUSE
  **************************************************************
BOOLEAN(LTURB); LTURB= T
  **************************************************************
    GROUP 1. Run title and other preliminaries
REAL(RINT,RIMP,RBAF,REXT,REYNO,PI,OMEGA,RLNK,DTHYD,TKEIN,EPSIN)
REAL(UIMP,HBOT,HIMP,HCOR,HRST,RCOR,RTMP,UTMP,COEF)
INTEGER(NX1,NY11,NY12,NY13,NY1,NZ1,NX2,NY21,NY22,NY2,NZ2,IBAF)
INTEGER(NZ11,NZ12,NZ13,NZ14,NZ21,NZ22,NZ23,NZ24,NY14,JC,KC,JJ1)
  ** Problem definition:
PI   = 3.14159
IF(LTURB) THEN
+ TEXT(MBFGE: 3D (Re=2000.; K-E) sliding grid.
+ REYNO= 2000.;  COEF= GRND3
ELSE
+ TEXT(MBFGE: 3D (Re=100.) sliding grid.
+ REYNO= 100.;  COEF= 1.0
ENDIF
REXT = 1.0;        RINT = 0.08*REXT;   RCOR= 0.15*REXT
RIMP = 0.52*REXT;  RLNK = 0.6*REXT;    RBAF= 0.84*REXT
HBOT = 0.42*REXT;  HIMP = 0.54*REXT;   HCOR= 0.84*REXT
HRST = 2.2*REXT;   OMEGA= 3.0;         UIMP= OMEGA*RIMP
TKEIN= 0.25*UIMP*UIMP*0.018;   EPSIN= TKEIN**1.5/RIMP
NX1 = 12;  NY11= 2;   NY12= 3;  NY13= 9;  NY14= 3
NY1 = NY11+NY12+NY13+NY14
NZ11=  8;  NZ12= 4;   NZ13= 5;  NZ14= 10
NZ1 = NZ11+NZ12+NZ13+NZ14
NX2 = 12;  NY21= 8;   NY22= 6;  NY2 = NY21+NY22
NZ21= NZ11;  NZ22= NZ12;  NZ23= NZ13;  NZ24= NZ14
NZ2 = NZ21+NZ22+NZ23+NZ24
IBAF= NX1/2; DTHYD= PI/2./NX1/OMEGA
    GROUP 2. Transience; time-step specification
STEADY= F; GRDPWR(T,NX1+1,1.0,1.0)
    GROUP 6. Body-fitted coordinates or grid distortion
BFC= T; GSET(D,NX1,NY1,NZ1,REXT,REXT,0.1)
GSET(P, P1,0.0, 0.0, 0.0);  GSET(P, P2,0.0,RINT, 0.0)
GSET(P, P3,0.0,RCOR, 0.0);  GSET(P, P4,0.0,RIMP, 0.0)
GSET(P, P5,0.0,RLNK, 0.0);  GSET(P, P6,0.0,RLNK,HBOT)
GSET(P, P7,0.0,RLNK,HIMP);  GSET(P, P8,0.0,RLNK,HCOR)
GSET(P, P9,0.0,RLNK,HRST);  GSET(P,P10,0.0,RIMP,HRST)
GSET(P,P11,0.0,RCOR,HRST);  GSET(P,P12,0.0,RINT,HRST)
GSET(P,P13,0.0, 0.0,HRST);  GSET(P,P14,0.0, 0.0,HCOR)
GSET(P,P15,0.0, 0.0,HIMP);  GSET(P,P16,0.0, 0.0,HBOT)
GSET(L, L1, P1, P2,NY11, 1.0); GSET(L, L2, P2, P3,NY12, 1.0)
GSET(L, L3, P3, P4,NY13, 1.0); GSET(L, L4, P4, P5,NY14, 1.0)
GSET(L, L5, P5, P6,NZ11,S1.5); GSET(L, L6, P6, P7,NZ12, 1.0)
GSET(L, L7, P7, P8,NZ13, 1.2); GSET(L, L8, P8, P9,NZ14, 1.2)
GSET(L, L9, P9,P10,NY14, 1.0); GSET(L,L10,P10,P11,NY13, 1.0)
GSET(L,L11,P11,P12,NY12, 1.0); GSET(L,L12,P12,P13,NY11, 1.0)
GSET(L,L13,P13,P14,NZ14,-1.2); GSET(L,L14,P14,P15,NZ13,-1.2)
GSET(L,L15,P15,P16,NZ12, 1.0); GSET(L,L16,P16, P1,NZ11,S1.5)
GSET(F,F1,P1,P2.P3.P4,P5,P6.P7.P8,P9,P10.P11.P12,P13,P14.P15.P16)
GSET(M,F1,+J+K,1,1,1)
GSET(C,I:NX1+1:,F,I1,1,NY1,1,NZ1,RZ,-PI/2.,0.0,0.0,INC,1.0)
 DUMPC(MBGR1)
BFC= T; GSET(D,NX2,NY2,NZ2,REXT,REXT,0.1)
GSET(P, P1,0.0,RLNK, 0.0);  GSET(P, P2,0.0,RBAF, 0.0)
GSET(P, P3,0.0,REXT, 0.0);  GSET(P, P4,0.0,REXT,HBOT)
GSET(P, P5,0.0,REXT,HIMP);  GSET(P, P6,0.0,REXT,HCOR)
GSET(P, P7,0.0,REXT,HRST);  GSET(P, P8,0.0,RBAF,HRST)
GSET(P, P9,0.0,RLNK,HRST);  GSET(P,P10,0.0,RLNK,HCOR)
GSET(P,P11,0.0,RLNK,HIMP);  GSET(P,P12,0.0,RLNK,HBOT)
GSET(L, L1, P1, P2,NY21,S1.5); GSET(L, L2, P2, P3,NY22,S1.5)
GSET(L, L3, P3, P4,NZ21,S1.5); GSET(L, L4, P4, P5,NZ22, 1.0)
GSET(L, L5, P5, P6,NZ23, 1.2); GSET(L, L6, P6, P7,NZ24, 1.2)
GSET(L, L7, P7, P8,NY22,S1.5); GSET(L, L8, P8, P9,NY21,S1.5)
GSET(L, L9, P9,P10,NZ24,-1.2); GSET(L,L10,P10,P11,NZ23,-1.2)
GSET(L,L11,P11,P12,NZ22, 1.0); GSET(L,L12,P12, P1,NZ21,S1.5)
GSET(F,F2,P1,P2,P3,P4.P5.P6,P7,P8,P9,P10.P11.P12)
GSET(M,F2,+J+K,1,1,1)
GSET(C,I:NX2+1:,F,I1,1,NY2,1,NZ2,RZ,-PI/2.,0.0,0.0,INC,1.0)
 DUMPC(MBGR2)
  ** Assemble blocks:
NUMBLK = 2; READCO(MBGR+Y);
VIEW
  ** Next is instead of MBLINK(1,SOUTH,2,NORTH)
     to set sliding LINK:
MPATCH(1,MBS1.2,NORTH,1,NX1,NY1,NY1,1,NZ1,1,LSTEP)
MPATCH(2,MBS2.1,SOUTH,1,NX2,  1,  1,1,NZ2,1,LSTEP)
    GROUP 7. Variables stored, solved & named
STORE(VPOR,EPOR); SOLVE(P1,U1,V1,W1)
IF(LTURB) THEN
+ TURMOD(KEMODL); KELIN= 1; STORE(ENUT,GEN1)
+ SOLUTN(KE,Y,Y,Y,N,P,P); SOLUTN(EP,Y,Y,Y,N,P,P)
+ FIINIT(KE)= TKEIN;  FIINIT(EP)= EPSIN
ENDIF
L($F150)
    GROUP 9. Properties of the medium (or media)
ENUL= RIMP*RIMP*OMEGA/REYNO;  RHO1= 1.189
    GROUP 11. Initialization of variable or porosity fields
INIADD= F
    *** Impeller:
    ***** Axis:
JC= NY11+NY12;  KC= NZ11+NZ12+NZ13
MPATCH(1,CORE,INIVAL,1,NX1,1,JC,NZ11+1,KC,1,1)
 COVAL(CORE,VPOR,0.0,0.0)
MPATCH(1,AXIS,INIVAL,1,NX1,1,NY11,KC+1,NZ1,1,1)
 COVAL(AXIS,VPOR,0.0,0.0)
    ***** Impeller Blade:
MPATCH(1,IMPEL,INIVAL,IBAF,IBAF,JC+1,JC+NY13,NZ11+1,NZ11+NZ12,1,1)
 COVAL(IMPEL,EPOR,0.0,0.0)
    ***** Baffle:
MPATCH(2,BAFF,INIVAL,IBAF,IBAF,NY21+1,NY2,1,NZ2,1,1)
 COVAL(BAFF,EPOR,0.0,0.0)
    GROUP 13. Boundary conditions and special sources
XCYCLE= T
    ** Coriolis forces for rotating domain:
ANGVEL= OMEGA
ROTAXA= 0.0;  ROTAYA= 0.0;  ROTAZA= 0.0
ROTAXB= 0.0;  ROTAYB= 0.0;  ROTAZB= 1.0
MPATCH(1,SLIDRT1,CELL,1,NX1,1,NY1,1,NZ1,1,LSTEP)
    ** Walls.
    ** Bottom stationary wall:
DO JJ= 1,NY1
+ JJ1=JJ+1; RTMP=(YC(1,JJ,1)+YC(2,JJ,1)+YC(1,JJ1,1)+YC(2,JJ1,1))/4.
+ UTMP= RTMP*OMEGA
+MPATCH(1,WLS:JJ:,LWALL,1,NX1,JJ,JJ,1,1,1,LSTEP)
+ COVAL(WLS:JJ:,UC1,COEF,UTMP); COVAL(WLS:JJ:,VC1,COEF,0.0)
+ COVAL(WLS:JJ:,WC1,COEF, 0.0)
+IF(LTURB) THEN
+ COVAL(WLS:JJ:,KE,GENLAW,GENLAW); COVAL(WLS:JJ:,EP,GENLAW,GENLAW)
+ENDIF
ENDDO
MPATCH(1,WS1,SWALL,1,NX1,JC+1,JC+1,NZ11+1,KC,1,LSTEP)
 COVAL(WS1,UC1,COEF,0.0); COVAL(WS1,VC1,COEF,0.0)
 COVAL(WS1,WC1,COEF,0.0)
MPATCH(1,WS2,SWALL,1,NX1,NY11+1,NY11+1,KC+1,NZ1,1,LSTEP)
 COVAL(WS2,UC1,COEF,0.0); COVAL(WS2,VC1,COEF,0.0)
 COVAL(WS2,WC1,COEF,0.0)
MPATCH(1, WH,HWALL,1,NX1,1,JC,NZ11,NZ11,1,LSTEP)
 COVAL(WH,UC1,COEF,0.0); COVAL(WH,VC1,COEF,0.0)
 COVAL(WH,WC1,COEF,0.0)
MPATCH(1, WL,LWALL,1,NX1,NY11+1,JC,KC+1,KC+1,1,LSTEP)
 COVAL(WL,UC1,COEF,0.0); COVAL(WL,VC1,COEF,0.0)
 COVAL(WL,WC1,COEF,0.0); KC= NZ11+NZ12
MPATCH(1,WIE,EWALL,IBAF,IBAF,JC+1,JC+NY13,NZ11+1,KC,1,LSTEP)
 COVAL(WIE,UC1,COEF,0.0); COVAL(WIE,VC1,COEF,0.0)
 COVAL(WIE,WC1,COEF,0.0)
MPATCH(1,WIW,WWALL,IBAF+1,IBAF+1,JC+1,JC+NY13,NZ11+1,KC,1,LSTEP)
 COVAL(WIW,UC1,COEF,0.0); COVAL(WIW,VC1,COEF,0.0)
 COVAL(WIW,WC1,COEF,0.0)
MPATCH(2, WN,NWALL,1,NX2,NY2,NY2,1,NZ2,1,LSTEP)
 COVAL(WN,UC1,COEF,0.0); COVAL(WN,VC1,COEF,0.0)
 COVAL(WN,WC1,COEF,0.0)
MPATCH(2,WLL,LWALL,1,NX2,1,NY2,1,1,1,LSTEP)
 COVAL(WLL,UC1,COEF,0.0); COVAL(WLL,VC1,COEF,0.0)
 COVAL(WLL,WC1,COEF,0.0)
MPATCH(2,WBE,EWALL,IBAF,IBAF,NY21+1,NY2,1,NZ2,1,LSTEP)
 COVAL(WBE,UC1,COEF,0.0); COVAL(WBE,VC1,COEF,0.0)
 COVAL(WBE,WC1,COEF,0.0)
MPATCH(2,WBW,WWALL,IBAF+1,IBAF+1,NY21+1,NY2,1,NZ2,1,LSTEP)
 COVAL(WBW,UC1,COEF,0.0); COVAL(WBW,VC1,COEF,0.0)
 COVAL(WBW,WC1,COEF,0.0)
    *** Free top surface:
MPATCH(1,POUT1,HIGH,1,NX1,1,NY1,NZ1,NZ1,1,LSTEP)
 COVAL(POUT1, P1,  FIXP, 0.0); COVAL(POUT1,UC1,ONLYMS,SAME)
 COVAL(POUT1,VC1,ONLYMS,SAME); COVAL(POUT1,WC1,ONLYMS, 0.0)
MPATCH(2,POUT2,HIGH,1,NX2,1,NY2,NZ2,NZ2,1,LSTEP)
 COVAL(POUT2, P1,  FIXP, 0.0); COVAL(POUT2,UC1,ONLYMS,SAME)
 COVAL(POUT2,VC1,ONLYMS,SAME); COVAL(POUT2,WC1,ONLYMS, 0.0)
IF(LTURB) THEN
+ COVAL(WS1,KE,GENLAW,GENLAW); COVAL(WS1,EP,GENLAW,GENLAW)
+ COVAL(WS2,KE,GENLAW,GENLAW); COVAL(WS2,EP,GENLAW,GENLAW)
+ COVAL( WH,KE,GENLAW,GENLAW); COVAL( WH,EP,GENLAW,GENLAW)
+ COVAL( WL,KE,GENLAW,GENLAW); COVAL( WL,EP,GENLAW,GENLAW)
+ COVAL(WIE,KE,GENLAW,GENLAW); COVAL(WIE,EP,GENLAW,GENLAW)
+ COVAL(WIW,KE,GENLAW,GENLAW); COVAL(WIW,EP,GENLAW,GENLAW)
+ COVAL( WN,KE,GENLAW,GENLAW); COVAL( WN,EP,GENLAW,GENLAW)
+ COVAL(WLL,KE,GENLAW,GENLAW); COVAL(WLL,EP,GENLAW,GENLAW)
+ COVAL(WBE,KE,GENLAW,GENLAW); COVAL(WBE,EP,GENLAW,GENLAW)
+ COVAL(WBW,KE,GENLAW,GENLAW); COVAL(WBW,EP,GENLAW,GENLAW)
+ COVAL(POUT1,KE,ONLYMS,SAME); COVAL(POUT1,EP,ONLYMS,SAME)
+ COVAL(POUT2,KE,ONLYMS,SAME); COVAL(POUT2,EP,ONLYMS,SAME)
ENDIF
    GROUP 15. Termination of sweeps
LSWEEP = 30; TSTSWP = -1
    GROUP 16. Termination of iterations
SELREF = T; RESFAC = 1.E-5
    GROUP 17. Under-relaxation devices
RELAX(P1,LINRLX,0.5)
    GROUP 19. Data communicated by satellite to GROUND
    * LSG3 = T, activates curvilinearity treatment;
LSG3= T
SPEDAT(SET,GXMONI,TRANSIENT,L,F)
    GROUP 22. Spot-value print-out
IXMON = NX1/2+1; IYMON = NY1/2+1; IZMON = NZ1/2+1
    GROUP 24. Dumps for restarts
    * User can get PHI-files for each step by uncommenting next:
    IDISPA= 1; CSG1  = A
    RESTRT(ALL)