Add meltpond coupling

[kieranricardo]

Add support for coupling melt ponds between UM and CICE6.

1. Export meltpond fractions from CICE. https://github.com/ACCESS-NRI/CICE/blob/04fbebc897fbd319a1de9b0359ca7bda3490d259/cicecore/drivers/access/cmeps/ice_import_export.F90#L2078
2. Couple fields through the mediator (already done): https://github.com/ACCESS-NRI/CMEPS/blob/d80e784c2208d4f6faeb72c23d0a7cf7ceac1a49/mediator/esmFldsExchange_access_mod.F90#L415
3. Import fields into UM: https://github.com/ACCESS-NRI/cm3-um/blob/b925487e9745af3584dcb9f995811edb2dd43d7d/src/control/top_level/um_cap_import_export.F90#L177
4. Turn on melt ponds for UM (and CICE): https://github.com/ACCESS-NRI/cm3-suite/blob/main/app/um/rose-app.conf
5. Process the import fields: https://github.com/ACCESS-NRI/cm3-um/blob/b925487e9745af3584dcb9f995811edb2dd43d7d/src/control/top_level/um_cap_import_export.F90#L987. Should match the UM oasis coupling here: https://github.com/ACCESS-NRI/cm3-um/blob/main/src/control/coupling/oasis3_geto2a.F90.

[anton-seaice]

Can we turn them off ? Is the UM expecting melt ponds ?

[kieranricardo]

@anton-seaice oh I was under the impression we wanted melt ponds, if not that makes life easier.

The UM isn't expecting melt ponds (but we could turn this on), and we're currently using the tr_pond_lvl scheme with the config here

[anton-seaice]

They should all be off, basically there aren't meltponds in the Antarctic, and previous attempts at using them have overrepresented their impacts in the Antarctic.

We changed this here for OM3:

https://github.com/ACCESS-NRI/access-om3-configs/commit/c45fb289997d82a719b4edc44ab9f2edd729d16e#diff-46cbcfe1a092910fbb4c812ff0cf91d611cd0c5f4201b86c82dfffb6a26ab74aR53

I suggest updating the whole config to the om3-configs, but in the interim, you could just remove tr_pond_lvl as the default is .false.

[blimlim]

Following further discussion, it sounds like we want to add the meltpond coupling. I've had a go at following the instructions in the issue description. Branches with the changes are available here:

• CICE: https://github.com/ACCESS-NRI/CICE/compare/cice-export...meltpond-coupling
• UM: https://github.com/ACCESS-NRI/cm3-um/compare/main...37-meltpond-coupling
• cm3-suite: https://github.com/ACCESS-NRI/cm3-suite/compare/main...meltpond-coupling

A lot of it is still a bit over my head, and so there is a chance I've missed something here. The above changes have no impact however. The following shows the surface temperature on the 31st day of simulations with and without the changes

I'm not sure whether anything is done by the UM with the pond_frac_cat and pond_depth_cat variables once they're imported. The only place I can see them outside of the variable definition routines is in rad_ctl.F90 https://github.com/ACCESS-NRI/cm3-um/blob/3eb408bc535ed183e32054d683c62541bced09d9/src/atmosphere/radiation_control/rad_ctl.F90#L2306-L2309 where they are used to set ainfo%pond_frac_cat_sicat and ainfo%pond_depth_cat_sicat, however these don't seem to be used anywhere else.

In /g/data/access/rose-meta/um/um-atmos/vn13.0/rose-meta.conf it has

[namelist:jules_sea_seaice=l_sice_meltponds_cice]
compulsory=true
description=Sea-ice albedo is affected by meltponds (from CICE meltponds scheme)
sort-key=35

I'm wondering whether the meltponds impacting albedo is part of JULES, and so not having any effect for us?

[anton-seaice]

I think JULES doesn't include Antarctica ?

If you have a look at the f_apond_ai field in the CICE output you should get a sense of if CICE is generating any meltponds? If you only ran for a month, it might just be the model needs to be run longer?

Also you might need to turn on restart_pond_lvl = .true. in ice_in

[blimlim]

Todo from meeting:

• [ ] Switch from tr_pond_lvl=.true. to tr_pond_topo=.true., also use restart_pond_topo=.true.
• [ ] Switch to exporting trcrn(:,:,nt_apnd_out,:,:) rather than apeffn for the melt pond fraction.

[blimlim]

With the current coupling and pond scheme, something is also incorrect with the pond depth trcrn(:,:,nt_hpnd,:,:) export to the mediator. The mediator receives zero everywhere while the pond depth in cice is not zero:

[anton-seaice]

Just confirming that there is mediator output on other fields ? We had an issue at one point where the mediator output didn't work at all?

