/[MITgcm]/MITgcm/model/src/calc_surf_dr.F
ViewVC logotype

Contents of /MITgcm/model/src/calc_surf_dr.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.1 - (show annotations) (download)
Mon Aug 27 18:50:10 2001 UTC (22 years, 9 months ago) by jmc
Branch: MAIN
New routines relative to NonLin-FreeSurf

1 C $Header: $
2 C $Name: $
3
4 #include "CPP_OPTIONS.h"
5
6 SUBROUTINE CALC_SURF_DR(etaFld,
7 I myTime, myIter, myThid )
8 C /==========================================================\
9 C | SUBROUTINE CALC_SURF_DR |
10 C | o Calculate the new surface level thickness according to |
11 C | the surface r-position (Non-Linear Free-Surf) |
12 C | o take decision if grid box becomes too thin or too thick|
13 C \==========================================================/
14 IMPLICIT NONE
15
16 C == Global variables
17 #include "SIZE.h"
18 #include "EEPARAMS.h"
19 #include "PARAMS.h"
20 c #include "DYNVARS.h"
21 #include "GRID.h"
22 #include "SURFACE.h"
23
24 C == Routine arguments ==
25 C myTime - Current time in simulation
26 C myIter - Current iteration number in simulation
27 C myThid - Thread number for this instance of the routine.
28 C etaFld - current eta field used to update the hFactor
29 _RL myTime
30 INTEGER myIter
31 INTEGER myThid
32 _RL etaFld(1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy)
33
34 #ifdef NONLIN_FRSURF
35
36 C Local variables
37 C i,j,k,bi,bj - loop counter
38 INTEGER i,j,k,bi,bj
39 INTEGER ks
40 _RL hFactmp
41 _RS hhm, hhp
42 CHARACTER*(MAX_LEN_MBUF) suff
43
44 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
45
46 C-- Initialise arrays :
47 IF (myIter.EQ.nIter0) THEN
48 DO bj=myByLo(myThid), myByHi(myThid)
49 DO bi=myBxLo(myThid), myBxHi(myThid)
50 DO j=1-Oly,sNy+Oly
51 DO i=1-Olx,sNx+Olx
52 hFac_surfC(i,j,bi,bj) = 0.
53 hFac_surfW(i,j,bi,bj) = 0.
54 hFac_surfS(i,j,bi,bj) = 0.
55 ENDDO
56 ENDDO
57 ENDDO
58 ENDDO
59 ENDIF
60
61 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
62
63 C-- Compute the new fractional thickness of surface level (ksurfC):
64 DO bj=myByLo(myThid), myByHi(myThid)
65 DO bi=myBxLo(myThid), myBxHi(myThid)
66 c DO j=1-Oly,sNy+Oly
67 c DO i=1-Olx,sNx+Olx
68 DO j=1,sNy
69 DO i=1,sNx
70 ks = ksurfC(i,j,bi,bj)
71 IF (ks.LE.Nr) THEN
72 hFac_surfC(i,j,bi,bj) = hFacC(i,j,ks,bi,bj)
73 hFactmp = ( Ro_surf(i,j,bi,bj)+etaFld(i,j,bi,bj)
74 & -MAX( rF(ks+1), R_low(i,j,bi,bj) )
75 & )*recip_drF(ks)
76 IF ( hFactmp.GE.hFacInf .AND. hFactmp.LE.hFacSup ) THEN
77 hFac_surfC(i,j,bi,bj) = maskC(i,j,ks,bi,bj)*hFactmp
78 ELSE
79 C---------- needs to do something :
80 IF (hFactmp.LT.hFacInf) THEN
81 write(0,'(2A,6I4,I10)') 'WARNING: hFacC < hFacInf at:',
82 & ' i,j,k,bi,bj,Thid,Iter=',i,j,ks,bi,bj,myThid,myIter
83 ENDIF
84 IF (hFactmp.GT.hFacSup) THEN
85 write(0,'(2A,6I4,I10)') 'WARNING: hFacC > hFacSup at:',
86 & ' i,j,k,bi,bj,Thid,Iter=',i,j,ks,bi,bj,myThid,myIter
87 ENDIF
88 write(0,'(A,2F10.6,1PE14.6)') 'hFac_n-1,hFac_n,eta =',
89 & hfacC(i,j,ks,bi,bj), hFactmp, etaFld(i,j,bi,bj)
90 write(0,'(2A)') 'STOP in CALC_SURF_DR :',
91 & ' too SMALL or too LARGE hFac !'
92 STOP 'ABNORMAL END: S/R CALC_SURF_DR'
93 C----------
94 ENDIF
95 ENDIF
96 ENDDO
97 ENDDO
98 ENDDO
99 ENDDO
100
101 _EXCH_XY_R4(hFac_surfC, myThid )
102
103 C-- Apply the mask after the exchange : <-- no need here
104 c DO bj=myByLo(myThid), myByHi(myThid)
105 c DO bi=myBxLo(myThid), myBxHi(myThid)
106 c DO j=1-Oly,sNy+Oly
107 c DO i=1-Olx,sNx+Olx
108 c ks = ksurfC(i,j,bi,bj)
109 c IF (ks.LE.Nr) THEN
110 c hFac_surfC(i,j,bi,bj) = hFac_surfC(i,j,bi,bj)
111 c & * maskC(i,j,ks,bi,bj)
112 c ELSE
113 c hFac_surfC(i,j,bi,bj) = 0.
114 c ENDIF
115 c ENDDO
116 c ENDDO
117 c ENDDO
118 c ENDDO
119
120 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
121
122 C-- Compute fractional thickness of surface level, for U & V point:
123
124 DO bj=myByLo(myThid), myByHi(myThid)
125 DO bi=myBxLo(myThid), myBxHi(myThid)
126 DO j=1,sNy
127 DO i=1,sNx+1
128 ks = ksurfW(i,j,bi,bj)
129 IF (ks.LE.Nr) THEN
130 hhm = hFacC(i-1,j,ks,bi,bj)
131 IF(ks.EQ.ksurfC(i-1,j,bi,bj)) hhm=hFac_surfC(i-1,j,bi,bj)
132 hhp = hFacC(i,j,ks,bi,bj)
133 IF(ks.EQ.ksurfC(i,j,bi,bj)) hhp=hFac_surfC(i,j,bi,bj)
134 hFac_surfW(i,j,bi,bj) = MIN(hhm,hhp)*maskW(i,j,ks,bi,bj)
135 ENDIF
136 ENDDO
137 ENDDO
138 DO j=1,sNy+1
139 DO i=1,sNx
140 ks = ksurfS(i,j,bi,bj)
141 IF (ks.LE.Nr) THEN
142 hhm = hFacC(i,j-1,ks,bi,bj)
143 IF(ks.EQ.ksurfC(i,j-1,bi,bj)) hhm=hFac_surfC(i,j-1,bi,bj)
144 hhp = hFacC(i,j,ks,bi,bj)
145 IF(ks.EQ.ksurfC(i,j,bi,bj)) hhp = hFac_surfC(i,j,bi,bj)
146 hFac_surfS(i,j,bi,bj) = MIN(hhm,hhp)*maskS(i,j,ks,bi,bj)
147 ENDIF
148 ENDDO
149 ENDDO
150 ENDDO
151 ENDDO
152
153 CALL EXCH_UV_XY_RS(hFac_surfW,hFac_surfS,.FALSE.,myThid)
154
155 C-- Apply the mask after the exchange : <-- no need here
156 c DO bj=myByLo(myThid), myByHi(myThid)
157 c DO bi=myBxLo(myThid), myBxHi(myThid)
158 c DO j=1-Oly,sNy+Oly
159 c DO i=1-Olx,sNx+Olx
160
161 c ks = ksurfW(i,j,bi,bj)
162 c IF (ks.LE.Nr) THEN
163 c hFac_surfW(i,j,bi,bj) = hFac_surfW(i,j,bi,bj)
164 c & * maskW(i,j,ks,bi,bj)
165 c ELSE
166 c hFac_surfW(i,j,bi,bj) = 0.
167 c ENDIF
168
169 c ks = ksurfS(i,j,bi,bj)
170 c IF (ks.LE.Nr) THEN
171 c hFac_surfS(i,j,bi,bj) = hFac_surfS(i,j,bi,bj)
172 c & * maskS(i,j,ks,bi,bj)
173 c ELSE
174 c hFac_surfS(i,j,bi,bj) = 0.
175 c ENDIF
176
177 c ENDDO
178 c ENDDO
179 c ENDDO
180 c ENDDO
181
182 WRITE(suff,'(I10.10)') myIter
183 c CALL WRITE_FLD_XY_RS('hFac_surfC.',suff,hFac_surfC,
184 c & myIter,myThid)
185
186 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
187 #endif /* NONLIN_FRSURF */
188
189 RETURN
190 END

  ViewVC Help
Powered by ViewVC 1.1.22