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

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

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


Revision 1.8 - (hide annotations) (download)
Wed Oct 10 15:03:10 2012 UTC (12 years, 9 months ago) by dgoldberg
Branch: MAIN
Changes since 1.7: +4 -4 lines
new parameters

1 dgoldberg 1.8 C $Header: /u/gcmpack/MITgcm_contrib/dgoldberg/streamice/streamice_init_varia.F,v 1.7 2012/10/04 15:40:16 dgoldberg Exp $
2 heimbach 1.1 C $Name: $
3    
4 dgoldberg 1.6 #include "COST_OPTIONS.h"
5 heimbach 1.1 #include "STREAMICE_OPTIONS.h"
6    
7     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
8    
9     CBOP
10     SUBROUTINE STREAMICE_INIT_VARIA( myThid )
11     C /============================================================\
12     C | SUBROUTINE STREAMICE_INIT_VARIA |
13     C | o Routine to initialize STREAMICE variables. |
14     C |============================================================|
15     C | Initialize STREAMICE parameters and variables. |
16     C \============================================================/
17     IMPLICIT NONE
18    
19     C === Global variables ===
20     #include "SIZE.h"
21     #include "GRID.h"
22 dgoldberg 1.4 #include "SET_GRID.h"
23 heimbach 1.1 #include "EEPARAMS.h"
24     #include "PARAMS.h"
25     #include "STREAMICE.h"
26     #include "STREAMICE_CG.h"
27     #include "STREAMICE_ADV.h"
28    
29     C === Routine arguments ===
30     C myThid - Number of this instance of STREAMICE_INIT_VARIA
31     INTEGER myThid
32     CEndOfInterface
33    
34     #ifdef ALLOW_STREAMICE
35     C === Local variables ===
36     C I,J,bi,bj - Loop counters
37     INTEGER i, j, k, bi, bj, Gi, Gj
38     INTEGER col_y, col_x
39 dgoldberg 1.4 _RL slope_pos, c1, x, y, lenx, leny
40 heimbach 1.1 CHARACTER*(MAX_LEN_MBUF) msgBuf
41     CEOP
42    
43     C ZERO OUT FLOATING POINT ARRAYS
44    
45     DO bj = myByLo(myThid), myByHi(myThid)
46     DO bi = myBxLo(myThid), myBxHi(myThid)
47     DO j=1-Oly,sNy+Oly
48     DO i=1-Olx,sNx+Olx
49     H_streamIce(i,j,bi,bj) = 0. _d 0
50     U_streamice(i,j,bi,bj) = 0. _d 0
51     V_streamice(i,j,bi,bj) = 0. _d 0
52     visc_streamice(i,j,bi,bj) = 0. _d 0
53     tau_beta_eff_streamice(i,j,bi,bj) = 0. _d 0
54     float_frac_streamice(i,j,bi,bj) = 0. _d 0
55     base_el_streamice(i,j,bi,bj) = 0. _d 0
56     surf_el_streamice(i,j,bi,bj) = 0. _d 0
57     area_shelf_streamice(i,j,bi,bj) = 0. _d 0
58     mass_ice_streamice(i,j,bi,bj) = 0. _d 0
59 dgoldberg 1.3 BDOT_streamice(i,j,bi,bj) = 0. _d 0
60 dgoldberg 1.4 ! C_basal_friction(i,j,bi,bj) = C_basal_fric_const
61 dgoldberg 1.3 A_glen(i,j,bi,bj) = A_glen_isothermal
62 dgoldberg 1.7 H_streamice_prev(i,j,bi,bj) = 0. _d 0
63 heimbach 1.1 #ifdef ALLOW_AUTODIFF_TAMC
64     ru_old_si(i,j,bi,bj) = 0. _d 0
65     rv_old_si(i,j,bi,bj) = 0. _d 0
66     zu_old_si(i,j,bi,bj) = 0. _d 0
67     zv_old_si(i,j,bi,bj) = 0. _d 0
68     h_after_uflux_SI(i,j,bi,bj) = 0. _d 0
69     #endif
70 dgoldberg 1.6 #ifdef USE_ALT_RLOW
71     R_low_si(i,j,bi,bj) = 0. _d 0
72     #endif
73    
74 dgoldberg 1.4 #ifdef STREAMICE_HYBRID_STRESS
75     do k=1,Nr
76     visc_streamice_full(i,j,k,bi,bj) =
77     & eps_glen_min**((1-n_glen)/n_glen)
78     enddo
79 heimbach 1.5 streamice_taubx (i,j,bi,bj) = 0. _d 0
80     streamice_tauby (i,j,bi,bj) = 0. _d 0
81 dgoldberg 1.4 #endif
82 heimbach 1.1 ENDDO
83     ENDDO
84 dgoldberg 1.6
85     #ifdef ALLOW_COST_TEST
86     cost_func1_streamice (bi,bj) = 0.0
87     #endif
88    
89 heimbach 1.1 ENDDO
90     ENDDO
91    
92     DO j = 1-oly, sNy+oly
93     DO i = 1-olx, sNx+olx
94     DO bj = myByLo(myThid), myByHi(myThid)
95     DO bi = myBxLo(myThid), myBxHi(myThid)
96     cc DO k=1,4
97     DO col_x=-1,1
98     DO col_y=-1,1
99     streamice_cg_A1(i,j,bi,bj,col_x,col_y)=0.0
100     streamice_cg_A2(i,j,bi,bj,col_x,col_y)=0.0
101     streamice_cg_A3(i,j,bi,bj,col_x,col_y)=0.0
102     streamice_cg_A4(i,j,bi,bj,col_x,col_y)=0.0
103     ENDDO
104     ENDDO
105     cc ENDDO
106     ENDDO
107     ENDDO
108     ENDDO
109     ENDDO
110    
111     C INIT. INTEGER ARRAYS
112    
113     DO bj = myByLo(myThid), myByHi(myThid)
114     DO bi = myBxLo(myThid), myBxHi(myThid)
115     DO j=1-Oly,sNy+Oly
116     DO i=1-Olx,sNx+Olx
117     STREAMICE_hmask(i,j,bi,bj) = -1.0
118     STREAMICE_umask(i,j,bi,bj) = 0.0
119     STREAMICE_vmask(i,j,bi,bj) = 0.0
120     STREAMICE_ufacemask(i,j,bi,bj) = 0.0
121     STREAMICE_vfacemask(i,j,bi,bj) = 0.0
122     STREAMICE_float_cond(i,j,bi,bj) = 0.0
123     ENDDO
124     ENDDO
125     ENDDO
126     ENDDO
127    
128 dgoldberg 1.6
129     #ifdef USE_ALT_RLOW
130     ! init alternate array for topog
131     IF ( bathyFile .NE. ' ' ) THEN
132     _BARRIER
133     C The 0 is the "iteration" argument. The ' ' is an empty suffix
134     CALL READ_FLD_XY_RS( bathyFile, '',
135     & R_low_si, 0, myThid )
136    
137     ELSE
138     WRITE(msgBuf,'(A)') 'STREAMICE TOPOG - FILENAME MISSING'
139     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
140     & SQUEEZE_RIGHT , 1)
141     ENDIF
142     #endif
143    
144 dgoldberg 1.4 ! initialize thickness
145 heimbach 1.1
146     IF ( STREAMICEthickInit.EQ.'PARAM' ) THEN
147    
148     WRITE(msgBuf,'(A)') 'initializing analytic thickness'
149     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
150     & SQUEEZE_RIGHT , 1)
151    
152     slope_pos = shelf_edge_pos - shelf_flat_width
153     c1 = 0.0
154     IF (shelf_slope_scale .GT. 0.0) THEN
155     c1 = 1.0 / shelf_slope_scale
156     ENDIF
157    
158     DO bj = myByLo(myThid), myByHi(myThid)
159     DO bi = myBxLo(myThid), myBxHi(myThid)
160     DO j=1,sNy
161     DO i=1,sNx
162     Gi = (myXGlobalLo-1)+(bi-1)*sNx+i
163     Gj = (myYGlobalLo-1)+(bj-1)*sNy+j
164    
165     IF ((Gi.lt.Nx).and.(Gj.lt.Ny)) THEN
166    
167     C IF (flow_dir .EQ. 2.0) THEN
168     IF (.TRUE.) THEN
169     IF (xC(i-1,j,bi,bj).GE.shelf_edge_pos) THEN
170     area_shelf_streamice(i,j,bi,bj) = 0. _d 0
171     STREAMICE_hmask(i,j,bi,bj) = 0. _d 0
172     ELSE
173    
174     IF (xC(i,j,bi,bj).GT.slope_pos) THEN
175     H_streamice (i,j,bi,bj) = shelf_min_draft
176     ELSE
177     H_streamice (i,j,bi,bj) = (shelf_min_draft +
178     & (shelf_max_draft - shelf_min_draft) *
179     & min (1.0, (c1*(slope_pos-xC(i,j,bi,bj)))**2))
180     ENDIF
181    
182     IF (xC(i,j,bi,bj).GT.shelf_edge_pos) THEN
183     area_shelf_streamice(i,j,bi,bj) = rA(i,j,bi,bj) *
184     & (shelf_edge_pos-xG(i,j,bi,bj)) /
185     & (xG(i+1,j,bi,bj)-xG(i,j,bi,bj))
186 heimbach 1.5 IF (area_shelf_streamice(i,j,bi,bj).gt. 0. _d 0) THEN
187 heimbach 1.1 STREAMICE_hmask(i,j,bi,bj) = 2.0
188     ELSE
189     STREAMICE_hmask(i,j,bi,bj) = 0.0
190     H_streamice(i,j,bi,bj) = 0.0
191     ENDIF
192     ELSE
193     area_shelf_streamice(i,j,bi,bj) = rA(i,j,bi,bj)
194     STREAMICE_hmask(i,j,bi,bj) = 1.0
195     ENDIF
196    
197    
198     ENDIF
199     ENDIF
200     ENDIF
201     ENDDO
202     ENDDO
203     ENDDO
204     ENDDO
205    
206     ELSE IF ( STREAMICEthickInit.EQ.'FILE' ) THEN
207    
208     IF ( STREAMICEthickFile .NE. ' ' ) THEN
209     _BARRIER
210     C The 0 is the "iteration" argument. The ' ' is an empty suffix
211     CALL READ_FLD_XY_RS( STREAMICEthickFile, ' ', H_streamice,
212     & 0, myThid )
213     DO bj = myByLo(myThid), myByHi(myThid)
214     DO bi = myBxLo(myThid), myBxHi(myThid)
215     DO j=1,sNy
216     DO i=1,sNx
217     Gi = (myXGlobalLo-1)+(bi-1)*sNx+i
218     Gj = (myYGlobalLo-1)+(bj-1)*sNy+j
219 dgoldberg 1.4 IF ((Gi.lt.Nx.OR.STREAMICE_EW_periodic).and.
220     & (Gj.lt.Ny.OR.STREAMICE_NS_periodic)) THEN
221 heimbach 1.1 IF (H_streamice(i,j,bi,bj).GT.0. _d 0) THEN
222     area_shelf_streamice(i,j,bi,bj) = rA(i,j,bi,bj)
223     STREAMICE_hmask(i,j,bi,bj) = 1.0
224     ELSE
225     area_shelf_streamice(i,j,bi,bj) = 0. _d 0
226     STREAMICE_hmask(i,j,bi,bj) = 0. _d 0
227     ENDIF
228     ENDIF
229     ENDDO
230     ENDDO
231     ENDDO
232     ENDDO
233     ELSE
234     WRITE(msgBuf,'(A)') 'INIT THICKNESS - FILENAME MISSING'
235     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
236     & SQUEEZE_RIGHT , 1)
237     ENDIF
238    
239     ELSE
240    
241     WRITE(msgBuf,'(A)') 'INIT THICKNESS - NOT IMPLENTED'
242     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
243     & SQUEEZE_RIGHT , 1)
244     ENDIF
245    
246 dgoldberg 1.4 ! finish initialize thickness
247    
248     ! initialize basal traction
249    
250     IF ( STREAMICEbasalTracConfig.EQ.'FILE' ) THEN
251    
252     IF ( STREAMICEbasalTracFile .NE. ' ' ) THEN
253     _BARRIER
254     C The 0 is the "iteration" argument. The ' ' is an empty suffix
255     CALL READ_FLD_XY_RS( STREAMICEbasalTracFile, ' ',
256     & C_basal_friction, 0, myThid )
257    
258     ELSE
259     WRITE(msgBuf,'(A)') 'INIT THICKNESS - FILENAME MISSING'
260     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
261     & SQUEEZE_RIGHT , 1)
262     ENDIF
263    
264     ELSE IF (STREAMICEbasalTracConfig.EQ.'UNIFORM' ) THEN
265    
266     DO bj = myByLo(myThid), myByHi(myThid)
267     DO bi = myBxLo(myThid), myBxHi(myThid)
268     DO j=1,sNy
269     DO i=1,sNx
270     C_basal_friction(i,j,bi,bj) = C_basal_fric_const
271     ENDDO
272     ENDDO
273     ENDDO
274     ENDDO
275    
276     ELSE IF (STREAMICEbasalTracConfig.EQ.'2DPERIODIC' ) THEN
277    
278     lenx = sNx*nSx*nPx*delX(1)
279     leny = sNy*nSy*nPy*delY(1)
280     print *, 'lenx', lenx
281     print *, 'leny', leny
282     DO bj = myByLo(myThid), myByHi(myThid)
283     DO bi = myBxLo(myThid), myBxHi(myThid)
284     DO j=1,sNy
285     DO i=1,sNx
286     x = xC(i,j,bi,bj)
287     y = yC(i,j,bi,bj)
288     C_basal_friction(i,j,bi,bj) =
289     & sqrt(C_basal_fric_const**2*
290     & (1+sin(2*streamice_kx_b_init*PI*x/lenx)*
291     & sin(2*streamice_ky_b_init*PI*y/leny)))
292     ENDDO
293     ENDDO
294     ENDDO
295     ENDDO
296    
297     ELSE IF (STREAMICEbasalTracConfig.EQ.'1DPERIODIC' ) THEN
298    
299     lenx = sNx*nSx*nPx*delX(1)
300     print *, 'lenx', lenx
301     DO bj = myByLo(myThid), myByHi(myThid)
302     DO bi = myBxLo(myThid), myBxHi(myThid)
303     DO j=1,sNy
304     DO i=1,sNx
305     x = xC(i,j,bi,bj)
306     y = yC(i,j,bi,bj)
307     C_basal_friction(i,j,bi,bj) =
308     & sqrt(C_basal_fric_const**2*(1+
309     & sin(2*streamice_kx_b_init*PI*x/lenx)))
310     ENDDO
311     ENDDO
312     ENDDO
313     ENDDO
314    
315     ELSE
316    
317 dgoldberg 1.8 WRITE(msgBuf,'(A)') 'INIT TRAC - NOT IMPLENTED'
318     CALL PRINT_ERROR( msgBuf, myThid)
319     STOP 'ABNORMAL END: S/R STREAMICE_INIT_VAR'
320 dgoldberg 1.4 ENDIF
321    
322     ! finish initialize basal traction
323    
324     CALL STREAMICE_UPD_FFRAC_UNCOUPLED ( myThid )
325 heimbach 1.1
326     _EXCH_XY_RL(H_streamice, myThid )
327     _EXCH_XY_RL(STREAMICE_hmask, myThid )
328     _EXCH_XY_RL(area_shelf_streamice, myThid )
329 dgoldberg 1.4 _EXCH_XY_RL(C_basal_friction, myThid )
330 dgoldberg 1.6 #ifdef USE_ALT_RLOW
331     _EXCH_XY_RL(R_low_si, myThid )
332     #endif
333 dgoldberg 1.4
334     #ifdef STREAMICE_HYBRID_STRESS
335    
336     CALL STREAMICE_VISC_BETA (myThid)
337    
338     #endif
339 dgoldberg 1.6
340     CALL WRITE_FLD_XY_RL ( "C_basal_fric", "",
341     & C_basal_friction, 0, myThid )
342 heimbach 1.1 CALL WRITE_FLD_XY_RL ( "H_streamice", "init",
343     & H_streamIce, 0, myThid )
344     CALL WRITE_FLD_XY_RL ( "area_shelf_streamice", "init",
345     & area_shelf_streamice, 0, myThid )
346     CALL WRITE_FLD_XY_RL ( "STREAMICE_hmask", "init",
347     & STREAMICE_hmask, 0, myThid )
348    
349 dgoldberg 1.2 ! CALL STREAMICE_VELMASK_UPD (myThid)
350     ! CALL STREAMICE_UPD_FFRAC_UNCOUPLED ( myThid )
351     ! CALL STREAMICE_VEL_SOLVE( myThid )
352 heimbach 1.1
353     CALL WRITE_FLD_XY_RL ( "U_init", "",
354 dgoldberg 1.4 & C_basal_friction, 0, myThid )
355 heimbach 1.1 CALL WRITE_FLD_XY_RL ( "V_init", "",
356     & V_streamice, 0, myThid )
357 dgoldberg 1.6 #ifdef USE_ALT_RLOW
358     CALL WRITE_FLD_XY_RL ( "R_low_si", "init",
359     & R_low_si, 0, myThid )
360     #endif
361 heimbach 1.1
362     ! CALL WRITE_FULLARRAY_RL ("H",H_streamice,1,0,0,1,0,myThid)
363     ! CALL WRITE_FULLARRAY_RL ("hmask",STREAMICE_hmask,1,0,0,1,0,myThid)
364     ! CALL WRITE_FULLARRAY_RL ("umask",STREAMICE_umask,1,0,0,1,0,myThid)
365    
366     #endif /* ALLOW_STREAMICE */
367    
368     RETURN
369     END
370    

  ViewVC Help
Powered by ViewVC 1.1.22