diff --git a/config_src/coupled_driver/MOM_surface_forcing_gfdl.F90 b/config_src/coupled_driver/MOM_surface_forcing_gfdl.F90 index 61d9c60d1d..766f2127c6 100644 --- a/config_src/coupled_driver/MOM_surface_forcing_gfdl.F90 +++ b/config_src/coupled_driver/MOM_surface_forcing_gfdl.F90 @@ -242,6 +242,8 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, sfc real :: delta_sss ! temporary storage for sss diff from restoring value [ppt] real :: delta_sst ! temporary storage for sst diff from restoring value [degC] + real :: kg_m2_s_conversion !< A combination of unit conversion factors for rescaling + !! mass fluxes [R Z s m2 kg-1 T-1 ~> 1]. real :: C_p ! heat capacity of seawater [J degC-1 kg-1] real :: sign_for_net_FW_bug ! Should be +1. but an old bug can be recovered by using -1. @@ -255,6 +257,7 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, sfc IsdB = G%IsdB ; IedB = G%IedB ; JsdB = G%JsdB ; JedB = G%JedB isr = is-isd+1 ; ier = ie-isd+1 ; jsr = js-jsd+1 ; jer = je-jsd+1 + kg_m2_s_conversion = US%kg_m3_to_R*US%m_to_Z*US%T_to_s C_p = fluxes%C_p open_ocn_mask(:,:) = 1.0 pme_adj(:,:) = 0.0 @@ -372,19 +375,21 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, sfc delta_sss = sfc_state%SSS(i,j) - data_restore(i,j) delta_sss = sign(1.0,delta_sss)*min(abs(delta_sss),CS%max_delta_srestore) fluxes%vprec(i,j) = (CS%basin_mask(i,j)*open_ocn_mask(i,j)*CS%srestore_mask(i,j))* & - (US%R_to_kg_m3*CS%Rho0*CS%Flux_const) * & + (US%m_to_Z*US%T_to_s*CS%Rho0*CS%Flux_const) * & delta_sss / (0.5*(sfc_state%SSS(i,j) + data_restore(i,j))) endif enddo ; enddo if (CS%adjust_net_srestore_to_zero) then if (CS%adjust_net_srestore_by_scaling) then - call adjust_area_mean_to_zero(fluxes%vprec, G, fluxes%vPrecGlobalScl) + call adjust_area_mean_to_zero(fluxes%vprec, G, fluxes%vPrecGlobalScl, & + unit_scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) fluxes%vPrecGlobalAdj = 0. else - work_sum(is:ie,js:je) = US%L_to_m**2*G%areaT(is:ie,js:je)*fluxes%vprec(is:ie,js:je) + work_sum(is:ie,js:je) = US%L_to_m**2*G%areaT(is:ie,js:je) * & + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%vprec(is:ie,js:je) fluxes%vPrecGlobalAdj = reproducing_sum(work_sum(:,:), isr, ier, jsr, jer) / CS%area_surf do j=js,je ; do i=is,ie - fluxes%vprec(i,j) = ( fluxes%vprec(i,j) - fluxes%vPrecGlobalAdj ) * G%mask2dT(i,j) + fluxes%vprec(i,j) = ( fluxes%vprec(i,j) - kg_m2_s_conversion*fluxes%vPrecGlobalAdj ) * G%mask2dT(i,j) enddo ; enddo endif endif @@ -408,31 +413,31 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, sfc do j=js,je ; do i=is,ie if (associated(IOB%lprec)) then - fluxes%lprec(i,j) = IOB%lprec(i-i0,j-j0) * G%mask2dT(i,j) + fluxes%lprec(i,j) = kg_m2_s_conversion * IOB%lprec(i-i0,j-j0) * G%mask2dT(i,j) if (CS%check_no_land_fluxes) & call check_mask_val_consistency(IOB%lprec(i-i0,j-j0), G%mask2dT(i,j), i, j, 'lprec', G) endif if (associated(IOB%fprec)) then - fluxes%fprec(i,j) = IOB%fprec(i-i0,j-j0) * G%mask2dT(i,j) + fluxes%fprec(i,j) = kg_m2_s_conversion * IOB%fprec(i-i0,j-j0) * G%mask2dT(i,j) if (CS%check_no_land_fluxes) & call check_mask_val_consistency(IOB%fprec(i-i0,j-j0), G%mask2dT(i,j), i, j, 'fprec', G) endif if (associated(IOB%q_flux)) then - fluxes%evap(i,j) = - IOB%q_flux(i-i0,j-j0) * G%mask2dT(i,j) + fluxes%evap(i,j) = - kg_m2_s_conversion * IOB%q_flux(i-i0,j-j0) * G%mask2dT(i,j) if (CS%check_no_land_fluxes) & call check_mask_val_consistency(IOB%q_flux(i-i0,j-j0), G%mask2dT(i,j), i, j, 'q_flux', G) endif if (associated(IOB%runoff)) then - fluxes%lrunoff(i,j) = IOB%runoff(i-i0,j-j0) * G%mask2dT(i,j) + fluxes%lrunoff(i,j) = kg_m2_s_conversion * IOB%runoff(i-i0,j-j0) * G%mask2dT(i,j) if (CS%check_no_land_fluxes) & call check_mask_val_consistency(IOB%runoff(i-i0,j-j0), G%mask2dT(i,j), i, j, 'runoff', G) endif if (associated(IOB%calving)) then - fluxes%frunoff(i,j) = IOB%calving(i-i0,j-j0) * G%mask2dT(i,j) + fluxes%frunoff(i,j) = kg_m2_s_conversion * IOB%calving(i-i0,j-j0) * G%mask2dT(i,j) if (CS%check_no_land_fluxes) & call check_mask_val_consistency(IOB%calving(i-i0,j-j0), G%mask2dT(i,j), i, j, 'calving', G) endif @@ -565,7 +570,8 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, sfc sign_for_net_FW_bug = 1. if (CS%use_net_FW_adjustment_sign_bug) sign_for_net_FW_bug = -1. do j=js,je ; do i=is,ie - net_FW(i,j) = (((fluxes%lprec(i,j) + fluxes%fprec(i,j)) + & + net_FW(i,j) = US%R_to_kg_m3*US%Z_to_m*US%s_to_T* & + (((fluxes%lprec(i,j) + fluxes%fprec(i,j)) + & (fluxes%lrunoff(i,j) + fluxes%frunoff(i,j))) + & (fluxes%evap(i,j) + fluxes%vprec(i,j)) ) * US%L_to_m**2*G%areaT(i,j) ! The following contribution appears to be calculating the volume flux of sea-ice @@ -583,13 +589,13 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, sfc if (CS%adjust_net_fresh_water_by_scaling) then call adjust_area_mean_to_zero(net_FW2, G, fluxes%netFWGlobalScl) do j=js,je ; do i=is,ie - fluxes%vprec(i,j) = fluxes%vprec(i,j) + & + fluxes%vprec(i,j) = fluxes%vprec(i,j) + US%kg_m3_to_R*US%m_to_Z*US%T_to_s * & (net_FW2(i,j) - net_FW(i,j)/(US%L_to_m**2*G%areaT(i,j))) * G%mask2dT(i,j) enddo ; enddo else fluxes%netFWGlobalAdj = reproducing_sum(net_FW(:,:), isr, ier, jsr, jer) / CS%area_surf do j=js,je ; do i=is,ie - fluxes%vprec(i,j) = ( fluxes%vprec(i,j) - fluxes%netFWGlobalAdj ) * G%mask2dT(i,j) + fluxes%vprec(i,j) = ( fluxes%vprec(i,j) - kg_m2_s_conversion * fluxes%netFWGlobalAdj ) * G%mask2dT(i,j) enddo ; enddo endif @@ -611,7 +617,7 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, sfc if (CS%allow_flux_adjustments) then ! Apply adjustments to fluxes - call apply_flux_adjustments(G, CS, Time, fluxes) + call apply_flux_adjustments(G, US, CS, Time, fluxes) endif ! Allow for user-written code to alter fluxes after all the above @@ -1086,8 +1092,9 @@ end subroutine extract_IOB_stresses !! - hflx_adj (Heat flux into the ocean [W m-2]) !! - sflx_adj (Salt flux into the ocean [kg salt m-2 s-1]) !! - prcme_adj (Fresh water flux into the ocean [kg m-2 s-1]) -subroutine apply_flux_adjustments(G, CS, Time, fluxes) +subroutine apply_flux_adjustments(G, US, CS, Time, fluxes) type(ocean_grid_type), intent(inout) :: G !< Ocean grid structure + type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(surface_forcing_CS), pointer :: CS !< Surface forcing control structure type(time_type), intent(in) :: Time !< Model time structure type(forcing), intent(inout) :: fluxes !< Surface fluxes structure @@ -1120,7 +1127,7 @@ subroutine apply_flux_adjustments(G, CS, Time, fluxes) call data_override('OCN', 'prcme_adj', temp_at_h(isc:iec,jsc:jec), Time, override=overrode_h) if (overrode_h) then ; do j=jsc,jec ; do i=isc,iec - fluxes%vprec(i,j) = fluxes%vprec(i,j) + temp_at_h(i,j)* G%mask2dT(i,j) + fluxes%vprec(i,j) = fluxes%vprec(i,j) + US%kg_m3_to_R*US%m_to_Z*US%T_to_s * temp_at_h(i,j)* G%mask2dT(i,j) enddo ; enddo ; endif ! Not needed? ! if (overrode_h) call pass_var(fluxes%vprec, G%Domain) end subroutine apply_flux_adjustments diff --git a/config_src/coupled_driver/ocean_model_MOM.F90 b/config_src/coupled_driver/ocean_model_MOM.F90 index c5d10c7aaf..9982754053 100644 --- a/config_src/coupled_driver/ocean_model_MOM.F90 +++ b/config_src/coupled_driver/ocean_model_MOM.F90 @@ -520,7 +520,7 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, time_start_upda if (OS%use_ice_shelf) & call shelf_calc_flux(OS%sfc_state, OS%fluxes, OS%Time, dt_coupling, OS%Ice_shelf_CSp) if (OS%icebergs_alter_ocean) & - call iceberg_fluxes(OS%grid, OS%fluxes, OS%use_ice_shelf, & + call iceberg_fluxes(OS%grid, OS%US, OS%fluxes, OS%use_ice_shelf, & OS%sfc_state, dt_coupling, OS%marine_ice_CSp) #ifdef _USE_GENERIC_TRACER @@ -541,7 +541,7 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, time_start_upda if (OS%use_ice_shelf) & call shelf_calc_flux(OS%sfc_state, OS%flux_tmp, OS%Time, dt_coupling, OS%Ice_shelf_CSp) if (OS%icebergs_alter_ocean) & - call iceberg_fluxes(OS%grid, OS%flux_tmp, OS%use_ice_shelf, & + call iceberg_fluxes(OS%grid, OS%US, OS%flux_tmp, OS%use_ice_shelf, & OS%sfc_state, dt_coupling, OS%marine_ice_CSp) call fluxes_accumulate(OS%flux_tmp, OS%fluxes, dt_coupling, OS%grid, weight) @@ -554,7 +554,7 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, time_start_upda ! The net mass forcing is not currently used in the MOM6 dynamics solvers, so this is may be unnecessary. if (do_dyn .and. associated(OS%forces%net_mass_src) .and. .not.OS%forces%net_mass_src_set) & - call get_net_mass_forcing(OS%fluxes, OS%grid, OS%forces%net_mass_src) + call get_net_mass_forcing(OS%fluxes, OS%grid, OS%US, OS%forces%net_mass_src) if (OS%use_waves .and. do_thermo) then ! For now, the waves are only updated on the thermodynamics steps, because that is where @@ -654,7 +654,7 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, time_start_upda if (OS%fluxes%fluxes_used .and. do_thermo) then call enable_averaging(OS%fluxes%dt_buoy_accum, OS%Time, OS%diag) call forcing_diagnostics(OS%fluxes, OS%sfc_state, OS%fluxes%dt_buoy_accum, & - OS%grid, OS%diag, OS%forcing_CSp%handles) + OS%grid, OS%US, OS%diag, OS%forcing_CSp%handles) call disable_averaging(OS%diag) endif diff --git a/config_src/ice_solo_driver/MOM_surface_forcing.F90 b/config_src/ice_solo_driver/MOM_surface_forcing.F90 index f86fc44101..ea3385e88e 100644 --- a/config_src/ice_solo_driver/MOM_surface_forcing.F90 +++ b/config_src/ice_solo_driver/MOM_surface_forcing.F90 @@ -275,7 +275,7 @@ subroutine set_forcing(sfc_state, forcing, fluxes, day_start, day_interval, G, U if ((CS%variable_buoyforce .or. CS%first_call_set_forcing) .and. & (.not.CS%adiabatic)) then - call set_net_mass_forcing(fluxes, forces, G) + call set_net_mass_forcing(fluxes, forces, G, US) endif CS%first_call_set_forcing = .false. @@ -670,7 +670,7 @@ subroutine buoyancy_forcing_from_files(sfc_state, fluxes, day, dt, G, CS) temp(:,:), G%Domain, timelevel=time_lev) do j=js,je ; do i=is,ie fluxes%latent(i,j) = -hlv*temp(i,j) - fluxes%evap(i,j) = -temp(i,j) + fluxes%evap(i,j) = -US%kg_m3_to_R*US%m_to_Z*US%T_to_s * temp(i,j) fluxes%latent_evap_diag(i,j) = fluxes%latent(i,j) enddo ; enddo @@ -688,20 +688,20 @@ subroutine buoyancy_forcing_from_files(sfc_state, fluxes, day, dt, G, CS) enddo ; enddo call MOM_read_data(trim(CS%inputdir)//trim(CS%snow_file), "snow", & - fluxes%fprec(:,:), G%Domain, timelevel=time_lev) + fluxes%fprec(:,:), G%Domain, timelevel=time_lev, scale=US%kg_m3_to_R*US%m_to_Z*US%T_to_s) call MOM_read_data(trim(CS%inputdir)//trim(CS%precip_file), "precip", & - fluxes%lprec(:,:), G%Domain, timelevel=time_lev) + fluxes%lprec(:,:), G%Domain, timelevel=time_lev, scale=US%kg_m3_to_R*US%m_to_Z*US%T_to_s) do j=js,je ; do i=is,ie fluxes%lprec(i,j) = fluxes%lprec(i,j) - fluxes%fprec(i,j) enddo ; enddo call MOM_read_data(trim(CS%inputdir)//trim(CS%freshdischarge_file), "disch_w", & - temp(:,:), G%Domain, timelevel=time_lev_monthly) + temp(:,:), G%Domain, timelevel=time_lev_monthly, scale=US%kg_m3_to_R*US%m_to_Z*US%T_to_s) do j=js,je ; do i=is,ie fluxes%lrunoff(i,j) = temp(i,j)*US%m_to_L**2*G%IareaT(i,j) enddo ; enddo call MOM_read_data(trim(CS%inputdir)//trim(CS%freshdischarge_file), "disch_s", & - temp(:,:), G%Domain, timelevel=time_lev_monthly) + temp(:,:), G%Domain, timelevel=time_lev_monthly, scale=US%kg_m3_to_R*US%m_to_Z*US%T_to_s) do j=js,je ; do i=is,ie fluxes%frunoff(i,j) = temp(i,j)*US%m_to_L**2*G%IareaT(i,j) enddo ; enddo @@ -731,10 +731,11 @@ subroutine buoyancy_forcing_from_files(sfc_state, fluxes, day, dt, G, CS) fluxes%sw(i,j) = fluxes%sw(i,j) * G%mask2dT(i,j) fluxes%latent(i,j) = fluxes%latent(i,j) * G%mask2dT(i,j) - fluxes%heat_content_lrunoff(i,j) = fluxes%C_p*fluxes%lrunoff(i,j)*sfc_state%SST(i,j) + fluxes%heat_content_lrunoff(i,j) = fluxes%C_p*US%R_to_kg_m3*US%Z_to_m*US%s_to_T * & + fluxes%lrunoff(i,j)*sfc_state%SST(i,j) fluxes%latent_evap_diag(i,j) = fluxes%latent_evap_diag(i,j) * G%mask2dT(i,j) - fluxes%latent_fprec_diag(i,j) = -fluxes%fprec(i,j)*hlf - fluxes%latent_frunoff_diag(i,j) = -fluxes%frunoff(i,j)*hlf + fluxes%latent_fprec_diag(i,j) = -US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%fprec(i,j)*hlf + fluxes%latent_frunoff_diag(i,j) = -US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%frunoff(i,j)*hlf enddo ; enddo endif ! time_lev /= CS%buoy_last_lev_read @@ -745,7 +746,7 @@ subroutine buoyancy_forcing_from_files(sfc_state, fluxes, day, dt, G, CS) if (G%mask2dT(i,j) > 0) then fluxes%heat_restore(i,j) = G%mask2dT(i,j) * & ((CS%T_Restore(i,j) - sfc_state%SST(i,j)) * rhoXcp * CS%Flux_const) - fluxes%vprec(i,j) = - ((US%R_to_kg_m3*CS%Rho0)*CS%Flux_const) * & + fluxes%vprec(i,j) = - ((US%m_to_Z*US%T_to_s*CS%Rho0)*CS%Flux_const) * & (CS%S_Restore(i,j) - sfc_state%SSS(i,j)) / & (0.5*(sfc_state%SSS(i,j) + CS%S_Restore(i,j))) else @@ -878,7 +879,7 @@ subroutine buoyancy_forcing_linear(sfc_state, fluxes, day, dt, G, CS) if (G%mask2dT(i,j) > 0) then fluxes%heat_restore(i,j) = G%mask2dT(i,j) * & ((T_Restore - sfc_state%SST(i,j)) * (((US%R_to_kg_m3*CS%Rho0) * fluxes%C_p) * CS%Flux_const)) - fluxes%vprec(i,j) = - ((US%R_to_kg_m3*CS%Rho0)*CS%Flux_const) * & + fluxes%vprec(i,j) = - ((US%m_to_Z*US%T_to_s*CS%Rho0)*CS%Flux_const) * & (S_Restore - sfc_state%SSS(i,j)) / & (0.5*(sfc_state%SSS(i,j) + S_Restore)) else diff --git a/config_src/ice_solo_driver/user_surface_forcing.F90 b/config_src/ice_solo_driver/user_surface_forcing.F90 index 53ed835af9..28d60c895a 100644 --- a/config_src/ice_solo_driver/user_surface_forcing.F90 +++ b/config_src/ice_solo_driver/user_surface_forcing.F90 @@ -79,11 +79,11 @@ module user_surface_forcing ! state variables. logical :: restorebuoy ! If true, use restoring surface buoyancy forcing. real :: Rho0 ! The density used in the Boussinesq - ! approximation [kg m-3]. + ! approximation [R ~> kg m-3]. real :: G_Earth ! The gravitational acceleration [L2 Z-1 T-2 ~> m s-2]. real :: Flux_const ! The restoring rate at the surface [m s-1]. real :: gust_const ! A constant unresolved background gustiness - ! that contributes to ustar [Pa]. + ! that contributes to ustar [R Z L T-1 ~> Pa]. type(diag_ctrl), pointer :: diag ! A structure that is used to regulate the ! timing of diagnostic output. @@ -140,9 +140,9 @@ subroutine USER_wind_forcing(sfc_state, forces, day, G, US, CS) ! Set the surface friction velocity [Z s-1 ~> m s-1]. ustar is always positive. if (associated(forces%ustar)) then ; do j=js,je ; do i=is,ie ! This expression can be changed if desired, but need not be. - forces%ustar(i,j) = US%m_to_Z*US%T_to_s * G%mask2dT(i,j) * sqrt(CS%gust_const/CS%Rho0 + & - US%R_to_kg_m3*US%L_T_to_m_s**2*US%Z_to_L*sqrt(0.5*(forces%taux(I-1,j)**2 + forces%taux(I,j)**2) + & - 0.5*(forces%tauy(i,J-1)**2 + forces%tauy(i,J)**2))/CS%Rho0) + forces%ustar(i,j) = G%mask2dT(i,j) * sqrt(US%L_to_Z * (CS%gust_const/CS%Rho0 + & + sqrt(0.5*(forces%taux(I-1,j)**2 + forces%taux(I,j)**2) + & + 0.5*(forces%tauy(i,J-1)**2 + forces%tauy(i,J)**2))/CS%Rho0)) enddo ; enddo ; endif end subroutine USER_wind_forcing @@ -174,7 +174,7 @@ subroutine USER_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS) ! (fprec, lrunoff and frunoff) left as arrays full of zeros. ! Evap is usually negative and precip is usually positive. All heat fluxes ! are in W m-2 and positive for heat going into the ocean. All fresh water -! fluxes are in kg m-2 s-1 and positive for water moving into the ocean. +! fluxes are in [R Z T-1 ~> kg m-2 s-1] and positive for water moving into the ocean. real :: Temp_restore ! The temperature that is being restored toward [C]. real :: Salin_restore ! The salinity that is being restored toward [ppt] @@ -250,7 +250,7 @@ subroutine USER_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS) call MOM_error(FATAL, "User_buoyancy_surface_forcing: " // & "Temperature and salinity restoring used without modification." ) - rhoXcp = CS%Rho0 * fluxes%C_p + rhoXcp = US%R_to_kg_m3*CS%Rho0 * fluxes%C_p do j=js,je ; do i=is,ie ! Set Temp_restore and Salin_restore to the temperature (in degC) and ! salinity (in ppt or PSU) that are being restored toward. @@ -259,9 +259,8 @@ subroutine USER_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS) fluxes%heat_added(i,j) = (G%mask2dT(i,j) * (rhoXcp * CS%Flux_const)) * & (Temp_restore - sfc_state%SST(i,j)) - fluxes%vprec(i,j) = - (G%mask2dT(i,j) * (CS%Rho0*CS%Flux_const)) * & - ((Salin_restore - sfc_state%SSS(i,j)) / & - (0.5 * (Salin_restore + sfc_state%SSS(i,j)))) + fluxes%vprec(i,j) = - (G%mask2dT(i,j) * (US%m_to_Z*US%T_to_s*CS%Rho0*CS%Flux_const)) * & + ((Salin_restore - sfc_state%SSS(i,j)) / (0.5 * (Salin_restore + sfc_state%SSS(i,j)))) enddo ; enddo else ! When modifying the code, comment out this error message. It is here @@ -270,7 +269,7 @@ subroutine USER_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS) "Buoyancy restoring used without modification." ) ! The -1 is because density has the opposite sign to buoyancy. - buoy_rest_const = -1.0 * (CS%G_Earth * US%m_to_Z*US%T_to_s*CS%Flux_const) / CS%Rho0 + buoy_rest_const = -1.0 * (CS%G_Earth * US%m_to_Z*US%T_to_s*CS%Flux_const) / (US%R_to_kg_m3*CS%Rho0) do j=js,je ; do i=is,ie ! Set density_restore to an expression for the surface potential ! density [kg m-3] that is being restored toward. @@ -320,10 +319,10 @@ subroutine USER_surface_forcing_init(Time, G, US, param_file, diag, CS) "calculate accelerations and the mass for conservation "//& "properties, or with BOUSSINSEQ false to convert some "//& "parameters from vertical units of m to kg m-2.", & - units="kg m-3", default=1035.0) + units="kg m-3", default=1035.0, scale=US%R_to_kg_m3) call get_param(param_file, mdl, "GUST_CONST", CS%gust_const, & - "The background gustiness in the winds.", units="Pa", & - default=0.02) + "The background gustiness in the winds.", & + units="Pa", default=0.02, scale=US%kg_m3_to_R*US%m_s_to_L_T**2*US%L_to_Z) call get_param(param_file, mdl, "RESTOREBUOY", CS%restorebuoy, & "If true, the buoyancy fluxes drive the model back "//& diff --git a/config_src/mct_driver/mom_ocean_model_mct.F90 b/config_src/mct_driver/mom_ocean_model_mct.F90 index 4f1c7d963a..7ae09cf615 100644 --- a/config_src/mct_driver/mom_ocean_model_mct.F90 +++ b/config_src/mct_driver/mom_ocean_model_mct.F90 @@ -529,10 +529,10 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, & endif if (OS%icebergs_alter_ocean) then if (do_dyn) & - call iceberg_forces(OS%grid, OS%forces, OS%use_ice_shelf, & + call iceberg_forces(OS%grid, OS%US, OS%forces, OS%use_ice_shelf, & OS%sfc_state, dt_coupling, OS%marine_ice_CSp) if (do_thermo) & - call iceberg_fluxes(OS%grid, OS%fluxes, OS%use_ice_shelf, & + call iceberg_fluxes(OS%grid, OS%US, OS%fluxes, OS%use_ice_shelf, & OS%sfc_state, dt_coupling, OS%marine_ice_CSp) endif @@ -583,7 +583,7 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, & endif call set_derived_forcing_fields(OS%forces, OS%fluxes, OS%grid, OS%US, OS%GV%Rho0) - call set_net_mass_forcing(OS%fluxes, OS%forces, OS%grid) + call set_net_mass_forcing(OS%fluxes, OS%forces, OS%grid, OS%US) if (OS%use_waves) then call Update_Surface_Waves(OS%grid, OS%GV, OS%US, OS%time, ocean_coupling_time_step, OS%waves) @@ -677,7 +677,7 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, & if (OS%fluxes%fluxes_used) then call enable_averaging(OS%fluxes%dt_buoy_accum, OS%Time, OS%diag) call forcing_diagnostics(OS%fluxes, OS%sfc_state, OS%fluxes%dt_buoy_accum, & - OS%grid, OS%diag, OS%forcing_CSp%handles) + OS%grid, OS%US, OS%diag, OS%forcing_CSp%handles) call disable_averaging(OS%diag) endif diff --git a/config_src/mct_driver/mom_surface_forcing_mct.F90 b/config_src/mct_driver/mom_surface_forcing_mct.F90 index 7072c406e8..f9489c8a42 100644 --- a/config_src/mct_driver/mom_surface_forcing_mct.F90 +++ b/config_src/mct_driver/mom_surface_forcing_mct.F90 @@ -237,9 +237,11 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, & !! is present, or false (no restoring) otherwise. logical :: restore_sst !< local copy of the argument restore_temp, if it !! is present, or false (no restoring) otherwise. - real :: delta_sss !< temporary storage for sss diff from restoring value + real :: delta_sss !< temporary storage for sss diff from restoring value2 real :: delta_sst !< temporary storage for sst diff from restoring value + real :: kg_m2_s_conversion !< A combination of unit conversion factors for rescaling + !! mass fluxes [R Z s m2 kg-1 T-1 ~> 1]. real :: C_p !< heat capacity of seawater ( J/(K kg) ) real :: sign_for_net_FW_bug !< Should be +1. but an old bug can be recovered by using -1. @@ -253,6 +255,7 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, & IsdB = G%IsdB ; IedB = G%IedB ; JsdB = G%JsdB ; JedB = G%JedB isr = is-isd+1 ; ier = ie-isd+1 ; jsr = js-jsd+1 ; jer = je-jsd+1 + kg_m2_s_conversion = US%kg_m3_to_R*US%m_to_Z*US%T_to_s C_p = fluxes%C_p open_ocn_mask(:,:) = 1.0 pme_adj(:,:) = 0.0 @@ -375,19 +378,21 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, & delta_sss = sfc_state%SSS(i,j) - data_restore(i,j) delta_sss = sign(1.0,delta_sss)*min(abs(delta_sss),CS%max_delta_srestore) fluxes%vprec(i,j) = (CS%basin_mask(i,j)*open_ocn_mask(i,j)*CS%srestore_mask(i,j))* & - (US%R_to_kg_m3*CS%Rho0*CS%Flux_const) * & + (US%m_to_Z*US%T_to_s * CS%Rho0*CS%Flux_const) * & delta_sss / (0.5*(sfc_state%SSS(i,j) + data_restore(i,j))) endif enddo; enddo if (CS%adjust_net_srestore_to_zero) then if (CS%adjust_net_srestore_by_scaling) then - call adjust_area_mean_to_zero(fluxes%vprec, G, fluxes%vPrecGlobalScl) + call adjust_area_mean_to_zero(fluxes%vprec, G, fluxes%vPrecGlobalScl, & + unit_scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) fluxes%vPrecGlobalAdj = 0. else - work_sum(is:ie,js:je) = US%L_to_m**2*G%areaT(is:ie,js:je)*fluxes%vprec(is:ie,js:je) + work_sum(is:ie,js:je) = US%L_to_m**2*G%areaT(is:ie,js:je) * & + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%vprec(is:ie,js:je) fluxes%vPrecGlobalAdj = reproducing_sum(work_sum(:,:), isr, ier, jsr, jer) / CS%area_surf do j=js,je ; do i=is,ie - fluxes%vprec(i,j) = ( fluxes%vprec(i,j) - fluxes%vPrecGlobalAdj ) * G%mask2dT(i,j) + fluxes%vprec(i,j) = ( fluxes%vprec(i,j) - kg_m2_s_conversion*fluxes%vPrecGlobalAdj ) * G%mask2dT(i,j) enddo; enddo endif endif @@ -410,28 +415,28 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, & do j=js,je ; do i=is,ie ! liquid precipitation (rain) if (associated(IOB%lprec)) & - fluxes%lprec(i,j) = IOB%lprec(i-i0,j-j0) * G%mask2dT(i,j) + fluxes%lprec(i,j) = kg_m2_s_conversion * IOB%lprec(i-i0,j-j0) * G%mask2dT(i,j) ! frozen precipitation (snow) if (associated(IOB%fprec)) & - fluxes%fprec(i,j) = IOB%fprec(i-i0,j-j0) * G%mask2dT(i,j) + fluxes%fprec(i,j) = kg_m2_s_conversion * IOB%fprec(i-i0,j-j0) * G%mask2dT(i,j) ! evaporation if (associated(IOB%q_flux)) & - fluxes%evap(i,j) = IOB%q_flux(i-i0,j-j0) * G%mask2dT(i,j) + fluxes%evap(i,j) = kg_m2_s_conversion * IOB%q_flux(i-i0,j-j0) * G%mask2dT(i,j) ! liquid runoff flux if (associated(IOB%rofl_flux)) then - fluxes%lrunoff(i,j) = IOB%rofl_flux(i-i0,j-j0) * G%mask2dT(i,j) + fluxes%lrunoff(i,j) = kg_m2_s_conversion * IOB%rofl_flux(i-i0,j-j0) * G%mask2dT(i,j) else if (associated(IOB%runoff)) then - fluxes%lrunoff(i,j) = IOB%runoff(i-i0,j-j0) * G%mask2dT(i,j) + fluxes%lrunoff(i,j) = kg_m2_s_conversion * IOB%runoff(i-i0,j-j0) * G%mask2dT(i,j) end if ! ice runoff flux if (associated(IOB%rofi_flux)) then - fluxes%frunoff(i,j) = IOB%rofi_flux(i-i0,j-j0) * G%mask2dT(i,j) + fluxes%frunoff(i,j) = kg_m2_s_conversion * IOB%rofi_flux(i-i0,j-j0) * G%mask2dT(i,j) else if (associated(IOB%calving)) then - fluxes%frunoff(i,j) = IOB%calving(i-i0,j-j0) * G%mask2dT(i,j) + fluxes%frunoff(i,j) = kg_m2_s_conversion * IOB%calving(i-i0,j-j0) * G%mask2dT(i,j) end if if (associated(IOB%ustar_berg)) & @@ -467,7 +472,7 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, & ! water flux due to sea ice and snow melt [kg/m2/s] if (associated(IOB%seaice_melt)) & - fluxes%seaice_melt(i,j) = G%mask2dT(i,j) * IOB%seaice_melt(i-i0,j-j0) + fluxes%seaice_melt(i,j) = G%mask2dT(i,j) * kg_m2_s_conversion * IOB%seaice_melt(i-i0,j-j0) ! latent heat flux (W/m^2) fluxes%latent(i,j) = 0.0 @@ -533,7 +538,8 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, & sign_for_net_FW_bug = 1. if (CS%use_net_FW_adjustment_sign_bug) sign_for_net_FW_bug = -1. do j=js,je ; do i=is,ie - net_FW(i,j) = (((fluxes%lprec(i,j) + fluxes%fprec(i,j) + fluxes%seaice_melt(i,j)) + & + net_FW(i,j) = US%R_to_kg_m3*US%Z_to_m*US%s_to_T * & + (((fluxes%lprec(i,j) + fluxes%fprec(i,j) + fluxes%seaice_melt(i,j)) + & (fluxes%lrunoff(i,j) + fluxes%frunoff(i,j))) + & (fluxes%evap(i,j) + fluxes%vprec(i,j)) ) * US%L_to_m**2*G%areaT(i,j) @@ -543,13 +549,13 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, & if (CS%adjust_net_fresh_water_by_scaling) then call adjust_area_mean_to_zero(net_FW2, G, fluxes%netFWGlobalScl) do j=js,je ; do i=is,ie - fluxes%vprec(i,j) = fluxes%vprec(i,j) + & + fluxes%vprec(i,j) = fluxes%vprec(i,j) + kg_m2_s_conversion * & (net_FW2(i,j) - net_FW(i,j)/(US%L_to_m**2*G%areaT(i,j))) * G%mask2dT(i,j) enddo; enddo else fluxes%netFWGlobalAdj = reproducing_sum(net_FW(:,:), isr, ier, jsr, jer) / CS%area_surf do j=js,je ; do i=is,ie - fluxes%vprec(i,j) = ( fluxes%vprec(i,j) - fluxes%netFWGlobalAdj ) * G%mask2dT(i,j) + fluxes%vprec(i,j) = ( fluxes%vprec(i,j) - kg_m2_s_conversion * fluxes%netFWGlobalAdj ) * G%mask2dT(i,j) enddo; enddo endif endif @@ -560,7 +566,7 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, & if (CS%allow_flux_adjustments) then ! Apply adjustments to fluxes - call apply_flux_adjustments(G, CS, Time, fluxes) + call apply_flux_adjustments(G, US, CS, Time, fluxes) endif ! Allow for user-written code to alter fluxes after all the above @@ -865,8 +871,9 @@ end subroutine convert_IOB_to_forces !! - hflx_adj (Heat flux into the ocean, in W m-2) !! - sflx_adj (Salt flux into the ocean, in kg salt m-2 s-1) !! - prcme_adj (Fresh water flux into the ocean, in kg m-2 s-1) -subroutine apply_flux_adjustments(G, CS, Time, fluxes) +subroutine apply_flux_adjustments(G, US, CS, Time, fluxes) type(ocean_grid_type), intent(inout) :: G !< Ocean grid structure + type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(surface_forcing_CS), pointer :: CS !< Surface forcing control structure type(time_type), intent(in) :: Time !< Model time structure type(forcing), intent(inout) :: fluxes !< Surface fluxes structure @@ -899,7 +906,7 @@ subroutine apply_flux_adjustments(G, CS, Time, fluxes) call data_override('OCN', 'prcme_adj', temp_at_h(isc:iec,jsc:jec), Time, override=overrode_h) if (overrode_h) then ; do j=jsc,jec ; do i=isc,iec - fluxes%vprec(i,j) = fluxes%vprec(i,j) + temp_at_h(i,j)* G%mask2dT(i,j) + fluxes%vprec(i,j) = fluxes%vprec(i,j) + US%kg_m3_to_R*US%m_to_Z*US%T_to_s * temp_at_h(i,j)* G%mask2dT(i,j) enddo ; enddo ; endif ! Not needed? ! if (overrode_h) call pass_var(fluxes%vprec, G%Domain) diff --git a/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 b/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 index e04064f672..726ad93ec0 100644 --- a/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 +++ b/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 @@ -526,7 +526,7 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, & call iceberg_forces(OS%grid, OS%forces, OS%use_ice_shelf, & OS%sfc_state, dt_coupling, OS%marine_ice_CSp) if (do_thermo) & - call iceberg_fluxes(OS%grid, OS%fluxes, OS%use_ice_shelf, & + call iceberg_fluxes(OS%grid, OS%US, OS%fluxes, OS%use_ice_shelf, & OS%sfc_state, dt_coupling, OS%marine_ice_CSp) endif @@ -557,7 +557,7 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, & call iceberg_forces(OS%grid, OS%forces, OS%use_ice_shelf, & OS%sfc_state, dt_coupling, OS%marine_ice_CSp) if (do_thermo) & - call iceberg_fluxes(OS%grid, OS%flux_tmp, OS%use_ice_shelf, & + call iceberg_fluxes(OS%grid, OS%US, OS%flux_tmp, OS%use_ice_shelf, & OS%sfc_state, dt_coupling, OS%marine_ice_CSp) endif @@ -571,7 +571,7 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, & #endif endif call set_derived_forcing_fields(OS%forces, OS%fluxes, OS%grid, OS%US, OS%GV%Rho0) - call set_net_mass_forcing(OS%fluxes, OS%forces, OS%grid) + call set_net_mass_forcing(OS%fluxes, OS%forces, OS%grid, OS%US) if (OS%use_waves) then call Update_Surface_Waves(OS%grid, OS%GV, OS%US, OS%time, ocean_coupling_time_step, OS%waves) @@ -664,7 +664,7 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, & if (OS%fluxes%fluxes_used) then call enable_averaging(OS%fluxes%dt_buoy_accum, OS%Time, OS%diag) call forcing_diagnostics(OS%fluxes, OS%sfc_state, OS%fluxes%dt_buoy_accum, & - OS%grid, OS%diag, OS%forcing_CSp%handles) + OS%grid, US%US, OS%diag, OS%forcing_CSp%handles) call disable_averaging(OS%diag) endif diff --git a/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90 b/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90 index 7e56780a36..f81ea561db 100644 --- a/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90 +++ b/config_src/nuopc_driver/mom_surface_forcing_nuopc.F90 @@ -244,6 +244,8 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, & !! is present, or false (no restoring) otherwise. real :: delta_sss !< temporary storage for sss diff from restoring value real :: delta_sst !< temporary storage for sst diff from restoring value + real :: kg_m2_s_conversion !< A combination of unit conversion factors for rescaling + !! mass fluxes [R Z s m2 kg-1 T-1 ~> 1]. real :: C_p !< heat capacity of seawater ( J/(K kg) ) real :: sign_for_net_FW_bug !< Should be +1. but an old bug can be recovered by using -1. @@ -258,6 +260,7 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, & IsdB = G%IsdB ; IedB = G%IedB ; JsdB = G%JsdB ; JedB = G%JedB isr = is-isd+1 ; ier = ie-isd+1 ; jsr = js-jsd+1 ; jer = je-jsd+1 + kg_m2_s_conversion = US%kg_m3_to_R*US%m_to_Z*US%T_to_s C_p = fluxes%C_p open_ocn_mask(:,:) = 1.0 pme_adj(:,:) = 0.0 @@ -387,13 +390,15 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, & enddo ; enddo if (CS%adjust_net_srestore_to_zero) then if (CS%adjust_net_srestore_by_scaling) then - call adjust_area_mean_to_zero(fluxes%vprec, G, fluxes%vPrecGlobalScl) + call adjust_area_mean_to_zero(fluxes%vprec, G, fluxes%vPrecGlobalScl, & + unit_scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) fluxes%vPrecGlobalAdj = 0. else - work_sum(is:ie,js:je) = US%L_to_m**2*G%areaT(is:ie,js:je)*fluxes%vprec(is:ie,js:je) + work_sum(is:ie,js:je) = US%L_to_m**2*G%areaT(is:ie,js:je) * & + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%vprec(is:ie,js:je) fluxes%vPrecGlobalAdj = reproducing_sum(work_sum(:,:), isr, ier, jsr, jer) / CS%area_surf do j=js,je ; do i=is,ie - fluxes%vprec(i,j) = ( fluxes%vprec(i,j) - fluxes%vPrecGlobalAdj ) * G%mask2dT(i,j) + fluxes%vprec(i,j) = ( fluxes%vprec(i,j) - kg_m2_s_conversion * fluxes%vPrecGlobalAdj ) * G%mask2dT(i,j) enddo ; enddo endif endif @@ -416,26 +421,26 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, & do j=js,je ; do i=is,ie if (associated(IOB%lprec)) & - fluxes%lprec(i,j) = IOB%lprec(i-i0,j-j0) * G%mask2dT(i,j) + fluxes%lprec(i,j) = kg_m2_s_conversion * IOB%lprec(i-i0,j-j0) * G%mask2dT(i,j) if (associated(IOB%fprec)) & - fluxes%fprec(i,j) = IOB%fprec(i-i0,j-j0) * G%mask2dT(i,j) + fluxes%fprec(i,j) = kg_m2_s_conversion * IOB%fprec(i-i0,j-j0) * G%mask2dT(i,j) if (associated(IOB%q_flux)) & - fluxes%evap(i,j) = IOB%q_flux(i-i0,j-j0) * G%mask2dT(i,j) + fluxes%evap(i,j) = kg_m2_s_conversion * IOB%q_flux(i-i0,j-j0) * G%mask2dT(i,j) ! liquid runoff flux if (associated(IOB%rofl_flux)) then - fluxes%lrunoff(i,j) = IOB%rofl_flux(i-i0,j-j0) * G%mask2dT(i,j) + fluxes%lrunoff(i,j) = kg_m2_s_conversion * IOB%rofl_flux(i-i0,j-j0) * G%mask2dT(i,j) else if (associated(IOB%runoff)) then - fluxes%lrunoff(i,j) = IOB%runoff(i-i0,j-j0) * G%mask2dT(i,j) + fluxes%lrunoff(i,j) = kg_m2_s_conversion * IOB%runoff(i-i0,j-j0) * G%mask2dT(i,j) endif ! ice runoff flux if (associated(IOB%rofi_flux)) then - fluxes%frunoff(i,j) = IOB%rofi_flux(i-i0,j-j0) * G%mask2dT(i,j) + fluxes%frunoff(i,j) = kg_m2_s_conversion * IOB%rofi_flux(i-i0,j-j0) * G%mask2dT(i,j) elseif (associated(IOB%calving)) then - fluxes%frunoff(i,j) = IOB%calving(i-i0,j-j0) * G%mask2dT(i,j) + fluxes%frunoff(i,j) = kg_m2_s_conversion * IOB%calving(i-i0,j-j0) * G%mask2dT(i,j) endif if (associated(IOB%ustar_berg)) & @@ -465,7 +470,7 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, & ! water flux due to sea ice and snow melt [kg/m2/s] if (associated(IOB%seaice_melt)) & - fluxes%seaice_melt(i,j) = G%mask2dT(i,j) * IOB%seaice_melt(i-i0,j-j0) + fluxes%seaice_melt(i,j) = kg_m2_s_conversion * G%mask2dT(i,j) * IOB%seaice_melt(i-i0,j-j0) fluxes%latent(i,j) = 0.0 if (associated(IOB%fprec)) then @@ -527,7 +532,8 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, & sign_for_net_FW_bug = 1. if (CS%use_net_FW_adjustment_sign_bug) sign_for_net_FW_bug = -1. do j=js,je ; do i=is,ie - net_FW(i,j) = (((fluxes%lprec(i,j) + fluxes%fprec(i,j) + fluxes%seaice_melt(i,j)) + & + net_FW(i,j) = US%R_to_kg_m3*US%Z_to_m*US%s_to_T * & + (((fluxes%lprec(i,j) + fluxes%fprec(i,j) + fluxes%seaice_melt(i,j)) + & (fluxes%lrunoff(i,j) + fluxes%frunoff(i,j))) + & (fluxes%evap(i,j) + fluxes%vprec(i,j)) ) * US%L_to_m**2*G%areaT(i,j) net_FW2(i,j) = net_FW(i,j) / (US%L_to_m**2*G%areaT(i,j)) @@ -536,13 +542,13 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, & if (CS%adjust_net_fresh_water_by_scaling) then call adjust_area_mean_to_zero(net_FW2, G, fluxes%netFWGlobalScl) do j=js,je ; do i=is,ie - fluxes%vprec(i,j) = fluxes%vprec(i,j) + & + fluxes%vprec(i,j) = fluxes%vprec(i,j) + US%kg_m3_to_R*US%m_to_Z*US%T_to_s * & (net_FW2(i,j) - net_FW(i,j)/(US%L_to_m**2*G%areaT(i,j))) * G%mask2dT(i,j) enddo ; enddo else fluxes%netFWGlobalAdj = reproducing_sum(net_FW(:,:), isr, ier, jsr, jer) / CS%area_surf do j=js,je ; do i=is,ie - fluxes%vprec(i,j) = ( fluxes%vprec(i,j) - fluxes%netFWGlobalAdj ) * G%mask2dT(i,j) + fluxes%vprec(i,j) = ( fluxes%vprec(i,j) - kg_m2_s_conversion * fluxes%netFWGlobalAdj ) * G%mask2dT(i,j) enddo ; enddo endif @@ -554,7 +560,7 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, G, US, CS, & if (CS%allow_flux_adjustments) then ! Apply adjustments to fluxes - call apply_flux_adjustments(G, CS, Time, fluxes) + call apply_flux_adjustments(G, US, CS, Time, fluxes) endif ! Allow for user-written code to alter fluxes after all the above @@ -862,8 +868,9 @@ end subroutine convert_IOB_to_forces !! - hflx_adj (Heat flux into the ocean, in W m-2) !! - sflx_adj (Salt flux into the ocean, in kg salt m-2 s-1) !! - prcme_adj (Fresh water flux into the ocean, in kg m-2 s-1) -subroutine apply_flux_adjustments(G, CS, Time, fluxes) +subroutine apply_flux_adjustments(G, US, CS, Time, fluxes) type(ocean_grid_type), intent(inout) :: G !< Ocean grid structure + type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(surface_forcing_CS), pointer :: CS !< Surface forcing control structure type(time_type), intent(in) :: Time !< Model time structure type(forcing), intent(inout) :: fluxes !< Surface fluxes structure @@ -896,7 +903,7 @@ subroutine apply_flux_adjustments(G, CS, Time, fluxes) call data_override('OCN', 'prcme_adj', temp_at_h(isc:iec,jsc:jec), Time, override=overrode_h) if (overrode_h) then ; do j=jsc,jec ; do i=isc,iec - fluxes%vprec(i,j) = fluxes%vprec(i,j) + temp_at_h(i,j)* G%mask2dT(i,j) + fluxes%vprec(i,j) = fluxes%vprec(i,j) + US%kg_m3_to_R*US%m_to_Z*US%T_to_s * temp_at_h(i,j)* G%mask2dT(i,j) enddo ; enddo ; endif ! Not needed? ! if (overrode_h) call pass_var(fluxes%vprec, G%Domain) end subroutine apply_flux_adjustments diff --git a/config_src/solo_driver/MESO_surface_forcing.F90 b/config_src/solo_driver/MESO_surface_forcing.F90 index ee3cd36b41..f828513dae 100644 --- a/config_src/solo_driver/MESO_surface_forcing.F90 +++ b/config_src/solo_driver/MESO_surface_forcing.F90 @@ -142,7 +142,7 @@ subroutine MESO_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS) ! Fluxes of fresh water through the surface are in units of [kg m-2 s-1] ! and are positive downward - i.e. evaporation should be negative. fluxes%evap(i,j) = -0.0 * G%mask2dT(i,j) - fluxes%lprec(i,j) = CS%PmE(i,j) * CS%Rho0 * G%mask2dT(i,j) + fluxes%lprec(i,j) = US%kg_m3_to_R*US%m_to_Z*US%T_to_s * CS%PmE(i,j) * CS%Rho0 * G%mask2dT(i,j) ! vprec will be set later, if it is needed for salinity restoring. fluxes%vprec(i,j) = 0.0 @@ -176,7 +176,7 @@ subroutine MESO_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS) if (G%mask2dT(i,j) > 0) then fluxes%heat_added(i,j) = G%mask2dT(i,j) * & ((CS%T_Restore(i,j) - sfc_state%SST(i,j)) * rhoXcp * CS%Flux_const) - fluxes%vprec(i,j) = - (CS%Rho0*CS%Flux_const) * & + fluxes%vprec(i,j) = - (US%kg_m3_to_R*US%m_to_Z*US%T_to_s*CS%Rho0*CS%Flux_const) * & (CS%S_Restore(i,j) - sfc_state%SSS(i,j)) / & (0.5*(sfc_state%SSS(i,j) + CS%S_Restore(i,j))) else diff --git a/config_src/solo_driver/MOM_driver.F90 b/config_src/solo_driver/MOM_driver.F90 index b057e06f9e..a6d6597c0e 100644 --- a/config_src/solo_driver/MOM_driver.F90 +++ b/config_src/solo_driver/MOM_driver.F90 @@ -580,7 +580,7 @@ program MOM_main if (.not. offline_tracer_mode) then if (fluxes%fluxes_used) then call enable_averaging(fluxes%dt_buoy_accum, Time, diag) - call forcing_diagnostics(fluxes, sfc_state, fluxes%dt_buoy_accum, grid, & + call forcing_diagnostics(fluxes, sfc_state, fluxes%dt_buoy_accum, grid, US, & diag, surface_forcing_CSp%handles) call disable_averaging(diag) else diff --git a/config_src/solo_driver/MOM_surface_forcing.F90 b/config_src/solo_driver/MOM_surface_forcing.F90 index 101956d283..7224d68d48 100644 --- a/config_src/solo_driver/MOM_surface_forcing.F90 +++ b/config_src/solo_driver/MOM_surface_forcing.F90 @@ -315,13 +315,13 @@ subroutine set_forcing(sfc_state, forces, fluxes, day_start, day_interval, G, US elseif (trim(CS%buoy_config) == "Neverland") then call Neverland_buoyancy_forcing(sfc_state, fluxes, day_center, dt, G, US, CS%Neverland_forcing_CSp) elseif (trim(CS%buoy_config) == "SCM_CVmix_tests") then - call SCM_CVmix_tests_buoyancy_forcing(sfc_state, fluxes, day_center, G, CS%SCM_CVmix_tests_CSp) + call SCM_CVmix_tests_buoyancy_forcing(sfc_state, fluxes, day_center, G, US, CS%SCM_CVmix_tests_CSp) elseif (trim(CS%buoy_config) == "USER") then call USER_buoyancy_forcing(sfc_state, fluxes, day_center, dt, G, US, CS%user_forcing_CSp) elseif (trim(CS%buoy_config) == "BFB") then call BFB_buoyancy_forcing(sfc_state, fluxes, day_center, dt, G, US, CS%BFB_forcing_CSp) elseif (trim(CS%buoy_config) == "dumbbell") then - call dumbbell_buoyancy_forcing(sfc_state, fluxes, day_center, dt, G, CS%dumbbell_forcing_CSp) + call dumbbell_buoyancy_forcing(sfc_state, fluxes, day_center, dt, G, US, CS%dumbbell_forcing_CSp) elseif (trim(CS%buoy_config) == "NONE") then call MOM_mesg("MOM_surface_forcing: buoyancy forcing has been set to omitted.") elseif (CS%variable_buoyforce .and. .not.CS%first_call_set_forcing) then @@ -348,7 +348,7 @@ subroutine set_forcing(sfc_state, forces, fluxes, day_start, day_interval, G, US if ((CS%variable_buoyforce .or. CS%first_call_set_forcing) .and. & (.not.CS%adiabatic)) then - call set_net_mass_forcing(fluxes, forces, G) + call set_net_mass_forcing(fluxes, forces, G, US) endif CS%first_call_set_forcing = .false. @@ -842,12 +842,12 @@ subroutine buoyancy_forcing_from_files(sfc_state, fluxes, day, dt, G, US, CS) G%Domain, timelevel=time_lev) do j=js,je ; do i=is,ie fluxes%latent(i,j) = -CS%latent_heat_vapor*temp(i,j) - fluxes%evap(i,j) = -temp(i,j) + fluxes%evap(i,j) = -US%kg_m3_to_R*US%m_to_Z*US%T_to_s*temp(i,j) fluxes%latent_evap_diag(i,j) = fluxes%latent(i,j) enddo ; enddo else call MOM_read_data(CS%evaporation_file, CS%evap_var, fluxes%evap(:,:), & - G%Domain, timelevel=time_lev) + G%Domain, timelevel=time_lev, scale=US%kg_m3_to_R*US%m_to_Z*US%T_to_s) endif CS%evap_last_lev = time_lev @@ -902,9 +902,9 @@ subroutine buoyancy_forcing_from_files(sfc_state, fluxes, day, dt, G, US, CS) case default ; time_lev = 1 end select call MOM_read_data(CS%snow_file, CS%snow_var, & - fluxes%fprec(:,:), G%Domain, timelevel=time_lev) + fluxes%fprec(:,:), G%Domain, timelevel=time_lev, scale=US%kg_m3_to_R*US%m_to_Z*US%T_to_s) call MOM_read_data(CS%rain_file, CS%rain_var, & - fluxes%lprec(:,:), G%Domain, timelevel=time_lev) + fluxes%lprec(:,:), G%Domain, timelevel=time_lev, scale=US%kg_m3_to_R*US%m_to_Z*US%T_to_s) if (CS%archaic_OMIP_file) then do j=js,je ; do i=is,ie fluxes%lprec(i,j) = fluxes%lprec(i,j) - fluxes%fprec(i,j) @@ -919,20 +919,20 @@ subroutine buoyancy_forcing_from_files(sfc_state, fluxes, day, dt, G, US, CS) end select if (CS%archaic_OMIP_file) then call MOM_read_data(CS%runoff_file, CS%lrunoff_var, temp(:,:), & - G%Domain, timelevel=time_lev) + G%Domain, timelevel=time_lev, scale=US%kg_m3_to_R*US%m_to_Z*US%T_to_s) do j=js,je ; do i=is,ie fluxes%lrunoff(i,j) = temp(i,j)*US%m_to_L**2*G%IareaT(i,j) enddo ; enddo call MOM_read_data(CS%runoff_file, CS%frunoff_var, temp(:,:), & - G%Domain, timelevel=time_lev) + G%Domain, timelevel=time_lev, scale=US%kg_m3_to_R*US%m_to_Z*US%T_to_s) do j=js,je ; do i=is,ie fluxes%frunoff(i,j) = temp(i,j)*US%m_to_L**2*G%IareaT(i,j) enddo ; enddo else call MOM_read_data(CS%runoff_file, CS%lrunoff_var, fluxes%lrunoff(:,:), & - G%Domain, timelevel=time_lev) + G%Domain, timelevel=time_lev, scale=US%kg_m3_to_R*US%m_to_Z*US%T_to_s) call MOM_read_data(CS%runoff_file, CS%frunoff_var, fluxes%frunoff(:,:), & - G%Domain, timelevel=time_lev) + G%Domain, timelevel=time_lev, scale=US%kg_m3_to_R*US%m_to_Z*US%T_to_s) endif CS%runoff_last_lev = time_lev @@ -976,8 +976,8 @@ subroutine buoyancy_forcing_from_files(sfc_state, fluxes, day, dt, G, US, CS) fluxes%latent(i,j) = fluxes%latent(i,j) * G%mask2dT(i,j) fluxes%latent_evap_diag(i,j) = fluxes%latent_evap_diag(i,j) * G%mask2dT(i,j) - fluxes%latent_fprec_diag(i,j) = -fluxes%fprec(i,j)*CS%latent_heat_fusion - fluxes%latent_frunoff_diag(i,j) = -fluxes%frunoff(i,j)*CS%latent_heat_fusion + fluxes%latent_fprec_diag(i,j) = -US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%fprec(i,j)*CS%latent_heat_fusion + fluxes%latent_frunoff_diag(i,j) = -US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%frunoff(i,j)*CS%latent_heat_fusion enddo ; enddo endif ! time_lev /= CS%buoy_last_lev_read @@ -991,12 +991,12 @@ subroutine buoyancy_forcing_from_files(sfc_state, fluxes, day, dt, G, US, CS) if (G%mask2dT(i,j) > 0) then fluxes%heat_added(i,j) = G%mask2dT(i,j) * & ((CS%T_Restore(i,j) - sfc_state%SST(i,j)) * rhoXcp * CS%Flux_const_T) - fluxes%vprec(i,j) = - (US%R_to_kg_m3*CS%Rho0*CS%Flux_const_S) * & + fluxes%vprec(i,j) = - (US%m_to_Z*US%T_to_s*CS%Rho0*CS%Flux_const_S) * & (CS%S_Restore(i,j) - sfc_state%SSS(i,j)) / & (0.5*(sfc_state%SSS(i,j) + CS%S_Restore(i,j))) else fluxes%heat_added(i,j) = 0.0 - fluxes%vprec(i,j) = 0.0 + fluxes%vprec(i,j) = 0.0 endif enddo ; enddo else @@ -1089,10 +1089,12 @@ subroutine buoyancy_forcing_from_data_override(sfc_state, fluxes, day, dt, G, US ! note the sign convention do j=js,je ; do i=is,ie - fluxes%evap(i,j) = -fluxes%evap(i,j) ! Normal convention is positive into the ocean - ! but evap is normally a positive quantity in the files - fluxes%latent(i,j) = CS%latent_heat_vapor*fluxes%evap(i,j) - fluxes%latent_evap_diag(i,j) = fluxes%latent(i,j) + ! This is dangerous because it is not clear whether the data files have been read! + fluxes%evap(i,j) = -fluxes%evap(i,j) ! Normal convention is positive into the ocean + ! but evap is normally a positive quantity in the files + fluxes%latent(i,j) = CS%latent_heat_vapor*fluxes%evap(i,j) + fluxes%latent_evap_diag(i,j) = fluxes%latent(i,j) + fluxes%evap(i,j) = US%kg_m3_to_R*US%m_to_Z*US%T_to_s*fluxes%evap(i,j) enddo ; enddo call data_override('OCN', 'sens', fluxes%sens(:,:), day, & @@ -1108,16 +1110,23 @@ subroutine buoyancy_forcing_from_data_override(sfc_state, fluxes, day, dt, G, US is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in) call data_override('OCN', 'snow', fluxes%fprec(:,:), day, & - is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in) + is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in) ! scale=US%kg_m3_to_R*US%m_to_Z*US%T_to_s call data_override('OCN', 'rain', fluxes%lprec(:,:), day, & - is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in) + is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in) ! scale=US%kg_m3_to_R*US%m_to_Z*US%T_to_s call data_override('OCN', 'runoff', fluxes%lrunoff(:,:), day, & - is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in) + is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in) ! scale=US%kg_m3_to_R*US%m_to_Z*US%T_to_s call data_override('OCN', 'calving', fluxes%frunoff(:,:), day, & - is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in) + is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in) ! scale=US%kg_m3_to_R*US%m_to_Z*US%T_to_s + + if (US%kg_m3_to_R*US%m_to_Z*US%T_to_s /= 1.0) then ; do j=js,je ; do i=is,ie + fluxes%lprec(i,j) = fluxes%lprec(i,j) * US%kg_m3_to_R*US%m_to_Z*US%T_to_s + fluxes%fprec(i,j) = fluxes%fprec(i,j) * US%kg_m3_to_R*US%m_to_Z*US%T_to_s + fluxes%lrunoff(i,j) = fluxes%lrunoff(i,j) * US%kg_m3_to_R*US%m_to_Z*US%T_to_s + fluxes%frunoff(i,j) = fluxes%frunoff(i,j) * US%kg_m3_to_R*US%m_to_Z*US%T_to_s + enddo ; enddo ; endif ! Read the SST and SSS fields for damping. if (CS%restorebuoy) then !#CTRL# .or. associated(CS%ctrl_forcing_CSp)) then @@ -1136,7 +1145,7 @@ subroutine buoyancy_forcing_from_data_override(sfc_state, fluxes, day, dt, G, US if (G%mask2dT(i,j) > 0) then fluxes%heat_added(i,j) = G%mask2dT(i,j) * & ((CS%T_Restore(i,j) - sfc_state%SST(i,j)) * rhoXcp * CS%Flux_const_T) - fluxes%vprec(i,j) = - (Rho0_mks*CS%Flux_const_S) * & + fluxes%vprec(i,j) = - (CS%Rho0*US%m_to_Z*US%T_to_s*CS%Flux_const_S) * & (CS%S_Restore(i,j) - sfc_state%SSS(i,j)) / & (0.5*(sfc_state%SSS(i,j) + CS%S_Restore(i,j))) else @@ -1180,8 +1189,8 @@ subroutine buoyancy_forcing_from_data_override(sfc_state, fluxes, day, dt, G, US fluxes%sw(i,j) = fluxes%sw(i,j) * G%mask2dT(i,j) fluxes%latent_evap_diag(i,j) = fluxes%latent_evap_diag(i,j) * G%mask2dT(i,j) - fluxes%latent_fprec_diag(i,j) = -fluxes%fprec(i,j)*CS%latent_heat_fusion - fluxes%latent_frunoff_diag(i,j) = -fluxes%frunoff(i,j)*CS%latent_heat_fusion + fluxes%latent_fprec_diag(i,j) = -US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%fprec(i,j)*CS%latent_heat_fusion + fluxes%latent_frunoff_diag(i,j) = -US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%frunoff(i,j)*CS%latent_heat_fusion enddo ; enddo @@ -1336,7 +1345,7 @@ subroutine buoyancy_forcing_linear(sfc_state, fluxes, day, dt, G, US, CS) if (G%mask2dT(i,j) > 0) then fluxes%heat_added(i,j) = G%mask2dT(i,j) * & ((T_Restore - sfc_state%SST(i,j)) * ((Rho0_mks * fluxes%C_p) * CS%Flux_const)) - fluxes%vprec(i,j) = - (Rho0_mks*CS%Flux_const) * & + fluxes%vprec(i,j) = - (US%m_to_Z*US%T_to_s*CS%Rho0*CS%Flux_const) * & (S_Restore - sfc_state%SSS(i,j)) / & (0.5*(sfc_state%SSS(i,j) + S_Restore)) else diff --git a/config_src/solo_driver/user_surface_forcing.F90 b/config_src/solo_driver/user_surface_forcing.F90 index 1afe999e51..1831503f1f 100644 --- a/config_src/solo_driver/user_surface_forcing.F90 +++ b/config_src/solo_driver/user_surface_forcing.F90 @@ -123,7 +123,7 @@ subroutine USER_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS) ! (fprec, lrunoff and frunoff) left as arrays full of zeros. ! Evap is usually negative and precip is usually positive. All heat fluxes ! are in W m-2 and positive for heat going into the ocean. All fresh water -! fluxes are in kg m-2 s-1 and positive for water moving into the ocean. +! fluxes are in [R Z T-1 ~> kg m-2 s-1] and positive for water moving into the ocean. ! Local variables real :: Temp_restore ! The temperature that is being restored toward [degC]. @@ -172,7 +172,7 @@ subroutine USER_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS) ! Set whichever fluxes are to be used here. Any fluxes that ! are always zero do not need to be changed here. do j=js,je ; do i=is,ie - ! Fluxes of fresh water through the surface are in units of [kg m-2 s-1] + ! Fluxes of fresh water through the surface are in units of [R Z T-1 ~> kg m-2 s-1] ! and are positive downward - i.e. evaporation should be negative. fluxes%evap(i,j) = -0.0 * G%mask2dT(i,j) fluxes%lprec(i,j) = 0.0 * G%mask2dT(i,j) @@ -211,9 +211,8 @@ subroutine USER_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS) fluxes%heat_added(i,j) = (G%mask2dT(i,j) * (rhoXcp * CS%Flux_const)) * & (Temp_restore - sfc_state%SST(i,j)) - fluxes%vprec(i,j) = - (G%mask2dT(i,j) * (Rho0_mks*CS%Flux_const)) * & - ((Salin_restore - sfc_state%SSS(i,j)) / & - (0.5 * (Salin_restore + sfc_state%SSS(i,j)))) + fluxes%vprec(i,j) = - (G%mask2dT(i,j) * (US%m_to_Z*US%T_to_s*CS%Rho0*CS%Flux_const)) * & + ((Salin_restore - sfc_state%SSS(i,j)) / (0.5 * (Salin_restore + sfc_state%SSS(i,j)))) enddo ; enddo else ! When modifying the code, comment out this error message. It is here diff --git a/src/core/MOM.F90 b/src/core/MOM.F90 index f4ef5a1376..b490311cf2 100644 --- a/src/core/MOM.F90 +++ b/src/core/MOM.F90 @@ -852,7 +852,7 @@ subroutine step_MOM(forces, fluxes, sfc_state, Time_start, time_interval, CS, & ! Accumulate the surface fluxes for assessing conservation if (do_thermo .and. fluxes%fluxes_used) & call accumulate_net_input(fluxes, sfc_state, CS%tv, fluxes%dt_buoy_accum, & - G, CS%sum_output_CSp) + G, US, CS%sum_output_CSp) if (MOM_state_is_synchronized(CS)) & call write_energy(CS%u, CS%v, CS%h, CS%tv, Time_local, CS%nstep_tot, & diff --git a/src/core/MOM_forcing_type.F90 b/src/core/MOM_forcing_type.F90 index 2b064a2834..a5e56b9ad1 100644 --- a/src/core/MOM_forcing_type.F90 +++ b/src/core/MOM_forcing_type.F90 @@ -82,13 +82,13 @@ module MOM_forcing_type ! water mass fluxes into the ocean [kg m-2 s-1]; these fluxes impact the ocean mass real, pointer, dimension(:,:) :: & - evap => NULL(), & !< (-1)*fresh water flux evaporated out of the ocean [kg m-2 s-1] - lprec => NULL(), & !< precipitating liquid water into the ocean [kg m-2 s-1] - fprec => NULL(), & !< precipitating frozen water into the ocean [kg m-2 s-1] - vprec => NULL(), & !< virtual liquid precip associated w/ SSS restoring [kg m-2 s-1] - lrunoff => NULL(), & !< liquid river runoff entering ocean [kg m-2 s-1] - frunoff => NULL(), & !< frozen river runoff (calving) entering ocean [kg m-2 s-1] - seaice_melt => NULL(), & !< snow/seaice melt (positive) or formation (negative) [kg m-2 s-1] + evap => NULL(), & !< (-1)*fresh water flux evaporated out of the ocean [R Z T-1 ~> kg m-2 s-1] + lprec => NULL(), & !< precipitating liquid water into the ocean [R Z T-1 ~> kg m-2 s-1] + fprec => NULL(), & !< precipitating frozen water into the ocean [R Z T-1 ~> kg m-2 s-1] + vprec => NULL(), & !< virtual liquid precip associated w/ SSS restoring [R Z T-1 ~> kg m-2 s-1] + lrunoff => NULL(), & !< liquid river runoff entering ocean [R Z T-1 ~> kg m-2 s-1] + frunoff => NULL(), & !< frozen river runoff (calving) entering ocean [R Z T-1 ~> kg m-2 s-1] + seaice_melt => NULL(), & !< snow/seaice melt (positive) or formation (negative) [R Z T-1 ~> kg m-2 s-1] netMassIn => NULL(), & !< Sum of water mass flux out of the ocean [kg m-2 s-1] netMassOut => NULL(), & !< Net water mass flux into of the ocean [kg m-2 s-1] netSalt => NULL() !< Net salt entering the ocean [kgSalt m-2 s-1] @@ -97,7 +97,7 @@ module MOM_forcing_type real, pointer, dimension(:,:) :: & heat_content_cond => NULL(), & !< heat content associated with condensating water [W m-2] heat_content_lprec => NULL(), & !< heat content associated with liquid >0 precip [W m-2] (diagnostic) - heat_content_icemelt => NULL(), & !< heat content associated with snow/seaice melt/formation [W/m^2] + heat_content_icemelt => NULL(), & !< heat content associated with snow/seaice melt/formation [W m-2] heat_content_fprec => NULL(), & !< heat content associated with frozen precip [W m-2] heat_content_vprec => NULL(), & !< heat content associated with virtual >0 precip [W m-2] heat_content_lrunoff => NULL(), & !< heat content associated with liquid runoff [W m-2] @@ -509,16 +509,18 @@ subroutine extractFluxes1d(G, GV, US, fluxes, optics, nsw, j, dt, endif ! net volume/mass of liquid and solid passing through surface boundary fluxes - netMassInOut(i) = dt * (scale * (((((( fluxes%lprec(i,j) & + netMassInOut(i) = dt * (scale * US%R_to_kg_m3*US%Z_to_m*US%s_to_T*& + (((((( fluxes%lprec(i,j) & + fluxes%fprec(i,j) ) & + fluxes%evap(i,j) ) & + fluxes%lrunoff(i,j) ) & + fluxes%vprec(i,j) ) & + fluxes%seaice_melt(i,j)) & - + fluxes%frunoff(i,j) )) + + fluxes%frunoff(i,j) )) if (do_NMIOr) then ! Repeat the above code w/ dt=1s for legacy reasons - netMassInOut_rate(i) = (scale * (((((( fluxes%lprec(i,j) & + netMassInOut_rate(i) = (scale * US%R_to_kg_m3*US%Z_to_m*US%s_to_T* & + (((((( fluxes%lprec(i,j) & + fluxes%fprec(i,j) ) & + fluxes%evap(i,j) ) & + fluxes%lrunoff(i,j) ) & @@ -545,25 +547,25 @@ subroutine extractFluxes1d(G, GV, US, fluxes, optics, nsw, j, dt, ! evap < 0 means evaporation of water from the ocean, in ! which case heat_content_evap is computed in MOM_diabatic_driver.F90 if (fluxes%evap(i,j) < 0.0) then - netMassOut(i) = netMassOut(i) + fluxes%evap(i,j) + netMassOut(i) = netMassOut(i) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%evap(i,j) ! if (associated(fluxes%heat_content_cond)) fluxes%heat_content_cond(i,j) = 0.0 !??? --AJA endif ! lprec < 0 means sea ice formation taking water from the ocean. ! smg: we should split the ice melt/formation from the lprec if (fluxes%lprec(i,j) < 0.0) then - netMassOut(i) = netMassOut(i) + fluxes%lprec(i,j) + netMassOut(i) = netMassOut(i) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%lprec(i,j) endif ! seaice_melt < 0 means sea ice formation taking water from the ocean. if (fluxes%seaice_melt(i,j) < 0.0) then - netMassOut(i) = netMassOut(i) + fluxes%seaice_melt(i,j) + netMassOut(i) = netMassOut(i) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%seaice_melt(i,j) endif ! vprec < 0 means virtual evaporation arising from surface salinity restoring, ! in which case heat_content_vprec is computed in MOM_diabatic_driver.F90. if (fluxes%vprec(i,j) < 0.0) then - netMassOut(i) = netMassOut(i) + fluxes%vprec(i,j) + netMassOut(i) = netMassOut(i) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%vprec(i,j) endif netMassOut(i) = dt * scale * netMassOut(i) @@ -603,15 +605,15 @@ subroutine extractFluxes1d(G, GV, US, fluxes, optics, nsw, j, dt, if (useRiverHeatContent) then ! remove lrunoff*SST here, to counteract its addition elsewhere net_heat(i) = (net_heat(i) + (scale*(dt*J_m2_to_H)) * fluxes%heat_content_lrunoff(i,j)) - & - (GV%kg_m2_to_H * (scale * dt)) * fluxes%lrunoff(i,j) * T(i,1) + (GV%kg_m2_to_H * (scale * dt)) * US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%lrunoff(i,j) * T(i,1) !BGR-Jul 5, 2017{ !Intentionally neglect the following contribution to rate for legacy reasons. !if (do_NHR) net_heat_rate(i) = (net_heat_rate(i) + (scale*(J_m2_to_H)) * fluxes%heat_content_lrunoff(i,j)) - & - ! (GV%kg_m2_to_H * (scale)) * fluxes%lrunoff(i,j) * T(i,1) + ! (GV%kg_m2_to_H * (scale)) * US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%lrunoff(i,j) * T(i,1) !}BGR if (calculate_diags .and. associated(tv%TempxPmE)) then tv%TempxPmE(i,j) = tv%TempxPmE(i,j) + (scale * dt) * & - (I_Cp*fluxes%heat_content_lrunoff(i,j) - fluxes%lrunoff(i,j)*T(i,1)) + (I_Cp*fluxes%heat_content_lrunoff(i,j) - US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%lrunoff(i,j)*T(i,1)) endif endif @@ -620,15 +622,15 @@ subroutine extractFluxes1d(G, GV, US, fluxes, optics, nsw, j, dt, if (useCalvingHeatContent) then ! remove frunoff*SST here, to counteract its addition elsewhere net_heat(i) = net_heat(i) + (scale*(dt*J_m2_to_H)) * fluxes%heat_content_frunoff(i,j) - & - (GV%kg_m2_to_H * (scale * dt)) * fluxes%frunoff(i,j) * T(i,1) + (GV%kg_m2_to_H * (scale * dt)) * US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%frunoff(i,j) * T(i,1) !BGR-Jul 5, 2017{ !Intentionally neglect the following contribution to rate for legacy reasons. ! if (do_NHR) net_heat_rate(i) = net_heat_rate(i) + (scale*(J_m2_to_H)) * fluxes%heat_content_frunoff(i,j) - & -! (GV%kg_m2_to_H * (scale)) * fluxes%frunoff(i,j) * T(i,1) +! (GV%kg_m2_to_H * (scale)) * US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%frunoff(i,j) * T(i,1) !}BGR if (calculate_diags .and. associated(tv%TempxPmE)) then tv%TempxPmE(i,j) = tv%TempxPmE(i,j) + (scale * dt) * & - (I_Cp*fluxes%heat_content_frunoff(i,j) - fluxes%frunoff(i,j)*T(i,1)) + (I_Cp*fluxes%heat_content_frunoff(i,j) - US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%frunoff(i,j)*T(i,1)) endif endif @@ -730,7 +732,7 @@ subroutine extractFluxes1d(G, GV, US, fluxes, optics, nsw, j, dt, ! wait until MOM_diabatic_driver.F90. if (associated(fluxes%heat_content_lprec)) then if (fluxes%lprec(i,j) > 0.0) then - fluxes%heat_content_lprec(i,j) = fluxes%C_p*fluxes%lprec(i,j)*T(i,1) + fluxes%heat_content_lprec(i,j) = fluxes%C_p*US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%lprec(i,j)*T(i,1) else fluxes%heat_content_lprec(i,j) = 0.0 endif @@ -741,7 +743,7 @@ subroutine extractFluxes1d(G, GV, US, fluxes, optics, nsw, j, dt, ! and until we do so fprec is treated like lprec and enters at SST. -AJA if (associated(fluxes%heat_content_fprec)) then if (fluxes%fprec(i,j) > 0.0) then - fluxes%heat_content_fprec(i,j) = fluxes%C_p*fluxes%fprec(i,j)*T(i,1) + fluxes%heat_content_fprec(i,j) = fluxes%C_p*US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%fprec(i,j)*T(i,1) else fluxes%heat_content_fprec(i,j) = 0.0 endif @@ -750,7 +752,7 @@ subroutine extractFluxes1d(G, GV, US, fluxes, optics, nsw, j, dt, ! Following lprec and fprec, water flux due to sea ice melt (seaice_melt) enters at SST - GMM if (associated(fluxes%heat_content_icemelt)) then if (fluxes%seaice_melt(i,j) > 0.0) then - fluxes%heat_content_icemelt(i,j) = fluxes%C_p*fluxes%seaice_melt(i,j)*T(i,1) + fluxes%heat_content_icemelt(i,j) = fluxes%C_p*US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%seaice_melt(i,j)*T(i,1) else fluxes%heat_content_icemelt(i,j) = 0.0 endif @@ -761,7 +763,7 @@ subroutine extractFluxes1d(G, GV, US, fluxes, optics, nsw, j, dt, ! vprec < 0 means remove water from ocean; set heat_content_vprec in MOM_diabatic_driver.F90 if (associated(fluxes%heat_content_vprec)) then if (fluxes%vprec(i,j) > 0.0) then - fluxes%heat_content_vprec(i,j) = fluxes%C_p*fluxes%vprec(i,j)*T(i,1) + fluxes%heat_content_vprec(i,j) = fluxes%C_p*US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%vprec(i,j)*T(i,1) else fluxes%heat_content_vprec(i,j) = 0.0 endif @@ -775,7 +777,7 @@ subroutine extractFluxes1d(G, GV, US, fluxes, optics, nsw, j, dt, ! Condensation is assumed to drop into the ocean at the SST, just like lprec. if (associated(fluxes%heat_content_cond)) then if (fluxes%evap(i,j) > 0.0) then - fluxes%heat_content_cond(i,j) = fluxes%C_p*fluxes%evap(i,j)*T(i,1) + fluxes%heat_content_cond(i,j) = fluxes%C_p*US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%evap(i,j)*T(i,1) else fluxes%heat_content_cond(i,j) = 0.0 endif @@ -784,14 +786,14 @@ subroutine extractFluxes1d(G, GV, US, fluxes, optics, nsw, j, dt, ! Liquid runoff enters ocean at SST if land model does not provide runoff heat content. if (.not. useRiverHeatContent) then if (associated(fluxes%lrunoff) .and. associated(fluxes%heat_content_lrunoff)) then - fluxes%heat_content_lrunoff(i,j) = fluxes%C_p*fluxes%lrunoff(i,j)*T(i,1) + fluxes%heat_content_lrunoff(i,j) = fluxes%C_p*US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%lrunoff(i,j)*T(i,1) endif endif ! Icebergs enter ocean at SST if land model does not provide calving heat content. if (.not. useCalvingHeatContent) then if (associated(fluxes%frunoff) .and. associated(fluxes%heat_content_frunoff)) then - fluxes%heat_content_frunoff(i,j) = fluxes%C_p*fluxes%frunoff(i,j)*T(i,1) + fluxes%heat_content_frunoff(i,j) = fluxes%C_p*US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%frunoff(i,j)*T(i,1) endif endif @@ -1045,15 +1047,16 @@ subroutine MOM_forcing_chksum(mesg, fluxes, G, US, haloshift) if (associated(fluxes%sens)) & call hchksum(fluxes%sens, mesg//" fluxes%sens",G%HI,haloshift=hshift) if (associated(fluxes%evap)) & - call hchksum(fluxes%evap, mesg//" fluxes%evap",G%HI,haloshift=hshift) + call hchksum(fluxes%evap, mesg//" fluxes%evap",G%HI,haloshift=hshift, scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) if (associated(fluxes%lprec)) & - call hchksum(fluxes%lprec, mesg//" fluxes%lprec",G%HI,haloshift=hshift) + call hchksum(fluxes%lprec, mesg//" fluxes%lprec",G%HI,haloshift=hshift, scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) if (associated(fluxes%fprec)) & - call hchksum(fluxes%fprec, mesg//" fluxes%fprec",G%HI,haloshift=hshift) + call hchksum(fluxes%fprec, mesg//" fluxes%fprec",G%HI,haloshift=hshift, scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) if (associated(fluxes%vprec)) & - call hchksum(fluxes%vprec, mesg//" fluxes%vprec",G%HI,haloshift=hshift) + call hchksum(fluxes%vprec, mesg//" fluxes%vprec",G%HI,haloshift=hshift, scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) if (associated(fluxes%seaice_melt)) & - call hchksum(fluxes%seaice_melt, mesg//" fluxes%seaice_melt",G%HI,haloshift=hshift) + call hchksum(fluxes%seaice_melt, mesg//" fluxes%seaice_melt",G%HI,haloshift=hshift, & + scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) if (associated(fluxes%seaice_melt_heat)) & call hchksum(fluxes%seaice_melt_heat, mesg//" fluxes%seaice_melt_heat",G%HI,haloshift=hshift) if (associated(fluxes%p_surf)) & @@ -1066,9 +1069,9 @@ subroutine MOM_forcing_chksum(mesg, fluxes, G, US, haloshift) if (associated(fluxes%ustar_tidal)) & call hchksum(fluxes%ustar_tidal, mesg//" fluxes%ustar_tidal",G%HI,haloshift=hshift, scale=US%Z_to_m*US%s_to_T) if (associated(fluxes%lrunoff)) & - call hchksum(fluxes%lrunoff, mesg//" fluxes%lrunoff",G%HI,haloshift=hshift) + call hchksum(fluxes%lrunoff, mesg//" fluxes%lrunoff",G%HI,haloshift=hshift, scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) if (associated(fluxes%frunoff)) & - call hchksum(fluxes%frunoff, mesg//" fluxes%frunoff",G%HI,haloshift=hshift) + call hchksum(fluxes%frunoff, mesg//" fluxes%frunoff",G%HI,haloshift=hshift, scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) if (associated(fluxes%heat_content_lrunoff)) & call hchksum(fluxes%heat_content_lrunoff, mesg//" fluxes%heat_content_lrunoff",G%HI,haloshift=hshift) if (associated(fluxes%heat_content_frunoff)) & @@ -1294,12 +1297,14 @@ subroutine register_forcing_type_diags(Time, diag, US, use_temperature, handles, 'Liquid + frozen precipitation into ocean', 'kg m-2 s-1') handles%id_fprec = register_diag_field('ocean_model', 'fprec', diag%axesT1, Time, & - 'Frozen precipitation into ocean', 'kg m-2 s-1', & + 'Frozen precipitation into ocean', & + units='kg m-2 s-1', conversion=US%R_to_kg_m3*US%Z_to_m*US%s_to_T, & standard_name='snowfall_flux', cmor_field_name='prsn', & cmor_standard_name='snowfall_flux', cmor_long_name='Snowfall Flux where Ice Free Ocean over Sea') handles%id_lprec = register_diag_field('ocean_model', 'lprec', diag%axesT1, Time, & - 'Liquid precipitation into ocean', 'kg m-2 s-1', & + 'Liquid precipitation into ocean', & + units='kg m-2 s-1', conversion=US%R_to_kg_m3*US%Z_to_m*US%s_to_T, & standard_name='rainfall_flux', & cmor_field_name='prlq', cmor_standard_name='rainfall_flux', & cmor_long_name='Rainfall Flux where Ice Free Ocean over Sea') @@ -2095,21 +2100,23 @@ end subroutine set_derived_forcing_fields !> This subroutine determines the net mass source to the ocean from !! a (thermodynamic) forcing type and stores it in a mech_forcing type. -subroutine set_net_mass_forcing(fluxes, forces, G) +subroutine set_net_mass_forcing(fluxes, forces, G, US) type(forcing), intent(in) :: fluxes !< A structure containing thermodynamic forcing fields type(mech_forcing), intent(inout) :: forces !< A structure with the driving mechanical forces + type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(ocean_grid_type), intent(in) :: G !< The ocean grid type if (associated(forces%net_mass_src)) & - call get_net_mass_forcing(fluxes, G, forces%net_mass_src) + call get_net_mass_forcing(fluxes, G, US, forces%net_mass_src) end subroutine set_net_mass_forcing !> This subroutine calculates determines the net mass source to the ocean from !! a (thermodynamic) forcing type and stores it in a provided array. -subroutine get_net_mass_forcing(fluxes, G, net_mass_src) +subroutine get_net_mass_forcing(fluxes, G, US, net_mass_src) type(forcing), intent(in) :: fluxes !< A structure containing thermodynamic forcing fields type(ocean_grid_type), intent(in) :: G !< The ocean grid type + type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type real, dimension(SZI_(G),SZJ_(G)), intent(out) :: net_mass_src !< The net mass flux of water into the ocean !! [kg m-2 s-1]. @@ -2118,25 +2125,25 @@ subroutine get_net_mass_forcing(fluxes, G, net_mass_src) net_mass_src(:,:) = 0.0 if (associated(fluxes%lprec)) then ; do j=js,je ; do i=is,ie - net_mass_src(i,j) = net_mass_src(i,j) + fluxes%lprec(i,j) + net_mass_src(i,j) = net_mass_src(i,j) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%lprec(i,j) enddo ; enddo ; endif if (associated(fluxes%fprec)) then ; do j=js,je ; do i=is,ie - net_mass_src(i,j) = net_mass_src(i,j) + fluxes%fprec(i,j) + net_mass_src(i,j) = net_mass_src(i,j) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%fprec(i,j) enddo ; enddo ; endif if (associated(fluxes%vprec)) then ; do j=js,je ; do i=is,ie - net_mass_src(i,j) = net_mass_src(i,j) + fluxes%vprec(i,j) + net_mass_src(i,j) = net_mass_src(i,j) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%vprec(i,j) enddo ; enddo ; endif if (associated(fluxes%lrunoff)) then ; do j=js,je ; do i=is,ie - net_mass_src(i,j) = net_mass_src(i,j) + fluxes%lrunoff(i,j) + net_mass_src(i,j) = net_mass_src(i,j) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%lrunoff(i,j) enddo ; enddo ; endif if (associated(fluxes%frunoff)) then ; do j=js,je ; do i=is,ie - net_mass_src(i,j) = net_mass_src(i,j) + fluxes%frunoff(i,j) + net_mass_src(i,j) = net_mass_src(i,j) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%frunoff(i,j) enddo ; enddo ; endif if (associated(fluxes%evap)) then ; do j=js,je ; do i=is,ie - net_mass_src(i,j) = net_mass_src(i,j) + fluxes%evap(i,j) + net_mass_src(i,j) = net_mass_src(i,j) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%evap(i,j) enddo ; enddo ; endif if (associated(fluxes%seaice_melt)) then ; do j=js,je ; do i=is,ie - net_mass_src(i,j) = net_mass_src(i,j) + fluxes%seaice_melt(i,j) + net_mass_src(i,j) = net_mass_src(i,j) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%seaice_melt(i,j) enddo ; enddo ; endif end subroutine get_net_mass_forcing @@ -2196,12 +2203,13 @@ end subroutine mech_forcing_diags !> Offer buoyancy forcing fields for diagnostics for those !! fields registered as part of register_forcing_type_diags. -subroutine forcing_diagnostics(fluxes, sfc_state, dt, G, diag, handles) +subroutine forcing_diagnostics(fluxes, sfc_state, dt, G, US, diag, handles) type(forcing), intent(in) :: fluxes !< A structure containing thermodynamic forcing fields type(surface), intent(in) :: sfc_state !< A structure containing fields that !! describe the surface state of the ocean. real, intent(in) :: dt !< time step type(ocean_grid_type), intent(in) :: G !< grid type + type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(diag_ctrl), intent(in) :: diag !< diagnostic regulator type(forcing_diags), intent(inout) :: handles !< diagnostic ids @@ -2228,14 +2236,15 @@ subroutine forcing_diagnostics(fluxes, sfc_state, dt, G, diag, handles) if (handles%id_prcme > 0 .or. handles%id_total_prcme > 0 .or. handles%id_prcme_ga > 0) then do j=js,je ; do i=is,ie res(i,j) = 0.0 - if (associated(fluxes%lprec)) res(i,j) = res(i,j)+fluxes%lprec(i,j) - if (associated(fluxes%fprec)) res(i,j) = res(i,j)+fluxes%fprec(i,j) + if (associated(fluxes%lprec)) res(i,j) = res(i,j)+US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%lprec(i,j) + if (associated(fluxes%fprec)) res(i,j) = res(i,j)+US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%fprec(i,j) ! fluxes%cond is not needed because it is derived from %evap > 0 - if (associated(fluxes%evap)) res(i,j) = res(i,j)+fluxes%evap(i,j) - if (associated(fluxes%lrunoff)) res(i,j) = res(i,j)+fluxes%lrunoff(i,j) - if (associated(fluxes%frunoff)) res(i,j) = res(i,j)+fluxes%frunoff(i,j) - if (associated(fluxes%vprec)) res(i,j) = res(i,j)+fluxes%vprec(i,j) - if (associated(fluxes%seaice_melt)) res(i,j) = res(i,j)+fluxes%seaice_melt(i,j) + if (associated(fluxes%evap)) res(i,j) = res(i,j)+US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%evap(i,j) + if (associated(fluxes%lrunoff)) res(i,j) = res(i,j)+US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%lrunoff(i,j) + if (associated(fluxes%frunoff)) res(i,j) = res(i,j)+US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%frunoff(i,j) + if (associated(fluxes%vprec)) res(i,j) = res(i,j)+US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%vprec(i,j) + if (associated(fluxes%seaice_melt)) res(i,j) = res(i,j) + & + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%seaice_melt(i,j) enddo ; enddo if (handles%id_prcme > 0) call post_data(handles%id_prcme, res, diag) if (handles%id_total_prcme > 0) then @@ -2252,17 +2261,17 @@ subroutine forcing_diagnostics(fluxes, sfc_state, dt, G, diag, handles) do j=js,je ; do i=is,ie res(i,j) = 0.0 if (associated(fluxes%lprec)) then - if (fluxes%lprec(i,j) < 0.0) res(i,j) = res(i,j) + fluxes%lprec(i,j) + if (fluxes%lprec(i,j) < 0.0) res(i,j) = res(i,j) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%lprec(i,j) endif if (associated(fluxes%vprec)) then - if (fluxes%vprec(i,j) < 0.0) res(i,j) = res(i,j) + fluxes%vprec(i,j) + if (fluxes%vprec(i,j) < 0.0) res(i,j) = res(i,j) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%vprec(i,j) endif if (associated(fluxes%evap)) then - if (fluxes%evap(i,j) < 0.0) res(i,j) = res(i,j) + fluxes%evap(i,j) + if (fluxes%evap(i,j) < 0.0) res(i,j) = res(i,j) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%evap(i,j) endif if (associated(fluxes%seaice_melt)) then if (fluxes%seaice_melt(i,j) < 0.0) & - res(i,j) = res(i,j) + fluxes%seaice_melt(i,j) + res(i,j) = res(i,j) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%seaice_melt(i,j) endif enddo ; enddo if (handles%id_net_massout > 0) call post_data(handles%id_net_massout, res, diag) @@ -2280,25 +2289,25 @@ subroutine forcing_diagnostics(fluxes, sfc_state, dt, G, diag, handles) res(i,j) = 0.0 if (associated(fluxes%fprec)) & - res(i,j) = res(i,j) + fluxes%fprec(i,j) + res(i,j) = res(i,j) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%fprec(i,j) if (associated(fluxes%lrunoff)) & - res(i,j) = res(i,j) + fluxes%lrunoff(i,j) + res(i,j) = res(i,j) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%lrunoff(i,j) if (associated(fluxes%frunoff)) & - res(i,j) = res(i,j) + fluxes%frunoff(i,j) + res(i,j) = res(i,j) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%frunoff(i,j) if (associated(fluxes%lprec)) then - if (fluxes%lprec(i,j) > 0.0) res(i,j) = res(i,j) + fluxes%lprec(i,j) + if (fluxes%lprec(i,j) > 0.0) res(i,j) = res(i,j) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%lprec(i,j) endif if (associated(fluxes%vprec)) then - if (fluxes%vprec(i,j) > 0.0) res(i,j) = res(i,j) + fluxes%vprec(i,j) + if (fluxes%vprec(i,j) > 0.0) res(i,j) = res(i,j) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%vprec(i,j) endif ! fluxes%cond is not needed because it is derived from %evap > 0 if (associated(fluxes%evap)) then - if (fluxes%evap(i,j) > 0.0) res(i,j) = res(i,j) + fluxes%evap(i,j) + if (fluxes%evap(i,j) > 0.0) res(i,j) = res(i,j) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%evap(i,j) endif if (associated(fluxes%seaice_melt)) then if (fluxes%seaice_melt(i,j) > 0.0) & - res(i,j) = res(i,j) + fluxes%seaice_melt(i,j) + res(i,j) = res(i,j) + US%R_to_kg_m3*US%Z_to_m*US%s_to_T*fluxes%seaice_melt(i,j) endif enddo ; enddo if (handles%id_net_massin > 0) call post_data(handles%id_net_massin, res, diag) @@ -2314,17 +2323,17 @@ subroutine forcing_diagnostics(fluxes, sfc_state, dt, G, diag, handles) if ((handles%id_evap > 0) .and. associated(fluxes%evap)) & call post_data(handles%id_evap, fluxes%evap, diag) if ((handles%id_total_evap > 0) .and. associated(fluxes%evap)) then - total_transport = global_area_integral(fluxes%evap,G) + total_transport = global_area_integral(fluxes%evap, G, scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) call post_data(handles%id_total_evap, total_transport, diag) endif if ((handles%id_evap_ga > 0) .and. associated(fluxes%evap)) then - ave_flux = global_area_mean(fluxes%evap,G) + ave_flux = global_area_mean(fluxes%evap, G, scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) call post_data(handles%id_evap_ga, ave_flux, diag) endif if (associated(fluxes%lprec) .and. associated(fluxes%fprec)) then do j=js,je ; do i=is,ie - res(i,j) = fluxes%lprec(i,j) + fluxes%fprec(i,j) + res(i,j) = US%R_to_kg_m3*US%Z_to_m*US%s_to_T* (fluxes%lprec(i,j) + fluxes%fprec(i,j)) enddo ; enddo if (handles%id_precip > 0) call post_data(handles%id_precip, res, diag) if (handles%id_total_precip > 0) then @@ -2340,11 +2349,11 @@ subroutine forcing_diagnostics(fluxes, sfc_state, dt, G, diag, handles) if (associated(fluxes%lprec)) then if (handles%id_lprec > 0) call post_data(handles%id_lprec, fluxes%lprec, diag) if (handles%id_total_lprec > 0) then - total_transport = global_area_integral(fluxes%lprec,G) + total_transport = global_area_integral(fluxes%lprec, G, scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) call post_data(handles%id_total_lprec, total_transport, diag) endif if (handles%id_lprec_ga > 0) then - ave_flux = global_area_mean(fluxes%lprec,G) + ave_flux = global_area_mean(fluxes%lprec, G, scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) call post_data(handles%id_lprec_ga, ave_flux, diag) endif endif @@ -2352,11 +2361,11 @@ subroutine forcing_diagnostics(fluxes, sfc_state, dt, G, diag, handles) if (associated(fluxes%fprec)) then if (handles%id_fprec > 0) call post_data(handles%id_fprec, fluxes%fprec, diag) if (handles%id_total_fprec > 0) then - total_transport = global_area_integral(fluxes%fprec,G) + total_transport = global_area_integral(fluxes%fprec ,G, scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) call post_data(handles%id_total_fprec, total_transport, diag) endif if (handles%id_fprec_ga > 0) then - ave_flux = global_area_mean(fluxes%fprec,G) + ave_flux = global_area_mean(fluxes%fprec, G, scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) call post_data(handles%id_fprec_ga, ave_flux, diag) endif endif @@ -2364,11 +2373,11 @@ subroutine forcing_diagnostics(fluxes, sfc_state, dt, G, diag, handles) if (associated(fluxes%vprec)) then if (handles%id_vprec > 0) call post_data(handles%id_vprec, fluxes%vprec, diag) if (handles%id_total_vprec > 0) then - total_transport = global_area_integral(fluxes%vprec,G) + total_transport = global_area_integral(fluxes%vprec, G, scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) call post_data(handles%id_total_vprec, total_transport, diag) endif if (handles%id_vprec_ga > 0) then - ave_flux = global_area_mean(fluxes%vprec,G) + ave_flux = global_area_mean(fluxes%vprec, G, scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) call post_data(handles%id_vprec_ga, ave_flux, diag) endif endif @@ -2376,7 +2385,7 @@ subroutine forcing_diagnostics(fluxes, sfc_state, dt, G, diag, handles) if (associated(fluxes%lrunoff)) then if (handles%id_lrunoff > 0) call post_data(handles%id_lrunoff, fluxes%lrunoff, diag) if (handles%id_total_lrunoff > 0) then - total_transport = global_area_integral(fluxes%lrunoff,G) + total_transport = global_area_integral(fluxes%lrunoff, G, scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) call post_data(handles%id_total_lrunoff, total_transport, diag) endif endif @@ -2384,7 +2393,7 @@ subroutine forcing_diagnostics(fluxes, sfc_state, dt, G, diag, handles) if (associated(fluxes%frunoff)) then if (handles%id_frunoff > 0) call post_data(handles%id_frunoff, fluxes%frunoff, diag) if (handles%id_total_frunoff > 0) then - total_transport = global_area_integral(fluxes%frunoff,G) + total_transport = global_area_integral(fluxes%frunoff, G, scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) call post_data(handles%id_total_frunoff, total_transport, diag) endif endif @@ -2392,7 +2401,7 @@ subroutine forcing_diagnostics(fluxes, sfc_state, dt, G, diag, handles) if (associated(fluxes%seaice_melt)) then if (handles%id_seaice_melt > 0) call post_data(handles%id_seaice_melt, fluxes%seaice_melt, diag) if (handles%id_total_seaice_melt > 0) then - total_transport = global_area_integral(fluxes%seaice_melt,G) + total_transport = global_area_integral(fluxes%seaice_melt, G, scale=US%R_to_kg_m3*US%Z_to_m*US%s_to_T) call post_data(handles%id_total_seaice_melt, total_transport, diag) endif endif diff --git a/src/diagnostics/MOM_sum_output.F90 b/src/diagnostics/MOM_sum_output.F90 index 1a8a9879b3..7bb8ba73e9 100644 --- a/src/diagnostics/MOM_sum_output.F90 +++ b/src/diagnostics/MOM_sum_output.F90 @@ -936,7 +936,7 @@ end subroutine write_energy !> This subroutine accumates the net input of volume, salt and heat, through !! the ocean surface for use in diagnosing conservation. -subroutine accumulate_net_input(fluxes, sfc_state, tv, dt, G, CS) +subroutine accumulate_net_input(fluxes, sfc_state, tv, dt, G, US, CS) type(forcing), intent(in) :: fluxes !< A structure containing pointers to any possible !! forcing fields. Unused fields are unallocated. type(surface), intent(in) :: sfc_state !< A structure containing fields that @@ -945,6 +945,7 @@ subroutine accumulate_net_input(fluxes, sfc_state, tv, dt, G, CS) !! thermodynamic variables. real, intent(in) :: dt !< The amount of time over which to average [s]. type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure. + type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(Sum_output_CS), pointer :: CS !< The control structure returned by a previous call !! to MOM_sum_output_init. ! Local variables @@ -977,7 +978,7 @@ subroutine accumulate_net_input(fluxes, sfc_state, tv, dt, G, CS) if (associated(fluxes%evap)) then if (associated(fluxes%lprec) .and. associated(fluxes%fprec)) then do j=js,je ; do i=is,ie - FW_in(i,j) = dt*G%US%L_to_m**2*G%areaT(i,j)*(fluxes%evap(i,j) + & + FW_in(i,j) = G%US%L_to_m**2*US%R_to_kg_m3*US%Z_to_m*US%s_to_T*dt*G%areaT(i,j)*(fluxes%evap(i,j) + & (((fluxes%lprec(i,j) + fluxes%vprec(i,j)) + fluxes%lrunoff(i,j)) + & (fluxes%fprec(i,j) + fluxes%frunoff(i,j)))) enddo ; enddo @@ -988,7 +989,8 @@ subroutine accumulate_net_input(fluxes, sfc_state, tv, dt, G, CS) endif if (associated(fluxes%seaice_melt)) then ; do j=js,je ; do i=is,ie - FW_in(i,j) = FW_in(i,j) + dt * G%US%L_to_m**2*G%areaT(i,j) * fluxes%seaice_melt(i,j) + FW_in(i,j) = FW_in(i,j) + G%US%L_to_m**2*US%R_to_kg_m3*US%Z_to_m*US%s_to_T*dt * & + G%areaT(i,j) * fluxes%seaice_melt(i,j) enddo ; enddo ; endif salt_in(:,:) = 0.0 ; heat_in(:,:) = 0.0 diff --git a/src/ice_shelf/MOM_ice_shelf.F90 b/src/ice_shelf/MOM_ice_shelf.F90 index ca8f3049ee..76f595ee06 100644 --- a/src/ice_shelf/MOM_ice_shelf.F90 +++ b/src/ice_shelf/MOM_ice_shelf.F90 @@ -967,10 +967,10 @@ subroutine add_shelf_flux(G, US, CS, state, fluxes) if (associated(fluxes%evap)) fluxes%evap(i,j) = 0.0 if (associated(fluxes%lprec)) then if (ISS%water_flux(i,j) > 0.0) then - fluxes%lprec(i,j) = frac_area*ISS%water_flux(i,j)*CS%flux_factor + fluxes%lprec(i,j) = US%kg_m3_to_R*US%m_to_Z*US%T_to_s*frac_area*ISS%water_flux(i,j)*CS%flux_factor else fluxes%lprec(i,j) = 0.0 - fluxes%evap(i,j) = frac_area*ISS%water_flux(i,j)*CS%flux_factor + fluxes%evap(i,j) = US%kg_m3_to_R*US%m_to_Z*US%T_to_s*frac_area*ISS%water_flux(i,j)*CS%flux_factor endif endif @@ -1061,6 +1061,8 @@ subroutine add_shelf_flux(G, US, CS, state, fluxes) fluxes%vprec(i,j) = -mean_melt_flux * CS%density_ice/1000. ! evap is negative fluxes%sens(i,j) = fluxes%vprec(i,j) * CS%Cp * CS%T0 ! W /m^2 fluxes%salt_flux(i,j) = fluxes%vprec(i,j) * CS%S0*1.0e-3 ! kg (salt)/(m^2 s) + ! Rescale fluxes%vprec to the proper units. + fluxes%vprec(i,j) = US%kg_m3_to_R*US%m_to_Z*US%T_to_s * fluxes%vprec(i,j) endif enddo ; enddo diff --git a/src/ice_shelf/MOM_marine_ice.F90 b/src/ice_shelf/MOM_marine_ice.F90 index 16b543387d..4042681803 100644 --- a/src/ice_shelf/MOM_marine_ice.F90 +++ b/src/ice_shelf/MOM_marine_ice.F90 @@ -15,6 +15,7 @@ module MOM_marine_ice use MOM_forcing_type, only : forcing, mech_forcing use MOM_grid, only : ocean_grid_type use MOM_time_manager, only : time_type +use MOM_unit_scaling, only : unit_scale_type use MOM_variables, only : surface implicit none ; private @@ -102,9 +103,10 @@ end subroutine iceberg_forces !> iceberg_fluxes adds ice-area-coverage and modifies various !! thermodynamic fluxes due to the presence of icebergs. -subroutine iceberg_fluxes(G, fluxes, use_ice_shelf, sfc_state, & +subroutine iceberg_fluxes(G, US, fluxes, use_ice_shelf, sfc_state, & time_step, CS) type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure + type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(forcing), intent(inout) :: fluxes !< A structure with pointers to themodynamic, !! tracer and mass exchange forcing fields type(surface), intent(inout) :: sfc_state !< A structure containing fields that @@ -113,8 +115,8 @@ subroutine iceberg_fluxes(G, fluxes, use_ice_shelf, sfc_state, & real, intent(in) :: time_step !< The coupling time step [s]. type(marine_ice_CS), pointer :: CS !< Pointer to the control structure for MOM_marine_ice - real :: fraz ! refreezing rate [kg m-2 s-1] - real :: I_dt_LHF ! The inverse of the timestep times the latent heat of fusion [kg J-1 s-1]. + real :: fraz ! refreezing rate [R Z T-1 ~> kg m-2 s-1] + real :: I_dt_LHF ! The inverse of the timestep times the latent heat of fusion [kg J-1 T-1 ~> kg J-1 s-1]. integer :: i, j, is, ie, js, je, isd, ied, jsd, jed is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec isd = G%isd ; jsd = G%jsd ; ied = G%ied ; jed = G%jed @@ -142,7 +144,7 @@ subroutine iceberg_fluxes(G, fluxes, use_ice_shelf, sfc_state, & !Zero'ing out other fluxes under the tabular icebergs if (CS%berg_area_threshold >= 0.) then - I_dt_LHF = 1.0 / (time_step * CS%latent_heat_fusion) + I_dt_LHF = 1.0 / (US%s_to_T*time_step * CS%latent_heat_fusion) do j=jsd,jed ; do i=isd,ied if (fluxes%frac_shelf_h(i,j) > CS%berg_area_threshold) then ! Only applying for ice shelf covering most of cell. @@ -153,13 +155,14 @@ subroutine iceberg_fluxes(G, fluxes, use_ice_shelf, sfc_state, & if (associated(fluxes%evap)) fluxes%evap(i,j) = 0.0 ! Add frazil formation diagnosed by the ocean model [J m-2] in the - ! form of surface layer evaporation [kg m-2 s-1]. Update lprec in the + ! form of surface layer evaporation [R Z T-1 ~> kg m-2 s-1]. Update lprec in the ! control structure for diagnostic purposes. if (associated(sfc_state%frazil)) then - fraz = sfc_state%frazil(i,j) * I_dt_LHF - if (associated(fluxes%evap)) fluxes%evap(i,j) = fluxes%evap(i,j) - fraz - !CS%lprec(i,j)=CS%lprec(i,j) - fraz + fraz = US%kg_m3_to_R*US%m_to_Z*sfc_state%frazil(i,j) * I_dt_LHF + if (associated(fluxes%evap)) & + fluxes%evap(i,j) = fluxes%evap(i,j) - fraz + ! fluxes%lprec(i,j) = fluxes%lprec(i,j) - fraz sfc_state%frazil(i,j) = 0.0 endif diff --git a/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 b/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 index aa101fb9f1..e09b21f251 100644 --- a/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 +++ b/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 @@ -514,10 +514,10 @@ subroutine bulkmixedlayer(h_3d, u_3d, v_3d, tv, fluxes, dt_in_T, ea, eb, G, GV, ! rivermix_depth = The prescribed depth over which to mix river inflow ! drho_ds = The gradient of density wrt salt at the ambient surface salinity. ! Sriver = 0 (i.e. rivers are assumed to be pure freshwater) - RmixConst = 0.5*CS%rivermix_depth * (GV%g_Earth*US%m_to_Z) * Irho0**2 + RmixConst = 0.5*CS%rivermix_depth * GV%g_Earth * Irho0**2 do i=is,ie TKE_river(i) = max(0.0, RmixConst*dR0_dS(i)* & - US%T_to_s*US%kg_m3_to_R*(fluxes%lrunoff(i,j) + fluxes%frunoff(i,j)) * S(i,1)) + (fluxes%lrunoff(i,j) + fluxes%frunoff(i,j)) * S(i,1)) enddo else do i=is,ie ; TKE_river(i) = 0.0 ; enddo diff --git a/src/parameterizations/vertical/MOM_diabatic_aux.F90 b/src/parameterizations/vertical/MOM_diabatic_aux.F90 index e21323f6b8..ad2f57f2d4 100644 --- a/src/parameterizations/vertical/MOM_diabatic_aux.F90 +++ b/src/parameterizations/vertical/MOM_diabatic_aux.F90 @@ -1130,12 +1130,12 @@ subroutine applyBoundaryFluxesInOut(CS, G, GV, US, dt, fluxes, optics, nsw, h, t ! drho_ds = The gradient of density wrt salt at the ambient surface salinity. ! Sriver = 0 (i.e. rivers are assumed to be pure freshwater) if (GV%Boussinesq) then - RivermixConst = -0.5*(CS%rivermix_depth*dt)*(US%m_to_Z) * ( US%L_to_Z**2*GV%g_Earth ) * GV%Rho0 + RivermixConst = -0.5*(CS%rivermix_depth*US%s_to_T*dt) * ( US%L_to_Z**2*GV%g_Earth ) * GV%Rho0 else - RivermixConst = -0.5*(CS%rivermix_depth*dt)*(US%m_to_Z) * GV%Rho0 * ( US%L_to_Z**2*GV%g_Earth ) + RivermixConst = -0.5*(CS%rivermix_depth*US%s_to_T*dt) * GV%Rho0 * ( US%L_to_Z**2*GV%g_Earth ) endif cTKE(i,j,k) = cTKE(i,j,k) + max(0.0, RivermixConst*dSV_dS(i,j,1) * & - US%kg_m3_to_R*(fluxes%lrunoff(i,j) + fluxes%frunoff(i,j)) * tv%S(i,j,1)) + (fluxes%lrunoff(i,j) + fluxes%frunoff(i,j)) * tv%S(i,j,1)) endif ! Update state diff --git a/src/tracer/MOM_generic_tracer.F90 b/src/tracer/MOM_generic_tracer.F90 index d12897038f..28f31c6fa1 100644 --- a/src/tracer/MOM_generic_tracer.F90 +++ b/src/tracer/MOM_generic_tracer.F90 @@ -457,7 +457,8 @@ subroutine MOM_generic_tracer_column_physics(h_old, h_new, ea, eb, fluxes, Hml, call g_tracer_get_pointer(g_tracer,g_tracer_name,'trunoff',trunoff_array) call g_tracer_get_pointer(g_tracer,g_tracer_name,'runoff_tracer_flux',runoff_tracer_flux_array) !nnz: Why is fluxes%river = 0? - runoff_tracer_flux_array = trunoff_array * fluxes%lrunoff + runoff_tracer_flux_array(:,:) = trunoff_array(:,:) * & + G%US%R_to_kg_m3*G%US%Z_to_m*G%US%s_to_T*fluxes%lrunoff(:,:) stf_array = stf_array + runoff_tracer_flux_array endif @@ -492,9 +493,10 @@ subroutine MOM_generic_tracer_column_physics(h_old, h_new, ea, eb, fluxes, Hml, ! call generic_tracer_source(tv%T,tv%S,rho_dzt,dzt,Hml,G%isd,G%jsd,1,dt,& - G%US%L_to_m**2*G%areaT, get_diag_time_end(CS%diag), & + G%US%L_to_m**2*G%areaT(:,:), get_diag_time_end(CS%diag), & optics%nbands, optics%max_wavelength_band, optics%sw_pen_band, optics%opacity_band, & - internal_heat=tv%internal_heat, frunoff=fluxes%frunoff, sosga=sosga) + internal_heat=tv%internal_heat, & + frunoff=G%US%R_to_kg_m3*G%US%Z_to_m*G%US%s_to_T*fluxes%frunoff(:,:), sosga=sosga) ! This uses applyTracerBoundaryFluxesInOut to handle the change in tracer due to freshwater fluxes ! usually in ALE mode diff --git a/src/user/BFB_surface_forcing.F90 b/src/user/BFB_surface_forcing.F90 index 558be86734..bce0698240 100644 --- a/src/user/BFB_surface_forcing.F90 +++ b/src/user/BFB_surface_forcing.F90 @@ -136,9 +136,8 @@ subroutine BFB_buoyancy_forcing(state, fluxes, day, dt, G, US, CS) fluxes%heat_added(i,j) = (G%mask2dT(i,j) * (rhoXcp * CS%Flux_const)) * & (Temp_restore - state%SST(i,j)) - fluxes%vprec(i,j) = - (G%mask2dT(i,j) * (CS%Rho0*CS%Flux_const)) * & - ((Salin_restore - state%SSS(i,j)) / & - (0.5 * (Salin_restore + state%SSS(i,j)))) + fluxes%vprec(i,j) = - (G%mask2dT(i,j) * (US%kg_m3_to_R*US%m_to_Z*US%T_to_s*CS%Rho0*CS%Flux_const)) * & + ((Salin_restore - state%SSS(i,j)) / (0.5 * (Salin_restore + state%SSS(i,j)))) enddo ; enddo else ! When modifying the code, comment out this error message. It is here diff --git a/src/user/SCM_CVMix_tests.F90 b/src/user/SCM_CVMix_tests.F90 index a61600fa56..960abd49ca 100644 --- a/src/user/SCM_CVMix_tests.F90 +++ b/src/user/SCM_CVMix_tests.F90 @@ -227,11 +227,12 @@ subroutine SCM_CVMix_tests_wind_forcing(state, forces, day, G, US, CS) end subroutine SCM_CVMix_tests_wind_forcing -subroutine SCM_CVMix_tests_buoyancy_forcing(state, fluxes, day, G, CS) +subroutine SCM_CVMix_tests_buoyancy_forcing(state, fluxes, day, G, US, CS) type(surface), intent(in) :: state !< Surface state structure type(forcing), intent(inout) :: fluxes !< Surface fluxes structure type(time_type), intent(in) :: day !< Current model time type(ocean_grid_type), intent(inout) :: G !< Grid structure + type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(SCM_CVMix_tests_CS), pointer :: CS !< Container for SCM parameters ! Local variables @@ -259,9 +260,9 @@ subroutine SCM_CVMix_tests_buoyancy_forcing(state, fluxes, day, G, CS) if (CS%UseEvaporation) then do J=Jsq,Jeq ; do i=is,ie ! Note CVMix test inputs give evaporation in [m s-1] - ! This therefore must be converted to mass flux - ! by multiplying by density - fluxes%evap(i,J) = CS%surf_evap * CS%Rho0 + ! This therefore must be converted to mass flux in [R Z T-1 ~> kg m-2 s-1] + ! by multiplying by density and some unit conversion factors. + fluxes%evap(i,J) = CS%surf_evap * US%kg_m3_to_R*US%m_to_Z*US%T_to_s * CS%Rho0 enddo ; enddo endif diff --git a/src/user/dumbbell_surface_forcing.F90 b/src/user/dumbbell_surface_forcing.F90 index d8b3ad269b..4f9483d7e5 100644 --- a/src/user/dumbbell_surface_forcing.F90 +++ b/src/user/dumbbell_surface_forcing.F90 @@ -47,7 +47,7 @@ module dumbbell_surface_forcing contains !> Surface buoyancy (heat and fresh water) fluxes for the dumbbell test case -subroutine dumbbell_buoyancy_forcing(state, fluxes, day, dt, G, CS) +subroutine dumbbell_buoyancy_forcing(state, fluxes, day, dt, G, US, CS) type(surface), intent(inout) :: state !< A structure containing fields that !! describe the surface state of the ocean. type(forcing), intent(inout) :: fluxes !< A structure containing pointers to any @@ -57,6 +57,7 @@ subroutine dumbbell_buoyancy_forcing(state, fluxes, day, dt, G, CS) real, intent(in) :: dt !< The amount of time over which !! the fluxes apply [s] type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure + type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(dumbbell_surface_forcing_CS), pointer :: CS !< A control structure returned by a previous !! call to dumbbell_surface_forcing_init ! Local variables @@ -123,9 +124,8 @@ subroutine dumbbell_buoyancy_forcing(state, fluxes, day, dt, G, CS) ! Set density_restore to an expression for the surface potential ! density [kg m-3] that is being restored toward. if (CS%forcing_mask(i,j)>0.) then - fluxes%vprec(i,j) = - (G%mask2dT(i,j) * (CS%Rho0*CS%Flux_const)) * & - ((CS%S_restore(i,j) - state%SSS(i,j)) / & - (0.5 * (CS%S_restore(i,j) + state%SSS(i,j)))) + fluxes%vprec(i,j) = - (G%mask2dT(i,j) * (US%kg_m3_to_R*US%m_to_Z*US%T_to_s*CS%Rho0*CS%Flux_const)) * & + ((CS%S_restore(i,j) - state%SSS(i,j)) / (0.5 * (CS%S_restore(i,j) + state%SSS(i,j)))) endif enddo ; enddo