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

Diff of /MITgcm_contrib/dgoldberg/streamice/streamice_advect_thickness.F

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

revision 1.3 by dgoldberg, Thu Jul 26 16:13:18 2012 UTC revision 1.9 by dgoldberg, Tue Jun 11 17:42:17 2013 UTC
# Line 7  C---+----1----+----2----+----3----+----4 Line 7  C---+----1----+----2----+----3----+----4
7    
8    
9  CBOP  CBOP
10        SUBROUTINE STREAMICE_ADVECT_THICKNESS ( myThid, time_step )        SUBROUTINE STREAMICE_ADVECT_THICKNESS ( myThid,myIter,time_step )
11    
12  C     /============================================================\  C     /============================================================\
13  C     | SUBROUTINE                                                 |    C     | SUBROUTINE                                                 |  
# Line 28  C     === Global variables === Line 28  C     === Global variables ===
28  # include "tamc.h"  # include "tamc.h"
29  #endif  #endif
30    
31        INTEGER myThid        INTEGER myThid, myIter
32        _RL time_step        _RL time_step
33    
34  #ifdef ALLOW_STREAMICE  #ifdef ALLOW_STREAMICE
35    
36        INTEGER i, j, bi, bj        INTEGER i, j, bi, bj
37        _RL thick_bd        _RL thick_bd, uflux, vflux, max_icfl, loc_icfl
38        _RL SLOPE_LIMITER        _RL time_step_full, time_step_rem
39        _RL sec_per_year, time_step_loc        _RL sec_per_year, time_step_loc, MR, SMB, TMB
40          _RL BCVALX(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
41          _RL BCVALY(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
42          _RS BCMASKX(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
43          _RS BCMASKY(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
44          _RL utrans(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
45          _RL vtrans(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
46          _RL h_after_uflux_SI(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
47          _RL h_after_vflux_SI(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
48          _RL hflux_x_SI(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
49          _RL hflux_y_SI(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
50          
51        CHARACTER*(MAX_LEN_MBUF) msgBuf        CHARACTER*(MAX_LEN_MBUF) msgBuf
       external SLOPE_LIMITER  
52    
53        sec_per_year = 365.*86400.        sec_per_year = 365.*86400.
54    
55        time_step_loc = time_step / sec_per_year        time_step_loc = time_step / sec_per_year
56          time_step_full = time_step_loc
57          time_step_rem = time_step_loc
58          PRINT *, "time_step_loc ", time_step_loc
59    
60  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
61  CADJ STORE streamice_hmask  = comlev1, key=ikey_dynamics  CADJ STORE streamice_hmask  = comlev1, key=ikey_dynamics
# Line 50  CADJ STORE streamice_hmask  = comlev1, k Line 63  CADJ STORE streamice_hmask  = comlev1, k
63    
64        DO bj=myByLo(myThid),myByHi(myThid)        DO bj=myByLo(myThid),myByHi(myThid)
65         DO bi=myBxLo(myThid),myBxHi(myThid)         DO bi=myBxLo(myThid),myBxHi(myThid)
66          DO j=1-OLy,sNy+OLy          DO j=1-3,sNy+3
67           DO i=1-OLx,sNx+OLx           DO i=1-3,sNx+3
68    
69              H_streamice_prev(i,j,bi,bj) =
70         &     H_streamice(i,j,bi,bj)
71    
72            hflux_x_SI (i,j,bi,bj) = 0. _d 0            hflux_x_SI (i,j,bi,bj) = 0. _d 0
73            hflux_y_SI (i,j,bi,bj) = 0. _d 0            hflux_y_SI (i,j,bi,bj) = 0. _d 0
74            hflux_x_SI2 (i,j,bi,bj) = 0. _d 0  
75            hflux_y_SI2 (i,j,bi,bj) = 0. _d 0  
76            IF (STREAMICE_hmask(i,j,bi,bj).eq.1.0) THEN            IF (STREAMICE_ufacemask(i,j,bi,bj).eq.3.0) THEN
77             h_after_uflux_SI (i,j,bi,bj) =             BCMASKX(i,j,bi,bj) = 3.0
78       &      H_streamice (i,j,bi,bj)             BCVALX(i,j,bi,bj) = h_ubdry_values_SI(i,j,bi,bj)
79               utrans(i,j,bi,bj) = .5 * (
80         &      u_streamice(i,j,bi,bj)+u_streamice(i,j+1,bi,bj))
81              ELSEIF (STREAMICE_ufacemask(i,j,bi,bj).eq.4.0) THEN
82               uflux = u_flux_bdry_SI(i,j,bi,bj)
83               BCMASKX(i,j,bi,bj) = 0.0
84               BCVALX(i,j,bi,bj) = uflux
85               utrans(i,j,bi,bj) = 1.0
86              ELSEIF (.not.(
87         &     STREAMICE_hmask(i,j,bi,bj).eq.1.0.OR.
88         &     STREAMICE_hmask(i-1,j,bi,bj).eq.1.0)) THEN
89               BCMASKX(i,j,bi,bj) = 0.0
90               BCVALX(i,j,bi,bj) = 0. _d 0
91               utrans(i,j,bi,bj) = 0. _d 0
92              ELSE
93               BCMASKX(i,j,bi,bj) = 0.0
94               BCVALX(i,j,bi,bj) = 0. _d 0
95               utrans(i,j,bi,bj) = .5 * (
96         &      u_streamice(i,j,bi,bj)+u_streamice(i,j+1,bi,bj))
97            ENDIF            ENDIF
98    
99            thick_bd = h_bdry_values_SI (i,j,bi,bj)            IF (STREAMICE_vfacemask(i,j,bi,bj).eq.3.0) THEN
100            IF (thick_bd .ne. 0. _d 0) THEN             BCMASKy(i,j,bi,bj) = 3.0
101              h_after_uflux_SI (i,j,bi,bj) = thick_bd             BCVALy(i,j,bi,bj) = h_vbdry_values_SI(i,j,bi,bj)
102               vtrans(i,j,bi,bj) = .5 * (
103         &      v_streamice(i,j,bi,bj)+v_streamice(i+1,j,bi,bj))
104              ELSEIF (STREAMICE_vfacemask(i,j,bi,bj).eq.4.0) THEN
105               vflux = v_flux_bdry_SI(i,j,bi,bj)
106               BCMASKy(i,j,bi,bj) = 0.0
107               BCVALy(i,j,bi,bj) = vflux
108               vtrans(i,j,bi,bj) = 1.0
109              ELSEIF (.not.(
110         &     STREAMICE_hmask(i,j,bi,bj).eq.1.0.OR.
111         &     STREAMICE_hmask(i,j-1,bi,bj).eq.1.0)) THEN
112               BCMASKY(i,j,bi,bj) = 0.0
113               BCVALY(i,j,bi,bj) = 0. _d 0
114               vtrans(i,j,bi,bj) = 0. _d 0
115              ELSE
116               BCMASKy(i,j,bi,bj) = 0.0
117               BCVALy(i,j,bi,bj) = 0. _d 0
118               vtrans(i,j,bi,bj) =  .5 * (
119         &      v_streamice(i,j,bi,bj)+v_streamice(i+1,j,bi,bj))
120              ENDIF
121    
122    
123             ENDDO
124            ENDDO
125           ENDDO
126          ENDDO
127    
128          _EXCH_XY_RL(utrans,myThid)
129          _EXCH_XY_RL(vtrans,myThid)
130          _EXCH_XY_RS(BCMASKx,myThid)
131          _EXCH_XY_RS(BCMASKy,myThid)
132          _EXCH_XY_RL(BCVALX,myThid)
133          _EXCH_XY_RL(BCVALY,myThid)
134    
135    #ifndef ALLOW_AUTODIFF_TAMC
136    
137          max_icfl = 1.e-20
138    
139          DO bj=myByLo(myThid),myByHi(myThid)
140           DO bi=myBxLo(myThid),myBxHi(myThid)
141            DO j=1,sNy
142             DO i=1,sNx
143              IF (streamice_hmask(i,j,bi,bj).eq.1.0) THEN
144               loc_icfl=max(abs(utrans(i,j,bi,bj)),
145         &                  abs(utrans(i+1,j,bi,bj))) / dxF(i,j,bi,bj)
146               loc_icfl=max(loc_icfl,max(abs(vtrans(i,j,bi,bj)),
147         &                  abs(vtrans(i,j+1,bi,bj))) / dyF(i,j,bi,bj))
148               if (loc_icfl.gt.max_icfl) then
149                max_icfl = loc_icfl
150               ENDIF
151            ENDIF            ENDIF
152           ENDDO           ENDDO
153          ENDDO          ENDDO
154         ENDDO         ENDDO
155        ENDDO        ENDDO
156    
157          CALL GLOBAL_MAX_R8 (max_icfl, myThid)
158    
159    #endif
160          
161    
162  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
 CADJ STORE h_after_uflux_si = comlev1, key=ikey_dynamics  
163  CADJ STORE streamice_hmask  = comlev1, key=ikey_dynamics  CADJ STORE streamice_hmask  = comlev1, key=ikey_dynamics
164    CADJ STORE H_streamice  = comlev1, key=ikey_dynamics
165    #endif
166    
167    #ifndef ALLOW_AUTODIFF_TAMC
168          do while (time_step_rem .gt. 1.e-15)
169           time_step_loc = min (
170         &  streamice_cfl_factor / max_icfl,
171         &  time_step_rem )
172           if (time_step_loc .lt. time_step_full) then
173            PRINT *, "TAKING PARTIAL TIME STEP", time_step_loc
174           endif
175  #endif  #endif
176    
177        CALL STREAMICE_ADVECT_THICKNESS_X ( myThid,        CALL STREAMICE_ADV_FLUX_FL_X ( myThid ,
178         I   utrans ,
179         I   H_streamice ,
180         I   BCMASKX,
181         I   BCVALX,
182       O   hflux_x_SI,       O   hflux_x_SI,
      O   h_after_uflux_SI,  
183       I   time_step_loc )       I   time_step_loc )
184    
185    
186    
187        DO bj=myByLo(myThid),myByHi(myThid)        DO bj=myByLo(myThid),myByHi(myThid)
188         DO bi=myBxLo(myThid),myBxHi(myThid)         DO bi=myBxLo(myThid),myBxHi(myThid)
189          DO j=1-OLy,sNy+OLy          DO j=1-3,sNy+3
190           DO i=1-OLx,sNx+OLx           DO i=1-1,sNx+1
191            h_after_vflux_SI (i,j,bi,bj) =            IF (STREAMICE_hmask(i,j,bi,bj).eq.1.0) THEN
192       &     h_after_uflux_SI (i,j,bi,bj)             h_after_uflux_SI (i,j,bi,bj) = H_streamice(i,j,bi,bj) -
193         &      (hflux_x_SI(i+1,j,bi,bj)*dyG(i+1,j,bi,bj) -
194         &        hflux_x_SI(i,j,bi,bj)*dyG(i,j,bi,bj))
195         &        * recip_rA (i,j,bi,bj) * time_step_loc
196               IF ( h_after_uflux_SI (i,j,bi,bj).le.0.0) THEN
197                PRINT *, "h neg after x", i,j,hflux_x_SI(i+1,j,bi,bj),
198         &        hflux_x_SI(i,j,bi,bj)
199               ENDIF
200              ENDIF
201           ENDDO           ENDDO
202          ENDDO          ENDDO
203         ENDDO         ENDDO
204        ENDDO        ENDDO
205    
206    
207    
208  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
 CADJ STORE h_after_vflux_si = comlev1, key=ikey_dynamics  
209  CADJ STORE streamice_hmask  = comlev1, key=ikey_dynamics  CADJ STORE streamice_hmask  = comlev1, key=ikey_dynamics
210  #endif  #endif
211    
212        CALL STREAMICE_ADVECT_THICKNESS_Y ( myThid,  !      CALL STREAMICE_ADVECT_THICKNESS_Y ( myThid,
213    !     O   hflux_y_SI,
214    !     O   h_after_vflux_SI,
215    !     I   time_step_loc )
216    
217          CALL STREAMICE_ADV_FLUX_FL_Y ( myThid ,
218         I   vtrans ,
219         I   h_after_uflux_si ,
220         I   BCMASKY,
221         I   BCVALY,
222       O   hflux_y_SI,       O   hflux_y_SI,
      O   h_after_vflux_SI,  
223       I   time_step_loc )       I   time_step_loc )
224    
225    
226          DO bj=myByLo(myThid),myByHi(myThid)
227           DO bi=myBxLo(myThid),myBxHi(myThid)
228            DO j=1-1,sNy+1
229             DO i=1-1,sNx+1
230              IF (STREAMICE_hmask(i,j,bi,bj).eq.1.0) THEN
231               h_after_vflux_SI (i,j,bi,bj) = h_after_uflux_SI(i,j,bi,bj) -
232         &      (hflux_y_SI(i,j+1,bi,bj)*dxG(i,j+1,bi,bj) -
233         &        hflux_y_SI(i,j,bi,bj)*dxG(i,j,bi,bj)) *
234         &        recip_rA (i,j,bi,bj) * time_step_loc
235               IF ( h_after_vflux_SI (i,j,bi,bj).le.0.0) THEN
236                PRINT *, "h neg after y", i,j,hflux_y_SI(i,j+1,bi,bj),
237         &        hflux_y_SI(i,j,bi,bj)
238               ENDIF
239    
240              ENDIF
241             ENDDO
242            ENDDO
243           ENDDO
244          ENDDO
245    
246        DO bj=myByLo(myThid),myByHi(myThid)        DO bj=myByLo(myThid),myByHi(myThid)
247         DO bi=myBxLo(myThid),myBxHi(myThid)         DO bi=myBxLo(myThid),myBxHi(myThid)
248          DO j=1-OLy,sNy+OLy          DO j=1,sNy
249           DO i=1-OLx,sNx+OLx           DO i=1,sNx
250            IF (STREAMICE_hmask(i,j,bi,bj).eq.1.0) THEN            IF (STREAMICE_hmask(i,j,bi,bj).eq.1.0) THEN
251             H_streamice (i,j,bi,bj) =             H_streamice (i,j,bi,bj) =
252       &      h_after_vflux_SI (i,j,bi,bj)       &      h_after_vflux_SI (i,j,bi,bj)
# Line 120  CADJ STORE streamice_hmask  = comlev1, k Line 256  CADJ STORE streamice_hmask  = comlev1, k
256         ENDDO         ENDDO
257        ENDDO        ENDDO
258    
259    !      CALL STREAMICE_ADV_FRONT (
260    !     &  myThid, time_step_loc,
261    !     &  hflux_x_SI, hflux_y_SI )
262    
263    
264    #ifdef ALLOW_STREAMICE_2DTRACER
265          CALL STREAMICE_ADVECT_2DTRACER(
266         &  myThid,
267         &  myIter,
268         &  time_step,
269         &  uTrans,
270         &  vTrans,
271         &  bcMaskx,
272         &  bcMasky )
273    #endif
274    
275    #ifndef ALLOW_AUTODIFF_TAMC
276          time_step_rem = time_step_rem - time_step_loc      
277          enddo
278    #endif
279    
280    
       CALL STREAMICE_ADV_FRONT ( myThid, time_step_loc )  
281    
282    ! NOW WE APPLY MELT RATES !!
283    ! THIS MAY BE MOVED TO A SEPARATE SUBROUTINE
284    
285          DO bj=myByLo(myThid),myByHi(myThid)
286           DO bi=myBxLo(myThid),myBxHi(myThid)
287            DO j=1,sNy
288             DO i=1,sNx
289               IF (STREAMICE_hmask(i,j,bi,bj).eq.1.0 .or.
290         &       STREAMICE_hmask(i,j,bi,bj).eq.2.0) THEN
291                 MR = (1.-float_frac_streamice(i,j,bi,bj)) *
292         &             BDOT_STREAMICE(i,j,bi,bj)
293                 SMB = ADOT_STREAMICE(i,j,bi,bj)
294                 TMB = SMB - MR
295                 IF ((TMB.lt.0.0) .and.
296         &         (MR * time_step_loc .gt.
297         &          H_streamice (i,j,bi,bj))) THEN
298                    H_streamice (i,j,bi,bj) = 0. _d 0
299                    STREAMICE_hmask(i,j,bi,bj) = 0.
300                 ELSE
301                   H_streamice (i,j,bi,bj) =
302         &          H_streamice (i,j,bi,bj) + TMB * time_step_loc
303                 ENDIF
304               ENDIF
305             ENDDO
306            ENDDO
307           ENDDO
308          ENDDO    
309    
310        _EXCH_XY_RL( H_streamice, myThid )        _EXCH_XY_RL (H_streamice,myThid)
       _EXCH_XY_RL( area_shelf_streamice, myThid )  
       _EXCH_XY_RL( STREAMICE_hmask, myThid )  
311                
312        WRITE(msgBuf,'(A)') 'END STREAMICE_ADVECT_THICKNESS'        WRITE(msgBuf,'(A)') 'END STREAMICE_ADVECT_THICKNESS'
313         CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,         CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
Legend:
Removed from v.1.3  
changed lines
  Added in v.1.9
  ViewVC Help
Powered by ViewVC 1.1.22