/[MITgcm]/MITgcm_contrib/dcarroll/highres_darwin/code/shelfice_thermodynamics.F
ViewVC logotype

Contents of /MITgcm_contrib/dcarroll/highres_darwin/code/shelfice_thermodynamics.F

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


Revision 1.1 - (show annotations) (download)
Sun Sep 22 21:23:47 2019 UTC (5 years, 10 months ago) by dcarroll
Branch: MAIN
CVS Tags: HEAD
Initial check in of high resolution Darwin simulation code

1 C $Header: /home/ubuntu/mnt/e9_copy/MITgcm/pkg/shelfice/shelfice_thermodynamics.F,v 1.47 2015/12/17 01:52:05 jmc Exp $
2 C $Name: $
3
4 #include "SHELFICE_OPTIONS.h"
5 #ifdef ALLOW_AUTODIFF
6 # include "AUTODIFF_OPTIONS.h"
7 #endif
8 #ifdef ALLOW_CTRL
9 # include "CTRL_OPTIONS.h"
10 #endif
11
12 CBOP
13 C !ROUTINE: SHELFICE_THERMODYNAMICS
14 C !INTERFACE:
15 SUBROUTINE SHELFICE_THERMODYNAMICS(
16 I myTime, myIter, myThid )
17 C !DESCRIPTION: \bv
18 C *=============================================================*
19 C | S/R SHELFICE_THERMODYNAMICS
20 C | o shelf-ice main routine.
21 C | compute temperature and (virtual) salt flux at the
22 C | shelf-ice ocean interface
23 C |
24 C | stresses at the ice/water interface are computed in separate
25 C | routines that are called from mom_fluxform/mom_vecinv
26
27 CIGF | ASSUMES
28 C--- | * SHELFICEconserve = true
29 C *=============================================================*
30 C \ev
31
32 C !USES:
33 IMPLICIT NONE
34
35 C === Global variables ===
36 #include "SIZE.h"
37 #include "EEPARAMS.h"
38 #include "PARAMS.h"
39 #include "GRID.h"
40 #include "DYNVARS.h"
41 #include "FFIELDS.h"
42 #include "SHELFICE.h"
43 #include "SHELFICE_COST.h"
44 #ifdef ALLOW_AUTODIFF
45 # include "CTRL_SIZE.h"
46 # include "ctrl.h"
47 # include "ctrl_dummy.h"
48 #endif /* ALLOW_AUTODIFF */
49 #ifdef ALLOW_AUTODIFF_TAMC
50 # ifdef SHI_ALLOW_GAMMAFRICT
51 # include "tamc.h"
52 # include "tamc_keys.h"
53 # endif /* SHI_ALLOW_GAMMAFRICT */
54 #endif /* ALLOW_AUTODIFF_TAMC */
55
56 C !INPUT/OUTPUT PARAMETERS:
57 C === Routine arguments ===
58 C myIter :: iteration counter for this thread
59 C myTime :: time counter for this thread
60 C myThid :: thread number for this instance of the routine.
61 _RL myTime
62 INTEGER myIter
63 INTEGER myThid
64
65 #ifdef ALLOW_SHELFICE
66 C !LOCAL VARIABLES :
67 C === Local variables ===
68 C I,J,K,Kp1,bi,bj :: loop counters
69 C tLoc, sLoc, pLoc :: local potential temperature, salinity, pressure
70 C theta/saltFreeze :: temperature and salinity of water at the
71 C ice-ocean interface (at the freezing point)
72 C freshWaterFlux :: local variable for fresh water melt flux due
73 C to melting in kg/m^2/s
74 C (negative density x melt rate)
75 C iceFrontCellThickness :: the ratio of the horizontal length
76 C of the ice front in each model grid cell
77 C divided by the grid cell area. The "thickness"
78 C of the colum perpendicular to the front
79 C iceFrontWidth :: the width of the ice front.
80
81 INTEGER I,J,K,Kp1
82 INTEGER bi,bj
83 INTEGER CURI, CURJ, FRONT_K
84
85 _RL tLoc
86 _RL sLoc
87 _RL pLoc
88 _RL uLoc
89 _RL vLoc
90 _RL wLoc
91 _RL speedLoc
92
93 C#ifndef SHI_USTAR_WETPOINT
94 C _RL uLoc(1-olx:snx+olx,1-oly:sny+oly)
95 C _RL vLoc(1-olx:snx+olx,1-oly:sny+oly)
96 C#endif
97 C _RL velSq(1-olx:snx+olx,1-oly:sny+oly)
98
99 _RL freshWaterFlux
100
101 _RL ice_bottom_Z_C, seafloor_N
102 _RL wet_top_Z_N, wet_bottom_Z_N
103 _RL iceFrontWetContact_Z_max, iceFrontContact_Z_min
104 _RL iceFrontContact_H
105 _RL iceFrontVertContactFrac, iceFrontCellThickness
106 _RL iceFrontWidth, iceFrontFaceArea
107 _RL thermalConductionDistance, thermalConductionTemp
108 _RL tmpHeatFlux, tmpFWFLX
109 _RL tmpForcingT, tmpForcingS
110 _RL tmpFac, icfgridareaFrac
111 INTEGER SI
112
113 #ifdef ALLOW_DIAGNOSTICS
114 _RL uStarDiag(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
115 #endif /* ALLOW_DIAGNOSTICS */
116
117 _RL epsilon_H
118
119 #ifdef ALLOW_SHIFWFLX_CONTROL
120 _RL xx_shifwflx_loc(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy)
121 #endif
122
123 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
124
125 C-- minimum fraction of a cell adjacent to an ice front that must be
126 C-- wet for exchange to happen
127 epsilon_H = 1. _d -03
128
129 C-- hard coded for now.
130 thermalConductionDistance = 100.0 _d 0
131 thermalConductionTemp = -20.0 _d 0
132 icfgridareaFrac = 1.0 _d 0
133
134 C heat flux into the ice shelf, default is diffusive flux
135 C (Holland and Jenkins, 1999, eq.21)
136
137 DO bj = myByLo(myThid), myByHi(myThid)
138 DO bi = myBxLo(myThid), myBxHi(myThid)
139 DO J = 1-OLy,sNy+OLy
140 DO I = 1-OLx,sNx+OLx
141 shelfIceHeatFlux (I,J,bi,bj) = 0. _d 0
142 shelfIceFreshWaterFlux(I,J,bi,bj) = 0. _d 0
143 SHIICFHeatFlux (I,J,bi,bj) = 0. _d 0
144 SHIICFFreshWaterFlux(I,J,bi,bj) = 0. _d 0
145 shelficeForcingT (I,J,bi,bj) = 0. _d 0
146 shelficeForcingS (I,J,bi,bj) = 0. _d 0
147 shelficeForcingTR (I,J,bi,bj) = 0. _d 0
148
149 DO K = 1, NR
150 iceFrontHeatFlux(I,J,K,bi,bj) = 0. _d 0
151 iceFrontFreshWaterFlux(I,J,K,bi,bj) = 0. _d 0
152 iceFrontForcingT(I,J,K,bi,bj) = 0. _d 0
153 iceFrontForcingS(I,J,K,bi,bj) = 0. _d 0
154 iceFrontForcingTR(I,J,K,bi,bj) = 0. _d 0
155 ENDDO /* K */
156
157 ENDDO /* I */
158 ENDDO /* J */
159
160 C 3-D velocity-dependent turbulent transfer coefficients
161
162 DO J = 1-OLy,sNy+OLy
163 DO I = 1-OLx,sNx+OLx
164 DO K = 1, NR
165
166 C make local copies of velocity
167 uLoc = 0.5*(uVel(I,J,K,bi,bj)+uVel(I+1,J,K,bi,bj))
168 vLoc = 0.5*(vVel(I,J,K,bi,bj)+vVel(I,J+1,K,bi,bj))
169 wLoc = 0.5*(wVel(I,J,K,bi,bj)+wVel(I,J,K+1,bi,bj))
170 speedLoc = SQRT(uLoc*uLoc + vLoc*vLoc + wLoc*wLoc)
171
172 #ifdef ALLOW_shiTransCoeff_3d
173 shiTransCoeffT(I,J,K,bi,bj) = 1.1 _d -3 *
174 & abs(speedLoc) + 1. _d -4
175
176 shiTransCoeffS(I,J,K,bi,bj) = 3.1 _d -5 *
177 & abs(speedLoc) + 5.05 _d -7
178 #endif
179 ENDDO /* K */
180 ENDDO /* I */
181 ENDDO /* J */
182
183 C-- First ice front then ice shelf. Loop through each i,j point
184 C-- process ice fronts in k, then process ice shelf.
185 DO J = 1, sNy
186 DO I = 1, sNx
187
188 C-- The K index where the ice front ends (0 if no ice front)
189 FRONT_K = K_icefront(I,J,bi,bj)
190
191 C-- If there is an ice front at this (I,J) continue
192 IF (FRONT_K .GT. 0) THEN
193
194 C-- Loop through all depths where the ice front is fround
195 DO K = 1, FRONT_K
196 C-- Loop around the four laterally neighboring cells of the ice front.
197 C-- If any neighboring points has wet volume in contact with the ice
198 C-- front at (I,J) then calculate ice-ocean exchanges.
199 C-- The four laterally neighboring point are at (CURI,CURJ)
200 DO SI = 1,4
201 IF (SI .EQ. 1) THEN
202 C-- Looking to right
203 CURI = I+1
204 CURJ = J
205
206 iceFrontWidth = dyG(I+1,J,bi,bj)
207
208 ELSEIF (SI .EQ. 2) THEN
209 C-- Looking to LEFT
210 CURI = I-1
211 CURJ = J
212
213 iceFrontWidth = dyG(I,J,bi,bj)
214 ELSEIF (SI .EQ. 3) THEN
215 C-- Looking to NORTH
216 CURI = I
217 CURJ = J+1
218
219 iceFrontWidth = dxG(I,J+1,bi,bj)
220 ELSEIF (SI .EQ. 4) THEN
221 C-- Looking to south
222 CURI = I
223 CURJ = J-1
224
225 iceFrontWidth = dxG(I,J,bi,bj)
226 endif
227
228 C-- cell depth describes the average distance
229 C-- perpendicular to the ice front fact
230
231 iceFrontCellThickness = 0. _d 0
232 IF(iceFrontWidth.NE.0. _d 0)
233 & iceFrontCellThickness = RA(CURI,CURJ,bi,bj)
234 & /iceFrontWidth
235 iceFrontFaceArea = DRF(K)*iceFrontWidth
236
237 C-- First, make sure the adjacent point has at least some water in it.
238 IF (_hFacC(CURI,CURJ,K,bi,bj) .GT. zeroRL) THEN
239
240 C-- we need to determine how much of the ice front is in contact with
241 C-- water in the neighboring grid cell at this depth level.
242
243 C-- 1. Determine the top depth with water in the current cell
244 C-- 2. Determine the top depth with water in the neighbor cell
245 C-- 3. Determine the depth where water gap between (1) and (2).
246 C-- 4. If there is a gap then ice front is in contact with water in
247 C-- the neighboring cell
248
249 C-- ice_bottom_Z_C: the depth (m) of the bottom of the ice in the
250 C-- current cell. Bounded between rF(K) and rF(K+1).
251 C-- * If the ice extends past the bottom of the cell then
252 C-- ice_bottom_Z_C = rF(K+1)
253 C-- [rF(k) >= ice_bottom_Z_C >= rF(K+1)] (rF is negative)
254 ice_bottom_Z_C = max(rF(K+1),
255 & min(Ro_surf(I,J, bi,bj), rF(K)))
256
257 C-- wet_top_Z_N: the depth (m) of the bottom of the ice in the
258 C-- neighboring grid. If the neighboring cell has ice in
259 C-- (in the form of a shelf or front) then wet_top_Z_N is
260 C-- the depth of this neighboring ice.
261 C--
262 C-- * If neighbor cell has no ice, then Ro_surf = 0 and
263 C-- wet_top_Z_N = rF(K)
264 C-- [rF(k) >= wet_top_Z_N >= rF(K+1)] (rF is negative)
265
266 wet_top_Z_N = max(rF(K+1),
267 & min(Ro_surf(CURI,CURJ, bi,bj), rF(K)))
268
269 C-- wet_bottom_Z_N: the depth (m) of the bottom of the wet part of the
270 C-- neighboring cell. If the seafloor reaches into
271 C-- the grid cell then the bottom of the wet part of the
272 C-- grid cell is at the seafloor.
273 C--
274 C-- * If the seafloor is deeper than this grid cell then
275 C-- wet_bottom_Z = rF(K+1)
276 C-- * If the seafloor is shallower than this grid cell then
277 C-- wet_bottom_Z = rF(K)
278 C-- * If the seafloor reaches partly into this grid cell
279 C-- then wet_bottom_Z = R_low
280
281 C-- [rF(k) >= wet_bottom_Z >= rF(K+1)] (rF is negative)
282
283 wet_bottom_Z_N = min(rF(K),
284 & max(R_low(CURI,CURJ, bi,bj), rF(K+1)))
285
286 C-- iceFrontWetContact_Z_max: The deepest point where the
287 C-- the ice front at (I,J) is in contact with water
288 C-- in the neighboring cell. The shallower of
289 C-- wet_bottom_Z_N (seafloor depth of neighboring point) and
290 C-- ice_bottom_Z_C (bottom of ice front in this center cell).
291
292 C-- * wet_bottom_Z_N if the seafloor of the neighboring
293 C-- cell is shallower than the ice draft at (I,J).
294 C-- * ice_bottom_Z_C if the ice draft at (I,J) is shallower
295 C-- than the seafloor of the neighboring cell.
296
297 IF (ice_bottom_Z_C .GT. wet_bottom_Z_N) THEN
298 iceFrontWetContact_Z_max = ice_bottom_Z_C
299 ELSE
300 iceFrontWetContact_Z_max = wet_bottom_Z_N
301 ENDIF
302
303 C-- The shallowest depth where the ice front at (I,J) is in contact
304 C-- with water in the neighboring cell. If the neighboring cell has
305 C-- no ice draft then wet_top_Z_N = rF(k), the top of the cell.
306 C-- Otherwise, the shallowest depth where the ice front at (I,J) can
307 C-- be in in contact with water (not ice) in (CURI, CURJ)
308 C-- is wet_top_Z_N.
309
310 C-- the fraction of the grid cell height that has ice draft in contact
311 C-- with water in the neighboring cell.
312 iceFrontVertContactFrac =
313 & (wet_top_Z_N - iceFrontWetContact_Z_max)/ DRF(K)
314
315
316 C-- Only proceed if iceFrontVertContactFrac is > 0, the
317 C-- ice draft at (I,J)
318 C-- is in contact with some water in the neighboring grid cell.
319 IF (iceFrontVertContactFrac .GT. epsilon_H) THEN
320 tLoc = theta(CURI,CURJ,K,bi,bj)
321 sLoc = MAX(salt(CURI,CURJ,K,bi,bj), zeroRL)
322
323 C-- use pressure at the halfway point between the top and bottom of
324 C-- points of the ice front where the ice front is in contact with
325 C-- open water.
326 pLoc = 0.5 _d 0 * ABS(wet_top_Z_N +
327 & iceFrontWetContact_Z_max)
328
329 CALL SHELFICE_SOLVE4FLUXES(
330 I tLoc, sLoc, pLoc,
331 #ifndef ALLOW_shiTransCoeff_3d
332 I shiTransCoeffT(CURI,CURJ,bi,bj),
333 I shiTransCoeffS(CURI,CURJ,bi,bj),
334 #else
335 I shiTransCoeffT(CURI,CURJ,K,bi,bj),
336 I shiTransCoeffS(CURI,CURJ,K,bi,bj),
337 #endif
338 I thermalConductionDistance,
339 I thermalConductionTemp,
340 O tmpHeatFlux, tmpFWFLX,
341 O tmpForcingT, tmpForcingS,
342 I bi, bj, myTime, myIter, myThid )
343
344 C-- fluxes and forcing must be scaled by iceFrontVertContactFract and
345 C-- iceFrontContactFrac some fraction of the heigth and width of the
346 C-- grid cell face may not ice in contact with water.
347
348 C tmpHeatFlux and tmpFWFLX come as W/m^2 and kg/m^2/s respectively
349 C-- but these rates only apply to the
350 C-- fraction of the grid cell that has ice in contact with seawater.
351 C-- we must scale by iceFrontVertContactFrac to get to the average
352 C-- fluxes in this grid cell.
353
354 C-- In units W/m^2
355 iceFrontHeatFlux(CURI,CURJ,K,bi,bj) =
356 & iceFrontHeatFlux(CURI,CURJ,K,bi,bj) +
357 & tmpHeatFlux*iceFrontVertContactFrac
358
359 C In units of kg/m^2/s
360 iceFrontFreshWaterFlux(CURI,CURJ,K,bi,bj) =
361 & iceFrontFreshWaterFlux(CURI,CURJ,K,bi,bj) +
362 & tmpFWFLX*iceFrontVertContactFrac
363
364 iceFrontForcingTR(CURI,CURJ,K,bi,bj) =
365 & iceFrontFreshWaterFlux(CURI,CURJ,K,bi,bj) *
366 & mass2rUnit
367
368 C ow - 06/29/2018
369 C ow - Verticallly sum up the 3D icefront heat and freshwater fluxes to
370 C ow - compute the total flux for the water column. The shelfice fluxes,
371 C ow - which are 2D, will be added later. NOTE that only
372 C ow - ice-front melts below shelf-ice are included to be consistent
373 C ow - with Rignot's data
374 if(k.GE.kTopC(I,J,bi,bj))then
375 if(RA(CURI,CURJ,bi,bj).NE.0. _d 0)then
376 icfgridareaFrac =
377 & iceFrontFaceArea/RA(CURI,CURJ,bi,bj)
378 SHIICFHeatFlux(CURI,CURJ,bi,bj) =
379 & SHIICFHeatFlux(CURI,CURJ,bi,bj) +
380 & iceFrontHeatFlux(CURI,CURJ,K,bi,bj)
381 & * icfgridareaFrac
382 SHIICFFreshWaterFlux(CURI,CURJ,bi,bj) =
383 & SHIICFFreshWaterFlux(CURI,CURJ,bi,bj) +
384 & iceFrontFreshWaterFlux(CURI,CURJ,K,bi,bj)
385 & * icfgridareaFrac
386 endif
387 endif
388 C iceFrontForcing[T,S] X m/s but these rates only apply to the
389 C-- fraction of the grid cell that has ice in contact with seawater.
390 C-- we must scale by iceFrontVertContactFrac to get to the average
391 C-- fluxes in this grid cell. We must also divide the by the length
392 C-- of the grid cell perpendicular to the face.
393
394 IF (iceFrontCellThickness .NE. 0. _d 0) THEN
395 C In units of K / s
396 iceFrontForcingT(CURI,CURJ,K,bi,bj) =
397 & iceFrontForcingT(CURI,CURJ,K,bi,bj) +
398 & tmpForcingT/iceFrontCellThickness*
399 & iceFrontVertContactFrac
400
401 C In units of psu /s
402 iceFrontForcingS(CURI,CURJ,K,bi,bj) =
403 & iceFrontForcingS(CURI,CURJ,K,bi,bj) +
404 & tmpForcingS/iceFrontCellThickness*
405 & iceFrontVertContactFrac
406
407 ENDIF /* iceFrontCellThickness */
408 C In units of kg /s
409 addMass(CURI,CURJ,K,bi,bj) =
410 & addMass(CURI,CURJ,K,bi,bj) -
411 & tmpFWFLX*iceFrontFaceArea*
412 & iceFrontVertContactFrac
413 ENDIF /* iceFrontVertContactFrac */
414 ENDIF /* hFacC(CURI,CURJ,K,bi,bj) */
415 ENDDO /* SI loop for adjacent cells */
416 ENDDO /* K LOOP */
417 ENDIF /* FRONT K */
418
419 C-- ice shelf
420 K = kTopC(I,J,bi,bj)
421
422 C-- If there is an ice front at this (I,J) continue
423 C-- I am assuming K is only .GT. when there is at least some
424 C-- nonzero wet point below the shelf in the grid cell.
425 IF (K .GT. 0) THEN
426 C-- Initialize these values to zero
427 pLoc = 0 _d 0
428 tLoc = 0 _d 0
429 sLoc = 0 _d 0
430
431 C-- make local copies of temperature, salinity and depth
432 C-- (pressure in deci-bar) underneath the ice
433 C-- for the ice shelf case we use hydrostatic pressure at the ice
434 C-- base of the ice shelf, top of the cavity.
435
436 pLoc = ABS(R_shelfIce(I,J,bi,bj))
437 tLoc = theta(I,J,K,bi,bj)
438 sLoc = MAX(salt(I,J,K,bi,bj), zeroRL)
439
440 CALL SHELFICE_SOLVE4FLUXES(
441 I tLoc, sLoc, pLoc,
442 #ifndef ALLOW_shiTransCoeff_3d
443 I shiTransCoeffT(I,J,bi,bj),
444 I shiTransCoeffS(I,J,bi,bj),
445 #else
446 I shiTransCoeffT(I,J,K,bi,bj),
447 I shiTransCoeffS(I,J,K,bi,bj),
448 #endif
449 I pLoc, thermalConductionTemp,
450 O tmpHeatFlux, tmpFWFLX,
451 O tmpForcingT, tmpForcingS,
452 I bi, bj, myTime, myIter, myThid )
453
454 C In units of W/m^2
455 shelficeHeatFlux(I,J,bi,bj) = tmpHeatFlux
456 C In units of kg/m^2/s
457 shelfIceFreshWaterFlux(I,J,bi,bj) = tmpFWFLX
458
459 shelficeForcingTR(I,J,bi,bj) =
460 & shelfIceFreshWaterFlux(I,J,bi,bj) * mass2rUnit
461
462 C ow - 06/29/2018
463 C ow - Now add shelfice heat and freshwater fluxes
464 SHIICFHeatFlux(i,j,bi,bj) =
465 & SHIICFHeatFlux(i,j,bi,bj) +
466 & shelficeHeatFlux(i,j,bi,bj)
467 SHIICFFreshWaterFlux(i,j,bi,bj) =
468 & SHIICFFreshWaterFlux(i,j,bi,bj) +
469 & shelfIceFreshWaterFlux(i,j,bi,bj)
470 C In units of K/s -- division by drF required first
471 shelficeForcingT(I,J,bi,bj) = tmpForcingT*
472 & recip_drF(K)* _recip_hFacC(i,j,K,bi,bj)
473 C In units of psu/s -- division by drF required first
474 shelficeForcingS(I,J,bi,bj) = tmpForcingS*
475 & recip_drF(K)* _recip_hFacC(i,j,K,bi,bj)
476 C In units of kg/s -- multiplication of area required first
477 addMass(I,J,K, bi,bj) = addMass(I,J,K, bi,bj) -
478 & tmpFWFLX*RA(I,J,bi,bj)
479 ENDIF /* SHELF K > 0 */
480 ENDDO /* i */
481 ENDDO /* j */
482 ENDDO /* bi */
483 ENDDO /* bj */
484
485
486 C-- Calculate new loading anomaly (in case the ice-shelf mass was updated)
487 #ifndef ALLOW_AUTODIFF
488 c IF ( SHELFICEloadAnomalyFile .EQ. ' ' ) THEN
489 DO bj = myByLo(myThid), myByHi(myThid)
490 DO bi = myBxLo(myThid), myBxHi(myThid)
491 DO j = 1-OLy, sNy+OLy
492 DO i = 1-OLx, sNx+OLx
493 shelficeLoadAnomaly(i,j,bi,bj) = gravity
494 & *( shelficeMass(i,j,bi,bj) + rhoConst*Ro_surf(i,j,bi,bj) )
495 ENDDO
496 ENDDO
497 ENDDO
498 ENDDO
499 c ENDIF
500 #endif /* ndef ALLOW_AUTODIFF */
501
502 #ifdef ALLOW_DIAGNOSTICS
503 IF ( useDiagnostics ) THEN
504 CALL DIAGNOSTICS_FILL_RS(shelfIceFreshWaterFlux,'SHIfwFlx',
505 & 0,1,0,1,1,myThid)
506 CALL DIAGNOSTICS_FILL_RS(shelfIceHeatFlux, 'SHIhtFlx',
507 & 0,1,0,1,1,myThid)
508
509 CALL DIAGNOSTICS_FILL_RS(SHIICFFreshWaterFlux,'SHIICFfwFlx',
510 & 0,1,0,1,1,myThid)
511 CALL DIAGNOSTICS_FILL_RS(SHIICFHeatFlux, 'SHIICFhtFlx',
512 & 0,1,0,1,1,myThid)
513
514 CALL DIAGNOSTICS_FILL(iceFrontFreshWaterFlux, 'ICFfwFlx',
515 & 0,Nr,0,1,1,myThid)
516 CALL DIAGNOSTICS_FILL(iceFrontHeatFlux, 'ICFhtFlx',
517 & 0,Nr,0,1,1,myThid)
518
519 CALL DIAGNOSTICS_FILL(shelfIceForcingTR, 'SHITR ',
520 & 0,Nr,0,1,1,myThid)
521
522 CALL DIAGNOSTICS_FILL(iceFrontForcingTR, 'ICFTR ',
523 & 0,Nr,0,1,1,myThid)
524
525 C SHIForcT (Ice shelf forcing for theta [W/m2], >0 increases theta)
526 tmpFac = HeatCapacity_Cp*rUnit2mass
527 CALL DIAGNOSTICS_SCALE_FILL(shelficeForcingT,tmpFac,1,
528 & 'SHIForcT',0,1,0,1,1,myThid)
529 C SHIForcS (Ice shelf forcing for salt [g/m2/s], >0 increases salt)
530 tmpFac = rUnit2mass
531 CALL DIAGNOSTICS_SCALE_FILL(shelficeForcingS,tmpFac,1,
532 & 'SHIForcS',0,1,0,1,1,myThid)
533
534 C ICFForcT (Ice front forcing for theta [W/m2], >0 increases theta)
535 tmpFac = HeatCapacity_Cp*rUnit2mass
536 CALL DIAGNOSTICS_SCALE_FILL(iceFrontForcingT,tmpFac,1,
537 & 'ICFForcT',0,Nr,0,1,1,myThid)
538 C ICFForcS (Ice front forcing for salt [g/m2/s], >0 increases salt)
539 tmpFac = rUnit2mass
540 CALL DIAGNOSTICS_SCALE_FILL(iceFrontForcingS,tmpFac,1,
541 & 'ICFForcS',0,Nr,0,1,1,myThid)
542
543 C Transfer coefficients
544 #ifndef ALLOW_shiTransCoeff_3d
545 CALL DIAGNOSTICS_FILL(shiTransCoeffT,'SHIgammT',
546 & 0,1,0,1,1,myThid)
547 CALL DIAGNOSTICS_FILL(shiTransCoeffS,'SHIgammS',
548 & 0,1,0,1,1,myThid)
549 #else
550 CALL DIAGNOSTICS_FILL(shiTransCoeffT,'SHIgammT',
551 & 0,Nr,0,1,1,myThid)
552 CALL DIAGNOSTICS_FILL(shiTransCoeffS,'SHIgammS',
553 & 0,Nr,0,1,1,myThid)
554 #endif
555 C Friction velocity
556 #ifdef SHI_ALLOW_GAMMAFRICT
557 IF ( SHELFICEuseGammaFrict )
558 & CALL DIAGNOSTICS_FILL(uStarDiag,'SHIuStar',0,1,0,1,1,myThid)
559 #endif /* SHI_ALLOW_GAMMAFRICT */
560 ENDIF
561 #endif
562
563 #endif /* ALLOW_SHELFICE */
564 RETURN
565 END

  ViewVC Help
Powered by ViewVC 1.1.22