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

aef96ed361 Jean*0001 #include "PACKAGES_CONFIG.h"
905f660335 Alis*0002 #include "CPP_OPTIONS.h"
                0003
aef96ed361 Jean*0004 CBOP
                0005 C     !ROUTINE: INI_PSURF
                0006 C     !INTERFACE:
905f660335 Alis*0007       SUBROUTINE INI_PSURF( myThid )
aef96ed361 Jean*0008
                0009 C     !DESCRIPTION: \bv
                0010 C     *==========================================================*
905f660335 Alis*0011 C     | SUBROUTINE INI_PSURF                                     |
                0012 C     | o Set model initial free-surface height/pressure.        |
aef96ed361 Jean*0013 C     *==========================================================*
905f660335 Alis*0014 C     | There are several options for setting the initial        |
                0015 C     | surface displacement (r unit) field.                     |
                0016 C     |  1. Inline code                                          |
                0017 C     |  2. Two-dimensional data from a file.                    |
aef96ed361 Jean*0018 C     *==========================================================*
                0019 C     \ev
905f660335 Alis*0020
aef96ed361 Jean*0021 C     !USES:
                0022       IMPLICIT NONE
905f660335 Alis*0023 C     === Global variables ===
                0024 #include "SIZE.h"
                0025 #include "EEPARAMS.h"
                0026 #include "PARAMS.h"
                0027 #include "GRID.h"
                0028 #include "DYNVARS.h"
aef96ed361 Jean*0029 #include "SURFACE.h"
                0030 #ifdef ALLOW_CD_CODE
                0031 # include "CD_CODE_VARS.h"
                0032 #endif
                0033 #ifdef ALLOW_SHELFICE
                0034 # include "SHELFICE.h"
                0035 #endif
905f660335 Alis*0036
aef96ed361 Jean*0037 C     !INPUT/OUTPUT PARAMETERS:
                0038 C     myThid :: my Thread Id number
905f660335 Alis*0039       INTEGER myThid
                0040
aef96ed361 Jean*0041 C     !LOCAL VARIABLES:
                0042 C     bi,bj  :: tiles indices
                0043 C     i, j   :: Loop counters
905f660335 Alis*0044       INTEGER bi, bj
aef96ed361 Jean*0045       INTEGER i, j
                0046       _RL omegaPrime, snFac, psFac
                0047 CEOP
905f660335 Alis*0048
                0049 C--   Initialise surface position anomaly to zero
aef96ed361 Jean*0050       omegaPrime = 80. _d 0 / rSphere
                0051 c     psFac = -(rSphere**2)*omegaPrime*(Omega+omegaPrime)
                0052 C     previous expression above was missing one "half" factor:
6b4f212f67 Jean*0053       psfac = -(rSphere*rSphere)*omegaPrime
                0054      &                          *( Omega + omegaPrime*0.5 _d 0 )
aef96ed361 Jean*0055       snFac = 1. _d 0 / (4. _d 0*Omega*Omega)
6b4f212f67 Jean*0056 C     recip_Bo is constant here, Z-coord: =1/gBaro
                0057 C                                P-coord: if uniformLin_PhiSurf=T, =rhoConst
905f660335 Alis*0058       DO bj = myByLo(myThid), myByHi(myThid)
                0059        DO bi = myBxLo(myThid), myBxHi(myThid)
aef96ed361 Jean*0060         DO j=1-OLy,sNy+OLy
                0061          DO i=1-OLx,sNx+OLx
                0062           etaN(i,j,bi,bj) = 0. _d 0
                0063      &     + psFac*( snFac*fCori(i,j,bi,bj)*fCori(i,j,bi,bj)
6b4f212f67 Jean*0064      &             - 1. _d 0 / 3. _d 0 )*recip_Bo(i,j,bi,bj)
905f660335 Alis*0065          ENDDO
                0066         ENDDO
                0067        ENDDO
                0068       ENDDO
                0069 C     Read an initial state
aef96ed361 Jean*0070       IF ( pSurfInitFile .NE. ' ' ) THEN
905f660335 Alis*0071        CALL READ_FLD_XY_RL( pSurfInitFile, ' ', etaN, 0, myThid )
aef96ed361 Jean*0072 C      fill the overlap (+ BARRIER)
                0073        _EXCH_XY_RL(etaN, myThid)
905f660335 Alis*0074       ENDIF
                0075
138482fdf6 Ed H*0076 #ifdef ALLOW_CD_CODE
aef96ed361 Jean*0077 C--   By default, initialize etaNm1 with etaN :
905f660335 Alis*0078       DO bj=myByLo(myThid),myByHi(myThid)
                0079        DO bi=myBxLo(myThid),myBxHi(myThid)
aef96ed361 Jean*0080         DO j=1-OLy,sNy+OLy
                0081          DO i=1-OLx,sNx+OLx
                0082           etaNm1(i,j,bi,bj) = etaN(i,j,bi,bj)
905f660335 Alis*0083          ENDDO
                0084         ENDDO
                0085        ENDDO
                0086       ENDDO
                0087 #endif
                0088
aef96ed361 Jean*0089 #ifdef EXACT_CONSERV
                0090 C--   By default, initialize etaH with etaN :
                0091       DO bj=myByLo(myThid),myByHi(myThid)
                0092        DO bi=myBxLo(myThid),myBxHi(myThid)
                0093         DO j=1-OLy,sNy+OLy
                0094          DO i=1-OLx,sNx+OLx
                0095           etaH(i,j,bi,bj) = etaN(i,j,bi,bj)
                0096           etaHnm1(i,j,bi,bj) = etaN(i,j,bi,bj)
                0097           dEtaHdt(i,j,bi,bj) = 0. _d 0
                0098          ENDDO
                0099         ENDDO
                0100        ENDDO
                0101       ENDDO
                0102 #endif /* EXACT_CONSERV */
                0103
                0104 #ifdef ALLOW_SHELFICE
                0105       IF ( useShelfIce .AND. usingZCoords ) THEN
                0106         DO bj=myByLo(myThid),myByHi(myThid)
                0107          DO bi=myBxLo(myThid),myBxHi(myThid)
                0108            DO j=1-OLy,sNy+OLy
                0109             DO i=1-OLx,sNx+OLx
                0110               phi0surf(i,j,bi,bj) = phi0surf(i,j,bi,bj)
                0111      &          + shelficeLoadAnomaly(i,j,bi,bj)*recip_rhoConst
                0112             ENDDO
                0113            ENDDO
                0114          ENDDO
                0115         ENDDO
                0116       ENDIF
                0117 #endif /* ALLOW_SHELFICE */
                0118
905f660335 Alis*0119       RETURN
                0120       END