/[MITgcm]/MITgcm_contrib/torge/itd/code/seaice_model.F
ViewVC logotype

Diff of /MITgcm_contrib/torge/itd/code/seaice_model.F

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

revision 1.1 by dimitri, Fri Apr 27 22:22:17 2012 UTC revision 1.11 by torge, Wed Mar 27 18:59:53 2013 UTC
# Line 55  CEndOfInterface Line 55  CEndOfInterface
55    
56  C !LOCAL VARIABLES: ====================================================  C !LOCAL VARIABLES: ====================================================
57  C     i,j,bi,bj :: Loop counters  C     i,j,bi,bj :: Loop counters
58  #if defined(SEAICE_GROWTH_LEGACY) || defined(ALLOW_AUTODIFF_TAMC)  #ifdef SEAICE_DEBUG
59    CToM<<<
60    C     msgBuf      :: Informational/error message buffer
61          CHARACTER*(MAX_LEN_MBUF) msgBuf
62          CHARACTER*10 HlimitMsgFormat
63    #endif
64    #ifdef SEAICE_ITD
65          INTEGER IT
66    #endif
67    #if defined(ALLOW_AUTODIFF_TAMC) || defined(SEAICE_ITD)
68    C>>>ToM
69        INTEGER i, j, bi, bj        INTEGER i, j, bi, bj
70  #endif  #endif
71  #ifdef ALLOW_SITRACER  #ifdef ALLOW_SITRACER
# Line 77  C--   Map thSice-variables to HEFF and A Line 87  C--   Map thSice-variables to HEFF and A
87        ENDIF        ENDIF
88  #endif /* ALLOW_THSICE */  #endif /* ALLOW_THSICE */
89    
 #ifdef SEAICE_GROWTH_LEGACY  
       IF ( .NOT.useThSice ) THEN  
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ STORE heff  = comlev1, key=ikey_dynamics, kind=isbyte  
 CADJ STORE heffm = comlev1, key=ikey_dynamics, kind=isbyte  
 CADJ STORE area  = comlev1, key=ikey_dynamics, kind=isbyte  
 CADJ STORE hsnow = comlev1, key=ikey_dynamics, kind=isbyte  
 CADJ STORE tice  = comlev1, key=ikey_dynamics, kind=isbyte  
 #ifdef SEAICE_VARIABLE_SALINITY  
 CADJ STORE hsalt = comlev1, key=ikey_dynamics, kind=isbyte  
 #endif  
 #endif  
       DO bj=myByLo(myThid),myByHi(myThid)  
        DO bi=myBxLo(myThid),myBxHi(myThid)  
         DO j=1-OLy,sNy+OLy  
          DO i=1-OLx,sNx+OLx  
           IF ( (heff(i,j,bi,bj).EQ.0.)  
      &     .OR.(area(i,j,bi,bj).EQ.0.)  
      &     ) THEN  
            HEFF(i,j,bi,bj) = 0. _d 0  
            AREA(i,j,bi,bj) = 0. _d 0  
            HSNOW(i,j,bi,bj) = 0. _d 0  
            TICE(i,j,bi,bj) = celsius2K  
 #ifdef SEAICE_VARIABLE_SALINITY  
            HSALT(i,j,bi,bj) = 0. _d 0  
 #endif  
           ENDIF  
          ENDDO  
         ENDDO  
        ENDDO  
       ENDDO  
       ENDIF  
 #endif  
   
90  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
91        DO bj=myByLo(myThid),myByHi(myThid)        DO bj=myByLo(myThid),myByHi(myThid)
92         DO bi=myBxLo(myThid),myBxHi(myThid)         DO bi=myBxLo(myThid),myBxHi(myThid)
93          DO j=1-OLy,sNy+OLy          DO j=1-OLy,sNy+OLy
94           DO i=1-OLx,sNx+OLx           DO i=1-OLx,sNx+OLx
 # ifdef SEAICE_GROWTH_LEGACY  
           areaNm1(i,j,bi,bj) = 0. _d 0  
           hEffNm1(i,j,bi,bj) = 0. _d 0  
 # endif  
