From 4e51baea178a1bec83017acc7d870fd992397d54 Mon Sep 17 00:00:00 2001
From: "denise.worthen" <Denise.Worthen@noaa.gov>
Date: Fri, 7 Aug 2020 10:48:38 -0400
Subject: [PATCH 01/12] Revert "Merge remote-tracking branch 'upstream/dev/emc'
 into dev/emc"

This reverts commit e4ca1dcb0fab0f6db08f23e3f15c43035a4b75c4, reversing
changes made to 40bfb4bfe375e66df5ea1127d9799656251dc7fd.
---
 config_src/nuopc_driver/mom_cap.F90         |  35 +---
 config_src/nuopc_driver/mom_cap_methods.F90 | 183 +-------------------
 2 files changed, 8 insertions(+), 210 deletions(-)

diff --git a/config_src/nuopc_driver/mom_cap.F90 b/config_src/nuopc_driver/mom_cap.F90
index d49f370a47..a8056129ff 100644
--- a/config_src/nuopc_driver/mom_cap.F90
+++ b/config_src/nuopc_driver/mom_cap.F90
@@ -26,7 +26,7 @@ module MOM_cap_mod
 use time_manager_mod,         only: fms_get_calendar_type => get_calendar_type
 use MOM_domains,              only: MOM_infra_init, num_pes, root_pe, pe_here
 use MOM_file_parser,          only: get_param, log_version, param_file_type, close_param_file
-use MOM_get_input,            only: get_MOM_input, directories
+use MOM_get_input,            only: Get_MOM_Input, directories
 use MOM_domains,              only: pass_var
 use MOM_error_handler,        only: MOM_error, FATAL, is_root_pe
 use MOM_ocean_model_nuopc,    only: ice_ocean_boundary_type
@@ -36,7 +36,7 @@ module MOM_cap_mod
 use MOM_ocean_model_nuopc,    only: ocean_model_init, update_ocean_model, ocean_model_end
 use MOM_ocean_model_nuopc,    only: get_ocean_grid, get_eps_omesh
 use MOM_cap_time,             only: AlarmInit
-use MOM_cap_methods,          only: mom_import, mom_export, mom_set_geomtype, state_diagnose
+use MOM_cap_methods,          only: mom_import, mom_export, mom_set_geomtype
 #ifdef CESMCOUPLED
 use shr_file_mod,             only: shr_file_setLogUnit, shr_file_getLogUnit
 #endif
@@ -124,7 +124,7 @@ module MOM_cap_mod
 integer              :: fldsFrOcn_num = 0
 type (fld_list_type) :: fldsFrOcn(fldsMax)
 
-integer              :: dbug = 0
+integer              :: debug = 0
 integer              :: import_slice = 1
 integer              :: export_slice = 1
 character(len=256)   :: tmpstr
@@ -273,14 +273,6 @@ subroutine InitializeP0(gcomp, importState, exportState, clock, rc)
   write(logmsg,*) grid_attach_area
   call ESMF_LogWrite('MOM_cap:GridAttachArea = '//trim(logmsg), ESMF_LOGMSG_INFO)
 
-  call NUOPC_CompAttributeGet(gcomp, name='dbug_flag', value=value, isPresent=isPresent, isSet=isSet, rc=rc)
-  if (ChkErr(rc,__LINE__,u_FILE_u)) return
-  if (isPresent .and. isSet) then
-   read(value,*) dbug
-  end if
-  write(logmsg,'(i6)') dbug
-  call ESMF_LogWrite('MOM_cap:dbug = '//trim(logmsg), ESMF_LOGMSG_INFO)
-
   scalar_field_name = ""
   call NUOPC_CompAttributeGet(gcomp, name="ScalarFieldName", value=value, &
        isPresent=isPresent, isSet=isSet, rc=rc)
@@ -366,7 +358,6 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc)
   type(ice_ocean_boundary_type), pointer :: Ice_ocean_boundary => NULL()
   type(ocean_internalstate_wrapper)      :: ocean_internalstate
   type(ocean_grid_type),         pointer :: ocean_grid => NULL()
-  type(directories)                      :: dirs
   type(time_type)                        :: Run_len      !< length of experiment
   type(time_type)                        :: time0        !< Start time of coupled model's calendar.
   type(time_type)                        :: time_start   !< The time at which to initialize the ocean model
@@ -529,13 +520,8 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc)
 
   restartfile = ""
   if (runtype == "initial") then
-    if (cesm_coupled) then
-      restartfile = "n"
-    else
-      call get_MOM_input(dirs=dirs)
-      restartfile = dirs%input_filename(1:1)
-    endif
-    call ESMF_LogWrite('MOM_cap:restartfile = '//trim(restartfile), ESMF_LOGMSG_INFO)
+
+     restartfile = "n"
 
   else if (runtype == "continue") then ! hybrid or branch or continuos runs
 
@@ -835,7 +821,7 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc)
   allocate(xb(ntiles),xe(ntiles),yb(ntiles),ye(ntiles),pe(ntiles))
   call mpp_get_compute_domains(ocean_public%domain, xbegin=xb, xend=xe, ybegin=yb, yend=ye)
   call mpp_get_pelist(ocean_public%domain, pe)
-  if (dbug > 1) then
+  if (debug > 0) then
      do n = 1,ntiles
         write(tmpstr,'(a,6i6)') subname//' tiles ',n,pe(n),xb(n),xe(n),yb(n),ye(n)
         call ESMF_LogWrite(trim(tmpstr), ESMF_LOGMSG_INFO)
@@ -1445,11 +1431,6 @@ subroutine ModelAdvance(gcomp, rc)
       enddo
      endif
 
