/[MITgcm]/MITgcm_contrib/ifenty/seaiceAdjointCode/seaice_growth.F
ViewVC logotype

Annotation of /MITgcm_contrib/ifenty/seaiceAdjointCode/seaice_growth.F

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


Revision 1.1 - (hide annotations) (download)
Fri Jun 29 18:54:03 2007 UTC (18 years ago) by gforget
Branch: MAIN
CVS Tags: HEAD
Ian Fenty sea-ice growth routines (adjointable, in developement)

1 gforget 1.1 C $Header: /u/gcmpack/MITgcm/pkg/seaice/seaice_growth.F,v 1.10 2007/01/09 13:33:49 mlosch Exp $
2     C $Name: $
3    
4     #include "SEAICE_OPTIONS.h"
5    
6     C StartOfInterface
7     SUBROUTINE SEAICE_GROWTH( myTime, myIter, myThid )
8     C /==========================================================\
9     C | SUBROUTINE seaice_growth |
10     C | o Updata ice thickness and snow depth |
11     C |==========================================================|
12     C \==========================================================/
13     IMPLICIT NONE
14    
15     C === Global variables ===
16     #include "SIZE.h"
17     #include "EEPARAMS.h"
18     #include "PARAMS.h"
19     #include "DYNVARS.h"
20     #include "GRID.h"
21     #include "FFIELDS.h"
22     #include "SEAICE_PARAMS.h"
23     #include "SEAICE.h"
24     #include "SEAICE_FFIELDS.h"
25    
26     #ifdef ALLOW_AUTODIFF_TAMC
27     # include "tamc.h"
28     #endif
29     C === Routine arguments ===
30     C myTime - Simulation time
31     C myIter - Simulation timestep number
32     C myThid - Thread no. that called this routine.
33     _RL myTime
34     INTEGER myIter, myThid
35     C EndOfInterface(global-font-lock-mode 1)
36    
37     C === Local variables ===
38     C i,j,bi,bj - Loop counters
39    
40     INTEGER i, j, bi, bj
41     C number of surface interface layer
42     INTEGER kSurface
43    
44     C constants
45     _RL TBC, salinity_ice, SDF, ICE2WATR, ICE2SNOW
46    
47     #ifdef ALLOW_SEAICE_FLOODING
48     _RL hDraft, hFlood
49     #endif /* ALLOW_SEAICE_FLOODING */
50    
51     C QNETI - net surface heat flux under ice in W/m^2
52     C QSWO - short wave heat flux over ocean in W/m^2
53     C QSWI - short wave heat flux under ice in W/m^2
54    
55     _RL QNETI (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
56     _RL QSWO (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
57     _RL QSWI (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
58    
59     _RL QSWO_IN_FIRST_LAYER
60     & (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
61     _RL QSWO_BELOW_FIRST_LAYER
62     & (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
63    
64     _RL QSW_absorb_in_ML (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
65     _RL QSW_absorb_below_ML (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
66    
67     C actual ice thickness with upper and lower limit
68     _RL HICE_ACTUAL (1-OLx:sNx+OLx, 1-OLy:sNy+OLy)
69    
70     C actual snow thickness
71     _RL HSNOW_ACTUAL(1-OLx:sNx+OLx, 1-OLy:sNy+OLy)
72    
73     C wind speed
74     _RL UG (1-OLx:sNx+OLx, 1-OLy:sNy+OLy)
75     _RL SPEED_SQ
76    
77     C IAN
78     _RL RHOI, RHOFW,CPW,LI,QI,QS,GAMMAT,GAMMA,RHOSW,RHOSN
79     _RL FL_C1,FL_C2,FL_C3,FL_C4,deltaHS,deltaHI
80    
81     _RL NetExistingIceGrowthRate (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
82     _RL IceGrowthRateUnderExistingIce (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
83     _RL IceGrowthRateFromSurface (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
84     _RL IceGrowthRateOpenWater (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
85     _RL IceGrowthRateMixedLayer (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
86    
87     _RL PredictedTemperatureChange
88     & (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
89     _RL PredictedTemperatureChangeFromQSW
90     & (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
91     _RL PredictedTemperatureChangeFromOA_MQNET
92     & (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
93     _RL PredictedTemperatureChangeFromFIA
94     & (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
95     _RL PredictedTemperatureChangeFromNewIceVol
96     & (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
97     _RL PredictedTemperatureChangeFromF_IA_NET
98     & (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
99     _RL PredictedTemperatureChangeFromF_IO_NET
100     & (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
101    
102     _RL ExpectedIceVolumeChange (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
103     _RL ExpectedSnowVolumeChange (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
104     _RL ActualNewTotalVolumeChange(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
105     _RL ActualNewTotalSnowMelt(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
106    
107     _RL EnergyInNewTotalIceVolume (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
108     _RL NetEnergyFluxOutOfSystem (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
109    
110     _RL ResidualHeatOutOfSystem (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
111    
112     _RL SnowAccumulationRateOverIce (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
113     _RL SnowAccumulationOverIce (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
114     _RL PrecipRateOverIceSurfaceToSea (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
115    
116     _RL PotSnowMeltRateFromSurf (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
117     _RL PotSnowMeltFromSurf (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
118     _RL SnowMeltFromSurface (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
119     _RL SnowMeltRateFromSurface (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
120    
121     _RL FreshwaterContribFromSnowMelt (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
122     _RL FreshwaterContribFromIce (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
123    
124     _RL SurfHeatFluxConvergToSnowMelt (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
125     _RL EnergyToMeltSnowAndIce (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
126     _RL EnergyToMeltSnowAndIce2 (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
127    
128     C dA/dt = S_a
129     _RL S_a (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
130     C dh/dt = S_h
131     _RL S_h (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
132     _RL S_hsnow (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
133     _RL HSNOW_ORIG (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
134    
135     C F_ia - heat flux from ice to atmosphere (W/m^2)
136     C >0 causes ice growth, <0 causes snow and sea ice melt
137     _RL F_ia (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
138     _RL F_ia_net (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
139     _RL F_ia_net_before_snow (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
140     _RL F_io_net (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
141    
142     C F_ao - heat flux from atmosphere to ocean (W/m^2)
143     _RL F_ao (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
144    
145     C F_mi - heat flux from mixed layer to ice (W/m^2)
146     _RL F_mi (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
147    
148     c INTEGER IMAX
149     c _RL FACTORM
150    
151     _RL FLUX_TO_DELTA_TEMP,ENERGY_TO_DELTA_TEMP
152    
153     if ( buoyancyRelation .eq. 'OCEANICP' ) then
154     kSurface = Nr
155     else
156     kSurface = 1
157     endif
158    
159     FLUX_TO_DELTA_TEMP = SEAICE_deltaTtherm*
160     & recip_Cp*recip_rhoConst * recip_drF(1)
161    
162     ENERGY_TO_DELTA_TEMP = recip_Cp*recip_rhoConst*recip_drF(1)
163    
164     C ICE SALINITY (g/kg)
165     salinity_ice = 4.0
166    
167     C FREEZING TEMP. OF SEA WATER (deg C)
168     TBC = SEAICE_freeze
169    
170     C IAN
171    
172     c Sea ice density (kg m^-3)
173     RHOI = 917.0
174    
175     c Seawater density (kg m^-3)
176     RHOSW = 1026.0
177    
178     c Freshwater density (KG M^-3)
179     RHOFW = 1000.0
180    
181     C Snow density
182     RHOSN = 330.0
183    
184     C Heat capacity of seawater (J m^-3 K^-1)
185     CPW = 4010.0
186    
187     c latent heat of fusion for ice (J kg^-1)
188     LI = 3.340e5
189     c conversion between Joules and m^3 of ice (m^3)
190     QI = 1/rhoi/Li
191     QS = 1/RHOSN/Li
192    
193     c FOR FLOODING
194     FL_C2 = RHOI/RHOSW
195     FL_C2 = RHOSN/RHOSW
196     FL_C3 = (RHOSW-RHOI)/RHOSN
197     FL_C4 = RHOSN/RHOI
198    
199     c Timescale for melting of ice from a warm ML (3 days in seconds)
200     c Damping term for mixed layer heat to melt existing ice
201     GAMMA = dRf(1)/SEAICE_gamma_t
202    
203     DO bj=myByLo(myThid),myByHi(myThid)
204     DO bi=myBxLo(myThid),myBxHi(myThid)
205     c
206     #ifdef ALLOW_AUTODIFF_TAMC
207     act1 = bi - myBxLo(myThid)
208     max1 = myBxHi(myThid) - myBxLo(myThid) + 1
209     act2 = bj - myByLo(myThid)
210     max2 = myByHi(myThid) - myByLo(myThid) + 1
211     act3 = myThid - 1
212     max3 = nTx*nTy
213     act4 = ikey_dynamics - 1
214     iicekey = (act1 + 1) + act2*max1
215     & + act3*max1*max2
216     & + act4*max1*max2*max3
217     #endif /* ALLOW_AUTODIFF_TAMC */
218     C
219     C initialise a few fields
220     C
221     #ifdef ALLOW_AUTODIFF_TAMC
222     CADJ STORE area(:,:,:,bi,bj) = comlev1_bibj,
223     CADJ & key = iicekey, byte = isbyte
224     CADJ STORE qnet(:,:,bi,bj) = comlev1_bibj,
225     CADJ & key = iicekey, byte = isbyte
226     CADJ STORE qsw(:,:,bi,bj) = comlev1_bibj,
227     CADJ & key = iicekey, byte = isbyte
228     #endif /* ALLOW_AUTODIFF_TAMC */
229     DO J=1,sNy
230     DO I=1,sNx
231     F_ia_net (I,J) = 0.0
232     F_ia_net_before_snow(I,J) = 0.0
233     F_io_net (I,J) = 0.0
234    
235     F_ia (I,J) = 0.0
236     F_ao (I,J) = 0.0
237     F_mi (I,J) = 0.0
238    
239     QNETI(I,J) = 0.0
240     QSWO (I,J) = 0.0
241     QSWI (I,J) = 0.0
242    
243     QSWO_BELOW_FIRST_LAYER (I,J) = 0.0
244     QSWO_IN_FIRST_LAYER (I,J) = 0.0
245    
246     S_a (I,J) = 0.0
247     S_h (I,J) = 0.0
248    
249     IceGrowthRateUnderExistingIce (I,J) = 0.0
250     IceGrowthRateFromSurface (I,J) = 0.0
251     NetExistingIceGrowthRate (I,J) = 0.0
252    
253     PredictedTemperatureChange (I,J) = 0.0
254     PredictedTemperatureChangeFromQSW (I,J) = 0.0
255     PredictedTemperatureChangeFromOA_MQNET (I,J) = 0.0
256     PredictedTemperatureChangeFromFIA (I,J) = 0.0
257     PredictedTemperatureChangeFromF_IA_NET (I,J) = 0.0
258     PredictedTemperatureChangeFromF_IO_NET (I,J) = 0.0
259     PredictedTemperatureChangeFromNewIceVol (I,J) = 0.0
260    
261     IceGrowthRateOpenWater (I,J) = 0.0
262     IceGrowthRateMixedLayer (I,J) = 0.0
263    
264     ExpectedIceVolumeChange (I,J) = 0.0
265     ExpectedSnowVolumeChange (I,J) = 0.0
266     ActualNewTotalVolumeChange (I,J) = 0.0
267     ActualNewTotalSnowMelt (I,J) = 0.0
268    
269     EnergyInNewTotalIceVolume (I,J) = 0.0
270     NetEnergyFluxOutOfSystem (I,J) = 0.0
271     ResidualHeatOutOfSystem (I,J) = 0.0
272     QSW_absorb_in_ML (I,J) = 0.0
273     QSW_absorb_below_ML (I,J) = 0.0
274    
275     SnowAccumulationRateOverIce (I,J) = 0.0
276     SnowAccumulationOverIce (I,J) = 0.0
277     PrecipRateOverIceSurfaceToSea (I,J) = 0.0
278    
279     PotSnowMeltRateFromSurf (I,J) = 0.0
280     PotSnowMeltFromSurf (I,J) = 0.0
281     SnowMeltFromSurface (I,J) = 0.0
282     SnowMeltRateFromSurface (I,J) = 0.0
283     SurfHeatFluxConvergToSnowMelt (I,J) = 0.0
284    
285     FreshwaterContribFromSnowMelt (I,J) = 0.0
286     FreshwaterContribFromIce (I,J) = 0.0
287    
288     ENDDO
289     ENDDO
290    
291    
292     #ifdef ALLOW_AUTODIFF_TAMC
293     CADJ STORE heff(:,:,:,bi,bj) = comlev1_bibj,
294     CADJ & key = iicekey, byte = isbyte
295     CADJ STORE hsnow(:,:,bi,bj) = comlev1_bibj,
296     CADJ & key = iicekey, byte = isbyte
297     #endif /* ALLOW_AUTODIFF_TAMC */
298    
299     DO J=1,sNy
300     DO I=1,sNx
301     IF (AREA(I,J,2,bi,bj) .GT. 0.) THEN
302     HICE_ACTUAL(I,J) =
303     & HEFF(I,J,2,bi,bj)/AREA(I,J,2,bi,bj)
304    
305     HSNOW_ACTUAL(I,J) = HSNOW(I,J,bi,bj)/
306     & AREA(I,J,2,bi,bj)
307    
308     c ACCUMULATE SNOW
309     c Is the ice/surface below freezing or at the freezing
310     c point (melting). If it is freezing the precip is
311     c felt as snow and will accumulate over the ice. Else,
312     c precip makes its way, like all things in time, to the sea.
313     IF (TICE(I,J,bi,bj) .LT. 273.15) THEN
314     c Snow falls onto freezing surface remaining as snow
315     SnowAccumulationRateOverIce(I,J) =
316     & PRECIP(I,J,bi,bj)*RHOFW/RHOSN
317    
318     c None of the precipitation falls into the sea
319     PrecipRateOverIceSurfaceToSea(I,J) = 0.0
320    
321     ELSE
322     c The snow melts on impact is is considered
323     c nothing more than rain. Since meltponds are
324     c not explicitly represented,this rain runs
325     c immediately into the sea
326    
327     SnowAccumulationRateOverIce(I,J) = 0.0
328     C The rate of rainfall over melting ice.
329     PrecipRateOverIceSurfaceToSea(I,J)=
330     & PRECIP(I,J,bi,bj)
331     ENDIF
332    
333     c In m of mean snow thickness.
334     SnowAccumulationOverIce(I,J) =
335     & SnowAccumulationRateOverIce(I,J)
336     & *SEAICE_deltaTtherm*Area(I,J,2,bi,bj)
337    
338     ELSE
339     HEFF(I,J,2,bi,bj) = 0.0
340     HICE_ACTUAL(I,J) = 0.0
341     HSNOW_ACTUAL(I,J) = 0.0
342     HSNOW(I,J,bi,bj) = 0.0
343     ENDIF
344     HSNOW_ORIG(I,J) = HSNOW(I,J,bi,bj)
345     ENDDO
346     ENDDO
347    
348     C FIND ATM. WIND SPEED
349     DO J=1,sNy
350     DO I=1,sNx
351     #ifdef SEAICE_EXTERNAL_FORCING
352     c USE EXF PACKAGE
353     UG(I,J) = MAX(SEAICE_EPS,wspeed(I,J,bi,bj))
354     #else
355     C CALCULATE IT HERE
356     SPEED_SQ = UWIND(I,J,bi,bj)**2 + VWIND(I,J,bi,bj)**2
357     IF ( SPEED_SQ .LE. SEAICE_EPS_SQ ) THEN
358     UG(I,J)=SEAICE_EPS
359     ELSE
360     UG(I,J)=SQRT(SPEED_SQ)
361     ENDIF
362     #endif /* SEAICE_EXTERNAL_FORCING */
363     ENDDO
364     ENDDO
365    
366     #ifdef ALLOW_AUTODIFF_TAMC
367     cphCADJ STORE heff = comlev1, key = ikey_dynamics
368     cphCADJ STORE hsnow = comlev1, key = ikey_dynamics
369     cphCADJ STORE uwind = comlev1, key = ikey_dynamics
370     cphCADJ STORE vwind = comlev1, key = ikey_dynamics
371     CADJ STORE tice = comlev1, key = ikey_dynamics
372     #endif /* ALLOW_AUTODIFF_TAMC */
373    
374     C SET LAYER TEMPERATURE IN KELVIN
375     DO J=1,sNy
376     DO I=1,sNx
377     TMIX(I,J,bi,bj)=
378     & theta(I,J,kSurface,bi,bj)+273.15 _d +00
379     ENDDO
380     ENDDO
381    
382     C NOW DO ICE
383     CALL SEAICE_BUDGET_ICE(
384     I UG, HICE_ACTUAL, HSNOW_ACTUAL,
385     U TICE,
386     O F_io_net,F_ia_net,F_ia, QSWI,
387     I bi, bj)
388    
389     C-- NET HEAT FLUX TO ICE FROM MIXED LAYER (POSITIVE MEANS NET OUT)
390     DO J=1,sNy
391     DO I=1,sNx
392     #ifdef SEAICE_DEBUG
393     print *,'I,J,F_ia,F_ia_net',I,J,F_ia(I,J),F_ia_net(I,J)
394     F_ia_net_before_snow(I,J) = F_ia_net(I,J)
395     #endif
396    
397    
398     IF (Area(I,J,2,bi,bj)*HEFF(I,J,2,bi,bj) .LE. 0.) THEN
399     IceGrowthRateUnderExistingIce(I,J) = 0.0
400     IceGrowthRateFromSurface(I,J) = 0.0
401     NetExistingIceGrowthRate(I,J) = 0.0
402     ELSE
403     c The growth rate under existing ice is given by the upward
404     c ocean-ice conductive flux, F_io_net, and QI, which converts
405     c Joules to meters of ice. This quantity has units of meters
406     c of ice per second.
407     IceGrowthRateUnderExistingIce(I,J)=F_io_net(I,J)*QI
408    
409     c Snow/Ice surface heat convergence is first used to melt
410     c snow. If all of this heat convergence went into melting
411     c snow, this is the rate at which it would do it
412     c F_ia_net must be negative, -> PSMRFW is positive for melting
413     PotSnowMeltRateFromSurf(I,J)= - F_ia_net(I,J)*QS
414    
415     c This is the depth of snow that would be melted at this rate
416     c and the seaice delta t. In meters of snow.
417     PotSnowMeltFromSurf(I,J) =
418     & PotSnowMeltRateFromSurf(I,J)* SEAICE_deltaTtherm
419    
420     c If we can melt MORE than is actually there, then we will
421     c reduce the melt rate so that only that which is there
422     c is melted in one time step. In this case not all of the
423     c heat flux convergence at the surface is used to melt snow,
424     c The leftover energy is going to melt ice.
425     c SurfHeatFluxConvergToSnowMelt is the part of the total heat
426     c flux convergence going to melt snow.
427    
428     IF (PotSnowMeltFromSurf(I,J) .GE.
429     & HSNOW_ACTUAL(I,J)) THEN
430     c Snow melt and melt rate in actual snow thickness.
431     SnowMeltFromSurface(I,J) = HSNOW_ACTUAL(I,J)
432    
433     SnowMeltRateFromSurface(I,J) =
434     & SnowMeltFromSurface(I,J)/ SEAICE_deltaTtherm
435    
436     c Since F_ia_net is focused only over ice, its reduction
437     c requires knowing how much snow is actually melted
438     SurfHeatFluxConvergToSnowMelt(I,J) =
439     & -HSNOW_ACTUAL(I,J)/QS/SEAICE_deltaTtherm
440     ELSE
441     c In this case there will be snow remaining after melting.
442     c All of the surface heat convergence will be redirected to
443     c this effort.
444     SnowMeltFromSurface(I,J)=PotSnowMeltFromSurf(I,J)
445    
446     SnowMeltRateFromSurface(I,J) =
447     & PotSnowMeltRateFromSurf(I,J)
448    
449     SurfHeatFluxConvergToSnowMelt(I,J) =F_ia_net(I,J)
450     ENDIF
451    
452     c Taken across entire cell, this is the amount of freshwater
453     c to add per square meter from melting snow. Positive with
454     c the addition of freshwater to the ocean
455     c FreshwaterContribFromSnowMelt(I,J) =
456     c & AREA(I,J,2,bi,bj)*SnowMeltFromSurface(I,J)
457     c & *RHOSN/RHOFW
458    
459     c Reduce the heat flux convergence available to melt surface
460     c ice by the amount used to melt snow
461     F_ia_net(I,J) =
462     & F_ia_net(I,J)-SurfHeatFluxConvergToSnowMelt(I,J)
463    
464     IceGrowthRateFromSurface(I,J) = F_ia_net(I,J)*QI
465    
466     NetExistingIceGrowthRate(I,J) =
467     & IceGrowthRateUnderExistingIce(I,J) +
468     & IceGrowthRateFromSurface(I,J)
469     ENDIF
470     ENDDO
471     ENDDO
472    
473     c HERE WE WILL MELT SNOW AND ADJUST NET EXISTING ICE GROWTH RATE
474     C TO REFLECT REDUCTION IN SEA ICE MELT.
475    
476     C NOW DETERMINE GROWTH RATES
477     C FIRST DO OPEN WATER
478     CALL SEAICE_BUDGET_OCEAN(
479     I UG,
480     U TMIX,
481     O F_ao, QSWO,
482     I bi, bj)
483    
484     C-- NET HEAT FLUX TO ICE FROM MIXED LAYER (POSITIVE MEANS NET OUT)
485     c-- not all of the sw radiation is absorbed in the first layer, only that
486     c which is absorbed melts ice. SWFRACB is calculated in seaice_init_vari.F
487     DO J=1,sNy
488     DO I=1,sNx
489     QSWO_BELOW_FIRST_LAYER(I,J)= QSWO(I,J)*SWFRACB
490     QSWO_IN_FIRST_LAYER(I,J) = QSWO(I,J)*(1.0 - SWFRACB)
491     c IceGrowthRateOpenWater(I,J)= F_ao(I,J)*QI
492     IceGrowthRateOpenWater(I,J)= QI*
493     & (F_ao(I,J) - QSWO(I,J) + QSWO_IN_FIRST_LAYER(I,J))
494     ENDDO
495     ENDDO
496    
497    
498     #ifdef ALLOW_AUTODIFF_TAMC
499     CADJ STORE theta(:,:,:,bi,bj)= comlev1_bibj,
500     CADJ & key = iicekey, byte = isbyte
501     CADJ STORE heff(:,:,:,bi,bj) = comlev1_bibj,
502     CADJ & key = iicekey, byte = isbyte
503     #endif /* ALLOW_AUTODIFF_TAMC */
504    
505    
506     C-- NET HEAT FLUX TO ICE FROM MIXED LAYER (POSITIVE MEANS FLUX INTO ICE
507     C AND MELTING)
508     DO J=1,sNy
509     DO I=1,sNx
510    
511     C FIND THE FREEZING POINT OF SEAWATER IN THIS CELL
512     #ifdef SEAICE_VARIABLE_FREEZING_POINT
513     TBC = -0.0575 _d 0*salt(I,J,kSurface,bi,bj) +
514     & 0.0901 _d 0
515     #endif /* SEAICE_VARIABLE_FREEZING_POINT */
516    
517     c example: theta(i,j,ksurf) = 0, tbc = -2,
518     c fmi = -gamm*rhocpw * (0-(-2)) = - 2 * gamm * rhocpw,
519     c a NEGATIVE number. Heat flux INTO ice.
520    
521     F_mi(I,J) = -GAMMA*RHOSW*CPW *
522     & (theta(I,J,kSurface,bi,bj) - TBC)
523    
524     IceGrowthRateMixedLayer(I,J) = F_mi(I,J)*QI;
525     ENDDO
526     ENDDO
527    
528     #ifdef ALLOW_AUTODIFF_TAMC
529     CADJ STORE S_h(:,:) = comlev1_bibj,
530     CADJ & key = iicekey, byte = isbyte
531     #endif /* ALLOW_AUTODIFF_TAMC */
532    
533     C CALCULATE THICKNESS DERIVATIVE (S_h)
534     DO J=1,sNy
535     DO I=1,sNx
536     S_h(I,J) =
537     & NetExistingIceGrowthRate(I,J)*AREA(I,J,2,bi,bj)+
538     & (1. -AREA(I,J,2,bi,bj))*
539     & IceGrowthRateOpenWater(I,J) +
540     & IceGrowthRateMixedLayer(I,J)
541    
542     c Both the accumulation and melt rates are in terms
543     c of actual snow thickness. As with ice, multiplying
544     c with area converts to mean snow thickness.
545     S_hsnow(I,J) = AREA(I,J,2,bi,bj)* (
546     & SnowAccumulationRateOverIce(I,J) -
547     & SnowMeltRateFromSurface(I,J) )
548     ENDDO
549     ENDDO
550    
551     #ifdef ALLOW_AUTODIFF_TAMC
552     CADJ STORE area(:,:,:,bi,bj) = comlev1_bibj,
553     CADJ & key = iicekey, byte = isbyte
554     CADJ STORE heff(:,:,:,bi,bj) = comlev1_bibj,
555     CADJ & key = iicekey, byte = isbyte
556     CADJ STORE hsnow(:,:,bi,bj) = comlev1_bibj,
557     CADJ & key = iicekey, byte = isbyte
558     CADJ STORE S_h(:,:) = comlev1_bibj,
559     CADJ & key = iicekey, byte = isbyte
560     CADJ STORE S_hsnow(:,:) = comlev1_bibj,
561     CADJ & key = iicekey, byte = isbyte
562     #endif /* ALLOW_AUTODIFF_TAMC */
563    
564     DO J=1,sNy
565     DO I=1,sNx
566     S_a(I,J) = 0.0
567     C IF THE OPEN WATER GROWTH RATE, F_ao, IS POSITIVE
568     C THEN EXTEND ICE AREAL COVER, S_a > 0
569     IF (F_ao(I,J) .GT. 0) THEN
570     S_a(I,J) = S_a(I,J) + (ONE - AREA(I,J,2,bi,bj))*
571     & IceGrowthRateOpenWater(I,J)/HO
572     ELSE
573     S_a(I,J) = S_a(I,J) + 0.0
574     ENDIF
575    
576    
577     C REDUCE THE ICE COVER IF ICE IS PRESENT
578     IF ( (S_h(I,J) .LT. 0.) .AND.
579     & (AREA(I,J,2,bi,bj).GT. 0.) .AND.
580     & (HEFF(I,J,2,bi,bj).NE. 0.) ) THEN
581    
582     S_a(I,J) = S_a(I,J)
583     & + AREA(I,J,2,bi,bj)/(2.0*HEFF(I,J,2,bi,bj))*
584     & IceGrowthRateOpenWater(I,J)*
585     & (1-AREA(I,J,2,bi,bj))
586     ELSE
587     S_a(I,J) = S_a(I,J) + 0.0
588     ENDIF
589    
590     C REDUCE THE ICE COVER IF ICE IS PRESENT
591     IF ( (IceGrowthRateMixedLayer(I,J) .LE. 0.) .AND.
592     & (AREA(I,J,2,bi,bj).GT. 0.) .AND.
593     & (HEFF(I,J,2,bi,bj).NE. 0.) ) THEN
594    
595     S_a(I,J) = S_a(I,J)
596     & + AREA(I,J,2,bi,bj)/(2.0*HEFF(I,J,2,bi,bj))*
597     & IceGrowthRateMixedLayer(I,J)
598    
599     ELSE
600     S_a(I,J) = S_a(I,J) + 0.0
601     ENDIF
602    
603     C REDUCE THE ICE COVER IF ICE IS PRESENT
604     IF ( (NetExistingIceGrowthRate(I,J) .LE. 0.) .AND.
605     & (AREA(I,J,2,bi,bj).GT. 0.) .AND.
606     & (HEFF(I,J,2,bi,bj).NE. 0.) ) THEN
607    
608     S_a(I,J) = S_a(I,J)
609     & + AREA(I,J,2,bi,bj)/(2.0*HEFF(I,J,2,bi,bj))*
610     & NetExistingIceGrowthRate(I,J)*AREA(I,J,2,bi,bj)
611    
612     ELSE
613     S_a(I,J) = S_a(I,J) + 0.0
614     ENDIF
615    
616     ENDDO
617     ENDDO
618    
619    
620     #ifdef ALLOW_AUTODIFF_TAMC
621     CADJ STORE heff(:,:,:,bi,bj) = comlev1_bibj,
622     CADJ & key = iicekey, byte = isbyte
623     CADJ STORE hsnow(:,:,bi,bj) = comlev1_bibj,
624     CADJ & key = iicekey, byte = isbyte
625     CADJ STORE S_a(:,:) = comlev1_bibj,
626     CADJ & key = iicekey, byte = isbyte
627     CADJ STORE S_h(:,:) = comlev1_bibj,
628     CADJ & key = iicekey, byte = isbyte
629     CADJ STORE f_ao(:,:) = comlev1_bibj,
630     CADJ & key = iicekey, byte = isbyte
631     CADJ STORE qswi(:,:) = comlev1_bibj,
632     CADJ & key = iicekey, byte = isbyte
633     CADJ STORE qswo(:,:) = comlev1_bibj,
634     CADJ & key = iicekey, byte = isbyte
635     CADJ STORE area(:,:,:,bi,bj) = comlev1_bibj,
636     CADJ & key = iicekey, byte = isbyte
637     #endif
638    
639     C ACTUALLY CHANGE THE AREA AND THICKNESS
640     DO J=1,sNy
641     DO I=1,sNx
642     AREA(I,J,1,bi,bj) = AREA(I,J,2,bi,bj) +
643     & SEAICE_deltaTtherm * S_a(I,J)
644     C SET LIMIT ON AREA
645     AREA(I,J,1,bi,bj) = MIN(1.,AREA(I,J,1,bi,bj))
646     AREA(I,J,1,bi,bj) = MAX(0.,AREA(I,J,1,bi,bj))
647    
648     ENDDO
649     ENDDO
650    
651     #ifdef ALLOW_AUTODIFF_TAMC
652     CADJ STORE area(:,:,:,bi,bj) = comlev1_bibj,
653     CADJ & key = iicekey, byte = isbyte
654     #endif
655     DO J=1,sNy
656     DO I=1,sNx
657     HEFF(I,J,1,bi,bj) = HEFF(I,J,2,bi,bj) +
658     & SEAICE_deltaTTherm * S_h(I,J)
659     HEFF(I,J,1,bi,bj) = MAX(0.0, HEFF(I,J,1,bi,bj))
660    
661     HSNOW(I,J,bi,bj) = HSNOW(I,J,bi,bj) +
662     & SEAICE_deltaTTherm * S_hsnow(I,J)
663     HSNOW(I,J,bi,bj) = MAX(0.0, HSNOW(I,J,bi,bj))
664    
665     IF (AREA(I,J,1,bi,bj) .GT. 0.0) THEN
666     HICE_ACTUAL(I,J) =
667     & HEFF(I,J,1,bi,bj)/AREA(I,J,1,bi,bj)
668     HSNOW_ACTUAL(I,J) =
669     & HSNOW(I,J,bi,bj)/AREA(I,J,1,bi,bj)
670     ELSE
671     HICE_ACTUAL(I,J) = 0.0
672     HSNOW_ACTUAL(I,J) = 0.0
673     ENDIF
674    
675     ENDDO
676     ENDDO
677    
678     #ifdef ALLOW_AUTODIFF_TAMC
679     CADJ STORE area(:,:,:,bi,bj) = comlev1_bibj,
680     CADJ & key = iicekey, byte = isbyte
681     CADJ STORE heff(:,:,:,bi,bj) = comlev1_bibj,
682     CADJ & key = iicekey, byte = isbyte
683     #endif /* ALLOW_AUTODIFF_TAMC */
684    
685     c constrain area is no thickness and vice versa.
686     DO J=1,sNy
687     DO I=1,sNx
688     IF (HEFF(I,J,1,bi,bj) .LE. 0.0 .OR.
689     & AREA(I,J,1,bi,bj) .LE. 0.0) THEN
690    
691     AREA(I,J,1,bi,bj) = 0.0
692     HEFF(I,J,1,bi,bj) = 0.0
693     HICE_ACTUAL(I,J) = 0.0
694     HSNOW(I,J,bi,bj) = 0.0
695     HSNOW_ACTUAL(I,J) = 0.0
696     ENDIF
697     ENDDO
698     ENDDO
699    
700     #ifdef ALLOW_AUTODIFF_TAMC
701     CADJ STORE area(:,:,:,bi,bj) = comlev1_bibj,
702     CADJ & key = iicekey, byte = isbyte
703     CADJ STORE heff(:,:,:,bi,bj) = comlev1_bibj,
704     CADJ & key = iicekey, byte = isbyte
705     #endif /* ALLOW_AUTODIFF_TAMC */
706    
707     DO J=1,sNy
708     DO I=1,sNx
709    
710     c The amount of new mean thickness we expect to grow
711     ExpectedIceVolumeChange(I,J) = S_h(I,J) *
712     & SEAICE_deltaTtherm
713    
714     ExpectedSnowVolumeChange(I,J) = S_hsnow(I,J)*
715     & SEAICE_deltaTtherm
716    
717     c THE EFFECTIVE SHORTWAVE HEATING RATE
718     QSW(I,J,bi,bj) =
719     & QSWI(I,J) * ( AREA(I,J,2,bi,bj)) +
720     & QSWO(I,J) * (1. - AREA(I,J,2,bi,bj))
721    
722     ActualNewTotalVolumeChange(I,J) =
723     & HEFF(I,J,1,bi,bj) - HEFF(I,J,2,bi,bj)
724    
725     ActualNewTotalSnowMelt(I,J) =
726     & HSNOW_ORIG(I,J) +
727     & SnowAccumulationOverIce(I,J) -
728     & HSNOW(I,J,bi,bj)
729    
730     c The latent heat of fusion of the new ice
731     EnergyInNewTotalIceVolume(I,J) =
732     & ActualNewTotalVolumeChange(I,J)/QI
733    
734     c THE EFFECTIVE SHORTWAVE HEATING RATE THE LIQUID OCEAN WILL FEEL
735     NetEnergyFluxOutOfSystem(I,J) = SEAICE_deltaTtherm *
736     & (AREA(I,J,2,bi,bj) *
737     & (F_ia_net(I,J) + F_io_net(I,J) + QSWI(I,J))
738     & + (1.0 - AREA(I,J,2,bi,bj)) * F_ao(I,J))
739    
740     c THE QUANTITY OF HEAT WHICH IS THE RESIDUAL TO THE QUANTITY OF
741     c ML temperature. If the net energy flux is exactly balanced by the
742     c latent energy of fusion in the new ice created then we won't
743     c change the ML temperature at all.
744    
745     ResidualHeatOutOfSystem(I,J) =
746     & NetEnergyFluxOutOfSystem(I,J) -
747     & EnergyInNewTotalIceVolume(I,J)
748    
749     c Like snow melt, if there is melting, this quantity is positive.
750     FreshwaterContribFromIce(I,J) =
751     & -ActualNewTotalVolumeChange(I,J)*RHOI/RHOFW
752    
753     c The freshwater contribution from snow comes only in the form of melt
754     c unlike ice, which takes freshwater upon growth and yields freshwater
755     c upon melt. Thus, the actual new total snow melt must be determined.
756     FreshwaterContribFromSnowMelt(I,J) =
757     & ActualNewTotalSnowMelt(I,J)*RHOSN/RHOFW
758    
759     c This seems to be in m/s, original time level 2 for area
760     EmPmR(I,J,bi,bj) = maskC(I,J,kSurface,bi,bj)*(
761     & ( EVAP(I,J,bi,bj)-PRECIP(I,J,bi,bj) )
762     & * ( ONE - AREA(I,J,2,bi,bj) )
763     & - PrecipRateOverIceSurfaceToSea(I,J)*
764     & AREA(I,J,2,bi,bj)
765     & - RUNOFF(I,J,bi,bj)
766     & - (FreshwaterContribFromIce(I,J) +
767     & FreshwaterContribFromSnowMelt(I,J))/
768     & SEAICE_deltaTtherm)
769    
770     C THE SHORTWAVE ENERGY FLUX ABSORBED IN THE SURFACE LAYER
771     QSW_absorb_in_ML(I,J) = QSW(I,J,bi,bj)*
772     & (1.0 - SWFRACB)
773    
774     C THE SHORTWAVE ENERGY FLUX PENETRATING BELOW THE SURFACE LAYER
775     QSW_absorb_below_ML(I,J) =
776     & QSW(I,J,bi,bj) - QSW_absorb_in_ML(I,J);
777    
778     PredictedTemperatureChangeFromQSW(I,J) =
779     & - QSW_absorb_in_ML(I,J) * FLUX_TO_DELTA_TEMP
780    
781     PredictedTemperatureChangeFromOA_MQNET(I,J) =
782     & -(QNET(I,J,bi,bj) - QSWO(I,J))*(1. -AREA(I,J,2,bi,bj))
783     & * FLUX_TO_DELTA_TEMP
784    
785     PredictedTemperatureChangeFromF_IO_NET(I,J) =
786     & -F_io_net(I,J)*AREA(I,J,2,bi,bj)*FLUX_TO_DELTA_TEMP
787    
788     PredictedTemperatureChangeFromF_IA_NET(I,J) =
789     & -F_ia_net(I,J)*AREA(I,J,2,bi,bj)*FLUX_TO_DELTA_TEMP
790    
791     PredictedTemperatureChangeFromNewIceVol(I,J) =
792     & EnergyInNewTotalIceVolume(I,J)*ENERGY_TO_DELTA_TEMP
793    
794     PredictedTemperatureChange(I,J) =
795     & PredictedTemperatureChangeFromQSW(I,J) +
796     & PredictedTemperatureChangeFromOA_MQNET(I,J) +
797     & PredictedTemperatureChangeFromF_IO_NET(I,J) +
798     & PredictedTemperatureChangeFromF_IA_NET(I,J) +
799     & PredictedTemperatureChangeFromNewIceVol(I,J)
800    
801     C NOW FORMULATE QNET, which time LEVEL, ORIG 2.
802     C THIS QNET WILL DETERMINE THE TEMPERATURE CHANGE OF THE MIXED LAYER
803     C QNET IS A DEPTH AVERAGED HEAT FLUX FOR THE OCEAN COLUMN
804     C BECAUSE OF THE
805     QNET(I,J,bi,bj) =
806     & ResidualHeatOutOfSystem(I,J) / SEAICE_deltaTtherm
807    
808     ENDDO
809     ENDDO
810    
811    
812     #ifdef ALLOW_AUTODIFF_TAMC
813     CADJ STORE area(:,:,:,bi,bj) = comlev1_bibj,
814     CADJ & key = iicekey, byte = isbyte
815     CADJ STORE hsnow(:,:,bi,bj) = comlev1_bibj,
816     CADJ & key = iicekey, byte = isbyte
817     CADJ STORE heff(:,:,:,bi,bj) = comlev1_bibj,
818     CADJ & key = iicekey, byte = isbyte
819     #endif /* ALLOW_AUTODIFF_TAMC */
820    
821     DO J=1,sNy
822     DO I=1,sNx
823     AREA(I,J,1,bi,bj) = AREA(I,J,1,bi,bj)*HEFFM(I,J,bi,bj)
824     HEFF(I,J,1,bi,bj) = HEFF(I,J,1,bi,bj)*HEFFM(I,J,bi,bj)
825     HSNOW(I,J,bi,bj) = HSNOW(I,J,bi,bj)*HEFFM(I,J,bi,bj)
826     ENDDO
827     ENDDO
828    
829     #ifdef ALLOW_AUTODIFF_TAMC
830     CADJ STORE area(:,:,:,bi,bj) = comlev1_bibj,
831     CADJ & key = iicekey, byte = isbyte
832     CADJ STORE hsnow(:,:,bi,bj) = comlev1_bibj,
833     CADJ & key = iicekey, byte = isbyte
834     CADJ STORE heff(:,:,:,bi,bj) = comlev1_bibj,
835     CADJ & key = iicekey, byte = isbyte
836     #endif /* ALLOW_AUTODIFF_TAMC */
837    
838    
839     #ifdef ALLOW_SEAICE_FLOODING
840     DO J = 1,sNy
841     DO I = 1,sNx
842     EnergyToMeltSnowAndIce(I,J) =
843     & HEFF(I,J,1,bi,bj)/QI +
844     & HSNOW(I,J,bi,bj)/QS
845    
846     deltaHS = FL_C2*( HSNOW_ACTUAL(I,J) -
847     & HICE_ACTUAL(I,J)*FL_C3 )
848    
849     IF (deltaHS .GT. 0.0) THEN
850     deltaHI = FL_C4*deltaHS
851    
852     HICE_ACTUAL(I,J) = HICE_ACTUAL(I,J)
853     & + deltaHI
854    
855     HSNOW_ACTUAL(I,J)= HSNOW_ACTUAL(I,J)
856     & - deltaHS
857    
858     HEFF(I,J,1,bi,bj)= HICE_ACTUAL(I,J) *
859     & AREA(I,J,1,bi,bj)
860    
861     HSNOW(I,J,bi,bj) = HSNOW_ACTUAL(I,J)*
862     & AREA(I,J,1,bi,bj)
863    
864     EnergyToMeltSnowAndIce2(I,J) =
865     & HEFF(I,J,1,bi,bj)/QI +
866     & HSNOW(I,J,bi,bj)/QS
867    
868     #ifdef SEAICE_DEBUG
869     print *,'Energy to melt snow+ice: pre,post,delta',
870     & EnergyToMeltSnowAndIce(I,J),
871     & EnergyToMeltSnowAndIce2(I,J),
872     & EnergyToMeltSnowAndIce(I,J) -
873     & EnergyToMeltSnowAndIce2(I,J)
874     #endif
875     ENDIF
876     ENDDO
877     ENDDO
878     #endif
879    
880     #ifdef ALLOW_AUTODIFF_TAMC
881     CADJ STORE area(:,:,:,bi,bj) = comlev1_bibj,
882     CADJ & key = iicekey, byte = isbyte
883     CADJ STORE hsnow(:,:,bi,bj) = comlev1_bibj,
884     CADJ & key = iicekey, byte = isbyte
885     CADJ STORE heff(:,:,:,bi,bj) = comlev1_bibj,
886     CADJ & key = iicekey, byte = isbyte
887     #endif /* ALLOW_AUTODIFF_TAMC */
888    
889    
890     #ifdef ATMOSPHERIC_LOADING
891     IF ( useRealFreshWaterFlux ) THEN
892     DO J=1,sNy
893     DO I=1,sNx
894     sIceLoad(i,j,bi,bj) = HEFF(I,J,1,bi,bj)*
895     & SEAICE_rhoIce + HSNOW(I,J,bi,bj)* 330. _d 0
896     ENDDO
897     ENDDO
898     ENDIF
899     #endif
900    
901     #ifdef SEAICE_DEBUG
902     DO j=1-OLy,sNy+OLy
903     DO i=1-OLx,sNx+OLx
904     IF ((i .EQ. SEAICE_debugPointX .and.
905     & j .EQ. SEAICE_debugPointY)) THEN
906    
907     print *,'ifsig: myTime,myIter:',myTime,myIter
908    
909     print '(A,2i4,2(1x,1P3E15.4))',
910     & 'ifice i j -------------- ',i,j
911    
912     print '(A,2i4,2(1x,1P3E15.4))',
913     & 'ifice i j IGR(ML OW ICE) ',i,j,
914     & IceGrowthRateMixedLayer(i,j),
915     & IceGrowthRateOpenWater(i,j),
916     & NetExistingIceGrowthRate(i,j),
917     & SEAICE_deltaTtherm
918    
919     print '(A,2i4,2(1x,1P3E15.4))',
920     & 'ifice i j F(mi ao) ',
921     & i,j,F_mi(i,j), F_ao(i,j)
922    
923     print '(A,2i4,2(1x,1P3E15.4))',
924     & 'ifice i j Fi(a,ant2/1 ont)',
925     & i,j,F_ia(i,j),
926     & F_ia_net_before_snow(i,j),
927     & F_ia_net(i,j),
928     & F_io_net(i,j)
929    
930     print '(A,2i4,2(1x,1P3E15.4))',
931     & 'ifice i j AREA2/1 HEFF2/1 ',i,j,
932     & AREA(i,j,2,bi,bj),
933     & AREA(i,j,1,bi,bj),
934     & HEFF(i,j,2,bi,bj),
935     & HEFF(i,j,1,bi,bj)
936    
937     print '(A,2i4,2(1x,1P3E15.4))',
938     & 'ifice i j HSNOW2/1 TMX TBC',i,j,
939     & HSNOW_ORIG(I,J),
940     & HSNOW(I,J,bi,bj),
941     & TMIX(i,j,bi,bj)-273.15,
942     & TBC
943    
944     print '(A,2i4,2(1x,1P3E15.4))',
945     & 'ifice i j TI ATP LWD ',i,j,
946     & TICE(i,j,bi,bj)-273.15,
947     & ATEMP(i,j,bi,bj)-273.15,
948     & LWDOWN(i,j,bi,bj)
949    
950    
951     print '(A,2i4,2(1x,1P3E15.4))',
952     & 'ifice i j S_a S_h S_hsnow ',i,j,
953     & S_a(i,j),
954     & S_h(i,j),
955     & S_hsnow(i,j)
956    
957     print '(A,2i4,2(1x,1P3E15.4))',
958     & 'ifice i j IVC(E A ENIN) ',i,j,
959     & ExpectedIceVolumeChange(i,j),
960     & ActualNewTotalVolumeChange(i,j),
961     & EnergyInNewTotalIceVolume(i,j)
962    
963     print '(A,2i4,2(1x,1P3E15.4))',
964     & 'ifice i j EF(NOS RE) QNET ',i,j,
965     & NetEnergyFluxOutOfSystem(i,j),
966     & ResidualHeatOutOfSystem(i,j),
967     & QNET(I,J,bi,bj)
968    
969     print '(A,2i4,3(1x,1P3E15.4))',
970     & 'ifice i j QSW QSWO QSWI ',i,j,
971     & QSW(i,j,bi,bj),
972     & QSWO(i,j),
973     & QSWI(i,j)
974    
975     print '(A,2i4,3(1x,1P3E15.4))',
976     & 'ifice i j SW(BML IML SW) ',i,j,
977     & QSW_absorb_below_ML(i,j),
978     & QSW_absorb_in_ML(i,j),
979     & SWFRACB
980    
981     print '(A,2i4,3(1x,1P3E15.4))',
982     & 'ifice i j ptc(to, qsw, oa)',i,j,
983     & PredictedTemperatureChange(i,j),
984     & PredictedTemperatureChangeFromQSW (i,j),
985     & PredictedTemperatureChangeFromOA_MQNET(i,j)
986    
987    
988     print '(A,2i4,3(1x,1P3E15.4))',
989     & 'ifice i j ptc(fion,ian,ia)',i,j,
990     & PredictedTemperatureChangeFromF_IO_NET(i,j),
991     & PredictedTemperatureChangeFromF_IA_NET(i,j),
992     & PredictedTemperatureChangeFromFIA(i,j)
993    
994     print '(A,2i4,3(1x,1P3E15.4))',
995     & 'ifice i j ptc(niv) ',i,j,
996     & PredictedTemperatureChangeFromNewIceVol(i,j)
997    
998    
999     print '(A,2i4,3(1x,1P3E15.4))',
1000     & 'ifice i j EmPmR EVP PRE RU',i,j,
1001     & EmPmR(I,J,bi,bj),
1002     & EVAP(I,J,bi,bj),
1003     & PRECIP(I,J,bi,bj),
1004     & RUNOFF(I,J,bi,bj)
1005    
1006     print '(A,2i4,3(1x,1P3E15.4))',
1007     & 'ifice i j PRROIS,SAOI(R .)',i,j,
1008     & PrecipRateOverIceSurfaceToSea(I,J),
1009     & SnowAccumulationRateOverIce(I,J),
1010     & SnowAccumulationOverIce(I,J)
1011    
1012     print '(A,2i4,4(1x,1P3E15.4))',
1013     & 'ifice i j SM(PM PMR . .R) ',i,j,
1014     & PotSnowMeltFromSurf(I,J),
1015     & PotSnowMeltRateFromSurf(I,J),
1016     & SnowMeltFromSurface(I,J),
1017     & SnowMeltRateFromSurface(I,J)
1018    
1019     print '(A,2i4,4(1x,1P3E15.4))',
1020     & 'ifice i j TotSnwMlt ExSnVC',i,j,
1021     & ActualNewTotalSnowMelt(I,J),
1022     & ExpectedSnowVolumeChange(I,J)
1023    
1024    
1025     print '(A,2i4,4(1x,1P3E15.4))',
1026     & 'ifice i j fw(CFICE, CFSM) ',i,j,
1027     & FreshwaterContribFromIce(I,J),
1028     & FreshwaterContribFromSnowMelt(I,J)
1029    
1030     print '(A,2i4,2(1x,1P3E15.4))',
1031     & 'ifice i j -------------- ',i,j
1032    
1033     ENDIF
1034     ENDDO
1035     ENDDO
1036     #endif /* SEAICE_DEBUG */
1037    
1038    
1039     C BI,BJ'S
1040     ENDDO
1041     ENDDO
1042    
1043     RETURN
1044     END

  ViewVC Help
Powered by ViewVC 1.1.22