Skip to content

Add GEOS coupling capability including heat and mass flux options in Icepack #1012

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
apcraig merged 13 commits into from
Jun 20, 2025
Merged

Add GEOS coupling capability including heat and mass flux options in Icepack #1012

Show file tree
Hide file tree
Changes from 12 commits
Commits
Show all changes
13 commits
Select commit Hold shift + click to select a range
21c9bcd
Add GEOS heatflux and massflux updates. This allows coupling to the …
apcraig Mar 16, 2025
c214f2b
Remove trailing blanks
apcraig Mar 16, 2025
6aacde6
updates/fixes to make it wiork with new GEOS PR
zhaobin74 Mar 26, 2025
25c9dd9
rearrange zero initialization; fixed a bug in restert filename_spec
zhaobin74 Mar 26, 2025
0466ca6
fixed indentation
zhaobin74 Mar 26, 2025
caa2cf6
allow A-grid currents and sea surface height tilt terms
zhaobin74 Apr 23, 2025
9341acf
Rebased to #95e581eb1d86cc06 on main and update implementation to
apcraig May 14, 2025
bc1b6e3
add missing variable import and definition; fixes due to icepack inte…
zhaobin74 May 15, 2025
8513eed
Merge pull request #122 from zhaobin74/geos
apcraig May 21, 2025
4bb7d15
Update CICE consistent with latest changes to Icepack for
apcraig May 28, 2025
d068a27
Merge branch 'main' into geos
apcraig May 28, 2025
9e1fd4a
Update Icepack to #ac57a1c48af244 includes
apcraig Jun 20, 2025
cd726e8
Update Icepack to #f61be9a1bd3578f89 with minor GEOS modifications
apcraig Jun 20, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
6 changes: 3 additions & 3 deletions cicecore/cicedyn/dynamics/ice_dyn_eap.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2090,8 +2090,7 @@ subroutine read_restart_eap()
use ice_boundary, only: ice_HaloUpdate_stress
use ice_constants, only: &
field_loc_center, field_type_scalar
use ice_domain, only: nblocks, halo_info
use ice_grid, only: grid_type
use ice_domain, only: nblocks, halo_info, ns_boundary_type
use ice_restart, only: read_restart_field

! local variables
Expand Down Expand Up @@ -2131,7 +2130,8 @@ subroutine read_restart_eap()
call read_restart_field(nu_restart_eap,0,a12_4,'ruf8', &
'a12_4',1,diag,field_loc_center,field_type_scalar) ! a12_4

if (trim(grid_type) == 'tripole') then
if (trim(ns_boundary_type) == 'tripole' .or. &
trim(ns_boundary_type) == 'tripoleT') then

call ice_HaloUpdate_stress(a11_1, a11_3, halo_info, &
field_loc_center, field_type_scalar)
Expand Down
9 changes: 5 additions & 4 deletions cicecore/cicedyn/dynamics/ice_dyn_evp.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -259,7 +259,8 @@ subroutine evp (dt)
use ice_boundary, only: ice_halo, ice_HaloMask, ice_HaloUpdate, &
ice_HaloDestroy, ice_HaloUpdate_stress
use ice_blocks, only: block, get_block, nx_block, ny_block, nghost
use ice_domain, only: nblocks, blocks_ice, halo_info, maskhalo_dyn
use ice_domain, only: nblocks, blocks_ice, halo_info, maskhalo_dyn, &
ns_boundary_type
use ice_domain_size, only: max_blocks, ncat
use ice_flux, only: rdg_conv, rdg_shear, strairxT, strairyT, &
strairxU, strairyU, uocn, vocn, ss_tltx, ss_tlty, fmU, &
Expand All @@ -280,8 +281,7 @@ subroutine evp (dt)
use ice_grid, only: tmask, umask, umaskCD, nmask, emask, uvm, epm, npm, &
dxE, dxN, dxT, dxU, dyE, dyN, dyT, dyU, &
tarear, uarear, earear, narear, grid_average_X2Y, uarea, &
grid_type, grid_ice, &
grid_atm_dynu, grid_atm_dynv, grid_ocn_dynu, grid_ocn_dynv
grid_ice, grid_atm_dynu, grid_atm_dynv, grid_ocn_dynu, grid_ocn_dynv
use ice_state, only: aice, aiU, vice, vsno, uvel, vvel, uvelN, vvelN, &
uvelE, vvelE, divu, shear, vort, &
aice_init, aice0, aicen, vicen, strength
Expand Down Expand Up @@ -1313,7 +1313,8 @@ subroutine evp (dt)
endif

