model/src/ini_masks_etc.F

C $Header: /u/gcmpack/models/MITgcmUV/model/src/ini_masks_etc.F,v 1.4 1998/08/22 17:51:08 cnh Exp $

#include "CPP_EEOPTIONS.h"

CStartOfInterface
      SUBROUTINE INI_MASKS_ETC( myThid )
C     /==========================================================\
C     | SUBROUTINE INI_MASKS_ETC                                 |
C     | o Initialise masks and topography factors                |
C     |==========================================================|
C     | These arrays are used throughout the code and describe   |
C     | the topography of the domain through masks (0s and 1s)   |
C     | and fractional height factors (0<hFac<1). The latter     |
C     | distinguish between the lopped-cell and full-step        |
C     | topographic representations.                             |
C     \==========================================================/

C     === Global variables ===
#include "SIZE.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "GRID.h"

C     == Routine arguments ==
C     myThid -  Number of this instance of INI_CARTESIAN_GRID
      INTEGER myThid
CEndOfInterface

C     == Local variables ==
C     bi,bj  - Loop counters
C     I,J,K
      INTEGER bi, bj
      INTEGER  I, J, K

C     Calculate quantities derived from XY depth map
      DO bj = myByLo(myThid), myByHi(myThid)
       DO bi = myBxLo(myThid), myBxHi(myThid)
        DO J=1,sNy
         DO I=1,sNx
C         Inverse of depth
          IF ( h(i,j,bi,bj) .EQ. 0. _d 0 ) THEN
           recip_H(i,j,bi,bj) = 0. _d 0
          ELSE
           recip_H(i,j,bi,bj) = 1. _d 0 /  abs( H(i,j,bi,bj) )
          ENDIF
         ENDDO
        ENDDO
       ENDDO
      ENDDO
      _EXCH_XY_R4(   recip_H, myThid )

C     Calculate lopping factor hFacC
      DO bj=myByLo(myThid), myByHi(myThid)
       DO bi=myBxLo(myThid), myBxHi(myThid)
        DO K=1, Nr
         DO J=1,sNy
          DO I=1,sNx
           IF     ( H(I,J,bi,bj) .GE. rF(K) ) THEN
C           Top of cell is below base of domain
            hFacC(I,J,K,bi,bj) = 0.
           ELSEIF ( H(I,J,bi,bj) .LE. rF(K+1) ) THEN
C           Base of domain is below bottom of this cell
            hFacC(I,J,K,bi,bj) = 1.
           ELSE
C           Base of domain is in this cell
C           Set hFac to the fraction of the cell that is open.
            hFacC(I,J,K,bi,bj) = (rF(K)-H(I,J,bi,bj))*recip_drF(K)
           ENDIF
C Impose minimum fraction
           IF (hFacC(I,J,K,bi,bj).LT.hFacMin) THEN
            IF (hFacC(I,J,K,bi,bj).LT.hFacMin*0.5) THEN
             hFacC(I,J,K,bi,bj)=0.
            ELSE
             hFacC(I,J,K,bi,bj)=hFacMin
            ENDIF
           ENDIF
C Impose minimum size (dimensional)
           IF (drF(k)*hFacC(I,J,K,bi,bj).LT.hFacMinDr) THEN
            IF (drF(k)*hFacC(I,J,K,bi,bj).LT.hFacMinDr*0.5) THEN
             hFacC(I,J,K,bi,bj)=0.
            ELSE
             hFacC(I,J,K,bi,bj)=hFacMinDr*recip_drF(k)
            ENDIF
           ENDIF
          ENDDO
         ENDDO
        ENDDO
       ENDDO
      ENDDO
      _EXCH_XYZ_R4(hFacC , myThid )

C     hFacW and hFacS (at U and V points)
      DO bj=myByLo(myThid), myByHi(myThid)
       DO bi=myBxLo(myThid), myBxHi(myThid)
        DO K=1, Nr
         DO J=1,sNy
          DO I=1,sNx
           hFacW(I,J,K,bi,bj)=
     &       MIN(hFacC(I,J,K,bi,bj),hFacC(I-1,J,K,bi,bj))
           hFacS(I,J,K,bi,bj)=
     &       MIN(hFacC(I,J,K,bi,bj),hFacC(I,J-1,K,bi,bj))
          ENDDO
         ENDDO
        ENDDO
       ENDDO
      ENDDO
      _EXCH_XYZ_R4(hFacW , myThid )
      _EXCH_XYZ_R4(hFacS , myThid )

