Tracking mass and energy out of PISM into a coupled GCM

[mankoff]

PISM conserves mass. But should it? Am I misunderstanding what this means? Isn't it important to be able to lose mass across various domain boundaries? I'm curious which variables do capture this, and across which cells.

If you prefer, I'm happy to break the long text below into multiple small Discussion posts.

Is it reasonable to think of 3 planes, top, bottom, and edges?

top is surface mass balance - I don't think I need to worry about surface mass balance because that is handled in our model, and we interface at the bottom of the firn. When we have surface melt in our model and no snow (e.g. blue ice in ablation zone) we track that and adjust the ice thickness by an appropriate amount for the next coupling step.

I'd like to separate bottom into grounded and floating ice.

For grounded ice, I see basal_mass_flux_grounded (units kg m-2 year-1). How is this different than tendency_of_subglacial_water_mass (units Gt year-1) other than units? I assume the latter is un-routed and the routed version(s) are tendency_of_subglacial_water_mass_at_grounded_margins ("subglacial water flux at grounded ice margins") or tendency_of_subglacial_water_mass_at_grounding_line ("subglacial water flux at grounding lines") and what is the difference between these latter two?

What terms go into this basal mass balance? I assume geothermal heat flux and frictional heating from sliding? What about viscous heat dissipation of surface-runoff routed to the bed?. Is surface runoff routed to the bed? I assume not, because that implies treatment of moulins and a level of detail I did not think existed in PISM.

For floating ice, do I subtract basal_mass_flux_grounded from tendency_of_ice_mass_due_to_basal_mass_flux' to get basal_mass_flux_floating`?

How is edge mass loss tracked? I see tendency_of_ice_amount_due_to_discharge which is "discharge flux (calving, frontal melt, forced retreat)", and I assume this is a combination of tendency_of_ice_mass_due_to_calving + tendency_of_ice_amount_due_to_frontal_melt + tendency_of_ice_mass_due_to_forced_retreat.

Constituent terms do not matter for mass conservation, but I think I need to work with constituent terms for energy conservation, because some of that is solid ice below the phase transition temperature (PTT) and some is liquid water I presume at the PTT. I need to track this properly into the appropriate GCM ocean cell(s).

[citibeth]

Ken, Conservation of mass does NOT mean the mass of the ice sheet always remains the same. It DOES mean that the mass continuity equation is satisfied over all time. I addressed this issue quite carefully, and presented an AGU poster on it as well. At the time I instrumented the code to track all mass going in / out of PISM, and was able to confirm that PISM does, in fact, conserve mass. I did the similar for energy and confirmed that PISM does NOT conserve energy. Are you aware of the place(s) in the code that perform this computation? I believe it's on the PISM side of code I wrote. I would start by searching for "epsilon" in the source code. "epsilon.mass" and "epsilon.enth" are the discrepancies from full conservation for both mass and enthalpy. As I mentioned before, in my experience *epsilon.mass* is numerically zero, whereas *epsilon.enth* is not. https://github.com/search?q=repo%3Acitibeth%2Ficebin%20epsilon&type=code I also had a place where I took the discrepancy of *epsilon.enth* and dumped it into the ocean, thereby keeping the entire coupled GCM / Ice Model conserving energy. In a separate email I will send you the 2014 AGU poster reporting on the code and results of evaluating conservation in PISM.

…

-- Elizabeth

