[MITgcm-devel] [MITgcm-cvs] MITgcm/pkg/seaice CVS Commit

[Martin Losch]

Hi Gael,

my preliminary conclusion about the "vintage bug" (I am waiting for Jinlun to confirm or contradict/explain) is, contrary to my previous statements, that is really a bug. And a hard one to catch, too!
On Oct 16, 2010, at 2:13 AM, Gael Forget wrote:

>> 
>> Obviously I could not convince you that the "vintage bug" really isn't one. How can you convince me? I guess we should ask Jinlun how to do it right. And then equivalent part for snow is treated differently ...
> It is very well possible that I did not get it right. But the more I circled around it, the more it felt 
> like a bug -- very much like the one you fixed in v1.34 of growth.F (Nov 2006) in fact. Maybe I just did not get the hsnow/heff treatment difference that you are referring to, and we will need a fix. I reported the questionable operation as a comment in seaice_growth.F  v1.89. Let's come back to this question in a few days, if that is ok with you.

I explain it now like this (forget about snow for a moment):
h: ice thickness
c: concentration
heff=h*c
qi,qo = heat flux out of ocean for ice part and for ocean part, so that qnet = qi*c + qo*(1-c), the sign convention is positive upward
qa = heat from atmosphere available to melt ice, we are interested in the case when qa is large enough to melt all ice and there is some residual heat left over to pass to the ocean.
lam = heat of fusion (conversion factor between ice thickness of ice and heat flux)

Then for the ice covered part we have
dh/dt = -qa/lam
but because we want h=0 after this step we compute the residual
qi = (h - dt*qa/lam)*lam/dt = h*lam/dt - qa
this is the part of the ice. If you now do everything with heff, then
dheff/dt = -qa/lam*c
and
qi*c = (heff-dt*qa/lam*c)*lam/dt = heff*lam/dt - qa*c
and clearly the residual heat already contains the fractional area, so that.

Fortunately this problem is not severe in practice: if you weight with c^2 instead of c, then errror maximally max(c-c^2) = .25 (for c=[0 1]), so 25% of qi at c=.5, but for c near 1 or 0 (more often the case than c=.5) the error is much smaller. I ran a comparison with my coases (cs32) coupled embm-ocean system and afer 350 years the streamfunctions are nearly indistiguishable. maximum ice thickness differences are < 4cm (although ice tends to be a little thicker on average, <1cm). MLD differs by less than 25m in a few places (remember the coarse vertical resolution of 15layers in this run). Differences might be much higher in high resolution runs, but I am not too concerned. I think that Dimitris should check this in one of his high-resolution experiments.

The missing "area" in evaportation (growth.F v1.34) was more severe, because the error was 100% underneath the ice.

Again: nice catch and I think we can get rid of the comment.

Martin