press_release/code/update_surf_dr.F

C $Header: /u/gcmpack/MITgcm_contrib/ksnow/press_release/code/update_surf_dr.F,v 1.1 2016/12/16 15:23:18 ksnow Exp $
C $Name:  $

#include "CPP_OPTIONS.h"

CBOP
C     !ROUTINE: UPDATE_SURF_DR
C     !INTERFACE:
      SUBROUTINE UPDATE_SURF_DR( useLatest, myTime, myIter, myThid )

C     !DESCRIPTION: \bv
C     *==========================================================*
C     | SUBROUTINE UPDATE_SURF_DR
C     | o Update the surface-level thickness fraction (hFacC,W,S)
C     |   according to the surface r-position = Non-Linear FrSurf
C     *==========================================================*
C     \ev

C     !USES:
      IMPLICIT NONE
C     == Global variables
#include "SIZE.h"
#include "DYNVARS.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "GRID.h"
#include "SURFACE.h"

C     !INPUT/OUTPUT PARAMETERS:
C     == Routine arguments ==
C     useLatest :: if true use hFac_surfC, else use hFac_surfNm1C
C     myTime    :: Current time in simulation
C     myIter    :: Current iteration number in simulation
C     myThid    :: Thread number for this instance of the routine.
      LOGICAL useLatest
      _RL myTime
      INTEGER myIter
      INTEGER myThid

C     !LOCAL VARIABLES:
#ifdef NONLIN_FRSURF
C     Local variables
C     i,j,bi,bj - loop counter
      INTEGER i,j,k,bi,bj
      INTEGER ks
CEOP

      DO bj=myByLo(myThid), myByHi(myThid)
       DO bi=myBxLo(myThid), myBxHi(myThid)

C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|

       IF (useLatest.AND.(nonlinFreeSurf.GT.0)) then

C-- Update the fractional thickness "hFacC" of the surface level kSurfC :
        DO j=1-Oly,sNy+Oly
         DO i=1-Olx,sNx+Olx
          ks = kSurfC(i,j,bi,bj)
          IF (ks.LE.Nr) THEN
           hFacC(i,j,ks,bi,bj) = hFac_surfC(i,j,bi,bj)
           recip_hFacC(i,j,ks,bi,bj)= 1. _d 0 / hFac_surfC(i,j,bi,bj)
          ENDIF
         ENDDO
        ENDDO

C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
C-- Update fractional thickness "hFacW" & "hFacS" (at U and V points)

        DO j=1-Oly,sNy+Oly
         DO i=2-Olx,sNx+Olx
          ks = kSurfW(i,j,bi,bj)
          IF (ks.LE.Nr) THEN
           hFacW(i,j,ks,bi,bj) = hFac_surfW(i,j,bi,bj)
           IF (hFac_surfW(i,j,bi,bj) .NE. 0.0) THEN
            recip_hFacW(i,j,ks,bi,bj)= 1. _d 0 / hFac_surfW(i,j,bi,bj)
           ELSE
            recip_hFacW(i,j,ks,bi,bj) = 0.0
           ENDIF
          ENDIF
         ENDDO
        ENDDO
        DO j=2-Oly,sNy+Oly
         DO i=1-Olx,sNx+Olx
          ks = kSurfS(i,j,bi,bj)
          IF (ks.LE.Nr) THEN
           hFacS(i,j,ks,bi,bj) = hFac_surfS(i,j,bi,bj)

           IF (hFac_surfS(i,j,bi,bj) .NE. 0.0) THEN
            recip_hFacS(i,j,ks,bi,bj)= 1. _d 0 / hFac_surfS(i,j,bi,bj)
           ELSE
            recip_hFacS(i,j,ks,bi,bj) = 0.0
           ENDIF

          ENDIF
         ENDDO
        ENDDO

       ELSEIF (nonlinFreeSurf.GT.0) THEN

C-- Update the fractional thickness "hFacC" of the surface level kSurfC
C:
        DO j=1-Oly,sNy+Oly
         DO i=1-Olx,sNx+Olx
          ks = kSurfC(i,j,bi,bj)
          IF (ks.LE.Nr) THEN
           hFacC(i,j,ks,bi,bj) = hFac_surfNm1C(i,j,bi,bj)
           recip_hFacC(i,j,ks,bi,bj)= 1. _d 0 / hFac_surfNm1C(i,j,bi,bj)
          ENDIF
         ENDDO
        ENDDO

