model/src/calc_3d_diffusivity.F

C $Header: /u/gcmpack/MITgcm/model/src/calc_3d_diffusivity.F,v 1.13 2007/11/05 18:52:21 jmc Exp $
C $Name:  $

#include "PACKAGES_CONFIG.h"
#include "CPP_OPTIONS.h"

CBOP
C     !ROUTINE: CALC_3D_DIFFUSIVITY
C     !INTERFACE:
      SUBROUTINE CALC_3D_DIFFUSIVITY(
     I        bi,bj,iMin,iMax,jMin,jMax,
     I        trIdentity, trUseGMRedi, trUseKPP,
     O        KappaRTr,
     I        myThid)

C     !DESCRIPTION: \bv
C     *==========================================================*
C     | SUBROUTINE CALC_3D_DIFFUSIVITY
C     | o Calculate net (3D) vertical diffusivity for 1 tracer
C     *==========================================================*
C     | Combines spatially varying diffusion coefficients from
C     | KPP and/or GM and/or convective stability test.
C     *==========================================================*
C     \ev

C     !USES:
      IMPLICIT NONE
C     == GLobal variables ==
#include "SIZE.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "DYNVARS.h"
#include "GRID.h"
#ifdef ALLOW_GENERIC_ADVDIFF
#include "GAD.h"
#endif
#ifdef ALLOW_PTRACERS
#include "PTRACERS_SIZE.h"
#include "PTRACERS_PARAMS.h"
#endif

C     !INPUT/OUTPUT PARAMETERS:
C     == Routine arguments ==
C     bi, bj     :: tile indices
C     iMin,iMax  :: Range of points for which calculation is performed.
C     jMin,jMax  :: Range of points for which calculation is performed.
C     trIdentity :: tracer identifier
C     trUseGMRedi:: this tracer use GM-Redi
C     trUseKPP   :: this tracer use KPP
C     myThid     :: Instance number for this innvocation of CALC_3D_DIFFUSIVITY
C     KappaRTr   :: Net diffusivity for this tracer (trIdentity)
      INTEGER bi,bj,iMin,iMax,jMin,jMax
      INTEGER trIdentity
      LOGICAL trUseGMRedi, trUseKPP
      _RL KappaRTr(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
      INTEGER myThid

#ifdef ALLOW_GENERIC_ADVDIFF
C     !LOCAL VARIABLES:
C     == Local variables ==
C     i, j, k    :: Loop counters
C     iTr        :: passive tracer index
C     msgBuf     :: message buffer
      INTEGER i,j,k
      _RL KbryanLewis79
#ifdef ALLOW_BL79_LAT_VARY
      _RL KbryanLewisEQ
#endif
      CHARACTER*(MAX_LEN_MBUF) msgBuf
#ifdef ALLOW_PTRACERS
      INTEGER iTr
#endif
CEOP

      IF ( .NOT. trUseKPP ) THEN
       DO k = 1,Nr
        KbryanLewis79=diffKrBL79surf+(diffKrBL79deep-diffKrBL79surf)
     &       *(atan(-(rF(k)-diffKrBL79Ho)/diffKrBL79scl)/PI+0.5 _d 0)
#ifdef ALLOW_BL79_LAT_VARY
        KbryanLewisEQ=diffKrBLEQsurf+(diffKrBLEQdeep-diffKrBLEQsurf)
     &       *(atan(-(rF(k)-diffKrBLEQHo)/diffKrBLEQscl)/PI+0.5 _d 0)
#endif
        DO j = 1-Oly,sNy+Oly
         DO i = 1-Olx,sNx+Olx
          KappaRTr(i,j,k) =
     &         IVDConvCount(i,j,k,bi,bj)*ivdc_kappa
     &         + KbryanLewis79
#ifdef ALLOW_BL79_LAT_VARY
     &         + (KbryanLewisEQ-KbryanLewis79)*BL79LatArray(i,j,bi,bj)
#endif
         ENDDO
        ENDDO
       ENDDO
       IF ( trIdentity.EQ.GAD_TEMPERATURE ) THEN
        DO k = 1,Nr
         DO j = 1-Oly,sNy+Oly
          DO i = 1-Olx,sNx+Olx
           KappaRTr(i,j,k) = KappaRTr(i,j,k)
#if (defined ALLOW_3D_DIFFKR || \
     (defined (ALLOW_AUTODIFF_TAMC) && defined (ALLOW_DIFFKR_CONTROL)))
     &          + diffKr(i,j,k,bi,bj)
