Skip to content

Commit

Permalink
+Add ALE_remap_interface_vals and ALE_remap_vertex_vals
Browse files Browse the repository at this point in the history
  Added ALE_remap_interface_vals and ALE_remap_vertex_vals to handle the
interpolation of variables that are at the interfaces atop tracer cells or above
the corners of the tracers cells from one grid to another.  Because these are
not yet used (but have been tested in calls that will be added with the next
commit, all answers are bitwise identical, but there are two new publicly
visible routines.
  • Loading branch information
Hallberg-NOAA committed Nov 9, 2022
1 parent 5cb7400 commit 27c13b9
Showing 1 changed file with 85 additions and 0 deletions.
85 changes: 85 additions & 0 deletions src/ALE/MOM_ALE.F90
Original file line number Diff line number Diff line change
Expand Up @@ -129,6 +129,8 @@ module MOM_ALE
public ALE_remap_scalar
public ALE_remap_tracers
public ALE_remap_velocities
public ALE_remap_interface_vals
public ALE_remap_vertex_vals
public ALE_PLM_edge_values
public TS_PLM_edge_values
public TS_PPM_edge_values
Expand Down Expand Up @@ -974,6 +976,89 @@ subroutine ALE_remap_velocities(CS, G, GV, h_old, h_new, u, v, OBC, dzInterface,

end subroutine ALE_remap_velocities

!> Interpolate to find an updated array of values at interfaces after remapping.
subroutine ALE_remap_interface_vals(CS, G, GV, h_old, h_new, int_val)
type(ALE_CS), intent(in) :: CS !< ALE control structure
type(ocean_grid_type), intent(in) :: G !< Ocean grid structure
type(verticalGrid_type), intent(in) :: GV !< Ocean vertical grid structure
real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h_old !< Thickness of source grid
!! [H ~> m or kg m-2]
real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h_new !< Thickness of destination grid
!! [H ~> m or kg m-2]
real, dimension(SZI_(G),SZJ_(G),SZK_(GV)+1), &
intent(inout) :: int_val !< The interface values to interpolate [A]

real :: val_src(GV%ke) ! A column of interface values on the source grid [A]
real :: val_tgt(GV%ke) ! A column of interface values on the target grid [A]
real :: h_src(GV%ke) ! A column of source grid layer thicknesses [H ~> m or kg m-2]
real :: h_tgt(GV%ke) ! A column of target grid layer thicknesses [H ~> m or kg m-2]
real :: h_neglect, h_neglect_edge ! Tiny thicknesses used in remapping [H ~> m or kg m-2]
logical :: show_call_tree
integer :: i, j, k, nz

do j=G%jsc,G%jec ; do i=G%isc,G%iec ; if (G%mask2dT(i,j)>0.) then
do k=1,nz
h_src(k) = h_old(i,j,k)
h_tgt(k) = h_new(i,j,k)
enddo

do K=1,nz+1
val_src(K) = int_val(i,j,K)
enddo

call interpolate_column(nz, h_src, val_src, nz, h_tgt, val_tgt, .false.)

do K=1,nz+1
int_val(i,j,K) = val_tgt(K)
enddo
endif ; enddo ; enddo

end subroutine ALE_remap_interface_vals

!> Interpolate to find an updated array of values at vertices of tracer cells after remapping.
subroutine ALE_remap_vertex_vals(CS, G, GV, h_old, h_new, vert_val)
type(ALE_CS), intent(in) :: CS !< ALE control structure
type(ocean_grid_type), intent(in) :: G !< Ocean grid structure
type(verticalGrid_type), intent(in) :: GV !< Ocean vertical grid structure
real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h_old !< Thickness of source grid
!! [H ~> m or kg m-2]
real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h_new !< Thickness of destination grid
!! [H ~> m or kg m-2]
real, dimension(SZIB_(G),SZJB_(G),SZK_(GV)+1), &
intent(inout) :: vert_val !< The interface values to interpolate [A]

real :: val_src(GV%ke) ! A column of interface values on the source grid [A]
real :: val_tgt(GV%ke) ! A column of interface values on the target grid [A]
real :: h_src(GV%ke) ! A column of source grid layer thicknesses [H ~> m or kg m-2]
real :: h_tgt(GV%ke) ! A column of target grid layer thicknesses [H ~> m or kg m-2]
real :: I_mask_sum ! The inverse of the tracer point masks surrounding a corner [nondim]
real :: h_neglect, h_neglect_edge ! Tiny thicknesses used in remapping [H ~> m or kg m-2]
logical :: show_call_tree
integer :: i, j, k, nz

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)))

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
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
enddo

do K=1,nz+1
val_src(K) = vert_val(I,J,K)
enddo

call interpolate_column(nz, h_src, val_src, nz, h_tgt, val_tgt, .false.)

do K=1,nz+1
vert_val(I,J,K) = val_tgt(K)
enddo
endif ; enddo ; enddo

end subroutine ALE_remap_vertex_vals

!> Mask out thicknesses to 0 when their running sum exceeds a specified value.
subroutine apply_partial_cell_mask(h1, h_mask)
Expand Down

0 comments on commit 27c13b9

Please sign in to comment.