-     if (dbug > 0) then
-       call state_diagnose(importState,subname//':IS ',rc=rc)
-       if (ChkErr(rc,__LINE__,u_FILE_u)) return
-     end if
-
      !---------------
      ! Get ocean grid
      !---------------
@@ -1478,10 +1459,6 @@ subroutine ModelAdvance(gcomp, rc)
      call mom_export(ocean_public, ocean_grid, ocean_state, exportState, clock, rc=rc)
      if (ChkErr(rc,__LINE__,u_FILE_u)) return
 
-     if (dbug > 0) then
-       call state_diagnose(exportState,subname//':ES ',rc=rc)
-       if (ChkErr(rc,__LINE__,u_FILE_u)) return
-     end if
   endif
 
   !---------------
diff --git a/config_src/nuopc_driver/mom_cap_methods.F90 b/config_src/nuopc_driver/mom_cap_methods.F90
index 0997fbc635..8aca45094f 100644
--- a/config_src/nuopc_driver/mom_cap_methods.F90
+++ b/config_src/nuopc_driver/mom_cap_methods.F90
@@ -5,7 +5,7 @@ module MOM_cap_methods
 use ESMF,                      only: ESMF_TimeInterval, ESMF_TimeIntervalGet
 use ESMF,                      only: ESMF_State, ESMF_StateGet
 use ESMF,                      only: ESMF_Field, ESMF_FieldGet, ESMF_FieldCreate
-use ESMF,                      only: ESMF_GridComp, ESMF_Mesh, ESMF_MeshGet, ESMF_Grid, ESMF_GridCreate
+use ESMF,                      only: ESMF_GridComp, ESMF_Mesh, ESMF_Grid, ESMF_GridCreate
 use ESMF,                      only: ESMF_DistGrid, ESMF_DistGridCreate
 use ESMF,                      only: ESMF_KIND_R8, ESMF_SUCCESS, ESMF_LogFoundError
 use ESMF,                      only: ESMF_LOGERR_PASSTHRU, ESMF_LOGMSG_INFO, ESMF_LOGWRITE
@@ -13,8 +13,7 @@ module MOM_cap_methods
 use ESMF,                      only: ESMF_StateItem_Flag, ESMF_STATEITEM_NOTFOUND
 use ESMF,                      only: ESMF_GEOMTYPE_FLAG, ESMF_GEOMTYPE_GRID, ESMF_GEOMTYPE_MESH
 use ESMF,                      only: ESMF_RC_VAL_OUTOFRANGE, ESMF_INDEX_DELOCAL, ESMF_MESHLOC_ELEMENT
-use ESMF,                      only: ESMF_TYPEKIND_R8, ESMF_FIELDSTATUS_COMPLETE
-use ESMF,                      only: ESMF_FieldStatus_Flag, ESMF_LOGMSG_ERROR, ESMF_FAILURE, ESMF_MAXSTR
+use ESMF,                      only: ESMF_TYPEKIND_R8
 use ESMF,                      only: operator(/=), operator(==)
 use MOM_ocean_model_nuopc,     only: ocean_public_type, ocean_state_type
 use MOM_surface_forcing_nuopc, only: ice_ocean_boundary_type
@@ -29,7 +28,6 @@ module MOM_cap_methods
 public :: mom_set_geomtype
 public :: mom_import
 public :: mom_export
-public :: state_diagnose
 
 private :: State_getImport
 private :: State_setExport
@@ -765,183 +763,6 @@ subroutine State_SetExport(state, fldname, isc, iec, jsc, jec, input, ocean_grid
 
 end subroutine State_SetExport
 
-subroutine state_diagnose(State, string, rc)
-
-  ! ----------------------------------------------
-  ! Diagnose status of State
-  ! ----------------------------------------------
-
-  type(ESMF_State), intent(in)  :: state
-  character(len=*), intent(in)  :: string
-  integer         , intent(out) :: rc
-
-  ! local variables
-  integer                         :: i,j,n
-  type(ESMf_Field)                :: lfield
-  integer                         :: fieldCount, lrank
-  character(ESMF_MAXSTR) ,pointer :: lfieldnamelist(:)
-  real(ESMF_KIND_R8), pointer     :: dataPtr1d(:)
-  real(ESMF_KIND_R8), pointer     :: dataPtr2d(:,:)
-  character(len=*),parameter      :: subname='(state_diagnose)'
-  character(len=ESMF_MAXSTR)      :: msgString
-  ! ----------------------------------------------
-
-  call ESMF_StateGet(state, itemCount=fieldCount, rc=rc)
-  if (chkerr(rc,__LINE__,u_FILE_u)) return
-  allocate(lfieldnamelist(fieldCount))
-
-  call ESMF_StateGet(state, itemNameList=lfieldnamelist, rc=rc)
-  if (chkerr(rc,__LINE__,u_FILE_u)) return
-
-  do n = 1, fieldCount
-
-     call ESMF_StateGet(state, itemName=lfieldnamelist(n), field=lfield, rc=rc)
-     if (chkerr(rc,__LINE__,u_FILE_u)) return
-
-     call field_getfldptr(lfield, fldptr1=dataPtr1d, fldptr2=dataPtr2d, rank=lrank, rc=rc)
-     if (chkerr(rc,__LINE__,u_FILE_u)) return
-
-     if (lrank == 0) then
-        ! no local data
-     elseif (lrank == 1) then
-        if (size(dataPtr1d) > 0) then
-           write(msgString,'(A,3g14.7,i8)') trim(string)//': '//trim(lfieldnamelist(n)), &
-                minval(dataPtr1d), maxval(dataPtr1d), sum(dataPtr1d), size(dataPtr1d)
-        else
-           write(msgString,'(A,a)') trim(string)//': '//trim(lfieldnamelist(n))," no data"
-        endif
-     elseif (lrank == 2) then
-        if (size(dataPtr2d) > 0) then
-           write(msgString,'(A,3g14.7,i8)') trim(string)//': '//trim(lfieldnamelist(n)), &
-                minval(dataPtr2d), maxval(dataPtr2d), sum(dataPtr2d), size(dataPtr2d)
-        else
-           write(msgString,'(A,a)') trim(string)//': '//trim(lfieldnamelist(n))," no data"
-        endif
-     else
-        call ESMF_LogWrite(trim(subname)//": ERROR rank not supported ", ESMF_LOGMSG_ERROR)
-        rc = ESMF_FAILURE
-        return
-     endif
-     call ESMF_LogWrite(trim(msgString), ESMF_LOGMSG_INFO)
-  enddo
-
-  deallocate(lfieldnamelist)
-
-end subroutine state_diagnose
-
-!===============================================================================
-
-subroutine field_getfldptr(field, fldptr1, fldptr2, rank, abort, rc)
-
-  ! ----------------------------------------------
-  ! for a field, determine rank and return fldptr1 or fldptr2
-  ! abort is true by default and will abort if fldptr is not yet allocated in field
-  ! rank returns 0, 1, or 2.  0 means fldptr not allocated and abort=false
-  ! ----------------------------------------------
-
-  ! input/output variables
-  type(ESMF_Field)  , intent(in)              :: field
-  real(ESMF_KIND_R8), pointer , intent(inout), optional :: fldptr1(:)
-  real(ESMF_KIND_R8), pointer , intent(inout), optional :: fldptr2(:,:)
-  integer           , intent(out)  , optional :: rank
-  logical           , intent(in)   , optional :: abort
-  integer           , intent(out)  , optional :: rc
-
-  ! local variables
-  type(ESMF_GeomType_Flag)    :: geomtype
-  type(ESMF_FieldStatus_Flag) :: status
-  type(ESMF_Mesh)             :: lmesh
-  integer                     :: lrank, nnodes, nelements
-  logical                     :: labort
-  character(len=*), parameter :: subname='(field_getfldptr)'
-  ! ----------------------------------------------
-
-  if (.not.present(rc)) then
-     call ESMF_LogWrite(trim(subname)//": ERROR rc not present ", &
-          ESMF_LOGMSG_ERROR, line=__LINE__, file=u_FILE_u)
-     rc = ESMF_FAILURE
-     return
-  endif
-
-  rc = ESMF_SUCCESS
-
-  labort = .true.
-  if (present(abort)) then
-     labort = abort
-  endif
-  lrank = -99
-
-  call ESMF_FieldGet(field, status=status, rc=rc)
-  if (chkerr(rc,__LINE__,u_FILE_u)) return
-
-  if (status /= ESMF_FIELDSTATUS_COMPLETE) then
-     lrank = 0
-     if (labort) then
-        call ESMF_LogWrite(trim(subname)//": ERROR data not allocated ", ESMF_LOGMSG_INFO, rc=rc)
-        rc = ESMF_FAILURE
-        return
-     else
-        call ESMF_LogWrite(trim(subname)//": WARNING data not allocated ", ESMF_LOGMSG_INFO, rc=rc)
-     endif
-  else
-
-     call ESMF_FieldGet(field, geomtype=geomtype, rc=rc)
-     if (chkerr(rc,__LINE__,u_FILE_u)) return
-
-     if (geomtype == ESMF_GEOMTYPE_GRID) then
-        call ESMF_FieldGet(field, rank=lrank, rc=rc)
-        if (chkerr(rc,__LINE__,u_FILE_u)) return
-     elseif (geomtype == ESMF_GEOMTYPE_MESH) then
-        call ESMF_FieldGet(field, rank=lrank, rc=rc)
-        if (chkerr(rc,__LINE__,u_FILE_u)) return
-        call ESMF_FieldGet(field, mesh=lmesh, rc=rc)
-        if (chkerr(rc,__LINE__,u_FILE_u)) return
-        call ESMF_MeshGet(lmesh, numOwnedNodes=nnodes, numOwnedElements=nelements, rc=rc)
-        if (chkerr(rc,__LINE__,u_FILE_u)) return
-        if (nnodes == 0 .and. nelements == 0) lrank = 0
-     else  
-        call ESMF_LogWrite(trim(subname)//": ERROR geomtype not supported ", &
-             ESMF_LOGMSG_INFO, rc=rc)
-        rc = ESMF_FAILURE
-        return
-     endif ! geomtype
-
-     if (lrank == 0) then
-        call ESMF_LogWrite(trim(subname)//": no local nodes or elements ", &
-             ESMF_LOGMSG_INFO)
-     elseif (lrank == 1) then
-        if (.not.present(fldptr1)) then
-           call ESMF_LogWrite(trim(subname)//": ERROR missing rank=1 array ", &
-                ESMF_LOGMSG_ERROR, line=__LINE__, file=u_FILE_u)
-           rc = ESMF_FAILURE
-           return
-        endif
-        call ESMF_FieldGet(field, farrayPtr=fldptr1, rc=rc)
-        if (chkerr(rc,__LINE__,u_FILE_u)) return
-     elseif (lrank == 2) then
-        if (.not.present(fldptr2)) then
-           call ESMF_LogWrite(trim(subname)//": ERROR missing rank=2 array ", &
-                ESMF_LOGMSG_ERROR, line=__LINE__, file=u_FILE_u)
-           rc = ESMF_FAILURE
-           return
-        endif
-        call ESMF_FieldGet(field, farrayPtr=fldptr2, rc=rc)
-        if (chkerr(rc,__LINE__,u_FILE_u)) return
-     else
-        call ESMF_LogWrite(trim(subname)//": ERROR in rank ", &
-             ESMF_LOGMSG_ERROR, line=__LINE__, file=u_FILE_u)
-        rc = ESMF_FAILURE
-        return
-     endif
-
-  endif  ! status
-
-  if (present(rank)) then
-     rank = lrank
-  endif
-
-end subroutine field_getfldptr
-
 logical function chkerr(rc, line, file)
   integer, intent(in) :: rc
   integer, intent(in) :: line

From 6164f63a9846b1eba00cac2b264dcf4d73b27013 Mon Sep 17 00:00:00 2001
From: "denise.worthen" <Denise.Worthen@noaa.gov>
Date: Tue, 1 Sep 2020 08:08:17 -0400
Subject: [PATCH 02/12] manual add of dev/emc changes (statediagnose)

---
 config_src/nuopc_driver/mom_cap.F90         |  32 +++-
 config_src/nuopc_driver/mom_cap_methods.F90 | 182 +++++++++++++++++++-
 2 files changed, 207 insertions(+), 7 deletions(-)

diff --git a/config_src/nuopc_driver/mom_cap.F90 b/config_src/nuopc_driver/mom_cap.F90
index a8056129ff..f70eddebdb 100644
--- a/config_src/nuopc_driver/mom_cap.F90
+++ b/config_src/nuopc_driver/mom_cap.F90
@@ -26,7 +26,7 @@ module MOM_cap_mod
 use time_manager_mod,         only: fms_get_calendar_type => get_calendar_type
 use MOM_domains,              only: MOM_infra_init, num_pes, root_pe, pe_here
 use MOM_file_parser,          only: get_param, log_version, param_file_type, close_param_file
-use MOM_get_input,            only: Get_MOM_Input, directories
+use MOM_get_input,            only: get_MOM_input, directories
 use MOM_domains,              only: pass_var
 use MOM_error_handler,        only: MOM_error, FATAL, is_root_pe
 use MOM_ocean_model_nuopc,    only: ice_ocean_boundary_type
@@ -36,7 +36,7 @@ module MOM_cap_mod
 use MOM_ocean_model_nuopc,    only: ocean_model_init, update_ocean_model, ocean_model_end
 use MOM_ocean_model_nuopc,    only: get_ocean_grid, get_eps_omesh
 use MOM_cap_time,             only: AlarmInit
-use MOM_cap_methods,          only: mom_import, mom_export, mom_set_geomtype
+use MOM_cap_methods,          only: mom_import, mom_export, mom_set_geomtype, state_diagnose
 #ifdef CESMCOUPLED
 use shr_file_mod,             only: shr_file_setLogUnit, shr_file_getLogUnit
 #endif
@@ -124,7 +124,7 @@ module MOM_cap_mod
 integer              :: fldsFrOcn_num = 0
 type (fld_list_type) :: fldsFrOcn(fldsMax)
 
-integer              :: debug = 0
+integer              :: dbug = 0
 integer              :: import_slice = 1
 integer              :: export_slice = 1
 character(len=256)   :: tmpstr
@@ -272,6 +272,13 @@ subroutine InitializeP0(gcomp, importState, exportState, clock, rc)
   if (isPresent .and. isSet) grid_attach_area=(trim(value)=="true")
   write(logmsg,*) grid_attach_area
   call ESMF_LogWrite('MOM_cap:GridAttachArea = '//trim(logmsg), ESMF_LOGMSG_INFO)
+  call NUOPC_CompAttributeGet(gcomp, name='dbug_flag', value=value, isPresent=isPresent, isSet=isSet, rc=rc)
+  if (ChkErr(rc,__LINE__,u_FILE_u)) return
+  if (isPresent .and. isSet) then
+   read(value,*) dbug
+  end if
+  write(logmsg,'(i6)') dbug
+  call ESMF_LogWrite('MOM_cap:dbug = '//trim(logmsg), ESMF_LOGMSG_INFO)
 
   scalar_field_name = ""
   call NUOPC_CompAttributeGet(gcomp, name="ScalarFieldName", value=value, &
@@ -358,6 +365,7 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc)
   type(ice_ocean_boundary_type), pointer :: Ice_ocean_boundary => NULL()
   type(ocean_internalstate_wrapper)      :: ocean_internalstate
   type(ocean_grid_type),         pointer :: ocean_grid => NULL()
+  type(directories)                      :: dirs
   type(time_type)                        :: Run_len      !< length of experiment
   type(time_type)                        :: time0        !< Start time of coupled model's calendar.
   type(time_type)                        :: time_start   !< The time at which to initialize the ocean model
@@ -521,7 +529,13 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc)
   restartfile = ""
   if (runtype == "initial") then
 
-     restartfile = "n"
+    if (cesm_coupled) then
+      restartfile = "n"
+    else
+      call get_MOM_input(dirs=dirs)
+      restartfile = dirs%input_filename(1:1)
+    endif
+    call ESMF_LogWrite('MOM_cap:restartfile = '//trim(restartfile), ESMF_LOGMSG_INFO)
 
   else if (runtype == "continue") then ! hybrid or branch or continuos runs
 
@@ -821,7 +835,7 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc)
   allocate(xb(ntiles),xe(ntiles),yb(ntiles),ye(ntiles),pe(ntiles))
   call mpp_get_compute_domains(ocean_public%domain, xbegin=xb, xend=xe, ybegin=yb, yend=ye)
   call mpp_get_pelist(ocean_public%domain, pe)
-  if (debug > 0) then
+  if (dbug > 1) then
      do n = 1,ntiles
         write(tmpstr,'(a,6i6)') subname//' tiles ',n,pe(n),xb(n),xe(n),yb(n),ye(n)
         call ESMF_LogWrite(trim(tmpstr), ESMF_LOGMSG_INFO)
@@ -1430,6 +1444,10 @@ subroutine ModelAdvance(gcomp, rc)
        endif
       enddo
      endif
+     if (dbug > 0) then
+       call state_diagnose(importState,subname//':IS ',rc=rc)
+       if (ChkErr(rc,__LINE__,u_FILE_u)) return
+     end if
 
      !---------------
      ! Get ocean grid
@@ -1459,6 +1477,10 @@ subroutine ModelAdvance(gcomp, rc)
      call mom_export(ocean_public, ocean_grid, ocean_state, exportState, clock, rc=rc)
      if (ChkErr(rc,__LINE__,u_FILE_u)) return
 
+     if (dbug > 0) then
+       call state_diagnose(exportState,subname//':ES ',rc=rc)
+       if (ChkErr(rc,__LINE__,u_FILE_u)) return
+     end if
   endif
 
   !---------------
diff --git a/config_src/nuopc_driver/mom_cap_methods.F90 b/config_src/nuopc_driver/mom_cap_methods.F90
index 8aca45094f..22938d91f4 100644
--- a/config_src/nuopc_driver/mom_cap_methods.F90
+++ b/config_src/nuopc_driver/mom_cap_methods.F90
@@ -5,7 +5,7 @@ module MOM_cap_methods
 use ESMF,                      only: ESMF_TimeInterval, ESMF_TimeIntervalGet
 use ESMF,                      only: ESMF_State, ESMF_StateGet
 use ESMF,                      only: ESMF_Field, ESMF_FieldGet, ESMF_FieldCreate
-use ESMF,                      only: ESMF_GridComp, ESMF_Mesh, ESMF_Grid, ESMF_GridCreate
+use ESMF,                      only: ESMF_GridComp, ESMF_Mesh, ESMF_MeshGet, ESMF_Grid, ESMF_GridCreate
 use ESMF,                      only: ESMF_DistGrid, ESMF_DistGridCreate
 use ESMF,                      only: ESMF_KIND_R8, ESMF_SUCCESS, ESMF_LogFoundError
 use ESMF,                      only: ESMF_LOGERR_PASSTHRU, ESMF_LOGMSG_INFO, ESMF_LOGWRITE
@@ -13,7 +13,8 @@ module MOM_cap_methods
 use ESMF,                      only: ESMF_StateItem_Flag, ESMF_STATEITEM_NOTFOUND
 use ESMF,                      only: ESMF_GEOMTYPE_FLAG, ESMF_GEOMTYPE_GRID, ESMF_GEOMTYPE_MESH
 use ESMF,                      only: ESMF_RC_VAL_OUTOFRANGE, ESMF_INDEX_DELOCAL, ESMF_MESHLOC_ELEMENT
-use ESMF,                      only: ESMF_TYPEKIND_R8
+use ESMF,                      only: ESMF_TYPEKIND_R8, ESMF_FIELDSTATUS_COMPLETE
+use ESMF,                      only: ESMF_FieldStatus_Flag, ESMF_LOGMSG_ERROR, ESMF_FAILURE, ESMF_MAXSTR
 use ESMF,                      only: operator(/=), operator(==)
 use MOM_ocean_model_nuopc,     only: ocean_public_type, ocean_state_type
 use MOM_surface_forcing_nuopc, only: ice_ocean_boundary_type
@@ -28,6 +29,7 @@ module MOM_cap_methods
 public :: mom_set_geomtype
 public :: mom_import
 public :: mom_export
+public :: state_diagnose
 
 private :: State_getImport
 private :: State_setExport
@@ -762,6 +764,182 @@ subroutine State_SetExport(state, fldname, isc, iec, jsc, jec, input, ocean_grid
   endif
 
 end subroutine State_SetExport
+subroutine state_diagnose(State, string, rc)
+
+  ! ----------------------------------------------
+  ! Diagnose status of State
+  ! ----------------------------------------------
+
+  type(ESMF_State), intent(in)  :: state
+  character(len=*), intent(in)  :: string
+  integer         , intent(out) :: rc
+
+  ! local variables
+  integer                         :: i,j,n
+  type(ESMf_Field)                :: lfield
+  integer                         :: fieldCount, lrank
+  character(ESMF_MAXSTR) ,pointer :: lfieldnamelist(:)
+  real(ESMF_KIND_R8), pointer     :: dataPtr1d(:)
+  real(ESMF_KIND_R8), pointer     :: dataPtr2d(:,:)
+  character(len=*),parameter      :: subname='(state_diagnose)'
+  character(len=ESMF_MAXSTR)      :: msgString
+  ! ----------------------------------------------
+
+  call ESMF_StateGet(state, itemCount=fieldCount, rc=rc)
+  if (chkerr(rc,__LINE__,u_FILE_u)) return
+  allocate(lfieldnamelist(fieldCount))
+
+  call ESMF_StateGet(state, itemNameList=lfieldnamelist, rc=rc)
+  if (chkerr(rc,__LINE__,u_FILE_u)) return
+
+  do n = 1, fieldCount
+
+     call ESMF_StateGet(state, itemName=lfieldnamelist(n), field=lfield, rc=rc)
+     if (chkerr(rc,__LINE__,u_FILE_u)) return
+
+     call field_getfldptr(lfield, fldptr1=dataPtr1d, fldptr2=dataPtr2d, rank=lrank, rc=rc)
+     if (chkerr(rc,__LINE__,u_FILE_u)) return
+
+     if (lrank == 0) then
+        ! no local data
+     elseif (lrank == 1) then
+        if (size(dataPtr1d) > 0) then
+           write(msgString,'(A,3g14.7,i8)') trim(string)//': '//trim(lfieldnamelist(n)), &
+                minval(dataPtr1d), maxval(dataPtr1d), sum(dataPtr1d), size(dataPtr1d)
+        else
+           write(msgString,'(A,a)') trim(string)//': '//trim(lfieldnamelist(n))," no data"
+        endif
+     elseif (lrank == 2) then
+        if (size(dataPtr2d) > 0) then
+           write(msgString,'(A,3g14.7,i8)') trim(string)//': '//trim(lfieldnamelist(n)), &
+                minval(dataPtr2d), maxval(dataPtr2d), sum(dataPtr2d), size(dataPtr2d)
+        else
+           write(msgString,'(A,a)') trim(string)//': '//trim(lfieldnamelist(n))," no data"
+        endif
+     else
+        call ESMF_LogWrite(trim(subname)//": ERROR rank not supported ", ESMF_LOGMSG_ERROR)
+        rc = ESMF_FAILURE
+        return
+     endif
+     call ESMF_LogWrite(trim(msgString), ESMF_LOGMSG_INFO)
+  enddo
+
+  deallocate(lfieldnamelist)
+
+end subroutine state_diagnose
+
+!===============================================================================
+
+subroutine field_getfldptr(field, fldptr1, fldptr2, rank, abort, rc)
+
+  ! ----------------------------------------------
+  ! for a field, determine rank and return fldptr1 or fldptr2
+  ! abort is true by default and will abort if fldptr is not yet allocated in field
+  ! rank returns 0, 1, or 2.  0 means fldptr not allocated and abort=false
+  ! ----------------------------------------------
+
+  ! input/output variables
+  type(ESMF_Field)  , intent(in)              :: field
+  real(ESMF_KIND_R8), pointer , intent(inout), optional :: fldptr1(:)
+  real(ESMF_KIND_R8), pointer , intent(inout), optional :: fldptr2(:,:)
+  integer           , intent(out)  , optional :: rank
+  logical           , intent(in)   , optional :: abort
+  integer           , intent(out)  , optional :: rc
+
+  ! local variables
+  type(ESMF_GeomType_Flag)    :: geomtype
+  type(ESMF_FieldStatus_Flag) :: status
+  type(ESMF_Mesh)             :: lmesh
+  integer                     :: lrank, nnodes, nelements
+  logical                     :: labort
+  character(len=*), parameter :: subname='(field_getfldptr)'
+  ! ----------------------------------------------
+
+  if (.not.present(rc)) then
+     call ESMF_LogWrite(trim(subname)//": ERROR rc not present ", &
+          ESMF_LOGMSG_ERROR, line=__LINE__, file=u_FILE_u)
+     rc = ESMF_FAILURE
+     return
+  endif
+
+  rc = ESMF_SUCCESS
+
+  labort = .true.
+  if (present(abort)) then
+     labort = abort
+  endif
+  lrank = -99
+
+  call ESMF_FieldGet(field, status=status, rc=rc)
+  if (chkerr(rc,__LINE__,u_FILE_u)) return
+
+  if (status /= ESMF_FIELDSTATUS_COMPLETE) then
+     lrank = 0
+     if (labort) then
+        call ESMF_LogWrite(trim(subname)//": ERROR data not allocated ", ESMF_LOGMSG_INFO, rc=rc)
+        rc = ESMF_FAILURE
+        return
+     else
+        call ESMF_LogWrite(trim(subname)//": WARNING data not allocated ", ESMF_LOGMSG_INFO, rc=rc)
+     endif
+  else
+
+     call ESMF_FieldGet(field, geomtype=geomtype, rc=rc)
+     if (chkerr(rc,__LINE__,u_FILE_u)) return
+
+     if (geomtype == ESMF_GEOMTYPE_GRID) then
+        call ESMF_FieldGet(field, rank=lrank, rc=rc)
+        if (chkerr(rc,__LINE__,u_FILE_u)) return
+     elseif (geomtype == ESMF_GEOMTYPE_MESH) then
+        call ESMF_FieldGet(field, rank=lrank, rc=rc)
+        if (chkerr(rc,__LINE__,u_FILE_u)) return
+        call ESMF_FieldGet(field, mesh=lmesh, rc=rc)
+        if (chkerr(rc,__LINE__,u_FILE_u)) return
+        call ESMF_MeshGet(lmesh, numOwnedNodes=nnodes, numOwnedElements=nelements, rc=rc)
+        if (chkerr(rc,__LINE__,u_FILE_u)) return
+        if (nnodes == 0 .and. nelements == 0) lrank = 0
+     else  
+        call ESMF_LogWrite(trim(subname)//": ERROR geomtype not supported ", &
+             ESMF_LOGMSG_INFO, rc=rc)
+        rc = ESMF_FAILURE
+        return
+     endif ! geomtype
+
+     if (lrank == 0) then
+        call ESMF_LogWrite(trim(subname)//": no local nodes or elements ", &
+             ESMF_LOGMSG_INFO)
+     elseif (lrank == 1) then
+        if (.not.present(fldptr1)) then
+           call ESMF_LogWrite(trim(subname)//": ERROR missing rank=1 array ", &
+                ESMF_LOGMSG_ERROR, line=__LINE__, file=u_FILE_u)
+           rc = ESMF_FAILURE
+           return
+        endif
+        call ESMF_FieldGet(field, farrayPtr=fldptr1, rc=rc)
+        if (chkerr(rc,__LINE__,u_FILE_u)) return
+     elseif (lrank == 2) then
+        if (.not.present(fldptr2)) then
+           call ESMF_LogWrite(trim(subname)//": ERROR missing rank=2 array ", &
+                ESMF_LOGMSG_ERROR, line=__LINE__, file=u_FILE_u)
+           rc = ESMF_FAILURE
+           return
+        endif
+        call ESMF_FieldGet(field, farrayPtr=fldptr2, rc=rc)
+        if (chkerr(rc,__LINE__,u_FILE_u)) return
+     else
+        call ESMF_LogWrite(trim(subname)//": ERROR in rank ", &
+             ESMF_LOGMSG_ERROR, line=__LINE__, file=u_FILE_u)
+        rc = ESMF_FAILURE
+        return
+     endif
+
+  endif  ! status
+
+  if (present(rank)) then
+     rank = lrank
+  endif
+
+end subroutine field_getfldptr
 
 logical function chkerr(rc, line, file)
   integer, intent(in) :: rc

From bfbd95aacb7f286b70f7f2c7d9c6d53fa6f15aec Mon Sep 17 00:00:00 2001
From: "denise.worthen" <Denise.Worthen@noaa.gov>
Date: Tue, 1 Sep 2020 08:18:34 -0400
Subject: [PATCH 03/12] white space changes

---
 config_src/nuopc_driver/mom_cap.F90         | 3 ++-
 config_src/nuopc_driver/mom_cap_methods.F90 | 1 +
 2 files changed, 3 insertions(+), 1 deletion(-)

diff --git a/config_src/nuopc_driver/mom_cap.F90 b/config_src/nuopc_driver/mom_cap.F90
index f70eddebdb..d49f370a47 100644
--- a/config_src/nuopc_driver/mom_cap.F90
+++ b/config_src/nuopc_driver/mom_cap.F90
@@ -272,6 +272,7 @@ subroutine InitializeP0(gcomp, importState, exportState, clock, rc)
   if (isPresent .and. isSet) grid_attach_area=(trim(value)=="true")
   write(logmsg,*) grid_attach_area
   call ESMF_LogWrite('MOM_cap:GridAttachArea = '//trim(logmsg), ESMF_LOGMSG_INFO)
+
   call NUOPC_CompAttributeGet(gcomp, name='dbug_flag', value=value, isPresent=isPresent, isSet=isSet, rc=rc)
   if (ChkErr(rc,__LINE__,u_FILE_u)) return
   if (isPresent .and. isSet) then
@@ -528,7 +529,6 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc)
 
   restartfile = ""
   if (runtype == "initial") then
-
     if (cesm_coupled) then
       restartfile = "n"
     else
@@ -1444,6 +1444,7 @@ subroutine ModelAdvance(gcomp, rc)
        endif
       enddo
      endif
+
      if (dbug > 0) then
        call state_diagnose(importState,subname//':IS ',rc=rc)
        if (ChkErr(rc,__LINE__,u_FILE_u)) return
diff --git a/config_src/nuopc_driver/mom_cap_methods.F90 b/config_src/nuopc_driver/mom_cap_methods.F90
index 22938d91f4..0997fbc635 100644
--- a/config_src/nuopc_driver/mom_cap_methods.F90
+++ b/config_src/nuopc_driver/mom_cap_methods.F90
@@ -764,6 +764,7 @@ subroutine State_SetExport(state, fldname, isc, iec, jsc, jec, input, ocean_grid
   endif
 
 end subroutine State_SetExport
+
 subroutine state_diagnose(State, string, rc)
 
   ! ----------------------------------------------

From d502516329fe3ebb28f5465d9ef10d4486f5c674 Mon Sep 17 00:00:00 2001
From: jiandewang <jiande.wang@noaa.gov>
Date: Sat, 19 Dec 2020 22:56:26 -0500
Subject: [PATCH 04/12] minor enhancement in MOM_wave_interface.F90 to aviod
 undocumented error

---
 src/user/MOM_wave_interface.F90 | 6 ++++++
 1 file changed, 6 insertions(+)

diff --git a/src/user/MOM_wave_interface.F90 b/src/user/MOM_wave_interface.F90
index 4ba1b779e3..e79f7b65ef 100644
--- a/src/user/MOM_wave_interface.F90
+++ b/src/user/MOM_wave_interface.F90
@@ -466,6 +466,12 @@ subroutine Update_Surface_Waves(G, GV, US, Day, dt, CS, forces)
     if (DataSource==DATAOVR) then
       call Surface_Bands_by_data_override(day_center, G, GV, US, CS)
     elseif (DataSource==Coupler) then
+      if (.not.present(FORCES)) then 
+        call MOM_error(FATAL,"The code cannot be run with the options  "//&
+             "SURFBAND_SOURCE = COUPLER for with this driver.  If you are using a "//&
+             "wave coupled driver then check the call to update_surface_waves, otherwise"//&
+             "select another option for SURFBAND_SOURCE.")
+      endif
       if (size(CS%WaveNum_Cen).ne.size(forces%stk_wavenumbers)) then
         call MOM_error(FATAL, "Number of wavenumber bands in WW3 does not match that in MOM6. "//&
              "Make sure that STK_BAND_COUPLER in MOM6 input is equal to the number of bands in "//&

From 971802a7d9b41aee75407e729e3bb17dafdcc0fd Mon Sep 17 00:00:00 2001
From: jiandewang <jiande.wang@noaa.gov>
Date: Mon, 4 Jan 2021 10:39:29 -0500
Subject: [PATCH 05/12] revert back to GFDL main for coupled_driver codes

---
 .../MOM_surface_forcing_gfdl.F90              | 25 +------------------
 config_src/coupled_driver/ocean_model_MOM.F90 |  2 +-
 2 files changed, 2 insertions(+), 25 deletions(-)

diff --git a/config_src/coupled_driver/MOM_surface_forcing_gfdl.F90 b/config_src/coupled_driver/MOM_surface_forcing_gfdl.F90
index 3c8f084b4a..7075fb7c10 100644
--- a/config_src/coupled_driver/MOM_surface_forcing_gfdl.F90
+++ b/config_src/coupled_driver/MOM_surface_forcing_gfdl.F90
@@ -189,12 +189,6 @@ module MOM_surface_forcing_gfdl
                                                             !! ice-shelves, expressed as a coefficient
                                                             !! for divergence damping, as determined
                                                             !! outside of the ocean model [m3 s-1]
-  real, pointer, dimension(:,:)   :: ustk0           => NULL() !<
-  real, pointer, dimension(:,:)   :: vstk0           => NULL() !<
-  real, pointer, dimension(:)     :: stk_wavenumbers => NULL() !<
-  real, pointer, dimension(:,:,:) :: ustkb           => NULL() !<
-  real, pointer, dimension(:,:,:) :: vstkb           => NULL() !<
-
   integer :: xtype                    !< The type of the exchange - REGRID, REDIST or DIRECT
   type(coupler_2d_bc_type) :: fluxes  !< A structure that may contain an array of named fields
                                       !! used for passive tracer fluxes.
@@ -202,7 +196,6 @@ module MOM_surface_forcing_gfdl
                                       !! This flag may be set by the flux-exchange code, based on what
                                       !! the sea-ice model is providing.  Otherwise, the value from
                                       !! the surface_forcing_CS is used.
-  integer :: num_stk_bands            !< Number of Stokes drift bands passed through the coupler
 end type ice_ocean_boundary_type
 
 integer :: id_clock_forcing !< A CPU time clock
@@ -684,7 +677,7 @@ subroutine convert_IOB_to_forces(IOB, forces, index_bounds, Time, G, US, CS, dt_
   real :: mass_eff      ! effective mass of sea ice for rigidity [R Z ~> kg m-2]
   real :: wt1, wt2      ! Relative weights of previous and current values of ustar, ND.
 
-  integer :: i, j, is, ie, js, je, Isq, Ieq, Jsq, Jeq, i0, j0, istk
+  integer :: i, j, is, ie, js, je, Isq, Ieq, Jsq, Jeq, i0, j0
   integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB, isr, ier, jsr, jer
   integer :: isc_bnd, iec_bnd, jsc_bnd, jec_bnd
 
@@ -725,9 +718,6 @@ subroutine convert_IOB_to_forces(IOB, forces, index_bounds, Time, G, US, CS, dt_
        (associated(IOB%mass_berg) .and. (.not. associated(forces%mass_berg))) ) &
     call allocate_mech_forcing(G, forces, iceberg=.true.)
 
-  if ( associated(IOB%ustk0) ) &
-    call allocate_mech_forcing(G, forces, waves=.true., num_stk_bands=IOB%num_stk_bands)
-
   if (associated(IOB%ice_rigidity)) then
     rigidity_at_h(:,:) = 0.0
     call safe_alloc_ptr(forces%rigidity_ice_u,IsdB,IedB,jsd,jed)
@@ -787,19 +777,6 @@ subroutine convert_IOB_to_forces(IOB, forces, index_bounds, Time, G, US, CS, dt_
       forces%ustar(i,j) = wt1*forces%ustar(i,j) + wt2*ustar_tmp(i,j)
     enddo ; enddo
   endif
-  forces%stk_wavenumbers(:) = IOB%stk_wavenumbers
-  do j=js,je; do i=is,ie
-     forces%ustk0(i,j) = IOB%ustk0(i-I0,j-J0) ! How to be careful here that the domains are right?
-     forces%vstk0(i,j) = IOB%vstk0(i-I0,j-J0)
-  enddo ; enddo
-  call pass_vector(forces%ustk0,forces%vstk0, G%domain )
-  do istk = 1,IOB%num_stk_bands
-    do j=js,je; do i=is,ie
-      forces%ustkb(i,j,istk) = IOB%ustkb(i-I0,j-J0,istk)
-      forces%vstkb(i,j,istk) = IOB%vstkb(i-I0,j-J0,istk)
-    enddo; enddo
-    call pass_vector(forces%ustkb(:,:,istk),forces%vstkb(:,:,istk), G%domain )
-  enddo
 
   ! Find the net mass source in the input forcing without other adjustments.
   if (CS%approx_net_mass_src .and. associated(forces%net_mass_src)) then
diff --git a/config_src/coupled_driver/ocean_model_MOM.F90 b/config_src/coupled_driver/ocean_model_MOM.F90
index 6cb358cdcb..082099158c 100644
--- a/config_src/coupled_driver/ocean_model_MOM.F90
+++ b/config_src/coupled_driver/ocean_model_MOM.F90
@@ -565,7 +565,7 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, time_start_upda
     ! For now, the waves are only updated on the thermodynamics steps, because that is where
     ! the wave intensities are actually used to drive mixing.  At some point, the wave updates
     ! might also need to become a part of the ocean dynamics, according to B. Reichl.
-    call Update_Surface_Waves(OS%grid, OS%GV, OS%US, OS%time, ocean_coupling_time_step, OS%waves, OS%forces)
+    call Update_Surface_Waves(OS%grid, OS%GV, OS%US, OS%time, ocean_coupling_time_step, OS%waves)
   endif
 
   if ((OS%nstep==0) .and. (OS%nstep_thermo==0)) then ! This is the first call to update_ocean_model.

From 5f5871c1df3e95201292039e3247439f2ff48c17 Mon Sep 17 00:00:00 2001
From: "Jessica.Meixner" <Jessica.Meixner@noaa.gov>
Date: Tue, 5 Jan 2021 20:20:02 +0000
Subject: [PATCH 06/12] Addressing reviewers comments: * document wave coupling
  variables * switch ii jj loops for computational efficiency * clarify error
 message * fix indexing and style

---
 .../mom_surface_forcing_nuopc.F90             | 14 +++--
 src/core/MOM_forcing_type.F90                 | 18 +++---
 src/user/MOM_wave_interface.F90               | 60 +++++++++----------
 3 files changed, 50 insertions(+), 42 deletions(-)

diff --git a/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90 b/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90
index 7b4e33a56a..585117c3ba 100644
--- a/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90
+++ b/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90
@@ -183,11 +183,15 @@ module MOM_surface_forcing_nuopc
                                                               !! ice-shelves, expressed as a coefficient
                                                               !! for divergence damping, as determined
                                                               !! outside of the ocean model in [m3/s]
-  real, pointer, dimension(:,:)   :: ustk0           => NULL() !<
-  real, pointer, dimension(:,:)   :: vstk0           => NULL() !<
-  real, pointer, dimension(:)     :: stk_wavenumbers => NULL() !<
-  real, pointer, dimension(:,:,:) :: ustkb           => NULL() !<
-  real, pointer, dimension(:,:,:) :: vstkb           => NULL() !<
+  real, pointer, dimension(:,:)   :: ustk0           => NULL() !< Surface Stokes drift, zonal [m/s]
+  real, pointer, dimension(:,:)   :: vstk0           => NULL() !< Surface Stokes drift, meridional [m/s]
+  real, pointer, dimension(:)     :: stk_wavenumbers => NULL() !< The central wave number of Stokes bands [rad/m]
+  real, pointer, dimension(:,:,:) :: ustkb           => NULL() !< Stokes Drift spectrum, zonal [m/s]
+                                                               !! Horizontal  - u points
+                                                               !! 3rd dimension - wavenumber
+  real, pointer, dimension(:,:,:) :: vstkb           => NULL() !< Stokes Drift spectrum, meridional [m/s]
+                                                               !! Horizontal  - v points
+                                                               !! 3rd dimension - wavenumber
   integer :: num_stk_bands            !< Number of Stokes drift bands passed through the coupler
   integer :: xtype                                            !< The type of the exchange - REGRID, REDIST or DIRECT
   type(coupler_2d_bc_type)      :: fluxes                     !< A structure that may contain an array of
diff --git a/src/core/MOM_forcing_type.F90 b/src/core/MOM_forcing_type.F90
index 6720414b2b..a135107025 100644
--- a/src/core/MOM_forcing_type.F90
+++ b/src/core/MOM_forcing_type.F90
@@ -250,13 +250,17 @@ module MOM_forcing_type
                                 !! reset to zero at the driver level when appropriate.
 
   real, pointer, dimension(:,:) :: &
-       ustk0 => NULL(), &
-       vstk0 => NULL()
+       ustk0 => NULL(), &       !< Surface Stokes drift, zonal [m/s]
+       vstk0 => NULL()          !< Surface Stokes drift, meridional [m/s]
   real, pointer, dimension(:) :: &
-       stk_wavenumbers => NULL()
+       stk_wavenumbers => NULL() !< The central wave number of Stokes bands [rad/m]
   real, pointer, dimension(:,:,:) :: &
-       ustkb => NULL(), &
-       vstkb => NULL()
+       ustkb => NULL(), &       !< Stokes Drift spectrum, zonal [m/s]
+                                !! Horizontal - u points
+                                !! 3rd dimension - wavenumber
+       vstkb => NULL()          !< Stokes Drift spectrum, meridional [m/s]
+                                !! Horizontal - v points
+                                !! 3rd dimension - wavenumber
 
   logical :: initialized = .false. !< This indicates whether the appropriate arrays have been initialized.
 end type mech_forcing
@@ -3041,12 +3045,12 @@ subroutine allocate_mech_forcing_by_group(G, forces, stress, ustar, shelf, &
     forces%stk_wavenumbers(:) = 0.0    
     allocate(forces%ustkb(isd:ied,jsd:jed,num_stk_bands)) 
     forces%ustkb(isd:ied,jsd:jed,:) = 0.0
-    endif; endif; endif
+  endif ; endif ; endif
 
   if (present(waves)) then; if (waves) then; if (.not.associated(forces%vstkb)) then    
     allocate(forces%vstkb(isd:ied,jsd:jed,num_stk_bands))
     forces%vstkb(isd:ied,jsd:jed,:) = 0.0
-  endif; endif; endif
+  endif ; endif ; endif
 
 end subroutine allocate_mech_forcing_by_group
 
diff --git a/src/user/MOM_wave_interface.F90 b/src/user/MOM_wave_interface.F90
index e79f7b65ef..6352774fa7 100644
--- a/src/user/MOM_wave_interface.F90
+++ b/src/user/MOM_wave_interface.F90
@@ -467,10 +467,10 @@ subroutine Update_Surface_Waves(G, GV, US, Day, dt, CS, forces)
       call Surface_Bands_by_data_override(day_center, G, GV, US, CS)
     elseif (DataSource==Coupler) then
       if (.not.present(FORCES)) then 
-        call MOM_error(FATAL,"The code cannot be run with the options  "//&
-             "SURFBAND_SOURCE = COUPLER for with this driver.  If you are using a "//&
-             "wave coupled driver then check the call to update_surface_waves, otherwise"//&
-             "select another option for SURFBAND_SOURCE.")
+        call MOM_error(FATAL,"The option SURFBAND = COUPLER can not be used with "//&
+             "this driver. If you are using a coupled driver with a wave model then "//& 
+             "check the arguments in the subroutine call to Update_Surface_Waves, "//&
+             "otherwise select another option for SURFBAND_SOURCE.")
       endif
       if (size(CS%WaveNum_Cen).ne.size(forces%stk_wavenumbers)) then
         call MOM_error(FATAL, "Number of wavenumber bands in WW3 does not match that in MOM6. "//&
@@ -481,13 +481,13 @@ subroutine Update_Surface_Waves(G, GV, US, Day, dt, CS, forces)
       do b=1,CS%NumBands
         CS%WaveNum_Cen(b) = forces%stk_wavenumbers(b)
         !Interpolate from a grid to c grid
-        do II=G%iscB,G%iecB
-          do jj=G%jsc,G%jec
+        do jj=G%jsc,G%jec
+          do II=G%iscB,G%iecB
             CS%STKx0(II,jj,b) = 0.5*(forces%UStkb(ii,jj,b)+forces%UStkb(ii+1,jj,b))
           enddo
         enddo
-        do ii=G%isc,G%iec
-          do JJ=G%jscB, G%jecB
+        do JJ=G%jscB, G%jecB
+          do ii=G%isc,G%iec
             CS%STKY0(ii,JJ,b) = 0.5*(forces%VStkb(ii,jj,b)+forces%VStkb(ii,jj+1,b))
           enddo
         enddo
@@ -495,13 +495,13 @@ subroutine Update_Surface_Waves(G, GV, US, Day, dt, CS, forces)
       enddo
     elseif (DataSource==Input) then
       do b=1,CS%NumBands
-        do II=G%isdB,G%iedB
-          do jj=G%jsd,G%jed
+        do jj=G%jsd,G%jed
+          do II=G%isdB,G%iedB
             CS%STKx0(II,jj,b) = CS%PrescribedSurfStkX(b)
           enddo
         enddo
-        do ii=G%isd,G%ied
-          do JJ=G%jsdB, G%jedB
+        do JJ=G%jsdB, G%jedB
+          do ii=G%isd,G%ied
             CS%STKY0(ii,JJ,b) = CS%PrescribedSurfStkY(b)
           enddo
         enddo
@@ -537,8 +537,8 @@ subroutine Update_Stokes_Drift(G, GV, US, CS, h, ustar)
   !    Computing mid-point value from surface value and decay wavelength
   if (WaveMethod==TESTPROF) then
     DecayScale = 4.*PI / TP_WVL !4pi
-    do II = G%isdB,G%iedB
-      do jj = G%jsd,G%jed
+    do jj = G%jsd,G%jed
+      do II = G%isdB,G%iedB
         IIm1 = max(1,II-1)
         Bottom = 0.0
         MidPoint = 0.0
@@ -550,8 +550,8 @@ subroutine Update_Stokes_Drift(G, GV, US, CS, h, ustar)
         enddo
       enddo
     enddo
-    do ii = G%isd,G%ied
-      do JJ = G%jsdB,G%jedB
+    do JJ = G%jsdB,G%jedB
+      do ii = G%isd,G%ied
         JJm1 = max(1,JJ-1)
         Bottom = 0.0
         MidPoint = 0.0
@@ -572,8 +572,8 @@ subroutine Update_Stokes_Drift(G, GV, US, CS, h, ustar)
     CS%Us0_x(:,:) = 0.0
     CS%Us0_y(:,:) = 0.0
     ! Computing X direction Stokes drift
-    do II = G%isdB,G%iedB
-      do jj = G%jsd,G%jed
+    do jj = G%jsd,G%jed
+      do II = G%isdB,G%iedB
         ! 1. First compute the surface Stokes drift
         !    by integrating over the partitionas.
         do b = 1,CS%NumBands
@@ -630,8 +630,8 @@ subroutine Update_Stokes_Drift(G, GV, US, CS, h, ustar)
       enddo
     enddo
     ! Computing Y direction Stokes drift
-    do ii = G%isd,G%ied
-      do JJ = G%jsdB,G%jedB
+    do JJ = G%jsdB,G%jedB
+      do ii = G%isd,G%ied
         ! Compute the surface values.
         do b = 1,CS%NumBands
           if (PartitionMode==0) then
@@ -688,8 +688,8 @@ subroutine Update_Stokes_Drift(G, GV, US, CS, h, ustar)
     enddo
   elseif (WaveMethod==DHH85) then
     if (.not.(StaticWaves .and. DHH85_is_set)) then
-      do II = G%isdB,G%iedB
-        do jj = G%jsd,G%jed
+      do jj = G%jsd,G%jed
+        do II = G%isdB,G%iedB
           bottom = 0.0
           do kk = 1,G%ke
             Top = Bottom
@@ -706,8 +706,8 @@ subroutine Update_Stokes_Drift(G, GV, US, CS, h, ustar)
           enddo
         enddo
       enddo
-      do ii = G%isd,G%ied
-        do JJ = G%jsdB,G%jedB
+      do JJ = G%jsdB,G%jedB
+        do ii = G%isd,G%ied
           Bottom = 0.0
           do kk=1, G%ke
             Top = Bottom
@@ -732,13 +732,13 @@ subroutine Update_Stokes_Drift(G, GV, US, CS, h, ustar)
     endif
   else! Keep this else, fallback to 0 Stokes drift
     do kk= 1,G%ke
-      do II = G%isdB,G%iedB
-        do jj = G%jsd,G%jed
+      do jj = G%jsd,G%jed
+        do II = G%isdB,G%iedB
           CS%Us_x(II,jj,kk) = 0.
         enddo
       enddo
-      do ii = G%isd,G%ied
-        do JJ = G%jsdB,G%jedB
+      do JJ = G%jsdB,G%jedB
+        do ii = G%isd,G%ied
           CS%Us_y(ii,JJ,kk) = 0.
         enddo
       enddo
@@ -748,8 +748,8 @@ subroutine Update_Stokes_Drift(G, GV, US, CS, h, ustar)
   ! Turbulent Langmuir number is computed here and available to use anywhere.
   ! SL Langmuir number requires mixing layer depth, and therefore is computed
   ! in the routine it is needed by (e.g. KPP or ePBL).
-  do ii = G%isc,G%iec
-    do jj = G%jsc, G%jec
+  do jj = G%jsc, G%jec
+    do ii = G%isc,G%iec
       Top = h(ii,jj,1)*GV%H_to_Z
       call get_Langmuir_Number( La, G, GV, US, Top, ustar(ii,jj), ii, jj, &
              H(ii,jj,:),Override_MA=.false.,WAVES=CS)

From 633dc9241b698c3512a9c09eb58f3f12897ced8a Mon Sep 17 00:00:00 2001
From: jiandewang <jiande.wang@noaa.gov>
Date: Wed, 6 Jan 2021 14:06:04 -0500
Subject: [PATCH 07/12] add documentation for mom_cap_methods remove duplicated
 declearation of function chkerr, change case to match

---
 config_src/nuopc_driver/mom_cap.F90         | 13 +---
 config_src/nuopc_driver/mom_cap_methods.F90 | 69 ++++++++++-----------
 config_src/nuopc_driver/mom_cap_time.F90    | 13 +---
 config_src/nuopc_driver/time_utils.F90      | 13 +---
 4 files changed, 36 insertions(+), 72 deletions(-)

diff --git a/config_src/nuopc_driver/mom_cap.F90 b/config_src/nuopc_driver/mom_cap.F90
index 4bfc9e1cb6..f0ce8720bb 100644
--- a/config_src/nuopc_driver/mom_cap.F90
+++ b/config_src/nuopc_driver/mom_cap.F90
@@ -36,6 +36,7 @@ module MOM_cap_mod
 use MOM_ocean_model_nuopc,    only: get_ocean_grid, get_eps_omesh
 use MOM_cap_time,             only: AlarmInit
 use MOM_cap_methods,          only: mom_import, mom_export, mom_set_geomtype, state_diagnose
+use MOM_cap_methods,          only: ChkErr
 #ifdef CESMCOUPLED
 use shr_file_mod,             only: shr_file_setLogUnit, shr_file_getLogUnit
 #endif
@@ -2049,18 +2050,6 @@ subroutine shr_file_getLogUnit(nunit)
 end subroutine shr_file_getLogUnit
 #endif
 
-logical function chkerr(rc, line, file)
-  integer, intent(in) :: rc
-  integer, intent(in) :: line
-  character(len=*), intent(in) :: file
-  integer :: lrc
-  chkerr = .false.
-  lrc = rc
-  if (ESMF_LogFoundError(rcToCheck=lrc, msg=ESMF_LOGERR_PASSTHRU, line=line, file=file)) then
-     chkerr = .true.
-  endif
-end function chkerr
-
 !>
 !! @page nuopc_cap NUOPC Cap
 !! @author Fei Liu (fei.liu@gmail.com)
diff --git a/config_src/nuopc_driver/mom_cap_methods.F90 b/config_src/nuopc_driver/mom_cap_methods.F90
index 1d51c1e6dd..39f450d453 100644
--- a/config_src/nuopc_driver/mom_cap_methods.F90
+++ b/config_src/nuopc_driver/mom_cap_methods.F90
@@ -30,6 +30,7 @@ module MOM_cap_methods
 public :: mom_import
 public :: mom_export
 public :: state_diagnose
+public :: ChkErr
 
 private :: State_getImport
 private :: State_setExport
@@ -251,9 +252,9 @@ subroutine mom_import(ocean_public, ocean_grid, importState, ice_ocean_boundary,
 
 
   !----
-  ! Partitioned Stokes Drift Components 
+  ! Partitioned Stokes Drift Components
   !----
-  if ( associated(ice_ocean_boundary%ustkb) ) then    
+  if ( associated(ice_ocean_boundary%ustkb) ) then
     allocate(stkx1(isc:iec,jsc:jec))
     allocate(stky1(isc:iec,jsc:jec))
     allocate(stkx2(isc:iec,jsc:jec))
@@ -765,15 +766,18 @@ subroutine State_SetExport(state, fldname, isc, iec, jsc, jec, input, ocean_grid
 
 end subroutine State_SetExport
 
+!> This subroutine writes the minimum, maximum and sum of each field
+!! contained within an ESMF state.
 subroutine state_diagnose(State, string, rc)
 
   ! ----------------------------------------------
   ! Diagnose status of State
   ! ----------------------------------------------
 
-  type(ESMF_State), intent(in)  :: state
-  character(len=*), intent(in)  :: string
-  integer         , intent(out) :: rc
+  type(ESMF_State), intent(in)    :: state  !< An ESMF State
+  character(len=*), intent(in)    :: string !< A string indicating whether the State is an
+                                            !! import or export State
+  integer         , intent(out)   :: rc     !< Return code
 
   ! local variables
   integer                         :: i,j,n
@@ -787,19 +791,19 @@ subroutine state_diagnose(State, string, rc)
   ! ----------------------------------------------
 
   call ESMF_StateGet(state, itemCount=fieldCount, rc=rc)
-  if (chkerr(rc,__LINE__,u_FILE_u)) return
+  if (ChkErr(rc,__LINE__,u_FILE_u)) return
   allocate(lfieldnamelist(fieldCount))
 
   call ESMF_StateGet(state, itemNameList=lfieldnamelist, rc=rc)
-  if (chkerr(rc,__LINE__,u_FILE_u)) return
+  if (ChkErr(rc,__LINE__,u_FILE_u)) return
 
   do n = 1, fieldCount
 
      call ESMF_StateGet(state, itemName=lfieldnamelist(n), field=lfield, rc=rc)
-     if (chkerr(rc,__LINE__,u_FILE_u)) return
+     if (ChkErr(rc,__LINE__,u_FILE_u)) return
 
      call field_getfldptr(lfield, fldptr1=dataPtr1d, fldptr2=dataPtr2d, rank=lrank, rc=rc)
-     if (chkerr(rc,__LINE__,u_FILE_u)) return
+     if (ChkErr(rc,__LINE__,u_FILE_u)) return
 
      if (lrank == 0) then
         ! no local data
@@ -829,23 +833,16 @@ subroutine state_diagnose(State, string, rc)
 
 end subroutine state_diagnose
 
-!===============================================================================
-
+!> Obtain a pointer to a rank 1 or rank 2 ESMF field
 subroutine field_getfldptr(field, fldptr1, fldptr2, rank, abort, rc)
 
-  ! ----------------------------------------------
-  ! for a field, determine rank and return fldptr1 or fldptr2
-  ! abort is true by default and will abort if fldptr is not yet allocated in field
-  ! rank returns 0, 1, or 2.  0 means fldptr not allocated and abort=false
-  ! ----------------------------------------------
-
   ! input/output variables
-  type(ESMF_Field)  , intent(in)              :: field
-  real(ESMF_KIND_R8), pointer , intent(inout), optional :: fldptr1(:)
-  real(ESMF_KIND_R8), pointer , intent(inout), optional :: fldptr2(:,:)
-  integer           , intent(out)  , optional :: rank
-  logical           , intent(in)   , optional :: abort
-  integer           , intent(out)  , optional :: rc
+  type(ESMF_Field)  , intent(in)                        :: field        !< An ESMF field
+  real(ESMF_KIND_R8), pointer , intent(inout), optional :: fldptr1(:)   !< A pointer to a rank 1 ESMF field
+  real(ESMF_KIND_R8), pointer , intent(inout), optional :: fldptr2(:,:) !< A pointer to a rank 2 ESMF field
+  integer           , intent(out)            , optional :: rank         !< Field rank
+  logical           , intent(in)             , optional :: abort        !< Abort code
+  integer           , intent(out)            , optional :: rc           !< Return code
 
   ! local variables
   type(ESMF_GeomType_Flag)    :: geomtype
@@ -872,7 +869,7 @@ subroutine field_getfldptr(field, fldptr1, fldptr2, rank, abort, rc)
   lrank = -99
 
   call ESMF_FieldGet(field, status=status, rc=rc)
-  if (chkerr(rc,__LINE__,u_FILE_u)) return
+  if (ChkErr(rc,__LINE__,u_FILE_u)) return
 
   if (status /= ESMF_FIELDSTATUS_COMPLETE) then
      lrank = 0
@@ -886,20 +883,20 @@ subroutine field_getfldptr(field, fldptr1, fldptr2, rank, abort, rc)
   else
 
      call ESMF_FieldGet(field, geomtype=geomtype, rc=rc)
-     if (chkerr(rc,__LINE__,u_FILE_u)) return
+     if (ChkErr(rc,__LINE__,u_FILE_u)) return
 
      if (geomtype == ESMF_GEOMTYPE_GRID) then
         call ESMF_FieldGet(field, rank=lrank, rc=rc)
-        if (chkerr(rc,__LINE__,u_FILE_u)) return
+        if (ChkErr(rc,__LINE__,u_FILE_u)) return
      elseif (geomtype == ESMF_GEOMTYPE_MESH) then
         call ESMF_FieldGet(field, rank=lrank, rc=rc)
-        if (chkerr(rc,__LINE__,u_FILE_u)) return
+        if (ChkErr(rc,__LINE__,u_FILE_u)) return
         call ESMF_FieldGet(field, mesh=lmesh, rc=rc)
-        if (chkerr(rc,__LINE__,u_FILE_u)) return
+        if (ChkErr(rc,__LINE__,u_FILE_u)) return
         call ESMF_MeshGet(lmesh, numOwnedNodes=nnodes, numOwnedElements=nelements, rc=rc)
-        if (chkerr(rc,__LINE__,u_FILE_u)) return
+        if (ChkErr(rc,__LINE__,u_FILE_u)) return
         if (nnodes == 0 .and. nelements == 0) lrank = 0
-     else  
+     else
         call ESMF_LogWrite(trim(subname)//": ERROR geomtype not supported ", &
              ESMF_LOGMSG_INFO)
         rc = ESMF_FAILURE
@@ -917,7 +914,7 @@ subroutine field_getfldptr(field, fldptr1, fldptr2, rank, abort, rc)
            return
         endif
         call ESMF_FieldGet(field, farrayPtr=fldptr1, rc=rc)
-        if (chkerr(rc,__LINE__,u_FILE_u)) return
+        if (ChkErr(rc,__LINE__,u_FILE_u)) return
      elseif (lrank == 2) then
         if (.not.present(fldptr2)) then
            call ESMF_LogWrite(trim(subname)//": ERROR missing rank=2 array ", &
@@ -926,7 +923,7 @@ subroutine field_getfldptr(field, fldptr1, fldptr2, rank, abort, rc)
            return
         endif
         call ESMF_FieldGet(field, farrayPtr=fldptr2, rc=rc)
-        if (chkerr(rc,__LINE__,u_FILE_u)) return
+        if (ChkErr(rc,__LINE__,u_FILE_u)) return
      else
         call ESMF_LogWrite(trim(subname)//": ERROR in rank ", &
              ESMF_LOGMSG_ERROR, line=__LINE__, file=u_FILE_u)
@@ -942,16 +939,16 @@ subroutine field_getfldptr(field, fldptr1, fldptr2, rank, abort, rc)
 
 end subroutine field_getfldptr
 
-logical function chkerr(rc, line, file)
+logical function ChkErr(rc, line, file)
   integer, intent(in) :: rc
   integer, intent(in) :: line
   character(len=*), intent(in) :: file
   integer :: lrc
-  chkerr = .false.
+  ChkErr = .false.
   lrc = rc
   if (ESMF_LogFoundError(rcToCheck=lrc, msg=ESMF_LOGERR_PASSTHRU, line=line, file=file)) then
-     chkerr = .true.
+     ChkErr = .true.
   endif
-end function chkerr
+end function ChkErr
 
 end module MOM_cap_methods
diff --git a/config_src/nuopc_driver/mom_cap_time.F90 b/config_src/nuopc_driver/mom_cap_time.F90
index aa1a6c7072..7f210bda71 100644
--- a/config_src/nuopc_driver/mom_cap_time.F90
+++ b/config_src/nuopc_driver/mom_cap_time.F90
@@ -16,6 +16,7 @@ module MOM_cap_time
 use ESMF                  , only : ESMF_RC_ARG_BAD
 use ESMF                  , only : operator(<), operator(/=), operator(+), operator(-), operator(*) , operator(>=)
 use ESMF                  , only : operator(<=), operator(>), operator(==)
+use MOM_cap_methods       , only : ChkErr
 
 implicit none; private
 
@@ -336,16 +337,4 @@ subroutine date2ymd (date, year, month, day)
 
 end subroutine date2ymd
 
-logical function chkerr(rc, line, file)
-  integer, intent(in) :: rc
-  integer, intent(in) :: line
-  character(len=*), intent(in) :: file
-  integer :: lrc
-  chkerr = .false.
-  lrc = rc
-  if (ESMF_LogFoundError(rcToCheck=lrc, msg=ESMF_LOGERR_PASSTHRU, line=line, file=file)) then
-     chkerr = .true.
-  endif
-end function chkerr
-
 end module MOM_cap_time
diff --git a/config_src/nuopc_driver/time_utils.F90 b/config_src/nuopc_driver/time_utils.F90
index db056f4bf5..81efcd2765 100644
--- a/config_src/nuopc_driver/time_utils.F90
+++ b/config_src/nuopc_driver/time_utils.F90
@@ -14,6 +14,7 @@ module time_utils_mod
 use ESMF,               only: ESMF_Time, ESMF_TimeGet, ESMF_LogFoundError
 use ESMF,               only: ESMF_LOGERR_PASSTHRU,ESMF_TimeInterval
 use ESMF,               only: ESMF_TimeIntervalGet, ESMF_TimeSet, ESMF_SUCCESS
+use MOM_cap_methods,    only: ChkErr
 
 implicit none; private
 
@@ -160,16 +161,4 @@ function string_to_date(string, rc)
 
 end function string_to_date
 
-logical function chkerr(rc, line, file)
-  integer, intent(in) :: rc
-  integer, intent(in) :: line
-  character(len=*), intent(in) :: file
-  integer :: lrc
-  chkerr = .false.
-  lrc = rc
-  if (ESMF_LogFoundError(rcToCheck=lrc, msg=ESMF_LOGERR_PASSTHRU, line=line, file=file)) then
-     chkerr = .true.
-  endif
-end function chkerr
-
 end module time_utils_mod

From 91282c10620abdb432987891f97e030f2d2f16f5 Mon Sep 17 00:00:00 2001
From: "denise.worthen" <Denise.Worthen@noaa.gov>
Date: Fri, 15 Jan 2021 11:41:46 -0500
Subject: [PATCH 08/12] add brandon's halo update fix for LI_2016

---
 config_src/nuopc_driver/mom_surface_forcing_nuopc.F90 | 2 ++
 1 file changed, 2 insertions(+)

diff --git a/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90 b/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90
index 7b4e33a56a..4c9c872a5a 100644
--- a/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90
+++ b/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90
@@ -790,6 +790,7 @@ subroutine convert_IOB_to_forces(IOB, forces, index_bounds, Time, G, US, CS)
       endif
       forces%ustar(i,j) = sqrt(gustiness*Irho0 + Irho0*tau_mag)
     enddo ; enddo
+    call pass_vector(forces%taux, forces%tauy, G%Domain, halo=1)
 
   elseif (wind_stagger == AGRID) then
     call pass_vector(taux_at_h, tauy_at_h, G%Domain, To_All+Omit_Corners, stagger=AGRID, halo=1)
@@ -816,6 +817,7 @@ subroutine convert_IOB_to_forces(IOB, forces, index_bounds, Time, G, US, CS)
       forces%ustar(i,j) = sqrt(gustiness*Irho0 + Irho0 * G%mask2dT(i,j) * &
                                sqrt(taux_at_h(i,j)**2 + tauy_at_h(i,j)**2))
     enddo ; enddo
+    call pass_vector(forces%taux, forces%tauy, G%Domain, halo=1)
 
   else ! C-grid wind stresses.
     if (G%symmetric) &

From 172025dd87ef1042d652286e0eccd0f56dc56fd9 Mon Sep 17 00:00:00 2001
From: jiandewang <jiande.wang@noaa.gov>
Date: Wed, 20 Jan 2021 12:24:14 -0500
Subject: [PATCH 09/12] add documentation for function ChkErr

---
 config_src/nuopc_driver/mom_cap_methods.F90 | 6 +++---
 1 file changed, 3 insertions(+), 3 deletions(-)

diff --git a/config_src/nuopc_driver/mom_cap_methods.F90 b/config_src/nuopc_driver/mom_cap_methods.F90
index 39f450d453..4946fc9927 100644
--- a/config_src/nuopc_driver/mom_cap_methods.F90
+++ b/config_src/nuopc_driver/mom_cap_methods.F90
@@ -940,9 +940,9 @@ subroutine field_getfldptr(field, fldptr1, fldptr2, rank, abort, rc)
 end subroutine field_getfldptr
 
 logical function ChkErr(rc, line, file)
-  integer, intent(in) :: rc
-  integer, intent(in) :: line
-  character(len=*), intent(in) :: file
+  integer, intent(in) :: rc            !< return code to check
+  integer, intent(in) :: line          !< Integer source line number
+  character(len=*), intent(in) :: file !< User-provided source file name
   integer :: lrc
   ChkErr = .false.
   lrc = rc

From c05b042d9e9e16efc6afb36a0e8d8a888d82bde2 Mon Sep 17 00:00:00 2001
From: Alistair Adcroft <Alistair.Adcroft@noaa.gov>
Date: Wed, 20 Jan 2021 14:54:46 -0500
Subject: [PATCH 10/12] Added API documentation for chkerr()

---
 config_src/nuopc_driver/mom_cap_methods.F90 | 1 +
 1 file changed, 1 insertion(+)

diff --git a/config_src/nuopc_driver/mom_cap_methods.F90 b/config_src/nuopc_driver/mom_cap_methods.F90
index 4946fc9927..78014f1c63 100644
--- a/config_src/nuopc_driver/mom_cap_methods.F90
+++ b/config_src/nuopc_driver/mom_cap_methods.F90
@@ -939,6 +939,7 @@ subroutine field_getfldptr(field, fldptr1, fldptr2, rank, abort, rc)
 
 end subroutine field_getfldptr
 
+!> Returns true if ESMF_LogFoundError() determines that rc is an error code. Otherwise false.
 logical function ChkErr(rc, line, file)
   integer, intent(in) :: rc            !< return code to check
   integer, intent(in) :: line          !< Integer source line number

From ea5e2458fbd8b39f07e7598195d39f54c352a78a Mon Sep 17 00:00:00 2001
From: Alistair Adcroft <Alistair.Adcroft@noaa.gov>
Date: Wed, 20 Jan 2021 18:02:52 -0500
Subject: [PATCH 11/12] Fixed white space and line length

---
 .../nuopc_driver/mom_surface_forcing_nuopc.F90   |  2 +-
 src/core/MOM_forcing_type.F90                    | 16 ++++++++--------
 src/tracer/MOM_lateral_boundary_diffusion.F90    | 12 ++++++------
 src/user/MOM_wave_interface.F90                  |  4 ++--
 4 files changed, 17 insertions(+), 17 deletions(-)

diff --git a/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90 b/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90
index 585117c3ba..689a9f0f4a 100644
--- a/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90
+++ b/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90
@@ -862,7 +862,7 @@ subroutine convert_IOB_to_forces(IOB, forces, index_bounds, Time, G, US, CS)
       enddo; enddo
       call pass_vector(forces%ustkb(:,:,istk),forces%vstkb(:,:,istk), G%domain )
     enddo
-  endif 
+  endif
 
   ! sea ice related dynamic fields
   if (associated(IOB%ice_rigidity)) then
diff --git a/src/core/MOM_forcing_type.F90 b/src/core/MOM_forcing_type.F90
index a135107025..f0cc8f553c 100644
--- a/src/core/MOM_forcing_type.F90
+++ b/src/core/MOM_forcing_type.F90
@@ -3007,8 +3007,8 @@ subroutine allocate_mech_forcing_by_group(G, forces, stress, ustar, shelf, &
   logical, optional,     intent(in) :: shelf   !< If present and true, allocate forces for ice-shelf
   logical, optional,     intent(in) :: press   !< If present and true, allocate p_surf and related fields
   logical, optional,     intent(in) :: iceberg !< If present and true, allocate forces for icebergs
-  logical, optional,     intent(in) :: waves   !< If present and true, allocate wave fields              
-  integer, optional,     intent(in) :: num_stk_bands !< Number of Stokes bands to allocate 
+  logical, optional,     intent(in) :: waves   !< If present and true, allocate wave fields
+  integer, optional,     intent(in) :: num_stk_bands !< Number of Stokes bands to allocate
 
   ! Local variables
   integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB
@@ -3038,16 +3038,16 @@ subroutine allocate_mech_forcing_by_group(G, forces, stress, ustar, shelf, &
   !These fields should only be allocated when waves
   call myAlloc(forces%ustk0,isd,ied,jsd,jed, waves)
   call myAlloc(forces%vstk0,isd,ied,jsd,jed, waves)
-  if (present(waves)) then; if (waves) then; if (.not.associated(forces%ustkb)) then       
+  if (present(waves)) then; if (waves) then; if (.not.associated(forces%ustkb)) then
     if (.not.(present(num_stk_bands))) call MOM_error(FATAL,"Requested to &
-      initialize with waves, but no waves are present.") 
-    allocate(forces%stk_wavenumbers(num_stk_bands)) 
-    forces%stk_wavenumbers(:) = 0.0    
-    allocate(forces%ustkb(isd:ied,jsd:jed,num_stk_bands)) 
+      initialize with waves, but no waves are present.")
+    allocate(forces%stk_wavenumbers(num_stk_bands))
+    forces%stk_wavenumbers(:) = 0.0
+    allocate(forces%ustkb(isd:ied,jsd:jed,num_stk_bands))
     forces%ustkb(isd:ied,jsd:jed,:) = 0.0
   endif ; endif ; endif
 
-  if (present(waves)) then; if (waves) then; if (.not.associated(forces%vstkb)) then    
+  if (present(waves)) then; if (waves) then; if (.not.associated(forces%vstkb)) then
     allocate(forces%vstkb(isd:ied,jsd:jed,num_stk_bands))
     forces%vstkb(isd:ied,jsd:jed,:) = 0.0
   endif ; endif ; endif
diff --git a/src/tracer/MOM_lateral_boundary_diffusion.F90 b/src/tracer/MOM_lateral_boundary_diffusion.F90
index 2b7a5646cc..01d8e4163b 100644
--- a/src/tracer/MOM_lateral_boundary_diffusion.F90
+++ b/src/tracer/MOM_lateral_boundary_diffusion.F90
@@ -579,13 +579,13 @@ subroutine fluxes_layer_method(boundary, ke, hbl_L, hbl_R, h_L, h_R, phi_L, phi_
   type(lbd_CS),          pointer       :: CS       !< Lateral diffusion control structure
                                                       !! the boundary layer
   ! Local variables
-  real, dimension(:), allocatable :: dz_top    !< The LBD z grid to be created                                    [L ~ m]
-  real, dimension(:), allocatable :: phi_L_z   !< Tracer values in the ztop grid (left)                            [conc]
-  real, dimension(:), allocatable :: phi_R_z   !< Tracer values in the ztop grid (right)                           [conc]
-  real, dimension(:), allocatable :: F_layer_z !< Diffusive flux at U- or V-point in the ztop grid [H L2 conc ~> m3 conc]
+  real, dimension(:), allocatable :: dz_top    !< The LBD z grid to be created                                [L ~ m]
+  real, dimension(:), allocatable :: phi_L_z   !< Tracer values in the ztop grid (left)                        [conc]
+  real, dimension(:), allocatable :: phi_R_z   !< Tracer values in the ztop grid (right)                       [conc]
+  real, dimension(:), allocatable :: F_layer_z !< Diffusive flux at U/V-point in the ztop grid [H L2 conc ~> m3 conc]
   real, dimension(ke)             :: h_vel     !< Thicknesses at u- and v-points in the native grid
-                                               !! The harmonic mean is used to avoid zero values       [H ~> m or kg m-2]
-  real    :: khtr_avg                !< Thickness-weighted diffusivity at the u-point                         [m^2 s^-1]
+                                               !! The harmonic mean is used to avoid zero values   [H ~> m or kg m-2]
+  real    :: khtr_avg                !< Thickness-weighted diffusivity at the u-point                      [m^2 s^-1]
                                      !! This is just to remind developers that khtr_avg should be
                                      !! computed once khtr is 3D.
   real    :: htot                    !< Total column thickness [H ~> m or kg m-2]
diff --git a/src/user/MOM_wave_interface.F90 b/src/user/MOM_wave_interface.F90
index 6352774fa7..0780ca5e3d 100644
--- a/src/user/MOM_wave_interface.F90
+++ b/src/user/MOM_wave_interface.F90
@@ -466,9 +466,9 @@ subroutine Update_Surface_Waves(G, GV, US, Day, dt, CS, forces)
     if (DataSource==DATAOVR) then
       call Surface_Bands_by_data_override(day_center, G, GV, US, CS)
     elseif (DataSource==Coupler) then
-      if (.not.present(FORCES)) then 
+      if (.not.present(FORCES)) then
         call MOM_error(FATAL,"The option SURFBAND = COUPLER can not be used with "//&
-             "this driver. If you are using a coupled driver with a wave model then "//& 
+             "this driver. If you are using a coupled driver with a wave model then "//&
              "check the arguments in the subroutine call to Update_Surface_Waves, "//&
              "otherwise select another option for SURFBAND_SOURCE.")
       endif

From e7d0976de4f8db42dab689115d65397792872676 Mon Sep 17 00:00:00 2001
From: Gustavo Marques <gmarques@ucar.edu>
Date: Thu, 21 Jan 2021 09:47:44 -0700
Subject: [PATCH 12/12] Fixes latent heat from fprec and frunoff

This patch fixes a sign bug, in both MCT and NUOPC, when
accounting for the latent heat from fprec and frunnoff.
Following MOM6's definition, both fprec
and frunoff are > 0 into the ocean. Therefore, the latent heat
associated with these terms should be negative.
---
 config_src/mct_driver/mom_surface_forcing_mct.F90    | 12 ++++++------
 .../nuopc_driver/mom_surface_forcing_nuopc.F90       | 10 ++++++----
 2 files changed, 12 insertions(+), 10 deletions(-)

diff --git a/config_src/mct_driver/mom_surface_forcing_mct.F90 b/config_src/mct_driver/mom_surface_forcing_mct.F90
index 82105e040e..ef0527dd1c 100644
--- a/config_src/mct_driver/mom_surface_forcing_mct.F90
+++ b/config_src/mct_driver/mom_surface_forcing_mct.F90
@@ -486,17 +486,17 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, valid_time, G,
 
     ! latent heat flux (W/m^2)
     fluxes%latent(i,j) = 0.0
-    ! contribution from frozen ppt
+    ! contribution from frozen ppt (notice minus sign since fprec is positive into the ocean)
     if (associated(IOB%fprec)) then
-      fluxes%latent(i,j)              = fluxes%latent(i,j) + &
+      fluxes%latent(i,j)              = fluxes%latent(i,j) - &
           IOB%fprec(i-i0,j-j0)*US%W_m2_to_QRZ_T*CS%latent_heat_fusion
-      fluxes%latent_fprec_diag(i,j)   = G%mask2dT(i,j) * IOB%fprec(i-i0,j-j0)*US%W_m2_to_QRZ_T*CS%latent_heat_fusion
+      fluxes%latent_fprec_diag(i,j)   = - G%mask2dT(i,j) * IOB%fprec(i-i0,j-j0)*US%W_m2_to_QRZ_T*CS%latent_heat_fusion
     endif
-    ! contribution from frozen runoff
+    ! contribution from frozen runoff (notice minus sign since rofi_flux is positive into the ocean)
     if (associated(fluxes%frunoff)) then
-      fluxes%latent(i,j)              = fluxes%latent(i,j) + &
+      fluxes%latent(i,j)              = fluxes%latent(i,j) - &
           IOB%rofi_flux(i-i0,j-j0)*US%W_m2_to_QRZ_T*CS%latent_heat_fusion
-      fluxes%latent_frunoff_diag(i,j) = G%mask2dT(i,j) * IOB%rofi_flux(i-i0,j-j0)*US%W_m2_to_QRZ_T*CS%latent_heat_fusion
+      fluxes%latent_frunoff_diag(i,j) = - G%mask2dT(i,j) * IOB%rofi_flux(i-i0,j-j0)*US%W_m2_to_QRZ_T*CS%latent_heat_fusion
     endif
     ! contribution from evaporation
     if (associated(IOB%q_flux)) then
diff --git a/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90 b/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90
index 7c5e8ee6d9..7168823fbc 100644
--- a/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90
+++ b/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90
@@ -497,15 +497,17 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, valid_time, G,
          fluxes%seaice_melt(i,j) = kg_m2_s_conversion * G%mask2dT(i,j) * IOB%seaice_melt(i-i0,j-j0)
 
     fluxes%latent(i,j) = 0.0
+    ! notice minus sign since fprec is positive into the ocean
     if (associated(IOB%fprec)) then
-       fluxes%latent(i,j)            = fluxes%latent(i,j) + &
+       fluxes%latent(i,j)            = fluxes%latent(i,j) - &
            IOB%fprec(i-i0,j-j0)*US%W_m2_to_QRZ_T*CS%latent_heat_fusion
-       fluxes%latent_fprec_diag(i,j) = G%mask2dT(i,j) * IOB%fprec(i-i0,j-j0)*US%W_m2_to_QRZ_T*CS%latent_heat_fusion
+       fluxes%latent_fprec_diag(i,j) = - G%mask2dT(i,j) * IOB%fprec(i-i0,j-j0)*US%W_m2_to_QRZ_T*CS%latent_heat_fusion
     endif
+    ! notice minus sign since frunoff is positive into the ocean
     if (associated(IOB%frunoff)) then
-       fluxes%latent(i,j)              = fluxes%latent(i,j) + &
+       fluxes%latent(i,j)              = fluxes%latent(i,j) - &
            IOB%frunoff(i-i0,j-j0) * US%W_m2_to_QRZ_T * CS%latent_heat_fusion
-       fluxes%latent_frunoff_diag(i,j) = G%mask2dT(i,j) * &
+       fluxes%latent_frunoff_diag(i,j) = - G%mask2dT(i,j) * &
            IOB%frunoff(i-i0,j-j0) * US%W_m2_to_QRZ_T * CS%latent_heat_fusion
     endif
     if (associated(IOB%q_flux)) then