Test_ModelVerticalGrid.pf Source File

#include "MAPL_TestErr.h"
#include "MAPL_Generic.h"

! We use ESMF_TestMethod rather than basic TestMethod just in case
! there are any implied barriers is the ESMF construction in these
! tests.  E.g., if we end up needing to create nested grid comps.
! Almost certainly, is unnecessary.

module Test_ModelVerticalGrid

   use mapl_ErrorHandling
   use mapl3g_VerticalDimSpec
   use mapl3g_ModelVerticalGrid
   use mapl3g_StateRegistry
   use mapl3g_VariableSpec
   use mapl3g_VirtualConnectionPt
   use mapl3g_ActualConnectionPt
   use mapl3g_StateItemSpec
   use mapl3g_StateItemExtension
   use mapl3g_GriddedComponentDriver
   use mapl3g_ComponentDriver
   use mapl3g_ComponentDriverVector
   use mapl3g_ComponentDriverPtrVector
   use mapl3g_MultiState
   use mapl3g_make_ItemSpec
   use mapl3g_geom_mgr
   use mapl3g_CouplerPhases, only: GENERIC_COUPLER_UPDATE
   use esmf
   ! testing framework
   use ESMF_TestMethod_mod
   use funit

   implicit none (type, external)

   integer, parameter :: IM=6, JM=7, LM=3

   ! Trying to avoid a complex test fixture
   type(StateRegistry), target :: r