On Mon, Dec 18, 2023 at 8:40 PM Ken Mankoff ***@***.***> wrote: PISM conserves mass. But should it? Am I misunderstanding what this means? Isn't it important to be able to lose mass across various domain boundaries? I'm curious which variables do capture this, and across which cells. If you prefer, I'm happy to break the long text below into multiple small Discussion posts. Is it reasonable to think of 3 planes, *top*, *bottom*, and *edges*? top is surface mass balance - I don't think I need to worry about surface mass balance because that is handled in our model, and we interface at the bottom of the firn. When we have surface melt in our model and no snow (e.g. blue ice in ablation zone) we track that and adjust the ice thickness by an appropriate amount for the next coupling step. I'd like to separate bottom into grounded and floating ice. For grounded ice, I see basal_mass_flux_grounded (units kg m-2 year-1). How is this different than tendency_of_subglacial_water_mass (units Gt year-1) other than units? I assume the latter is un-routed and the routed version(s) are tendency_of_subglacial_water_mass_at_grounded_margins ("subglacial water flux at grounded ice margins") or tendency_of_subglacial_water_mass_at_grounding_line ("subglacial water flux at grounding lines") and what is the difference between these latter two? What terms go into this basal mass balance? I assume geothermal heat flux and frictional heating from sliding? What about viscous heat dissipation of surface-runoff routed to the bed?. Is surface runoff routed to the bed? I assume not, because that implies treatment of moulins and a level of detail I did not think existed in PISM. For floating ice, do I subtract basal_mass_flux_grounded from tendency_of_ice_mass_due_to_basal_mass_flux' to get basal_mass_flux_floating`? How is edge mass loss tracked? I see tendency_of_ice_amount_due_to_discharge which is "discharge flux (calving, frontal melt, forced retreat)", and I assume this is a combination of tendency_of_ice_mass_due_to_calving + tendency_of_ice_amount_due_to_frontal_melt + tendency_of_ice_mass_due_to_forced_retreat. Constituent terms do not matter for mass conservation, but I think I need to work with constituent terms for energy conservation, because some of that is solid ice below the phase transition temperature (PTT) and some is liquid water I presume at the PTT. I need to track this properly into the appropriate GCM ocean cell(s). — Reply to this email directly, view it on GitHub <#527>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/AAOVY57XIBSWEDHTTLJD6BLYKESFPAVCNFSM6AAAAABA2TKSOKVHI2DSMVQWIX3LMV43ERDJONRXK43TNFXW4OZVHE3TONJTGA> . You are receiving this because you are subscribed to this thread.Message ID: ***@***.***>

-- Elizabeth Ashley Fischer, she/her University of Alaska, Geophysical Institute cell: 617-308-0436 UA is an AA/EO employer and educational institution and prohibits illegal discrimination against any individual: www.alaska.edu/nondiscrimination. UAF is located at Troth Yeddha', on the traditional homelands of the Lower Tanana Dené People and UAF therefore benefits from the long stewardship of these lands.

[mankoff]

Thanks. This addresses my first paragraph. I'm looking forward to hearing from PISM developers about the many other questions I asked. I also thought from discussions with you that icebin is for atmospheric coupling, and you said it had nothing to do with and should not be used for ocean coupling. Did I misunderstand that?

I note that there seems to be easier coupling methods than icebin. Kreuzer 2021 appear to use PISM as-is treating it like a black box but claim conservation of mass and energy.

[bueler]

Confirming what Elizabeth said, the meaning of "conservation of X", where X=mass,energy,momentum, is that the X quantity plus X flux "books" are balance-able. A general mathematical statement of conservation, typically an integral equation, should be read as: within some region, enclosed by surface S, the quantity X only changes by the identified source mechanisms within the region and by the identified fluxes across S.

You are correct that you should be able to sum various tendency_of_... PISM diagnostic outputs and come up with zero. For the PISM operational view of this see the User Manual page
https://www.pism.io/docs/manual/practical-usage/mass-conservation.html
Note that the goal is that the post-hoc computed quantity ice_mass_accounting_error described on that page should be zero within rounding error.

From the mathematical side, the fact that a time-stepping ice sheet or glacial model has to track a lateral moving glacier margin turns out to be a barrier to exact mass conservation. (This barrier is independent of spatial discretization.) In the PISM operational view, this means that PISM programmers cannot make tendency_of_ice_mass_due_to_conservation_error identically zero even by perfect numerical techniques and programming. This mathematical barrier is described in
https://doi.org/10.1137/20M135217X

[mankoff]

From that mass conservation page, what is the difference between tendency_of_subglacial_water_mass_at_grounded_margins and tendency_of_subglacial_water_mass_at_grounding_line?

[bueler]

The code that computes these is in method Hydrology::enforce_bounds() in file src/hydrology/Hydrology.cc, which is called on both the water in the till and the transportable subglacial water. (See src/hydrology/{Routing.cc|Distributed.cc} for where it is called.) The upshot is that the tendency_of_subglacial_water_mass_at_grounded_margins diagnostic is reporting water that leaves the subglacial system from a grounded glacier margin onto land above sea level, i.e. into rivers and lakes. (Which PISM mostly does not track, but there is https://doi.org/10.5194/tc-16-941-2022) And ..._at_grounding_line is water leaving the subglacial system into the ocean, i.e. below sea level. When sea level rises or falls there are additional complications.

[citibeth]

From the mathematical side, the fact that a time-stepping ice sheet or glacial model has to track a lateral moving glacier margin turns out to be a barrier to exact mass conservation. (This barrier is independent of spatial discretization.) In the PISM operational view, this means that PISM programmers cannot make tendency_of_ice_mass_due_to_conservation_error identically zero even by perfect numerical techniques and programming. This mathematical barrier is described in https://doi.org/10.1137/20M135217X

Math aside, my past experiments with PISM seemed to confirm exact convservation of mass. Or at least exact to some number of digits of precision that was satisfactory to me. Ed, do you know what gives here?

…

-- Elizabeth

[ckhroulev]

Ed's analysis (as I understand it) states that exact mass conservation cannot be achieved in general, which does not mean that PISM (and other models) never conserve mass. In particular, you have to have an ablation margin where changes in geometry due to ice flow are "balanced" by ablation (negative mass balance terms) for conservation issues to manifest themselves.

And even in that case you can perform accounting that "balances the books", in a way. If PISM is given the surface mass balance rate that is negative at some locations it can report how much mass it managed to remove in ablation areas before it ran out of ice to remove during a given time step. There will likely be a difference between "provided SMB" and "applied SMB" which constitutes a conservation error: after all, you told PISM to remove a certain amount of mass per time step and it failed to do so.

So, we can balance the books if we treat "applied SMB" as "what SMB should have been like" and ignore the difference mentioned above.

In addition to this, you got a bit lucky when you performed your experiments since you did not run into conservation issues described in [1]. These show up when bed topography is "rough" -- and this rarely happens at low grid resolutions. (Note that this issue is resolved in the current PISM version.)

[1]: A. H. Jarosch, C. G. Schoof, and F. S. Anslow, “Restoring mass conservation to shallow ice flow models over complex terrain,” The Cryosphere, vol. 7, no. 1, pp. 229–240, Feb. 2013, doi: 10.5194/tc-7-229-2013.

[citibeth]

I also thought from discussions with you that icebin is for atmospheric coupling, and you said it had nothing to do with and should not be used for ocean coupling. Did I misunderstand that?

That is true. But one would think that whatever techniques I used in icebin to check conservation could be re-used for ocean coupling.

…

-- Elizabeth

[crodehacke]

Conservation in a coupled earth system-ice sheet/shelf model system shall be archived if the following applies:

• The grids of the interacting models are aligned / identical.
• The surface elevations (interface atmosphere-ground/ice) are identical.

In common cases, the atmospheric model has a coarser resolution than the ice sheet, which leads to a difference in surface elevation if we interpolate from the coarse atmosphere to the finer resolved ice sheet grid. Under common cases, height corrections between these grids, such as the lapse-rate correction of the air temperature, are applied, violating strict energy conservation. Here, I do not take the smoothing of the surface elevation in the atmosphere model into account to reduce the complexity; as a consequence of the smoothing, the strong topographic gradients at the ice sheet margins are not well resolved in the atmosphere -- a topic of its own.

The coupling between the atmosphere and ice sheet surface requires the computation of the surface mass balance (SMB). Since we have different grids, we apply the given atmospheric conditions to obtain the SMB at different heights than the original grid. At the same time, we commonly do not really adjust the fluxes between the atmosphere and snow/ice surface.

Let's assume the atmosphere has delivered too much energy; we could feed it back by upward longwave radiation as well as sensible and latent heat flux. The first may be partly backscattered by clouds, atmospheric humidity, or tracer gases. The sensible heat warms the overlying atmosphere, which may dampen the common too-weak gravitational-driven downward winds (catabatic winds). And latent heat has its own cascade of effects counteracting partly the correction.

In the other case, where the surface has an energy deficit due to the height correction, we could artificially lower the albedo (even if it should not be lowered) to take up more solar radiation -- at least in the summer half year. Otherwise, we somehow cool the atmosphere.

In both cases, shall we remove/provide the positive/negative energy deficit at the height of the ice sheet or the atmosphere's bottom layer? - I drop the case where the ice sheet has a lower elevation than the surface elevation in the atmosphere.

Considering all these feedback loops and dirty fixes, are height corrections a viable solution if we require strict energy conservation while we do not want to interfere with the atmosphere conditions? On the other side, height corrections have been successful in providing SMB fields meeting our expectionations and observations estimates.

What does it mean for energy conservation? I do not know, but it might be hard.