...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
| Table of Contents | ||||||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|
|
Overview
Ice sheet coupling with the rest of the model depends on which ice sheet is being modeled, Antarctica or Greenland.
...
ocean coupling: (glc_c2_ocn): Sg_icemask, Sg_blit, Sg_blis, Sg_lithop, Sg_icemask_grounded
compset (how tested) | glc_c2_lnd | lnd_c2_glc | ocn_c2_glcshelf | glcshelf_c2_ocn | glcshelf_c2_ice | ocn_c2_glctf | glc_c2_ocn | glc_c2_ice |
|---|---|---|---|---|---|---|---|---|
MPAS_LISIO_JRA1p (not tested without mod below) | T | |||||||
MPAS_LISIO_JRA1p with mod mpaso/ocn_glcshelf ais8to30, MALI%SIA | T | T | T | T | ||||
IGELM_MLI gis20 MALI | T | T | ||||||
BGWCYCL1850 gis20 MALI | T | T | T |
There is no ocn_c2_glc, ice_c2_glcshelf., ice_c2_glc
Attribute strings
Fields passed to/from the coupler related to the ice sheet.
Several fields are calculated per elevation class and so have two-digit numbers appended to them.
land to glc adds the following coupling fields
l2x states: Sl_tsrf00 through Sl_tsrf10 and Sl_topo00 through Sl_topo10
l2x states to glc: Sl_tsrf00 to Sl_tsrf10 AND Sl_topo00 to Sl_topo10
Sl_tsrf: Surface temperature of glacier (one for each elevation class)
Sl_topo: surface height of glacier (one for each elevation class)
l2x fluxes: Flgl_qice00 through Flgl_qice10
l2x_fluxes_to_glc: Flgl_qice00 to Flgl_qice10
Flgl_qice: New glacier ice flux
glc to lnd adds the following coupling fields:
...
Assuming glc_present and glc_prognostic are TRUE so not including.
init
prep_lnd_init (glc_c2_lnd)
if(glc_c2_lnd) init mapper_Sg2l, mapper_Fg2l; derive field lists glc2lnd_non_ec_fields and glc2lnd_ec_extra_fields
prep_ocn_init (glc_c2_ocn, glcshelf_c2_ocn)
if (glc_c2_ocn) init mapper_Rg2o_liq and mapper_Rg2o_ice
if(glcshelf_c2_ocn) init mapper_Sg2o and mapper_Fg2o
prep_ice_init (glc_c2_ice, glcshelf_c2_ice)
if (glc_c2_ice) init mapper_Rg2i
if(glcshelf_c2_ice) init mapper_Sg2i,mapper_Fg2i
prep_glc_init( lnd_c2_glc, ocn_c2_glctf, ocn_c2_glcshelf)
if ( lnd_c2_glc .or. do_hist_l2x1yrg) initialize l2gacc_lx with seq_flds_l2x_fields_to_glc
if ( lnd_c2_glc) initialize l2x_gx with seq_flds_x2g_fields
if (lnd_c2_glc) init mapper_Sl2g, mapper_Fl2g, mapper_Fg2l; form g2x_lx_fields
if (ocn_c2_glctf OR ocn_c2_glcshelf) init o2x_gx with seq_flds_o2x_fields; init x2gacc_gx
also init mapper_So2g, mapper_Fo2g and arrays for compute_melt_fluxes
if (ocn_c2_glctf) init mapper_So2g_tf
if(ocn_c2_glcshelf) init mapper_So2g_shelf
if (glc_c2_ocn) do mapping Rg2o_liq, Rg2o_ice; g2x_gx to g2x_ox; seq_flds_g2o_liq/ice_fluxes,
if (glcshelf_c2_ocn) do mapping Sg2o, Fg2o; g2x_gx to g2x_ox MIGHT OVERWRITE ABOVE
if(glc_c2_ice) do mapping Rg2i; g2x_gx to g2x_ix ‘Fixx_rofi’ only
if(glcshelf_c2_ice) do mapping_Sg2i; g2x_gx, g2x_ix, Sg_icemask_coupled_fluxes
if (glc_c2_lnd) prep_lnd_calc_g2x_lx do mapping with Fg2l, from g2x_gx to g2x_lx, only glc2lnd_non_ec_fields (Sg_icemask, Sg_icemask_coupled_fluxes)
loop over elevation classes: call map_glc2lnd_ec do special mapping of elevation class fields with from glc to lnd elevation classes using mapper_Fg2l, glc2lnd_ec_extra_fields
...
call cime_run_glc_setup_send
if (lnd_c2_glc OR ocn_c2_glctf OR ocn_c2_glcshelf)
if (glcrun_avg_alarm) accumulate for lnd and ocn
if (lnd_c2_glc) call prep_glc_calc_l2x_gx
...
if (smb_renormalize) prep_glc_renormalize_smb;(calls map_glc2lnd_ec)
call prep_glc_map_one_state_field_lnd2glc: special bilinear mapping one state at a time: also call map_lnd2glc
if (ocn_c2_glctf) call prep_glc_mrg_ocn
else prep_glc_zero_flds
...
call cime_run_ocnglc_coupling()
if (ocn_c2_glctf OR (glcshelf_to_ocn AND ocn_to_glcshelf))
call prep_glc_calc_o2x_gx:
if (ocn_c2_gtctf) call mapper_So2g_tf; o2x_ox to o2x_gx seq_flds_x2g_tf_states_from_ocn
if (ocn_c2_glcshelf) call mapper_So2g_shelf; o2x_ox to o2x_gx seq_flds_x2g_shelf_states_from_ocn
if (ocn_c2_glcshelf AND glcshelf_c2_ocn)
call prep_glc_calculate_subshelf_boundary_fluxes
call mapper_So2g_shelf; o2x_ox, x2g_gx; seq_flds_x2g_shelf_states_from_ocn
compute_melt_fluxes
call prep_ocn_shelf_calc_g2x_ox: map_Sg2o, map_Fg2o; g2x_gx, g2x_ox
call prep_glc_accum_ocn (accumulate x2g_g fields)
if ( glcshelf_c2_ice) call mapper map_Sg2i; g2x_gx, g2x_ix; Sg_icemask_coupled_fluxes
Special functions in land ice coupling
prep_glc_accum_avg
called by prepcime_lndrun_calcglc_g2xsetup_lxsend
| Code Block |
|---|
call map_glc2lnd_ec( & g2x_g = g2x_gx, & frac_field = glc_frac_field, &call prep_glc_accum_avg(timer='CPL:glcprep_avg', & topo_field = glc_topo_field, & icemask_field = glc_icemask_field, & extra_fields = glc2lnd_ec_extra_fields, & mapper = mapper_Fg2l, & g2x_l = g2x_lx(egi) |
glc_frac_field = ‘Sg_ice_covered’
glc_topo_field = ‘Sg_topo’
glc_icemask_field = 'Sg_icemask'
glc2lnd_ec_extra_fields = Flgg_hflx
At each grid point on land ice mesh, Flgg_hflx is defined for several elevation classes.
Given topographic height at a point, find which elevation class its in.
...
lnd2glc_averaged_now=lnd2glc_averaged_now) |
accumulation Avs: l2gacc_lx and x2gacc_gx are local save private Avs in prep_glc_mod
if (l2gacc_lx_cnt > 1) call mct_avect_avg(l2gacc_lx(eli), l2gacc_lx_cnt)
if (x2gacc_gx_cnt > 1) then
call mct_avect_avg(x2gacc_gx(egi), x2gacc_gx_cnt)
call mct_avect_copy(x2gacc_gx(egi), x2g_gx)
prep_glc_calc_l2x_gx
called by cime_run_glc_setup_send
call prep_glc_calc_l2x_gx(fractions_lx, timer='CPL:glcprep_lnd2glc')
loop over 1, number of fields in(seq_flds_x2g_fluxes_from_lnd) which is 1: Flgl_qice
call prep_glc_map_qice_conservative_lnd2glc
end loop
loop over 1, number of fields in seq_flds_x2g_states_from_lnd which is 1: Sl_tsrf
call prep_glc_map_one_state_field_lnd2glc
end loop
prep_glc_map_qice_conservate_lnd2glc
called by prep_glc_calc_l2x_gx which is called by cime_run_glc_setup_send
| Code Block |
|---|
! Use a bilinear (Sl2g) mapper, as for states.
! The Fg2l mapper is needed to map some glc fields to the land grid
! for purposes of conservation.
call prep_glc_map_qice_conservative_lnd2glc(egi=egi, eli=eli, &
fractions_lx = fractions_lx(efi), &
mapper_Sl2g = mapper_Sl2g, &
mapper_Fg2l = mapper_Fg2l) |
using the accumulated l2gacc_lx, map only Flgl_qice (the SMB) using map_lnd2glc with mapper_Sl2g and l2x_gx
export Flgl_qice from l2x_gx to an array qice_g
multiply qice_g by area_g/aream_g
if (smb_renormalize) call prep_glc_renormalize_smb using g2x_gx with mapper_Fg2l. Also modifies qice_g
import qice_g back in to l2x_gx in the Flgl_qice attribute.
prep_glc_map_one_state_field_lnd2glc
called only by prep_glc_calc_l2x_gx which is called by cime_run_glc_setup_send
called only for Sl_tsrf
calls map_lnd2glc for Sl_tsrf mapping local l2gacc_lx to l2x_gx. Also uses fields from g2x_gx
map_glc2lnd_ec
called by prep_lnd_calc_g2x_lx and prep_glc_renormalize_smb
| Code Block |
|---|
call map_glc2lnd_ec( &
g2x_g = g2x_gx, &
frac_field = glc_frac_field, &
topo_field = glc_topo_field, &
icemask_field = glc_icemask_field, &
extra_fields = glc2lnd_ec_extra_fields, &
mapper = mapper_Fg2l, &
g2x_l = g2x_lx(egi) |
glc_frac_field = ‘Sg_ice_covered’
glc_topo_field = ‘Sg_topo’
glc_icemask_field = 'Sg_icemask'
glc2lnd_ec_extra_fields = Flgg_hflx
This will map Sg_ice_covered, Sg_topo and Flgg_hflx from glc mesh to land model mesh with elevation classes: Sg_ice_covered00, Sg_ice_covered01, etc. Sg_topo00, Sg_topo01, etc., Flgg_hflx00, Flgg_hflx01, etc. There is a loop over elevations classes.
The call from prep_glc_renormalize_smb does not include extra_fields.
map_lnd2glc
called by both prep_glc_map_qice_conservative_lnd2glc and prep_glc_map_one_state_field_lnd2glc
maps one field from land with elevation classes to “glacier”
Note: “fieldname” for the “conservative” is Flql_qice which is calculated on multiple elevation classes in the land model.
export Sg_ice_covered and Sg_topo into arrays.
call get_glc_elevation_classes(glc_ice_covered, glc_topo, glc_elevclass) which returns the integer array glc_elevclass that has the elevation class index that glc_topo is contained in. glc_ice_covered is either 0 (bare land) or 1 (fully ice).
call map_bare_land (l2x_l, landfrac_l, “fieldname”, mapper, array) creates an glacier av with fieldname_bare_land = fieldname // elevclass_as_string as the only attribute and maps the one field (which is on elev class 0) to that glacier av from l2x_l. Then pulls out the field into an array and passes it back.
That field is copied to data_g
call map_ice_covered(l2x_l, landfrac_l, “fieldname”, glc_topo, mapper, array) see below. returns an array callled data_g_ice_covered
Form the final data_g field with
| Code Block |
|---|
where (glc_elevclass /= 0)
data_g = data_g_ice_covered
end where |
map_ice_covered
called only by map_lnd2glc
create a list of of land model attributes with elevation classes except for class 0. Result will be something like 'Flgl_qice01:Flgl_qice02: ... :Flgl_qice10:Sl_topo01:Sl_topo02: ... :Sltopo10'
Create a glacer-sized aV with above attributes.
map l2x_l to that aV with normalization and weights lfrac
loop over elevation classes and export the Flgl_qice and Sl_topo to a 2D array for each elevation class and point.
Create the field on the glacier mesh by looking at the elevation classes.
do i= 1, size_of_glacier_mesh
if topo(i) lt mapped_topo(i,1) then
the actual topography is lower then the height of the lowest elevation class so assign the field value of the lowest elevation class: data(i) = mapped_data(i, 1)
if topo(i) gt mapped_topo(i,nEC) then
the actual topography is higher then highest mapped elevation class so assign the field value of the highest elevation class: data(i) = mapped_data(i,nEC)
else need to do an interpolation with the 2 elevation classes nearest the actual topo.
do j = 2, nEC
if (topo(i) < mapped_topo(i,j) then
elev_lower = mapped_topo(i,j-1)
elev_upper = mapped_topo(i,j)
ediff = elev_upper - elev_lower
data(i) = mapped_data(i,j-1) * (elev-upper - topo(i))/ediff + mapped_data(i,j) * (topo(i)- elev_lower)/ediff
enddo
endif
enddo
return the data array.
prep_glc_renormalize_smb
called only by prep_glc_map_qice_conservate_lnd2glc
| Code Block |
|---|
! Renormalizes surface mass balance (smb, here named qice_g) so that the global
! integral on the glc grid is equal to the global integral on the land grid. |
qice_g is passed in as an array.
map Sg_icemask to land grid with normalization. (BUG note from Bill Sacks: I think we actually want norm = .false. here, but this requires some more thought. See https://github.com/ESMCI/cime/issues/1516)
Then pull out mapped field in to an array.
map Sg_ice_covered and Sg_topo to land grid and elevation classes with map_glc2lnd_ec. No extra fields as in other call to map_glc2lnd_ec.
loop over elevation classes
export the mapped Sg_ice_coverd and local time averaged Flgl_qice to 2D arrays sized by (land grid, elevation classes)
Find the local sum of negative and positive Flgl_qice with a double loop over land points and elevation classes.
loop over land points
effective_area = min(lfrac(n),Sg_icemask_l(n)) * aream_l(n) where lfrac is from fractions_lx
loop over elevation classes
if (qice >= 0) local_accum_on_land_grid = local_accum_on_land_grid + effective_area * frac_l(n,ec) * qice_l(n,ec) frac_l is the mapped Sg_ice_covered.
else local_ablat_on_land_grid = local_ablat_on_land_grid + effective_area * frac_l(n,ec) * qice_l(n,ec
end double loop
Use shr_mpi_sum to find the global sum global_accum_on_land_grid and global_ablat_on_land_grid from the local sums computed above NOT BFB with PROC COUNT. Then bcast global sum back to all procs
Loop over glacier points and find the local weighted sum of qice_g (the input array)
if(qice >=0) local_accum_on_glc_grid = local_accum_on_glc_grid + Sg_icemask_g(n) * aream_g(n) * qice_g(n)
else local_ablat_on_glc_grid = local_ablat_on_glc_grid + Sg_icemask_g(n) * aream_g(n) * qice_g(n)
end loop
use shr_mpi_sum again to find global sums global_accum_on_glc_grid, global_ablat_on_glc_grid NOT BFB with PROC COUNT. Then bcast global sum back to all procs.
Find the renormalization factor
| Code Block |
|---|
if (global_ablat_on_glc_grid < 0.0_r8) then ! negative by definition
ablat_renorm_factor = global_ablat_on_land_grid / global_ablat_on_glc_grid
if (global_accum_on_glc_grid > 0.0_r8) then
accum_renorm_factor = global_accum_on_land_grid / global_accum_on_glc_grid |
Apply renormalization
| Code Block |
|---|
do n = 1, lsize_g
if (qice_g(n) >= 0.0_r8) then
qice_g(n) = qice_g(n) * accum_renorm_factor
else
qice_g(n) = qice_g(n) * ablat_renorm_factor
endif
enddo |
return.
prep_glc_calculate_subshelf_boundary_fluxes
called by cime_run_ocnglc_coupling
On the ice sheet grid, calculate shelf boundary fluxes
Pull out input arrays from attribute vectors all on glacier grid
So_blt, So_bls, So_htv, So_stv. Sg_tbot, Sg_dztbot, Sg_lithop, Sg_icemask_floating
call compute_melt_fluxes
take output from above and place in attributes: Sg_blis, Sg_blit, Fogx_qiceho, Fogx_gicelo, Fogx_qicehi, Fogx_qiceli
compute_melt_fluxes
a physics calculation that works only with Fortran arrays (no Avs). called only by prep_glc_calculate_subshelf_boundary_fluxes
(Text in routine from Xylar) This routine computes melt fluxes (melt rate, temperature fluxes into the ice and the ocean, and salt flux) as well as the interface temperature and salinity. This routine expects an ice temperature in the bottom layer of ice and ocean temperature and salinity in the top ocean layer as well as the pressure at the ice/ocean interface.
The ocean heat and salt transfer velocities are determined based on observations of turbulent mixing rates in the under-ice boundary layer. They should be the product of the friction velocity and a (possibly spatially variable) non-dimensional transfer coefficient.
The iceTemperatureDistance (Sg_dztbot) is the distance between the location where the iceTemperature is supplied and the ice-ocean interface, used to compute a temperature gradient. The ice thermal conductivity, SHR_CONST_KAPPA_LAND_ICE, is zero for the freezing solution from Holland and Jenkins (1999) in which the ice is purely insulating.