95            uIceNm1(i,j,bi,bj) = 0. _d 0            uIceNm1(i,j,bi,bj) = 0. _d 0
96            vIceNm1(i,j,bi,bj) = 0. _d 0            vIceNm1(i,j,bi,bj) = 0. _d 0
97  # ifdef ALLOW_SITRACER  # ifdef ALLOW_SITRACER
# Line 189  C--   Apply ice velocity open boundary c Line 161  C--   Apply ice velocity open boundary c
161  #endif /* ALLOW_OBCS */  #endif /* ALLOW_OBCS */
162    
163  #ifdef ALLOW_THSICE  #ifdef ALLOW_THSICE
164        IF ( .NOT.useThSice ) THEN        IF ( useThSice ) THEN
165    #ifndef OLD_THSICE_CALL_SEQUENCE
166    #ifdef ALLOW_DEBUG
167            IF (debugMode) CALL DEBUG_CALL( 'THSICE_DO_ADVECT', myThid )
168    #endif
169            CALL THSICE_DO_ADVECT( 0, 0, myTime, myIter, myThid )
170    #endif /* OLD_THSICE_CALL_SEQUENCE */
171          ELSE
172  #endif  #endif
173  C--   Only call advection of heff, area, snow, and salt and  C--   Only call advection of heff, area, snow, and salt and
174  C--   growth for the generic 0-layer thermodynamics of seaice  C--   growth for the generic 0-layer thermodynamics of seaice
# Line 199  C--   (called from DO_OCEANIC_PHYSICS) t Line 178  C--   (called from DO_OCEANIC_PHYSICS) t
178  C NOW DO ADVECTION and DIFFUSION  C NOW DO ADVECTION and DIFFUSION
179        IF ( SEAICEadvHeff .OR. SEAICEadvArea .OR. SEAICEadvSnow        IF ( SEAICEadvHeff .OR. SEAICEadvArea .OR. SEAICEadvSnow
180       &        .OR. SEAICEadvSalt ) THEN       &        .OR. SEAICEadvSalt ) THEN
181    CToM<<<
182    #ifdef SEAICE_ITD
183    #ifdef SEAICE_DEBUG
184    C     ToM: generate some test output
185           WRITE(HlimitMsgFormat,'(A,I2,A)') '(A,',nITD,'F8.4)'
186           DO bj=myByLo(myThid),myByHi(myThid)
187            DO bi=myBxLo(myThid),myBxHi(myThid)
188               WRITE(msgBuf,'(A,F8.4,x,F8.4)')
189         &       ' SEAICE_MODEL: AREA and HEFF before advection: ',
190         &       AREA(1,1,bi,bj), HEFF(1,1,bi,bj)
191    c     &       ' SEAICE_MODEL: AREA and HEFF/AREA before advection: ',
192    c     &       AREA(1,1,bi,bj), HEFF(1,1,bi,bj)/AREA(1,1,bi,bj)
193               CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
194         &       SQUEEZE_RIGHT , myThid)
195               WRITE(msgBuf,HlimitMsgFormat)
196         &       ' SEAICE_MODEL: HEFFITD       before advection: ',
197         &       HEFFITD(1,1,:,bi,bj)
198    c     &       ' SEAICE_MODEL: HEFFITD/AREAITD before advection: ',
199    c     &       HEFFITD(1,1,:,bi,bj) / AREAITD(1,1,:,bi,bj)
200               CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
201         &       SQUEEZE_RIGHT , myThid)
202               WRITE(msgBuf,HlimitMsgFormat)
203         &       ' SEAICE_MODEL: AREAITD       before advection: ',
204         &       AREAITD(1,1,:,bi,bj)
205               CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
206         &       SQUEEZE_RIGHT , myThid)
207            ENDDO
208           ENDDO
209    #endif
210    #endif
211    C>>>ToM
212  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
213         IF (debugMode) CALL DEBUG_CALL( 'SEAICE_ADVDIFF', myThid )         IF (debugMode) CALL DEBUG_CALL( 'SEAICE_ADVDIFF', myThid )
214  #endif  #endif
215         CALL SEAICE_ADVDIFF( myTime, myIter, myThid )         CALL SEAICE_ADVDIFF( myTime, myIter, myThid )
216  #ifdef SEAICE_GROWTH_LEGACY  CToM<<<
217        ELSE  #ifdef SEAICE_ITD
218    #ifdef SEAICE_DEBUG
219    C     ToM: generate some test output
220           WRITE(HlimitMsgFormat,'(A,I2,A)') '(A,',nITD,'F8.4)'
221         DO bj=myByLo(myThid),myByHi(myThid)         DO bj=myByLo(myThid),myByHi(myThid)
222          DO bi=myBxLo(myThid),myBxHi(myThid)          DO bi=myBxLo(myThid),myBxHi(myThid)
223           DO j=1-OLy,sNy+OLy             WRITE(msgBuf,HlimitMsgFormat)
224            DO i=1-OLx,sNx+OLx       &       ' SEAICE_MODEL: HEFFITD        after advection: ',
225             areaNm1(i,j,bi,bj) = AREA(i,j,bi,bj)       &       HEFFITD(1,1,:,bi,bj)
226             hEffNm1(i,j,bi,bj) = HEFF(i,j,bi,bj)             CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
227            ENDDO       &       SQUEEZE_RIGHT , myThid)
228           ENDDO             WRITE(msgBuf,HlimitMsgFormat)
229         &       ' SEAICE_MODEL: AREAITD        after advection: ',
230         &       AREAITD(1,1,:,bi,bj)
231               CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
232         &       SQUEEZE_RIGHT , myThid)
233                WRITE(msgBuf,'(A)')
234         &        ' --------------------------------------------- '
235                CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
236         &        SQUEEZE_RIGHT , myThid)
237            ENDDO
238           ENDDO
239    #endif
240    C
241    C     check that all ice thickness categories meet their limits
242    C     (includes Hibler-type ridging)
243    #ifdef ALLOW_DEBUG
244            IF (debugMode) CALL DEBUG_CALL( 'SEAICE_ITD_REDIST', myThid )
245    #endif
246           DO bj=myByLo(myThid),myByHi(myThid)
247            DO bi=myBxLo(myThid),myBxHi(myThid)
248             CALL SEAICE_ITD_REDIST(bi, bj, myTime, myIter, myThid)
249          ENDDO          ENDDO
250         ENDDO         ENDDO
251  #endif /* SEAICE_GROWTH_LEGACY */  C     update mean ice thickness HEFF and total ice concentration AREA
252    C     to match single category values
253    #ifdef ALLOW_DEBUG
254            IF (debugMode) CALL DEBUG_CALL( 'SEAICE_ITD_SUM', myThid )
255    #endif
256           DO bj=myByLo(myThid),myByHi(myThid)
257            DO bi=myBxLo(myThid),myBxHi(myThid)
258             CALL SEAICE_ITD_SUM(bi, bj, myTime, myIter, myThid)
259            ENDDO
260           ENDDO
261    #endif
262    C>>>ToM
263        ENDIF        ENDIF
264  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
265  CADJ STORE heffm  = comlev1, key=ikey_dynamics, kind=isbyte  CADJ STORE heffm  = comlev1, key=ikey_dynamics, kind=isbyte
266  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
267    
268    #ifdef SEAICE_ITD
269    #ifdef SEAICE_DEBUG
270    C     ToM: generate some test output
271           WRITE(HlimitMsgFormat,'(A,I2,A)') '(A,',nITD,'F8.4)'
272           DO bj=myByLo(myThid),myByHi(myThid)
273            DO bi=myBxLo(myThid),myBxHi(myThid)
274               WRITE(msgBuf,HlimitMsgFormat)
275         &       ' SEAICE_MODEL: HEFFITD          before growth: ',
276         &       HEFFITD(1,1,:,bi,bj)
277               CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
278         &       SQUEEZE_RIGHT , myThid)
279               WRITE(msgBuf,HlimitMsgFormat)
280         &       ' SEAICE_MODEL: AREAITD          before growth: ',
281         &       AREAITD(1,1,:,bi,bj)
282               CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
283         &       SQUEEZE_RIGHT , myThid)
284               WRITE(msgBuf,HlimitMsgFormat)
285         &       ' SEAICE_MODEL: HSNOWITD         before growth: ',
286         &       HSNOWITD(1,1,:,bi,bj)
287               CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
288         &       SQUEEZE_RIGHT , myThid)
289            ENDDO
290           ENDDO
291    #endif
292    #endif
293    
294  #ifndef DISABLE_SEAICE_GROWTH  #ifndef DISABLE_SEAICE_GROWTH
295  C     thermodynamics growth  C     thermodynamics growth
296  C     must call growth after calling advection  C     must call growth after calling advection
# Line 230  C     because of ugly time level busines Line 300  C     because of ugly time level busines
300          IF (debugMode) CALL DEBUG_CALL( 'SEAICE_GROWTH', myThid )          IF (debugMode) CALL DEBUG_CALL( 'SEAICE_GROWTH', myThid )
301  #endif  #endif
302          CALL SEAICE_GROWTH( myTime, myIter, myThid )          CALL SEAICE_GROWTH( myTime, myIter, myThid )
303    CToM<<<
304    #ifdef SEAICE_ITD
305    #ifdef SEAICE_DEBUG
306    C     ToM: generate some test output
307            WRITE(HlimitMsgFormat,'(A,I2,A)') '(A,',nITD,'F8.4)'
308            DO bj=myByLo(myThid),myByHi(myThid)
309             DO bi=myBxLo(myThid),myBxHi(myThid)
310                WRITE(msgBuf,HlimitMsgFormat)
311         &        ' SEAICE_MODEL: HEFFITD           after growth: ',
312         &        HEFFITD(1,1,:,bi,bj)
313                CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
314         &        SQUEEZE_RIGHT , myThid)
315                WRITE(msgBuf,HlimitMsgFormat)
316         &        ' SEAICE_MODEL: AREAITD           after growth: ',
317         &        AREAITD(1,1,:,bi,bj)
318                CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
319         &        SQUEEZE_RIGHT , myThid)
320                WRITE(msgBuf,HlimitMsgFormat)
321         &        ' SEAICE_MODEL: HSNOWITD          after growth: ',
322         &        HSNOWITD(1,1,:,bi,bj)
323                CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
324         &        SQUEEZE_RIGHT , myThid)
325                WRITE(msgBuf,'(A)')
326         &        ' --------------------------------------------- '
327                CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
328         &        SQUEEZE_RIGHT , myThid)
329             ENDDO
330            ENDDO
331    #endif
332    C
333    C     redistribute sea ice into proper sea ice category after growth/melt
334    C     in case model runs with ice thickness distribution
335    C---+-|--1----+----2----+----3----+----4----+----5----+----6----+----7-|
336    #ifdef ALLOW_DEBUG
337            IF (debugMode) CALL DEBUG_CALL( 'SEAICE_ITD_REDIST', myThid )
338    #endif
339           DO bj=myByLo(myThid),myByHi(myThid)
340            DO bi=myBxLo(myThid),myBxHi(myThid)
341             CALL SEAICE_ITD_REDIST(bi, bj, myTime, myIter, myThid)
342            ENDDO
343           ENDDO
344    C     store the mean ice thickness in HEFF (for dynamic solver and diagnostics)
345           DO bj=myByLo(myThid),myByHi(myThid)
346            DO bi=myBxLo(myThid),myBxHi(myThid)
347             CALL SEAICE_ITD_SUM(bi, bj, myTime, myIter, myThid)
348            ENDDO
349           ENDDO
350    
351    #ifdef SEAICE_DEBUG
352    C     ToM: generate some test output
353            WRITE(HlimitMsgFormat,'(A,I2,A)') '(A,',nITD,'F8.4)'
354            DO bj=myByLo(myThid),myByHi(myThid)
355             DO bi=myBxLo(myThid),myBxHi(myThid)
356                WRITE(msgBuf,HlimitMsgFormat)
357         &        ' SEAICE_MODEL: HEFFITD    after final sorting: ',
358         &        HEFFITD(1,1,:,bi,bj)
359                CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
360         &        SQUEEZE_RIGHT , myThid)
361                WRITE(msgBuf,HlimitMsgFormat)
362         &        ' SEAICE_MODEL: AREAITD    after final sorting: ',
363         &        AREAITD(1,1,:,bi,bj)
364                CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
365         &        SQUEEZE_RIGHT , myThid)
366                WRITE(msgBuf,HlimitMsgFormat)
367         &        ' SEAICE_MODEL: HSNOWITD   after final sorting: ',
368         &        HSNOWITD(1,1,:,bi,bj)
369                CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
370         &        SQUEEZE_RIGHT , myThid)
371                WRITE(msgBuf,'(A)')
372         &        ' ============================================= '
373                CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
374         &        SQUEEZE_RIGHT , myThid)
375             ENDDO
376            ENDDO
377    #endif
378    #endif
379    C
380    C>>>ToM
381         ENDIF         ENDIF
382  #endif /* DISABLE_SEAICE_GROWTH */  #endif /* DISABLE_SEAICE_GROWTH */
383    
# Line 304  C     the seaice model Line 452  C     the seaice model
452          CALL DIAGNOSTICS_FILL_RS(EmPmR,'SIempmr ',0,1 ,0,1,1,myThid)          CALL DIAGNOSTICS_FILL_RS(EmPmR,'SIempmr ',0,1 ,0,1,1,myThid)
453          CALL DIAGNOSTICS_FILL_RS(Qnet ,'SIqnet  ',0,1 ,0,1,1,myThid)          CALL DIAGNOSTICS_FILL_RS(Qnet ,'SIqnet  ',0,1 ,0,1,1,myThid)
454          CALL DIAGNOSTICS_FILL_RS(Qsw  ,'SIqsw   ',0,1 ,0,1,1,myThid)          CALL DIAGNOSTICS_FILL_RS(Qsw  ,'SIqsw   ',0,1 ,0,1,1,myThid)
455    #ifdef SEAICE_ITD
456            CALL DIAGNOSTICS_FILL(HEFFITD ,'SIheffN ',0,nITD,0,1,1,myThid)
457            CALL DIAGNOSTICS_FILL(AREAITD ,'SIareaN ',0,nITD,0,1,1,myThid)
458    #endif
459         ENDIF         ENDIF
460  #endif /* ALLOW_DIAGNOSTICS */  #endif /* ALLOW_DIAGNOSTICS */
461    

Legend:
Removed from v.1.1  
changed lines
  Added in v.1.11

  ViewVC Help
Powered by ViewVC 1.1.22