Also - wondering if trcrn(:,:,nt_hpnd,:,:) returns a 4d array ? Maybe try writing it to a temporary array ?

[blimlim]

Sorry about this – turns out the issue was me not understanding how to build properly, and the code changes didn't make their way into the runs... After doing that properly and swapping to the tr_pond_topo scheme we get the following pond fraction and pond depth imports/exports for ice category 1 in the mediator:

(I've rolled the ice export data so that it better lines up with the atmosphere in the plots)

The pond fraction to the mediator very closely matches apondn in the ice history:

But the pond depth sent to the mediator is greater than the hpondn history variable:

I think because the hpondn output weights the pond depth by the pond fraction to get the average pond depth over the sea ice.

Need to try and understand what the UM requires when the l_sice_meltponds_cice is switched on. E.g. The raw pond area tracer trcrn(:,:,nt_apnd,1:ncat_hist,iblk) has a different meaning for tr_pond_topo vs tr_pond_level. For tr_pond_topo it's the fraction of ice covered by ponds, and for tr_pond_level, it's the fraction of level ice covered by ponds (useful info here).

[ofa001]

HI @blimlim I will have a look next week, I thought the pond fraction looked quite high, for cat 1 and 2 mor like what total ice conc might be. Might still be something not quite right but progress over the levl ice set up.

[blimlim]

Thanks @ofa001. I've had a look into the JULES jules_ssi_albedo routine and am putting together some notes with the details.=. My understanding is that we should send both the pond fraction and pond depth without any weighting when using the topo scheme - which would correspond to a_pnd and h_pnd on this page, and I think these would be trcrn(:,:,nt_apnd,:,:) and trcrn(:,:,nt_hpnd,:,:) in the code.

[ofa001]

HI @bimlim It will be best to line it up with what we sent in ACCESS-CM2, I will get you the relevent lines of code in a couple of hours after my next set of meetings, and will check if/when weighting where applied on the CICE or UM side of the coupling. Glad you are looking at some of the JULES routines!

[blimlim]

Some todos/questions to resolve from the meeting:

• [ ] Run for a couple of years and investigate ponds further from the initialisation – do we still get large pond fractions?
• [ ] Check pond areas are conserved through the coupling
• [ ] Abort simulations if tr_pond_lvl is selected.
• [ ] Incorporate ice lids into pond areas
• [ ] Multiple ice lid parameterisations available – which are we using/should we use?
• [ ] Send effective pond fraction apeffn rather than raw pond fraction apndn to incorporate the automatic effects of ice lids

[ofa001]

If we are using tr_pond_topo, which is what I would recommend then icepack_shortwave has this section of code already set up. So if the lid is

Calculate effective pond area for HadGEM

     if (tr_pond_topo) then
        do n = 1, ncat
           apeffn(n) = c0
           if (aicen(n) > puny) then
           ! Lid effective if thicker than hp1
             if (apndn(n)*aicen(n) > puny .and. ipndn(n) < hp1) then
                apeffn(n) = apndn(n)
             else
                apeffn(n) = c0
             endif
             if (apndn(n) < puny) apeffn(n) = c0
           endif
        enddo  ! ncat

[ofa001]

The sentence above should read. "so if the lid is > hp (=0.01) then the system already sets the apeffn=0. So you don't need to change". Apologies, I had a system problem and need to reboot, after several trues to edit.

[ofa001]

The other lid scheme is the Stefan approximation is part of icepack_melpond_lvl scheme but we are going to set a switch to abort if that option is chosen, so I think the scheme should already include lids if they are forming in the topo scheme but it would be good to check that it is in practise.

[blimlim]

Hi @ofa001, thanks for looking into this! It turns out I'm currently sending the raw pond fraction apndn rather than the effective pond fraction apeffn, meaning that the effects of the ice lids aren't showing up in our results. This was due to a mistake I made earlier on. I thought apeffn was only set when calc_Tsfc=.true. and hence swapped to sending the raw pond fraction apndn, but this was incorrect! Apologies about this!

Swapping to apeffn should help a lot with the large pond fractions appearing on the thinner ice. The following shows the raw pond fraction 'apndn, the effective fractionapeffn`, and their difference on Jan 30 of the first year:

I've added fixing this to my todo list above.

[ofa001]

Hi @blimlim, Thats good news its a big difference, I thought I would have seen values like this before in CM2.

I am going to look into the ice thickness issue next, again it may not be as big as your plot showed, as it will only effect the lowest ice thickness category, so around the ice edge mostly, and I am wondering if I am mixing up our old CSIRO model where our multi layer model did become 'zero layer' as it got thinner, same with snow models (in pre-CABLE set up).