add Cgrid-related fixes for nuopc/cmeps

[DeniseWorthen]

For detailed information about submitting Pull Requests (PRs) to the CICE-Consortium,
please refer to: https://github.com/CICE-Consortium/About-Us/wiki/Resource-Index#information-for-developers

PR checklist

• Short (1 sentence) summary of your PR:
  Updates drivers/nuopc/cmeps for C-grid code changes. See C-grid coupled validation #721
• Developer(s):
  @DeniseWorthen
  @apcraig
  @dabail10
• Suggest PR reviewers from list in the column to the right. @apcraig, @dabail10
• Please copy the PR test results link or provide a summary of testing completed below.
  See below
• How much do the PR code changes differ from the unmodified code?
  • bit for bit
  • different at roundoff level
  • more substantial
• Does this PR create or have dependencies on Icepack or any other models?
  • Yes
  • No
• Does this PR add any new test cases?
  • Yes
  • No
• Is the documentation being updated? ("Documentation" includes information on the wiki or in the .rst files from doc/source/, which are used to create the online technical docs at https://readthedocs.org/projects/cice-consortium-cice/. A test build of the technical docs will be performed as part of the PR testing.)
  • Yes
  • No, does the documentation need to be updated at a later time?
    • Yes
    • No
• Please provide any additional information or relevant details below:

This PR will change baselines even when the C-grid option is not exercised because of the change in cicecore/cicedynB/dynamics/ice_dyn_evp.F90:

      call grid_average_X2Y('S', uocn     , grid_ocn_dynu, uocnU   , 'U')
      call grid_average_X2Y('S', vocn     , grid_ocn_dynv, vocnU   , 'U')
      call grid_average_X2Y('S', ss_tltx  , grid_ocn_dynu, ss_tltxU, 'U')
      call grid_average_X2Y('S', ss_tlty  , grid_ocn_dynv, ss_tltyU, 'U')

These lines determine whether the variables are averaged as though they were fluxes ("F") vs state variables ("S"). The option "S" is correct for these terms (ocean velocities and surface tilts). The pre-Cgrid code is equivalent to averaging the terms as if they were fluxes ("F"). Using the coupled application of UWM, retaining these four lines as "F" gave B4B answers.

A change is also made to cicecore/cicedynB/dynamics/ice_dyn_shared.F90 to initialize E/N variables.

[apcraig]

There seems to be a bug in the implementation which github Actions is picking up.

/Users/runner/cice/cicecore/cicedynB/dynamics/ice_dyn_shared.F90:263:32:
  263 |            iceEmask(i,j,iblk) = c0
      |                                1
Error: Cannot convert REAL(8) to LOGICAL(4) at (1)
/Users/runner/cice/cicecore/cicedynB/dynamics/ice_dyn_shared.F90:279:32:
  279 |            iceNmask(i,j,iblk) = c0
      |                                1
Error: Cannot convert REAL(8) to LOGICAL(4) at (1)
make: *** [/Users/runner/cice/testsuite.macos-latest/conda_macos_smoke_gx3_1x2_run2day.macos-latest/Makefile:178: ice_dyn_shared.o] Error 1

Can you fix this and retest please. thanks!

[apcraig]

Also, it would be great if some of the initialization were moved to other parts of the code. The C-grid forcing/history terms should be initialized in cicecore/cicedynB/general/ice_flux.F90:init_history_dyn. An "if (grid_ice == "CD" .or. grid_ice == "C")" can be added there. The iceEmask and iceNmask can stay in ice_dyn_shared.F90 but need to be initialized to false there.

[dabail10]

I checkout out Denise's code and ran an ice-only (D compset) test with grid_ice == 'B'. It is still hanging with threading. I am going to try turning off threading.

[DeniseWorthen]

@apcraig Sorry about the sloppy setting of the mask variables as c0. I will fix that.

Also, I did originally have the changes in the ice_flux routine; I only moved them because I saw that uvelE/N are initialized in ice_dyn_shared. So, other than the mask variables (logical, .false. on initialization) all of the others should be moved to the history initialization ?

[apcraig]

@DeniseWorthen, no problems. I think the velocity can stay where it is as well. I think we should have the C grid variables initialized, as much as we can, in the same places where the B-grid equivalents are initialized. So init_dyn should have velocities (already there) and iceemask/icenmask near iceumask. I think the other variables go in ice_flux.F90. Let me know if there are questions. Just trying to keep some consistency. thanks!

[dabail10]

This version works with threading off for me! So, there is something likely bad with an OMP loop in the cap. Still hunting.

[apcraig]

@dabail10, I'll wait for your approval to merge.

[dabail10]

I still have some threading issues to sort out, but I don't think we should hold this up.

[dabail10]

I did find a threading bug in ice_import_export.F90 where inst_pres_heigh_lowest needs to be private. I don't think this is the problem though. I wonder though, we don't have the "SCHEDULE" thing in the cap. Does this matter?

else if (State_fldChk(importState, 'inst_pres_height_lowest')) then
   !$OMP PARALLEL DO PRIVATE(iblk,i,j,inst_pres_height_lowest)
   do iblk = 1, nblocks
      do j = 1,ny_block
         do i = 1,nx_block
            inst_pres_height_lowest = aflds(i,j,6,iblk)
            if (inst_pres_height_lowest > 0.0_ESMF_KIND_R8) then
               potT (i,j,iblk) = Tair(i,j,iblk) * (100000._ESMF_KIND_R8/inst_pres_height_lowest)**0.286_ESMF_KIND_R8
            else
               potT (i,j,iblk) = 0.0_ESMF_KIND_R8
            end if
            if (Tair(i,j,iblk) /= 0._ESMF_KIND_R8) then
               rhoa(i,j,iblk) = inst_pres_height_lowest / &
                    (287.058_ESMF_KIND_R8*(1._ESMF_KIND_R8+0.608_ESMF_KIND_R8*Qa(i,j,iblk))*Tair(i,j,iblk))
            else
               rhoa(i,j,iblk) = 1.2_ESMF_KIND_R8
            endif
         end do !i
      end do !j
   end do !iblk
   !$OMP END PARALLEL DO
end if

[apcraig]

The SCHEDULE only affects performance, it won't impact usage.

It was a big pain to debug the OpenMP when I went thru the exercise earlier this year in CICE. One thing I often did is turn off all OpenMP loops in a suspicious subroutine or file. Then I would turn them on one or a few at a time and see if I could identify which one(s) were causing problems. Then I would focus debugging that OpenMP loop. Let me know if I can help.

I will merge this now.