From 8c09c14404650f9f31c3f0dcc56812c3e585bf1f Mon Sep 17 00:00:00 2001 From: "brandon.reichl" Date: Thu, 26 Aug 2021 11:06:20 -0400 Subject: [PATCH] Adds option to homogenize IOB fields in FMS_cap/ocean_model_MOM - With new option HOMOGENIZE_IOB_FORCINGS = True, a copy of the IOB structure is created with horizontally homogenized fields across all processor elements. - This new feature is primarily for running in single column like configurations with the coupler, which requires perfectly equal forcing across all cells. --- .../drivers/FMS_cap/ocean_model_MOM.F90 | 214 +++++++++++++++++- 1 file changed, 205 insertions(+), 9 deletions(-) diff --git a/config_src/drivers/FMS_cap/ocean_model_MOM.F90 b/config_src/drivers/FMS_cap/ocean_model_MOM.F90 index c3e13329f2..be80a973b3 100644 --- a/config_src/drivers/FMS_cap/ocean_model_MOM.F90 +++ b/config_src/drivers/FMS_cap/ocean_model_MOM.F90 @@ -25,6 +25,7 @@ module ocean_model_mod use MOM_diag_mediator, only : diag_mediator_close_registration, diag_mediator_end use MOM_domains, only : MOM_domain_type, domain2d, clone_MOM_domain, get_domain_extent use MOM_domains, only : pass_var, pass_vector, AGRID, BGRID_NE, CGRID_NE, TO_ALL, Omit_Corners +use MOM_domains, only : sum_across_PEs use MOM_error_handler, only : MOM_error, MOM_mesg, FATAL, WARNING, is_root_pe use MOM_error_handler, only : callTree_enter, callTree_leave use MOM_EOS, only : gsw_sp_from_sr, gsw_pt_from_ct @@ -164,7 +165,8 @@ module ocean_model_mod !! fields necessary to integrate only the tracer advection !! and diffusion equation read in from files stored from !! a previous integration of the prognostic model. - + logical :: use_homogenized_IOB !< Creates a homogenized version of the IOB structure + !! that yields homogeneous forcing/fluxes. logical :: single_step_call !< If true, advance the state of MOM with a single !! step including both dynamics and thermodynamics. !! If false, the two phases are advanced with @@ -358,6 +360,10 @@ subroutine ocean_model_init(Ocean_sfc, OS, Time_init, Time_in, wind_stagger, gas use_melt_pot=.false. endif + call get_param(param_file, mdl, "HOMOGENIZE_IOB_FORCINGS", OS%use_homogenized_IOB, & + "If true, creates a copy of all IOB fluxes and averages them "//& + "horizontally across the entire domain.", default=.false.) + !allocate(OS%sfc_state) call allocate_surface_state(OS%sfc_state, OS%grid, use_temperature, do_integrals=.true., & gas_fields_ocn=gas_fields_ocn, use_meltpot=use_melt_pot) @@ -446,6 +452,8 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, time_start_upda real, optional, intent(in) :: cycle_length !< The duration of a coupled time stepping cycle [s]. ! Local variables + type(ice_ocean_boundary_type) :: Ice_ocean_boundary_homogenized ! A copy of the structure containing the + ! forcing fields coming from the ice and atmosphere, that is horizontally averaged. type(time_type) :: Time_seg_start ! Stores the dynamic or thermodynamic ocean model time at the ! start of this call to allow step_MOM to temporarily change the time ! as seen by internal modules. @@ -504,9 +512,17 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, time_start_upda ! Translate Ice_ocean_boundary into fluxes and forces. call get_domain_extent(Ocean_sfc%Domain, index_bnds(1), index_bnds(2), index_bnds(3), index_bnds(4)) + if (OS%use_homogenized_IOB) & + call homogenize_IOB(Ice_ocean_boundary, Ice_ocean_boundary_homogenized, OS%grid, index_bnds) + if (do_dyn) then - call convert_IOB_to_forces(Ice_ocean_boundary, OS%forces, index_bnds, OS%Time_dyn, OS%grid, OS%US, & - OS%forcing_CSp, dt_forcing=dt_coupling, reset_avg=OS%fluxes%fluxes_used) + if (OS%use_homogenized_IOB) then + call convert_IOB_to_forces(Ice_ocean_boundary_homogenized, OS%forces, index_bnds, OS%Time_dyn, OS%grid, OS%US, & + OS%forcing_CSp, dt_forcing=dt_coupling, reset_avg=OS%fluxes%fluxes_used) + else + call convert_IOB_to_forces(Ice_ocean_boundary, OS%forces, index_bnds, OS%Time_dyn, OS%grid, OS%US, & + OS%forcing_CSp, dt_forcing=dt_coupling, reset_avg=OS%fluxes%fluxes_used) + endif if (OS%use_ice_shelf) & call add_shelf_forces(OS%grid, OS%US, OS%Ice_shelf_CSp, OS%forces) if (OS%icebergs_alter_ocean) & @@ -516,9 +532,13 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, time_start_upda if (do_thermo) then if (OS%fluxes%fluxes_used) then - call convert_IOB_to_fluxes(Ice_ocean_boundary, OS%fluxes, index_bnds, OS%Time, dt_coupling, & - OS%grid, OS%US, OS%forcing_CSp, OS%sfc_state) - + if (OS%use_homogenized_IOB) then + call convert_IOB_to_fluxes(Ice_ocean_boundary_homogenized, OS%fluxes, index_bnds, OS%Time, dt_coupling, & + OS%grid, OS%US, OS%forcing_CSp, OS%sfc_state) + else + call convert_IOB_to_fluxes(Ice_ocean_boundary, OS%fluxes, index_bnds, OS%Time, dt_coupling, & + OS%grid, OS%US, OS%forcing_CSp, OS%sfc_state) + endif ! Add ice shelf fluxes if (OS%use_ice_shelf) & call shelf_calc_flux(OS%sfc_state, OS%fluxes, OS%Time, dt_coupling, OS%Ice_shelf_CSp) @@ -534,9 +554,15 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, time_start_upda else ! The previous fluxes have not been used yet, so translate the input fluxes ! into a temporary type and then accumulate them in about 20 lines. - OS%flux_tmp%C_p = OS%fluxes%C_p - call convert_IOB_to_fluxes(Ice_ocean_boundary, OS%flux_tmp, index_bnds, OS%Time, dt_coupling, & - OS%grid, OS%US, OS%forcing_CSp, OS%sfc_state) + if (OS%use_homogenized_IOB) then + OS%flux_tmp%C_p = OS%fluxes%C_p + call convert_IOB_to_fluxes(Ice_ocean_boundary_homogenized, OS%flux_tmp, index_bnds, OS%Time, dt_coupling, & + OS%grid, OS%US, OS%forcing_CSp, OS%sfc_state) + else + OS%flux_tmp%C_p = OS%fluxes%C_p + call convert_IOB_to_fluxes(Ice_ocean_boundary, OS%flux_tmp, index_bnds, OS%Time, dt_coupling, & + OS%grid, OS%US, OS%forcing_CSp, OS%sfc_state) + endif if (OS%use_ice_shelf) & call shelf_calc_flux(OS%sfc_state, OS%flux_tmp, OS%Time, dt_coupling, OS%Ice_shelf_CSp) @@ -1194,4 +1220,174 @@ subroutine ocean_model_get_UV_surf(OS, Ocean, name, array2D, isc, jsc) end subroutine ocean_model_get_UV_surf +!> Homogenizes all associated IOB fields. +subroutine homogenize_IOB(IOB,IOB_h, G, index_bounds) + type(ice_ocean_boundary_type), & + target, intent(in) :: IOB !< An ice-ocean boundary type with fluxes to drive + !! the ocean in a coupled model + type(ice_ocean_boundary_type), & + target, intent(out) :: IOB_h !< Similar to IOB, but with homogenized fields. + type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure + integer, dimension(4), intent(in) :: index_bounds !< The i- and j- size of the arrays in IOB. + + integer :: is, ie, js, je, i0, j0, isc_bnd, jsc_bnd, iec_bnd, jec_bnd, i, j + real :: tmp_sum, tmp_mean, N + + is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec + isc_bnd = index_bounds(1) ; iec_bnd = index_bounds(2) + jsc_bnd = index_bounds(3) ; jec_bnd = index_bounds(4) + + i0 = is - isc_bnd ; j0 = js - jsc_bnd + + if (associated(IOB%u_flux)) then + allocate(IOB_h%u_flux (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%u_flux = 0.0 + call homogenize_field(IOB%u_flux,IOB_h%u_flux, G, index_bounds) + endif + + if (associated(IOB%v_flux)) then + allocate(IOB_h%v_flux (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%v_flux = 0.0 + call homogenize_field(IOB%v_flux,IOB_h%v_flux, G, index_bounds) + endif + + if (associated(IOB%t_flux)) then + allocate(IOB_h%t_flux (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%t_flux = 0.0 + call homogenize_field(IOB%t_flux,IOB_h%t_flux, G, index_bounds) + endif + + if (associated(IOB%q_flux)) then + allocate(IOB_h%q_flux (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%q_flux = 0.0 + call homogenize_field(IOB%q_flux,IOB_h%q_flux, G, index_bounds) + endif + + if (associated(IOB%salt_flux)) then + allocate(IOB_h%salt_flux (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%salt_flux = 0.0 + call homogenize_field(IOB%salt_flux,IOB_h%salt_flux, G, index_bounds) + endif + + + if (associated(IOB%lw_flux)) then + allocate(IOB_h%lw_flux (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%lw_flux = 0.0 + call homogenize_field(IOB%lw_flux,IOB_h%lw_flux, G, index_bounds) + endif + + if (associated(IOB%sw_flux_vis_dir)) then + allocate(IOB_h%sw_flux_vis_dir (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%sw_flux_vis_dir = 0.0 + call homogenize_field(IOB%sw_flux_vis_dir,IOB_h%sw_flux_vis_dir, G, index_bounds) + endif + + if (associated(IOB%sw_flux_vis_dif)) then + allocate(IOB_h%sw_flux_vis_dif (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%sw_flux_vis_dif = 0.0 + call homogenize_field(IOB%sw_flux_vis_dif,IOB_h%sw_flux_vis_dif, G, index_bounds) + endif + + if (associated(IOB%sw_flux_nir_dir)) then + allocate(IOB_h%sw_flux_nir_dir (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%sw_flux_nir_dir = 0.0 + call homogenize_field(IOB%sw_flux_nir_dir,IOB_h%sw_flux_nir_dir, G, index_bounds) + endif + + if (associated(IOB%sw_flux_nir_dif)) then + allocate(IOB_h%sw_flux_nir_dif (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%sw_flux_nir_dif = 0.0 + call homogenize_field(IOB%sw_flux_nir_dif,IOB_h%sw_flux_nir_dif, G, index_bounds) + endif + + if (associated(IOB%lprec)) then + allocate(IOB_h%lprec (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%lprec = 0.0 + call homogenize_field(IOB%lprec,IOB_h%lprec, G, index_bounds) + endif + + if (associated(IOB%fprec)) then + allocate(IOB_h%fprec (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%fprec = 0.0 + call homogenize_field(IOB%fprec,IOB_h%fprec, G, index_bounds) + endif + + if (associated(IOB%runoff)) then + allocate(IOB_h%runoff (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%runoff = 0.0 + call homogenize_field(IOB%runoff,IOB_h%runoff, G, index_bounds) + endif + + if (associated(IOB%calving)) then + allocate(IOB_h%calving (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%calving = 0.0 + call homogenize_field(IOB%calving,IOB_h%calving, G, index_bounds) + endif + + if (associated(IOB%stress_mag)) then + allocate(IOB_h%stress_mag (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%stress_mag = 0.0 + call homogenize_field(IOB%stress_mag,IOB_h%stress_mag, G, index_bounds) + endif + + if (associated(IOB%ustar_berg)) then + allocate(IOB_h%ustar_berg (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%ustar_berg = 0.0 + call homogenize_field(IOB%ustar_berg,IOB_h%ustar_berg, G, index_bounds) + endif + + if (associated(IOB%area_berg)) then + allocate(IOB_h%area_berg (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%area_berg = 0.0 + call homogenize_field(IOB%area_berg,IOB_h%area_berg, G, index_bounds) + endif + + if (associated(IOB%mass_berg)) then + allocate(IOB_h%mass_berg (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%mass_berg = 0.0 + call homogenize_field(IOB%mass_berg,IOB_h%mass_berg, G, index_bounds) + endif + + if (associated(IOB%runoff_hflx)) then + allocate(IOB_h%runoff_hflx (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%runoff_hflx = 0.0 + call homogenize_field(IOB%runoff_hflx,IOB_h%runoff_hflx, G, index_bounds) + endif + + if (associated(IOB%calving_hflx)) then + allocate(IOB_h%calving_hflx (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%calving_hflx = 0.0 + call homogenize_field(IOB%calving_hflx,IOB_h%calving_hflx, G, index_bounds) + endif + + if (associated(IOB%p)) then + allocate(IOB_h%p (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%p = 0.0 + call homogenize_field(IOB%p,IOB_h%p, G, index_bounds) + endif + + if (associated(IOB%mi)) then + allocate(IOB_h%mi (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%mi = 0.0 + call homogenize_field(IOB%mi,IOB_h%mi, G, index_bounds) + endif + + if (associated(IOB%ice_rigidity)) then + allocate(IOB_h%ice_rigidity (isc_bnd:iec_bnd,jsc_bnd:jec_bnd) ) ; IOB_h%ice_rigidity = 0.0 + call homogenize_field(IOB%ice_rigidity,IOB_h%ice_rigidity, G, index_bounds) + endif + +end subroutine homogenize_IOB + +!> Homogenize fields for homogenize_IOB +subroutine homogenize_field(field_in,field_out,G, index_bounds) + real, dimension(:,:), intent(in) :: field_in !< Input 2d array meant to be homogenized + real, dimension(:,:), intent(out) :: field_out !< Output 2d array that is homogenized + type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure + integer, dimension(4), intent(in) :: index_bounds !< The i- and j- size of the arrays in IOB. + + real :: tmp_sum, N, tmp_mean + integer :: is, ie, js, je, isc_bnd, jsc_bnd + integer :: i0, j0, i, j + + is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec + isc_bnd = index_bounds(1) + jsc_bnd = index_bounds(3) + + i0 = is - isc_bnd ; j0 = js - jsc_bnd + + tmp_sum = 0.0 ; N = 0.0 + do j=js,je ; do i=is,ie + tmp_sum = tmp_sum + field_in(i-i0,j-j0) * G%mask2dT(i,j) + N = N + G%mask2dT(i,j) + enddo; enddo + + call sum_across_PEs(tmp_sum) + call sum_across_PEs(N) + + tmp_mean = tmp_sum/N + do j=js,je ; do i=is,ie + field_out(i-i0,j-j0) = tmp_mean + enddo; enddo + +end subroutine homogenize_field + end module ocean_model_mod