[MITgcm-devel] potential problem in kpp_routines.F, S/R wscale

[Martin Losch]

Hi there, 

just for later reference: There may be a bug in S/R wscale that has been found by my FESOM colleague Dimitry Sidorenko and has been simultaneously discussed in a mom4 users group, e.g. here: <https://groups.google.com/forum/#!msg/mom-users/XGyu_EWvhPA/1s3FpEChDQAJ> 
All models share essentially the same KPP code.

Apparently, the statement
 zdiff = zehat - zmin
can lead to negative zdiff (in fact, I checked and this happens very often). zmin = -4e-7 by default. zdiff is then used to interpolate between values of a lookup table (wmt and wst), which is defined only between zmin = -4.e-7 and zmax = 0. For zdiff < 0, this leads to a linear extrapolation beyond the lower limit of the lookup table. In most cases this is fine, because zdiff is small (i.e. > -1.e-7), but for extreme values of zehat (i.e. extreme values of negative buoyance forcing bfsfc), this can lead to numerical problems and even floating point exceptions.

Dima’s fix to this is to replace the linear extrapolation to a nearest neighbor interpolation by setting zdiff = max(0, zehat - zmin). I will not implement this fix, because it will change the results, but we might want to keep it in mind, in case someone stumbles of this problem. Or what do you think?

Martin