/[MITgcm]/MITgcm/pkg/seaice/seaice_init_varia.F
ViewVC logotype

Annotation of /MITgcm/pkg/seaice/seaice_init_varia.F

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


Revision 1.2 - (hide annotations) (download)
Wed Dec 20 12:25:15 2006 UTC (17 years, 5 months ago) by mlosch
Branch: MAIN
CVS Tags: checkpoint58u_post, checkpoint58w_post, checkpoint58x_post, checkpoint58t_post, checkpoint58v_post
Changes since 1.1: +3 -3 lines
o fix multi-category seaice:
 - change cpp flag SEAICE_MULTILEVEL to more meaningful name:
   SEAICE_MULTICATEGORY
 - fix short wave heat flux
o replace field areaLoc by scalar variable

1 mlosch 1.2 C $Header: /u/gcmpack/MITgcm/pkg/seaice/seaice_init_varia.F,v 1.1 2006/10/22 01:11:02 heimbach Exp $
2 heimbach 1.1 C $Name: $
3    
4     #include "SEAICE_OPTIONS.h"
5    
6     CStartOfInterface
7     SUBROUTINE SEAICE_INIT_VARIA( myThid )
8     C /==========================================================\
9     C | SUBROUTINE SEAICE_INIT_VARIA |
10     C | o Initialization of sea ice model. |
11     C |==========================================================|
12     C \==========================================================/
13     IMPLICIT NONE
14    
15     C === Global variables ===
16     #include "SIZE.h"
17     #include "EEPARAMS.h"
18     #include "PARAMS.h"
19     #include "GRID.h"
20     #include "SEAICE.h"
21     CML#include "SEAICE_GRID.h"
22     #include "SEAICE_DIAGS.h"
23     #include "SEAICE_PARAMS.h"
24     #include "FFIELDS.h"
25     #ifdef ALLOW_EXCH2
26     #include "W2_EXCH2_TOPOLOGY.h"
27     #include "W2_EXCH2_PARAMS.h"
28     #endif /* ALLOW_EXCH2 */
29    
30     C === Routine arguments ===
31     C myThid - Thread no. that called this routine.
32     INTEGER myThid
33     CEndOfInterface
34    
35     C === Local variables ===
36     C i,j,k,bi,bj - Loop counters
37    
38     INTEGER i, j, k, bi, bj
39     _RS mask_uice
40     INTEGER myIter, myTile
41    
42     cph(
43     cph make sure TAF sees proper initialisation
44     cph to avoid partial recomputation issues
45     DO bj=myByLo(myThid),myByHi(myThid)
46     DO bi=myBxLo(myThid),myBxHi(myThid)
47     c
48     DO K=1,3
49     DO J=1-OLy,sNy+OLy
50     DO I=1-OLx,sNx+OLx
51     HEFF(I,J,k,bi,bj)=ZERO
52     AREA(I,J,k,bi,bj)=ZERO
53     UICE(I,J,k,bi,bj)=ZERO
54     VICE(I,J,k,bi,bj)=ZERO
55     ENDDO
56     ENDDO
57     ENDDO
58     c
59     DO J=1-OLy,sNy+OLy
60     DO I=1-OLx,sNx+OLx
61     HSNOW(I,J,bi,bj)=ZERO
62     ZETA(I,J,bi,bj)=ZERO
63     ENDDO
64     ENDDO
65     c
66     ENDDO
67     ENDDO
68     cph)
69    
70     C-- Initialize grid info
71     DO bj=myByLo(myThid),myByHi(myThid)
72     DO bi=myBxLo(myThid),myBxHi(myThid)
73     DO j=1-OLy,sNy+OLy
74     DO i=1-OLx,sNx+OLx
75     HEFFM(i,j,bi,bj)=ONE
76     IF (_hFacC(i,j,1,bi,bj).eq.0.) HEFFM(i,j,bi,bj)=ZERO
77     ENDDO
78     ENDDO
79     DO J=1-OLy+1,sNy+OLy
80     DO I=1-OLx+1,sNx+OLx
81     #ifndef SEAICE_CGRID
82     UVM(i,j,bi,bj)=ZERO
83     mask_uice=HEFFM(I,J, bi,bj)+HEFFM(I-1,J-1,bi,bj)
84     & +HEFFM(I,J-1,bi,bj)+HEFFM(I-1,J, bi,bj)
85     IF(mask_uice.GT.3.5) UVM(I,J,bi,bj)=ONE
86     #else
87     seaiceMaskU(I,J,bi,bj)= 0.0 _d 0
88     seaiceMaskV(I,J,bi,bj)= 0.0 _d 0
89     mask_uice=HEFFM(I,J,bi,bj)+HEFFM(I-1,J ,bi,bj)
90     IF(mask_uice.GT.1.5) seaiceMaskU(I,J,bi,bj)=ONE
91     mask_uice=HEFFM(I,J,bi,bj)+HEFFM(I ,J-1,bi,bj)
92     IF(mask_uice.GT.1.5) seaiceMaskV(I,J,bi,bj)=ONE
93     #endif /* not SEAICE_CGRID */
94     ENDDO
95     ENDDO
96    
97     #ifdef ALLOW_EXCH2
98     #ifdef SEAICE_CGRID
99     #else
100     C-- Special stuff for cubed sphere: assume grid is rectangular and
101     C set UV mask to zero except for Arctic and Antarctic cube faces.
102     IF (useCubedSphereExchange) THEN
103     myTile = W2_myTileList(bi)
104     IF ( exch2_myFace(myTile) .EQ. 1 .OR.
105     & exch2_myFace(myTile) .EQ. 2 .OR.
106     & exch2_myFace(myTile) .EQ. 4 .OR.
107     & exch2_myFace(myTile) .EQ. 5 ) THEN
108     DO J=1-OLy,sNy+OLy
109     DO I=1-OLx,sNx+OLx
110     UVM(i,j,bi,bj)=ZERO
111     ENDDO
112     ENDDO
113     ELSEIF ( exch2_isWedge(myTile) .EQ. 1 ) THEN
114     I=1
115     DO J=1-OLy,sNy+OLy
116     UVM(i,j,bi,bj)=ZERO
117     ENDDO
118     ELSEIF ( exch2_isSedge(myTile) .EQ. 1 ) THEN
119     J=1
120     DO I=1-OLx,sNx+OLx
121     UVM(i,j,bi,bj)=ZERO
122     ENDDO
123     ENDIF
124     ENDIF
125     #endif
126     #endif /* ALLOW_EXCH2 */
127    
128     DO j=1-OLy,sNy+OLy
129     DO i=1-OLx,sNx+OLx
130     TICE(I,J,bi,bj)=273.0 _d 0
131 mlosch 1.2 #ifdef SEAICE_MULTICATEGORY
132 heimbach 1.1 DO k=1,MULTDIM
133     TICES(I,J,k,bi,bj)=273.0 _d 0
134     ENDDO
135 mlosch 1.2 #endif /* SEAICE_MULTICATEGORY */
136 heimbach 1.1 UICEC (I,J,bi,bj)=ZERO
137     VICEC (I,J,bi,bj)=ZERO
138     DAIRN (I,J,bi,bj)=ZERO
139     DWATN (I,J,bi,bj)=ZERO
140     #ifndef SEAICE_CGRID
141     AMASS (I,J,bi,bj)=1000.0 _d 0
142     #else
143     seaiceMassC(I,J,bi,bj)=1000.0 _d 0
144     seaiceMassU(I,J,bi,bj)=1000.0 _d 0
145     seaiceMassV(I,J,bi,bj)=1000.0 _d 0
146     #endif
147     GWATX (I,J,bi,bj)=ZERO
148     GWATY (I,J,bi,bj)=ZERO
149     ENDDO
150     ENDDO
151    
152     C-- Choose a proxy level for geostrophic velocity,
153     DO j=1-OLy,sNy+OLy
154     DO i=1-OLx,sNx+OLx
155     #ifdef SEAICE_TEST_ICE_STRESS_1
156     KGEO(I,J,bi,bj) = 1
157     #else /* SEAICE_TEST_ICE_STRESS_1 */
158     IF (klowc(i,j,bi,bj) .LT. 2) THEN
159     KGEO(I,J,bi,bj) = 1
160     ELSE
161     KGEO(I,J,bi,bj) = 2
162     DO WHILE ( abs(rC(KGEO(I,J,bi,bj))) .LT. 50.0 .AND.
163     & KGEO(I,J,bi,bj) .LT. (klowc(i,j,bi,bj)-1) )
164     KGEO(I,J,bi,bj) = KGEO(I,J,bi,bj) + 1
165     ENDDO
166     ENDIF
167     #endif /* SEAICE_TEST_ICE_STRESS_1 */
168     ENDDO
169     ENDDO
170    
171     ENDDO
172     ENDDO
173    
174     C-- Update overlap regions
175     #ifdef SEAICE_CGRID
176     CALL EXCH_UV_XY_RL(seaiceMaskU,seaiceMaskV,.FALSE.,myThid)
177     #else
178     _EXCH_XY_R8(UVM, myThid)
179     #endif
180    
181     C-- Now lets look at all these beasts
182     IF ( debugLevel .GE. debLevB ) THEN
183     myIter=0
184     CALL PLOT_FIELD_XYRL( HEFFM , 'Current HEFFM ' ,
185     & myIter, myThid )
186     #ifdef SEAICE_CGRID
187     CALL PLOT_FIELD_XYRL( seaiceMaskU, 'Current seaiceMaskU',
188     & myIter, myThid )
189     CALL PLOT_FIELD_XYRL( seaiceMaskV, 'Current seaiceMaskV',
190     & myIter, myThid )
191     #else
192     CALL PLOT_FIELD_XYRL( UVM , 'Current UVM ' ,
193     & myIter, myThid )
194     #endif
195     ENDIF
196    
197     #if (defined (SEAICE_CGRID) && defined (SEAICE_ALLOW_EVP))
198     DO bj=myByLo(myThid),myByHi(myThid)
199     DO bi=myBxLo(myThid),myBxHi(myThid)
200     DO j=1-OLy,sNy+OLy
201     DO i=1-OLx,sNx+OLx
202     seaice_sigma1 (I,J,bi,bj) = 0. _d 0
203     seaice_sigma2 (I,J,bi,bj) = 0. _d 0
204     seaice_sigma12(I,J,bi,bj) = 0. _d 0
205     seaice_div (I,J,bi,bj) = 0. _d 0
206     seaice_tension(I,J,bi,bj) = 0. _d 0
207     seaice_shear (I,J,bi,bj) = 0. _d 0
208     ENDDO
209     ENDDO
210     ENDDO
211     ENDDO
212     #endif /* SEAICE_ALLOW_EVP and SEAICE_CGRID */
213    
214     C-- Set model variables to initial/restart conditions
215     IF ( nIter0 .NE. 0 ) THEN
216    
217     CALL SEAICE_READ_PICKUP ( myThid )
218    
219     ELSE
220    
221     DO bj=myByLo(myThid),myByHi(myThid)
222     DO bi=myBxLo(myThid),myBxHi(myThid)
223     DO j=1-OLy,sNy+OLy
224     DO i=1-OLx,sNx+OLx
225     HSNOW(I,J,bi,bj)=0.2*HEFFM(i,j,bi,bj)
226     YNEG(I,J,bi,bj)=ZERO
227     TMIX(I,J,bi,bj)=TICE(I,J,bi,bj)
228     DO k=1,3
229     HEFF(I,J,k,bi,bj)=SEAICE_initialHEFF*HEFFM(i,j,bi,bj)
230     UICE(I,J,k,bi,bj)=ZERO
231     VICE(I,J,k,bi,bj)=ZERO
232     ENDDO
233     ENDDO
234     ENDDO
235     ENDDO
236     ENDDO
237    
238     C-- Read initial sea-ice thickness from file if available.
239     IF ( HeffFile .NE. ' ' ) THEN
240     _BEGIN_MASTER( myThid )
241     CALL READ_FLD_XY_RL( HeffFile, ' ', ZETA, 0, myThid )
242     _END_MASTER(myThid)
243     _EXCH_XY_R8(ZETA,myThid)
244     DO bj=myByLo(myThid),myByHi(myThid)
245     DO bi=myBxLo(myThid),myBxHi(myThid)
246     DO j=1-OLy,sNy+OLy
247     DO i=1-OLx,sNx+OLx
248     DO k=1,3
249     HEFF(I,J,k,bi,bj) = MAX(ZETA(i,j,bi,bj),ZERO)
250     ENDDO
251     ENDDO
252     ENDDO
253     ENDDO
254     ENDDO
255     ENDIF
256    
257     DO bj=myByLo(myThid),myByHi(myThid)
258     DO bi=myBxLo(myThid),myBxHi(myThid)
259     DO j=1-OLy,sNy+OLy
260     DO i=1-OLx,sNx+OLx
261     DO k=1,3
262     IF(HEFF(I,J,k,bi,bj).GT.ZERO) AREA(I,J,k,bi,bj)=ONE
263     ENDDO
264     ENDDO
265     ENDDO
266     ENDDO
267     ENDDO
268    
269     ENDIF
270    
271     C--- Complete initialization
272     DO bj=myByLo(myThid),myByHi(myThid)
273     DO bi=myBxLo(myThid),myBxHi(myThid)
274     DO j=1-OLy,sNy+OLy
275     DO i=1-OLx,sNx+OLx
276     ZETA(I,J,bi,bj)=HEFF(I,J,1,bi,bj)*(1.0 _d 11)
277     ETA(I,J,bi,bj)=ZETA(I,J,bi,bj)/4.0 _d 0
278     ENDDO
279     ENDDO
280     #ifdef ATMOSPHERIC_LOADING
281     IF ( useRealFreshWaterFlux .AND. .NOT.useThSIce ) THEN
282     DO j=1-OLy,sNy+OLy
283     DO i=1-OLx,sNx+OLx
284     sIceLoad(i,j,bi,bj) = HEFF(I,J,1,bi,bj)*SEAICE_rhoIce
285     & + HSNOW(I,J,bi,bj)* 330. _d 0
286    
287     ENDDO
288     ENDDO
289     ENDIF
290     #endif
291     ENDDO
292     ENDDO
293    
294     RETURN
295     END

  ViewVC Help
Powered by ViewVC 1.1.22