model/src/external_forcing.F

C $Header: /u/gcmpack/MITgcm/model/src/external_forcing.F,v 1.38 2006/06/07 01:55:12 heimbach Exp $
C $Name:  $

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

CBOP
C     !ROUTINE: EXTERNAL_FORCING_U
C     !INTERFACE:
      SUBROUTINE EXTERNAL_FORCING_U(
     I           iMin,iMax, jMin,jMax, bi,bj, kLev,
     I           myTime, myThid )
C     !DESCRIPTION: \bv
C     *==========================================================*
C     | S/R EXTERNAL_FORCING_U
C     | o Contains problem specific forcing for zonal velocity.
C     *==========================================================*
C     | Adds terms to gU for forcing by external sources
C     | e.g. wind stress, bottom friction etc ...
C     *==========================================================*
C     \ev

C     !USES:
      IMPLICIT NONE
C     == Global data ==
#include "SIZE.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "GRID.h"
#include "DYNVARS.h"
#include "FFIELDS.h"

C     !INPUT/OUTPUT PARAMETERS:
C     == Routine arguments ==
C     iMin,iMax :: Working range of x-index for applying forcing.
C     jMin,jMax :: Working range of y-index for applying forcing.
C     bi,bj     :: Current tile indices
C     kLev      :: Current vertical level index
C     myTime    :: Current time in simulation
C     myThid    :: Thread Id number
      INTEGER iMin, iMax, jMin, jMax, kLev, bi, bj
      _RL myTime
      INTEGER myThid

C     !LOCAL VARIABLES:
C     == Local variables ==
C     i,j       :: Loop counters
C     kSurface  :: index of surface layer
      INTEGER i, j
      INTEGER kSurface
CEOP

      IF ( fluidIsAir ) THEN
       kSurface = 0
      ELSEIF ( usingPCoords ) THEN
       kSurface = Nr
      ELSE
       kSurface = 1
      ENDIF

C--   Forcing term
#ifdef ALLOW_AIM
      IF ( useAIM ) CALL AIM_TENDENCY_APPLY_U(
     &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
     &                      myTime, myThid )
#endif /* ALLOW_AIM */

#ifdef ALLOW_FIZHI
      IF ( useFIZHI ) CALL FIZHI_TENDENCY_APPLY_U(
     &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
     &                      myTime, myThid )
#endif /* ALLOW_FIZHI */

#ifdef ALLOW_MYPACKAGE
      IF ( useMYPACKAGE ) CALL MYPACKAGE_TENDENCY_APPLY_U(
     &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
     &                      myTime, myThid )
#endif /* ALLOW_MYPACKAGE */

C     Add windstress momentum impulse into the top-layer
      IF ( kLev .EQ. kSurface ) THEN
c      DO j=1,sNy
C-jmc: Without CD-scheme, this is OK ; but with CD-scheme, needs to cover [0:sNy+1]
       DO j=0,sNy+1
        DO i=1,sNx+1
         gU(i,j,kLev,bi,bj) = gU(i,j,kLev,bi,bj)
     &   +foFacMom*surfaceForcingU(i,j,bi,bj)
     &   *recip_drF(kLev)*_recip_hFacW(i,j,kLev,bi,bj)
        ENDDO
       ENDDO
      ENDIF

#if (defined (ALLOW_TAU_EDDY))
       CALL TAUEDDY_EXTERNAL_FORCING_U(
     I           iMin,iMax, jMin,jMax, bi,bj, kLev,
     I           myTime, myThid )
#endif

#ifdef ALLOW_OBCS
      IF (useOBCS) THEN
       CALL OBCS_SPONGE_U(
     I           iMin,iMax, jMin,jMax, bi,bj, kLev,
     I           myTime, myThid )
      ENDIF
#endif

      RETURN
      END

C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
CBOP
C     !ROUTINE: EXTERNAL_FORCING_V
C     !INTERFACE:
      SUBROUTINE EXTERNAL_FORCING_V(
     I           iMin,iMax, jMin,jMax, bi,bj, kLev,
     I           myTime, myThid )
C     !DESCRIPTION: \bv
C     *==========================================================*
C     | S/R EXTERNAL_FORCING_V
C     | o Contains problem specific forcing for merid velocity.
C     *==========================================================*
C     | Adds terms to gV for forcing by external sources
C     | e.g. wind stress, bottom friction etc ...
C     *==========================================================*
C     \ev

C     !USES:
      IMPLICIT NONE
C     == Global data ==
#include "SIZE.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "GRID.h"
#include "DYNVARS.h"
#include "FFIELDS.h"

C     !INPUT/OUTPUT PARAMETERS:
C     == Routine arguments ==
C     iMin,iMax :: Working range of x-index for applying forcing.
C     jMin,jMax :: Working range of y-index for applying forcing.
C     bi,bj     :: Current tile indices
C     kLev      :: Current vertical level index
C     myTime    :: Current time in simulation
C     myThid    :: Thread Id number
      INTEGER iMin, iMax, jMin, jMax, kLev, bi, bj
      _RL myTime
      INTEGER myThid

