diff --git a/CHANGELOG.md b/CHANGELOG.md index 29e90034c156..4199be2efd80 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -36,6 +36,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 - Added capability for HistoryCollectionGridComp to extract field names from expressions - Added ability for HistoryCollectionGridComp to extract multiple field names from expressions - Added vertical and ungridded dimensions to output for History3G +- Create rank-agnostic representation of `ESMF_Field` objects as rank-3 array pointers. ### Changed diff --git a/GeomIO/SharedIO.F90 b/GeomIO/SharedIO.F90 index e2d75441a8d1..77c1774d93f0 100644 --- a/GeomIO/SharedIO.F90 +++ b/GeomIO/SharedIO.F90 @@ -8,7 +8,7 @@ module mapl3g_SharedIO use MAPL_BaseMod use mapl3g_UngriddedDims use mapl3g_UngriddedDim - use mapl3g_output_info + use mapl3g_FieldDimensionInfo implicit none diff --git a/base/CMakeLists.txt b/base/CMakeLists.txt index 4a8120b9ced3..8da90b1e4cb4 100644 --- a/base/CMakeLists.txt +++ b/base/CMakeLists.txt @@ -56,7 +56,6 @@ set (srcs MAPL_XYGridFactory.F90 MAPL_NetCDF.F90 Plain_netCDF_Time.F90 MAPL_DateTime_Parsing_ESMF.F90 MAPL_ObsUtil.F90 - MAPL_ESMF_InfoKeys.F90 # Orphaned program: should not be in this library. # tstqsat.F90 ) diff --git a/esmf_utils/CMakeLists.txt b/esmf_utils/CMakeLists.txt index 362155ea897f..fdb11f971418 100644 --- a/esmf_utils/CMakeLists.txt +++ b/esmf_utils/CMakeLists.txt @@ -1,7 +1,7 @@ esma_set_this (OVERRIDE MAPL.esmf_utils) set(srcs - OutputInfo.F90 + FieldDimensionInfo.F90 UngriddedDim.F90 UngriddedDims.F90 UngriddedDimVector.F90 @@ -10,7 +10,7 @@ set(srcs esma_add_library(${this} SRCS ${srcs} - DEPENDENCIES MAPL.shared MAPL.base + DEPENDENCIES MAPL.shared TYPE SHARED ) @@ -18,3 +18,6 @@ target_include_directories (${this} PUBLIC $) target_link_libraries (${this} PUBLIC ESMF::ESMF) +if (PFUNIT_FOUND) + add_subdirectory(tests) +endif () diff --git a/esmf_utils/OutputInfo.F90 b/esmf_utils/FieldDimensionInfo.F90 similarity index 89% rename from esmf_utils/OutputInfo.F90 rename to esmf_utils/FieldDimensionInfo.F90 index 43248e648204..941005341b34 100644 --- a/esmf_utils/OutputInfo.F90 +++ b/esmf_utils/FieldDimensionInfo.F90 @@ -1,5 +1,5 @@ #include "MAPL_Generic.h" -module mapl3g_output_info +module mapl3g_FieldDimensionInfo use mapl3g_UngriddedDim use mapl3g_UngriddedDimVector @@ -10,8 +10,8 @@ module mapl3g_output_info use esmf, only: ESMF_Info, ESMF_InfoIsPresent use esmf, only: ESMF_InfoDestroy, ESMF_InfoCreate use esmf, only: ESMF_InfoGet, ESMF_InfoGetFromHost - use esmf, only: ESMF_InfoGetAlloc, ESMF_InfoGetCharAlloc - use esmf, only: ESMF_InfoPrint + use esmf, only: ESMF_InfoGetAlloc, ESMF_InfoPrint + use esmf, only: ESMF_MAXSTR, ESMF_SUCCESS use Mapl_ErrorHandling implicit none @@ -93,12 +93,11 @@ integer function get_num_levels_info(info, rc) result(num) type(ESMF_Info), intent(in) :: info integer, optional, intent(out) :: rc integer :: status - logical :: is_none + character(len=:), allocatable :: spec_name num = 0 - is_none = VERT_DIM_NONE == get_vertical_dim_spec_info(info, _RC) - _RETURN_IF(is_none) - + spec_name = get_vertical_dim_spec_info(info, _RC) + _RETURN_IF(spec_name == VERT_DIM_NONE) call ESMF_InfoGet(info, key=KEY_NUM_LEVELS, value=num, _RC) _RETURN(_SUCCESS) @@ -123,12 +122,12 @@ function get_vertical_dim_spec_names_bundle_info(info, rc) result(names) integer, optional, intent(out) :: rc integer :: status integer :: i - character(len=:), allocatable :: name + character(len=:), allocatable :: spec_name names = StringVector() do i=1, size(info) - name = get_vertical_dim_spec_info(info(i), _RC) - if(find_index(names, name) == 0) call names%push_back(name) + spec_name = get_vertical_dim_spec_info(info(i), _RC) + if(find_index(names, spec_name) == 0) call names%push_back(spec_name) end do _RETURN(_SUCCESS) @@ -152,8 +151,14 @@ function get_vertical_dim_spec_info(info, rc) result(spec_name) type(ESMF_Info), intent(in) :: info integer, optional, intent(out) :: rc integer :: status + logical :: isPresent + character(len=ESMF_MAXSTR) :: raw + + isPresent = ESMF_InfoIsPresent(info, key=KEY_VLOC, _RC) + _ASSERT(isPresent, 'Failed to get vertical dim spec name.') + call ESMF_InfoGet(info, key=KEY_VLOC, value=raw, _RC) + spec_name = trim(adjustl(raw)) - call ESMF_InfoGetCharAlloc(info, key=KEY_VLOC, value=spec_name, _RC) _RETURN(_SUCCESS) end function get_vertical_dim_spec_info @@ -225,8 +230,9 @@ function make_ungridded_dim(info, n, rc) result(ungridded_dim) type(ESMF_Info), intent(in) :: info integer, optional, intent(out) :: rc integer :: status - character(len=:), allocatable :: key type(ESMF_Info) :: dim_info + character(len=ESMF_MAXSTR) :: raw + character(len=:), allocatable :: key character(len=:), allocatable :: name character(len=:), allocatable :: units real, allocatable :: coordinates(:) @@ -237,11 +243,13 @@ function make_ungridded_dim(info, n, rc) result(ungridded_dim) call ESMF_InfoGet(info, key=key, isPresent=is_present, _RC) if(.not. is_present) then call ESMF_InfoPrint(info, unit=json_repr, _RC) + _FAIL('Key ' // trim(key) // ' not found in ' // trim(json_repr)) end if - _ASSERT(is_present, 'Key ' // key // ' not found in ' // trim(json_repr)) dim_info = ESMF_InfoCreate(info, key=key, _RC) - call ESMF_InfoGetCharAlloc(dim_info, key=KEY_UNGRIDDED_NAME, value=name, _RC) - call ESMF_InfoGetCharAlloc(dim_info, key=KEY_UNGRIDDED_UNITS, value=units, _RC) + call ESMF_InfoGet(dim_info, key=KEY_UNGRIDDED_NAME, value=raw, _RC) + name = trim(adjustl(raw)) + call ESMF_InfoGet(dim_info, key=KEY_UNGRIDDED_UNITS, value=raw, _RC) + units = trim(adjustl(raw)) call ESMF_InfoGetAlloc(dim_info, key=KEY_UNGRIDDED_COORD, values=coordinates, _RC) call ESMF_InfoDestroy(dim_info, _RC) ungridded_dim = UngriddedDim(coordinates, name=name, units=units) @@ -284,7 +292,6 @@ subroutine check_duplicate(vec, udim, rc) class(UngriddedDimVector), intent(in) :: vec class(UngriddedDim), intent(in) :: udim integer, optional, intent(out) :: rc - integer :: status type(UngriddedDimVectorIterator) :: iter type(UngriddedDim) :: vdim @@ -300,20 +307,12 @@ subroutine check_duplicate(vec, udim, rc) end subroutine check_duplicate - logical function is_vertical_dim_none(s) - character(len=*), intent(in) :: s - - is_vertical_dim_none = s == 'VERTICAL_DIM_NONE' - - end function is_vertical_dim_none - function create_bundle_info(bundle, rc) result(bundle_info) type(ESMF_Info), allocatable :: bundle_info(:) type(ESMF_FieldBundle), intent(in) :: bundle integer, optional, intent(out) :: rc integer :: status integer :: field_count, i - type(ESMF_Field) :: field type(ESMF_Field), allocatable :: fields(:) type(ESMF_Info) :: info @@ -343,4 +342,4 @@ subroutine destroy_bundle_info(bundle_info, rc) end subroutine destroy_bundle_info -end module mapl3g_output_info +end module mapl3g_FieldDimensionInfo diff --git a/esmf_utils/tests/CMakeLists.txt b/esmf_utils/tests/CMakeLists.txt new file mode 100644 index 000000000000..4dbe5299ae66 --- /dev/null +++ b/esmf_utils/tests/CMakeLists.txt @@ -0,0 +1,25 @@ +set(MODULE_DIRECTORY "${esma_include}/MAPL.esmf_utils.tests") + +set (test_srcs + Test_FieldDimensionInfo.pf + ) + +add_pfunit_ctest(MAPL.esmf_utils.tests + TEST_SOURCES ${test_srcs} + LINK_LIBRARIES MAPL.esmf_utils MAPL.pfunit + EXTRA_INITIALIZE Initialize + EXTRA_USE MAPL_pFUnit_Initialize + WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR} + MAX_PES 1 + ) +set_target_properties(MAPL.esmf_utils.tests PROPERTIES Fortran_MODULE_DIRECTORY ${MODULE_DIRECTORY}) +set_tests_properties(MAPL.esmf_utils.tests PROPERTIES LABELS "ESSENTIAL") + +if (APPLE) + set(LD_PATH "DYLD_LIBRARY_PATH") +else() + set(LD_PATH "LD_LIBRARY_PATH") +endif () +set_property(TEST MAPL.esmf_utils.tests PROPERTY ENVIRONMENT "${LD_PATH}=${CMAKE_CURRENT_BINARY_DIR}/esmf_utils:$ENV{${LD_PATH}}") + +add_dependencies(build-tests MAPL.esmf_utils.tests) diff --git a/gridcomps/History3G/tests/Test_OutputInfo.pf b/esmf_utils/tests/Test_FieldDimensionInfo.pf similarity index 98% rename from gridcomps/History3G/tests/Test_OutputInfo.pf rename to esmf_utils/tests/Test_FieldDimensionInfo.pf index 3e8ca30b8fcc..54110565fac2 100644 --- a/gridcomps/History3G/tests/Test_OutputInfo.pf +++ b/esmf_utils/tests/Test_FieldDimensionInfo.pf @@ -5,8 +5,8 @@ #define _SUCCESS 0 #define _FAILURE _SUCCESS-1 #include "MAPL_TestErr.h" -module Test_OutputInfo - use mapl3g_output_info +module Test_FieldDimensionInfo + use mapl3g_FieldDimensionInfo use mapl3g_esmf_info_keys use mapl3g_UngriddedDim use mapl3g_UngriddedDimVector @@ -16,10 +16,8 @@ module Test_OutputInfo implicit none - integer, parameter :: NUM_FIELDS_DEFAULT = 2 integer, parameter :: NUM_LEVELS_DEFAULT = 3 character(len=*), parameter :: VLOC_DEFAULT = 'VERTICAL_DIM_CENTER' - integer, parameter :: NUM_UNGRIDDED_DEFAULT = 3 character(len=*), parameter :: NAME_DEFAULT = 'A1' character(len=*), parameter :: UNITS_DEFAULT = 'stones' real, parameter :: COORDINATES_DEFAULT(*) = [2.0, 2.4, 2.5] @@ -250,4 +248,4 @@ contains if(allocated(info)) call deallocate_destroy(info) end subroutine safe_dealloc -end module Test_OutputInfo +end module Test_FieldDimensionInfo diff --git a/field_utils/CMakeLists.txt b/field_utils/CMakeLists.txt index 7fec50a25cf0..fec2a17ccc3e 100644 --- a/field_utils/CMakeLists.txt +++ b/field_utils/CMakeLists.txt @@ -8,6 +8,8 @@ set(srcs FieldUnaryFunctions.F90 FieldBinaryOperations.F90 FieldUnits.F90 + FieldCondensedArray.F90 + FieldCondensedArray_private.F90 ) # To use extended udunits2 procedures, udunits2.c must be built and linked. @@ -24,9 +26,10 @@ endif () esma_add_library(${this} SRCS ${srcs} - DEPENDENCIES MAPL.shared PFLOGGER::pflogger udunits2f + DEPENDENCIES MAPL.shared MAPL.esmf_utils PFLOGGER::pflogger udunits2f TYPE SHARED ) + #DEPENDENCIES MAPL.shared PFLOGGER::pflogger udunits2f #add_subdirectory(specs) #add_subdirectory(registry) diff --git a/field_utils/FieldCondensedArray.F90 b/field_utils/FieldCondensedArray.F90 new file mode 100644 index 000000000000..7d59ab717017 --- /dev/null +++ b/field_utils/FieldCondensedArray.F90 @@ -0,0 +1,72 @@ +#include "MAPL_Generic.h" +module mapl3g_FieldCondensedArray + use, intrinsic :: iso_c_binding, only: c_ptr, c_f_pointer + use mapl3g_FieldCondensedArray_private, only: ARRAY_RANK, get_shape => get_fptr_shape + use mapl3g_FieldDimensionInfo, only: get_vertical_dim_spec_name + use MAPL_FieldPointerUtilities, only: FieldGetLocalElementCount, assign_fptr + use MAPL_ExceptionHandling + use ESMF, only: ESMF_Field, ESMF_FieldGet + use ESMF, only: ESMF_KIND_R4, ESMF_KIND_R8, ESMF_KIND_I8 + + implicit none + private + public :: assign_fptr_condensed_array + + interface assign_fptr_condensed_array + module procedure :: assign_fptr_condensed_array_r4 + module procedure :: assign_fptr_condensed_array_r8 + end interface assign_fptr_condensed_array + +contains + + subroutine assign_fptr_condensed_array_r4(x, fptr, rc) + type(ESMF_Field), intent(inout) :: x + real(kind=ESMF_KIND_R4), pointer, intent(out) :: fptr(:,:,:) + integer, optional, intent(out) :: rc + integer(ESMF_KIND_I8) :: fp_shape(ARRAY_RANK) + integer :: status + + fp_shape = get_fptr_shape(x, _RC) + call assign_fptr(x, fp_shape, fptr, _RC) + _RETURN(_SUCCESS) + + end subroutine assign_fptr_condensed_array_r4 + + subroutine assign_fptr_condensed_array_r8(x, fptr, rc) + type(ESMF_Field), intent(inout) :: x + real(kind=ESMF_KIND_R8), pointer, intent(out) :: fptr(:,:,:) + integer, optional, intent(out) :: rc + integer(ESMF_KIND_I8) :: fp_shape(ARRAY_RANK) + integer :: status + + fp_shape = get_fptr_shape(x, _RC) + call assign_fptr(x, fp_shape, fptr, _RC) + _RETURN(_SUCCESS) + + end subroutine assign_fptr_condensed_array_r8 + + function get_fptr_shape(f, rc) result(fptr_shape) + integer :: fptr_shape(ARRAY_RANK) + type(ESMF_Field), intent(inout) :: f + integer, optional, intent(out) :: rc + integer :: status + integer :: rank + integer, allocatable :: gridToFieldMap(:) + integer, allocatable :: localElementCount(:) + logical :: has_vertical + character(len=:), allocatable :: spec_name + character(len=*), parameter :: VERTICAL_DIM_NONE_NAME = 'VERTICAL_DIM_NONE' + + call ESMF_FieldGet(f, gridToFieldMap=gridToFieldMap, _RC) + call ESMF_FieldGet(f, rank=rank, _RC) + allocate(localElementCount(rank)) + ! Due to an ESMF bug, getting the localElementCount must use the module function. + ! See FieldGetLocalElementCount (specific function) comments. + localElementCount = FieldGetLocalElementCount(f, _RC) + spec_name = get_vertical_dim_spec_name(f, _RC) + has_vertical = spec_name /= VERTICAL_DIM_NONE_NAME + fptr_shape = get_shape(gridToFieldMap, localElementCount, has_vertical, _RC) + + end function get_fptr_shape + +end module mapl3g_FieldCondensedArray diff --git a/field_utils/FieldCondensedArray_private.F90 b/field_utils/FieldCondensedArray_private.F90 new file mode 100644 index 000000000000..acc6db269038 --- /dev/null +++ b/field_utils/FieldCondensedArray_private.F90 @@ -0,0 +1,44 @@ +#include "MAPL_Generic.h" +module mapl3g_FieldCondensedArray_private + + use MAPL_ExceptionHandling + implicit none + + private + public :: get_fptr_shape, ARRAY_RANK + + integer, parameter :: ARRAY_RANK = 3 + +contains + + function get_fptr_shape(gridToFieldMap, localElementCount, has_vertical, rc) & + &result(fptr_shape) + integer :: fptr_shape(ARRAY_RANK) + integer, intent(in) :: gridToFieldMap(:) + integer, intent(in) :: localElementCount(:) + logical, intent(in) :: has_vertical + integer, optional, intent(out) :: rc + integer :: rank, i + integer, allocatable :: grid_dims(:) + integer, allocatable :: ungridded_dims(:) + integer :: horz_size, vert_size, ungridded_size + integer :: vert_dim + + vert_dim = 0 + vert_size = 1 + + rank = size(localElementCount) + grid_dims = pack(gridToFieldMap, gridToFieldMap /= 0) + _ASSERT(all(grid_dims <= size(grid_dims)), 'MAPL expects geom dims before ungridded.') + if(has_vertical) vert_dim = 1 + if(size(grid_dims) > 0) vert_dim = maxval(grid_dims) + vert_dim + ungridded_dims = pack([(i,i=1,rank)], [(all([vert_dim, grid_dims] /= i), i=1, rank)]) + horz_size = product([(localElementCount(grid_dims(i)), i=1, size(grid_dims))]) + if(has_vertical) vert_size = localElementCount(vert_dim) + ungridded_size = product([(localElementCount(ungridded_dims(i)), i=1, size(ungridded_dims))]) + fptr_shape = [horz_size, vert_size, ungridded_size] + _RETURN(_SUCCESS) + + end function get_fptr_shape + +end module mapl3g_FieldCondensedArray_private diff --git a/field_utils/FieldPointerUtilities.F90 b/field_utils/FieldPointerUtilities.F90 index 52a0f75e5eff..238b8ba24f9b 100644 --- a/field_utils/FieldPointerUtilities.F90 +++ b/field_utils/FieldPointerUtilities.F90 @@ -31,6 +31,8 @@ module MAPL_FieldPointerUtilities module procedure assign_fptr_r8_rank1 module procedure assign_fptr_r4_rank2 module procedure assign_fptr_r8_rank2 + module procedure assign_fptr_r4_rank3 + module procedure assign_fptr_r8_rank3 end interface assign_fptr interface FieldGetCptr @@ -78,8 +80,8 @@ module MAPL_FieldPointerUtilities interface MAPL_FieldDestroy procedure destroy end interface -contains +contains subroutine assign_fptr_r4_rank1(x, fptr, rc) type(ESMF_Field), intent(inout) :: x @@ -129,6 +131,7 @@ subroutine assign_fptr_r4_rank2(x, fp_shape, fptr, rc) type(c_ptr) :: cptr integer :: status + _ASSERT(size(fp_shape) == rank(fptr), 'Shape size must match pointer rank.') call FieldGetCptr(x, cptr, _RC) call c_f_pointer(cptr, fptr, fp_shape) @@ -145,12 +148,47 @@ subroutine assign_fptr_r8_rank2(x, fp_shape, fptr, rc) type(c_ptr) :: cptr integer :: status + _ASSERT(size(fp_shape) == rank(fptr), 'Shape size must match pointer rank.') call FieldGetCptr(x, cptr, _RC) call c_f_pointer(cptr, fptr, fp_shape) _RETURN(_SUCCESS) end subroutine assign_fptr_r8_rank2 + subroutine assign_fptr_r4_rank3(x, fp_shape, fptr, rc) + type(ESMF_Field), intent(inout) :: x + integer(ESMF_KIND_I8), intent(in) :: fp_shape(:) + real(kind=ESMF_KIND_R4), pointer, intent(out) :: fptr(:,:,:) + integer, optional, intent(out) :: rc + + ! local declarations + type(c_ptr) :: cptr + integer :: status + + _ASSERT(size(fp_shape) == rank(fptr), 'Shape size must match pointer rank.') + call FieldGetCptr(x, cptr, _RC) + call c_f_pointer(cptr, fptr, fp_shape) + + _RETURN(_SUCCESS) + end subroutine assign_fptr_r4_rank3 + + subroutine assign_fptr_r8_rank3(x, fp_shape, fptr, rc) + type(ESMF_Field), intent(inout) :: x + integer(ESMF_KIND_I8), intent(in) :: fp_shape(:) + real(kind=ESMF_KIND_R8), pointer, intent(out) :: fptr(:,:,:) + integer, optional, intent(out) :: rc + + ! local declarations + type(c_ptr) :: cptr + integer :: status + + _ASSERT(size(fp_shape) == rank(fptr), 'Shape size must match pointer rank.') + call FieldGetCptr(x, cptr, _RC) + call c_f_pointer(cptr, fptr, fp_shape) + + _RETURN(_SUCCESS) + end subroutine assign_fptr_r8_rank3 + subroutine get_cptr(x, cptr, rc) type(ESMF_Field), intent(inout) :: x type(c_ptr), intent(out) :: cptr @@ -525,7 +563,7 @@ logical function are_same_type_kind(x, y, rc) result(same_tk) _RETURN(_SUCCESS) end function are_same_type_kind - subroutine verify_typekind_scalar(x, expected_tk, rc) + subroutine verify_typekind_scalar(x, expected_tk, rc) type(ESMF_Field), intent(inout) :: x type(ESMF_TypeKind_Flag), intent(in) :: expected_tk integer, optional, intent(out) :: rc @@ -720,7 +758,7 @@ subroutine copy(x, y, rc) call FieldGetCptr(y, cptr_y, _RC) call ESMF_FieldGet(y, typekind = tk_y, _RC) - !wdb fixme convert between precisions ? get rid of extra cases + !wdb fixme convert between precisions ? get rid of extra cases y_is_double = (tk_y == ESMF_TYPEKIND_R8) _ASSERT(y_is_double .or. (tk_y == ESMF_TYPEKIND_R4), UNSUPPORTED_TK//'y.') @@ -796,159 +834,160 @@ subroutine copy_r8_r8(cptr_x, cptr_y, n) y_ptr=x_ptr end subroutine copy_r8_r8 -! this procedure must go away as soon as ESMF Fixes their bug - - subroutine MAPL_FieldGetLocalElementCount(field,local_count,rc) - type(ESMF_Field), intent(inout) :: field - integer, allocatable, intent(out) :: local_count(:) - integer, optional, intent(out) :: rc - - integer :: status, rank - type(ESMF_TypeKind_Flag) :: tk - - real(kind=ESMF_KIND_R4), pointer :: r4_1d(:),r4_2d(:,:),r4_3d(:,:,:),r4_4d(:,:,:,:) - real(kind=ESMF_KIND_R8), pointer :: r8_1d(:),r8_2d(:,:),r8_3d(:,:,:),r8_4d(:,:,:,:) - - call ESMF_FieldGet(field,rank=rank,typekind=tk,_RC) - if (tk == ESMF_TypeKind_R4) then - if (rank==1) then - call ESMF_FieldGet(field,0,farrayptr=r4_1d,_RC) - local_count = shape(r4_1d) - else if (rank ==2) then - call ESMF_FieldGet(field,0,farrayptr=r4_2d,_RC) - local_count = shape(r4_2d) - else if (rank ==3) then - call ESMF_FieldGet(field,0,farrayptr=r4_3d,_RC) - local_count = shape(r4_3d) - else if (rank ==4) then - call ESMF_FieldGet(field,0,farrayptr=r4_4d,_RC) - local_count = shape(r4_4d) - else - _FAIL("Unsupported rank") - end if - else if (tk == ESMF_TypeKind_R8) then - if (rank==1) then - call ESMF_FieldGet(field,0,farrayptr=r8_1d,_RC) - local_count = shape(r8_1d) - else if (rank ==2) then - call ESMF_FieldGet(field,0,farrayptr=r8_2d,_RC) - local_count = shape(r8_2d) - else if (rank ==3) then - call ESMF_FieldGet(field,0,farrayptr=r8_3d,_RC) - local_count = shape(r8_3d) - else if (rank ==4) then - call ESMF_FieldGet(field,0,farrayptr=r8_4d,_RC) - local_count = shape(r8_4d) - else - _FAIL("Unsupported rank") - end if - else - _FAIL("Unsupported type") - end if - _RETURN(_SUCCESS) - end subroutine MAPL_FieldGetLocalElementCount - - function FieldsHaveUndef(fields,rc) result(all_have_undef) - logical :: all_have_undef - type(ESMF_Field), intent(inout) :: fields(:) - integer, optional, intent(out) :: rc - - integer :: status, i - logical :: isPresent - type(ESMF_Info) :: infoh - - all_have_undef = .true. - do i =1,size(fields) - call ESMF_InfoGetFromHost(fields(i),infoh,_RC) - isPresent = ESMF_InfoIsPresent(infoh,"missing_value",_RC) - all_have_undef = (all_have_undef .and. isPresent) - enddo - _RETURN(_SUCCESS) - end function - - subroutine GetFieldsUndef_r4(fields,undef_values,rc) - type(ESMF_Field), intent(inout) :: fields(:) - real(kind=ESMF_KIND_R4), allocatable,intent(inout) :: undef_values(:) - integer, optional, intent(out) :: rc - - integer :: status, i - logical :: isPresent - type(ESMF_Info) :: infoh - - allocate(undef_values(size(fields))) - do i =1,size(fields) - call ESMF_InfoGetFromHost(fields(i),infoh,_RC) - isPresent = ESMF_InfoIsPresent(infoh,"missing_value",_RC) - _ASSERT(isPresent,"missing undef value") - call ESMF_InfoGet(infoh,value=undef_values(i),key="missing_value",_RC) - enddo - _RETURN(_SUCCESS) - end subroutine GetFieldsUndef_r4 - - subroutine GetFieldsUndef_r8(fields,undef_values,rc) - type(ESMF_Field), intent(inout) :: fields(:) - real(kind=ESMF_KIND_R8), allocatable,intent(inout) :: undef_values(:) - integer, optional, intent(out) :: rc - - integer :: status, i - logical :: isPresent - type(ESMF_Info) :: infoh - - allocate(undef_values(size(fields))) - do i =1,size(fields) - call ESMF_InfoGetFromHost(fields(i),infoh,_RC) - isPresent = ESMF_InfoIsPresent(infoh,"missing_value",_RC) - _ASSERT(isPresent,"missing undef value") - call ESMF_InfoGet(infoh,value=undef_values(i),key="missing_value",_RC) - enddo - _RETURN(_SUCCESS) - end subroutine GetFieldsUndef_r8 - -subroutine Destroy(Field,RC) - type(ESMF_Field), intent(INOUT) :: Field - integer, optional, intent(OUT ) :: RC - - integer :: STATUS - - real(kind=ESMF_KIND_R4), pointer :: VR4_1D(:), VR4_2D(:,:), VR4_3D(:,:,:), VR4_4D(:,:,:,:) - real(kind=ESMF_KIND_R8), pointer :: VR8_1D(:), VR8_2D(:,:), VR8_3D(:,:,:), VR8_4D(:,:,:,:) - integer :: rank - type(ESMF_TypeKind_Flag) :: tk - logical :: esmf_allocated - - call ESMF_FieldGet(Field,typekind=tk,dimCount=rank,isESMFAllocated=esmf_allocated,_RC) - if (.not. esmf_allocated) then - if (tk == ESMF_TYPEKIND_R4 .and. rank == 1) then - call ESMF_FieldGet(Field,0,VR4_1d,_RC) - deallocate(VR4_1d,_STAT) - else if (tk == ESMF_TYPEKIND_R8 .and. rank == 1) then - call ESMF_FieldGet(Field,0,VR8_1d,_RC) - deallocate(VR8_1d,_STAT) - else if (tk == ESMF_TYPEKIND_R4 .and. rank == 2) then - call ESMF_FieldGet(Field,0,VR4_2d,_RC) - deallocate(VR4_2d,_STAT) - else if (tk == ESMF_TYPEKIND_R8 .and. rank == 2) then - call ESMF_FieldGet(Field,0,VR8_2d,_RC) - deallocate(VR8_2d,_STAT) - else if (tk == ESMF_TYPEKIND_R4 .and. rank == 3) then - call ESMF_FieldGet(Field,0,VR4_3D,_RC) - deallocate(VR4_3d,_STAT) - else if (tk == ESMF_TYPEKIND_R8 .and. rank == 3) then - call ESMF_FieldGet(Field,0,VR8_3D,_RC) - deallocate(VR8_3d,_STAT) - else if (tk == ESMF_TYPEKIND_R4 .and. rank == 4) then - call ESMF_FieldGet(Field,0,VR4_4D,_RC) - deallocate(VR4_3d,_STAT) - else if (tk == ESMF_TYPEKIND_R8 .and. rank == 4) then - call ESMF_FieldGet(Field,0,VR8_4D,_RC) - deallocate(VR8_3d,_STAT) - else - _FAIL( 'unsupported typekind+rank') - end if - end if - call ESMF_FieldDestroy(Field,noGarbage = .true., rc=status) - _VERIFY(STATUS) - _RETURN(ESMF_SUCCESS) - - end subroutine Destroy -end module + ! this procedure must go away as soon as ESMF Fixes their bug + + subroutine MAPL_FieldGetLocalElementCount(field,local_count,rc) + type(ESMF_Field), intent(inout) :: field + integer, allocatable, intent(out) :: local_count(:) + integer, optional, intent(out) :: rc + + integer :: status, rank + type(ESMF_TypeKind_Flag) :: tk + + real(kind=ESMF_KIND_R4), pointer :: r4_1d(:),r4_2d(:,:),r4_3d(:,:,:),r4_4d(:,:,:,:) + real(kind=ESMF_KIND_R8), pointer :: r8_1d(:),r8_2d(:,:),r8_3d(:,:,:),r8_4d(:,:,:,:) + + call ESMF_FieldGet(field,rank=rank,typekind=tk,_RC) + if (tk == ESMF_TypeKind_R4) then + if (rank==1) then + call ESMF_FieldGet(field,0,farrayptr=r4_1d,_RC) + local_count = shape(r4_1d) + else if (rank ==2) then + call ESMF_FieldGet(field,0,farrayptr=r4_2d,_RC) + local_count = shape(r4_2d) + else if (rank ==3) then + call ESMF_FieldGet(field,0,farrayptr=r4_3d,_RC) + local_count = shape(r4_3d) + else if (rank ==4) then + call ESMF_FieldGet(field,0,farrayptr=r4_4d,_RC) + local_count = shape(r4_4d) + else + _FAIL("Unsupported rank") + end if + else if (tk == ESMF_TypeKind_R8) then + if (rank==1) then + call ESMF_FieldGet(field,0,farrayptr=r8_1d,_RC) + local_count = shape(r8_1d) + else if (rank ==2) then + call ESMF_FieldGet(field,0,farrayptr=r8_2d,_RC) + local_count = shape(r8_2d) + else if (rank ==3) then + call ESMF_FieldGet(field,0,farrayptr=r8_3d,_RC) + local_count = shape(r8_3d) + else if (rank ==4) then + call ESMF_FieldGet(field,0,farrayptr=r8_4d,_RC) + local_count = shape(r8_4d) + else + _FAIL("Unsupported rank") + end if + else + _FAIL("Unsupported type") + end if + _RETURN(_SUCCESS) + end subroutine MAPL_FieldGetLocalElementCount + + function FieldsHaveUndef(fields,rc) result(all_have_undef) + logical :: all_have_undef + type(ESMF_Field), intent(inout) :: fields(:) + integer, optional, intent(out) :: rc + + integer :: status, i + logical :: isPresent + type(ESMF_Info) :: infoh + + all_have_undef = .true. + do i =1,size(fields) + call ESMF_InfoGetFromHost(fields(i),infoh,_RC) + isPresent = ESMF_InfoIsPresent(infoh,"missing_value",_RC) + all_have_undef = (all_have_undef .and. isPresent) + enddo + _RETURN(_SUCCESS) + end function + + subroutine GetFieldsUndef_r4(fields,undef_values,rc) + type(ESMF_Field), intent(inout) :: fields(:) + real(kind=ESMF_KIND_R4), allocatable,intent(inout) :: undef_values(:) + integer, optional, intent(out) :: rc + + integer :: status, i + logical :: isPresent + type(ESMF_Info) :: infoh + + allocate(undef_values(size(fields))) + do i =1,size(fields) + call ESMF_InfoGetFromHost(fields(i),infoh,_RC) + isPresent = ESMF_InfoIsPresent(infoh,"missing_value",_RC) + _ASSERT(isPresent,"missing undef value") + call ESMF_InfoGet(infoh,value=undef_values(i),key="missing_value",_RC) + enddo + _RETURN(_SUCCESS) + end subroutine GetFieldsUndef_r4 + + subroutine GetFieldsUndef_r8(fields,undef_values,rc) + type(ESMF_Field), intent(inout) :: fields(:) + real(kind=ESMF_KIND_R8), allocatable,intent(inout) :: undef_values(:) + integer, optional, intent(out) :: rc + + integer :: status, i + logical :: isPresent + type(ESMF_Info) :: infoh + + allocate(undef_values(size(fields))) + do i =1,size(fields) + call ESMF_InfoGetFromHost(fields(i),infoh,_RC) + isPresent = ESMF_InfoIsPresent(infoh,"missing_value",_RC) + _ASSERT(isPresent,"missing undef value") + call ESMF_InfoGet(infoh,value=undef_values(i),key="missing_value",_RC) + enddo + _RETURN(_SUCCESS) + end subroutine GetFieldsUndef_r8 + + subroutine Destroy(Field,RC) + type(ESMF_Field), intent(INOUT) :: Field + integer, optional, intent(OUT ) :: RC + + integer :: STATUS + + real(kind=ESMF_KIND_R4), pointer :: VR4_1D(:), VR4_2D(:,:), VR4_3D(:,:,:), VR4_4D(:,:,:,:) + real(kind=ESMF_KIND_R8), pointer :: VR8_1D(:), VR8_2D(:,:), VR8_3D(:,:,:), VR8_4D(:,:,:,:) + integer :: rank + type(ESMF_TypeKind_Flag) :: tk + logical :: esmf_allocated + + call ESMF_FieldGet(Field,typekind=tk,dimCount=rank,isESMFAllocated=esmf_allocated,_RC) + if (.not. esmf_allocated) then + if (tk == ESMF_TYPEKIND_R4 .and. rank == 1) then + call ESMF_FieldGet(Field,0,VR4_1d,_RC) + deallocate(VR4_1d,_STAT) + else if (tk == ESMF_TYPEKIND_R8 .and. rank == 1) then + call ESMF_FieldGet(Field,0,VR8_1d,_RC) + deallocate(VR8_1d,_STAT) + else if (tk == ESMF_TYPEKIND_R4 .and. rank == 2) then + call ESMF_FieldGet(Field,0,VR4_2d,_RC) + deallocate(VR4_2d,_STAT) + else if (tk == ESMF_TYPEKIND_R8 .and. rank == 2) then + call ESMF_FieldGet(Field,0,VR8_2d,_RC) + deallocate(VR8_2d,_STAT) + else if (tk == ESMF_TYPEKIND_R4 .and. rank == 3) then + call ESMF_FieldGet(Field,0,VR4_3D,_RC) + deallocate(VR4_3d,_STAT) + else if (tk == ESMF_TYPEKIND_R8 .and. rank == 3) then + call ESMF_FieldGet(Field,0,VR8_3D,_RC) + deallocate(VR8_3d,_STAT) + else if (tk == ESMF_TYPEKIND_R4 .and. rank == 4) then + call ESMF_FieldGet(Field,0,VR4_4D,_RC) + deallocate(VR4_3d,_STAT) + else if (tk == ESMF_TYPEKIND_R8 .and. rank == 4) then + call ESMF_FieldGet(Field,0,VR8_4D,_RC) + deallocate(VR8_3d,_STAT) + else + _FAIL( 'unsupported typekind+rank') + end if + end if + call ESMF_FieldDestroy(Field,noGarbage = .true., rc=status) + _VERIFY(STATUS) + _RETURN(ESMF_SUCCESS) + + end subroutine Destroy + +end module MAPL_FieldPointerUtilities diff --git a/field_utils/tests/CMakeLists.txt b/field_utils/tests/CMakeLists.txt index 8ff68dd04668..880af840fc07 100644 --- a/field_utils/tests/CMakeLists.txt +++ b/field_utils/tests/CMakeLists.txt @@ -4,6 +4,7 @@ set(MODULE_DIRECTORY "${esma_include}/MAPL.field_utils.tests") set (test_srcs Test_FieldBLAS.pf Test_FieldArithmetic.pf + Test_FieldCondensedArray_private.pf ) diff --git a/field_utils/tests/Test_FieldCondensedArray_private.pf b/field_utils/tests/Test_FieldCondensedArray_private.pf new file mode 100644 index 000000000000..76078952d61f --- /dev/null +++ b/field_utils/tests/Test_FieldCondensedArray_private.pf @@ -0,0 +1,161 @@ +#include "MAPL_TestErr.h" +module Test_FieldCondensedArray_private + + use MAPL_ExceptionHandling + use pfunit + use mapl3g_FieldCondensedArray_private + implicit none + + character, parameter :: GENERIC_MESSAGE = 'actual does not match expected.' + +contains + + @Test + subroutine test_get_fptr_shape_3D() + integer :: expected(ARRAY_RANK), actual(ARRAY_RANK) + integer, allocatable :: gridToFieldMap(:) + integer, allocatable :: localElementCount(:) + logical :: has_vertical + + has_vertical = .TRUE. + gridToFieldMap = [1, 2] + localElementCount = [3, 5, 7] + expected = [product(localElementCount(1:2)), localElementCount(3), 1] + actual = get_fptr_shape(gridToFieldMap, localElementCount, has_vertical) + @assert_that(GENERIC_MESSAGE, actual, is(equal_to(expected))) + + end subroutine test_get_fptr_shape_3D + + @Test + subroutine test_get_fptr_shape_2D() + integer :: expected(ARRAY_RANK), actual(ARRAY_RANK) + integer, allocatable :: gridToFieldMap(:) + integer, allocatable :: localElementCount(:) + logical :: has_vertical + + has_vertical = .FALSE. + gridToFieldMap = [1, 2] + localElementCount = [3, 5] + expected = [product(localElementCount), 1, 1] + actual = get_fptr_shape(gridToFieldMap, localElementCount, has_vertical) + @assert_that(GENERIC_MESSAGE, actual, is(equal_to(expected))) + + end subroutine test_get_fptr_shape_2D + + @Test + subroutine test_get_fptr_shape_general() + integer :: expected(ARRAY_RANK), actual(ARRAY_RANK) + integer, allocatable :: gridToFieldMap(:) + integer, allocatable :: localElementCount(:) + logical :: has_vertical + + has_vertical = .TRUE. + gridToFieldMap = [1, 2] + localElementCount = [2, 3, 5, 7, 11] + expected = [product(localElementCount(1:2)), localElementCount(3), product(localElementCount(4:))] + actual = get_fptr_shape(gridToFieldMap, localElementCount, has_vertical) + @assert_that(GENERIC_MESSAGE, actual, is(equal_to(expected))) + + end subroutine test_get_fptr_shape_general + + @Test + subroutine test_get_fptr_shape_noz() + integer :: expected(ARRAY_RANK), actual(ARRAY_RANK) + integer, allocatable :: gridToFieldMap(:) + integer, allocatable :: localElementCount(:) + logical :: has_vertical + + has_vertical = .FALSE. + + gridToFieldMap = [1, 2] + localElementCount = [2, 3, 5, 7] + expected = [product(localElementCount(1:2)), 1, product(localElementCount(3:))] + actual = get_fptr_shape(gridToFieldMap, localElementCount, has_vertical) + @assert_that(GENERIC_MESSAGE, actual, is(equal_to(expected))) + + end subroutine test_get_fptr_shape_noz + + @Test + subroutine test_get_fptr_shape_0D() + integer :: expected(ARRAY_RANK), actual(ARRAY_RANK) + integer, allocatable :: gridToFieldMap(:) + integer, allocatable :: localElementCount(:) + logical :: has_vertical + + has_vertical = .FALSE. + gridToFieldMap = [0, 0] + localElementCount = [5, 7, 11] + expected = [1, 1, product(localElementCount)] + actual = get_fptr_shape(gridToFieldMap, localElementCount, has_vertical) + @assert_that(GENERIC_MESSAGE, actual, is(equal_to(expected))) + + end subroutine test_get_fptr_shape_0D + + @Test + subroutine test_get_fptr_shape_vert_only() + integer :: expected(ARRAY_RANK), actual(ARRAY_RANK) + integer, allocatable :: gridToFieldMap(:) + integer, allocatable :: localElementCount(:) + logical :: has_vertical + + has_vertical = .TRUE. + gridToFieldMap = [0, 0] + localElementCount = [3] + expected = [1, localElementCount(1), 1] + actual = get_fptr_shape(gridToFieldMap, localElementCount, has_vertical) + @assert_that(GENERIC_MESSAGE, actual, is(equal_to(expected))) + + end subroutine test_get_fptr_shape_vert_only + + @Test + subroutine test_get_fptr_shape_vert_ungrid() + integer :: expected(ARRAY_RANK), actual(ARRAY_RANK) + integer, allocatable :: gridToFieldMap(:) + integer, allocatable :: localElementCount(:) + logical :: has_vertical + + gridToFieldMap = [0, 0] + has_vertical = .TRUE. + localElementCount = [3, 5, 7] + expected = [1, localElementCount(1), product(localElementCount(2:))] + actual = get_fptr_shape(gridToFieldMap, localElementCount, has_vertical) + @assert_that(GENERIC_MESSAGE, actual, is(equal_to(expected))) + + end subroutine test_get_fptr_shape_vert_ungrid + + @Test + subroutine test_get_fptr_shape_2D_ungrid() + integer :: expected(ARRAY_RANK), actual(ARRAY_RANK) + integer, allocatable :: gridToFieldMap(:) + integer, allocatable :: localElementCount(:) + logical :: has_vertical + + has_vertical = .FALSE. + gridToFieldMap = [1, 2] + localElementCount = [3, 5, 7, 11] + expected = [product(localElementCount(1:2)), 1, product(localElementCount(3:))] + actual = get_fptr_shape(gridToFieldMap, localElementCount, has_vertical) + @assert_that(GENERIC_MESSAGE, actual, is(equal_to(expected))) + + end subroutine test_get_fptr_shape_2D_ungrid + + @Test + subroutine test_get_fptr_shape_wrong_order_raise_exception() + integer :: expected(ARRAY_RANK), actual(ARRAY_RANK) + integer, allocatable :: gridToFieldMap(:) + integer, allocatable :: localElementCount(:) + logical :: has_vertical + integer :: status + + gridToFieldMap = [4, 5] + has_vertical = .TRUE. + localElementCount = [2, 3, 5, 7, 11] + expected = [product(localElementCount(4:5)), localElementCount(3), product(localElementCount(1:2))] + ! This tests throws an Exception for improper input arguments. + ! In other words, the improper input arguments ARE the point. + actual = get_fptr_shape(gridToFieldMap, localElementCount, has_vertical, rc=status) + @assertExceptionRaised() + + end subroutine test_get_fptr_shape_wrong_order_raise_exception + +end module Test_FieldCondensedArray_private diff --git a/generic3g/Generic3g.F90 b/generic3g/Generic3g.F90 index 46fa1f9f5482..79527a2934ef 100644 --- a/generic3g/Generic3g.F90 +++ b/generic3g/Generic3g.F90 @@ -10,5 +10,5 @@ module Generic3g use mapl3g_GriddedComponentDriver use mapl3g_UserSetServices use mapl3g_ESMF_HConfigUtilities, only: MAPL_HConfigMatch - use mapl3g_output_info + use mapl3g_FieldDimensionInfo end module Generic3g diff --git a/gridcomps/History3G/HistoryCollectionGridComp_private.F90 b/gridcomps/History3G/HistoryCollectionGridComp_private.F90 index 25d89ff53079..90177190e2b5 100644 --- a/gridcomps/History3G/HistoryCollectionGridComp_private.F90 +++ b/gridcomps/History3G/HistoryCollectionGridComp_private.F90 @@ -10,8 +10,8 @@ module mapl3g_HistoryCollectionGridComp_private use MAPL_NewArthParserMod, only: parser_variables_in_expression use MAPL_TimeStringConversion use MAPL_BaseMod, only: MAPL_UnpackTime - use mapl3g_output_info, only: get_num_levels, get_vertical_dim_spec_names - use mapl3g_output_info, only: get_vertical_dim_spec_name, get_ungridded_dims + use mapl3g_FieldDimensionInfo, only: get_num_levels, get_vertical_dim_spec_names + use mapl3g_FieldDimensionInfo, only: get_vertical_dim_spec_name, get_ungridded_dims use mapl3g_UngriddedDims use gFTL2_StringSet diff --git a/gridcomps/History3G/tests/CMakeLists.txt b/gridcomps/History3G/tests/CMakeLists.txt index 431cdc92d582..1a298effd79c 100644 --- a/gridcomps/History3G/tests/CMakeLists.txt +++ b/gridcomps/History3G/tests/CMakeLists.txt @@ -3,7 +3,6 @@ set(MODULE_DIRECTORY "${esma_include}/MAPL.history3g.tests") set (test_srcs Test_HistoryGridComp.pf Test_HistoryCollectionGridComp.pf - Test_OutputInfo.pf ) add_pfunit_ctest(MAPL.history3g.tests diff --git a/shared/CMakeLists.txt b/shared/CMakeLists.txt index 3668b6d60808..34baf28f4e11 100644 --- a/shared/CMakeLists.txt +++ b/shared/CMakeLists.txt @@ -29,6 +29,7 @@ set (srcs ShaveMantissa.c MAPL_Sleep.F90 MAPL_CF_Time.F90 + MAPL_ESMF_InfoKeys.F90 # Fortran submodules Interp/Interp.F90 Interp/Interp_implementation.F90 Shmem/Shmem.F90 Shmem/Shmem_implementation.F90 diff --git a/base/MAPL_ESMF_InfoKeys.F90 b/shared/MAPL_ESMF_InfoKeys.F90 similarity index 100% rename from base/MAPL_ESMF_InfoKeys.F90 rename to shared/MAPL_ESMF_InfoKeys.F90