Merge pull request #21 from apcraig/icep230317

Merge #1eae173fafa9b29 from CICE-Consortium

.gitignore

@@ -11,3 +11,9 @@ doc/build
 # Ignore testsuite file/directories
 testsuite*
 caselist*
+
+# Ignore compiled .mod files
+*.mod
+
+# Ignore test case directories (no consistent name so we'll just ignore conda)
+conda_*

columnphysics/icepack_fsd.F90

@@ -120,18 +120,10 @@ subroutine icepack_init_fsd_bounds(nfsd, &
       real (kind=dbl_kind), dimension(:), allocatable :: &
          lims
 
-      logical (kind=log_kind) :: &
-         l_write_diags  ! local write diags
-
       character(len=8) :: c_fsd1,c_fsd2
       character(len=2) :: c_nf
       character(len=*), parameter :: subname='(icepack_init_fsd_bounds)'
 
-      l_write_diags = .true.
-      if (present(write_diags)) then
-         l_write_diags = write_diags
-      endif
-
       if (nfsd.eq.24) then
 
          allocate(lims(24+1))
@@ -230,7 +222,8 @@ subroutine icepack_init_fsd_bounds(nfsd, &
          c_fsd_range(n)=c_fsd1//'m < fsd Cat '//c_nf//' < '//c_fsd2//'m'
       enddo
 
-      if (l_write_diags) then
+      if (present(write_diags)) then
+      if (write_diags) then
          write(warnstr,*) ' '
          call icepack_warnings_add(warnstr)
          write(warnstr,*) subname
@@ -244,6 +237,7 @@ subroutine icepack_init_fsd_bounds(nfsd, &
          write(warnstr,*) ' '
          call icepack_warnings_add(warnstr)
       endif
+      endif
 
       end subroutine icepack_init_fsd_bounds
 
@@ -705,7 +699,10 @@ subroutine fsd_add_new_ice (ncat, n,    nfsd,          &
          DO WHILE (elapsed_t.lt.dt)
 
              nsubt = nsubt + 1
-             if (nsubt.gt.100) print *, 'latg not converging'
+             if (nsubt.gt.100) then
+                write(warnstr,*) subname,'latg not converging'
+                call icepack_warnings_add(warnstr)
+             endif
 
              ! finite differences
              df_flx(:) = c0 ! NB could stay zero if all in largest FS cat
@@ -717,8 +714,6 @@ subroutine fsd_add_new_ice (ncat, n,    nfsd,          &
                 df_flx(k) = f_flx(k+1) - f_flx(k)
              end do
 
-!         if (abs(sum(df_flx)) > puny) print*,'fsd_add_new ERROR df_flx /= 0'
-
              dafsd_tmp(:) = c0
              do k = 1, nfsd
                 dafsd_tmp(k) = (-df_flx(k) + c2 * G_radial * afsdn_latg(k,n) &
@@ -937,7 +932,6 @@ subroutine fsd_weld_thermo (ncat,  nfsd,   &
          gain, loss     ! welding tendencies
 
       real(kind=dbl_kind) :: &
-         prefac     , & ! multiplies kernel
          kern       , & ! kernel
          subdt      , & ! subcycling time step for stability (s)
          elapsed_t      ! elapsed subcycling time
@@ -948,7 +942,6 @@ subroutine fsd_weld_thermo (ncat,  nfsd,   &
       afsdn  (:,:) = c0
       afsd_init(:) = c0
       stability    = c0
-      prefac       = p5
 
       do n = 1, ncat
 
@@ -992,8 +985,7 @@ subroutine fsd_weld_thermo (ncat,  nfsd,   &
                    if (k > i) then
                        kern = c_weld * floe_area_c(i) * aicen(n)
                        loss(i) = loss(i) + kern*afsd_tmp(i)*afsd_tmp(j)
-                       if (i.eq.j) prefac = c1 ! otherwise 0.5
-                       gain(k) = gain(k) + prefac*kern*afsd_tmp(i)*afsd_tmp(j)
+                       gain(k) = gain(k) + kern*afsd_tmp(i)*afsd_tmp(j)
                    end if
                end do
                end do
@@ -1017,11 +1009,11 @@ subroutine fsd_weld_thermo (ncat,  nfsd,   &
 
             END DO ! time
 
+            afsdn(:,n) = afsd_tmp(:)
             call icepack_cleanup_fsdn (nfsd, afsdn(:,n))
             if (icepack_warnings_aborted(subname)) return
 
             do k = 1, nfsd
-               afsdn(k,n) = afsd_tmp(k)
                trcrn(nt_fsd+k-1,n) = afsdn(k,n)
                ! history/diagnostics
                d_afsdn_weld(k,n) = afsdn(k,n) - afsd_init(k)

columnphysics/icepack_itd.F90

@@ -1006,10 +1006,12 @@ subroutine cleanup_itd (dt,          ntrcr,      &
            faero_ocn(it) = faero_ocn(it) + dfaero_ocn(it)
          enddo
       endif
-      if (tr_iso) then
-         do it = 1, n_iso
-           fiso_ocn(it) = fiso_ocn(it) + dfiso_ocn(it)
-         enddo
+      if (present(fiso_ocn)) then
+         if (tr_iso) then
+            do it = 1, n_iso
+              fiso_ocn(it) = fiso_ocn(it) + dfiso_ocn(it)
+            enddo
+         endif
       endif
       if (present(flux_bio)) then
          do it = 1, nbtrcr

columnphysics/icepack_mechred.F90

@@ -58,10 +58,7 @@ module icepack_mechred
       implicit none
 
       private
-      public :: ridge_ice, &
-                asum_ridging, &
-                ridge_itd, &
-                icepack_ice_strength, &
+      public :: icepack_ice_strength, &
                 icepack_step_ridge
 
       real (kind=dbl_kind), parameter :: &
@@ -113,7 +110,7 @@ subroutine ridge_ice (dt,          ndtd,       &
                             dardg1ndt,   dardg2ndt,  &
                             dvirdgndt,   Tf,         &
                             araftn,      vraftn,     &
-                            closing_flag,closing )
+                            closing )
 
       integer (kind=int_kind), intent(in) :: &
          ndtd       , & ! number of dynamics subcycles
@@ -164,7 +161,6 @@ subroutine ridge_ice (dt,          ndtd,       &
          krdg_redist    ! selects redistribution function
 
       logical (kind=log_kind), intent(in) :: &
-         closing_flag, &! flag if closing is valid
          tr_brine       ! if .true., brine height differs from ice thickness
 
       ! optional history fields
@@ -299,7 +295,7 @@ subroutine ridge_ice (dt,          ndtd,       &
       ! Compute the area opening and closing.
       !-----------------------------------------------------------------
 
-      if (closing_flag) then
+      if (present(opening) .and. present(closing)) then
 
          opning = opening
          closing_net = closing
@@ -600,11 +596,13 @@ subroutine ridge_ice (dt,          ndtd,       &
             faero_ocn(it) = faero_ocn(it) + maero(it)*dti
          enddo
       endif
-      if (tr_iso) then
-         ! check size fiso_ocn vs n_iso ???
-         do it = 1, n_iso
-            fiso_ocn(it) = fiso_ocn(it) + miso(it)*dti
-         enddo
+      if (present(fiso_ocn)) then
+         if (tr_iso) then
+            ! check size fiso_ocn vs n_iso ???
+            do it = 1, n_iso
+               fiso_ocn(it) = fiso_ocn(it) + miso(it)*dti
+            enddo
+         endif
       endif
       if (present(fpond)) then
          fpond = fpond - mpond ! units change later
@@ -1837,12 +1835,6 @@ subroutine icepack_step_ridge (dt,           ndtd,         &
       real (kind=dbl_kind) :: &
          dtt          ! thermo time step
 
-      real (kind=dbl_kind) :: &
-         l_closing      ! local rate of closing due to divergence/shear (1/s)
-
-      logical (kind=log_kind) :: &
-         l_closing_flag ! flag if closing is passed
-
       logical (kind=log_kind), save :: &
          first_call = .true.   ! first call flag
 
@@ -1870,14 +1862,6 @@ subroutine icepack_step_ridge (dt,           ndtd,         &
       !        it may be out of whack, which the ridging helps fix).-ECH
       !-----------------------------------------------------------------
 
-      if (present(closing)) then
-         l_closing_flag = .true.
-         l_closing = closing
-      else
-         l_closing_flag = .false.
-         l_closing = c0
-      endif
-
       call ridge_ice (dt,           ndtd,           &
                       ncat,         n_aero,         &
                       nilyr,        nslyr,          &
@@ -1903,8 +1887,7 @@ subroutine icepack_step_ridge (dt,           ndtd,         &
                       dardg1ndt,    dardg2ndt,      &
                       dvirdgndt,    Tf,             &
                       araftn,       vraftn,         &
-                      l_closing_flag,               &
-                      l_closing )
+                      closing )
       if (icepack_warnings_aborted(subname)) return
 
       !-----------------------------------------------------------------

columnphysics/icepack_orbital.F90

@@ -177,11 +177,24 @@ subroutine compute_coszen (tlat,          tlon,        &
 
       real (kind=dbl_kind) :: ydayp1 ! day of year plus one time step
 
+      logical (kind=log_kind), save :: &
+         first_call = .true.   ! first call flag
+
       character(len=*),parameter :: subname='(compute_coszen)'
 
 ! Solar declination for next time step
 
 #ifdef CESMCOUPLED
+      if (icepack_chkoptargflag(first_call)) then
+         if (.not.(present(days_per_year) .and. &
+                   present(nextsw_cday)   .and. &
+                   present(calendar_type))) then
+            call icepack_warnings_add(subname//' error in CESMCOUPLED args')
+            call icepack_warnings_setabort(.true.,__FILE__,__LINE__)
+            return
+         endif
+      endif
+
       if (calendar_type == "GREGORIAN") then
          ydayp1 = min(nextsw_cday, real(days_per_year,kind=dbl_kind))
       else
@@ -206,6 +219,8 @@ subroutine compute_coszen (tlat,          tlon,        &
              *cos((sec/secday-p5)*c2*pi + tlon) !cos(hour angle)
 #endif
 
+      first_call = .false.
+
       end subroutine compute_coszen
 
 !===============================================================================

columnphysics/icepack_shortwave.F90

@@ -111,8 +111,8 @@ module icepack_shortwave
       implicit none
 
       private
-      public :: icepack_init_radiation, &
-                icepack_prep_radiation, &
+      public :: icepack_prep_radiation, &
+                icepack_init_radiation, &
                 icepack_step_radiation
 
       real (kind=dbl_kind), parameter :: &