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

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

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


Revision 1.2 - (hide annotations) (download)
Thu May 3 15:39:22 2012 UTC (13 years, 2 months ago) by heimbach
Branch: MAIN
Changes since 1.1: +106 -30 lines
Add changes to repo

1 heimbach 1.1 C $Header: /u/gcmpack/MITgcm/pkg/streamice/streamice_init_varia.F,v 1.6 2011/06/29 16:24:10 dng Exp $
2     C $Name: $
3    
4     #include "STREAMICE_OPTIONS.h"
5    
6     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
7    
8     CBOP
9     SUBROUTINE STREAMICE_ADV_FRONT ( myThid, time_step )
10    
11     C /============================================================\
12     C | SUBROUTINE |
13     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     #include "STREAMICE_ADV.h"
26 heimbach 1.2 #ifdef ALLOW_AUTODIFF_TAMC
27     # include "tamc.h"
28     #endif
29 heimbach 1.1
30     INTEGER myThid
31     _RL time_step
32    
33     #ifdef ALLOW_STREAMICE
34    
35 heimbach 1.2 INTEGER i, j, bi, bj, k, iter_count
36 heimbach 1.1 INTEGER Gi, Gj
37     INTEGER new_partial(4)
38 heimbach 1.2 INTEGER ikey_front, ikey_1
39     _RL iter_flag
40     _RL n_flux_1, n_flux_2
41 heimbach 1.1 _RL href, rho, partial_vol, tot_flux, hpot
42    
43     rho = streamice_density
44 heimbach 1.2 cph iter_count = 0
45     iter_flag = 1. _d 0
46    
47     DO iter_count = 0, 3
48    
49     #ifdef ALLOW_AUTODIFF_TAMC
50     ikey_front = (ikey_dynamics-1)*4 + iter_count + 1
51     CADJ STORE area_shelf_streamice
52     CADJ & = comlev1_stream_front, key = ikey_front
53     CADJ STORE h_streamice
54     CADJ & = comlev1_stream_front, key = ikey_front
55     CADJ STORE hflux_x_si
56     CADJ & = comlev1_stream_front, key = ikey_front
57     CADJ STORE hflux_x_si2
58     CADJ & = comlev1_stream_front, key = ikey_front
59     CADJ STORE hflux_y_si
60     CADJ & = comlev1_stream_front, key = ikey_front
61     CADJ STORE hflux_y_si2
62     CADJ & = comlev1_stream_front, key = ikey_front
63     CADJ STORE streamice_hmask
64     CADJ & = comlev1_stream_front, key = ikey_front
65     CADJ STORE iter_flag
66     CADJ & = comlev1_stream_front, key = ikey_front
67     #endif
68 heimbach 1.1
69 heimbach 1.2 IF ( iter_flag .GT. 0. ) THEN
70    
71     iter_flag = 0. _d 0
72 heimbach 1.1
73     IF (iter_count .gt. 0) then
74     DO bj=myByLo(myThid),myByHi(myThid)
75     DO bi=myBxLo(myThid),myBxHi(myThid)
76     DO j=1-OLy,sNy+OLy
77     DO i=1-OLx,sNx+OLx
78     hflux_x_SI(i,j,bi,bj)=hflux_x_SI2(i,j,bi,bj)
79     hflux_y_SI(i,j,bi,bj)=hflux_y_SI2(i,j,bi,bj)
80     hflux_x_SI2(i,j,bi,bj) = 0. _d 0
81     hflux_y_SI2(i,j,bi,bj) = 0. _d 0
82     ENDDO
83     ENDDO
84     ENDDO
85     ENDDO
86     ENDIF
87    
88 heimbach 1.2 ! iter_count = iter_count + 1
89 heimbach 1.1
90     DO bj=myByLo(myThid),myByHi(myThid)
91     DO bi=myBxLo(myThid),myBxHi(myThid)
92 heimbach 1.2
93 heimbach 1.1 DO j=1-1,sNy+1
94     Gj = (myYGlobalLo-1)+(bj-1)*sNy+j
95 heimbach 1.2 cph IF ((Gj .ge. 1) .and. (Gj .le. Ny)) THEN
96 heimbach 1.1 DO i=1-1,sNx+1
97     Gi = (myXGlobalLo-1)+(bi-1)*sNx+i
98 heimbach 1.2
99     #ifdef ALLOW_AUTODIFF_TAMC
100     act1 = bi - myBxLo(myThid)
101     max1 = myBxHi(myThid) - myBxLo(myThid) + 1
102     act2 = bj - myByLo(myThid)
103     max2 = myByHi(myThid) - myByLo(myThid) + 1
104     act3 = myThid - 1
105     max3 = nTx*nTy
106     act4 = ikey_front - 1
107     ikey_1 = i
108     & + sNx*(j-1)
109     & + sNx*sNy*act1
110     & + sNx*sNy*max1*act2
111     & + sNx*sNy*max1*max2*act3
112     & + sNx*sNy*max1*max2*max3*act4
113     CADJ STORE area_shelf_streamice(i,j,bi,bj)
114     CADJ & = comlev1_stream_ij, key = ikey_1
115     CADJ STORE h_streamice(i,j,bi,bj)
116     CADJ & = comlev1_stream_ij, key = ikey_1
117     CADJ STORE hflux_x_si(i,j,bi,bj)
118     CADJ & = comlev1_stream_ij, key = ikey_1
119     CADJ STORE hflux_y_si(i,j,bi,bj)
120     CADJ & = comlev1_stream_ij, key = ikey_1
121     CADJ STORE streamice_hmask(i,j,bi,bj)
122     CADJ & = comlev1_stream_ij, key = ikey_1
123     #endif
124    
125     cph IF ((Gi .ge. 1) .and. (Gi .le. Nx) .and.
126     IF ((STREAMICE_Hmask(i,j,bi,bj).eq.0.0 .or.
127 heimbach 1.1 & STREAMICE_Hmask(i,j,bi,bj).eq.2.0)) THEN
128 heimbach 1.2 n_flux_1 = 0. _d 0
129 heimbach 1.1 href = 0. _d 0
130     tot_flux = 0. _d 0
131    
132 heimbach 1.2 #ifdef ALLOW_AUTODIFF_TAMC
133     CADJ STORE hflux_x_SI(i,j,bi,bj)
134     CADJ & = comlev1_stream_ij, key = ikey_1
135     #endif
136 heimbach 1.1 IF (hflux_x_SI(i,j,bi,bj).gt. 0. _d 0) THEN
137 heimbach 1.2 n_flux_1 = n_flux_1 + 1. _d 0
138 heimbach 1.1 href = href + H_streamice(i-1,j,bi,bj)
139     tot_flux = tot_flux + hflux_x_SI(i,j,bi,bj) *
140     & dxG(i,j,bi,bj) * time_step
141     hflux_x_SI(i,j,bi,bj) = 0. _d 0
142     ENDIF
143    
144 heimbach 1.2 #ifdef ALLOW_AUTODIFF_TAMC
145     CADJ STORE hflux_x_SI(i,j,bi,bj)
146     CADJ & = comlev1_stream_ij, key = ikey_1
147     #endif
148 heimbach 1.1 IF (hflux_x_SI(i+1,j,bi,bj).lt. 0. _d 0) THEN
149 heimbach 1.2 n_flux_1 = n_flux_1 + 1. _d 0
150 heimbach 1.1 href = href + H_streamice(i+1,j,bi,bj)
151     tot_flux = tot_flux - hflux_x_SI(i+1,j,bi,bj) *
152     & dxG(i+1,j,bi,bj) * time_step
153     hflux_x_SI(i+1,j,bi,bj) = 0. _d 0
154     ENDIF
155    
156 heimbach 1.2 #ifdef ALLOW_AUTODIFF_TAMC
157     CADJ STORE hflux_y_SI(i,j,bi,bj)
158     CADJ & = comlev1_stream_ij, key = ikey_1
159     #endif
160 heimbach 1.1 IF (hflux_y_SI(i,j,bi,bj).gt. 0. _d 0) THEN
161 heimbach 1.2 n_flux_1 = n_flux_1 + 1. _d 0
162 heimbach 1.1 href = href + H_streamice(i,j-1,bi,bj)
163     tot_flux = tot_flux + hflux_y_SI(i,j,bi,bj) *
164     & dyG(i,j,bi,bj) * time_step
165     hflux_y_SI(i,j,bi,bj) = 0. _d 0
166     ENDIF
167    
168 heimbach 1.2 #ifdef ALLOW_AUTODIFF_TAMC
169     CADJ STORE hflux_y_SI(i,j,bi,bj)
170     CADJ & = comlev1_stream_ij, key = ikey_1
171     #endif
172 heimbach 1.1 IF (hflux_y_SI(i,j+1,bi,bj).lt. 0. _d 0) THEN
173 heimbach 1.2 n_flux_1 = n_flux_1 + 1. _d 0
174 heimbach 1.1 href = href + H_streamice(i,j+1,bi,bj)
175     tot_flux = tot_flux - hflux_y_SI(i,j+1,bi,bj) *
176     & dyG(i,j+1,bi,bj) * time_step
177     hflux_y_SI(i,j+1,bi,bj) = 0. _d 0
178     ENDIF
179    
180 heimbach 1.2 IF (n_flux_1 .gt. 0.) THEN
181 heimbach 1.1
182 heimbach 1.2 href = href / n_flux_1
183 heimbach 1.1 partial_vol = H_streamice (i,j,bi,bj) *
184     & area_shelf_streamice (i,j,bi,bj) + tot_flux
185     hpot = partial_vol * recip_rA(i,j,bi,bj)
186    
187     IF (hpot .eq. href) THEN ! cell is exactly covered, no overflow
188     STREAMICE_hmask (i,j,bi,bj) = 1.0
189     H_streamice (i,j,bi,bj) = href
190     area_shelf_streamice(i,j,bi,bj) =
191     & rA(i,j,bi,bj)
192     ELSEIF (hpot .lt. href) THEN ! cell still unfilled
193    
194     ! PRINT *, "PARTIAL CELL INCREASED", tot_flux, i,
195     ! & area_shelf_streamice (i,j,bi,bj),
196     ! & H_streamice (i,j,bi,bj)
197    
198     STREAMICE_hmask (i,j,bi,bj) = 2.0
199     area_shelf_streamice (i,j,bi,bj) = partial_vol / href
200     H_streamice (i,j,bi,bj) = href
201     ELSE ! cell is filled - do overflow
202    
203     ! PRINT *, "CELL FILLED"
204    
205     STREAMICE_hmask (i,j,bi,bj) = 1.0
206     area_shelf_streamice(i,j,bi,bj) =
207     & rA(i,j,bi,bj)
208    
209    
210     partial_vol = partial_vol - href * rA(i,j,bi,bj)
211    
212 heimbach 1.2 iter_flag = 1. _d 0
213 heimbach 1.1
214 heimbach 1.2 n_flux_2 = 0. _d 0 ;
215 heimbach 1.1 DO k=1,4
216     new_partial (:) = 0
217     ENDDO
218    
219     DO k=1,2
220     IF (STREAMICE_ufacemask(i-1+k,j,bi,bj).eq.2.0) THEN ! at a permanent calving boundary - no advance allowed
221 heimbach 1.2 n_flux_2 = n_flux_2 + 1. _d 0
222 heimbach 1.1 ELSEIF (STREAMICE_hmask(i+2*k-3,j,bi,bj).eq.0 _d 0) THEN ! adjacent cell is completely ice free
223 heimbach 1.2 n_flux_2 = n_flux_2 + 1. _d 0
224 heimbach 1.1 new_partial (k) = 1
225     ENDIF
226     ENDDO
227     DO k=1,2
228     IF (STREAMICE_vfacemask (i,j-1+k,bi,bj).eq.2.0) THEN
229 heimbach 1.2 n_flux_2 = n_flux_2 + 1. _d 0
230 heimbach 1.1 ELSEIF (STREAMICE_hmask(i,j+2*k-3,bi,bj).eq.0 _d 0) THEN
231 heimbach 1.2 n_flux_2 = n_flux_2 + 1. _d 0
232 heimbach 1.1 new_partial (k+2) = 1
233     ENDIF
234     ENDDO
235    
236 heimbach 1.2 IF (n_flux_2 .eq. 0.) THEN ! there is nowhere to put the extra ice!
237 heimbach 1.1 H_streamice(i,j,bi,bj) = href + partial_vol *
238     & recip_rA(i,j,bi,bj)
239     ELSE
240     H_streamice(i,j,bi,bj) = href
241    
242     DO k=1,2
243     IF (new_partial(k) .eq. 1) THEN
244     hflux_x_SI2(i-1+k,j,bi,bj) =
245 heimbach 1.2 & partial_vol/time_step/n_flux_2/
246 heimbach 1.1 & dxG(i-1+k,j,bi,bj)
247     ENDIF
248     ENDDO
249    
250     DO k=1,2
251     IF (new_partial(k+2) .eq. 1) THEN
252     hflux_y_SI2(i,j-1+k,bi,bj) =
253 heimbach 1.2 & partial_vol/time_step/n_flux_2/
254 heimbach 1.1 & dxG(i,j-1+k,bi,bj)
255     ENDIF
256     ENDDO
257    
258     ENDIF
259     ENDIF
260     ENDIF
261 heimbach 1.2
262 heimbach 1.1 ENDIF
263     ENDDO
264 heimbach 1.2 cph ENDIF
265 heimbach 1.1 ENDDO
266 heimbach 1.2 c
267 heimbach 1.1 ENDDO
268     ENDDO
269 heimbach 1.2 c
270     ENDIF
271 heimbach 1.1 ENDDO
272    
273 heimbach 1.2 cph IF (iter_count.gt.1) THEN
274     cph PRINT *, "FRONT ADVANCE: ", iter_count, " ITERATIONS"
275     cph ENDIF
276 heimbach 1.1
277    
278    
279     #endif
280     RETURN
281     END

  ViewVC Help
Powered by ViewVC 1.1.22