Merge from cesm_cice branch into trunk, r992:1074.

BUG FIX:  changed shortwave calculation to depend on the net shortwave sum 
instead of cosine of the zenith angle.  This changes the answers.

In addition to the different shortwave calculation, albedo output in history
is different.

Did not adopt CESM's naming convention for history files (commented out), 
since that potentially breaks the post-processing software of other users.


M       mpi/ice_boundary.F90
M       mpi/ice_global_reductions.F90
 - remove save attribute

MM      source/ice_init.F90
 - write nslyr to diagnostics file

M       source/ice_domain.F90
 - add 1 to latitude distribution

M       source/ice_calendar.F90
 - add CCSMCOUPLED wrapper

M       source/ice_restart_shared.F90
 - lengthen character length for restart_file

M       source/ice_itd.F90
 - add CCSMCOUPLED wrapper
 - remove (:,:) from atrcr in compute_tracers call

M       source/ice_history_pond.F90
 - use flag f_apond_ai for variable apond_ai

M       source/ice_shortwave.F90
 - add new variables snowfracn
 - compute netsw instead of using coszen for history averaging

M       source/ice_aerosol.F90
 - remove scaling by aice from dhs, dhi terms

M       source/ice_diagnostics.F90
 - use work1 instead of coszen

M       source/ice_step_mod.F90
 - add snowfracn

M       source/ice_flux.F90
 - add CCSMCOUPLED values for l_winter, l_spring
 - add snowfrac

MM      source/ice_therm_bl99.F90
 - add CCSMCOUPLED values for frac, dTemp

M       source/ice_history_shared.F90
 - change VGRD to VGRDi for 3D variables
 - add snowfrac, snowfracn, fswup, alvdr_ai etc
 - cosmetic

M       source/ice_history.F90
M       input_templates/gx1/ice_in
M       input_templates/gx3/ice_in
M       input_templates/col/ice_in
 - add snowfrac, snowfracn, fswup, alvdr_ai etc

M       drivers/cice/CICE_RunMod.F90
M       drivers/cice/CICE_RunMod.F90_debug
M       drivers/hadgem3/CICE_RunMod.F90
 - add snowfrac, snowfracn
 - compute netsw instead of using coszen for history averaging

M       drivers/cesm/CICE_RunMod.F90
 - add snowfrac, snowfracn, l_mpond_fresh
 - skip dynamics call if kdyn=0
 - compute netsw instead of using coszen for history averaging

M       drivers/cesm/ice_comp_esmf.F90
 - change idate

M       drivers/cesm/ice_comp_mct.F90
 - change idate
 - initialize data

M       io_pio/ice_history_write.F90
 - add histfreqn, time_period_freq, iodesc4ds,
 - change nzlyr to nzilyr
 - add nzslyr
 - add 4d fields (categories, vertical ice)

cice/drivers/cesm/CICE_RunMod.F90

@@ -197,7 +197,7 @@ subroutine ice_step
       ! dynamics, transport, ridging
       !-----------------------------------------------------------------
 
-         if (.not.prescribed_ice) then
+         if (.not.prescribed_ice .and. kdyn>0) then
          do k = 1, ndtd
             call step_dynamics (dt_dyn, ndtd)
          enddo
@@ -281,17 +281,18 @@ subroutine coupling_prep (iblk)
       use ice_constants, only: c0, c1, puny, rhofresh
       use ice_domain_size, only: ncat
       use ice_flux, only: alvdf, alidf, alvdr, alidr, albice, albsno, &
-          albpnd, albcnt, apeff_ai, coszen, fpond, fresh, &
+          albpnd, albcnt, apeff_ai, coszen, fpond, fresh, l_mpond_fresh, &
           alvdf_ai, alidf_ai, alvdr_ai, alidr_ai, fhocn_ai, &
           fresh_ai, fsalt_ai, fsalt, &
           fswthru_ai, fhocn, fswthru, scale_factor, &
           swvdr, swidr, swvdf, swidf, Tf, Tair, Qa, strairxT, strairyt, &
           fsens, flat, fswabs, flwout, evap, Tref, Qref, Uref, faero_ocn, &
-          fsurfn_f, flatn_f, scale_fluxes, frzmlt_init, frzmlt, wind
+          fsurfn_f, flatn_f, scale_fluxes, frzmlt_init, frzmlt, wind, &
+          snowfrac
       use ice_grid, only: tmask
       use ice_ocean, only: oceanmixed_ice, ocean_mixed_layer
       use ice_shortwave, only: alvdfn, alidfn, alvdrn, alidrn, &