! Force symmetry across the tripole seam
if (trim(grid_type) == 'tripole') then
if (trim(ns_boundary_type) == 'tripole' .or. &
trim(ns_boundary_type) == 'tripoleT') then
! TODO: C/CD-grid
if (maskhalo_dyn) then
!-------------------------------------------------------
Expand Down
10 changes: 9 additions & 1 deletion cicecore/cicedyn/dynamics/ice_dyn_shared.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -240,7 +240,8 @@ subroutine init_dyn_shared (dt)

use ice_blocks, only: block, get_block
use ice_boundary, only: ice_halo, ice_haloUpdate
use ice_domain, only: nblocks, halo_dynbundle, blocks_ice, halo_info
use ice_domain, only: nblocks, halo_dynbundle, blocks_ice, halo_info, &
ns_boundary_type
use ice_domain_size, only: max_blocks
use ice_flux, only: &
stressp_1, stressp_2, stressp_3, stressp_4, &
Expand Down Expand Up @@ -268,6 +269,13 @@ subroutine init_dyn_shared (dt)

character(len=*), parameter :: subname = '(init_dyn_shared)'

! checks
if (kdyn == 1 .and. evp_algorithm == 'shared_mem_1d' .and. &
(ns_boundary_type == 'tripole' .or. ns_boundary_type == 'tripoleT')) then
call abort_ice(subname//' ERROR: evp_alg shared mem 1d not supported with tripole', &
file=__FILE__, line=__LINE__)
endif

call set_evp_parameters (dt)
! allocate dyn shared (init_uvel,init_vvel)
call alloc_dyn_shared
Expand Down
8 changes: 5 additions & 3 deletions cicecore/cicedyn/dynamics/ice_dyn_vp.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -155,7 +155,8 @@ subroutine implicit_solver (dt)
use ice_boundary, only: ice_HaloMask, ice_HaloUpdate, &
ice_HaloDestroy, ice_HaloUpdate_stress
use ice_blocks, only: block, get_block, nx_block, ny_block
use ice_domain, only: blocks_ice, halo_info, maskhalo_dyn
use ice_domain, only: blocks_ice, halo_info, maskhalo_dyn, &
ns_boundary_type
use ice_domain_size, only: max_blocks, ncat
use ice_dyn_shared, only: deformations, iceTmask, iceUmask, &
cxp, cyp, cxm, cym
Expand All @@ -168,7 +169,7 @@ subroutine implicit_solver (dt)
stressm_1, stressm_2, stressm_3, stressm_4, &
stress12_1, stress12_2, stress12_3, stress12_4
use ice_grid, only: tmask, umask, dxT, dyT, dxU, dyU, &
tarear, grid_type, grid_average_X2Y, &
tarear, grid_average_X2Y, &
grid_atm_dynu, grid_atm_dynv, grid_ocn_dynu, grid_ocn_dynv
use ice_state, only: aice, aiU, vice, vsno, uvel, vvel, divu, shear, vort, &
aice_init, aice0, aicen, vicen, strength
Expand Down Expand Up @@ -546,7 +547,8 @@ subroutine implicit_solver (dt)
endif

! Force symmetry across the tripole seam
if (trim(grid_type) == 'tripole') then
if (trim(ns_boundary_type) == 'tripole' .or. &
trim(ns_boundary_type) == 'tripoleT') then
if (maskhalo_dyn) then
!-------------------------------------------------------
! set halomask to zero because ice_HaloMask always keeps
Expand Down
68 changes: 54 additions & 14 deletions cicecore/cicedyn/general/ice_flux.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -114,8 +114,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), &
dimension (:,:,:,:), allocatable, public :: &
real (kind=dbl_kind), dimension (:,:,:,:), allocatable, public :: &
! ridging diagnostics in categories
dardg1ndt, & ! rate of area loss by ridging ice (1/s)
dardg2ndt, & ! rate of area gain by new ridges (1/s)
Expand Down Expand Up @@ -177,13 +176,26 @@ 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), &
dimension (:,:,:,:), allocatable, public :: &
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)
fsensn_f , & ! sensible heat flux (W m-2)
flatn_f ! latent heat flux (W m-2)

! in from atmosphere
! required for coupling in GEOS

real (kind=dbl_kind), dimension (:,:,:,:), allocatable, public :: &
evapn_f, & ! evaporation/sublimation (kg m-2 s-1)
dflatndTs_f, & ! derivative of latent flux w.r.t. Ts
dfsurfndTs_f ! derivative of surface flux w.r.t. Ts

