The surface water flux is included in the temperature equation via

   DiffTup = -T(nlev)*wflux-hflux/(rho0*cp)

where wflux is the water flux.

just above this line, I read the comment

!  Note that this is the surface temperature flux for rigid-lid models like GOTM.
!  For a free surface model the surface temperature flux must be - -hflux/(rho_0*cp)

I see that this wflux term was added in a commit to include the plume dynamics 1e65f75

As mentioned in the comment I am not sure this term is relevant for the open ocean. Would it be possible to add a flag to choose if we keep this term in the temperature equation?

there is already a flux_impact flag in the precip section of the .yaml file.

I'd suggest to add another flag like free_surface_fresh_water (set by default to true -- i.e. revert to pre-plume commit). but any other name is fine, and default to false is also fine as long as it is well documented.