solid-body.cs-32x32x1/code/ini_vel.F

905f660335 Alis*0001 #include "CPP_OPTIONS.h"
                0002
644a731b64 Jean*0003 CBOP
                0004 C     !ROUTINE: INI_VEL
                0005 C     !INTERFACE:
905f660335 Alis*0006       SUBROUTINE INI_VEL( myThid )
                0007
644a731b64 Jean*0008 C     !DESCRIPTION: \bv
                0009 C     *=================================================================
                0010 C     | SUBROUTINE INI_VEL
                0011 C     | o Initialize flow field (either to zero or from input files)
                0012 C     *=================================================================
                0013 C     \ev
905f660335 Alis*0014
644a731b64 Jean*0015 C     !USES:
                0016       IMPLICIT NONE
905f660335 Alis*0017 C     === Global variables ===
                0018 #include "SIZE.h"
                0019 #include "EEPARAMS.h"
                0020 #include "PARAMS.h"
                0021 #include "GRID.h"
                0022 #include "DYNVARS.h"
                0023
644a731b64 Jean*0024 C     !INPUT/OUTPUT PARAMETERS:
aef96ed361 Jean*0025 C     myThid :: my Thread Id number
905f660335 Alis*0026       INTEGER myThid
                0027
aef96ed361 Jean*0028 C     !LOCAL VARIABLES:
                0029 C     i,j,k  :: Loop counters
                0030 C     bi,bj  :: tile indices
                0031       INTEGER i, j, k
                0032       INTEGER bi, bj
                0033       INTEGER ip1, jp1
                0034       _RL psi, omegaprime, fac
                0035 CEOP
905f660335 Alis*0036
aef96ed361 Jean*0037 C     internal function:
                0038       psi(i,j,bi,bj) = fac*fCoriG(i,j,bi,bj)
905f660335 Alis*0039
                0040 C--   Initialise velocity fields to zero
                0041 C
                0042 C     If you want to specify an analytic initial state for the flow
                0043 C     field then customize the following section of code.
                0044 C     It is, however, often easier to generate initial conditions
                0045 C     off-line and read them from input files...
                0046 C
aef96ed361 Jean*0047 C     this relative rotation corresponds to a 5.day rotation period:
                0048       omegaprime = 80. _d 0 / rSphere
                0049       fac = -(rSphere*rSphere)*omegaprime/(2. _d 0*Omega)
905f660335 Alis*0050       DO bj = myByLo(myThid), myByHi(myThid)
                0051        DO bi = myBxLo(myThid), myBxHi(myThid)
                0052         DO k=1,Nr
aef96ed361 Jean*0053          DO j=1-OLy,sNy+OLy
                0054           jp1 = MIN(j+1,sNy+OLy)
                0055           DO i=1-OLx,sNx+OLx
                0056            ip1 = MIN(i+1,sNx+OLx)
                0057            uVel(i,j,k,bi,bj) = 0. _d 0
                0058      &      + (psi(i,j,bi,bj)-psi(i,jp1,bi,bj))*recip_dyG(i,j,bi,bj)
                0059            vVel(i,j,k,bi,bj) = 0. _d 0
                0060      &      + (psi(ip1,j,bi,bj)-psi(i,j,bi,bj))*recip_dxG(i,j,bi,bj)
905f660335 Alis*0061           ENDDO
                0062          ENDDO
                0063         ENDDO
                0064        ENDDO
                0065       ENDDO
                0066
aef96ed361 Jean*0067 c     IF ( uVelInitFile .NE. ' ' .OR. vVelInitFile .NE. ' ' ) THEN
905f660335 Alis*0068 C     Read an initial state for each component if required
aef96ed361 Jean*0069        IF ( uVelInitFile .NE. ' ' )
                0070      & CALL READ_FLD_XYZ_RL( uVelInitFile, ' ', uVel, 0, myThid )
905f660335 Alis*0071
aef96ed361 Jean*0072        IF ( vVelInitFile .NE. ' ' )
                0073      & CALL READ_FLD_XYZ_RL( vVelInitFile, ' ', vVel, 0, myThid )
905f660335 Alis*0074
644a731b64 Jean*0075        CALL EXCH_UV_XYZ_RL(uVel,vVel,.TRUE.,myThid)
                0076 c     ENDIF
                0077
905f660335 Alis*0078       DO bj = myByLo(myThid), myByHi(myThid)
                0079        DO bi = myBxLo(myThid), myBxHi(myThid)
                0080         DO k=1,Nr
aef96ed361 Jean*0081          DO j=1-OLy,sNy+OLy
                0082           DO i=1-OLx,sNx+OLx
                0083            uVel(i,j,k,bi,bj) = uVel(i,j,k,bi,bj)*_maskW(i,j,k,bi,bj)
                0084            vVel(i,j,k,bi,bj) = vVel(i,j,k,bi,bj)*_maskS(i,j,k,bi,bj)
905f660335 Alis*0085           ENDDO
                0086          ENDDO
                0087         ENDDO
                0088        ENDDO
                0089       ENDDO
                0090
                0091       RETURN
                0092       END