real (kind=dbl_kind), dimension (:,:,:), allocatable, public :: &
swuvrdr , & ! vis uvr flux, direct (W m-2)
swuvrdf , & ! vis uvr flux, diffuse (W m-2)
swpardr , & ! vis par flux, direct (W m-2)
swpardf ! vis par flux, diffuse (W m-2)

! in from atmosphere

real (kind=dbl_kind), dimension (:,:,:), allocatable, public :: &
Expand Down Expand Up @@ -258,7 +270,11 @@ module ice_flux
fswthru_vdr , & ! vis dir shortwave penetrating to ocean (W/m^2)
fswthru_vdf , & ! vis dif shortwave penetrating to ocean (W/m^2)
fswthru_idr , & ! nir dir shortwave penetrating to ocean (W/m^2)
fswthru_idf ! nir dif shortwave penetrating to ocean (W/m^2)
fswthru_idf , & ! nir dif shortwave penetrating to ocean (W/m^2)
fswthru_uvrdr,& ! vis dir uvr SW penetrating to ocean (W/m^2)
fswthru_uvrdf,& ! vis dif uvr SW penetrating to ocean (W/m^2)
fswthru_pardr,& ! nir dir par SW penetrating to ocean (W/m^2)
fswthru_pardf ! nir dif par SW penetrating to ocean (W/m^2)

! internal

Expand Down Expand Up @@ -326,16 +342,14 @@ 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), &
dimension (:,:,:,:), allocatable, public :: &
real (kind=dbl_kind), dimension (:,:,:,:), allocatable, public :: &
fsurfn, & ! category fsurf
fcondtopn,& ! category fcondtop
fcondbotn,& ! category fcondbot
fsensn, & ! category sensible heat flux
flatn ! category latent heat flux

real (kind=dbl_kind), &
dimension (:,:,:,:), allocatable, public :: &
real (kind=dbl_kind), dimension (:,:,:,:), allocatable, public :: &
snwcnt ! counter for presence of snow