#else
     &          + diffKrNrT(k)
#endif
          ENDDO
         ENDDO
        ENDDO
       ELSEIF ( trIdentity.EQ.GAD_SALINITY) THEN
        DO k = 1,Nr
         DO j = 1-Oly, sNy+Oly
          DO i = 1-Olx, sNx+Olx
           KappaRTr(i,j,k) = KappaRTr(i,j,k)
#if (defined ALLOW_3D_DIFFKR || \
     (defined (ALLOW_AUTODIFF_TAMC) && defined (ALLOW_DIFFKR_CONTROL)))
     &          + diffKr(i,j,k,bi,bj)
#else
     &          + diffKrNrS(k)
#endif
          ENDDO
         ENDDO
        ENDDO
#ifdef ALLOW_PTRACERS
       ELSEIF ( trIdentity.GE.GAD_TR1
     &   .AND. trIdentity.LT.GAD_TR1+PTRACERS_numInUse) THEN

        iTr = trIdentity - GAD_TR1 + 1
        DO k = 1,Nr
         DO j = 1-Oly, sNy+Oly
          DO i = 1-Olx, sNx+Olx
           KappaRTr(i,j,k) = KappaRTr(i,j,k)
#if (defined ALLOW_3D_DIFFKR || \
     (defined (ALLOW_AUTODIFF_TAMC) && defined (ALLOW_DIFFKR_CONTROL)))
     &          + diffKr(i,j,k,bi,bj)
#else
     &          + PTRACERS_diffKrNr(k,iTr)
#endif
          ENDDO
         ENDDO
        ENDDO
#endif /* ALLOW_PTRACERS */
       ELSE
        WRITE(msgBuf,'(A,I4)')
     &       ' CALC_3D_DIFFUSIVITY: Invalid tracer Id: ',trIdentity
        CALL PRINT_ERROR(msgBuf, myThid)
        STOP 'ABNORMAL END: S/R CALC_3D_DIFFUSIVITY'
       ENDIF
      ENDIF

C--   Add physical pacakge contributions:

#ifdef ALLOW_KPP
      IF (trUseKPP) THEN
C--   Set vertical diffusivity contribution from KPP
       IF (trIdentity.EQ.GAD_TEMPERATURE) THEN
         CALL KPP_CALC_DIFF_T(
     I        bi,bj,iMin,iMax,jMin,jMax,0,Nr,
     O        KappaRTr,
     I        myThid)
       ELSE
         CALL KPP_CALC_DIFF_S(
     I        bi,bj,iMin,iMax,jMin,jMax,0,Nr,
     O        KappaRTr,
     I        myThid)
       ENDIF
#if (defined ALLOW_PTRACERS && ! (defined ALLOW_3D_DIFFKR || \
     (defined (ALLOW_AUTODIFF_TAMC) && defined (ALLOW_DIFFKR_CONTROL))))
       IF ( trIdentity.GE.GAD_TR1
     &   .AND. trIdentity.LT.GAD_TR1+PTRACERS_numInUse) THEN
        iTr = trIdentity - GAD_TR1 + 1
        DO k = 1,Nr
         DO j = 1-Oly, sNy+Oly
          DO i = 1-Olx, sNx+Olx
           KappaRTr(i,j,k) = KappaRTr(i,j,k)
     &          - diffKrNrS(k) + PTRACERS_diffKrNr(k,iTr)
          ENDDO
         ENDDO
        ENDDO
       ENDIF
#endif
      ENDIF
#endif /* ALLOW_KPP */

#ifdef ALLOW_GMREDI
      IF (trUseGMRedi) THEN
         CALL GMREDI_CALC_DIFF(
     I        bi,bj,iMin,iMax,jMin,jMax,0,Nr,
     U        KappaRTr,
     I        trIdentity,myThid)
      ENDIF
#endif

#ifdef ALLOW_PP81
      IF (usePP81) THEN
         CALL PP81_CALC_DIFF(
     I        bi,bj,iMin,iMax,jMin,jMax,0,Nr,
     U        KappaRTr,
     I        myThid)
      ENDIF
#endif