C     !LOCAL VARIABLES:
C     == Local variables ==
C     i,j       :: Loop counters
C     kSurface  :: index of surface layer
      INTEGER i, j
      INTEGER kSurface
CEOP

      IF ( fluidIsAir ) THEN
       kSurface = 0
      ELSEIF ( usingPCoords ) THEN
       kSurface = Nr
      ELSE
       kSurface = 1
      ENDIF

C--   Forcing term
#ifdef ALLOW_AIM
      IF ( useAIM ) CALL AIM_TENDENCY_APPLY_V(
     &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
     &                      myTime, myThid )
#endif /* ALLOW_AIM */

#ifdef ALLOW_FIZHI
      IF ( useFIZHI ) CALL FIZHI_TENDENCY_APPLY_V(
     &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
     &                      myTime, myThid )
#endif /* ALLOW_FIZHI */

#ifdef ALLOW_MYPACKAGE
      IF ( useMYPACKAGE ) CALL MYPACKAGE_TENDENCY_APPLY_V(
     &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
     &                      myTime, myThid )
#endif /* ALLOW_MYPACKAGE */

C     Add windstress momentum impulse into the top-layer
      IF ( kLev .EQ. kSurface ) THEN
       DO j=1,sNy+1
c       DO i=1,sNx
C-jmc: Without CD-scheme, this is OK ; but with CD-scheme, needs to cover [0:sNx+1]
        DO i=0,sNx+1
         gV(i,j,kLev,bi,bj) = gV(i,j,kLev,bi,bj)
     &   +foFacMom*surfaceForcingV(i,j,bi,bj)
     &   *recip_drF(kLev)*_recip_hFacS(i,j,kLev,bi,bj)
        ENDDO
       ENDDO
      ENDIF

#if (defined (ALLOW_TAU_EDDY))
       CALL TAUEDDY_EXTERNAL_FORCING_V(
     I           iMin,iMax, jMin,jMax, bi,bj, kLev,
     I           myTime, myThid )
#endif

#ifdef ALLOW_OBCS
      IF (useOBCS) THEN
       CALL OBCS_SPONGE_V(
     I           iMin,iMax, jMin,jMax, bi,bj, kLev,
     I           myTime, myThid )
      ENDIF
#endif

      RETURN
      END

C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
CBOP
C     !ROUTINE: EXTERNAL_FORCING_T
C     !INTERFACE:
      SUBROUTINE EXTERNAL_FORCING_T(
     I           iMin,iMax, jMin,jMax, bi,bj, kLev,
     I           myTime, myThid )
C     !DESCRIPTION: \bv
C     *==========================================================*
C     | S/R EXTERNAL_FORCING_T
C     | o Contains problem specific forcing for temperature.
C     *==========================================================*
C     | Adds terms to gT for forcing by external sources
C     | e.g. heat flux, climatalogical relaxation, etc ...
C     *==========================================================*
C     \ev

C     !USES:
      IMPLICIT NONE
C     == Global data ==
#include "SIZE.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "GRID.h"
#include "DYNVARS.h"
#include "FFIELDS.h"

C     !INPUT/OUTPUT PARAMETERS:
C     == Routine arguments ==
C     iMin,iMax :: Working range of x-index for applying forcing.
C     jMin,jMax :: Working range of y-index for applying forcing.
C     bi,bj     :: Current tile indices
C     kLev      :: Current vertical level index
C     myTime    :: Current time in simulation
C     myThid    :: Thread Id number
      INTEGER iMin, iMax, jMin, jMax, kLev, bi, bj
      _RL myTime
      INTEGER myThid

C     !LOCAL VARIABLES:
C     == Local variables ==
C     i,j       :: Loop counters
C     kSurface  :: index of surface layer
      INTEGER i, j
      INTEGER kSurface
CEOP
#ifdef SHORTWAVE_HEATING
      integer two
      _RL minusone
      parameter (two=2,minusone=-1.)
      _RL swfracb(two)
      INTEGER kp1
#endif

      IF ( fluidIsAir ) THEN
       kSurface = 0
      ELSEIF ( usingPCoords ) THEN
       kSurface = Nr
      ELSE
       kSurface = 1
      ENDIF

C--   Forcing term
#ifdef ALLOW_AIM
      IF ( useAIM ) CALL AIM_TENDENCY_APPLY_T(
     &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
     &                      myTime, myThid )
#endif /* ALLOW_AIM */

#ifdef ALLOW_FIZHI
      IF ( useFIZHI ) CALL FIZHI_TENDENCY_APPLY_T(
     &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
     &                      myTime, myThid )
#endif /* ALLOW_FIZHI */

#ifdef ALLOW_MYPACKAGE
      IF ( useMYPACKAGE ) CALL MYPACKAGE_TENDENCY_APPLY_T(
     &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
     &                      myTime, myThid )
