Hi everyone,
I would like to report a diagnostic-output bug and a verified fix for CROCO v.2.1.3 (released April 2026), which we encountered while running a regional configuration of the Arabian Gulf using the BULK_FLUX, RAIN_FLUX and AVERAGES cpp keys. It is found that when RAIN_FLUX is active, the averaged swflux (surface freshwater flux) output variable is written as zero in the output files, even though the freshwater forcing is correctly applied internally.
Identified cause: In bulk_flux.F, under #ifdef RAIN_FLUX, the net freshwater flux E-P is computed and stored in the EmP array, then applied as a real surface volume flux. In that branch, stflx(i,j,isalt) is never assigned the line that sets it, stflx(isalt) = (evap-prate)*S, sits in the #else branch and is compiled out. The averaging routine set_avg.F accumulates stflx_avg(:,:,isalt), which wrt_avg.F writes to swflux. Since stflx(isalt) remains at its initialized zero under RAIN_FLUX, the output is identically zero. The EmP array itself has no output plumbing in set_avg.F.
Verified Fix (in set_avg.F): accumulate stf_cff*EmP under RAIN_FLUX instead of the salinity-weighted salt flux in set_avg.F. EmP is already E-P in m/s; the non-RAIN_FLUX branch divides stflx(isalt) by surface salinity to recover the same quantity, so both paths yield identical units (cm/day via stf_cff). The non-RAIN_FLUX path is unchanged. Kindly look at the following snippet.
After the fix, the output swflux now shows a physically sensible E−P field positive (net evaporation) over the basin with the expected spatial structure.