C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
C-- Update fractional thickness "hFacW" & "hFacS" (at U and V points)

        DO j=1-Oly,sNy+Oly
         DO i=2-Olx,sNx+Olx
          ks = kSurfW(i,j,bi,bj)
          IF (ks.LE.Nr) THEN
           hFacW(i,j,ks,bi,bj) = hFac_surfNm1W(i,j,bi,bj)           
           IF (hFac_surfNm1W(i,j,bi,bj) .NE. 0.0) THEN
            recip_hFacW(i,j,ks,bi,bj)= 1. _d 0 / 
     &       hFac_surfNm1W(i,j,bi,bj)
           ELSE
            recip_hFacW(i,j,ks,bi,bj) = 0.0
           ENDIF
          ENDIF
         ENDDO
        ENDDO
        DO j=2-Oly,sNy+Oly
         DO i=1-Olx,sNx+Olx
          ks = kSurfS(i,j,bi,bj)
          IF (ks.LE.Nr) THEN
           hFacS(i,j,ks,bi,bj) = hFac_surfNm1S(i,j,bi,bj)
           IF (hFac_surfNm1S(i,j,bi,bj) .NE. 0.0) THEN
            recip_hFacS(i,j,ks,bi,bj)= 1. _d 0 / 
     &       hFac_surfNm1S(i,j,bi,bj)
           ELSE
            recip_hFacS(i,j,ks,bi,bj) = 0.0
           ENDIF
          ENDIF
         ENDDO
        ENDDO

       ELSE

        DO k=1,Nr
         DO j=1-Oly,sNy+Oly
          DO i=1-Olx,sNx+Olx
           hFacC(i,j,k,bi,bj)=h0FacC(i,j,k,bi,bj)
           IF (h0FacC(i,j,k,bi,bj) .NE. 0. ) THEN
            recip_hFacC(i,j,k,bi,bj) = 1. _d 0 / h0FacC(i,j,k,bi,bj)
           ELSE
            recip_hFacC(i,j,k,bi,bj) = 0.
           ENDIF
          ENDDO
         ENDDO
        ENDDO

       ENDIF

C- end bi,bj loop
       ENDDO
      ENDDO

#ifdef ALLOW_PRESSURE_RELEASE_CODE
      DO bj=myByLo(myThid), myByHi(myThid)
       DO bi=myBxLo(myThid), myBxHi(myThid)
        DO J=1-OLy,sNy+OLy
         DO I=1-OLx,sNx+OLx

          depthColW(i,j,bi,bj) = 0.0
          depthColS(i,j,bi,bj) = 0.0

         ENDDO
        ENDDO
       ENDDO
      ENDDO
      DO K=1, Nr
       DO bj=myByLo(myThid), myByHi(myThid)
        DO bi=myBxLo(myThid), myBxHi(myThid)
         DO J=1-OLy,sNy+OLy
          DO I=1-OLx,sNx+OLx

           depthColW(i,j,bi,bj) = 
     &      depthColW(i,j,bi,bj) + 
     &      drF(k) * hFacW(i,j,k,bi,bj) 
           depthColS(i,j,bi,bj) = 
     &      depthColS(i,j,bi,bj) + 
     &      drF(k) * hFacS(i,j,k,bi,bj)

          ENDDO
         ENDDO
        ENDDO
       ENDDO
      ENDDO

      _EXCH_XY_RL (depthColW, myThid)
      _EXCH_XY_RL (depthColS, myThid)
#endif

C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
#endif /* NONLIN_FRSURF */

       

      RETURN
      END
