Skip to content

Add HYBRID_3D and HYBRID_MAP options to ALE_COORDINATE_CONFIG #965

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
Hallberg-NOAA merged 5 commits into from
Sep 8, 2025
Merged

Add HYBRID_3D and HYBRID_MAP options to ALE_COORDINATE_CONFIG #965

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
b6dbc37
Add ALE_COORDINATE_CONFIG HYBRID_3D and HYBRID_MAP
awallcraft Jun 12, 2025
9f77d55
Set the halos of new 3D arrays
awallcraft Jun 15, 2025
2b87ecf
Add SHALLOW_ALE_RESOLUTION
awallcraft Jun 17, 2025
3294a2d
Non-integer HYBRID_MAP values
awallcraft Aug 27, 2025
0deb6ce
indent continuations, source code <= 100 chars
awallcraft Sep 6, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
55 changes: 32 additions & 23 deletions src/ALE/MOM_ALE.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -165,8 +165,9 @@ module MOM_ALE
!! before the main time integration loop to initialize the regridding stuff.
!! We read the MOM_input file to register the values of different
!! regridding/remapping parameters.
subroutine ALE_init( param_file, GV, US, max_depth, CS)
subroutine ALE_init( param_file, G, GV, US, max_depth, CS)
type(param_file_type), intent(in) :: param_file !< Parameter file
type(ocean_grid_type), intent(in) :: G !< Grid structure
type(verticalGrid_type), intent(in) :: GV !< Ocean vertical grid structure
type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
real, intent(in) :: max_depth !< The maximum depth of the ocean [Z ~> m].
Expand Down Expand Up @@ -205,8 +206,9 @@ subroutine ALE_init( param_file, GV, US, max_depth, CS)
default=.false.)

! Initialize and configure regridding
call ALE_initRegridding(GV, US, max_depth, param_file, mdl, CS%regridCS)
call regridding_preadjust_reqs(CS%regridCS, CS%do_conv_adj, CS%use_hybgen_unmix, hybgen_CS=hybgen_regridCS)
call ALE_initRegridding(G, GV, US, max_depth, param_file, mdl, CS%regridCS)
call regridding_preadjust_reqs(CS%regridCS, CS%do_conv_adj, CS%use_hybgen_unmix, &
hybgen_CS=hybgen_regridCS)

! Initialize and configure remapping that is orchestrated by ALE.
call get_param(param_file, mdl, "REMAPPING_SCHEME", string, &
Expand Down Expand Up @@ -321,12 +323,12 @@ subroutine ALE_init( param_file, GV, US, max_depth, CS)
default=-0.001, units="m", scale=GV%m_to_H)
call get_param(param_file, mdl, "REMAP_VEL_MASK_H_THIN", CS%h_vel_mask, &
"A thickness at velocity points below which near-bottom layers are zeroed out "//&
"after remapping, following practice with Hybgen remapping, or a negative value "//&
"to avoid such filtering altogether.", &
"after remapping, following practice with Hybgen remapping, "//&
"or a negative value to avoid such filtering altogether.", &
default=1.0e-6, units="m", scale=GV%m_to_H, do_not_log=(CS%BBL_h_vel_mask<=0.0))

if (CS%use_hybgen_unmix) &
call init_hybgen_unmix(CS%hybgen_unmixCS, GV, US, param_file, hybgen_regridCS)
call init_hybgen_unmix(CS%hybgen_unmixCS, GV, US, param_file, hybgen_regridCS)