-                               albicen, albsnon, albpndn, apeffn
+                               albicen, albsnon, albpndn, apeffn, snowfracn
       use ice_state, only: aicen, aice, aice_init, nbtrcr
       use ice_therm_shared, only: calc_Tsfc
       use ice_timers, only: timer_couple, ice_timer_start, ice_timer_stop
@@ -307,8 +308,10 @@ subroutine coupling_prep (iblk)
          i,j         , & ! horizontal indices
          k               ! tracer index
 
-      real (kind=dbl_kind) :: cszn ! counter for history averaging
-
+      real (kind=dbl_kind) :: &
+         cszn        , & ! counter for history averaging
+         netsw           ! flag for shortwave radiation presence
+
       !-----------------------------------------------------------------
       ! Save current value of frzmlt for diagnostics.
       ! Update mixed layer with heat and radiation from ice.
@@ -340,10 +343,12 @@ subroutine coupling_prep (iblk)
             albsno(i,j,iblk) = c0
             albpnd(i,j,iblk) = c0
             apeff_ai(i,j,iblk) = c0
+            snowfrac(i,j,iblk) = c0
 
             ! for history averaging
             cszn = c0
-            if (coszen(i,j,iblk) > puny) cszn = c1
+            netsw = swvdr(i,j,iblk)+swidr(i,j,iblk)+swvdf(i,j,iblk)+swidf(i,j,iblk)
+            if (netsw > puny) cszn = c1
             do n = 1, nstreams
                albcnt(i,j,iblk,n) = albcnt(i,j,iblk,n) + cszn
             enddo
@@ -372,6 +377,8 @@ subroutine coupling_prep (iblk)
 
             apeff_ai(i,j,iblk) = apeff_ai(i,j,iblk) &       ! for history
                + apeffn(i,j,n,iblk)*aicen(i,j,n,iblk)
+            snowfrac(i,j,iblk) = snowfrac(i,j,iblk) &       ! for history
+               + snowfracn(i,j,n,iblk)*aicen(i,j,n,iblk)
          enddo
          enddo
          enddo
@@ -383,8 +390,10 @@ subroutine coupling_prep (iblk)
       ! reduce fresh by fpond for coupling
       !-----------------------------------------------------------------
 
-            fpond(i,j,iblk) = fpond(i,j,iblk) * rhofresh/dt
-            fresh(i,j,iblk) = fresh(i,j,iblk) - fpond(i,j,iblk)
+            if (l_mpond_fresh) then
+               fpond(i,j,iblk) = fpond(i,j,iblk) * rhofresh/dt
+               fresh(i,j,iblk) = fresh(i,j,iblk) - fpond(i,j,iblk)
+            endif
 
       !----------------------------------------------------------------
       ! Store grid box mean albedos and fluxes before scaling by aice

cice/drivers/cesm/ice_comp_esmf.F90

@@ -48,7 +48,7 @@ module ice_comp_esmf
   use ice_constants,   only : c0, c1, spval_dbl, rad_to_deg, radius, secday
   use ice_communicate, only : my_task, master_task, MPI_COMM_ICE
   use ice_calendar,    only : istep, istep1, force_restart_now, write_ic,&
-                              idate, mday, time, month, daycal,          &
+                              idate, idate0, mday, time, month, daycal,  &
 		              sec, dt, dt_dyn, calendar,                 &
                               calendar_type, nextsw_cday, days_per_year, &
                               nyr, new_year, time2sec, year_init
@@ -354,7 +354,9 @@ subroutine ice_init_esmf(comp, import_state, export_state, EClock, rc)
           write(nu_diag,*) trim(subname)//' resetting idate to match sync clock'
        end if
 
-       idate = curr_ymd - (year_init*10000)      ! adjust for year_init
+       idate0 = curr_ymd
+       idate = curr_ymd
+
        if (idate < 0) then
           write(nu_diag,*) trim(subname),' ERROR curr_ymd,year_init =',curr_ymd,year_init
           write(nu_diag,*) trim(subname),' ERROR idate lt zero',idate

cice/drivers/cesm/ice_comp_mct.F90

@@ -45,7 +45,7 @@ module ice_comp_mct
   use ice_constants,   only : c0, c1, spval_dbl, rad_to_deg, radius, secday
   use ice_communicate, only : my_task, master_task, MPI_COMM_ICE
   use ice_calendar,    only : istep, istep1, force_restart_now, write_ic,&
-                              idate, mday, time, month, daycal,          &
+                              idate, idate0, mday, time, month, daycal,  &
 		              sec, dt, dt_dyn, calendar,                 &
                               calendar_type, nextsw_cday, days_per_year, &
                               nyr, new_year, time2sec, year_init
