Skip to content

Commit

Permalink
add Cgrid-related fixes for nuopc/cmeps (CICE-Consortium#728)
Browse files Browse the repository at this point in the history 
* replace save_init with step_prep in CICE_RunMod

* fixes for cgrid repro

* remove added haloupdates

* baselines pass with these extra halo updates removed

* change F->S for ocean velocities and tilts

* fix debug failure when grid_ice=C

* compiling in debug mode using -init=snan,arrays requires
initialization of variables

* respond to review comments

* remove inserted whitespace for uvelE,N and vvelE,N

Co-authored-by: apcraig <[email protected]>
Co-authored-by: David Bailey <[email protected]>
Co-authored-by: Elizabeth Hunke <[email protected]>
Co-authored-by: Mariana Vertenstein <[email protected]>
Co-authored-by: Minsuk Ji <[email protected]>
Co-authored-by: Tony Craig <[email protected]>
Co-authored-by: Philippe Blain <[email protected]>
  • Loading branch information
@DeniseWorthen @apcraig
@dabail10 @eclare108213 @MinsukJi-NOAA @phil-blain
8 people authored Jun 23, 2022
1 parent 7705e13 commit 471c010
Show file tree
Hide file tree
Showing 5 changed files with 71 additions and 51 deletions.
7 changes: 5 additions & 2 deletions cicecore/cicedynB/dynamics/ice_dyn_shared.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -192,7 +192,7 @@ subroutine init_dyn (dt)
use ice_blocks, only: nx_block, ny_block
use ice_domain, only: nblocks, halo_dynbundle
use ice_domain_size, only: max_blocks
use ice_flux, only: rdg_conv, rdg_shear, iceumask, &
use ice_flux, only: rdg_conv, rdg_shear, iceumask, iceemask, icenmask, &
stressp_1, stressp_2, stressp_3, stressp_4, &
stressm_1, stressm_2, stressm_3, stressm_4, &
stress12_1, stress12_2, stress12_3, stress12_4, &
Expand Down Expand Up @@ -311,7 +311,10 @@ subroutine init_dyn (dt)

! ice extent mask on velocity points
iceumask(i,j,iblk) = .false.

if (grid_ice == 'CD' .or. grid_ice == 'C') then
iceemask(i,j,iblk) = .false.
icenmask(i,j,iblk) = .false.
end if
enddo ! i
enddo ! j
enddo ! iblk
Expand Down
91 changes: 59 additions & 32 deletions cicecore/cicedynB/general/ice_flux.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -39,7 +39,7 @@ module ice_flux

real (kind=dbl_kind), dimension (:,:,:), allocatable, public :: &

! in from atmos (if .not.calc_strair)
! in from atmos (if .not.calc_strair)
strax , & ! wind stress components (N/m^2), on grid_atm_dynu
stray , & ! on grid_atm_dynv

Expand All @@ -48,7 +48,7 @@ module ice_flux
vocn , & ! ocean current, y-direction (m/s), on grid_ocn_dynv
ss_tltx , & ! sea surface slope, x-direction (m/m), on grid_ocn_dynu
ss_tlty , & ! sea surface slope, y-direction, on grid_ocn_dynv
hwater , & ! water depth for seabed stress calc (landfast ice)
hwater , & ! water depth for seabed stress calc (landfast ice)

! out to atmosphere
strairxT, & ! stress on ice by air, x-direction at T points, computed in icepack
Expand Down Expand Up @@ -103,7 +103,7 @@ module ice_flux
dvirdgdt, & ! rate of ice volume ridged (m/s)
opening ! rate of opening due to divergence/shear (1/s)

real (kind=dbl_kind), &
real (kind=dbl_kind), &
dimension (:,:,:,:), allocatable, public :: &
! ridging diagnostics in categories
dardg1ndt, & ! rate of area loss by ridging ice (1/s)
Expand All @@ -114,7 +114,7 @@ module ice_flux
ardgn, & ! fractional area of ridged ice
vrdgn, & ! volume of ridged ice
araftn, & ! rafting ice area
vraftn, & ! rafting ice volume
vraftn, & ! rafting ice volume
aredistn, & ! redistribution function: fraction of new ridge area
vredistn ! redistribution function: fraction of new ridge volume

Expand Down Expand Up @@ -178,7 +178,7 @@ module ice_flux
! NOTE: when in CICE_IN_NEMO mode, these are gridbox mean fields,
! not per ice area. When in standalone mode, these are per ice area.

real (kind=dbl_kind), &
real (kind=dbl_kind), &
dimension (:,:,:,:), allocatable, public :: &
fsurfn_f , & ! net flux to top surface, excluding fcondtop
fcondtopn_f, & ! downward cond flux at top surface (W m-2)
Expand All @@ -201,7 +201,7 @@ module ice_flux
Tf , & ! freezing temperature (C)
qdp , & ! deep ocean heat flux (W/m^2), negative upward
hmix , & ! mixed layer depth (m)
daice_da ! data assimilation concentration increment rate
daice_da ! data assimilation concentration increment rate
! (concentration s-1)(only used in hadgem drivers)

! out to atmosphere (if calc_Tsfc)
Expand Down Expand Up @@ -247,8 +247,8 @@ module ice_flux
dimension(:,:,:,:), allocatable, public :: &
albcnt ! counter for zenith angle

! out to ocean
! (Note CICE_IN_NEMO does not use these for coupling.
! out to ocean
! (Note CICE_IN_NEMO does not use these for coupling.
! It uses fresh_ai,fsalt_ai,fhocn_ai and fswthru_ai)
real (kind=dbl_kind), dimension (:,:,:), allocatable, public :: &
fpond , & ! fresh water flux to ponds (kg/m^2/s)
Expand Down Expand Up @@ -280,7 +280,7 @@ module ice_flux
snoicen ! snow-ice formation in category n (m)

real (kind=dbl_kind), dimension (:,:,:,:), allocatable, public :: &
keffn_top ! effective thermal conductivity of the top ice layer
keffn_top ! effective thermal conductivity of the top ice layer
! on categories (W/m^2/K)

! quantities passed from ocean mixed layer to atmosphere
Expand Down Expand Up @@ -324,7 +324,7 @@ module ice_flux
frz_onset, &! day of year that freezing begins (congel or frazil)
frazil_diag ! frazil ice growth diagnostic (m/step-->cm/day)

real (kind=dbl_kind), &
real (kind=dbl_kind), &
dimension (:,:,:,:), allocatable, public :: &
fsurfn, & ! category fsurf
fcondtopn,& ! category fcondtop
Expand All @@ -339,7 +339,7 @@ module ice_flux
! As above but these remain grid box mean values i.e. they are not
! divided by aice at end of ice_dynamics. These are used in
! CICE_IN_NEMO for coupling and also for generating
! ice diagnostics and history files as these are more accurate.
! ice diagnostics and history files as these are more accurate.
! (The others suffer from problem of incorrect values at grid boxes
! that change from an ice free state to an icy state.)

Expand Down Expand Up @@ -369,12 +369,12 @@ module ice_flux
rside , & ! fraction of ice that melts laterally
fside , & ! lateral heat flux (W/m^2)
fsw , & ! incoming shortwave radiation (W/m^2)
coszen , & ! cosine solar zenith angle, < 0 for sun below horizon
coszen , & ! cosine solar zenith angle, < 0 for sun below horizon
rdg_conv, & ! convergence term for ridging (1/s)
rdg_shear ! shear term for ridging (1/s)

real (kind=dbl_kind), dimension(:,:,:,:), allocatable, public :: &
salinz ,& ! initial salinity profile (ppt)
salinz ,& ! initial salinity profile (ppt)
Tmltz ! initial melting temperature (^oC)

!=======================================================================
Expand All @@ -383,7 +383,7 @@ module ice_flux

!=======================================================================
!
! Allocate space for all variables
! Allocate space for all variables
!
subroutine alloc_flux

Expand All @@ -393,12 +393,12 @@ subroutine alloc_flux

allocate( &
strax (nx_block,ny_block,max_blocks), & ! wind stress components (N/m^2)
stray (nx_block,ny_block,max_blocks), & !
stray (nx_block,ny_block,max_blocks), & !
uocn (nx_block,ny_block,max_blocks), & ! ocean current, x-direction (m/s)
vocn (nx_block,ny_block,max_blocks), & ! ocean current, y-direction (m/s)
ss_tltx (nx_block,ny_block,max_blocks), & ! sea surface slope, x-direction (m/m)
ss_tlty (nx_block,ny_block,max_blocks), & ! sea surface slope, y-direction
hwater (nx_block,ny_block,max_blocks), & ! water depth for seabed stress calc (landfast ice)
hwater (nx_block,ny_block,max_blocks), & ! water depth for seabed stress calc (landfast ice)
strairxT (nx_block,ny_block,max_blocks), & ! stress on ice by air, x-direction
strairyT (nx_block,ny_block,max_blocks), & ! stress on ice by air, y-direction
strocnxT (nx_block,ny_block,max_blocks), & ! ice-ocean stress, x-direction
Expand Down Expand Up @@ -544,7 +544,7 @@ subroutine alloc_flux
rside (nx_block,ny_block,max_blocks), & ! fraction of ice that melts laterally
fside (nx_block,ny_block,max_blocks), & ! lateral melt rate (W/m^2)
fsw (nx_block,ny_block,max_blocks), & ! incoming shortwave radiation (W/m^2)
coszen (nx_block,ny_block,max_blocks), & ! cosine solar zenith angle, < 0 for sun below horizon
coszen (nx_block,ny_block,max_blocks), & ! cosine solar zenith angle, < 0 for sun below horizon
rdg_conv (nx_block,ny_block,max_blocks), & ! convergence term for ridging (1/s)
rdg_shear (nx_block,ny_block,max_blocks), & ! shear term for ridging (1/s)
dardg1ndt (nx_block,ny_block,ncat,max_blocks), & ! rate of area loss by ridging ice (1/s)
Expand All @@ -555,7 +555,7 @@ subroutine alloc_flux
ardgn (nx_block,ny_block,ncat,max_blocks), & ! fractional area of ridged ice
vrdgn (nx_block,ny_block,ncat,max_blocks), & ! volume of ridged ice
araftn (nx_block,ny_block,ncat,max_blocks), & ! rafting ice area
vraftn (nx_block,ny_block,ncat,max_blocks), & ! rafting ice volume
vraftn (nx_block,ny_block,ncat,max_blocks), & ! rafting ice volume
aredistn (nx_block,ny_block,ncat,max_blocks), & ! redistribution function: fraction of new ridge area
vredistn (nx_block,ny_block,ncat,max_blocks), & ! redistribution function: fraction of new ridge volume
fsurfn_f (nx_block,ny_block,ncat,max_blocks), & ! net flux to top surface, excluding fcondtop
Expand All @@ -575,7 +575,7 @@ subroutine alloc_flux
flatn (nx_block,ny_block,ncat,max_blocks), & ! category latent heat flux
albcnt (nx_block,ny_block,max_blocks,max_nstrm), & ! counter for zenith angle
snwcnt (nx_block,ny_block,max_blocks,max_nstrm), & ! counter for snow
salinz (nx_block,ny_block,nilyr+1,max_blocks), & ! initial salinity profile (ppt)
salinz (nx_block,ny_block,nilyr+1,max_blocks), & ! initial salinity profile (ppt)
Tmltz (nx_block,ny_block,nilyr+1,max_blocks), & ! initial melting temperature (^oC)
stat=ierr)
if (ierr/=0) call abort_ice('(alloc_flux): Out of memory')
Expand Down Expand Up @@ -719,7 +719,7 @@ subroutine init_coupler_flux
fcondtopn_f(:,:,:,:) = c0 ! conductive heat flux (W/m^2)
flatn_f (:,:,:,:) = -1.0_dbl_kind ! latent heat flux (W/m^2)
fsensn_f (:,:,:,:) = c0 ! sensible heat flux (W/m^2)
endif !
endif !

fiso_atm (:,:,:,:) = c0 ! isotope deposition rate (kg/m2/s)
faero_atm (:,:,:,:) = c0 ! aerosol deposition rate (kg/m2/s)
Expand Down Expand Up @@ -762,7 +762,7 @@ subroutine init_coupler_flux
flat (:,:,:) = c0
fswabs (:,:,:) = c0
fswint_ai(:,:,:) = c0
flwout (:,:,:) = -stefan_boltzmann*Tffresh**4
flwout (:,:,:) = -stefan_boltzmann*Tffresh**4
! in case atm model diagnoses Tsfc from flwout
evap (:,:,:) = c0
evaps (:,:,:) = c0
Expand Down Expand Up @@ -816,7 +816,7 @@ subroutine init_coupler_flux

coszen (:,:,:) = c0 ! Cosine of the zenith angle
fsw (:,:,:) = c0 ! shortwave radiation (W/m^2)
scale_factor(:,:,:) = c1 ! shortwave scaling factor
scale_factor(:,:,:) = c1 ! shortwave scaling factor
wind (:,:,:) = sqrt(uatm(:,:,:)**2 &
+ vatm(:,:,:)**2) ! wind speed, (m/s)
Cdn_atm(:,:,:) = (vonkar/log(zref/iceruf)) &
Expand Down Expand Up @@ -986,8 +986,8 @@ subroutine init_history_therm
snowfrac (:,:,:) = c0
frazil_diag (:,:,:) = c0

! drag coefficients are computed prior to the atmo_boundary call,
! during the thermodynamics section
! drag coefficients are computed prior to the atmo_boundary call,
! during the thermodynamics section
Cdn_ocn(:,:,:) = dragio
Cdn_atm(:,:,:) = (vonkar/log(zref/iceruf)) &
* (vonkar/log(zref/iceruf)) ! atmo drag for RASM
Expand Down Expand Up @@ -1023,6 +1023,7 @@ end subroutine init_history_therm
subroutine init_history_dyn

use ice_state, only: aice, vice, trcr, strength
use ice_grid, only: grid_ice

logical (kind=log_kind) :: &
tr_iage
Expand Down Expand Up @@ -1074,6 +1075,32 @@ subroutine init_history_dyn
aredistn (:,:,:,:) = c0
vredistn (:,:,:,:) = c0

if (grid_ice == "CD" .or. grid_ice == "C") then
taubxE (:,:,:) = c0
taubyE (:,:,:) = c0
strocnxE (:,:,:) = c0
strocnyE (:,:,:) = c0
strairxE (:,:,:) = c0
strairyE (:,:,:) = c0
strtltxE (:,:,:) = c0
strtltyE (:,:,:) = c0
strintxE (:,:,:) = c0
strintyE (:,:,:) = c0
fmE (:,:,:) = c0
TbE (:,:,:) = c0
taubxN (:,:,:) = c0
taubyN (:,:,:) = c0
strocnxN (:,:,:) = c0
strocnyN (:,:,:) = c0
strairxN (:,:,:) = c0
strairyN (:,:,:) = c0
strtltxN (:,:,:) = c0
strtltyN (:,:,:) = c0
strintxN (:,:,:) = c0
strintyN (:,:,:) = c0
fmN (:,:,:) = c0
TbN (:,:,:) = c0
end if
end subroutine init_history_dyn

!=======================================================================
Expand Down Expand Up @@ -1166,8 +1193,8 @@ subroutine scale_fluxes (nx_block, ny_block, &

! zsalinity fluxes
real (kind=dbl_kind), dimension(nx_block,ny_block), intent(inout) :: &
fzsal , & ! salt flux to ocean with prognositic salinity (kg/m2/s)
fzsal_g ! Gravity drainage salt flux to ocean (kg/m2/s)
fzsal , & ! salt flux to ocean with prognositic salinity (kg/m2/s)
fzsal_g ! Gravity drainage salt flux to ocean (kg/m2/s)

! isotopes
real (kind=dbl_kind), dimension(nx_block,ny_block,icepack_max_iso), &
Expand Down Expand Up @@ -1221,8 +1248,8 @@ subroutine scale_fluxes (nx_block, ny_block, &
alidr (i,j) = alidr (i,j) * ar
alvdf (i,j) = alvdf (i,j) * ar
alidf (i,j) = alidf (i,j) * ar
fzsal (i,j) = fzsal (i,j) * ar
fzsal_g (i,j) = fzsal_g (i,j) * ar
fzsal (i,j) = fzsal (i,j) * ar
fzsal_g (i,j) = fzsal_g (i,j) * ar
flux_bio (i,j,:) = flux_bio (i,j,:) * ar
faero_ocn(i,j,:) = faero_ocn(i,j,:) * ar
if (present(Qref_iso )) Qref_iso (i,j,:) = Qref_iso (i,j,:) * ar
Expand Down Expand Up @@ -1251,10 +1278,10 @@ subroutine scale_fluxes (nx_block, ny_block, &
fswthru_idf (i,j) = c0
alvdr (i,j) = c0 ! zero out albedo where ice is absent
alidr (i,j) = c0
alvdf (i,j) = c0
alvdf (i,j) = c0
alidf (i,j) = c0
fzsal (i,j) = c0
fzsal_g (i,j) = c0
fzsal (i,j) = c0
fzsal_g (i,j) = c0
flux_bio (i,j,:) = c0
faero_ocn(i,j,:) = c0
if (present(Qref_iso )) Qref_iso (i,j,:) = c0
Expand All @@ -1265,7 +1292,7 @@ subroutine scale_fluxes (nx_block, ny_block, &
enddo ! j

! Scale fluxes for history output
if (present(fsurf) .and. present(fcondtop) ) then
if (present(fsurf) .and. present(fcondtop) ) then

do j = 1, ny_block
do i = 1, nx_block
Expand Down
2 changes: 1 addition & 1 deletion cicecore/drivers/nuopc/cmeps/CICE_RunMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -129,7 +129,7 @@ subroutine ice_step
use ice_restoring, only: restore_ice, ice_HaloRestore
use ice_step_mod, only: prep_radiation, step_therm1, step_therm2, &
update_state, step_dyn_horiz, step_dyn_ridge, step_radiation, &
biogeochemistry, save_init, step_dyn_wave, step_snow
biogeochemistry, step_prep, step_dyn_wave, step_snow
use ice_timers, only: ice_timer_start, ice_timer_stop, &
timer_diags, timer_column, timer_thermo, timer_bound, &
timer_hist, timer_readwrite
Expand Down
1 change: 0 additions & 1 deletion cicecore/drivers/nuopc/cmeps/ice_comp_nuopc.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -549,7 +549,6 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc)

scol_valid = (scol_mask == 1)
if (.not. scol_valid) then
write(6,*)'DEBUG: i am here'
! Advertise fields
call ice_advertise_fields(gcomp, importState, exportState, flds_scalar_name, rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
Expand Down
21 changes: 6 additions & 15 deletions cicecore/drivers/nuopc/cmeps/ice_import_export.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -187,9 +187,9 @@ subroutine ice_advertise_fields(gcomp, importState, exportState, flds_scalar_nam
! in the cmeps esmFldsExchange_xxx_mod.F90 that is model specific
! from atm - black carbon deposition fluxes (3)
call fldlist_add(fldsToIce_num, fldsToIce, 'Faxa_bcph', ungridded_lbound=1, ungridded_ubound=3)
! from atm - wet dust deposition frluxes (4 sizes)
! from atm - wet dust deposition fluxes (4 sizes)
call fldlist_add(fldsToIce_num, fldsToIce, 'Faxa_dstwet', ungridded_lbound=1, ungridded_ubound=4)
! from - atm dry dust deposition frluxes (4 sizes)
! from atm - dry dust deposition fluxes (4 sizes)
call fldlist_add(fldsToIce_num, fldsToIce, 'Faxa_dstdry', ungridded_lbound=1, ungridded_ubound=4)

do n = 1,fldsToIce_num
Expand Down Expand Up @@ -800,19 +800,10 @@ subroutine ice_import( importState, rc )

if (.not.prescribed_ice) then
call t_startf ('cice_imp_t2u')
call ice_HaloUpdate(uocn, halo_info, field_loc_center, field_type_scalar)
call ice_HaloUpdate(vocn, halo_info, field_loc_center, field_type_scalar)
call ice_HaloUpdate(ss_tltx, halo_info, field_loc_center, field_type_scalar)
call ice_HaloUpdate(ss_tlty, halo_info, field_loc_center, field_type_scalar)
! tcraig, moved to dynamics for consistency
!work = uocn
!call grid_average_X2Y('F',work,'T',uocn,'U')
!work = vocn
!call grid_average_X2Y('F',work,'T',vocn,'U')
!work = ss_tltx
!call grid_average_X2Y('F',work,'T',ss_tltx,'U')
!work = ss_tlty
!call grid_average_X2Y('F',work,'T',ss_tlty,'U')
call ice_HaloUpdate(uocn, halo_info, field_loc_center, field_type_vector)
call ice_HaloUpdate(vocn, halo_info, field_loc_center, field_type_vector)
call ice_HaloUpdate(ss_tltx, halo_info, field_loc_center, field_type_vector)
call ice_HaloUpdate(ss_tlty, halo_info, field_loc_center, field_type_vector)
call t_stopf ('cice_imp_t2u')
end if

Expand Down

0 comments on commit 471c010

Please sign in to comment.