diff --git a/src/diagnostics/MOM_sum_output.F90 b/src/diagnostics/MOM_sum_output.F90 index 03204e4322..8b24463af9 100644 --- a/src/diagnostics/MOM_sum_output.F90 +++ b/src/diagnostics/MOM_sum_output.F90 @@ -4,7 +4,7 @@ module MOM_sum_output ! This file is part of MOM6. See LICENSE.md for the license. use iso_fortran_env, only : int64 -use MOM_coms, only : sum_across_PEs, PE_here, root_PE, num_PEs, max_across_PEs +use MOM_coms, only : sum_across_PEs, PE_here, root_PE, num_PEs, max_across_PEs, field_chksum use MOM_coms, only : reproducing_sum, reproducing_sum_EFP, EFP_to_real, real_to_EFP use MOM_coms, only : EFP_type, operator(+), operator(-), assignment(=), EFP_sum_across_PEs use MOM_error_handler, only : MOM_error, FATAL, WARNING, is_root_pe, MOM_mesg @@ -25,7 +25,6 @@ module MOM_sum_output use MOM_unit_scaling, only : unit_scale_type use MOM_variables, only : surface, thermo_var_ptrs use MOM_verticalGrid, only : verticalGrid_type -use mpp_mod, only : mpp_chksum use netcdf @@ -1511,13 +1510,13 @@ subroutine get_depth_list_checksums(G, depth_chksum, area_chksum) do j=G%jsc,G%jec ; do i=G%isc,G%iec field(i,j) = G%bathyT(i,j) enddo ; enddo - write(depth_chksum, '(Z16)') mpp_chksum(field(:,:)) + write(depth_chksum, '(Z16)') field_chksum(field(:,:)) ! Area checksum do j=G%jsc,G%jec ; do i=G%isc,G%iec field(i,j) = G%mask2dT(i,j) * G%US%L_to_m**2*G%areaT(i,j) enddo ; enddo - write(area_chksum, '(Z16)') mpp_chksum(field(:,:)) + write(area_chksum, '(Z16)') field_chksum(field(:,:)) deallocate(field) end subroutine get_depth_list_checksums diff --git a/src/framework/MOM_checksums.F90 b/src/framework/MOM_checksums.F90 index c3174dbe7b..5c503836f0 100644 --- a/src/framework/MOM_checksums.F90 +++ b/src/framework/MOM_checksums.F90 @@ -3,10 +3,11 @@ module MOM_checksums ! This file is part of MOM6. See LICENSE.md for the license. -use MOM_array_transform, only: rotate_array, rotate_array_pair, rotate_vector +use MOM_array_transform, only : rotate_array, rotate_array_pair, rotate_vector +use MOM_array_transform, only : allocate_rotated_array use MOM_coms, only : PE_here, root_PE, num_PEs, sum_across_PEs use MOM_coms, only : min_across_PEs, max_across_PEs -use MOM_coms, only : reproducing_sum +use MOM_coms, only : reproducing_sum, field_chksum use MOM_error_handler, only : MOM_error, FATAL, is_root_pe use MOM_file_parser, only : log_version, param_file_type use MOM_hor_index, only : hor_index_type, rotate_hor_index @@ -15,7 +16,7 @@ module MOM_checksums implicit none ; private -public :: chksum0, zchksum +public :: chksum0, zchksum, rotated_field_chksum public :: hchksum, Bchksum, uchksum, vchksum, qchksum, is_NaN, chksum public :: hchksum_pair, uvchksum, Bchksum_pair public :: MOM_checksums_init @@ -75,6 +76,15 @@ module MOM_checksums module procedure is_NaN_0d, is_NaN_1d, is_NaN_2d, is_NaN_3d end interface +!> Rotate and compute the checksum of a field +interface rotated_field_chksum + module procedure rotated_field_chksum_real_0d + module procedure rotated_field_chksum_real_1d + module procedure rotated_field_chksum_real_2d + module procedure rotated_field_chksum_real_3d + module procedure rotated_field_chksum_real_4d +end interface rotated_field_chksum + integer, parameter :: bc_modulus = 1000000000 !< Modulus of checksum bitcount integer, parameter :: default_shift=0 !< The default array shift logical :: calculateStatistics=.true. !< If true, report min, max and mean. @@ -2021,16 +2031,16 @@ function is_NaN_1d(x, skip_mpp) logical :: is_NaN_1d integer :: i, n - logical :: call_mpp + logical :: global_check n = 0 do i = LBOUND(x,1), UBOUND(x,1) if (is_NaN_0d(x(i))) n = n + 1 enddo - call_mpp = .true. - if (present(skip_mpp)) call_mpp = .not.skip_mpp + global_check = .true. + if (present(skip_mpp)) global_check = .not.skip_mpp - if (call_mpp) call sum_across_PEs(n) + if (global_check) call sum_across_PEs(n) is_NaN_1d = .false. if (n>0) is_NaN_1d = .true. @@ -2072,6 +2082,121 @@ function is_NaN_3d(x) end function is_NaN_3d +! The following set of routines do a checksum across the computational domain of +! a field, with the potential for rotation of this field and masking. + +!> Compute the field checksum of a scalar. +function rotated_field_chksum_real_0d(field, pelist, mask_val, turns) & + result(chksum) + real, intent(in) :: field !< Input scalar + integer, optional, intent(in) :: pelist(:) !< PE list of ranks to checksum + real, optional, intent(in) :: mask_val !< FMS mask value + integer, optional, intent(in) :: turns !< Number of quarter turns + integer :: chksum !< checksum of scalar + + if (present(turns)) call MOM_error(FATAL, "Rotation not supported for 0d fields.") + + chksum = field_chksum(field, pelist=pelist, mask_val=mask_val) +end function rotated_field_chksum_real_0d + + +!> Compute the field checksum of a 1d field. +function rotated_field_chksum_real_1d(field, pelist, mask_val, turns) & + result(chksum) + real, dimension(:), intent(in) :: field !< Input array + integer, optional, intent(in) :: pelist(:) !< PE list of ranks to checksum + real, optional, intent(in) :: mask_val !< FMS mask value + integer, optional, intent(in) :: turns !< Number of quarter turns + integer :: chksum !< checksum of array + + if (present(turns)) call MOM_error(FATAL, "Rotation not supported for 1d fields.") + + chksum = field_chksum(field, pelist=pelist, mask_val=mask_val) +end function rotated_field_chksum_real_1d + + +!> Compute the field checksum of a rotated 2d field. +function rotated_field_chksum_real_2d(field, pelist, mask_val, turns) & + result(chksum) + real, dimension(:,:), intent(in) :: field !< Unrotated input field + integer, optional, intent(in) :: pelist(:) !< PE list of ranks to checksum + real, optional, intent(in) :: mask_val !< FMS mask value + integer, optional, intent(in) :: turns !< Number of quarter turns + integer :: chksum !< checksum of array + + ! Local variables + real, allocatable :: field_rot(:,:) ! A rotated version of field, with the same units + integer :: qturns ! The number of quarter turns through which to rotate field + + qturns = 0 + if (present(turns)) & + qturns = modulo(turns, 4) + + if (qturns == 0) then + chksum = field_chksum(field, pelist=pelist, mask_val=mask_val) + else + call allocate_rotated_array(field, [1,1], qturns, field_rot) + call rotate_array(field, qturns, field_rot) + chksum = field_chksum(field_rot, pelist=pelist, mask_val=mask_val) + deallocate(field_rot) + endif +end function rotated_field_chksum_real_2d + +!> Compute the field checksum of a rotated 3d field. +function rotated_field_chksum_real_3d(field, pelist, mask_val, turns) & + result(chksum) + real, dimension(:,:,:), intent(in) :: field !< Unrotated input field + integer, optional, intent(in) :: pelist(:) !< PE list of ranks to checksum + real, optional, intent(in) :: mask_val !< FMS mask value + integer, optional, intent(in) :: turns !< Number of quarter turns + integer :: chksum !< checksum of array + + ! Local variables + real, allocatable :: field_rot(:,:,:) ! A rotated version of field, with the same units + integer :: qturns ! The number of quarter turns through which to rotate field + + qturns = 0 + if (present(turns)) & + qturns = modulo(turns, 4) + + if (qturns == 0) then + chksum = field_chksum(field, pelist=pelist, mask_val=mask_val) + else + call allocate_rotated_array(field, [1,1,1], qturns, field_rot) + call rotate_array(field, qturns, field_rot) + chksum = field_chksum(field_rot, pelist=pelist, mask_val=mask_val) + deallocate(field_rot) + endif +end function rotated_field_chksum_real_3d + +!> Compute the field checksum of a rotated 4d field. +function rotated_field_chksum_real_4d(field, pelist, mask_val, turns) & + result(chksum) + real, dimension(:,:,:,:), intent(in) :: field !< Unrotated input field + integer, optional, intent(in) :: pelist(:) !< PE list of ranks to checksum + real, optional, intent(in) :: mask_val !< FMS mask value + integer, optional, intent(in) :: turns !< Number of quarter turns + integer :: chksum !< checksum of array + + ! Local variables + real, allocatable :: field_rot(:,:,:,:) ! A rotated version of field, with the same units + integer :: qturns ! The number of quarter turns through which to rotate field + + qturns = 0 + if (present(turns)) & + qturns = modulo(turns, 4) + + if (qturns == 0) then + chksum = field_chksum(field, pelist=pelist, mask_val=mask_val) + else + call allocate_rotated_array(field, [1,1,1,1], qturns, field_rot) + call rotate_array(field, qturns, field_rot) + chksum = field_chksum(field_rot, pelist=pelist, mask_val=mask_val) + deallocate(field_rot) + endif +end function rotated_field_chksum_real_4d + + !> Write a message including the checksum of the non-shifted array subroutine chk_sum_msg1(fmsg, bc0, mesg, iounit) character(len=*), intent(in) :: fmsg !< A checksum code-location specific preamble diff --git a/src/framework/MOM_coms.F90 b/src/framework/MOM_coms.F90 index 0c6b948980..04ed46ad22 100644 --- a/src/framework/MOM_coms.F90 +++ b/src/framework/MOM_coms.F90 @@ -9,13 +9,13 @@ module MOM_coms use memutils_mod, only : print_memuse_stats use mpp_mod, only : PE_here => mpp_pe, root_PE => mpp_root_pe, num_PEs => mpp_npes use mpp_mod, only : Set_PElist => mpp_set_current_pelist, Get_PElist => mpp_get_current_pelist -use mpp_mod, only : broadcast => mpp_broadcast +use mpp_mod, only : broadcast => mpp_broadcast, field_chksum => mpp_chksum use mpp_mod, only : sum_across_PEs => mpp_sum, max_across_PEs => mpp_max, min_across_PEs => mpp_min implicit none ; private public :: PE_here, root_PE, num_PEs, MOM_infra_init, MOM_infra_end -public :: broadcast, sum_across_PEs, min_across_PEs, max_across_PEs +public :: broadcast, sum_across_PEs, min_across_PEs, max_across_PEs, field_chksum public :: reproducing_sum, reproducing_sum_EFP, EFP_sum_across_PEs, EFP_list_sum_across_PEs public :: EFP_plus, EFP_minus, EFP_to_real, real_to_EFP, EFP_real_diff public :: operator(+), operator(-), assignment(=) diff --git a/src/framework/MOM_domain_init.F90 b/src/framework/MOM_domain_init.F90 new file mode 100644 index 0000000000..25064cf24e --- /dev/null +++ b/src/framework/MOM_domain_init.F90 @@ -0,0 +1,330 @@ +!> Describes the decomposed MOM domain and has routines for communications across PEs +module MOM_domain_init + +! This file is part of MOM6. See LICENSE.md for the license. + +use MOM_coms, only : num_PEs +use MOM_domains, only : MOM_domain_type, create_MOM_domain, MOM_define_layout +use MOM_domains, only : MOM_thread_affinity_set, set_MOM_thread_affinity +use MOM_error_handler, only : MOM_error, MOM_mesg, NOTE, WARNING, FATAL +use MOM_file_parser, only : get_param, log_param, log_version, param_file_type +use MOM_io, only : file_exists +use MOM_string_functions, only : slasher + +implicit none ; private + +public :: MOM_domains_init, MOM_domain_type + +contains + +!> MOM_domains_init initializes a MOM_domain_type variable, based on the information +!! read in from a param_file_type, and optionally returns data describing various' +!! properties of the domain type. +subroutine MOM_domains_init(MOM_dom, param_file, symmetric, static_memory, & + NIHALO, NJHALO, NIGLOBAL, NJGLOBAL, NIPROC, NJPROC, & + min_halo, domain_name, include_name, param_suffix) + type(MOM_domain_type), pointer :: MOM_dom !< A pointer to the MOM_domain_type + !! being defined here. + type(param_file_type), intent(in) :: param_file !< A structure to parse for + !! run-time parameters + logical, optional, intent(in) :: symmetric !< If present, this specifies + !! whether this domain is symmetric, regardless of + !! whether the macro SYMMETRIC_MEMORY_ is defined. + logical, optional, intent(in) :: static_memory !< If present and true, this + !! domain type is set up for static memory and + !! error checking of various input values is + !! performed against those in the input file. + integer, optional, intent(in) :: NIHALO !< Default halo sizes, required + !! with static memory. + integer, optional, intent(in) :: NJHALO !< Default halo sizes, required + !! with static memory. + integer, optional, intent(in) :: NIGLOBAL !< Total domain sizes, required + !! with static memory. + integer, optional, intent(in) :: NJGLOBAL !< Total domain sizes, required + !! with static memory. + integer, optional, intent(in) :: NIPROC !< Processor counts, required with + !! static memory. + integer, optional, intent(in) :: NJPROC !< Processor counts, required with + !! static memory. + integer, dimension(2), optional, intent(inout) :: min_halo !< If present, this sets the + !! minimum halo size for this domain in the i- and j- + !! directions, and returns the actual halo size used. + character(len=*), optional, intent(in) :: domain_name !< A name for this domain, "MOM" + !! if missing. + character(len=*), optional, intent(in) :: include_name !< A name for model's include file, + !! "MOM_memory.h" if missing. + character(len=*), optional, intent(in) :: param_suffix !< A suffix to apply to + !! layout-specific parameters. + + ! Local variables + integer, dimension(2) :: layout ! The number of logical processors in the i- and j- directions + integer, dimension(2) :: io_layout ! The layout of logical processors for input and output + !$ integer :: ocean_nthreads ! Number of openMP threads + !$ logical :: ocean_omp_hyper_thread ! If true use openMP hyper-threads + integer, dimension(2) :: n_global ! The number of i- and j- points in the global computational domain + integer, dimension(2) :: n_halo ! The number of i- and j- points in the halos + integer :: nihalo_dflt, njhalo_dflt ! The default halo sizes + integer :: PEs_used ! The number of processors used + logical, dimension(2) :: reentrant ! True if the x- and y- directions are periodic. + logical, dimension(2,2) :: tripolar ! A set of flag indicating whether there is tripolar + ! connectivity for any of the four logical edges of the grid. + ! Currently only tripolar_N is implemented. + logical :: is_static ! If true, static memory is being used for this domain. + logical :: is_symmetric ! True if the domain being set up will use symmetric memory. + logical :: nonblocking ! If true, nonblocking halo updates will be used. + logical :: thin_halos ! If true, If true, optional arguments may be used to specify the + ! width of the halos that are updated with each call. + logical :: mask_table_exists ! True if there is a mask table file + character(len=128) :: inputdir ! The directory in which to find the diag table + character(len=200) :: mask_table ! The file name and later the full path to the diag table + character(len=64) :: inc_nm ! The name of the memory include file + character(len=200) :: mesg ! A string to use for error messages + + integer :: nip_parsed, njp_parsed + character(len=8) :: char_xsiz, char_ysiz, char_niglobal, char_njglobal + character(len=40) :: nihalo_nm, njhalo_nm, layout_nm, io_layout_nm, masktable_nm + character(len=40) :: niproc_nm, njproc_nm + ! This include declares and sets the variable "version". +# include "version_variable.h" + character(len=40) :: mdl ! This module's name. + + PEs_used = num_PEs() + + mdl = "MOM_domains" !### Change this to "MOM_domain_init" + + is_symmetric = .true. ; if (present(symmetric)) is_symmetric = symmetric + if (present(min_halo)) mdl = trim(mdl)//" min_halo" + + inc_nm = "MOM_memory.h" ; if (present(include_name)) inc_nm = trim(include_name) + + nihalo_nm = "NIHALO" ; njhalo_nm = "NJHALO" + layout_nm = "LAYOUT" ; io_layout_nm = "IO_LAYOUT" ; masktable_nm = "MASKTABLE" + niproc_nm = "NIPROC" ; njproc_nm = "NJPROC" + if (present(param_suffix)) then ; if (len(trim(adjustl(param_suffix))) > 0) then + nihalo_nm = "NIHALO"//(trim(adjustl(param_suffix))) + njhalo_nm = "NJHALO"//(trim(adjustl(param_suffix))) + layout_nm = "LAYOUT"//(trim(adjustl(param_suffix))) + io_layout_nm = "IO_LAYOUT"//(trim(adjustl(param_suffix))) + masktable_nm = "MASKTABLE"//(trim(adjustl(param_suffix))) + niproc_nm = "NIPROC"//(trim(adjustl(param_suffix))) + njproc_nm = "NJPROC"//(trim(adjustl(param_suffix))) + endif ; endif + + is_static = .false. ; if (present(static_memory)) is_static = static_memory + if (is_static) then + if (.not.present(NIHALO)) call MOM_error(FATAL, "NIHALO must be "// & + "present in the call to MOM_domains_init with static memory.") + if (.not.present(NJHALO)) call MOM_error(FATAL, "NJHALO must be "// & + "present in the call to MOM_domains_init with static memory.") + if (.not.present(NIGLOBAL)) call MOM_error(FATAL, "NIGLOBAL must be "// & + "present in the call to MOM_domains_init with static memory.") + if (.not.present(NJGLOBAL)) call MOM_error(FATAL, "NJGLOBAL must be "// & + "present in the call to MOM_domains_init with static memory.") + if (.not.present(NIPROC)) call MOM_error(FATAL, "NIPROC must be "// & + "present in the call to MOM_domains_init with static memory.") + if (.not.present(NJPROC)) call MOM_error(FATAL, "NJPROC must be "// & + "present in the call to MOM_domains_init with static memory.") + endif + + ! Read all relevant parameters and write them to the model log. + call log_version(param_file, mdl, version, "", log_to_all=.true., layout=.true.) + call get_param(param_file, mdl, "REENTRANT_X", reentrant(1), & + "If true, the domain is zonally reentrant.", default=.true.) + call get_param(param_file, mdl, "REENTRANT_Y", reentrant(2), & + "If true, the domain is meridionally reentrant.", & + default=.false.) + tripolar(1:2,1:2) = .false. + call get_param(param_file, mdl, "TRIPOLAR_N", tripolar(2,2), & + "Use tripolar connectivity at the northern edge of the "//& + "domain. With TRIPOLAR_N, NIGLOBAL must be even.", & + default=.false.) + +# ifndef NOT_SET_AFFINITY + !$ if (.not.MOM_thread_affinity_set()) then + !$ call get_param(param_file, mdl, "OCEAN_OMP_THREADS", ocean_nthreads, & + !$ "The number of OpenMP threads that MOM6 will use.", & + !$ default = 1, layoutParam=.true.) + !$ call get_param(param_file, mdl, "OCEAN_OMP_HYPER_THREAD", ocean_omp_hyper_thread, & + !$ "If True, use hyper-threading.", default = .false., layoutParam=.true.) + !$ call set_MOM_thread_affinity(ocean_nthreads, ocean_omp_hyper_thread) + !$ endif +# endif + + call log_param(param_file, mdl, "!SYMMETRIC_MEMORY_", is_symmetric, & + "If defined, the velocity point data domain includes every face of the "//& + "thickness points. In other words, some arrays are larger than others, "//& + "depending on where they are on the staggered grid. Also, the starting "//& + "index of the velocity-point arrays is usually 0, not 1. "//& + "This can only be set at compile time.",& + layoutParam=.true.) + call get_param(param_file, mdl, "NONBLOCKING_UPDATES", nonblocking, & + "If true, non-blocking halo updates may be used.", & + default=.false., layoutParam=.true.) + !### Note the duplicated "the the" in the following description, which should be fixed as a part + ! of a larger commit that also changes other MOM_parameter_doc file messages, but for now + ! reproduces the existing output files. + call get_param(param_file, mdl, "THIN_HALO_UPDATES", thin_halos, & + "If true, optional arguments may be used to specify the the width of the "//& + "halos that are updated with each call.", & + default=.true., layoutParam=.true.) + + nihalo_dflt = 4 ; njhalo_dflt = 4 + if (present(NIHALO)) nihalo_dflt = NIHALO + if (present(NJHALO)) njhalo_dflt = NJHALO + + call log_param(param_file, mdl, "!STATIC_MEMORY_", is_static, & + "If STATIC_MEMORY_ is defined, the principle variables will have sizes that "//& + "are statically determined at compile time. Otherwise the sizes are not "//& + "determined until run time. The STATIC option is substantially faster, but "//& + "does not allow the PE count to be changed at run time. This can only be "//& + "set at compile time.", layoutParam=.true.) + + if (is_static) then + call get_param(param_file, mdl, "NIGLOBAL", n_global(1), & + "The total number of thickness grid points in the x-direction in the physical "//& + "domain. With STATIC_MEMORY_ this is set in "//trim(inc_nm)//" at compile time.", & + static_value=NIGLOBAL) + call get_param(param_file, mdl, "NJGLOBAL", n_global(2), & + "The total number of thickness grid points in the y-direction in the physical "//& + "domain. With STATIC_MEMORY_ this is set in "//trim(inc_nm)//" at compile time.", & + static_value=NJGLOBAL) + if (n_global(1) /= NIGLOBAL) call MOM_error(FATAL,"MOM_domains_init: " // & + "static mismatch for NIGLOBAL_ domain size. Header file does not match input namelist") + if (n_global(2) /= NJGLOBAL) call MOM_error(FATAL,"MOM_domains_init: " // & + "static mismatch for NJGLOBAL_ domain size. Header file does not match input namelist") + + ! Check the requirement of equal sized compute domains when STATIC_MEMORY_ is used. + if ((MOD(NIGLOBAL, NIPROC) /= 0) .OR. (MOD(NJGLOBAL, NJPROC) /= 0)) then + write( char_xsiz, '(i4)' ) NIPROC + write( char_ysiz, '(i4)' ) NJPROC + write( char_niglobal, '(i4)' ) NIGLOBAL + write( char_njglobal, '(i4)' ) NJGLOBAL + call MOM_error(WARNING, 'MOM_domains: Processor decomposition (NIPROC_,NJPROC_) = ('//& + trim(char_xsiz)//','//trim(char_ysiz)//') does not evenly divide size '//& + 'set by preprocessor macro ('//trim(char_niglobal)//','//trim(char_njglobal)//').') + call MOM_error(FATAL,'MOM_domains: #undef STATIC_MEMORY_ in '//trim(inc_nm)//' to use '//& + 'dynamic allocation, or change processor decomposition to evenly divide the domain.') + endif + else + call get_param(param_file, mdl, "NIGLOBAL", n_global(1), & + "The total number of thickness grid points in the x-direction in the physical "//& + "domain. With STATIC_MEMORY_ this is set in "//trim(inc_nm)//" at compile time.", & + fail_if_missing=.true.) + call get_param(param_file, mdl, "NJGLOBAL", n_global(2), & + "The total number of thickness grid points in the y-direction in the physical "//& + "domain. With STATIC_MEMORY_ this is set in "//trim(inc_nm)//" at compile time.", & + fail_if_missing=.true.) + endif + + call get_param(param_file, mdl, trim(nihalo_nm), n_halo(1), & + "The number of halo points on each side in the x-direction. How this is set "//& + "varies with the calling component and static or dynamic memory configuration.", & + default=nihalo_dflt, static_value=nihalo_dflt) + call get_param(param_file, mdl, trim(njhalo_nm), n_halo(2), & + "The number of halo points on each side in the y-direction. How this is set "//& + "varies with the calling component and static or dynamic memory configuration.", & + default=njhalo_dflt, static_value=njhalo_dflt) + if (present(min_halo)) then + n_halo(1) = max(n_halo(1), min_halo(1)) + min_halo(1) = n_halo(1) + n_halo(2) = max(n_halo(2), min_halo(2)) + min_halo(2) = n_halo(2) + ! These are generally used only with static memory, so they are considerd layout params. + call log_param(param_file, mdl, "!NIHALO min_halo", n_halo(1), layoutParam=.true.) + call log_param(param_file, mdl, "!NJHALO min_halo", n_halo(2), layoutParam=.true.) + endif + if (is_static .and. .not.present(min_halo)) then + if (n_halo(1) /= NIHALO) call MOM_error(FATAL,"MOM_domains_init: " // & + "static mismatch for "//trim(nihalo_nm)//" domain size") + if (n_halo(2) /= NJHALO) call MOM_error(FATAL,"MOM_domains_init: " // & + "static mismatch for "//trim(njhalo_nm)//" domain size") + endif + + call get_param(param_file, mdl, "INPUTDIR", inputdir, do_not_log=.true., default=".") + inputdir = slasher(inputdir) + + call get_param(param_file, mdl, trim(masktable_nm), mask_table, & + "A text file to specify n_mask, layout and mask_list. This feature masks out "//& + "processors that contain only land points. The first line of mask_table is the "//& + "number of regions to be masked out. The second line is the layout of the "//& + "model and must be consistent with the actual model layout. The following "//& + "(n_mask) lines give the logical positions of the processors that are masked "//& + "out. The mask_table can be created by tools like check_mask. The following "//& + "example of mask_table masks out 2 processors, (1,2) and (3,6), out of the 24 "//& + "in a 4x6 layout: \n 2\n 4,6\n 1,2\n 3,6\n", default="MOM_mask_table", & + layoutParam=.true.) + mask_table = trim(inputdir)//trim(mask_table) + mask_table_exists = file_exists(mask_table) + + if (is_static) then + layout(1) = NIPROC ; layout(2) = NJPROC + else + call get_param(param_file, mdl, trim(layout_nm), layout, & + "The processor layout to be used, or 0, 0 to automatically set the layout "//& + "based on the number of processors.", default=0, do_not_log=.true.) + call get_param(param_file, mdl, trim(niproc_nm), nip_parsed, & + "The number of processors in the x-direction.", default=-1, do_not_log=.true.) + call get_param(param_file, mdl, trim(njproc_nm), njp_parsed, & + "The number of processors in the y-direction.", default=-1, do_not_log=.true.) + if (nip_parsed > -1) then + if ((layout(1) > 0) .and. (layout(1) /= nip_parsed)) & + call MOM_error(FATAL, trim(layout_nm)//" and "//trim(niproc_nm)//" set inconsistently. "//& + "Only LAYOUT should be used.") + layout(1) = nip_parsed + call MOM_mesg(trim(niproc_nm)//" used to set "//trim(layout_nm)//" in dynamic mode. "//& + "Shift to using "//trim(layout_nm)//" instead.") + endif + if (njp_parsed > -1) then + if ((layout(2) > 0) .and. (layout(2) /= njp_parsed)) & + call MOM_error(FATAL, trim(layout_nm)//" and "//trim(njproc_nm)//" set inconsistently. "//& + "Only "//trim(layout_nm)//" should be used.") + layout(2) = njp_parsed + call MOM_mesg(trim(njproc_nm)//" used to set "//trim(layout_nm)//" in dynamic mode. "//& + "Shift to using "//trim(layout_nm)//" instead.") + endif + + if ( (layout(1) == 0) .and. (layout(2) == 0) ) & + call MOM_define_layout( (/ 1, n_global(1), 1, n_global(2) /), PEs_used, layout) + if ( (layout(1) /= 0) .and. (layout(2) == 0) ) layout(2) = PEs_used / layout(1) + if ( (layout(1) == 0) .and. (layout(2) /= 0) ) layout(1) = PEs_used / layout(2) + + if (layout(1)*layout(2) /= PEs_used .and. (.not. mask_table_exists) ) then + write(mesg,'("MOM_domains_init: The product of the two components of layout, ", & + & 2i4,", is not the number of PEs used, ",i5,".")') & + layout(1), layout(2), PEs_used + call MOM_error(FATAL, mesg) + endif + endif + call log_param(param_file, mdl, trim(niproc_nm), layout(1), & + "The number of processors in the x-direction. With STATIC_MEMORY_ this "//& + "is set in "//trim(inc_nm)//" at compile time.", layoutParam=.true.) + call log_param(param_file, mdl, trim(njproc_nm), layout(2), & + "The number of processors in the y-direction. With STATIC_MEMORY_ this "//& + "is set in "//trim(inc_nm)//" at compile time.", layoutParam=.true.) + call log_param(param_file, mdl, trim(layout_nm), layout, & + "The processor layout that was actually used.", layoutParam=.true.) + + ! Idiot check that fewer PEs than columns have been requested + if (layout(1)*layout(2) > n_global(1)*n_global(2)) then + write(mesg,'(a,2(i5,x,a))') 'You requested to use',layout(1)*layout(2), & + 'PEs but there are only', n_global(1)*n_global(2), 'columns in the model' + call MOM_error(FATAL, mesg) + endif + + if (mask_table_exists) & + call MOM_error(NOTE, 'MOM_domains_init: reading maskmap information from '//trim(mask_table)) + + ! Set up the I/O layout, it will be checked later that it uses an even multiple of the number of + ! PEs in each direction. + io_layout(:) = (/ 1, 1 /) + call get_param(param_file, mdl, trim(io_layout_nm), io_layout, & + "The processor layout to be used, or 0,0 to automatically set the io_layout "//& + "to be the same as the layout.", default=1, layoutParam=.true.) + + call create_MOM_domain(MOM_dom, n_global, n_halo, reentrant, tripolar, layout, & + io_layout=io_layout, domain_name=domain_name, mask_table=mask_table, & + symmetric=symmetric, thin_halos=thin_halos, nonblocking=nonblocking) + +end subroutine MOM_domains_init + +end module MOM_domain_init diff --git a/src/framework/MOM_restart.F90 b/src/framework/MOM_restart.F90 index 73c5cc94e1..8034e19048 100644 --- a/src/framework/MOM_restart.F90 +++ b/src/framework/MOM_restart.F90 @@ -3,10 +3,10 @@ module MOM_restart ! This file is part of MOM6. See LICENSE.md for the license. +use MOM_checksums, only : chksum => rotated_field_chksum use MOM_domains, only : PE_here, num_PEs use MOM_error_handler, only : MOM_error, FATAL, WARNING, NOTE, is_root_pe use MOM_file_parser, only : get_param, log_param, log_version, param_file_type -use MOM_string_functions, only : lowercase use MOM_grid, only : ocean_grid_type use MOM_io, only : create_file, fieldtype, file_exists, open_file, close_file use MOM_io, only : MOM_read_data, read_data, MOM_write_field, read_field_chksum @@ -14,9 +14,9 @@ module MOM_restart use MOM_io, only : vardesc, var_desc, query_vardesc, modify_vardesc, get_filename_appendix use MOM_io, only : MULTIPLE, NETCDF_FILE, READONLY_FILE, SINGLE_FILE use MOM_io, only : CENTER, CORNER, NORTH_FACE, EAST_FACE +use MOM_string_functions, only : lowercase use MOM_time_manager, only : time_type, time_type_to_real, real_to_time use MOM_time_manager, only : days_in_month, get_date, set_date -use MOM_transform_FMS, only : chksum => rotated_mpp_chksum use MOM_verticalGrid, only : verticalGrid_type implicit none ; private @@ -874,7 +874,7 @@ subroutine save_restart(directory, time, G, CS, time_stamped, filename, GV, num_ ! this should be 2 Gb or less. integer :: start_var, next_var ! The starting variables of the ! current and next files. - integer :: unit ! The mpp unit of the open file. + integer :: unit ! The I/O unit of the open file. integer :: m, nz, num_files, var_periods integer :: seconds, days, year, month, hour, minute character(len=8) :: hor_grid, z_grid, t_grid ! Variable grid info. @@ -1086,7 +1086,7 @@ subroutine restore_state(filename, directory, day, G, CS) integer :: sizes(7) integer :: ndim, nvar, natt, ntime, pos - integer :: unit(CS%max_fields) ! The mpp unit of all open files. + integer :: unit(CS%max_fields) ! The I/O units of all open files. character(len=200) :: unit_path(CS%max_fields) ! The file names. logical :: unit_is_global(CS%max_fields) ! True if the file is global. @@ -1363,7 +1363,7 @@ function open_restart_units(filename, directory, G, CS, units, file_paths, & type(MOM_restart_CS), pointer :: CS !< The control structure returned by a previous !! call to restart_init. integer, dimension(:), & - optional, intent(out) :: units !< The mpp units of all opened files. + optional, intent(out) :: units !< The I/O units of all opened files. character(len=*), dimension(:), & optional, intent(out) :: file_paths !< The full paths to open files. logical, dimension(:), & diff --git a/src/framework/MOM_transform_FMS.F90 b/src/framework/MOM_transform_FMS.F90 index a6c28717b3..a4a3f7c2c4 100644 --- a/src/framework/MOM_transform_FMS.F90 +++ b/src/framework/MOM_transform_FMS.F90 @@ -3,28 +3,16 @@ module MOM_transform_FMS -use horiz_interp_mod, only : horiz_interp_type +use MOM_array_transform, only : allocate_rotated_array, rotate_array use MOM_error_handler, only : MOM_error, FATAL -use mpp_mod, only : mpp_chksum +use horiz_interp_mod, only : horiz_interp_type use time_manager_mod, only : time_type use time_interp_external_mod, only : time_interp_external -use MOM_array_transform, only : allocate_rotated_array, rotate_array - implicit none ; private -public rotated_mpp_chksum public rotated_time_interp_external -!> Rotate and compute the FMS (mpp) checksum of a field -interface rotated_mpp_chksum - module procedure rotated_mpp_chksum_real_0d - module procedure rotated_mpp_chksum_real_1d - module procedure rotated_mpp_chksum_real_2d - module procedure rotated_mpp_chksum_real_3d - module procedure rotated_mpp_chksum_real_4d -end interface rotated_mpp_chksum - !> Read a field based on model time, and rotate to the model domain interface rotated_time_interp_external module procedure rotated_time_interp_external_0d @@ -37,122 +25,6 @@ module MOM_transform_FMS ! NOTE: No transformations are applied to the 0d and 1d field implementations, ! but are provided to maintain compatibility with the FMS interfaces. - -!> Compute the FMS (mpp) checksum of a scalar. -!! This function is provided to support the full FMS mpp_chksum interface. -function rotated_mpp_chksum_real_0d(field, pelist, mask_val, turns) & - result(chksum) - real, intent(in) :: field !> Input scalar - integer, optional, intent(in) :: pelist(:) !> PE list of ranks to checksum - real, optional, intent(in) :: mask_val !> FMS mask value - integer, optional, intent(in) :: turns !> Number of quarter turns - integer :: chksum !> FMS checksum of scalar - - if (present(turns)) & - call MOM_error(FATAL, "Rotation not supported for 0d fields.") - - chksum = mpp_chksum(field, pelist=pelist, mask_val=mask_val) -end function rotated_mpp_chksum_real_0d - - -!> Compute the FMS (mpp) checksum of a 1d field. -!! This function is provided to support the full FMS mpp_chksum interface. -function rotated_mpp_chksum_real_1d(field, pelist, mask_val, turns) & - result(chksum) - real, intent(in) :: field(:) !> Input field - integer, optional, intent(in) :: pelist(:) !> PE list of ranks to checksum - real, optional, intent(in) :: mask_val !> FMS mask value - integer, optional, intent(in) :: turns !> Number of quarter-turns - integer :: chksum !> FMS checksum of field - - if (present(turns)) & - call MOM_error(FATAL, "Rotation not supported for 1d fields.") - - chksum = mpp_chksum(field, pelist=pelist, mask_val=mask_val) -end function rotated_mpp_chksum_real_1d - - -!> Compute the FMS (mpp) checksum of a rotated 2d field. -function rotated_mpp_chksum_real_2d(field, pelist, mask_val, turns) & - result(chksum) - real, intent(in) :: field(:,:) !> Unrotated input field - integer, optional, intent(in) :: pelist(:) !> PE list of ranks to checksum - real, optional, intent(in) :: mask_val !> FMS mask value - integer, optional, intent(in) :: turns !> Number of quarter-turns - integer :: chksum !> FMS checksum of field - - real, allocatable :: field_rot(:,:) - integer :: qturns - - qturns = 0 - if (present(turns)) & - qturns = modulo(turns, 4) - - if (qturns == 0) then - chksum = mpp_chksum(field, pelist=pelist, mask_val=mask_val) - else - call allocate_rotated_array(field, [1,1], qturns, field_rot) - call rotate_array(field, qturns, field_rot) - chksum = mpp_chksum(field_rot, pelist=pelist, mask_val=mask_val) - deallocate(field_rot) - endif -end function rotated_mpp_chksum_real_2d - - -!> Compute the FMS (mpp) checksum of a rotated 3d field. -function rotated_mpp_chksum_real_3d(field, pelist, mask_val, turns) & - result(chksum) - real, intent(in) :: field(:,:,:) !> Unrotated input field - integer, optional, intent(in) :: pelist(:) !> PE list of ranks to checksum - real, optional, intent(in) :: mask_val !> FMS mask value - integer, optional, intent(in) :: turns !> Number of quarter-turns - integer :: chksum !> FMS checksum of field - - real, allocatable :: field_rot(:,:,:) - integer :: qturns - - qturns = 0 - if (present(turns)) & - qturns = modulo(turns, 4) - - if (qturns == 0) then - chksum = mpp_chksum(field, pelist=pelist, mask_val=mask_val) - else - call allocate_rotated_array(field, [1,1,1], qturns, field_rot) - call rotate_array(field, qturns, field_rot) - chksum = mpp_chksum(field_rot, pelist=pelist, mask_val=mask_val) - deallocate(field_rot) - endif -end function rotated_mpp_chksum_real_3d - - -!> Compute the FMS (mpp) checksum of a rotated 4d field. -function rotated_mpp_chksum_real_4d(field, pelist, mask_val, turns) & - result(chksum) - real, intent(in) :: field(:,:,:,:) !> Unrotated input field - integer, optional, intent(in) :: pelist(:) !> PE list of ranks to checksum - real, optional, intent(in) :: mask_val !> FMS mask value - integer, optional, intent(in) :: turns !> Number of quarter-turns - integer :: chksum !> FMS checksum of field - - real, allocatable :: field_rot(:,:,:,:) - integer :: qturns - - qturns = 0 - if (present(turns)) & - qturns = modulo(turns, 4) - - if (qturns == 0) then - chksum = mpp_chksum(field, pelist=pelist, mask_val=mask_val) - else - call allocate_rotated_array(field, [1,1,1,1], qturns, field_rot) - call rotate_array(field, qturns, field_rot) - chksum = mpp_chksum(field_rot, pelist=pelist, mask_val=mask_val) - deallocate(field_rot) - endif -end function rotated_mpp_chksum_real_4d - - !> Read a scalar field based on model time !! This function is provided to support the full FMS time_interp_external !! interface.