/[MITgcm]/MITgcm/pkg/kpp/kpp_calc.F
ViewVC logotype

Diff of /MITgcm/pkg/kpp/kpp_calc.F

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

revision 1.41 by mlosch, Thu May 3 14:51:05 2007 UTC revision 1.42 by jmc, Thu May 3 21:37:34 2007 UTC
# Line 7  CBOP Line 7  CBOP
7  C !ROUTINE: KPP_CALC  C !ROUTINE: KPP_CALC
8    
9  C !INTERFACE: ==========================================================  C !INTERFACE: ==========================================================
10        subroutine KPP_CALC(        SUBROUTINE KPP_CALC(
11       I     bi, bj, myTime, myThid )       I     bi, bj, myTime, myIter, myThid )
12    
13  C !DESCRIPTION: \bv  C !DESCRIPTION: \bv
14  C     /==========================================================\  C     *==========================================================*
15  C     | SUBROUTINE KPP_CALC                                      |  C     | SUBROUTINE KPP_CALC                                      |
16  C     | o Compute all KPP fields defined in KPP.h                |  C     | o Compute all KPP fields defined in KPP.h                |
17  C     |==========================================================|  C     *==========================================================*
18  C     | This subroutine serves as an interface between MITGCMUV  |  C     | This subroutine serves as an interface between MITGCMUV  |
19  C     | code and NCOM 1-D routines in kpp_routines.F             |  C     | code and NCOM 1-D routines in kpp_routines.F             |
20  C     \==========================================================/  C     *==========================================================*
21        IMPLICIT NONE        IMPLICIT NONE
22    
23  c=======================================================================  c=======================================================================
# Line 52  c         within the boundary layer, abo Line 52  c         within the boundary layer, abo
52  c         determined by turbulent surface fluxes, and interior mixing at  c         determined by turbulent surface fluxes, and interior mixing at
53  c         the lower boundary, i.e. at hbl.  c         the lower boundary, i.e. at hbl.
54  c      c    
55  c     this subroutine provides the interface between the MIT GCM UV and the  c     this subroutine provides the interface between the MITGCM and
56  c     subroutine "kppmix", where boundary layer depth, vertical  c     the routine "kppmix", where boundary layer depth, vertical
57  c     viscosity, vertical diffusivity, and counter gradient term (ghat)  c     viscosity, vertical diffusivity, and counter gradient term (ghat)
58  c     are computed slabwise.  c     are computed slabwise.
59  c     note: subroutine "kppmix" uses m-k-s units.  c     note: subroutine "kppmix" uses m-k-s units.
# Line 120  C !USES: =============================== Line 120  C !USES: ===============================
120    
121  C !INPUT PARAMETERS: ===================================================  C !INPUT PARAMETERS: ===================================================
122  c Routine arguments  c Routine arguments
123  c     bi, bj - array indices on which to apply calculations  c     bi, bj :: Current tile indices
124  c     myTime - Current time in simulation  c     myTime :: Current time in simulation
125    c     myIter :: Current iteration number in simulation
126    c     myThid :: My Thread Id. number
127    
128        INTEGER bi, bj        INTEGER bi, bj
       INTEGER myThid  
