pkg/gmredi/gmredi_rtransport.F

C $Header: /u/gcmpack/MITgcm/pkg/gmredi/gmredi_rtransport.F,v 1.9 2002/11/14 22:43:49 heimbach Exp $
C $Name:  $

#include "GMREDI_OPTIONS.h"

      subroutine GMREDI_RTRANSPORT(
     I     iMin,iMax,jMin,jMax,bi,bj,K,
     I     Tracer,tracerIdentity,
     U     df,
     I     myThid)
C     /==========================================================\
C     | o SUBROUTINE GMREDI_RTRANSPORT                           |
C     |   Add vertical transport terms from GM/Redi              |
C     |   parameterization.                                      |
C     |==========================================================|
C     \==========================================================/
      IMPLICIT NONE

C     == GLobal variables ==
#include "SIZE.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "GRID.h"
#include "GMREDI.h"

C     == Routine arguments ==
C     iMin,iMax,jMin,  - Range of points for which calculation
C     jMax,bi,bj,k       results will be set.
C     xA               - Area of X face
C     Tracer           - 3D Tracer field
C     df               - Diffusive flux component work array.
      INTEGER iMin,iMax,jMin,jMax,bi,bj,k
      _RL Tracer(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
      INTEGER tracerIdentity
      _RL df    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
      INTEGER myThid

#ifdef ALLOW_GMREDI

C     == Local variables ==
C     I, J - Loop counters
      INTEGER I, J
      _RL dTdx,dTdy
      _RL rTrans

C     Surface flux is zero
      IF (useGMRedi .AND. K.GT.1) THEN

      DO j=jMin,jMax
       DO i=iMin,iMax

C-      Horizontal gradients interpolated to W points
        dTdx = 0.5d0*(
     &   +0.5d0*(_maskW(i+1,j,k,bi,bj)
     &         *_recip_dxC(i+1,j,bi,bj)*
     &           (Tracer(i+1,j,k,bi,bj)-Tracer(i,j,k,bi,bj))
     &        +_maskW(i,j,k,bi,bj)
     &         *_recip_dxC(i,j,bi,bj)*
     &           (Tracer(i,j,k,bi,bj)-Tracer(i-1,j,k,bi,bj)))
     &   +0.5d0*(_maskW(i+1,j,k-1,bi,bj)
     &         *_recip_dxC(i+1,j,bi,bj)*
     &           (Tracer(i+1,j,k-1,bi,bj)-Tracer(i,j,k-1,bi,bj))
     &        +_maskW(i,j,k-1,bi,bj)
     &         *_recip_dxC(i,j,bi,bj)*
     &           (Tracer(i,j,k-1,bi,bj)-Tracer(i-1,j,k-1,bi,bj)))
     &       )

        dTdy = 0.5d0*(
     &   +0.5d0*(_maskS(i,j,k,bi,bj)
     &         *_recip_dyC(i,j,bi,bj)*
     &           (Tracer(i,j,k,bi,bj)-Tracer(i,j-1,k,bi,bj))
     &        +_maskS(i,j+1,k,bi,bj)
     &         *_recip_dyC(i,j+1,bi,bj)*
     &           (Tracer(i,j+1,k,bi,bj)-Tracer(i,j,k,bi,bj)))
     &   +0.5d0*(_maskS(i,j,k-1,bi,bj)
     &         *_recip_dyC(i,j,bi,bj)*
     &           (Tracer(i,j,k-1,bi,bj)-Tracer(i,j-1,k-1,bi,bj))
     &        +_maskS(i,j+1,k-1,bi,bj)
     &         *_recip_dyC(i,j+1,bi,bj)*
     &           (Tracer(i,j+1,k-1,bi,bj)-Tracer(i,j,k-1,bi,bj)))
     &       )

C-      Off-diagonal components of vertical flux
        df(i,j) = df(i,j)
     &      - _rA(i,j,bi,bj)
     &        *( Kwx(i,j,k,bi,bj)*dTdx +Kwy(i,j,k,bi,bj)*dTdy )
cph     &        *GM_adjointRescale

       ENDDO
      ENDDO

#ifdef GM_BOLUS_ADVEC
      IF (GM_AdvForm .AND. GM_AdvSeparate) THEN
       DO j=jMin,jMax
        DO i=iMin,iMax
         rTrans = 
     &      dyG(i+1,j,bi,bj)*GM_PsiX(i+1,j,k,bi,bj)
     &     -dyG( i ,j,bi,bj)*GM_PsiX( i ,j,k,bi,bj)
     &     +dxG(i,j+1,bi,bj)*GM_PsiY(i,j+1,k,bi,bj)
     &     -dxG(i, j ,bi,bj)*GM_PsiY(i, j ,k,bi,bj)
         df(i,j) = df(i,j)
     &    +rTrans*0.5d0
     &           *(Tracer(i,j,k,bi,bj)+Tracer(i,j,k-1,bi,bj))
        ENDDO
       ENDDO
      ENDIF
#endif /* GM_BOLUS_ADVEC */  

c     IF (.NOT.implicitDiffusion) THEN
c
c This vertical diffusion term is currently implemented
c by adding the VisbeckK*Kwz diffusivity to KappaRT/S
c See calc_diffusivity.F and calc_gt.F (calc_gs.F)
c
c      DO j=jMin,jMax
c       DO i=iMin,iMax
c        df(i,j) = df(i,j) - _rA(i,j,bi,bj)
c    &    *maskUp(i,j)*VisbeckK(i,j,bi,bj)*Kwz(i,j,k,bi,bj)
c    &    *recip_drC(k)*rkfac
c    &    *(Tracer(i,j,k-1,bi,bj)-Tracer(i,j,k,bi,bj))
c       ENDDO
c      ENDDO
c     ENDIF

      ENDIF
#endif /* ALLOW_GMREDI */

      RETURN
      END
1	C $Header: /u/gcmpack/MITgcm/pkg/gmredi/gmredi_rtransport.F,v 1.9 2002/11/14 22:43:49 heimbach Exp $
2	C $Name: $
3
4	#include "GMREDI_OPTIONS.h"
5
6	subroutine GMREDI_RTRANSPORT(
7	I iMin,iMax,jMin,jMax,bi,bj,K,
8	I Tracer,tracerIdentity,
9	U df,
10	I myThid)
11	C /==========================================================\
12	C \| o SUBROUTINE GMREDI_RTRANSPORT \|
13	C \| Add vertical transport terms from GM/Redi \|
14	C \| parameterization. \|
15	C \|==========================================================\|
16	C \==========================================================/
17	IMPLICIT NONE
18
19	C == GLobal variables ==
20	#include "SIZE.h"
21	#include "EEPARAMS.h"
22	#include "PARAMS.h"
23	#include "GRID.h"
24	#include "GMREDI.h"
25
26	C == Routine arguments ==
27	C iMin,iMax,jMin, - Range of points for which calculation
28	C jMax,bi,bj,k results will be set.
29	C xA - Area of X face
30	C Tracer - 3D Tracer field
31	C df - Diffusive flux component work array.
32	INTEGER iMin,iMax,jMin,jMax,bi,bj,k
33	_RL Tracer(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
34	INTEGER tracerIdentity
35	_RL df (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
36	INTEGER myThid
37
38	#ifdef ALLOW_GMREDI
39
40	C == Local variables ==
41	C I, J - Loop counters
42	INTEGER I, J
43	_RL dTdx,dTdy
44	_RL rTrans
45
46	C Surface flux is zero
47	IF (useGMRedi .AND. K.GT.1) THEN
48
49	DO j=jMin,jMax
50	DO i=iMin,iMax
51
52	C- Horizontal gradients interpolated to W points
53	dTdx = 0.5d0*(
54	& +0.5d0*(_maskW(i+1,j,k,bi,bj)
55	& _recip_dxC(i+1,j,bi,bj)
56	& (Tracer(i+1,j,k,bi,bj)-Tracer(i,j,k,bi,bj))
57	& +_maskW(i,j,k,bi,bj)
58	& _recip_dxC(i,j,bi,bj)
59	& (Tracer(i,j,k,bi,bj)-Tracer(i-1,j,k,bi,bj)))
60	& +0.5d0*(_maskW(i+1,j,k-1,bi,bj)
61	& _recip_dxC(i+1,j,bi,bj)
62	& (Tracer(i+1,j,k-1,bi,bj)-Tracer(i,j,k-1,bi,bj))
63	& +_maskW(i,j,k-1,bi,bj)
64	& _recip_dxC(i,j,bi,bj)
65	& (Tracer(i,j,k-1,bi,bj)-Tracer(i-1,j,k-1,bi,bj)))
66	& )
67
68	dTdy = 0.5d0*(
69	& +0.5d0*(_maskS(i,j,k,bi,bj)
70	& _recip_dyC(i,j,bi,bj)
71	& (Tracer(i,j,k,bi,bj)-Tracer(i,j-1,k,bi,bj))
72	& +_maskS(i,j+1,k,bi,bj)
73	& _recip_dyC(i,j+1,bi,bj)
74	& (Tracer(i,j+1,k,bi,bj)-Tracer(i,j,k,bi,bj)))
75	& +0.5d0*(_maskS(i,j,k-1,bi,bj)
76	& _recip_dyC(i,j,bi,bj)
77	& (Tracer(i,j,k-1,bi,bj)-Tracer(i,j-1,k-1,bi,bj))
78	& +_maskS(i,j+1,k-1,bi,bj)
79	& _recip_dyC(i,j+1,bi,bj)
80	& (Tracer(i,j+1,k-1,bi,bj)-Tracer(i,j,k-1,bi,bj)))
81	& )
82
83	C- Off-diagonal components of vertical flux
84	df(i,j) = df(i,j)
85	& - _rA(i,j,bi,bj)
86	& ( Kwx(i,j,k,bi,bj)dTdx +Kwy(i,j,k,bi,bj)*dTdy )
87	cph & *GM_adjointRescale
88
89	ENDDO
90	ENDDO
91
92	#ifdef GM_BOLUS_ADVEC
93	IF (GM_AdvForm .AND. GM_AdvSeparate) THEN
94	DO j=jMin,jMax
95	DO i=iMin,iMax
96	rTrans =
97	& dyG(i+1,j,bi,bj)*GM_PsiX(i+1,j,k,bi,bj)
98	& -dyG( i ,j,bi,bj)*GM_PsiX( i ,j,k,bi,bj)
99	& +dxG(i,j+1,bi,bj)*GM_PsiY(i,j+1,k,bi,bj)
100	& -dxG(i, j ,bi,bj)*GM_PsiY(i, j ,k,bi,bj)
101	df(i,j) = df(i,j)
102	& +rTrans*0.5d0
103	& *(Tracer(i,j,k,bi,bj)+Tracer(i,j,k-1,bi,bj))
104	ENDDO
105	ENDDO
106	ENDIF
107	#endif /* GM_BOLUS_ADVEC */
108
109	c IF (.NOT.implicitDiffusion) THEN
110	c
111	c This vertical diffusion term is currently implemented
112	c by adding the VisbeckK*Kwz diffusivity to KappaRT/S
113	c See calc_diffusivity.F and calc_gt.F (calc_gs.F)
114	c
115	c DO j=jMin,jMax
116	c DO i=iMin,iMax
117	c df(i,j) = df(i,j) - _rA(i,j,bi,bj)
118	c & maskUp(i,j)VisbeckK(i,j,bi,bj)*Kwz(i,j,k,bi,bj)
119	c & recip_drC(k)rkfac
120	c & *(Tracer(i,j,k-1,bi,bj)-Tracer(i,j,k,bi,bj))
121	c ENDDO
122	c ENDDO
123	c ENDIF
124
125	ENDIF
126	#endif /* ALLOW_GMREDI */
127
128	RETURN
129	END