contains

   subroutine setup(var_name, vgrid, rc)
      character(*), intent(in) :: var_name
      type(ModelVerticalGrid), intent(out) :: vgrid
      integer, intent(out) :: rc

      type(VerticalDimSpec) :: vertical_dim_spec
      type(ESMF_Geom) :: geom
      type(VirtualConnectionPt) :: v_pt
      type(VariableSpec) :: var_spec
      class(StateItemSpec), allocatable :: fld_spec
      type(StateItemExtension), pointer :: extension
      class(StateItemSpec), pointer :: spec
      integer :: status

      select case (var_name)
      case ("PLE")
         vertical_dim_spec = VERTICAL_DIM_EDGE
      case ("PL")
         vertical_dim_spec = VERTICAL_DIM_CENTER
      case default
         _FAIL("var_name should be one of PLE/PL, not" // trim(var_name))
      end select

      rc = 0
      ! Inside user "set_geom" phase.
      geom = make_geom(_RC)
      vgrid = ModelVerticalGrid(standard_name="air_pressure", units="hPa", num_levels=LM)
      call vgrid%add_short_name(edge="PLE", center="PL")

      ! inside OuterMeta
      r = StateRegistry("dyn")
      call vgrid%set_registry(r) ! MAPL_SetVerticalGrid(...)

      v_pt = VirtualConnectionPt(state_intent="export", short_name=var_name)
      var_spec = VariableSpec(&
           short_name=var_name, &
           state_intent=ESMF_STATEINTENT_EXPORT, &
           standard_name="air_pressure", &
           units="hPa", &
           vertical_dim_spec=vertical_dim_spec, &
           default_value=3.)
      allocate(fld_spec, source=make_itemSpec(var_spec, r, rc=status))
      _VERIFY(status)
      call fld_spec%set_geometry(geom=geom, vertical_grid=vgrid, _RC)

      call r%add_primary_spec(v_pt, fld_spec)

      extension => r%get_primary_extension(v_pt, _RC)
      spec => extension%get_spec()
      call spec%set_active()
      call spec%create(_RC)
      call spec%allocate(_RC)
   end subroutine setup

   function make_geom(rc) result(geom)
      integer, intent(out) :: rc
      type(ESMF_Geom) :: geom
      integer :: status
      type(ESMF_HConfig) :: hconfig
      type(GeomManager), pointer :: geom_mgr
      type(MaplGeom), pointer :: mapl_geom

      rc = 0
      geom_mgr => get_geom_manager()
      hconfig = ESMF_HConfigCreate(content="{class: latlon, im_world: 6, jm_world: 7, pole: PC, dateline: DC}", _RC)
      mapl_geom => geom_mgr%get_mapl_geom(hconfig, _RC)
      geom = mapl_geom%get_geom()
   end function make_geom

   @test
   subroutine test_num_levels()
      type(ModelVerticalGrid) :: vgrid
      integer :: num_levels

      num_levels = 10
      vgrid = ModelVerticalGrid(standard_name="height", units="m", num_levels=num_levels)
      @assert_that(vgrid%get_num_levels(), is(num_levels))
   end subroutine test_num_levels

   @test(type=ESMF_TestMethod, npes=[1])
   subroutine test_created_fields_have_num_levels(this)
      class(ESMF_TestMethod), intent(inout) :: this
      type(ModelVerticalGrid) :: vgrid
      integer :: rank
      integer, allocatable :: localElementCount(:)
      type(VirtualConnectionPt) :: ple_pt
      class(StateItemSpec), pointer :: spec
      type(MultiState) :: multi_state
      type(StateItemExtension), pointer :: extension
      type(ESMF_Field) :: ple
      integer :: rc, status

      call setup("PLE", vgrid, _RC)

      ple_pt = VirtualConnectionPt(state_intent="export", short_name="PLE")
      extension => r%get_primary_extension(ple_pt, _RC)
      spec => extension%get_spec()

      multi_state = MultiState()
      call spec%add_to_state(multi_state, ActualConnectionPt(ple_pt), _RC)
      call ESMF_StateGet(multi_state%exportState, itemName="PLE", field=ple, _RC)
      call ESMF_FieldGet(ple, rank=rank, _RC)
      allocate(localElementCount(rank))
      call ESMF_FieldGet(ple, localElementCount=localElementCount, _RC)
      @assert_that(localElementCount, is(equal_to([IM,JM,LM+1])))
      _UNUSED_DUMMY(this)
   end subroutine test_created_fields_have_num_levels

   @test(type=ESMF_TestMethod, npes=[1])
   ! Request the specific coordinate corresponding particular geom/unit.
   ! In this case we start with one that already exists.  A later test
   ! should force extensions.
   subroutine test_get_coordinate_field_simple(this)
      class(ESMF_TestMethod), intent(inout) :: this
      type(ModelVerticalGrid) :: vgrid
      type(GriddedComponentDriver), pointer :: coupler
      type(ESMF_Field) :: vcoord
      type(ESMF_Geom) :: geom
      integer :: rc, status
      real(ESMF_KIND_R4), pointer :: a(:,:,:)

      call setup("PLE", vgrid, _RC)
      geom = make_geom(_RC)

      call vgrid%get_coordinate_field( &
           vcoord, coupler, &
           standard_name="air_pressure", &
           geom=geom, &
           typekind=ESMF_TYPEKIND_R4, &
           units="hPa", &
           vertical_dim_spec=VERTICAL_DIM_EDGE, &
           _RC)
      @assert_that(associated(coupler), is(false()))

      call ESMF_FieldGet(vcoord, fArrayPtr=a, _RC)
      @assert_that(a, every_item(is(equal_to(3.))))
      _UNUSED_DUMMY(this)
   end subroutine test_get_coordinate_field_simple

   @test(type=ESMF_TestMethod, npes=[1])
   ! Request the specific coordinate corresponding particular geom/unit.
   ! Here we request different units which should return a coordinate
   ! scaled by 100 (hPa = 100 Pa)
   subroutine test_get_coordinate_field_change_units_edge(this)
      class(ESMF_TestMethod), intent(inout) :: this
      type(ModelVerticalGrid) :: vgrid
      type(ESMF_Field) :: vcoord
      type(ESMF_Geom) :: geom
      integer :: status
      real(ESMF_KIND_R4), pointer :: a(:,:,:)
      type(ComponentDriverPtrVector) :: couplers
      type(ComponentDriverPtr) :: driver
      type(GriddedComponentDriver), pointer :: coupler
      integer :: i, rc

      call setup("PLE", vgrid, _RC)
      geom = make_geom(_RC)

      call vgrid%get_coordinate_field( &
           vcoord, coupler, &
           standard_name="air_pressure", &
           geom=geom, &
           typekind=ESMF_TYPEKIND_R4, &
           units="Pa", &
           vertical_dim_spec=VERTICAL_DIM_EDGE, &
           _RC)
      @assert_that(associated(coupler), is(true()))

      call r%allocate(_RC)

      call ESMF_FieldGet(vcoord, fArrayPtr=a, _RC)
      ! usually update is called on imports, but here we don't have an import handy,
      ! so we force updates on all export couplers in registry r.
      couplers = r%get_export_couplers()
      do i = 1, couplers%size()
         driver = couplers%of(i)
         call driver%ptr%initialize(_RC)
         call driver%ptr%run(phase_idx=GENERIC_COUPLER_UPDATE, _RC)
      end do
      @assert_that(shape(a), is(equal_to([IM, JM, LM+1])))
      @assert_that(a, every_item(is(equal_to(300.))))
      _UNUSED_DUMMY(this)
   end subroutine test_get_coordinate_field_change_units_edge

   @test(type=ESMF_TestMethod, npes=[1])
   ! Request the specific coordinate corresponding particular geom/unit.
   ! Here we request different units which should return a coordinate
   ! scaled by 100 (hPa = 100 Pa)
   subroutine test_get_coordinate_field_change_units_center(this)
      class(ESMF_TestMethod), intent(inout) :: this
      type(ModelVerticalGrid) :: vgrid
      type(ESMF_Field) :: vcoord
      type(ESMF_Geom) :: geom
      integer :: status
      real(ESMF_KIND_R4), pointer :: a(:,:,:)
      type(ComponentDriverPtrVector) :: couplers
      type(ComponentDriverPtr) :: driver
      type(GriddedComponentDriver), pointer :: coupler
      integer :: i, rc

      call setup("PL", vgrid, _RC)
      geom = make_geom(_RC)

      call vgrid%get_coordinate_field( &
           vcoord, coupler, &
           standard_name="air_pressure", &
           geom=geom, &
           typekind=ESMF_TYPEKIND_R4, units="Pa", &
           vertical_dim_spec=VERTICAL_DIM_CENTER, &
           _RC)
      @assert_that(associated(coupler), is(true()))

      call r%allocate(_RC)

      call ESMF_FieldGet(vcoord, fArrayPtr=a, _RC)
      ! usually update is called on imports, but here we don't have an import handy,
      ! so we force updates on all export couplers in registry r.
      couplers = r%get_export_couplers()
      do i = 1, couplers%size()
         driver = couplers%of(i)
         call driver%ptr%initialize(_RC)
         call driver%ptr%run(phase_idx=GENERIC_COUPLER_UPDATE, _RC)
      end do
      @assert_that(shape(a), is(equal_to([IM, JM, LM])))
      @assert_that(a, every_item(is(equal_to(300.))))
      _UNUSED_DUMMY(this)
   end subroutine test_get_coordinate_field_change_units_center

end module Test_ModelVerticalGrid