#endif /* ALLOW_MYPACKAGE */

C     Add heat in top-layer
      IF ( kLev .EQ. kSurface ) THEN
       DO j=1,sNy
        DO i=1,sNx
         gT(i,j,kLev,bi,bj)=gT(i,j,kLev,bi,bj)
     &     +surfaceForcingT(i,j,bi,bj)
     &     *recip_drF(kLev)*_recip_hFacC(i,j,kLev,bi,bj)
        ENDDO
       ENDDO
      ENDIF

#ifdef ALLOW_SHELFICE
      IF ( useShelfIce )
     &     CALL SHELFICE_FORCING_T(
     I     iMin,iMax, jMin,jMax, bi,bj, kLev,
     I     myTime, myThid )
#endif /* ALLOW_SHELFICE */

#ifdef SHORTWAVE_HEATING
C Penetrating SW radiation
c     IF ( usePenetratingSW ) THEN
       swfracb(1)=abs(rF(klev))
       swfracb(2)=abs(rF(klev+1))
       CALL SWFRAC(
     I     two,minusone,
     I     myTime,myThid,
     U     swfracb)
       kp1 = klev+1
       IF (klev.EQ.Nr) THEN
        kp1 = klev
        swfracb(2)=0. _d 0
       ENDIF
       DO j=1,sNy
        DO i=1,sNx
         gT(i,j,klev,bi,bj) = gT(i,j,klev,bi,bj)
     &   -Qsw(i,j,bi,bj)*(swfracb(1)*maskC(i,j,klev,bi,bj)
     &                   -swfracb(2)*maskC(i,j,kp1, bi,bj))
     &    *recip_Cp*recip_rhoConst
     &    *recip_drF(klev)*_recip_hFacC(i,j,kLev,bi,bj)
        ENDDO
       ENDDO
c     ENDIF
#endif

#ifdef ALLOW_RBCS
       if (useRBCS) then
          call RBCS_ADD_TENDENCY(bi,bj,klev, 1,
     &                            myTime, myThid )
       endif
#endif

#ifdef ALLOW_OBCS
      IF (useOBCS) THEN
       CALL OBCS_SPONGE_T(
     I           iMin,iMax, jMin,jMax, bi,bj, kLev,
     I           myTime, myThid )
      ENDIF
#endif

      RETURN
      END

C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
CBOP
C     !ROUTINE: EXTERNAL_FORCING_S
C     !INTERFACE:
      SUBROUTINE EXTERNAL_FORCING_S(
     I           iMin,iMax, jMin,jMax, bi,bj, kLev,
     I           myTime, myThid )

C     !DESCRIPTION: \bv
C     *==========================================================*
C     | S/R EXTERNAL_FORCING_S
C     | o Contains problem specific forcing for merid velocity.
C     *==========================================================*
C     | Adds terms to gS for forcing by external sources
C     | e.g. fresh-water flux, climatalogical relaxation, etc ...
C     *==========================================================*
C     \ev

C     !USES:
      IMPLICIT NONE
C     == Global data ==
#include "SIZE.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "GRID.h"
#include "DYNVARS.h"
#include "FFIELDS.h"

C     !INPUT/OUTPUT PARAMETERS:
C     == Routine arguments ==
C     iMin,iMax :: Working range of x-index for applying forcing.
C     jMin,jMax :: Working range of y-index for applying forcing.
C     bi,bj     :: Current tile indices
C     kLev      :: Current vertical level index
C     myTime    :: Current time in simulation
C     myThid    :: Thread Id number
      INTEGER iMin, iMax, jMin, jMax, kLev, bi, bj
      _RL myTime
      INTEGER myThid

C     !LOCAL VARIABLES:
C     == Local variables ==
C     i,j       :: Loop counters
C     kSurface  :: index of surface layer
      INTEGER i, j
      INTEGER kSurface
CEOP

      IF ( fluidIsAir ) THEN
       kSurface = 0
      ELSEIF ( usingPCoords ) THEN
       kSurface = Nr
      ELSE
       kSurface = 1
      ENDIF

C--   Forcing term
#ifdef ALLOW_AIM
      IF ( useAIM ) CALL AIM_TENDENCY_APPLY_S(
     &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
     &                      myTime, myThid )
#endif /* ALLOW_AIM */

#ifdef ALLOW_FIZHI
      IF ( useFIZHI ) CALL FIZHI_TENDENCY_APPLY_S(
     &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
     &                      myTime, myThid )
#endif /* ALLOW_FIZHI */

#ifdef ALLOW_MYPACKAGE
      IF ( useMYPACKAGE ) CALL MYPACKAGE_TENDENCY_APPLY_S(
     &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
     &                      myTime, myThid )
#endif /* ALLOW_MYPACKAGE */