C     Masks and reciprocals of hFac[CWS]
      DO bj = myByLo(myThid), myByHi(myThid)
       DO bi = myBxLo(myThid), myBxHi(myThid)
        DO K=1,Nr
         DO J=1,sNy
          DO I=1,sNx
           IF (HFacC(I,J,K,bi,bj) .NE. 0. D0 ) THEN
            recip_HFacC(I,J,K,bi,bj) = 1. D0 / HFacC(I,J,K,bi,bj)
           ELSE
            recip_HFacC(I,J,K,bi,bj) = 0. D0
           ENDIF
           IF (HFacW(I,J,K,bi,bj) .NE. 0. D0 ) THEN
            recip_HFacW(I,J,K,bi,bj) = 1. D0 / HFacW(I,J,K,bi,bj)
            maskW(I,J,K,bi,bj) = 1. D0
           ELSE
            recip_HFacW(I,J,K,bi,bj) = 0. D0
            maskW(I,J,K,bi,bj) = 0.0 D0
           ENDIF
           IF (HFacS(I,J,K,bi,bj) .NE. 0. D0 ) THEN
            recip_HFacS(I,J,K,bi,bj) = 1. D0 / HFacS(I,J,K,bi,bj)
            maskS(I,J,K,bi,bj) = 1. D0
           ELSE
            recip_HFacS(I,J,K,bi,bj) = 0. D0
            maskS(I,J,K,bi,bj) = 0. D0
           ENDIF
          ENDDO
         ENDDO
        ENDDO
       ENDDO
      ENDDO
      _EXCH_XYZ_R4(recip_HFacC    , myThid )
      _EXCH_XYZ_R4(recip_HFacW    , myThid )
      _EXCH_XYZ_R4(recip_HFacS    , myThid )
      _EXCH_XYZ_R4(maskW    , myThid )
      _EXCH_XYZ_R4(maskS    , myThid )

C     Calculate recipricols grid lengths
      DO bj = myByLo(myThid), myByHi(myThid)
       DO bi = myBxLo(myThid), myBxHi(myThid)
        DO J=1,sNy
         DO I=1,sNx
          recip_dxG(I,J,bi,bj)=1.d0/dxG(I,J,bi,bj)
          recip_dyG(I,J,bi,bj)=1.d0/dyG(I,J,bi,bj)
          recip_dxC(I,J,bi,bj)=1.d0/dxC(I,J,bi,bj)
          recip_dyC(I,J,bi,bj)=1.d0/dyC(I,J,bi,bj)
          recip_dxF(I,J,bi,bj)=1.d0/dxF(I,J,bi,bj)
          recip_dyF(I,J,bi,bj)=1.d0/dyF(I,J,bi,bj)
          recip_dxV(I,J,bi,bj)=1.d0/dxV(I,J,bi,bj)
          recip_dyU(I,J,bi,bj)=1.d0/dyU(I,J,bi,bj)
         ENDDO
        ENDDO
       ENDDO
      ENDDO
      _EXCH_XY_R4(recip_dxG, myThid )
      _EXCH_XY_R4(recip_dyG, myThid )
      _EXCH_XY_R4(recip_dxC, myThid )
      _EXCH_XY_R4(recip_dyC, myThid )
      _EXCH_XY_R4(recip_dxF, myThid )
      _EXCH_XY_R4(recip_dyF, myThid )
      _EXCH_XY_R4(recip_dxV, myThid )
      _EXCH_XY_R4(recip_dyU, myThid )

C
      RETURN
      END