@@ -291,7 +291,9 @@ subroutine ice_init_mct( EClock, cdata_i, x2i_i, i2x_i, NLFilename )
                trim(subname),' resetting idate to match sync clock'
        end if
 
-       idate = curr_ymd - (year_init*10000)      ! adjust for year_init
+       idate0 = curr_ymd
+       idate = curr_ymd
+
        if (idate < 0) then
           write(nu_diag,*) trim(subname),' ERROR curr_ymd,year_init =',curr_ymd,year_init
           write(nu_diag,*) trim(subname),' ERROR idate lt zero',idate
@@ -789,6 +791,7 @@ subroutine ice_domain_mct( lsize, gsMap_i, dom_i )
     allocate(work_dom(lsize)) 
     work_dom(:) = 0.0_dbl_kind
 
+    data(:) = -9999.0_R8 
     n=0
     do iblk = 1, nblocks
        this_block = get_block(blocks_ice(iblk),iblk)         
@@ -806,6 +809,7 @@ subroutine ice_domain_mct( lsize, gsMap_i, dom_i )
     enddo       !iblk
     call mct_gGrid_importRattr(dom_i,"lon",data,lsize) 
 
+    data(:) = -9999.0_R8 
     n=0
     do iblk = 1, nblocks
        this_block = get_block(blocks_ice(iblk),iblk)         
@@ -823,6 +827,7 @@ subroutine ice_domain_mct( lsize, gsMap_i, dom_i )
     enddo      !iblk
     call mct_gGrid_importRattr(dom_i,"lat",data,lsize) 
 
+    data(:) = -9999.0_R8 
     n=0
     do iblk = 1, nblocks
        this_block = get_block(blocks_ice(iblk),iblk)         
@@ -840,6 +845,7 @@ subroutine ice_domain_mct( lsize, gsMap_i, dom_i )
     enddo      !iblk
     call mct_gGrid_importRattr(dom_i,"area",data,lsize) 
 
+    data(:) = 0.0_R8
     n=0
     do iblk = 1, nblocks
        this_block = get_block(blocks_ice(iblk),iblk)         
@@ -857,6 +863,7 @@ subroutine ice_domain_mct( lsize, gsMap_i, dom_i )
     enddo      !iblk
     call mct_gGrid_importRattr(dom_i,"mask",data,lsize) 
 
+    data(:) = 0.0_R8
     n=0
     do iblk = 1, nblocks
        this_block = get_block(blocks_ice(iblk),iblk)         

cice/drivers/cice/CICE_RunMod.F90

@@ -279,11 +279,11 @@ subroutine coupling_prep (iblk)
           fswthru_ai, fhocn, fswthru, scale_factor, &
           swvdr, swidr, swvdf, swidf, Tf, Tair, Qa, strairxT, strairyt, &
           fsens, flat, fswabs, flwout, evap, Tref, Qref, faero_ocn, &
-          fsurfn_f, flatn_f, scale_fluxes, frzmlt_init, frzmlt
+          fsurfn_f, flatn_f, scale_fluxes, frzmlt_init, frzmlt, snowfrac
       use ice_grid, only: tmask
       use ice_ocean, only: oceanmixed_ice, ocean_mixed_layer
       use ice_shortwave, only: alvdfn, alidfn, alvdrn, alidrn, &
-                               albicen, albsnon, albpndn, apeffn
+                               albicen, albsnon, albpndn, apeffn, snowfracn
       use ice_state, only: aicen, aice, aice_init, nbtrcr
       use ice_therm_shared, only: calc_Tsfc
       use ice_timers, only: timer_couple, ice_timer_start, ice_timer_stop
@@ -299,7 +299,9 @@ subroutine coupling_prep (iblk)
          i,j         , & ! horizontal indices
          k               ! tracer index
 
-      real (kind=dbl_kind) :: cszn ! counter for history averaging
+      real (kind=dbl_kind) :: &
+         cszn        , & ! counter for history averaging
+         netsw           ! flag for shortwave radiation presence
 
       !-----------------------------------------------------------------
       ! Save current value of frzmlt for diagnostics.
@@ -332,10 +334,12 @@ subroutine coupling_prep (iblk)
             albsno(i,j,iblk) = c0
             albpnd(i,j,iblk) = c0
             apeff_ai(i,j,iblk) = c0
+            snowfrac(i,j,iblk) = c0
 
             ! for history averaging
             cszn = c0
