Fix for sigma_i going to zero in ITD landfast ice.

[kshedstrom]

• We can't divide by zero, so don't let sigma_i be zero.

[Hallberg-NOAA]

I do not understand why this change is necessary. When I look at the code as a whole, I see that the expression for sigma_i about 20 lines above this change, sigma_i = max(sqrt(log(1.0 + v_i/m_i**2)), CS%puny) at line 1959 of this same file, already includes a max with CS%puny, and I am not seeing anythng between these lines that changes the value of sigma_i. What am I missing here?

[kshedstrom]

Yes, you are right. I don't remember how I thought sigma_i was getting a zero, except that it was thin new ice.

[kshedstrom]

Now I remember - thin new ice lead to v_i being zero.

[kshedstrom]

Sorry, @Hallberg-NOAA, getting my branches confused. How about this one?

[Hallberg-NOAA]

This change makes much more sense than the previous version!

[marshallward]

@kshedstrom Could you do something like this?

if (sigma_i > 0) then
   g_k(:) = exp(-(log(x_k(:)/CS%onemeter) - mu_i) ** 2 / (2.0 * sigma_i ** 2)) / &
                 (x_k(:) * sigma_i * sqrt(2.0 * pi))
                 
   ...
   Tbt(:) = ...
else
   ! g_k(:) = 0., since exp(-a**2/n**2)/n -> 0 as n->0
   ! P_x(:) = 0.
   ! tb_tmp(:) = 0.

   Tbt(:) = 0.
endif

It would avoid the need for CS%puny here, and might also reduce the computation.

[kshedstrom]

This isn't the exact form I tested, given that I need the new SIS_restart code. It's in the queue for a longer test.

[marshallward]

Thanks for testing it out, let's see what happens 😅.

However, I think the better test would be to remove puny_m2, rather than introducing it alongside the sigma_i > 0 test. sigma_i is only zero if v_i is zero, so if you bound v_i then I would assume that it's not even possible for sigma_i to be zero.

I was thinking that the sigma_i > 0 test would avoid the need to introduce another "puny" parameter, but I suppose that either approach gets you to the same place.

[kshedstrom]

Try that one?

[marshallward]

That's what I had in mind, I guess the question is whether it fixed your problem without creating any new ones?

[kshedstrom]

My short job ran fine, a longer one is stuck in the queue.

[kshedstrom]

Yes, the year-long job finished successfully!

[marshallward]

Gaea regression: https://gitlab.gfdl.noaa.gov/ogrp/SIS2/-/pipelines/19489 ✔️

[marshallward]

@kshedstrom Can I squash-merge this one?

[kshedstrom]

Yes, of course.