1	ksnow	1.2	C $Header: /u/gcmpack/MITgcm_contrib/ksnow/press_release/code/update_surf_dr.F,v 1.1 2016/12/16 15:23:18 ksnow Exp $
2	ksnow	1.1	C $Name: $
3
4			#include "CPP_OPTIONS.h"
5
6			CBOP
7			C !ROUTINE: UPDATE_SURF_DR
8			C !INTERFACE:
9			SUBROUTINE UPDATE_SURF_DR( useLatest, myTime, myIter, myThid )
10
11			C !DESCRIPTION: \bv
12			C ==========================================================
13			C \| SUBROUTINE UPDATE_SURF_DR
14			C \| o Update the surface-level thickness fraction (hFacC,W,S)
15			C \| according to the surface r-position = Non-Linear FrSurf
16			C ==========================================================
17			C \ev
18
19			C !USES:
20			IMPLICIT NONE
21			C == Global variables
22			#include "SIZE.h"
23	ksnow	1.2	#include "DYNVARS.h"
24	ksnow	1.1	#include "EEPARAMS.h"
25			#include "PARAMS.h"
26			#include "GRID.h"
27			#include "SURFACE.h"
28
29			C !INPUT/OUTPUT PARAMETERS:
30			C == Routine arguments ==
31			C useLatest :: if true use hFac_surfC, else use hFac_surfNm1C
32			C myTime :: Current time in simulation
33			C myIter :: Current iteration number in simulation
34			C myThid :: Thread number for this instance of the routine.
35			LOGICAL useLatest
36			_RL myTime
37			INTEGER myIter
38			INTEGER myThid
39
40			C !LOCAL VARIABLES:
41			#ifdef NONLIN_FRSURF
42			C Local variables
43			C i,j,bi,bj - loop counter
44			INTEGER i,j,k,bi,bj
45			INTEGER ks
46			CEOP
47
48			DO bj=myByLo(myThid), myByHi(myThid)
49			DO bi=myBxLo(myThid), myBxHi(myThid)
50
51			C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-\|--+----\|
52
53			IF (useLatest.AND.(nonlinFreeSurf.GT.0)) then
54
55			C-- Update the fractional thickness "hFacC" of the surface level kSurfC :
56			DO j=1-Oly,sNy+Oly
57			DO i=1-Olx,sNx+Olx
58			ks = kSurfC(i,j,bi,bj)
59			IF (ks.LE.Nr) THEN
60			hFacC(i,j,ks,bi,bj) = hFac_surfC(i,j,bi,bj)
61			recip_hFacC(i,j,ks,bi,bj)= 1. _d 0 / hFac_surfC(i,j,bi,bj)
62			ENDIF
63			ENDDO
64			ENDDO
65
66			C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-\|--+----\|
67			C-- Update fractional thickness "hFacW" & "hFacS" (at U and V points)
68
69			DO j=1-Oly,sNy+Oly
70			DO i=2-Olx,sNx+Olx
71			ks = kSurfW(i,j,bi,bj)
72			IF (ks.LE.Nr) THEN
73			hFacW(i,j,ks,bi,bj) = hFac_surfW(i,j,bi,bj)
74			IF (hFac_surfW(i,j,bi,bj) .NE. 0.0) THEN
75			recip_hFacW(i,j,ks,bi,bj)= 1. _d 0 / hFac_surfW(i,j,bi,bj)
76			ELSE
77			recip_hFacW(i,j,ks,bi,bj) = 0.0
78			ENDIF
79			ENDIF
80			ENDDO
81			ENDDO
82			DO j=2-Oly,sNy+Oly
83			DO i=1-Olx,sNx+Olx
84			ks = kSurfS(i,j,bi,bj)
85			IF (ks.LE.Nr) THEN
86			hFacS(i,j,ks,bi,bj) = hFac_surfS(i,j,bi,bj)
87	ksnow	1.2
88	ksnow	1.1	IF (hFac_surfS(i,j,bi,bj) .NE. 0.0) THEN
89			recip_hFacS(i,j,ks,bi,bj)= 1. _d 0 / hFac_surfS(i,j,bi,bj)
90			ELSE
91			recip_hFacS(i,j,ks,bi,bj) = 0.0
92			ENDIF
93	ksnow	1.2
94	ksnow	1.1	ENDIF
95			ENDDO
96			ENDDO
97
98			ELSEIF (nonlinFreeSurf.GT.0) THEN
99
100			C-- Update the fractional thickness "hFacC" of the surface level kSurfC
101			C:
102			DO j=1-Oly,sNy+Oly
103			DO i=1-Olx,sNx+Olx
104			ks = kSurfC(i,j,bi,bj)
105			IF (ks.LE.Nr) THEN
106			hFacC(i,j,ks,bi,bj) = hFac_surfNm1C(i,j,bi,bj)
107			recip_hFacC(i,j,ks,bi,bj)= 1. _d 0 / hFac_surfNm1C(i,j,bi,bj)
108			ENDIF
109			ENDDO
110			ENDDO
111
112			C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-\|--+----\|
113			C-- Update fractional thickness "hFacW" & "hFacS" (at U and V points)
114
115			DO j=1-Oly,sNy+Oly
116			DO i=2-Olx,sNx+Olx
117			ks = kSurfW(i,j,bi,bj)
118			IF (ks.LE.Nr) THEN
119			hFacW(i,j,ks,bi,bj) = hFac_surfNm1W(i,j,bi,bj)
120			IF (hFac_surfNm1W(i,j,bi,bj) .NE. 0.0) THEN
121			recip_hFacW(i,j,ks,bi,bj)= 1. _d 0 /
122			& hFac_surfNm1W(i,j,bi,bj)
123			ELSE
124			recip_hFacW(i,j,ks,bi,bj) = 0.0
125			ENDIF
126			ENDIF
127			ENDDO
128			ENDDO
129			DO j=2-Oly,sNy+Oly
130			DO i=1-Olx,sNx+Olx
131			ks = kSurfS(i,j,bi,bj)
132			IF (ks.LE.Nr) THEN
133			hFacS(i,j,ks,bi,bj) = hFac_surfNm1S(i,j,bi,bj)
134			IF (hFac_surfNm1S(i,j,bi,bj) .NE. 0.0) THEN
135			recip_hFacS(i,j,ks,bi,bj)= 1. _d 0 /
136			& hFac_surfNm1S(i,j,bi,bj)
137			ELSE
138			recip_hFacS(i,j,ks,bi,bj) = 0.0
139			ENDIF
140			ENDIF
141			ENDDO
142			ENDDO
143
144			ELSE
145
146			DO k=1,Nr
147			DO j=1-Oly,sNy+Oly
148			DO i=1-Olx,sNx+Olx
149			hFacC(i,j,k,bi,bj)=h0FacC(i,j,k,bi,bj)
150			IF (h0FacC(i,j,k,bi,bj) .NE. 0. ) THEN
151			recip_hFacC(i,j,k,bi,bj) = 1. _d 0 / h0FacC(i,j,k,bi,bj)
152			ELSE
153			recip_hFacC(i,j,k,bi,bj) = 0.
154			ENDIF
155			ENDDO
156			ENDDO
157			ENDDO
158
159			ENDIF
160
161			C- end bi,bj loop
162			ENDDO
163			ENDDO
164
165			#ifdef ALLOW_PRESSURE_RELEASE_CODE
166			DO bj=myByLo(myThid), myByHi(myThid)
167			DO bi=myBxLo(myThid), myBxHi(myThid)
168			DO J=1-OLy,sNy+OLy
169			DO I=1-OLx,sNx+OLx
170
171			depthColW(i,j,bi,bj) = 0.0
172			depthColS(i,j,bi,bj) = 0.0
173
174			ENDDO
175			ENDDO
176			ENDDO
177			ENDDO
178			DO K=1, Nr
179			DO bj=myByLo(myThid), myByHi(myThid)
180			DO bi=myBxLo(myThid), myBxHi(myThid)
181			DO J=1-OLy,sNy+OLy
182			DO I=1-OLx,sNx+OLx
183
184			depthColW(i,j,bi,bj) =
185			& depthColW(i,j,bi,bj) +
186			& drF(k) * hFacW(i,j,k,bi,bj)
187			depthColS(i,j,bi,bj) =
188			& depthColS(i,j,bi,bj) +
189			& drF(k) * hFacS(i,j,k,bi,bj)
190
191			ENDDO
192			ENDDO
193			ENDDO
194			ENDDO
195			ENDDO
196
197			_EXCH_XY_RL (depthColW, myThid)
198			_EXCH_XY_RL (depthColS, myThid)
199			#endif
200
201			C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-\|--+----\|
202			#endif /* NONLIN_FRSURF */
203
204
205
206			RETURN
207			END