diff --git a/src/core_ocean/shared/mpas_ocn_diagnostics.F b/src/core_ocean/shared/mpas_ocn_diagnostics.F
index c0516260f7..ad7729bab6 100644
--- a/src/core_ocean/shared/mpas_ocn_diagnostics.F
+++ b/src/core_ocean/shared/mpas_ocn_diagnostics.F
@@ -1367,7 +1367,7 @@ subroutine ocn_compute_KPP_input_fields(statePool, forcingPool, meshPool, diagno
       real (kind=RKIND), dimension(:), pointer :: dcEdge, dvEdge, areaCell
       real (kind=RKIND), dimension(:), pointer :: penetrativeTemperatureFlux, surfaceThicknessFlux, &
            surfaceBuoyancyForcing, surfaceFrictionVelocity, penetrativeTemperatureFluxOBL, &
-           normalVelocitySurfaceLayer, surfaceThicknessFluxRunoff
+           normalVelocitySurfaceLayer, surfaceThicknessFluxRunoff, rainFlux, evaporationFlux
       real (kind=RKIND), pointer :: config_flux_attenuation_coefficient, config_flux_attenuation_coefficient_runoff
 
       real (kind=RKIND), dimension(:), pointer :: surfaceStress, surfaceStressMagnitude
@@ -1441,6 +1441,8 @@ subroutine ocn_compute_KPP_input_fields(statePool, forcingPool, meshPool, diagno
       call mpas_pool_get_array(forcingPool, 'penetrativeTemperatureFlux', penetrativeTemperatureFlux)
       call mpas_pool_get_array(forcingPool, 'surfaceStress', surfaceStress)
       call mpas_pool_get_array(forcingPool, 'surfaceStressMagnitude', surfaceStressMagnitude)
+      call mpas_pool_get_array(forcingPool, 'rainFlux', rainFlux)
+      call mpas_pool_get_array(forcingPool, 'evaporationFlux', evaporationFlux)
 
       call mpas_pool_get_array(tracersPool, 'activeTracers', activeTracers, 1)
       call mpas_pool_get_array(tracersSurfaceFluxPool, 'activeTracersSurfaceFlux', activeTracersSurfaceFlux)
@@ -1496,12 +1498,14 @@ subroutine ocn_compute_KPP_input_fields(statePool, forcingPool, meshPool, diagno
        ! transport code.  This includes tracer forcing due to thickness
 
        nonLocalSurfaceTracerFlux(indexTempFlux, iCell) = activeTracersSurfaceFlux(indexTempFlux,iCell) &
-               + penetrativeTemperatureFlux(iCell) - penetrativeTemperatureFluxOBL(iCell) - fracAbsorbed * &
-               surfaceThicknessFlux(iCell) * activeTracers(indexTempFlux,1,iCell) + &
-               activeTracersSurfaceFluxRunoff(indexTempFlux,iCell) * fracAbsorbedRunoff
+               + penetrativeTemperatureFlux(iCell) - penetrativeTemperatureFluxOBL(iCell)  &
+               - fracAbsorbed * (rainFlux(iCell) + evaporationFlux(iCell)) * activeTracers(indexTempFlux,1,iCell)/rho_sw  &
+               - fracAbsorbedRunoff * surfaceThicknessFluxRunoff(iCell) &
+                * min(activeTracers(indexTempFlux,1,iCell),0.0_RKIND)/rho_sw
 
        nonLocalSurfaceTracerFlux(indexSaltFlux,iCell) = activeTracersSurfaceFlux(indexSaltFlux,iCell) &
-               - fracAbsorbed * surfaceThicknessFlux(iCell) * activeTracers(indexSaltFlux,1,iCell)
+               - fracAbsorbed * surfaceThicknessFlux(iCell) * activeTracers(indexSaltFlux,1,iCell) &
+               - fracAbsorbedRunoff * surfaceThicknessFluxRunoff(iCell) * activeTracers(indexSaltFlux,1,iCell)
 
        surfaceBuoyancyForcing(iCell) =  thermalExpansionCoeff (1,iCell) &
               * nonLocalSurfaceTracerFlux(indexTempFlux,iCell) &
diff --git a/testing_and_setup/compass/ocean/global_ocean/template_forward.xml b/testing_and_setup/compass/ocean/global_ocean/template_forward.xml
index 36d21a9a47..0af5703b74 100644
--- a/testing_and_setup/compass/ocean/global_ocean/template_forward.xml
+++ b/testing_and_setup/compass/ocean/global_ocean/template_forward.xml
@@ -14,6 +14,7 @@
 		<option name="config_eos_type">'jm'</option>
 		<option name="config_implicit_bottom_drag_coeff">1.0e-3</option>
 		<option name="config_use_bulk_wind_stress">.true.</option>
+		<option name="config_use_bulk_thickness_flux">.true.</option>
 		<option name="config_use_activeTracers_surface_restoring">.true.</option>
 	</namelist>
 </template>