#ifdef ALLOW_MY82
      IF (useMY82) THEN
         CALL MY82_CALC_DIFF(
     I        bi,bj,iMin,iMax,jMin,jMax,0,Nr,
     U        KappaRTr,
     I        myThid)
      ENDIF
#endif

#ifdef ALLOW_GGL90
      IF (useGGL90) THEN
         CALL GGL90_CALC_DIFF(
     I        bi,bj,iMin,iMax,jMin,jMax,0,Nr,
     O        KappaRTr,
     I        myThid)
      ENDIF
#endif

C-    Apply mask to vertical diffusivity
C jmc: don't have the impression that masking is needed
C      but could be removed later if it's the case.
c     DO j = 1-Oly, sNy+Oly
c      DO i = 1-Olx, sNx+Olx
c       KappaRTr(i,j,k) = maskUp(i,j)*KappaRTr(i,j,k)
c      ENDDO
c     ENDDO

#endif /* ALLOW_GENERIC_ADVDIFF */

      RETURN
      END
1	dfer	1.14	C $Header: /u/gcmpack/MITgcm/model/src/calc_3d_diffusivity.F,v 1.13 2007/11/05 18:52:21 jmc Exp $
2	jmc	1.1	C $Name: $
3
4			#include "PACKAGES_CONFIG.h"
5			#include "CPP_OPTIONS.h"
6
7			CBOP
8			C !ROUTINE: CALC_3D_DIFFUSIVITY
9	jmc	1.11	C !INTERFACE:
10			SUBROUTINE CALC_3D_DIFFUSIVITY(
11	jmc	1.1	I bi,bj,iMin,iMax,jMin,jMax,
12			I trIdentity, trUseGMRedi, trUseKPP,
13			O KappaRTr,
14			I myThid)
15
16			C !DESCRIPTION: \bv
17			C ==========================================================
18			C \| SUBROUTINE CALC_3D_DIFFUSIVITY
19			C \| o Calculate net (3D) vertical diffusivity for 1 tracer
20			C ==========================================================
21			C \| Combines spatially varying diffusion coefficients from
22			C \| KPP and/or GM and/or convective stability test.
23			C ==========================================================
24			C \ev
25
26			C !USES:
27			IMPLICIT NONE
28			C == GLobal variables ==
29			#include "SIZE.h"
30			#include "EEPARAMS.h"
31			#include "PARAMS.h"
32			#include "DYNVARS.h"
33			#include "GRID.h"
34	jmc	1.3	#ifdef ALLOW_GENERIC_ADVDIFF
35	jmc	1.1	#include "GAD.h"
36	jmc	1.3	#endif
37	jmc	1.1	#ifdef ALLOW_PTRACERS
38			#include "PTRACERS_SIZE.h"
39	jmc	1.13	#include "PTRACERS_PARAMS.h"
40	jmc	1.1	#endif
41
42			C !INPUT/OUTPUT PARAMETERS:
43			C == Routine arguments ==
44			C bi, bj :: tile indices
45			C iMin,iMax :: Range of points for which calculation is performed.
46			C jMin,jMax :: Range of points for which calculation is performed.
47			C trIdentity :: tracer identifier
48			C trUseGMRedi:: this tracer use GM-Redi
49			C trUseKPP :: this tracer use KPP
50			C myThid :: Instance number for this innvocation of CALC_3D_DIFFUSIVITY
51			C KappaRTr :: Net diffusivity for this tracer (trIdentity)
52			INTEGER bi,bj,iMin,iMax,jMin,jMax
53			INTEGER trIdentity
54			LOGICAL trUseGMRedi, trUseKPP
55			_RL KappaRTr(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
56			INTEGER myThid
57
58	jmc	1.3	#ifdef ALLOW_GENERIC_ADVDIFF
59	jmc	1.1	C !LOCAL VARIABLES:
60			C == Local variables ==
61			C i, j, k :: Loop counters
62			C iTr :: passive tracer index
63			C msgBuf :: message buffer
64			INTEGER i,j,k
65	jmc	1.11	_RL KbryanLewis79
66			#ifdef ALLOW_BL79_LAT_VARY
67			_RL KbryanLewisEQ
68			#endif
69	jmc	1.1	CHARACTER*(MAX_LEN_MBUF) msgBuf
70	jmc	1.4	#ifdef ALLOW_PTRACERS
71			INTEGER iTr
72			#endif
73	jmc	1.1	CEOP
74
75	dimitri	1.10	IF ( .NOT. trUseKPP ) THEN
76	jmc	1.1	DO k = 1,Nr
77	jmc	1.8	KbryanLewis79=diffKrBL79surf+(diffKrBL79deep-diffKrBL79surf)
78	dimitri	1.10	& *(atan(-(rF(k)-diffKrBL79Ho)/diffKrBL79scl)/PI+0.5 _d 0)
79	dimitri	1.9	#ifdef ALLOW_BL79_LAT_VARY
80			KbryanLewisEQ=diffKrBLEQsurf+(diffKrBLEQdeep-diffKrBLEQsurf)
81	dimitri	1.10	& *(atan(-(rF(k)-diffKrBLEQHo)/diffKrBLEQscl)/PI+0.5 _d 0)
82	dimitri	1.9	#endif
83	jmc	1.1	DO j = 1-Oly,sNy+Oly
84			DO i = 1-Olx,sNx+Olx
85			KappaRTr(i,j,k) =
86			& IVDConvCount(i,j,k,bi,bj)*ivdc_kappa
87	dimitri	1.10	& + KbryanLewis79
88			#ifdef ALLOW_BL79_LAT_VARY
89			& + (KbryanLewisEQ-KbryanLewis79)*BL79LatArray(i,j,bi,bj)
90			#endif
91			ENDDO
92			ENDDO
93			ENDDO
94			IF ( trIdentity.EQ.GAD_TEMPERATURE ) THEN
95			DO k = 1,Nr
96			DO j = 1-Oly,sNy+Oly
97			DO i = 1-Olx,sNx+Olx
98	jmc	1.11	KappaRTr(i,j,k) = KappaRTr(i,j,k)
99	dimitri	1.6	#if (defined ALLOW_3D_DIFFKR \|\| \
100			(defined (ALLOW_AUTODIFF_TAMC) && defined (ALLOW_DIFFKR_CONTROL)))
101	dimitri	1.10	& + diffKr(i,j,k,bi,bj)
102	jmc	1.1	#else
103	dimitri	1.10	& + diffKrNrT(k)
104	dimitri	1.9	#endif
105	dimitri	1.10	ENDDO
106	jmc	1.1	ENDDO
107			ENDDO
108	dimitri	1.10	ELSEIF ( trIdentity.EQ.GAD_SALINITY) THEN
109			DO k = 1,Nr
110			DO j = 1-Oly, sNy+Oly
111			DO i = 1-Olx, sNx+Olx
112			KappaRTr(i,j,k) = KappaRTr(i,j,k)
113	dimitri	1.6	#if (defined ALLOW_3D_DIFFKR \|\| \
114			(defined (ALLOW_AUTODIFF_TAMC) && defined (ALLOW_DIFFKR_CONTROL)))
115	dimitri	1.10	& + diffKr(i,j,k,bi,bj)
116	jmc	1.1	#else
117	dimitri	1.10	& + diffKrNrS(k)
118	dimitri	1.9	#endif
119	dimitri	1.10	ENDDO
120	jmc	1.1	ENDDO
121			ENDDO
122			#ifdef ALLOW_PTRACERS
123	dimitri	1.10	ELSEIF ( trIdentity.GE.GAD_TR1
124	jmc	1.1	& .AND. trIdentity.LT.GAD_TR1+PTRACERS_numInUse) THEN
125
126	dimitri	1.10	iTr = trIdentity - GAD_TR1 + 1
127			DO k = 1,Nr
128			DO j = 1-Oly, sNy+Oly
129			DO i = 1-Olx, sNx+Olx
130			KappaRTr(i,j,k) = KappaRTr(i,j,k)
131	dimitri	1.6	#if (defined ALLOW_3D_DIFFKR \|\| \
132			(defined (ALLOW_AUTODIFF_TAMC) && defined (ALLOW_DIFFKR_CONTROL)))
133	dimitri	1.10	& + diffKr(i,j,k,bi,bj)
134	jmc	1.1	#else
135	dimitri	1.10	& + PTRACERS_diffKrNr(k,iTr)
136	dimitri	1.9	#endif
137	dimitri	1.10	ENDDO
138	jmc	1.1	ENDDO
139			ENDDO
140			#endif /* ALLOW_PTRACERS */
141	dimitri	1.10	ELSE
142	jmc	1.1	WRITE(msgBuf,'(A,I4)')
143	dimitri	1.10	& ' CALC_3D_DIFFUSIVITY: Invalid tracer Id: ',trIdentity
144	jmc	1.1	CALL PRINT_ERROR(msgBuf, myThid)
145			STOP 'ABNORMAL END: S/R CALC_3D_DIFFUSIVITY'
146	dimitri	1.10	ENDIF
147	jmc	1.1	ENDIF
148
149			C-- Add physical pacakge contributions:
150
151			#ifdef ALLOW_KPP
152			IF (trUseKPP) THEN
153	dimitri	1.12	C-- Set vertical diffusivity contribution from KPP
154	jmc	1.1	IF (trIdentity.EQ.GAD_TEMPERATURE) THEN
155			CALL KPP_CALC_DIFF_T(
156	jmc	1.2	I bi,bj,iMin,iMax,jMin,jMax,0,Nr,
157	dimitri	1.12	O KappaRTr,
158	jmc	1.1	I myThid)
159			ELSE
160			CALL KPP_CALC_DIFF_S(
161	jmc	1.2	I bi,bj,iMin,iMax,jMin,jMax,0,Nr,
162	dimitri	1.12	O KappaRTr,
163	jmc	1.1	I myThid)
164			ENDIF
165	dimitri	1.10	#if (defined ALLOW_PTRACERS && ! (defined ALLOW_3D_DIFFKR \|\| \
166			(defined (ALLOW_AUTODIFF_TAMC) && defined (ALLOW_DIFFKR_CONTROL))))
167			IF ( trIdentity.GE.GAD_TR1
168			& .AND. trIdentity.LT.GAD_TR1+PTRACERS_numInUse) THEN
169			iTr = trIdentity - GAD_TR1 + 1
170			DO k = 1,Nr
171			DO j = 1-Oly, sNy+Oly
172			DO i = 1-Olx, sNx+Olx
173			KappaRTr(i,j,k) = KappaRTr(i,j,k)
174			& - diffKrNrS(k) + PTRACERS_diffKrNr(k,iTr)
175			ENDDO
176			ENDDO
177			ENDDO
178			ENDIF
179			#endif
180			ENDIF
181			#endif /* ALLOW_KPP */
182
183			#ifdef ALLOW_GMREDI
184	jmc	1.11	IF (trUseGMRedi) THEN
185	dimitri	1.10	CALL GMREDI_CALC_DIFF(
186			I bi,bj,iMin,iMax,jMin,jMax,0,Nr,
187			U KappaRTr,
188	dfer	1.14	I trIdentity,myThid)
189	jmc	1.1	ENDIF
190			#endif
191
192			#ifdef ALLOW_PP81
193			IF (usePP81) THEN
194			CALL PP81_CALC_DIFF(
195	jmc	1.2	I bi,bj,iMin,iMax,jMin,jMax,0,Nr,
196	jmc	1.1	U KappaRTr,
197			I myThid)
198			ENDIF
199			#endif
200
201			#ifdef ALLOW_MY82
202			IF (useMY82) THEN
203			CALL MY82_CALC_DIFF(
204	jmc	1.2	I bi,bj,iMin,iMax,jMin,jMax,0,Nr,
205	jmc	1.1	U KappaRTr,
206			I myThid)
207			ENDIF
208			#endif
209	jmc	1.11
210	jmc	1.1	#ifdef ALLOW_GGL90
211			IF (useGGL90) THEN
212			CALL GGL90_CALC_DIFF(
213	jmc	1.2	I bi,bj,iMin,iMax,jMin,jMax,0,Nr,
214	jmc	1.1	O KappaRTr,
215			I myThid)
216			ENDIF
217			#endif
218	jmc	1.11
219			C- Apply mask to vertical diffusivity
220			C jmc: don't have the impression that masking is needed
221	jmc	1.1	C but could be removed later if it's the case.
222			c DO j = 1-Oly, sNy+Oly
223			c DO i = 1-Olx, sNx+Olx
224			c KappaRTr(i,j,k) = maskUp(i,j)*KappaRTr(i,j,k)
225			c ENDDO
226			c ENDDO
227
228	jmc	1.3	#endif /* ALLOW_GENERIC_ADVDIFF */
229
230	jmc	1.1	RETURN
231			END