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

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

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


Revision 1.1 - (hide 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 jmc 1.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