/[MITgcm]/MITgcm_contrib/dgoldberg/streamice/streamice_velmask_upd.F
ViewVC logotype

Annotation of /MITgcm_contrib/dgoldberg/streamice/streamice_velmask_upd.F

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


Revision 1.7 - (hide annotations) (download)
Wed Aug 27 19:29:14 2014 UTC (10 years, 10 months ago) by dgoldberg
Branch: MAIN
CVS Tags: HEAD
Changes since 1.6: +66 -120 lines
updating contrib streamice repo with latest files, and separated out convergence checks; and parameterised maximum iteration counts and interface w shelfice for coupling

1 dgoldberg 1.7 C $Header: /u/gcmpack/MITgcm/pkg/streamice/streamice_velmask_upd.F,v 1.7 2014/04/30 12:18:32 dgoldberg Exp $
2 heimbach 1.1 C $Name: $
3    
4     #include "STREAMICE_OPTIONS.h"
5    
6     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
7    
8     CBOP
9     SUBROUTINE STREAMICE_VELMASK_UPD ( myThid )
10    
11     C /============================================================\
12 dgoldberg 1.7 C | SUBROUTINE |
13 heimbach 1.1 C | o |
14     C |============================================================|
15     C | |
16     C \============================================================/
17     IMPLICIT NONE
18    
19     C === Global variables ===
20     #include "SIZE.h"
21     #include "GRID.h"
22     #include "EEPARAMS.h"
23     #include "PARAMS.h"
24     #include "STREAMICE.h"
25 dgoldberg 1.7 !#ifdef ALLOW_PETSC
26     !# ifdef ALLOW_USE_MPI
27     !# include "EESUPPORT.h"
28     !# endif
29     !#endif
30     ! #include "STREAMICE_ADV.h"
31    
32     #ifdef ALLOW_AUTODIFF_TAMC
33     # include "tamc.h"
34 dgoldberg 1.4 #endif
35 heimbach 1.1
36     INTEGER myThid
37    
38     #ifdef ALLOW_STREAMICE
39    
40     INTEGER i, j, bi, bj, ki, kj
41 dgoldberg 1.7 INTEGER maskFlag, myThidCopy
42 dgoldberg 1.4 CHARACTER*(MAX_LEN_MBUF) msgBuf
43     #ifdef ALLOW_USE_MPI
44     integer mpiRC, mpiMyWid
45     #endif
46     #ifdef ALLOW_PETSC
47     _RS DoFCount
48     integer n_dofs_proc_loc (0:nPx*nPy-1)
49     integer n_dofs_cum_sum (0:nPx*nPy-1)
50     #endif
51 heimbach 1.1
52 dgoldberg 1.2 _EXCH_XY_RL( H_streamice, myThid )
53     _EXCH_XY_RL( area_shelf_streamice, myThid )
54 dgoldberg 1.7 _EXCH_XY_RS( STREAMICE_hmask, myThid )
55 dgoldberg 1.2
56 heimbach 1.1 DO bj=myByLo(myThid),myByHi(myThid)
57     DO bi=myBxLo(myThid),myBxHi(myThid)
58     DO j=1-OLy,sNy+OLy
59     DO i=1-OLx,sNx+OLx
60 dgoldberg 1.7 STREAMICE_umask(i,j,bi,bj) = -1. _d 0
61     STREAMICE_vmask(i,j,bi,bj) = -1. _d 0
62 heimbach 1.1 STREAMICE_ufacemask(i,j,bi,bj) = 0. _d 0
63     STREAMICE_vfacemask(i,j,bi,bj) = 0. _d 0
64     ENDDO
65     ENDDO
66     ENDDO
67     ENDDO
68    
69     DO bj=myByLo(myThid),myByHi(myThid)
70     DO bi=myBxLo(myThid),myBxHi(myThid)
71     DO j=0,sNy+1
72     DO i=0,sNx+1
73 dgoldberg 1.7 IF (STREAMICE_hmask(i,j,bi,bj) .eq. 1.0) THEN
74 heimbach 1.1
75     DO kj=0,1
76     DO ki=0,1
77 dgoldberg 1.7 if (STREAMICE_umask(i+ki,j+kj,bi,bj).eq.-1.0) then
78     STREAMICE_umask (i+ki,j+kj,bi,bj) = 1.0
79     endif
80     if (STREAMICE_vmask(i+ki,j+kj,bi,bj).eq.-1.0) then
81     STREAMICE_vmask (i+ki,j+kj,bi,bj) = 1.0
82     endif
83 heimbach 1.1 ENDDO
84     ENDDO
85    
86     DO ki=0,1
87     maskFlag=INT(STREAMICE_ufacemask_bdry(i+ki,j,bi,bj))
88     IF (maskFlag.EQ.3) THEN
89     DO kj=0,1
90 dgoldberg 1.7 if (STREAMICE_umask(i+ki,j+kj,bi,bj).ne.0.0) then
91     STREAMICE_umask(i+ki,j+kj,bi,bj) = 3.0
92     endif
93     if(STREAMICE_vmask(i+ki,j+kj,bi,bj).ne.0.0) then
94     STREAMICE_vmask(i+ki,j+kj,bi,bj) = 3.0
95     endif
96 heimbach 1.1 ENDDO
97     STREAMICE_ufacemask(i+ki,j,bi,bj) = 3.0
98     ELSE IF (maskFlag.EQ.2) THEN
99 dgoldberg 1.3 !DO kj=0,1
100     STREAMICE_ufacemask(i+ki,j,bi,bj) = 2.0
101     !ENDDO
102 heimbach 1.1 ELSE IF (maskFlag.EQ.4) THEN
103     DO kj=0,1
104     STREAMICE_umask(i+ki,j+kj,bi,bj) = 0.0
105 dgoldberg 1.7 STREAMICE_vmask(i+ki,j+kj,bi,bj) = 0.0
106 heimbach 1.1 ENDDO
107     STREAMICE_ufacemask(i+ki,j,bi,bj) = 4.0
108     ELSE IF (maskFlag.EQ.0) THEN
109     DO kj=0,1
110     STREAMICE_umask(i+ki,j+kj,bi,bj) = 0.0
111 dgoldberg 1.7 STREAMICE_vmask(i+ki,j+kj,bi,bj) = 0.0
112 heimbach 1.1 ENDDO
113     STREAMICE_ufacemask(i+ki,j,bi,bj) = 0.0
114     ELSE IF (maskFlag.EQ.1) THEN
115     DO kj=0,1
116     STREAMICE_umask(i+ki,j+kj,bi,bj) = 0.0
117     ENDDO
118     END IF
119     ENDDO
120    
121     DO kj=0,1
122     maskFlag=INT(STREAMICE_vfacemask_bdry(i,j+kj,bi,bj))
123     IF (maskFlag.EQ.3) THEN
124     DO ki=0,1
125 dgoldberg 1.7 if(STREAMICE_vmask(i+ki,j+kj,bi,bj).ne.0.0) then
126     STREAMICE_vmask(i+ki,j+kj,bi,bj) = 3.0
127     endif
128     if(STREAMICE_umask(i+ki,j+kj,bi,bj).ne.0.0) then
129     STREAMICE_umask(i+ki,j+kj,bi,bj) = 3.0
130     endif
131 heimbach 1.1 ENDDO
132 dgoldberg 1.3 STREAMICE_vfacemask(i,j+kj,bi,bj) = 3.0
133 heimbach 1.1 ELSE IF (maskFlag.EQ.2) THEN
134 dgoldberg 1.3 !DO ki=0,1
135     STREAMICE_vfacemask(i,j+kj,bi,bj) = 2.0
136     !ENDDO
137 heimbach 1.1 ELSE IF (maskFlag.EQ.4) THEN
138     DO ki=0,1
139     STREAMICE_umask(i+ki,j+kj,bi,bj) = 0.0
140 dgoldberg 1.7 STREAMICE_vmask(i+ki,j+kj,bi,bj) = 0.0
141 heimbach 1.1 ENDDO
142 dgoldberg 1.3 STREAMICE_vfacemask(i,j+kj,bi,bj) = 4.0
143 heimbach 1.1 ELSE IF (maskFlag.EQ.0) THEN
144     DO ki=0,1
145     STREAMICE_umask(i+ki,j+kj,bi,bj) = 0.0
146 dgoldberg 1.7 STREAMICE_vmask(i+ki,j+kj,bi,bj) = 0.0
147 heimbach 1.1 ENDDO
148     STREAMICE_vfacemask(i+ki,j,bi,bj) = 0.0
149     ELSE IF (maskFlag.EQ.1) THEN
150     DO ki=0,1
151     STREAMICE_vmask(i+ki,j+kj,bi,bj) = 0.0
152     ENDDO
153     ENDIF
154     ENDDO
155    
156     IF (i .lt. sNx+OLx) THEN
157     IF ((STREAMICE_hmask(i+1,j,bi,bj) .eq. 0.0) .OR.
158     & (STREAMICE_hmask(i+1,j,bi,bj) .eq. 2.0)) THEN
159     !right boundary or adjacent to unfilled cell
160     STREAMICE_ufacemask(i+1,j,bi,bj) = 2.0
161     ENDIF
162     ENDIF
163    
164     IF (i .gt. 1-OLx) THEN
165     IF ((STREAMICE_hmask(i-1,j,bi,bj) .eq. 0.0) .OR.
166     & (STREAMICE_hmask(i-1,j,bi,bj) .eq. 2.0)) THEN
167     !left boundary or adjacent to unfilled cell
168     STREAMICE_ufacemask(i,j,bi,bj) = 2
169     ENDIF
170     ENDIF
171    
172     IF (j .lt. sNy+OLy) THEN
173     IF ((STREAMICE_hmask(i,j+1,bi,bj) .eq. 0.0) .OR.
174     & (STREAMICE_hmask(i,j+1,bi,bj) .eq. 2.0)) THEN
175     !top boundary or adjacent to unfilled cell
176     STREAMICE_vfacemask(i,j+1,bi,bj) = 2
177     ENDIF
178     ENDIF
179    
180     IF (j .gt. 1-OLy) THEN
181     IF ((STREAMICE_hmask(i,j-1,bi,bj) .eq. 0.0) .OR.
182     & (STREAMICE_hmask(i,j-1,bi,bj) .eq. 2.0)) THEN
183     !bot boundary or adjacent to unfilled cell
184     STREAMICE_vfacemask(i,j,bi,bj) = 2.0
185     ENDIF
186     ENDIF
187    
188     ENDIF
189     ENDDO
190     ENDDO
191     ENDDO
192     ENDDO
193    
194 dgoldberg 1.7 !$TAF STORE streamice_umask = comlev1, key=ikey_dynamics
195     !$TAF STORE streamice_vmask = comlev1, key=ikey_dynamics
196 heimbach 1.1
197 dgoldberg 1.7 DO bj=myByLo(myThid),myByHi(myThid)
198     DO bi=myBxLo(myThid),myBxHi(myThid)
199 dgoldberg 1.4 DO j=1,sNy
200     DO i=1,sNx
201 dgoldberg 1.7 IF(streamice_umask(i,j,bi,bj).eq.-1.0) THEN
202     streamice_umask(i,j,bi,bj)=0.0
203 dgoldberg 1.4 ENDIF
204 dgoldberg 1.7 IF(streamice_vmask(i,j,bi,bj).eq.-1.0) THEN
205     streamice_vmask(i,j,bi,bj)=0.0
206 dgoldberg 1.4 ENDIF
207     ENDDO
208     ENDDO
209     ENDDO
210     ENDDO
211    
212 dgoldberg 1.7
213    
214     _EXCH_XY_RS( STREAMICE_ufacemask, myThid )
215     _EXCH_XY_RS( STREAMICE_vfacemask, myThid )
216     _EXCH_XY_RS( STREAMICE_umask, myThid )
217     _EXCH_XY_RS( STREAMICE_vmask, myThid )
218 dgoldberg 1.4
219 dgoldberg 1.7 ! CALL WRITE_FULLARRAY_RL ("umask",STREAMICE_umask,
220     ! c 1,0,0,1,0,myThid)
221     ! CALL WRITE_FLD_XY_RL ("umask","",STREAMICE_umask,0,myThid)
222     ! CALL WRITE_FLD_XY_RL ("vmask","",STREAMICE_vmask,0,myThid)
223     ! CALL WRITE_FLD_XY_RL ("ufacemask","",STREAMICE_ufacemask,0,myThid)
224     ! CALL WRITE_FLD_XY_RL ("vfacemask","",STREAMICE_vfacemask,0,myThid)
225     ! CALL WRITE_FLD_XY_RL ("ufacemaskbdry","",
226     ! & STREAMICE_ufacemask_bdry,0,myThid)
227     ! CALL WRITE_FLD_XY_RL ("vfacemaskbdry","",
228     ! & STREAMICE_vfacemask_bdry,0,myThid)
229 dgoldberg 1.4
230 dgoldberg 1.7 #ifdef ALLOW_PETSC
231     myThidCopy = myThid
232     call streamice_petsc_numerate (myThidCopy)
233     #endif
234 heimbach 1.1
235    
236     #endif
237     RETURN
238     END

  ViewVC Help
Powered by ViewVC 1.1.22