/[MITgcm]/MITgcm_contrib/atnguyen/code_21Dec2012_saltplume/salt_plume_apply.F
ViewVC logotype

Contents of /MITgcm_contrib/atnguyen/code_21Dec2012_saltplume/salt_plume_apply.F

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


Revision 1.7 - (show annotations) (download)
Sun May 4 09:25:16 2014 UTC (11 years, 3 months ago) by atn
Branch: MAIN
CVS Tags: HEAD
Changes since 1.6: +14 -7 lines
bug fix, making sure SPforcing[T,S] are calculated.

1 C $Header: /u/gcmpack/MITgcm_contrib/atnguyen/code_21Dec2012_saltplume/salt_plume_apply.F,v 1.6 2014/05/02 06:09:11 atn Exp $
2 C $Name: $
3
4 #include "SALT_PLUME_OPTIONS.h"
5
6 CBOP
7 C !ROUTINE: SALT_PLUME_APPLY
8 C !INTERFACE:
9 SUBROUTINE SALT_PLUME_APPLY(
10 I trIdentity, bi, bj,
11 I recip_hFac_arg,
12 I tracer,trApplyFlag,
13 U trStar,
14 I myTime, myIter, myThid )
15
16 C !DESCRIPTION: \bv
17 C *==========================================================*
18 C | SUBROUTINE SALT_PLUME_APPLY
19 C | o Apply the salt_pume-transport to tracer field
20 C *==========================================================*
21 C \ev
22
23 C !USES:
24 IMPLICIT NONE
25
26 C === Global variables ===
27 #include "SIZE.h"
28 #include "GRID.h"
29 #include "EEPARAMS.h"
30 #include "PARAMS.h"
31 #include "DYNVARS.h"
32 #include "SALT_PLUME.h"
33 #ifdef ALLOW_GENERIC_ADVDIFF
34 # include "GAD.h"
35 #endif
36
37 C !INPUT/OUTPUT PARAMETERS:
38 C === Routine arguments ===
39 C trIdentity :: tracer identification number
40 C bi,bj :: Tile indices
41 C recip_drF :: Reciprol of cell thickness
42 C recip_hFac_arg :: Reciprol of cell open-depth factor
43 C tracer :: tracer field at current time (input)
44 C trApplyFlag:: [0]=update saltplume forcing T/S terms
45 C :: [1]=update gTr tendency
46 C trStar :: future tracer field (modified)
47 C myTime :: Current time in simulation
48 C myIter :: Current time-step number
49 C myThid :: my Thread Id. number
50
51 INTEGER trIdentity, trApplyFlag
52 INTEGER bi, bj
53 _RS recip_hFac_arg(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
54 _RL tracer (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
55 _RL trStar (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
56
57 _RL myTime
58 INTEGER myIter, myThid
59
60 #ifdef ALLOW_SALT_PLUME
61 #ifdef SALT_PLUME_VOLUME
62
63 C !LOCAL VARIABLES:
64 C === Local variables ===
65 C msgBuf :: Informational/error message buffer
66 C plumetend :: forcing terms [W/m2 or kg/m2/s*psu]
67 C work :: working array
68 C CHARACTER*(MAX_LEN_MBUF) msgBuf
69 INTEGER i, j, k
70 INTEGER upward
71 LOGICAL onOffFlag
72 _RL gTr_Surf2kLev, gTr_Below2kLev, gTr_kLev2Above,
73 & dSPvolBelow2kLev, gTr_totSurf2Below,
74 & SurfVal, SEAICE_Tfrz, ConvertFac, recip_ConvertFac
75 integer kp1, Nrp1
76 _RL plumetend(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
77 _RL work(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr+1)
78
79 #ifdef ALLOW_DIAGNOSTICS
80 CHARACTER*8 diagName
81 CHARACTER*4 diagSufx
82 LOGICAL doDiagSPtend
83 C- Functions:
84 LOGICAL DIAGNOSTICS_IS_ON
85 EXTERNAL DIAGNOSTICS_IS_ON
86 #ifdef ALLOW_GENERIC_ADVDIFF
87 CHARACTER*5 GAD_DIAG_SUFX
88 EXTERNAL GAD_DIAG_SUFX
89 #endif /* ALLOW_GENERIC_ADVDIFF */
90 #endif /* ALLOW_DIAGNOSTICS */
91
92 CEOP
93
94 IF ( trApplyFlag.LT.0 .OR. trApplyFlag.GT.1) THEN
95 STOP 'S/R SALT_PLUME_APPLY: incorrect setting of trApplyFlag!'
96 ELSE
97
98 SEAICE_Tfrz = -1.96 _d 0
99
100 onOffFlag = .FALSE.
101 #ifdef ALLOW_GENERIC_ADVDIFF
102 IF ( trIdentity.EQ.GAD_TEMPERATURE ) onOffFlag = .TRUE.
103 IF ( trIdentity.EQ.GAD_SALINITY ) onOffFlag = .TRUE.
104 #endif
105 IF ( onOffFlag ) THEN
106 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
107
108 c upward = rkSign*NINT(-gravitySign)
109 upward = 1
110 IF (usingZCoords) upward = -1
111
112 IF ( trIdentity.EQ.1 ) THEN
113 SurfVal = SEAICE_Tfrz
114 ConvertFac = HeatCapacity_Cp*rhoConst
115 recip_ConvertFac = recip_Cp*mass2rUnit
116 #ifdef ALLOW_DIAGNOSTICS
117 IF ( useDiagnostics ) diagSufx = '_TH '
118 #endif /* ALLOW_DIAGNOSTICS */
119 ENDIF
120 IF ( trIdentity.EQ.2 ) THEN
121 SurfVal = SPbrineSconst
122 ConvertFac = rhoConst
123 recip_ConvertFac = mass2rUnit
124 #ifdef ALLOW_DIAGNOSTICS
125 IF ( useDiagnostics ) diagSufx = '_SLT'
126 #endif /* ALLOW_DIAGNOSTICS */
127 ENDIF
128
129 #ifdef ALLOW_DIAGNOSTICS
130 doDiagSPtend = .FALSE.
131 diagName = 'SPtd'
132 IF ( useDiagnostics ) THEN
133 C-- Set diagnostic suffix for the current tracer
134 #ifdef ALLOW_GENERIC_ADVDIFF
135 diagSufx = GAD_DIAG_SUFX( trIdentity, myThid )
136 #endif /* ALLOW_GENERIC_ADVDIFF */
137 diagName = 'SPtd'//diagSufx
138 doDiagSPtend = DIAGNOSTICS_IS_ON(diagName,myThid)
139 C WRITE(msgBuf,'(3A)') 'diagSufx,diagName: ',
140 C & diagSufx,diagName
141 C CALL PRINT_MESSAGE(msgBuf,standardMessageUnit,
142 C & SQUEEZE_RIGHT, myThid )
143 ENDIF
144 #endif /* ALLOW_DIAGNOSTICS */
145
146 C-- initializing:
147 Nrp1=Nr+1
148 DO k=1,Nr
149 DO j=1,OLy
150 DO i=1,OLx
151 plumetend(i,j,k) = 0. _d 0
152 work(i,j,k) = tracer(i,j,k,bi,bj)
153 ENDDO
154 ENDDO
155 ENDDO
156 DO j=1,OLy
157 DO i=1,OLx
158 work(i,j,Nrp1) = 0. _d 0
159 ENDDO
160 ENDDO
161 C-----------------
162
163 Catn: After discussion with JM, it's cleaner to remove the negative
164 c buoyancy associated with saltplumeflux from the surface lev
165 c here instead of inside salt_plume_forcing_surf.F and
166 c kpp_forcing_surf.F:
167
168 DO k=Nr,1,-1
169 kp1=k+1
170 DO j=1-OLy,sNy+OLy
171 DO i=1-OLx,sNx+OLx
172 C IF(trIdentity.EQ.GAD_SALINITY) SurfVal=SPbrineSalt(i,j,bi,bj)
173 Catn: m/s*[degC,psu]
174 gTr_totSurf2Below = SPbrineVolFlux(i,j,bi,bj)*SurfVal
175 C
176 dSPvolBelow2kLev = -dSPvolkLev2Above(i,j,kp1,bi,bj)
177 gTr_Surf2kLev = dSPvolSurf2kLev(i,j,k,bi,bj) * SurfVal
178 gTr_Below2kLev= dSPvolBelow2kLev * work(i,j,kp1)
179 Catn: gTr_kLev2Above works even for kLev=1 because this is how much
180 C volume of original [salinity,heat] associated with [SSS,SST]
181 C was replaced by same volume of brine [salt,heat(from SEAICE_Tfrz)].
182 C Note: by design, dSPvolkLev2Above already is negative
183 gTr_kLev2Above= dSPvolkLev2Above(i,j,k,bi,bj) * work(i,j,k)
184
185 C salt: [m/s * psu * kg/m3] = [kg/s/m2 psu] = unit of saltPlumeFlux
186 C theta:[m/s * kg/m3 * J/kg/degC * degC] = [W/m2]
187 plumetend(i,j,k) = ConvertFac *
188 & ( gTr_Surf2kLev + gTr_Below2kLev + gTr_kLev2Above )
189 IF(k.EQ.1) THEN
190 plumetend(i,j,k) = plumetend(i,j,k)
191 & - ConvertFac*gTr_totSurf2Below
192 ENDIF
193
194 IF(trApplyFlag.EQ.0) THEN
195 Catn: report T/S SPforcing[T,S] related to saltplumeflux for kpp
196 C and return zero to do_oceanic_phys.F; unit: [g/m2/s or W/m2]
197 trStar(i,j,k,bi,bj) = 0. _d 0
198 IF (trIdentity.EQ.GAD_SALINITY) THEN
199 SPforcingS(i,j,k,bi,bj)=plumetend(i,j,k)
200 ENDIF
201 IF (trIdentity.EQ.GAD_TEMPERATURE) THEN
202 SPforcingT(i,j,k,bi,bj)=plumetend(i,j,k)
203 ENDIF
204 ELSE
205 Catn: updating tendency gTr (gS,gT); unit: [psu/s or degC/s]
206 trStar(i,j,k,bi,bj)=trStar(i,j,k,bi,bj)+plumetend(i,j,k)
207 & *recip_drF(k)*recip_hFac_arg(i,j,k,bi,bj)
208 & *recip_ConvertFac
209 ENDIF
210 ENDDO
211 ENDDO
212 ENDDO
213
214 #ifdef ALLOW_DIAGNOSTICS
215 IF ( doDiagSPtend )
216 & CALL DIAGNOSTICS_FILL(plumetend, diagName, 0,Nr,2,bi,bj,myThid)
217 #endif /* ALLOW_DIAGNOSTICS */
218
219 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
220 C-- end if on-off-flag
221 ENDIF
222 C-- end trApplyFlag
223 ENDIF
224
225 #endif /* SALT_PLUME_VOLUME */
226 #endif /* ALLOW_SALT_PLUME */
227
228 RETURN
229 END

  ViewVC Help
Powered by ViewVC 1.1.22