From 0b8dc7beaecf4b3cb33bc36b01e237a3889125d4 Mon Sep 17 00:00:00 2001 From: Alistair Adcroft Date: Tue, 14 Jun 2016 20:01:14 -0400 Subject: [PATCH] New function open_boundary_impose_normal_slope() - initialize_masks() in MOM_grid_initialization.F90 used to read the open-boundary parameters and then adjust the topography on the edge of the domain only. - There is no mention of "Flather" in MOM_grid_initialization.F90 any more. - MOM_fixed_initialization.F90 now calls this new routine prior to calling initialize_masks(). - I had to move a pass_var(G%bathyT) out of initialize_masks() and up to MOM_initialize_fixed() to do this. - This is part of a larger re-factor of OBC code. - No answer changes. --- src/core/MOM_open_boundary.F90 | 35 ++++++++++++ .../MOM_fixed_initialization.F90 | 11 +++- src/initialization/MOM_grid_initialize.F90 | 56 +------------------ 3 files changed, 44 insertions(+), 58 deletions(-) diff --git a/src/core/MOM_open_boundary.F90 b/src/core/MOM_open_boundary.F90 index 027a5b3474..5968bfd924 100644 --- a/src/core/MOM_open_boundary.F90 +++ b/src/core/MOM_open_boundary.F90 @@ -21,6 +21,7 @@ module MOM_open_boundary public open_boundary_init public open_boundary_query public open_boundary_end +public open_boundary_impose_normal_slope public Radiation_Open_Bdry_Conds public set_Flather_positions public set_Flather_data @@ -132,6 +133,14 @@ subroutine open_boundary_config(G, param_file, OBC) "Apply a Flather open boundary condition on the southern\n"//& "side of the global domain", & default=.false.) + + ! Safety check + if ((OBC%apply_OBC_u_flather_west .or. OBC%apply_OBC_v_flather_south) .and. & + .not.G%symmetric ) call MOM_error(FATAL, & + "MOM_open_boundary, open_boundary_config: "//& + "Symmetric memory must be used when APPLY_OBC_U_FLATHER_WEST "//& + "or APPLY_OBC_V_FLATHER_SOUTH is true.") + if (.not.(OBC%apply_OBC_u .or. OBC%apply_OBC_v .or. & OBC%apply_OBC_v_flather_north .or. OBC%apply_OBC_v_flather_south .or. & OBC%apply_OBC_u_flather_east .or. OBC%apply_OBC_u_flather_west)) then @@ -199,6 +208,32 @@ subroutine open_boundary_end(OBC) deallocate(OBC) end subroutine open_boundary_end +!> Sets the slope of bathymetry normal to an open bounndary to zero. +subroutine open_boundary_impose_normal_slope(OBC, G, depth) + type(ocean_OBC_type), pointer :: OBC !< Open boundary control structure + type(ocean_grid_type), intent(inout) :: G !< Ocean grid structure + real, dimension(SZI_(G),SZJ_(G)), intent(inout) :: depth !< Bathymetry at h-points + ! Local variables + integer :: i, j + + if (.not.associated(OBC)) return + + if (associated(OBC%OBC_kind_u)) then + do j=G%jsd,G%jed ; do I=G%isd,G%ied-1 + if (OBC%OBC_kind_u(I,j) == OBC_FLATHER_E) depth(i+1,j) = depth(i,j) + if (OBC%OBC_kind_u(I,j) == OBC_FLATHER_W) depth(i,j) = depth(i+1,j) + enddo ; enddo + endif + + if (associated(OBC%OBC_kind_v)) then + do J=G%jsd,G%jed-1 ; do i=G%isd,G%ied + if (OBC%OBC_kind_v(i,J) == OBC_FLATHER_N) depth(i,j+1) = depth(i,j) + if (OBC%OBC_kind_v(i,J) == OBC_FLATHER_S) depth(i,j) = depth(i,j+1) + enddo ; enddo + endif + +end subroutine open_boundary_impose_normal_slope + !> Diagnose radiation conditions at open boundaries subroutine Radiation_Open_Bdry_Conds(OBC, u_new, u_old, v_new, v_old, & h_new, h_old, G) diff --git a/src/initialization/MOM_fixed_initialization.F90 b/src/initialization/MOM_fixed_initialization.F90 index e9cd32e0ad..fd31128063 100644 --- a/src/initialization/MOM_fixed_initialization.F90 +++ b/src/initialization/MOM_fixed_initialization.F90 @@ -20,7 +20,7 @@ module MOM_fixed_initialization use MOM_grid_initialize, only : initialize_masks, set_grid_metrics use MOM_open_boundary, only : ocean_OBC_type use MOM_open_boundary, only : open_boundary_config, open_boundary_query -use MOM_open_boundary, only : set_Flather_positions +use MOM_open_boundary, only : set_Flather_positions, open_boundary_impose_normal_slope use MOM_string_functions, only : uppercase use user_initialization, only : user_initialize_topography, USER_set_OBC_positions use DOME_initialization, only : DOME_initialize_topography, DOME_set_OBC_positions @@ -104,8 +104,13 @@ subroutine MOM_initialize_fixed(G, OBC, PF, write_geom, output_dir) call set_Flather_positions(G, OBC) endif -! This call sets seamasks that prohibit flow over any point with ! -! a bottom that is shallower than min_depth from PF. ! + ! To initialize masks, the bathymetry in halo regions must be filled in + call pass_var(G%bathyT, G%Domain) + + ! Make bathymetry consistent with open boundaries + call open_boundary_impose_normal_slope(OBC, G, G%bathyT) + + ! This call sets masks that prohibit flow over any point interpreted as land call initialize_masks(G, PF) if (debug) then call hchksum(G%bathyT, 'MOM_initialize_fixed: depth ', G%HI, haloshift=1) diff --git a/src/initialization/MOM_grid_initialize.F90 b/src/initialization/MOM_grid_initialize.F90 index 03332babf1..6bb268cd3a 100644 --- a/src/initialization/MOM_grid_initialize.F90 +++ b/src/initialization/MOM_grid_initialize.F90 @@ -1343,8 +1343,6 @@ subroutine initialize_masks(G, PF) real :: Dmin, min_depth, mask_depth integer :: i, j - logical :: apply_OBC_u_flather_east, apply_OBC_u_flather_west - logical :: apply_OBC_v_flather_north, apply_OBC_v_flather_south character(len=40) :: mod = "MOM_grid_init initialize_masks" call callTree_enter("initialize_masks(), MOM_grid_initialize.F90") @@ -1358,65 +1356,13 @@ subroutine initialize_masks(G, PF) "The depth below which to mask points as land points, for which all\n"//& "fluxes are zeroed out. MASKING_DEPTH is ignored if negative.", & units="m", default=-9999.0) - call get_param(PF, mod, "APPLY_OBC_U_FLATHER_EAST", apply_OBC_u_flather_east,& - "Apply a Flather open boundary condition on the eastern \n"//& - "side of the global domain", default=.false.) - call get_param(PF, mod, "APPLY_OBC_U_FLATHER_WEST", apply_OBC_u_flather_west,& - "Apply a Flather open boundary condition on the western \n"//& - "side of the global domain", default=.false.) - call get_param(PF, mod, "APPLY_OBC_V_FLATHER_NORTH", apply_OBC_v_flather_north,& - "Apply a Flather open boundary condition on the northern \n"//& - "side of the global domain", default=.false.) - call get_param(PF, mod, "APPLY_OBC_V_FLATHER_SOUTH", apply_OBC_v_flather_south,& - "Apply a Flather open boundary condition on the southern \n"//& - "side of the global domain", default=.false.) - - if ((apply_OBC_u_flather_west .or. apply_OBC_v_flather_south) .and. & - .not.G%symmetric ) & - call MOM_error(FATAL, "Symmetric memory must be used when "//& - "APPLY_OBC_U_FLATHER_WEST or APPLY_OBC_V_FLATHER_SOUTH is true.") Dmin = min_depth if (mask_depth>=0.) Dmin = mask_depth - call pass_var(G%bathyT, G%Domain) G%mask2dCu(:,:) = 0.0 ; G%mask2dCv(:,:) = 0.0 ; G%mask2dBu(:,:) = 0.0 - ! Extrapolate the bottom depths at any points that are subject to Flather - ! open boundary conditions. This should be generalized for Flather OBCs - ! that are not necessarily at the edges of the domain. - if (apply_OBC_u_flather_west) then - do j=G%jsd,G%jed ; do I=G%isd+1,G%ied - if ((I+G%idg_offset) == G%isg) then - G%bathyT(i-1,j) = G%bathyT(i,j) - endif - enddo; enddo - endif - - if (apply_OBC_u_flather_east) then - do j=G%jsd,G%jed ; do I=G%isd,G%ied-1 - if ((i+G%idg_offset) == G%ieg) then - G%bathyT(i+1,j) = G%bathyT(i,j) - endif - enddo; enddo - endif - - if (apply_OBC_v_flather_north) then - do J=G%jsd,G%jed-1 ; do i=G%isd,G%ied - if ((j+G%jdg_offset) == G%jeg) then - G%bathyT(i,j+1) = G%bathyT(i,j) - endif - enddo; enddo - endif - - if (apply_OBC_v_flather_south) then - do J=G%jsd+1,G%jed ; do i=G%isd,G%ied - if ((J+G%jdg_offset) == G%jsg) then - G%bathyT(i,j-1) = G%bathyT(i,j) - endif - enddo; enddo - endif - + ! Construct the h-point or T-point mask do j=G%jsd,G%jed ; do i=G%isd,G%ied if (G%bathyT(i,j) <= Dmin) then G%mask2dT(i,j) = 0.0