Floating point overflow linked to snow thickness in SI3

[skeeley]

Context

We are having floating point exceptions in rare occasions in icethd.F90.

Analysis

The floating point exceptions happens in the following:

     WHERE( h_s_1d(1:npti)>0._wp ) e_s_1d(1:npti,jk) = e_s_1d(1:npti,jk) / (h_s_1d(1:npti) * a_i_1d(1:npti)) * nlay_s

where h_s_1d was extremely small (~1e-300 but >0.0) but a_i_1d identical to 0.0 causing a FPE.

Recommendation

While we can add a hack to check for very small values here to get the code running (and a similar hack in icedyn_adv_umx.F90) to make the code run then we suspect that this may be fixed more elegantly elsewhere by someone who knows the code better...
Happy to test potential fixes!

[clem]

I have modified one thermodynamic routine(attached here), so that snow thickness becomes 0 when ice completely disappears (I also made sure that we still conserve mass and heat). I do not think other routines should be changed.
However I cannot reproduce your FPE error. Could you try the routine before I commit it?
Thanks

[skeeley]

Thanks Clem!
We have the code in for testing now - will let you know.

[clem]

In 10785:

fix tickets #2256 and #2257

[clem]

In 10786:

fix tickets #2256 and #2257