call get_param(param_file, mdl, "REMAP_VEL_CONSERVE_KE", CS%conserve_ke, &
"If true, a correction is applied to the baroclinic component of velocity "//&
Expand Down Expand Up @@ -640,7 +642,8 @@ end subroutine ALE_offline_inputs

!> For a state-based coordinate, accelerate the process of regridding by
!! repeatedly applying the grid calculation algorithm
subroutine ALE_regrid_accelerated(CS, G, GV, US, h, tv, n_itt, u, v, OBC, Reg, dt, dzRegrid, initial)
subroutine ALE_regrid_accelerated(CS, G, GV, US, h, tv, n_itt, u, v, OBC, Reg, dt, &
dzRegrid, initial)
type(ALE_CS), pointer :: CS !< ALE control structure
type(ocean_grid_type), intent(inout) :: G !< Ocean grid
type(verticalGrid_type), intent(in) :: GV !< Vertical grid
Expand Down Expand Up @@ -689,7 +692,7 @@ subroutine ALE_regrid_accelerated(CS, G, GV, US, h, tv, n_itt, u, v, OBC, Reg, d

! initial total interface displacement due to successive regridding
if (CS%remap_uv_using_old_alg) &
dzIntTotal(:,:,:) = 0.
dzIntTotal(:,:,:) = 0.

call create_group_pass(pass_T_S_h, T, G%domain)
call create_group_pass(pass_T_S_h, S, G%domain)
Expand All @@ -708,7 +711,7 @@ subroutine ALE_regrid_accelerated(CS, G, GV, US, h, tv, n_itt, u, v, OBC, Reg, d

! Apply timescale to regridding (for e.g. filtered_grid_motion)
if (present(dt)) &
call ALE_update_regrid_weights(dt, CS)
call ALE_update_regrid_weights(dt, CS)

do itt = 1, n_itt

Expand All @@ -722,12 +725,14 @@ subroutine ALE_regrid_accelerated(CS, G, GV, US, h, tv, n_itt, u, v, OBC, Reg, d

call regridding_main(CS%remapCS, CS%regridCS, G, GV, US, h_loc, tv_local, h, dzInterface)
if (CS%remap_uv_using_old_alg) &
dzIntTotal(:,:,:) = dzIntTotal(:,:,:) + dzInterface(:,:,:)
dzIntTotal(:,:,:) = dzIntTotal(:,:,:) + dzInterface(:,:,:)

! remap from original grid onto new grid
do j = G%jsc-1,G%jec+1 ; do i = G%isc-1,G%iec+1
call remapping_core_h(CS%remapCS, nz, h_orig(i,j,:), tv%S(i,j,:), nz, h(i,j,:), tv_local%S(i,j,:))
call remapping_core_h(CS%remapCS, nz, h_orig(i,j,:), tv%T(i,j,:), nz, h(i,j,:), tv_local%T(i,j,:))
call remapping_core_h(CS%remapCS, nz, h_orig(i,j,:), tv%S(i,j,:), nz, h(i,j,:), &
tv_local%S(i,j,:))
call remapping_core_h(CS%remapCS, nz, h_orig(i,j,:), tv%T(i,j,:), nz, h(i,j,:), &
tv_local%T(i,j,:))
enddo ; enddo

! starting grid for next iteration
Expand Down Expand Up @@ -1146,7 +1151,7 @@ subroutine ALE_remap_velocities(CS, G, GV, h_old_u, h_old_v, h_new_u, h_new_v, u
if (CS%id_remap_delta_integ_v2>0) dv2h_tot(:,:) = 0.

if (((CS%id_remap_delta_integ_u2>0) .or. (CS%id_remap_delta_integ_v2>0)) .and. .not.present(dt))&
call MOM_error(FATAL, "ALE KE diagnostics requires passing dt into ALE_remap_velocities")
call MOM_error(FATAL, "ALE KE diagnostics requires passing dt into ALE_remap_velocities")

nz = GV%ke

Expand Down Expand Up @@ -1212,7 +1217,7 @@ subroutine ALE_remap_velocities(CS, G, GV, h_old_u, h_old_v, h_new_u, h_new_v, u
endif

if ((CS%BBL_h_vel_mask > 0.0) .and. (CS%h_vel_mask > 0.0)) &
call mask_near_bottom_vel(u_tgt, h2, CS%BBL_h_vel_mask, CS%h_vel_mask, nz)
call mask_near_bottom_vel(u_tgt, h2, CS%BBL_h_vel_mask, CS%h_vel_mask, nz)

! Copy the column of new velocities back to the 3-d array
do k=1,nz
Expand Down Expand Up @@ -1361,13 +1366,14 @@ subroutine ALE_remap_vertex_vals(CS, G, GV, h_old, h_new, vert_val)

do J=G%JscB,G%JecB ; do I=G%IscB,G%IecB
if ((G%mask2dT(i,j) + G%mask2dT(i+1,j+1)) + (G%mask2dT(i+1,j) + G%mask2dT(i,j+1)) > 0.0 ) then
I_mask_sum = 1.0 / ((G%mask2dT(i,j) + G%mask2dT(i+1,j+1)) + (G%mask2dT(i+1,j) + G%mask2dT(i,j+1)))
I_mask_sum = 1.0 / ((G%mask2dT(i,j) + G%mask2dT(i+1,j+1)) + &
(G%mask2dT(i+1,j) + G%mask2dT(i,j+1)))

do k=1,nz
h_src(k) = ((G%mask2dT(i,j) * h_old(i,j,k) + G%mask2dT(i+1,j+1) * h_old(i+1,j+1,k)) + &
(G%mask2dT(i+1,j) * h_old(i+1,j,k) + G%mask2dT(i,j+1) * h_old(i,j+1,k)) ) * I_mask_sum
(G%mask2dT(i+1,j) * h_old(i+1,j,k) + G%mask2dT(i,j+1) * h_old(i,j+1,k)) ) * I_mask_sum
h_tgt(k) = ((G%mask2dT(i,j) * h_new(i,j,k) + G%mask2dT(i+1,j+1) * h_new(i+1,j+1,k)) + &
(G%mask2dT(i+1,j) * h_new(i+1,j,k) + G%mask2dT(i,j+1) * h_new(i,j+1,k)) ) * I_mask_sum
(G%mask2dT(i+1,j) * h_new(i+1,j,k) + G%mask2dT(i,j+1) * h_new(i,j+1,k)) ) * I_mask_sum
enddo

do K=1,nz+1
Expand Down Expand Up @@ -1549,7 +1555,8 @@ subroutine ALE_PLM_edge_values( CS, G, GV, h, Q, bdry_extrap, Q_t, Q_b )
do j = G%jsc-1,G%jec+1 ; do i = G%isc-1,G%iec+1
slp(1) = 0.
do k = 2, GV%ke-1
slp(k) = PLM_slope_wa(h(i,j,k-1), h(i,j,k), h(i,j,k+1), h_neglect, Q(i,j,k-1), Q(i,j,k), Q(i,j,k+1))
slp(k) = PLM_slope_wa(h(i,j,k-1), h(i,j,k), h(i,j,k+1), h_neglect, &
Q(i,j,k-1), Q(i,j,k), Q(i,j,k+1))
enddo
slp(GV%ke) = 0.

Expand All @@ -1562,7 +1569,8 @@ subroutine ALE_PLM_edge_values( CS, G, GV, h, Q, bdry_extrap, Q_t, Q_b )
mslp = - PLM_extrapolate_slope(h(i,j,2), h(i,j,1), h_neglect, Q(i,j,2), Q(i,j,1))
Q_t(i,j,1) = Q(i,j,1) - 0.5 * mslp
Q_b(i,j,1) = Q(i,j,1) + 0.5 * mslp
mslp = PLM_extrapolate_slope(h(i,j,GV%ke-1), h(i,j,GV%ke), h_neglect, Q(i,j,GV%ke-1), Q(i,j,GV%ke))
mslp = PLM_extrapolate_slope(h(i,j,GV%ke-1), h(i,j,GV%ke), h_neglect, &
Q(i,j,GV%ke-1), Q(i,j,GV%ke))
Q_t(i,j,GV%ke) = Q(i,j,GV%ke) - 0.5 * mslp
Q_b(i,j,GV%ke) = Q(i,j,GV%ke) + 0.5 * mslp
else
Expand Down Expand Up @@ -1630,7 +1638,7 @@ subroutine TS_PPM_edge_values( CS, S_t, S_b, T_t, T_b, G, GV, tv, h, bdry_extrap
call PPM_reconstruction( GV%ke, hTmp, tmp, ppol_E, ppol_coefs, h_neglect, &
answer_date=CS%answer_date )
if (bdry_extrap) &
call PPM_boundary_extrapolation( GV%ke, hTmp, tmp, ppol_E, ppol_coefs, h_neglect )
call PPM_boundary_extrapolation( GV%ke, hTmp, tmp, ppol_E, ppol_coefs, h_neglect )

do k = 1,GV%ke
S_t(i,j,k) = ppol_E(k,1)
Expand All @@ -1651,7 +1659,7 @@ subroutine TS_PPM_edge_values( CS, S_t, S_b, T_t, T_b, G, GV, tv, h, bdry_extrap
call PPM_reconstruction( GV%ke, hTmp, tmp, ppol_E, ppol_coefs, h_neglect, &
answer_date=CS%answer_date )
if (bdry_extrap) &
call PPM_boundary_extrapolation(GV%ke, hTmp, tmp, ppol_E, ppol_coefs, h_neglect )
call PPM_boundary_extrapolation(GV%ke, hTmp, tmp, ppol_E, ppol_coefs, h_neglect )

do k = 1,GV%ke
T_t(i,j,k) = ppol_E(k,1)
Expand All @@ -1664,7 +1672,8 @@ end subroutine TS_PPM_edge_values


!> Initializes regridding for the main ALE algorithm
subroutine ALE_initRegridding(GV, US, max_depth, param_file, mdl, regridCS)
subroutine ALE_initRegridding(G, GV, US, max_depth, param_file, mdl, regridCS)
type(ocean_grid_type), intent(in) :: G !< Grid structure
type(verticalGrid_type), intent(in) :: GV !< Ocean vertical grid structure
type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
real, intent(in) :: max_depth !< The maximum depth of the ocean [Z ~> m].
Expand All @@ -1680,7 +1689,7 @@ subroutine ALE_initRegridding(GV, US, max_depth, param_file, mdl, regridCS)
trim(regriddingCoordinateModeDoc), &
default=DEFAULT_COORDINATE_MODE, fail_if_missing=.true.)

call initialize_regridding(regridCS, GV, US, max_depth, param_file, mdl, coord_mode, '', '')
call initialize_regridding(regridCS, G, GV, US, max_depth, param_file, mdl, coord_mode, '', '')

end subroutine ALE_initRegridding

Expand Down
Loading