! As above but these remain grid box mean values i.e. they are not
Expand Down Expand Up @@ -453,6 +467,10 @@ subroutine alloc_flux
swvdf (nx_block,ny_block,max_blocks), & ! sw down, visible, diffuse (W/m^2)
swidr (nx_block,ny_block,max_blocks), & ! sw down, near IR, direct (W/m^2)
swidf (nx_block,ny_block,max_blocks), & ! sw down, near IR, diffuse (W/m^2)
swuvrdr (nx_block,ny_block,max_blocks), & ! vis uvr flux, direct (W m-2)
swuvrdf (nx_block,ny_block,max_blocks), & ! vis uvr flux, diffuse (W m-2)
swpardr (nx_block,ny_block,max_blocks), & ! vis par flux, direct (W m-2)
swpardf (nx_block,ny_block,max_blocks), & ! vis par flux, diffuse (W m-2)
flw (nx_block,ny_block,max_blocks), & ! incoming longwave radiation (W/m^2)
frain (nx_block,ny_block,max_blocks), & ! rainfall rate (kg/m^2 s)
fsnow (nx_block,ny_block,max_blocks), & ! snowfall rate (kg/m^2 s)
Expand Down Expand Up @@ -499,11 +517,15 @@ subroutine alloc_flux
fhocn (nx_block,ny_block,max_blocks), & ! net heat flux to ocean (W/m^2)
fsloss (nx_block,ny_block,max_blocks), & ! rate of snow loss to leads (kg/m^2/s)
fswthru (nx_block,ny_block,max_blocks), & ! shortwave penetrating to ocean (W/m^2)
fswthru_vdr (nx_block,ny_block,max_blocks), & ! vis dir shortwave penetrating to ocean (W/m^2)
fswthru_vdf (nx_block,ny_block,max_blocks), & ! vis dif shortwave penetrating to ocean (W/m^2)
fswthru_idr (nx_block,ny_block,max_blocks), & ! nir dir shortwave penetrating to ocean (W/m^2)
fswthru_idf (nx_block,ny_block,max_blocks), & ! nir dif shortwave penetrating to ocean (W/m^2)
scale_factor (nx_block,ny_block,max_blocks), & ! scaling factor for shortwave components
fswthru_vdr(nx_block,ny_block,max_blocks), & ! vis dir shortwave penetrating to ocean (W/m^2)
fswthru_vdf(nx_block,ny_block,max_blocks), & ! vis dif shortwave penetrating to ocean (W/m^2)
fswthru_idr(nx_block,ny_block,max_blocks), & ! nir dir shortwave penetrating to ocean (W/m^2)
fswthru_idf(nx_block,ny_block,max_blocks), & ! nir dif shortwave penetrating to ocean (W/m^2)
fswthru_uvrdr (nx_block,ny_block,max_blocks), & ! vis dir uvr SW penetrating to ocean (W/m^2)
fswthru_uvrdf (nx_block,ny_block,max_blocks), & ! vis dir uvr SW penetrating to ocean (W/m^2)
fswthru_pardr (nx_block,ny_block,max_blocks), & ! vis dir par SW penetrating to ocean (W/m^2)
fswthru_pardf (nx_block,ny_block,max_blocks), & ! vis dir par SW penetrating to ocean (W/m^2)
scale_factor (nx_block,ny_block,max_blocks), & ! scaling factor for shortwave components
strairx_ocn(nx_block,ny_block,max_blocks), & ! stress on ocean by air, x-direction
strairy_ocn(nx_block,ny_block,max_blocks), & ! stress on ocean by air, y-direction
fsens_ocn (nx_block,ny_block,max_blocks), & ! sensible heat flux (W/m^2)
Expand Down Expand Up @@ -566,6 +588,9 @@ subroutine alloc_flux
fcondtopn_f(nx_block,ny_block,ncat,max_blocks), & ! downward cond flux at top surface (W m-2)
fsensn_f (nx_block,ny_block,ncat,max_blocks), & ! sensible heat flux (W m-2)
flatn_f (nx_block,ny_block,ncat,max_blocks), & ! latent heat flux (W m-2)
evapn_f (nx_block,ny_block,ncat,max_blocks), & ! evaporative water flux (kg/m^2/s) by atmosphere model
dflatndTs_f(nx_block,ny_block,ncat,max_blocks), & ! derivative of flatn with respect to Ts
dfsurfndTs_f(nx_block,ny_block,ncat,max_blocks), & ! derivative of fsurfn with respect to Ts
Comment thread
eclare108213 marked this conversation as resolved.
Outdated
meltsn (nx_block,ny_block,ncat,max_blocks), & ! snow melt in category n (m)
melttn (nx_block,ny_block,ncat,max_blocks), & ! top melt in category n (m)
meltbn (nx_block,ny_block,ncat,max_blocks), & ! bottom melt in category n (m)
Expand All @@ -584,6 +609,11 @@ subroutine alloc_flux
stat=ierr)
if (ierr/=0) call abort_ice('(alloc_flux): Out of memory')

swuvrdr(:,:,:) = c0
swuvrdf(:,:,:) = c0
swpardr(:,:,:) = c0
swpardf(:,:,:) = c0

if (grid_ice == "CD" .or. grid_ice == "C") &
allocate( &
taubxN (nx_block,ny_block,max_blocks), & ! seabed stress (x) at N points (N/m^2)
Expand Down Expand Up @@ -792,6 +822,10 @@ subroutine init_coupler_flux
fswthru_vdf (:,:,:) = c0
fswthru_idr (:,:,:) = c0
fswthru_idf (:,:,:) = c0
fswthru_uvrdr (:,:,:) = c0
fswthru_uvrdf (:,:,:) = c0
fswthru_pardr (:,:,:) = c0
fswthru_pardf (:,:,:) = c0
fresh_da(:,:,:) = c0 ! data assimilation
fsalt_da(:,:,:) = c0
flux_bio (:,:,:,:) = c0 ! bgc
Expand Down Expand Up @@ -853,6 +887,8 @@ subroutine init_flux_atm
! strairxT(:,:,:) = 0.15_dbl_kind
! strairyT(:,:,:) = 0.15_dbl_kind

fsurf (:,:,:) = c0
fcondtop(:,:,:) = c0
fsens (:,:,:) = c0
flat (:,:,:) = c0
fswabs (:,:,:) = c0
Expand Down Expand Up @@ -898,6 +934,10 @@ subroutine init_flux_ocn
fswthru_vdf (:,:,:) = c0
fswthru_idr (:,:,:) = c0
fswthru_idf (:,:,:) = c0
fswthru_uvrdr(:,:,:) = c0
fswthru_uvrdf(:,:,:) = c0
fswthru_pardr(:,:,:) = c0
fswthru_pardf(:,:,:) = c0

faero_ocn (:,:,:,:) = c0
fiso_ocn (:,:,:,:) = c0
Expand Down
Loading