| 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 |