Temp/salinity and thickness bug when initializing with an ice shelf

[claireyung]

Hello, I've come across strange features in the initialization of temperature/salinity and layer thickness/interface levels when an ice shelf is present. I initialized the model with a simple triangular ice shelf and TS_CONFIG = "ISOMIP" and THICKNESS_CONFIG = "ISOMIP" (though the same happens with other config options, e.g. TS=linear or Thickness=uniform) which is meant to produce a linear stratification according to given temperature and salinity bounds. However, although it looks mostly linear, there are strange periodic deviations from linearity in the initial conditions file which are highlighted below in this figure, showing the anomaly of the salt initial conditions from the theoretical linear profile.

Sea level/layer thicknesses also have very slight perturbations from their expected smooth distributions (for sigma_shelf_zstar coordinates) when TRIM_IC_FOR_P_SURF is not used (turning it on seems to fix those discontinuities but not the T/S ones). When I run the model, horizontal velocities emerge very quickly with the same spatial frequency as the salt anomalies suggesting that the anomalous horizontal density gradients drive this unwanted flow.

A temporary fix: I can get around this issue by using netcdf files as input for the thicknesses (with ice shelf accounted for), and by asking the model to skip the calc_sfc_displacement step in MOM_state_initialization.F90. This somehow results in a linear T/S profile and very low (order 10^(-8) m/s) velocities for the simple triangular ice shelf, however requires running the model twice. If it is helpful, I have linked a github repo that includes the files/MOM_input/override used in my test cases.