1	cnh	1.5	C $Header: /u/gcmpack/models/MITgcmUV/model/src/ini_masks_etc.F,v 1.4 1998/08/22 17:51:08 cnh Exp $
2	adcroft	1.1
3			#include "CPP_EEOPTIONS.h"
4
5			CStartOfInterface
6			SUBROUTINE INI_MASKS_ETC( myThid )
7			C /==========================================================\
8			C \| SUBROUTINE INI_MASKS_ETC \|
9			C \| o Initialise masks and topography factors \|
10			C \|==========================================================\|
11			C \| These arrays are used throughout the code and describe \|
12			C \| the topography of the domain through masks (0s and 1s) \|
13			C \| and fractional height factors (0<hFac<1). The latter \|
14			C \| distinguish between the lopped-cell and full-step \|
15			C \| topographic representations. \|
16			C \==========================================================/
17
18			C === Global variables ===
19			#include "SIZE.h"
20			#include "EEPARAMS.h"
21			#include "PARAMS.h"
22			#include "GRID.h"
23
24			C == Routine arguments ==
25			C myThid - Number of this instance of INI_CARTESIAN_GRID
26			INTEGER myThid
27			CEndOfInterface
28
29			C == Local variables ==
30			C bi,bj - Loop counters
31			C I,J,K
32			INTEGER bi, bj
33			INTEGER I, J, K
34
35			C Calculate quantities derived from XY depth map
36			DO bj = myByLo(myThid), myByHi(myThid)
37			DO bi = myBxLo(myThid), myBxHi(myThid)
38			DO J=1,sNy
39			DO I=1,sNx
40			C Inverse of depth
41			IF ( h(i,j,bi,bj) .EQ. 0. _d 0 ) THEN
42	cnh	1.4	recip_H(i,j,bi,bj) = 0. _d 0
43	adcroft	1.1	ELSE
44	cnh	1.4	recip_H(i,j,bi,bj) = 1. _d 0 / abs( H(i,j,bi,bj) )
45	adcroft	1.1	ENDIF
46			ENDDO
47			ENDDO
48			ENDDO
49			ENDDO
50	cnh	1.4	_EXCH_XY_R4( recip_H, myThid )
51	adcroft	1.2
52			C Calculate lopping factor hFacC
53			DO bj=myByLo(myThid), myByHi(myThid)
54			DO bi=myBxLo(myThid), myBxHi(myThid)
55	cnh	1.4	DO K=1, Nr
56	adcroft	1.2	DO J=1,sNy
57			DO I=1,sNx
58	cnh	1.4	IF ( H(I,J,bi,bj) .GE. rF(K) ) THEN
59	adcroft	1.3	C Top of cell is below base of domain
60	adcroft	1.2	hFacC(I,J,K,bi,bj) = 0.
61	cnh	1.4	ELSEIF ( H(I,J,bi,bj) .LE. rF(K+1) ) THEN
62	adcroft	1.3	C Base of domain is below bottom of this cell
63	adcroft	1.2	hFacC(I,J,K,bi,bj) = 1.
64			ELSE
65			C Base of domain is in this cell
66	adcroft	1.3	C Set hFac to the fraction of the cell that is open.
67	cnh	1.4	hFacC(I,J,K,bi,bj) = (rF(K)-H(I,J,bi,bj))*recip_drF(K)
68	adcroft	1.3	ENDIF
69			C Impose minimum fraction
70			IF (hFacC(I,J,K,bi,bj).LT.hFacMin) THEN
71			IF (hFacC(I,J,K,bi,bj).LT.hFacMin*0.5) THEN
72			hFacC(I,J,K,bi,bj)=0.
73			ELSE
74			hFacC(I,J,K,bi,bj)=hFacMin
75			ENDIF
76			ENDIF
77			C Impose minimum size (dimensional)
78	cnh	1.5	IF (drF(k)*hFacC(I,J,K,bi,bj).LT.hFacMinDr) THEN
79			IF (drF(k)hFacC(I,J,K,bi,bj).LT.hFacMinDr0.5) THEN
80	adcroft	1.3	hFacC(I,J,K,bi,bj)=0.
81			ELSE
82	cnh	1.5	hFacC(I,J,K,bi,bj)=hFacMinDr*recip_drF(k)
83	adcroft	1.3	ENDIF
84	adcroft	1.2	ENDIF
85			ENDDO
86			ENDDO
87			ENDDO
88			ENDDO
89			ENDDO
90			_EXCH_XYZ_R4(hFacC , myThid )
91	adcroft	1.1
92			C hFacW and hFacS (at U and V points)
93			DO bj=myByLo(myThid), myByHi(myThid)
94			DO bi=myBxLo(myThid), myBxHi(myThid)
95	cnh	1.4	DO K=1, Nr
96	adcroft	1.1	DO J=1,sNy
97			DO I=1,sNx
98			hFacW(I,J,K,bi,bj)=
99			& MIN(hFacC(I,J,K,bi,bj),hFacC(I-1,J,K,bi,bj))
100			hFacS(I,J,K,bi,bj)=
101			& MIN(hFacC(I,J,K,bi,bj),hFacC(I,J-1,K,bi,bj))
102			ENDDO
103			ENDDO
104			ENDDO
105			ENDDO
106			ENDDO
107			_EXCH_XYZ_R4(hFacW , myThid )
108			_EXCH_XYZ_R4(hFacS , myThid )
109
110			C Masks and reciprocals of hFac[CWS]
111			DO bj = myByLo(myThid), myByHi(myThid)
112			DO bi = myBxLo(myThid), myBxHi(myThid)
113	cnh	1.4	DO K=1,Nr
114	adcroft	1.1	DO J=1,sNy
115			DO I=1,sNx
116			IF (HFacC(I,J,K,bi,bj) .NE. 0. D0 ) THEN
117	cnh	1.4	recip_HFacC(I,J,K,bi,bj) = 1. D0 / HFacC(I,J,K,bi,bj)
118	adcroft	1.1	ELSE
119	cnh	1.4	recip_HFacC(I,J,K,bi,bj) = 0. D0
120	adcroft	1.1	ENDIF
121			IF (HFacW(I,J,K,bi,bj) .NE. 0. D0 ) THEN
122	cnh	1.4	recip_HFacW(I,J,K,bi,bj) = 1. D0 / HFacW(I,J,K,bi,bj)
123	adcroft	1.1	maskW(I,J,K,bi,bj) = 1. D0
124			ELSE
125	cnh	1.4	recip_HFacW(I,J,K,bi,bj) = 0. D0
126	adcroft	1.1	maskW(I,J,K,bi,bj) = 0.0 D0
127			ENDIF
128			IF (HFacS(I,J,K,bi,bj) .NE. 0. D0 ) THEN
129	cnh	1.4	recip_HFacS(I,J,K,bi,bj) = 1. D0 / HFacS(I,J,K,bi,bj)
130	adcroft	1.1	maskS(I,J,K,bi,bj) = 1. D0
131			ELSE
132	cnh	1.4	recip_HFacS(I,J,K,bi,bj) = 0. D0
133	adcroft	1.1	maskS(I,J,K,bi,bj) = 0. D0
134			ENDIF
135			ENDDO
136			ENDDO
137			ENDDO
138			ENDDO
139			ENDDO
140	cnh	1.4	_EXCH_XYZ_R4(recip_HFacC , myThid )
141			_EXCH_XYZ_R4(recip_HFacW , myThid )
142			_EXCH_XYZ_R4(recip_HFacS , myThid )
143	adcroft	1.1	_EXCH_XYZ_R4(maskW , myThid )
144			_EXCH_XYZ_R4(maskS , myThid )
145
146			C Calculate recipricols grid lengths
147			DO bj = myByLo(myThid), myByHi(myThid)
148			DO bi = myBxLo(myThid), myBxHi(myThid)
149			DO J=1,sNy
150			DO I=1,sNx
151	cnh	1.4	recip_dxG(I,J,bi,bj)=1.d0/dxG(I,J,bi,bj)
152			recip_dyG(I,J,bi,bj)=1.d0/dyG(I,J,bi,bj)
153			recip_dxC(I,J,bi,bj)=1.d0/dxC(I,J,bi,bj)
154			recip_dyC(I,J,bi,bj)=1.d0/dyC(I,J,bi,bj)
155			recip_dxF(I,J,bi,bj)=1.d0/dxF(I,J,bi,bj)
156			recip_dyF(I,J,bi,bj)=1.d0/dyF(I,J,bi,bj)
157			recip_dxV(I,J,bi,bj)=1.d0/dxV(I,J,bi,bj)
158			recip_dyU(I,J,bi,bj)=1.d0/dyU(I,J,bi,bj)
159	adcroft	1.1	ENDDO
160			ENDDO
161			ENDDO
162			ENDDO
163	cnh	1.4	_EXCH_XY_R4(recip_dxG, myThid )
164			_EXCH_XY_R4(recip_dyG, myThid )
165			_EXCH_XY_R4(recip_dxC, myThid )
166			_EXCH_XY_R4(recip_dyC, myThid )
167			_EXCH_XY_R4(recip_dxF, myThid )
168			_EXCH_XY_R4(recip_dyF, myThid )
169			_EXCH_XY_R4(recip_dxV, myThid )
170			_EXCH_XY_R4(recip_dyU, myThid )
171	adcroft	1.1
172			C
173			RETURN
174			END