-            if (coszen(i,j,iblk) > puny) cszn = c1
+            netsw = swvdr(i,j,iblk)+swidr(i,j,iblk)+swvdf(i,j,iblk)+swidf(i,j,iblk)
+            if (netsw > puny) cszn = c1
             do n = 1, nstreams
                albcnt(i,j,iblk,n) = albcnt(i,j,iblk,n) + cszn
             enddo
@@ -364,6 +368,8 @@ subroutine coupling_prep (iblk)
 
             apeff_ai(i,j,iblk) = apeff_ai(i,j,iblk) &       ! for history
                + apeffn(i,j,n,iblk)*aicen(i,j,n,iblk)
+            snowfrac(i,j,iblk) = snowfrac(i,j,iblk) &       ! for history
+               + snowfracn(i,j,n,iblk)*aicen(i,j,n,iblk)
          enddo
          enddo
          enddo

cice/drivers/cice/CICE_RunMod.F90_debug

@@ -312,11 +312,11 @@
           fswthru_ai, fhocn, fswthru, scale_factor, &
           swvdr, swidr, swvdf, swidf, Tf, Tair, Qa, strairxT, strairyt, &
           fsens, flat, fswabs, flwout, evap, Tref, Qref, faero_ocn, &
-          fsurfn_f, flatn_f, scale_fluxes, frzmlt_init, frzmlt
+          fsurfn_f, flatn_f, scale_fluxes, frzmlt_init, frzmlt, snowfrac
       use ice_grid, only: tmask
       use ice_ocean, only: oceanmixed_ice, ocean_mixed_layer
       use ice_shortwave, only: alvdfn, alidfn, alvdrn, alidrn, &
-                               albicen, albsnon, albpndn, apeffn
+                               albicen, albsnon, albpndn, apeffn, snowfracn
       use ice_state, only: aicen, aice, aice_init, nbtrcr
       use ice_therm_shared, only: calc_Tsfc
       use ice_timers, only: timer_couple, ice_timer_start, ice_timer_stop
@@ -332,7 +332,9 @@
          i,j         , & ! horizontal indices
          k               ! tracer index
 
-      real (kind=dbl_kind) :: cszn ! counter for history averaging
+      real (kind=dbl_kind) :: &
+         cszn        , & ! counter for history averaging
+         netsw           ! flag for shortwave radiation presence
 
       !-----------------------------------------------------------------
       ! Save current value of frzmlt for diagnostics.
@@ -365,10 +367,12 @@
             albsno(i,j,iblk) = c0
             albpnd(i,j,iblk) = c0
             apeff_ai(i,j,iblk) = c0
+            snowfrac(i,j,iblk) = c0
 
             ! for history averaging
             cszn = c0
-            if (coszen(i,j,iblk) > puny) cszn = c1
+            netsw = swvdr(i,j,iblk)+swidr(i,j,iblk)+swvdf(i,j,iblk)+swidf(i,j,iblk)
+            if (netsw > puny) cszn = c1
             do n = 1, nstreams
                albcnt(i,j,iblk,n) = albcnt(i,j,iblk,n) + cszn
             enddo
@@ -397,6 +401,8 @@
 
             apeff_ai(i,j,iblk) = apeff_ai(i,j,iblk) &       ! for history
                + apeffn(i,j,n,iblk)*aicen(i,j,n,iblk)
+            snowfrac(i,j,iblk) = snowfrac(i,j,iblk) &       ! for history
+               + snowfracn(i,j,n,iblk)*aicen(i,j,n,iblk)
          enddo
          enddo
          enddo

cice/drivers/hadgem3/CICE_RunMod.F90

@@ -286,11 +286,11 @@ subroutine coupling_prep (iblk)
           fswthru_ai, fhocn, fswthru, scale_factor, &
           swvdr, swidr, swvdf, swidf, Tf, Tair, Qa, strairxT, strairyt, &
           fsens, flat, fswabs, flwout, evap, Tref, Qref, faero_ocn, &
-          fsurfn_f, flatn_f, scale_fluxes, frzmlt_init, frzmlt
+          fsurfn_f, flatn_f, scale_fluxes, frzmlt_init, frzmlt, snowfrac
       use ice_grid, only: tmask
       use ice_ocean, only: oceanmixed_ice, ocean_mixed_layer
       use ice_shortwave, only: alvdfn, alidfn, alvdrn, alidrn, &
-                               albicen, albsnon, albpndn, apeffn
+                               albicen, albsnon, albpndn, apeffn, snowfracn
       use ice_state, only: aicen, aice, aice_init, nbtrcr
       use ice_therm_shared, only: calc_Tsfc
       use ice_timers, only: timer_couple, ice_timer_start, ice_timer_stop
