Skip to content

Commit

Permalink
Prestore axes_data in set_up_ALE_sponge_field
Browse files Browse the repository at this point in the history
fix style errors

Modify axes_data to make it allocatable

fix style
  • Loading branch information
yichengt900 authored and marshallward committed Jun 3, 2024
1 parent 9059671 commit f0badc6
Show file tree
Hide file tree
Showing 2 changed files with 21 additions and 8 deletions.
15 changes: 13 additions & 2 deletions src/framework/MOM_horizontal_regridding.F90
Original file line number Diff line number Diff line change
Expand Up @@ -627,7 +627,8 @@ end subroutine horiz_interp_and_extrap_tracer_record
subroutine horiz_interp_and_extrap_tracer_fms_id(field, Time, G, tr_z, mask_z, &
z_in, z_edges_in, missing_value, scale, &
homogenize, spongeOngrid, m_to_Z, &
answers_2018, tr_iter_tol, answer_date)
answers_2018, tr_iter_tol, answer_date, &
axes)

type(external_field), intent(in) :: field !< Handle for the time interpolated field
type(time_type), intent(in) :: Time !< A FMS time type
Expand Down Expand Up @@ -663,6 +664,7 @@ subroutine horiz_interp_and_extrap_tracer_fms_id(field, Time, G, tr_z, mask_z, &
!! Dates before 20190101 give the same answers
!! as the code did in late 2018, while later versions
!! add parentheses for rotational symmetry.
type(axis_info), allocatable, dimension(:), optional, intent(inout) :: axes !< Axis types for the input data

! Local variables
! In the following comments, [A] is used to indicate the arbitrary, possibly rescaled units of the
Expand Down Expand Up @@ -742,7 +744,16 @@ subroutine horiz_interp_and_extrap_tracer_fms_id(field, Time, G, tr_z, mask_z, &

call cpu_clock_begin(id_clock_read)

call get_external_field_info(field, size=fld_sz, axes=axes_data, missing=missing_val_in)
if (present(axes) .and. allocated(axes)) then
call get_external_field_info(field, size=fld_sz, missing=missing_val_in)
axes_data = axes
else
call get_external_field_info(field, size=fld_sz, axes=axes_data, missing=missing_val_in)
if (present(axes)) then
allocate(axes(4))
axes = axes_data
endif
endif
missing_value = scale*missing_val_in

verbosity = MOM_get_verbosity()
Expand Down
14 changes: 8 additions & 6 deletions src/parameterizations/vertical/MOM_ALE_sponge.F90
Original file line number Diff line number Diff line change
Expand Up @@ -23,6 +23,7 @@ module MOM_ALE_sponge
use MOM_horizontal_regridding, only : horiz_interp_and_extrap_tracer
use MOM_interpolate, only : init_external_field, get_external_field_info, time_interp_external_init
use MOM_interpolate, only : external_field
use MOM_io, only : axis_info
use MOM_remapping, only : remapping_cs, remapping_core_h, initialize_remapping
use MOM_spatial_means, only : global_i_mean
use MOM_time_manager, only : time_type
Expand Down Expand Up @@ -86,6 +87,7 @@ module MOM_ALE_sponge
character(len=:), allocatable :: name !< The name of the input field
character(len=:), allocatable :: long_name !< The long name of the input field
character(len=:), allocatable :: unit !< The unit of the input field
type(axis_info), allocatable :: axes_data(:) !< Axis types for the input field
end type p2d

!> ALE sponge control structure
Expand Down Expand Up @@ -770,7 +772,7 @@ subroutine set_up_ALE_sponge_field_varying(filename, fieldname, Time, G, GV, US,
CS%Ref_val(CS%fldno)%long_name = long_name
CS%Ref_val(CS%fldno)%unit = unit
fld_sz(1:4) = -1
call get_external_field_info(CS%Ref_val(CS%fldno)%field, size=fld_sz)
call get_external_field_info(CS%Ref_val(CS%fldno)%field, size=fld_sz, axes=CS%Ref_val(CS%fldno)%axes_data)
nz_data = fld_sz(3)
CS%Ref_val(CS%fldno)%nz_data = nz_data !< individual sponge fields may reside on a different vertical grid
CS%Ref_val(CS%fldno)%num_tlevs = fld_sz(4)
Expand Down Expand Up @@ -868,7 +870,7 @@ subroutine set_up_ALE_sponge_vel_field_varying(filename_u, fieldname_u, filename
CS%Ref_val_u%field = init_external_field(filename_u, fieldname_u)
endif
fld_sz(1:4) = -1
call get_external_field_info(CS%Ref_val_u%field, size=fld_sz)
call get_external_field_info(CS%Ref_val_u%field, size=fld_sz, axes=CS%Ref_val_u%axes_data)
CS%Ref_val_u%nz_data = fld_sz(3)
CS%Ref_val_u%num_tlevs = fld_sz(4)
CS%Ref_val_u%scale = US%m_s_to_L_T ; if (present(scale)) CS%Ref_val_u%scale = scale
Expand All @@ -879,7 +881,7 @@ subroutine set_up_ALE_sponge_vel_field_varying(filename_u, fieldname_u, filename
CS%Ref_val_v%field = init_external_field(filename_v, fieldname_v)
endif
fld_sz(1:4) = -1
call get_external_field_info(CS%Ref_val_v%field, size=fld_sz)
call get_external_field_info(CS%Ref_val_v%field, size=fld_sz, axes=CS%Ref_val_v%axes_data)
CS%Ref_val_v%nz_data = fld_sz(3)
CS%Ref_val_v%num_tlevs = fld_sz(4)
CS%Ref_val_v%scale = US%m_s_to_L_T ; if (present(scale)) CS%Ref_val_v%scale = scale
Expand Down Expand Up @@ -963,7 +965,7 @@ subroutine apply_ALE_sponge(h, tv, dt, G, GV, US, CS, Time)
call horiz_interp_and_extrap_tracer(CS%Ref_val(m)%field, Time, G, sp_val, &
mask_z, z_in, z_edges_in, missing_value, &
scale=CS%Ref_val(m)%scale, spongeOnGrid=CS%SpongeDataOngrid, m_to_Z=US%m_to_Z, &
answer_date=CS%hor_regrid_answer_date)
answer_date=CS%hor_regrid_answer_date, axes=CS%Ref_val(m)%axes_data)
allocate( dz_src(nz_data) )
allocate( tmpT1d(nz_data) )
do c=1,CS%num_col
Expand Down Expand Up @@ -1053,7 +1055,7 @@ subroutine apply_ALE_sponge(h, tv, dt, G, GV, US, CS, Time)
call horiz_interp_and_extrap_tracer(CS%Ref_val_u%field, Time, G, sp_val, &
mask_z, z_in, z_edges_in, missing_value, &
scale=CS%Ref_val_u%scale, spongeOnGrid=CS%SpongeDataOngrid, m_to_Z=US%m_to_Z, &
answer_date=CS%hor_regrid_answer_date)
answer_date=CS%hor_regrid_answer_date, axes=CS%Ref_val_u%axes_data)

! Initialize mask_z halos to zero before pass_var, in case of no update
mask_z(G%isc-1, G%jsc:G%jec, :) = 0.
Expand Down Expand Up @@ -1101,7 +1103,7 @@ subroutine apply_ALE_sponge(h, tv, dt, G, GV, US, CS, Time)
call horiz_interp_and_extrap_tracer(CS%Ref_val_v%field, Time, G, sp_val, &
mask_z, z_in, z_edges_in, missing_value, &
scale=CS%Ref_val_v%scale, spongeOnGrid=CS%SpongeDataOngrid, m_to_Z=US%m_to_Z,&
answer_date=CS%hor_regrid_answer_date)
answer_date=CS%hor_regrid_answer_date, axes=CS%Ref_val_v%axes_data)
! Initialize mask_z halos to zero before pass_var, in case of no update
mask_z(G%isc:G%iec, G%jsc-1, :) = 0.
mask_z(G%isc:G%iec, G%jec+1, :) = 0.
Expand Down

0 comments on commit f0badc6

Please sign in to comment.