129        _RL     myTime        _RL     myTime
130          INTEGER myIter
131          INTEGER myThid
132    
133  #ifdef ALLOW_KPP  #ifdef ALLOW_KPP
134    
# Line 189  cph) Line 192  cph)
192        _RL uRef  ( 1-OLx:sNx+OLx, 1-OLy:sNy+OLy                )        _RL uRef  ( 1-OLx:sNx+OLx, 1-OLy:sNy+OLy                )
193        _RL vRef  ( 1-OLx:sNx+OLx, 1-OLy:sNy+OLy                )        _RL vRef  ( 1-OLx:sNx+OLx, 1-OLy:sNy+OLy                )
194  #endif /* KPP_ESTIMATE_UREF */  #endif /* KPP_ESTIMATE_UREF */
195          
196        _RL tempvar2        _RL tempvar2
197        integer i, j, k, kp1, km1, im1, ip1, jm1, jp1        integer i, j, k, kp1, km1, im1, ip1, jm1, jp1
198    
# Line 214  CEOP Line 217  CEOP
217  c     Check to see if new vertical mixing coefficient should be computed now?  c     Check to see if new vertical mixing coefficient should be computed now?
218        IF ( DIFFERENT_MULTIPLE(kpp_freq,myTime,deltaTClock)        IF ( DIFFERENT_MULTIPLE(kpp_freq,myTime,deltaTClock)
219       1     .OR. myTime .EQ. startTime ) THEN       1     .OR. myTime .EQ. startTime ) THEN
220            
221  c-----------------------------------------------------------------------  c-----------------------------------------------------------------------
222  c     prepare input arrays for subroutine "kppmix" to compute  c     prepare input arrays for subroutine "kppmix" to compute
223  c     viscosity and diffusivity and ghat.  c     viscosity and diffusivity and ghat.
# Line 273  c     levels therefore k+1 mask must be Line 276  c     levels therefore k+1 mask must be
276        DO k = 1, Nr-1        DO k = 1, Nr-1
277           CALL SMOOTH_HORIZ (           CALL SMOOTH_HORIZ (
278       I        k+1, bi, bj,       I        k+1, bi, bj,
279       U        ghat (1-OLx,1-OLy,k),       U        ghat (1-OLx,1-OLy,k),
280       I        myThid )       I        myThid )
281        ENDDO        ENDDO
282    
# Line 283  c     levels therefore k+1 mask must be Line 286  c     levels therefore k+1 mask must be
286  c     horizontally smooth density related quantities with 121 filters  c     horizontally smooth density related quantities with 121 filters
287        CALL SMOOTH_HORIZ (        CALL SMOOTH_HORIZ (
288       I     1, bi, bj,       I     1, bi, bj,
289       U     work2,       U     work2,
290       I     myThid )       I     myThid )
291        DO k = 1, Nr        DO k = 1, Nr
292           CALL SMOOTH_HORIZ (           CALL SMOOTH_HORIZ (
293       I        k+1, bi, bj,       I        k+1, bi, bj,
294       U        dbloc (1-OLx,1-OLy,k),       U        dbloc (1-OLx,1-OLy,k),
295       I        myThid )       I        myThid )
296           CALL SMOOTH_HORIZ (           CALL SMOOTH_HORIZ (
297       I        k, bi, bj,       I        k, bi, bj,
298       U        Ritop (1-OLx,1-OLy,k),       U        Ritop (1-OLx,1-OLy,k),
299       I        myThid )       I        myThid )
300           CALL SMOOTH_HORIZ (           CALL SMOOTH_HORIZ (
301       I        k, bi, bj,       I        k, bi, bj,
302       U        TTALPHA(1-OLx,1-OLy,k),       U        TTALPHA(1-OLx,1-OLy,k),
303       I        myThid )       I        myThid )
304           CALL SMOOTH_HORIZ (           CALL SMOOTH_HORIZ (
305       I        k, bi, bj,       I        k, bi, bj,
306       U        SSBETA(1-OLx,1-OLy,k),       U        SSBETA(1-OLx,1-OLy,k),
307       I        myThid )       I        myThid )
308        ENDDO        ENDDO
309  #endif /* KPP_SMOOTH_DENS */  #endif /* KPP_SMOOTH_DENS */
# Line 329  c     note: land and ocean bottom values Line 332  c     note: land and ocean bottom values
332  c     so that the subroutine "bldepth" works correctly  c     so that the subroutine "bldepth" works correctly
333                 Ritop(i,j,k) = (zgrid(1)-zgrid(k)) * Ritop(i,j,k)                 Ritop(i,j,k) = (zgrid(1)-zgrid(k)) * Ritop(i,j,k)
334    
335              END DO              ENDDO
336           END DO           ENDDO
337        END DO        ENDDO
338    
339  cph(  cph(
340  cph  this avoids a single or double recomp./call of statekpp  cph  this avoids a single or double recomp./call of statekpp
# Line 357  CML        surfForcT = surfaceForcingT(i Line 360  CML        surfForcT = surfaceForcingT(i
360  CML     &       + shelficeForcingT(i,j,bi,bj) * shelfIceFac  CML     &       + shelficeForcingT(i,j,bi,bj) * shelfIceFac
361  CML        surfForcS = surfaceForcingS(i,j,bi,bj)  CML        surfForcS = surfaceForcingS(i,j,bi,bj)
362  CML     &       + shelficeForcingS(i,j,bi,bj) * shelfIceFac  CML     &       + shelficeForcingS(i,j,bi,bj) * shelfIceFac
363  CML       END DO  CML       ENDDO
364  CML      END DO  CML      ENDDO
365  CML#endif /* ALLOW_SHELFICE */  CML#endif /* ALLOW_SHELFICE */
366    
367  c------------------------------------------------------------------------  c------------------------------------------------------------------------
# Line 382  c     note: Vref can depend on the surfa Line 385  c     note: Vref can depend on the surfa
385  c     dVsq in the subroutine that does the surface related stuff  c     dVsq in the subroutine that does the surface related stuff
386  c     (admittedly this is a bit messy)  c     (admittedly this is a bit messy)
387  c------------------------------------------------------------------------  c------------------------------------------------------------------------
388          
389        CALL KPP_FORCING_SURF(        CALL KPP_FORCING_SURF(
390       I     work2, surfaceForcingU, surfaceForcingV,       I     work2, surfaceForcingU, surfaceForcingV,
391       I     surfaceForcingT, surfaceForcingS, surfaceForcingTice,       I     surfaceForcingT, surfaceForcingS, surfaceForcingTice,
392       I     Qsw, ttalpha, ssbeta,         I     Qsw, ttalpha, ssbeta,
393       O     ustar, bo, bosol, dVsq,       O     ustar, bo, bosol, dVsq,
394       I     ikppkey, iMin, iMax, jMin, jMax, bi, bj, myTime, myThid )       I     ikppkey, iMin, iMax, jMin, jMax, bi, bj, myTime, myThid )
395    
# Line 399  c initialize arrays to zero Line 402  c initialize arrays to zero
402         DO j = 1-OLy, sNy+OLy         DO j = 1-OLy, sNy+OLy
403          DO i = 1-OLx, sNx+OLx          DO i = 1-OLx, sNx+OLx
404           shsq(i,j,k) = p0           shsq(i,j,k) = p0
405          END DO          ENDDO
406         END DO         ENDDO
407        END DO        ENDDO
408    
409  c     shsq computation  c     shsq computation
410        DO k = 1, Nrm1        DO k = 1, Nrm1
# Line 413  c     shsq computation Line 416  c     shsq computation
416           im1 = i - 1           im1 = i - 1
417           ip1 = i + 1           ip1 = i + 1
418           shsq(i,j,k) = p5 * (           shsq(i,j,k) = p5 * (
419       $        (uVel(i,  j,  k,bi,bj)-uVel(i,  j,  kp1,bi,bj)) *       &        (uVel(i,  j,  k,bi,bj)-uVel(i,  j,  kp1,bi,bj)) *
420       $        (uVel(i,  j,  k,bi,bj)-uVel(i,  j,  kp1,bi,bj)) +       &        (uVel(i,  j,  k,bi,bj)-uVel(i,  j,  kp1,bi,bj)) +
421       $        (uVel(ip1,j,  k,bi,bj)-uVel(ip1,j,  kp1,bi,bj)) *       &        (uVel(ip1,j,  k,bi,bj)-uVel(ip1,j,  kp1,bi,bj)) *
422       $        (uVel(ip1,j,  k,bi,bj)-uVel(ip1,j,  kp1,bi,bj)) +       &        (uVel(ip1,j,  k,bi,bj)-uVel(ip1,j,  kp1,bi,bj)) +
423       $        (vVel(i,  j,  k,bi,bj)-vVel(i,  j,  kp1,bi,bj)) *       &        (vVel(i,  j,  k,bi,bj)-vVel(i,  j,  kp1,bi,bj)) *
424       $        (vVel(i,  j,  k,bi,bj)-vVel(i,  j,  kp1,bi,bj)) +       &        (vVel(i,  j,  k,bi,bj)-vVel(i,  j,  kp1,bi,bj)) +
425       $        (vVel(i,  jp1,k,bi,bj)-vVel(i,  jp1,kp1,bi,bj)) *       &        (vVel(i,  jp1,k,bi,bj)-vVel(i,  jp1,kp1,bi,bj)) *
426       $        (vVel(i,  jp1,k,bi,bj)-vVel(i,  jp1,kp1,bi,bj)) )       &        (vVel(i,  jp1,k,bi,bj)-vVel(i,  jp1,kp1,bi,bj)) )
427  #ifdef KPP_SMOOTH_SHSQ  #ifdef KPP_SMOOTH_SHSQ
428           shsq(i,j,k) = p5 * shsq(i,j,k) + p125 * (           shsq(i,j,k) = p5 * shsq(i,j,k) + p125 * (
429       $        (uVel(i,  jm1,k,bi,bj)-uVel(i,  jm1,kp1,bi,bj)) *       &        (uVel(i,  jm1,k,bi,bj)-uVel(i,  jm1,kp1,bi,bj)) *
430       $        (uVel(i,  jm1,k,bi,bj)-uVel(i,  jm1,kp1,bi,bj)) +       &        (uVel(i,  jm1,k,bi,bj)-uVel(i,  jm1,kp1,bi,bj)) +
431       $        (uVel(ip1,jm1,k,bi,bj)-uVel(ip1,jm1,kp1,bi,bj)) *       &        (uVel(ip1,jm1,k,bi,bj)-uVel(ip1,jm1,kp1,bi,bj)) *
432       $        (uVel(ip1,jm1,k,bi,bj)-uVel(ip1,jm1,kp1,bi,bj)) +       &        (uVel(ip1,jm1,k,bi,bj)-uVel(ip1,jm1,kp1,bi,bj)) +
433       $        (uVel(i,  jp1,k,bi,bj)-uVel(i,  jp1,kp1,bi,bj)) *       &        (uVel(i,  jp1,k,bi,bj)-uVel(i,  jp1,kp1,bi,bj)) *
434       $        (uVel(i,  jp1,k,bi,bj)-uVel(i,  jp1,kp1,bi,bj)) +       &        (uVel(i,  jp1,k,bi,bj)-uVel(i,  jp1,kp1,bi,bj)) +
435       $        (uVel(ip1,jp1,k,bi,bj)-uVel(ip1,jp1,kp1,bi,bj)) *       &        (uVel(ip1,jp1,k,bi,bj)-uVel(ip1,jp1,kp1,bi,bj)) *
436       $        (uVel(ip1,jp1,k,bi,bj)-uVel(ip1,jp1,kp1,bi,bj)) +       &        (uVel(ip1,jp1,k,bi,bj)-uVel(ip1,jp1,kp1,bi,bj)) +
437       $        (vVel(im1,j,  k,bi,bj)-vVel(im1,j,  kp1,bi,bj)) *       &        (vVel(im1,j,  k,bi,bj)-vVel(im1,j,  kp1,bi,bj)) *
438       $        (vVel(im1,j,  k,bi,bj)-vVel(im1,j,  kp1,bi,bj)) +       &        (vVel(im1,j,  k,bi,bj)-vVel(im1,j,  kp1,bi,bj)) +
439       $        (vVel(im1,jp1,k,bi,bj)-vVel(im1,jp1,kp1,bi,bj)) *       &        (vVel(im1,jp1,k,bi,bj)-vVel(im1,jp1,kp1,bi,bj)) *
440       $        (vVel(im1,jp1,k,bi,bj)-vVel(im1,jp1,kp1,bi,bj)) +       &        (vVel(im1,jp1,k,bi,bj)-vVel(im1,jp1,kp1,bi,bj)) +
441       $        (vVel(ip1,j,  k,bi,bj)-vVel(ip1,j,  kp1,bi,bj)) *       &        (vVel(ip1,j,  k,bi,bj)-vVel(ip1,j,  kp1,bi,bj)) *
442       $        (vVel(ip1,j,  k,bi,bj)-vVel(ip1,j,  kp1,bi,bj)) +       &        (vVel(ip1,j,  k,bi,bj)-vVel(ip1,j,  kp1,bi,bj)) +
443       $        (vVel(ip1,jp1,k,bi,bj)-vVel(ip1,jp1,kp1,bi,bj)) *       &        (vVel(ip1,jp1,k,bi,bj)-vVel(ip1,jp1,kp1,bi,bj)) *
444       $        (vVel(ip1,jp1,k,bi,bj)-vVel(ip1,jp1,kp1,bi,bj)) )       &        (vVel(ip1,jp1,k,bi,bj)-vVel(ip1,jp1,kp1,bi,bj)) )
445  #endif  #endif
446          END DO          ENDDO
447         END DO         ENDDO
448        END DO        ENDDO
449    
450  cph(  cph(
451  #ifdef KPP_AUTODIFF_EXCESSIVE_STORE  #ifdef KPP_AUTODIFF_EXCESSIVE_STORE
# Line 471  c     matching diffusivities at bottom o Line 474  c     matching diffusivities at bottom o
474           DO i = 1-OLx, sNx+OLx           DO i = 1-OLx, sNx+OLx
475              work1(i,j) = nzmax(i,j,bi,bj)              work1(i,j) = nzmax(i,j,bi,bj)
476              work2(i,j) = Fcori(i,j,bi,bj)              work2(i,j) = Fcori(i,j,bi,bj)
477           END DO           ENDDO
478        END DO        ENDDO
479        CALL TIMER_START('KPPMIX [KPP_CALC]', myThid)        CALL TIMER_START('KPPMIX [KPP_CALC]', myThid)
480        CALL KPPMIX (        CALL KPPMIX (
481       I       work1, shsq, dVsq, ustar       I       work1, shsq, dVsq, ustar
# Line 503  c--------------------------------------- Line 506  c---------------------------------------
506       &        * maskC(i,j,km1,bi,bj)       &        * maskC(i,j,km1,bi,bj)
507           KPPghat(i,j,k,bi,bj)   = ghat(i,j,k)       * maskC(i,j,k,bi,bj)           KPPghat(i,j,k,bi,bj)   = ghat(i,j,k)       * maskC(i,j,k,bi,bj)
508       &        * maskC(i,j,km1,bi,bj)       &        * maskC(i,j,km1,bi,bj)
509          END DO          ENDDO
510          k = 1          k = 1
511  #ifdef ALLOW_SHELFICE  #ifdef ALLOW_SHELFICE
512          if ( useShelfIce ) k = kTopC(i,j,bi,bj)          if ( useShelfIce ) k = kTopC(i,j,bi,bj)
513  #endif /* ALLOW_SHELFICE */  #endif /* ALLOW_SHELFICE */
514          KPPhbl(i,j,bi,bj) = hbl(i,j) * maskC(i,j,k,bi,bj)          KPPhbl(i,j,bi,bj) = hbl(i,j) * maskC(i,j,k,bi,bj)
515    
516         END DO         ENDDO
517        END DO        ENDDO
518    
519  #ifdef KPP_SMOOTH_VISC  #ifdef KPP_SMOOTH_VISC
520  c     horizontal smoothing of vertical viscosity  c     horizontal smoothing of vertical viscosity
521        DO k = 1, Nr        DO k = 1, Nr
522           CALL SMOOTH_HORIZ (           CALL SMOOTH_HORIZ (
523       I        k, bi, bj,       I        k, bi, bj,
524       U        KPPviscAz(1-OLx,1-OLy,k,bi,bj),       U        KPPviscAz(1-OLx,1-OLy,k,bi,bj),
525       I        myThid )       I        myThid )
526        END DO        ENDDO
527        _EXCH_XYZ_R8(KPPviscAz  , myThid )  C jmc: No EXCH inside bi,bj loop !!!
528    c     _EXCH_XYZ_R8(KPPviscAz  , myThid )
529  #endif /* KPP_SMOOTH_VISC */  #endif /* KPP_SMOOTH_VISC */
530    
531  #ifdef KPP_SMOOTH_DIFF  #ifdef KPP_SMOOTH_DIFF
# Line 529  c     horizontal smoothing of vertical d Line 533  c     horizontal smoothing of vertical d
533        DO k = 1, Nr        DO k = 1, Nr
534           CALL SMOOTH_HORIZ (           CALL SMOOTH_HORIZ (
535       I        k, bi, bj,       I        k, bi, bj,
536       U        KPPdiffKzS(1-OLx,1-OLy,k,bi,bj),       U        KPPdiffKzS(1-OLx,1-OLy,k,bi,bj),
537       I        myThid )       I        myThid )
538           CALL SMOOTH_HORIZ (           CALL SMOOTH_HORIZ (
539       I        k, bi, bj,       I        k, bi, bj,
540       U        KPPdiffKzT(1-OLx,1-OLy,k,bi,bj),       U        KPPdiffKzT(1-OLx,1-OLy,k,bi,bj),
541       I        myThid )       I        myThid )
542        END DO        ENDDO
       _EXCH_XYZ_R8(KPPdiffKzS , myThid )  
       _EXCH_XYZ_R8(KPPdiffKzT , myThid )  
543  #endif /* KPP_SMOOTH_DIFF */  #endif /* KPP_SMOOTH_DIFF */
544    
545  cph(  cph(
# Line 554  C     the bottom of the mixing layer. Line 556  C     the bottom of the mixing layer.
556        ENDDO        ENDDO
557        CALL SWFRAC(        CALL SWFRAC(
558       I     (sNx+2*OLx)*(sNy+2*OLy), minusone,       I     (sNx+2*OLx)*(sNy+2*OLy), minusone,
559       I     mytime, mythid,       U     worka,
560       U     worka )       I     myTime, myIter, myThid )
561        DO j=1-OLy,sNy+OLy        DO j=1-OLy,sNy+OLy
562           DO i=1-OLx,sNx+OLx           DO i=1-OLx,sNx+OLx
563              KPPfrac(i,j,bi,bj) = worka(i,j)              KPPfrac(i,j,bi,bj) = worka(i,j)
# Line 569  C     the bottom of the mixing layer. Line 571  C     the bottom of the mixing layer.
571        RETURN        RETURN
572        END        END
573    
574        subroutine KPP_CALC_DUMMY(        SUBROUTINE KPP_CALC_DUMMY(
575       I     bi, bj, myTime, myThid )       I     bi, bj, myTime, myIter, myThid )
576  C     /==========================================================\  C     *==========================================================*
577  C     | SUBROUTINE KPP_CALC_DUMMY                                |  C     | SUBROUTINE KPP_CALC_DUMMY                                |
578  C     | o Compute all KPP fields defined in KPP.h                |  C     | o Compute all KPP fields defined in KPP.h                |
579  C     | o Dummy routine for TAMC  C     | o Dummy routine for TAMC
580  C     |==========================================================|  C     *==========================================================*
581  C     | This subroutine serves as an interface between MITGCMUV  |  C     | This subroutine serves as an interface between MITGCMUV  |
582  C     | code and NCOM 1-D routines in kpp_routines.F             |  C     | code and NCOM 1-D routines in kpp_routines.F             |
583  C     \==========================================================/  C     *==========================================================*
584        IMPLICIT NONE        IMPLICIT NONE
585    
586  #include "SIZE.h"  #include "SIZE.h"
# Line 590  C     \================================= Line 592  C     \=================================
592  #include "GAD.h"  #include "GAD.h"
593    
594  c Routine arguments  c Routine arguments
595  c     bi, bj - array indices on which to apply calculations  c     bi, bj :: Current tile indices
596  c     myTime - Current time in simulation  c     myTime :: Current time in simulation
597    c     myIter :: Current iteration number in simulation
598    c     myThid :: My Thread Id. number
599    
600        INTEGER bi, bj        INTEGER bi, bj
       INTEGER myThid  
601        _RL     myTime        _RL     myTime
602          INTEGER myIter
603          INTEGER myThid
604    
605  #ifdef ALLOW_KPP  #ifdef ALLOW_KPP
606    
Legend:
Removed from v.1.41  
changed lines
  Added in v.1.42
  ViewVC Help
Powered by ViewVC 1.1.22