C     Add fresh-water in top-layer
      IF ( kLev .EQ. kSurface ) THEN
       DO j=1,sNy
        DO i=1,sNx
         gS(i,j,kLev,bi,bj)=gS(i,j,kLev,bi,bj)
     &     +surfaceForcingS(i,j,bi,bj)
     &     *recip_drF(kLev)*_recip_hFacC(i,j,kLev,bi,bj)
        ENDDO
       ENDDO
      ENDIF

#ifdef ALLOW_SHELFICE
      IF ( useShelfIce )
     &     CALL SHELFICE_FORCING_S(
     I     iMin,iMax, jMin,jMax, bi,bj, kLev,
     I     myTime, myThid )
#endif /* ALLOW_SHELFICE */

#ifdef ALLOW_RBCS
       if (useRBCS) then
          call RBCS_ADD_TENDENCY(bi,bj,klev, 2,
     &                            myTime, myThid )
       endif
#endif

#ifdef ALLOW_OBCS
      IF (useOBCS) THEN
       CALL OBCS_SPONGE_S(
     I           iMin,iMax, jMin,jMax, bi,bj, kLev,
     I           myTime, myThid )
      ENDIF
#endif

      RETURN
      END
1	heimbach	1.39	C $Header: /u/gcmpack/MITgcm/model/src/external_forcing.F,v 1.38 2006/06/07 01:55:12 heimbach Exp $
2	adcroft	1.12	C $Name: $
3	cnh	1.1
4	edhill	1.20	#include "PACKAGES_CONFIG.h"
5	cnh	1.1	#include "CPP_OPTIONS.h"
6
7	cnh	1.13	CBOP
8			C !ROUTINE: EXTERNAL_FORCING_U
9			C !INTERFACE:
10	cnh	1.1	SUBROUTINE EXTERNAL_FORCING_U(
11	jmc	1.31	I iMin,iMax, jMin,jMax, bi,bj, kLev,
12			I myTime, myThid )
13	cnh	1.13	C !DESCRIPTION: \bv
14			C ==========================================================
15	jmc	1.31	C \| S/R EXTERNAL_FORCING_U
16			C \| o Contains problem specific forcing for zonal velocity.
17	cnh	1.13	C ==========================================================
18	jmc	1.31	C \| Adds terms to gU for forcing by external sources
19			C \| e.g. wind stress, bottom friction etc ...
20	cnh	1.13	C ==========================================================
21			C \ev
22
23			C !USES:
24	cnh	1.2	IMPLICIT NONE
25	cnh	1.1	C == Global data ==
26			#include "SIZE.h"
27			#include "EEPARAMS.h"
28			#include "PARAMS.h"
29			#include "GRID.h"
30			#include "DYNVARS.h"
31	cnh	1.2	#include "FFIELDS.h"
32	cnh	1.13
33			C !INPUT/OUTPUT PARAMETERS:
34	cnh	1.1	C == Routine arguments ==
35	jmc	1.31	C iMin,iMax :: Working range of x-index for applying forcing.
36			C jMin,jMax :: Working range of y-index for applying forcing.
37			C bi,bj :: Current tile indices
38			C kLev :: Current vertical level index
39			C myTime :: Current time in simulation
40			C myThid :: Thread Id number
41	cnh	1.1	INTEGER iMin, iMax, jMin, jMax, kLev, bi, bj
42	jmc	1.31	_RL myTime
43	adcroft	1.4	INTEGER myThid
44	cnh	1.1
45	cnh	1.13	C !LOCAL VARIABLES:
46	cnh	1.2	C == Local variables ==
47	jmc	1.31	C i,j :: Loop counters
48			C kSurface :: index of surface layer
49			INTEGER i, j
50	mlosch	1.17	INTEGER kSurface
51	cnh	1.13	CEOP
52	cnh	1.2
53	jmc	1.28	IF ( fluidIsAir ) THEN
54	jmc	1.21	kSurface = 0
55	jmc	1.28	ELSEIF ( usingPCoords ) THEN
56	mlosch	1.17	kSurface = Nr
57	jmc	1.28	ELSE
58	mlosch	1.17	kSurface = 1
59	jmc	1.28	ENDIF
60	mlosch	1.17
61	cnh	1.2	C-- Forcing term
62	jmc	1.21	#ifdef ALLOW_AIM
63			IF ( useAIM ) CALL AIM_TENDENCY_APPLY_U(
64			& iMin,iMax, jMin,jMax, bi,bj, kLev,
65	jmc	1.31	& myTime, myThid )
66	jmc	1.21	#endif /* ALLOW_AIM */
67	jmc	1.31
68	molod	1.23	#ifdef ALLOW_FIZHI
69			IF ( useFIZHI ) CALL FIZHI_TENDENCY_APPLY_U(
70			& iMin,iMax, jMin,jMax, bi,bj, kLev,
71	jmc	1.31	& myTime, myThid )
72	molod	1.23	#endif /* ALLOW_FIZHI */
73	jmc	1.21
74	heimbach	1.39	#ifdef ALLOW_MYPACKAGE
75			IF ( useMYPACKAGE ) CALL MYPACKAGE_TENDENCY_APPLY_U(
76			& iMin,iMax, jMin,jMax, bi,bj, kLev,
77			& myTime, myThid )
78			#endif /* ALLOW_MYPACKAGE */
79
80	cnh	1.2	C Add windstress momentum impulse into the top-layer
81	mlosch	1.17	IF ( kLev .EQ. kSurface ) THEN
82	jmc	1.32	c DO j=1,sNy
83			C-jmc: Without CD-scheme, this is OK ; but with CD-scheme, needs to cover [0:sNy+1]
84			DO j=0,sNy+1
85	jmc	1.31	DO i=1,sNx+1
86	cnh	1.2	gU(i,j,kLev,bi,bj) = gU(i,j,kLev,bi,bj)
87	jmc	1.26	& +foFacMom*surfaceForcingU(i,j,bi,bj)
88	heimbach	1.38	& recip_drF(kLev)_recip_hFacW(i,j,kLev,bi,bj)
89	cnh	1.2	ENDDO
90			ENDDO
91			ENDIF
92
93	heimbach	1.30	#if (defined (ALLOW_TAU_EDDY))
94	heimbach	1.29	CALL TAUEDDY_EXTERNAL_FORCING_U(
95	jmc	1.31	I iMin,iMax, jMin,jMax, bi,bj, kLev,
96			I myTime, myThid )
97	heimbach	1.29	#endif
98
99	jmc	1.31	#ifdef ALLOW_OBCS
100	heimbach	1.16	IF (useOBCS) THEN
101			CALL OBCS_SPONGE_U(
102	jmc	1.31	I iMin,iMax, jMin,jMax, bi,bj, kLev,
103			I myTime, myThid )
104	heimbach	1.16	ENDIF
105	heimbach	1.14	#endif
106
107	cnh	1.1	RETURN
108			END
109	jmc	1.31
110			C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-\|--+----\|
111	cnh	1.13	CBOP
112			C !ROUTINE: EXTERNAL_FORCING_V
113			C !INTERFACE:
114	cnh	1.1	SUBROUTINE EXTERNAL_FORCING_V(
115	jmc	1.31	I iMin,iMax, jMin,jMax, bi,bj, kLev,
116			I myTime, myThid )
117	cnh	1.13	C !DESCRIPTION: \bv
118			C ==========================================================
119	jmc	1.31	C \| S/R EXTERNAL_FORCING_V
120			C \| o Contains problem specific forcing for merid velocity.
121	cnh	1.13	C ==========================================================
122	jmc	1.31	C \| Adds terms to gV for forcing by external sources
123			C \| e.g. wind stress, bottom friction etc ...
124	cnh	1.13	C ==========================================================
125			C \ev
126
127			C !USES:
128	cnh	1.2	IMPLICIT NONE
129	cnh	1.1	C == Global data ==
130			#include "SIZE.h"
131			#include "EEPARAMS.h"
132			#include "PARAMS.h"
133			#include "GRID.h"
134			#include "DYNVARS.h"
135	cnh	1.2	#include "FFIELDS.h"
136
137	cnh	1.13	C !INPUT/OUTPUT PARAMETERS:
138	cnh	1.1	C == Routine arguments ==
139	jmc	1.31	C iMin,iMax :: Working range of x-index for applying forcing.
140			C jMin,jMax :: Working range of y-index for applying forcing.
141			C bi,bj :: Current tile indices
142			C kLev :: Current vertical level index
143			C myTime :: Current time in simulation
144			C myThid :: Thread Id number
145	cnh	1.1	INTEGER iMin, iMax, jMin, jMax, kLev, bi, bj
146	jmc	1.31	_RL myTime
147	adcroft	1.4	INTEGER myThid
148	cnh	1.13
149			C !LOCAL VARIABLES:
150	cnh	1.2	C == Local variables ==
151	jmc	1.31	C i,j :: Loop counters
152			C kSurface :: index of surface layer
153			INTEGER i, j
154	mlosch	1.17	INTEGER kSurface
155	cnh	1.13	CEOP
156	cnh	1.2
157	jmc	1.28	IF ( fluidIsAir ) THEN
158	jmc	1.21	kSurface = 0
159	jmc	1.28	ELSEIF ( usingPCoords ) THEN
160	mlosch	1.17	kSurface = Nr
161	jmc	1.28	ELSE
162	mlosch	1.17	kSurface = 1
163	jmc	1.28	ENDIF
164	mlosch	1.17
165	cnh	1.2	C-- Forcing term
166	jmc	1.21	#ifdef ALLOW_AIM
167			IF ( useAIM ) CALL AIM_TENDENCY_APPLY_V(
168			& iMin,iMax, jMin,jMax, bi,bj, kLev,
169	jmc	1.31	& myTime, myThid )
170	jmc	1.21	#endif /* ALLOW_AIM */
171
172	molod	1.23	#ifdef ALLOW_FIZHI
173			IF ( useFIZHI ) CALL FIZHI_TENDENCY_APPLY_V(
174			& iMin,iMax, jMin,jMax, bi,bj, kLev,
175	jmc	1.31	& myTime, myThid )
176	molod	1.23	#endif /* ALLOW_FIZHI */
177	jmc	1.31
178	heimbach	1.39	#ifdef ALLOW_MYPACKAGE
179			IF ( useMYPACKAGE ) CALL MYPACKAGE_TENDENCY_APPLY_V(
180			& iMin,iMax, jMin,jMax, bi,bj, kLev,
181			& myTime, myThid )
182			#endif /* ALLOW_MYPACKAGE */
183
184	cnh	1.2	C Add windstress momentum impulse into the top-layer
185	mlosch	1.17	IF ( kLev .EQ. kSurface ) THEN
186	jmc	1.31	DO j=1,sNy+1
187	jmc	1.32	c DO i=1,sNx
188			C-jmc: Without CD-scheme, this is OK ; but with CD-scheme, needs to cover [0:sNx+1]
189			DO i=0,sNx+1
190	cnh	1.2	gV(i,j,kLev,bi,bj) = gV(i,j,kLev,bi,bj)
191	jmc	1.26	& +foFacMom*surfaceForcingV(i,j,bi,bj)
192	heimbach	1.38	& recip_drF(kLev)_recip_hFacS(i,j,kLev,bi,bj)
193	cnh	1.2	ENDDO
194			ENDDO
195			ENDIF
196	cnh	1.1
197	heimbach	1.30	#if (defined (ALLOW_TAU_EDDY))
198	heimbach	1.29	CALL TAUEDDY_EXTERNAL_FORCING_V(
199	jmc	1.31	I iMin,iMax, jMin,jMax, bi,bj, kLev,
200			I myTime, myThid )
201	heimbach	1.29	#endif
202
203	jmc	1.31	#ifdef ALLOW_OBCS
204	heimbach	1.16	IF (useOBCS) THEN
205			CALL OBCS_SPONGE_V(
206	jmc	1.31	I iMin,iMax, jMin,jMax, bi,bj, kLev,
207			I myTime, myThid )
208	heimbach	1.16	ENDIF
209	heimbach	1.14	#endif
210
211	cnh	1.1	RETURN
212			END
213	jmc	1.31
214			C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-\|--+----\|
215	cnh	1.13	CBOP
216			C !ROUTINE: EXTERNAL_FORCING_T
217			C !INTERFACE:
218	cnh	1.1	SUBROUTINE EXTERNAL_FORCING_T(
219	jmc	1.31	I iMin,iMax, jMin,jMax, bi,bj, kLev,
220			I myTime, myThid )
221	cnh	1.13	C !DESCRIPTION: \bv
222			C ==========================================================
223	jmc	1.31	C \| S/R EXTERNAL_FORCING_T
224			C \| o Contains problem specific forcing for temperature.
225	cnh	1.13	C ==========================================================
226	jmc	1.31	C \| Adds terms to gT for forcing by external sources
227			C \| e.g. heat flux, climatalogical relaxation, etc ...
228	cnh	1.13	C ==========================================================
229			C \ev
230
231			C !USES:
232	cnh	1.2	IMPLICIT NONE
233	cnh	1.1	C == Global data ==
234			#include "SIZE.h"
235			#include "EEPARAMS.h"
236			#include "PARAMS.h"
237			#include "GRID.h"
238			#include "DYNVARS.h"
239			#include "FFIELDS.h"
240	heimbach	1.7
241	cnh	1.13	C !INPUT/OUTPUT PARAMETERS:
242	cnh	1.1	C == Routine arguments ==
243	jmc	1.31	C iMin,iMax :: Working range of x-index for applying forcing.
244			C jMin,jMax :: Working range of y-index for applying forcing.
245			C bi,bj :: Current tile indices
246			C kLev :: Current vertical level index
247			C myTime :: Current time in simulation
248			C myThid :: Thread Id number
249	cnh	1.1	INTEGER iMin, iMax, jMin, jMax, kLev, bi, bj
250	jmc	1.31	_RL myTime
251	adcroft	1.4	INTEGER myThid
252	cnh	1.1
253	cnh	1.13	C !LOCAL VARIABLES:
254	cnh	1.2	C == Local variables ==
255	jmc	1.31	C i,j :: Loop counters
256			C kSurface :: index of surface layer
257			INTEGER i, j
258	mlosch	1.17	INTEGER kSurface
259	jmc	1.31	CEOP
260	jmc	1.24	#ifdef SHORTWAVE_HEATING
261			integer two
262			_RL minusone
263			parameter (two=2,minusone=-1.)
264			_RL swfracb(two)
265			INTEGER kp1
266			#endif
267	cnh	1.2
268	jmc	1.28	IF ( fluidIsAir ) THEN
269	jmc	1.21	kSurface = 0
270	jmc	1.28	ELSEIF ( usingPCoords ) THEN
271	mlosch	1.17	kSurface = Nr
272	jmc	1.28	ELSE
273	mlosch	1.17	kSurface = 1
274	jmc	1.28	ENDIF
275	mlosch	1.17
276	cnh	1.2	C-- Forcing term
277	jmc	1.21	#ifdef ALLOW_AIM
278			IF ( useAIM ) CALL AIM_TENDENCY_APPLY_T(
279			& iMin,iMax, jMin,jMax, bi,bj, kLev,
280	jmc	1.31	& myTime, myThid )
281	jmc	1.21	#endif /* ALLOW_AIM */
282
283	molod	1.23	#ifdef ALLOW_FIZHI
284			IF ( useFIZHI ) CALL FIZHI_TENDENCY_APPLY_T(
285			& iMin,iMax, jMin,jMax, bi,bj, kLev,
286	jmc	1.31	& myTime, myThid )
287	molod	1.23	#endif /* ALLOW_FIZHI */
288	heimbach	1.25
289	heimbach	1.39	#ifdef ALLOW_MYPACKAGE
290			IF ( useMYPACKAGE ) CALL MYPACKAGE_TENDENCY_APPLY_T(
291			& iMin,iMax, jMin,jMax, bi,bj, kLev,
292			& myTime, myThid )
293			#endif /* ALLOW_MYPACKAGE */
294
295	cnh	1.2	C Add heat in top-layer
296	mlosch	1.17	IF ( kLev .EQ. kSurface ) THEN
297	jmc	1.31	DO j=1,sNy
298			DO i=1,sNx
299	cnh	1.2	gT(i,j,kLev,bi,bj)=gT(i,j,kLev,bi,bj)
300	jmc	1.26	& +surfaceForcingT(i,j,bi,bj)
301	heimbach	1.38	& recip_drF(kLev)_recip_hFacC(i,j,kLev,bi,bj)
302	cnh	1.2	ENDDO
303			ENDDO
304			ENDIF
305	adcroft	1.5
306	mlosch	1.37	#ifdef ALLOW_SHELFICE
307			IF ( useShelfIce )
308			& CALL SHELFICE_FORCING_T(
309			I iMin,iMax, jMin,jMax, bi,bj, kLev,
310			I myTime, myThid )
311			#endif /* ALLOW_SHELFICE */
312
313	adcroft	1.5	#ifdef SHORTWAVE_HEATING
314			C Penetrating SW radiation
315	jmc	1.31	c IF ( usePenetratingSW ) THEN
316			swfracb(1)=abs(rF(klev))
317			swfracb(2)=abs(rF(klev+1))
318			CALL SWFRAC(
319	heimbach	1.8	I two,minusone,
320	jmc	1.31	I myTime,myThid,
321	dimitri	1.18	U swfracb)
322	jmc	1.31	kp1 = klev+1
323			IF (klev.EQ.Nr) THEN
324	jmc	1.24	kp1 = klev
325			swfracb(2)=0. _d 0
326	jmc	1.31	ENDIF
327			DO j=1,sNy
328			DO i=1,sNx
329			gT(i,j,klev,bi,bj) = gT(i,j,klev,bi,bj)
330	jmc	1.24	& -Qsw(i,j,bi,bj)(swfracb(1)maskC(i,j,klev,bi,bj)
331			& -swfracb(2)*maskC(i,j,kp1, bi,bj))
332	jmc	1.27	& recip_Cprecip_rhoConst
333	heimbach	1.38	& recip_drF(klev)_recip_hFacC(i,j,kLev,bi,bj)
334	jmc	1.31	ENDDO
335	adcroft	1.5	ENDDO
336	jmc	1.31	c ENDIF
337	adcroft	1.5	#endif
338	heimbach	1.14
339	stephd	1.35	#ifdef ALLOW_RBCS
340			if (useRBCS) then
341			call RBCS_ADD_TENDENCY(bi,bj,klev, 1,
342			& myTime, myThid )
343			endif
344			#endif
345
346	jmc	1.31	#ifdef ALLOW_OBCS
347	heimbach	1.16	IF (useOBCS) THEN
348			CALL OBCS_SPONGE_T(
349	jmc	1.31	I iMin,iMax, jMin,jMax, bi,bj, kLev,
350			I myTime, myThid )
351	heimbach	1.16	ENDIF
352	heimbach	1.14	#endif
353
354	cnh	1.1	RETURN
355			END
356	jmc	1.31
357			C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-\|--+----\|
358	cnh	1.13	CBOP
359			C !ROUTINE: EXTERNAL_FORCING_S
360			C !INTERFACE:
361	cnh	1.1	SUBROUTINE EXTERNAL_FORCING_S(
362	jmc	1.31	I iMin,iMax, jMin,jMax, bi,bj, kLev,
363			I myTime, myThid )
364	cnh	1.13
365			C !DESCRIPTION: \bv
366			C ==========================================================
367	jmc	1.31	C \| S/R EXTERNAL_FORCING_S
368			C \| o Contains problem specific forcing for merid velocity.
369	cnh	1.13	C ==========================================================
370	jmc	1.31	C \| Adds terms to gS for forcing by external sources
371			C \| e.g. fresh-water flux, climatalogical relaxation, etc ...
372	cnh	1.13	C ==========================================================
373			C \ev
374
375			C !USES:
376	cnh	1.2	IMPLICIT NONE
377	cnh	1.1	C == Global data ==
378			#include "SIZE.h"
379			#include "EEPARAMS.h"
380			#include "PARAMS.h"
381			#include "GRID.h"
382			#include "DYNVARS.h"
383	cnh	1.2	#include "FFIELDS.h"
384	cnh	1.1
385	cnh	1.13	C !INPUT/OUTPUT PARAMETERS:
386	cnh	1.1	C == Routine arguments ==
387	jmc	1.31	C iMin,iMax :: Working range of x-index for applying forcing.
388			C jMin,jMax :: Working range of y-index for applying forcing.
389			C bi,bj :: Current tile indices
390			C kLev :: Current vertical level index
391			C myTime :: Current time in simulation
392			C myThid :: Thread Id number
393	cnh	1.1	INTEGER iMin, iMax, jMin, jMax, kLev, bi, bj
394	jmc	1.31	_RL myTime
395	adcroft	1.4	INTEGER myThid
396	cnh	1.2
397	cnh	1.13	C !LOCAL VARIABLES:
398	cnh	1.2	C == Local variables ==
399	jmc	1.31	C i,j :: Loop counters
400			C kSurface :: index of surface layer
401			INTEGER i, j
402	mlosch	1.17	INTEGER kSurface
403	cnh	1.13	CEOP
404	cnh	1.2
405	jmc	1.28	IF ( fluidIsAir ) THEN
406	jmc	1.21	kSurface = 0
407	jmc	1.28	ELSEIF ( usingPCoords ) THEN
408	mlosch	1.17	kSurface = Nr
409	jmc	1.28	ELSE
410	mlosch	1.17	kSurface = 1
411	jmc	1.28	ENDIF
412	mlosch	1.17
413	cnh	1.2	C-- Forcing term
414	jmc	1.21	#ifdef ALLOW_AIM
415			IF ( useAIM ) CALL AIM_TENDENCY_APPLY_S(
416			& iMin,iMax, jMin,jMax, bi,bj, kLev,
417	jmc	1.31	& myTime, myThid )
418	jmc	1.21	#endif /* ALLOW_AIM */
419
420	molod	1.23	#ifdef ALLOW_FIZHI
421			IF ( useFIZHI ) CALL FIZHI_TENDENCY_APPLY_S(
422			& iMin,iMax, jMin,jMax, bi,bj, kLev,
423	jmc	1.31	& myTime, myThid )
424	molod	1.23	#endif /* ALLOW_FIZHI */
425	heimbach	1.25
426	heimbach	1.39	#ifdef ALLOW_MYPACKAGE
427			IF ( useMYPACKAGE ) CALL MYPACKAGE_TENDENCY_APPLY_S(
428			& iMin,iMax, jMin,jMax, bi,bj, kLev,
429			& myTime, myThid )
430			#endif /* ALLOW_MYPACKAGE */
431
432	cnh	1.2	C Add fresh-water in top-layer
433	mlosch	1.17	IF ( kLev .EQ. kSurface ) THEN
434	jmc	1.31	DO j=1,sNy
435			DO i=1,sNx
436	cnh	1.2	gS(i,j,kLev,bi,bj)=gS(i,j,kLev,bi,bj)
437	jmc	1.26	& +surfaceForcingS(i,j,bi,bj)
438	heimbach	1.38	& recip_drF(kLev)_recip_hFacC(i,j,kLev,bi,bj)
439	cnh	1.2	ENDDO
440			ENDDO
441			ENDIF
442	heimbach	1.14
443	mlosch	1.37	#ifdef ALLOW_SHELFICE
444			IF ( useShelfIce )
445			& CALL SHELFICE_FORCING_S(
446			I iMin,iMax, jMin,jMax, bi,bj, kLev,
447			I myTime, myThid )
448			#endif /* ALLOW_SHELFICE */
449
450	stephd	1.35	#ifdef ALLOW_RBCS
451			if (useRBCS) then
452			call RBCS_ADD_TENDENCY(bi,bj,klev, 2,
453			& myTime, myThid )
454			endif
455			#endif
456
457	jmc	1.31	#ifdef ALLOW_OBCS
458	heimbach	1.16	IF (useOBCS) THEN
459			CALL OBCS_SPONGE_S(
460	jmc	1.31	I iMin,iMax, jMin,jMax, bi,bj, kLev,
461			I myTime, myThid )
462	heimbach	1.16	ENDIF
463	heimbach	1.14	#endif
464	cnh	1.1
465			RETURN
466			END