@@ -309,7 +309,9 @@ subroutine coupling_prep (iblk)
          k           , & ! tracer index
          ilo,ihi,jlo,jhi ! beginning and end of physical domain
 
-      real (kind=dbl_kind) :: cszn ! counter for history averaging
+      real (kind=dbl_kind) :: &
+         cszn        , & ! counter for history averaging
+         netsw           ! flag for shortwave radiation presence
 
       !-----------------------------------------------------------------
       ! Save current value of frzmlt for diagnostics.
@@ -342,10 +344,12 @@ subroutine coupling_prep (iblk)
             albsno(i,j,iblk) = c0
             albpnd(i,j,iblk) = c0
             apeff_ai(i,j,iblk) = c0
+            snowfrac(i,j,iblk) = c0
 
             ! for history averaging
             cszn = c0
-            if (coszen(i,j,iblk) > puny) cszn = c1
+            netsw = swvdr(i,j,iblk)+swidr(i,j,iblk)+swvdf(i,j,iblk)+swidf(i,j,iblk)
+            if (netsw > puny) cszn = c1
             do n = 1, nstreams
                albcnt(i,j,iblk,n) = albcnt(i,j,iblk,n) + cszn
             enddo
@@ -374,6 +378,8 @@ subroutine coupling_prep (iblk)
 
             apeff_ai(i,j,iblk) = apeff_ai(i,j,iblk) &       ! for history
                + apeffn(i,j,n,iblk)*aicen(i,j,n,iblk)
+            snowfrac(i,j,iblk) = snowfrac(i,j,iblk) &       ! for history
+               + snowfracn(i,j,n,iblk)*aicen(i,j,n,iblk)
          enddo
          enddo
          enddo

cice/input_templates/col/ice_in

@@ -207,6 +207,7 @@
   , f_vatm         = 'm' 
   , f_fswdn        = 'm' 
   , f_flwdn        = 'm'
+  , f_snowfrac     = 'x'
   , f_snow         = 'x' 
   , f_snow_ai      = 'm' 
   , f_rain         = 'x' 
@@ -225,6 +226,10 @@
   , f_alidr        = 'x'
   , f_alvdf        = 'x'
   , f_alidf        = 'x'
+  , f_alvdr_ai     = 'x'
+  , f_alidr_ai     = 'x'
+  , f_alvdf_ai     = 'x'
+  , f_alidf_ai     = 'x'
   , f_albice       = 'x'
   , f_albsno       = 'x'
   , f_albpnd       = 'x'
@@ -233,6 +238,7 @@
   , f_flat_ai      = 'm' 
   , f_fsens        = 'x' 
   , f_fsens_ai     = 'm' 
+  , f_fswup        = 'x' 
   , f_flwup        = 'x' 
   , f_flwup_ai     = 'm' 
   , f_evap         = 'x' 
@@ -291,6 +297,7 @@
   , f_aicen        = 'm'
   , f_vicen        = 'm'
   , f_vsnon        = 'x'
+  , f_snowfracn    = 'x'
   , f_keffn_top    = 'x'
   , f_Tinz         = 'm'
   , f_Sinz         = 'm'

cice/input_templates/gx1/ice_in

@@ -207,6 +207,7 @@
   , f_vatm         = 'm' 
   , f_fswdn        = 'm' 
   , f_flwdn        = 'm'
+  , f_snowfrac     = 'x'
   , f_snow         = 'x' 
   , f_snow_ai      = 'm' 
   , f_rain         = 'x' 
@@ -225,6 +226,10 @@
   , f_alidr        = 'x'
   , f_alvdf        = 'x'
   , f_alidf        = 'x'
+  , f_alvdr_ai     = 'x'
+  , f_alidr_ai     = 'x'
+  , f_alvdf_ai     = 'x'
+  , f_alidf_ai     = 'x'
   , f_albice       = 'x'
   , f_albsno       = 'x'
   , f_albpnd       = 'x'
@@ -233,6 +238,7 @@
   , f_flat_ai      = 'm' 
   , f_fsens        = 'x' 
   , f_fsens_ai     = 'm' 
+  , f_fswup        = 'x' 
   , f_flwup        = 'x' 
   , f_flwup_ai     = 'm' 
   , f_evap         = 'x' 
@@ -291,6 +297,7 @@
   , f_aicen        = 'x'
   , f_vicen        = 'x'
   , f_vsnon        = 'x'
+  , f_snowfracn    = 'x'
   , f_keffn_top    = 'x'
   , f_Tinz         = 'x'
   , f_Sinz         = 'x'