From 23da80a9b5fcf2b0097de6c7e77e525594a18171 Mon Sep 17 00:00:00 2001 From: Purnendu Chakraborty Date: Tue, 26 Nov 2024 13:20:33 -0500 Subject: [PATCH] Added VerticalGrid::is_identical_to to check if the dst grid is identical to this --- generic3g/specs/FieldSpec.F90 | 45 +------------- generic3g/vertical/BasicVerticalGrid.F90 | 20 +++---- .../vertical/FixedLevelsVerticalGrid.F90 | 30 ++++++++++ generic3g/vertical/MirrorVerticalGrid.F90 | 13 +++++ generic3g/vertical/ModelVerticalGrid.F90 | 58 ++++++++++++++++++- generic3g/vertical/VerticalGrid.F90 | 8 +++ 6 files changed, 118 insertions(+), 56 deletions(-) diff --git a/generic3g/specs/FieldSpec.F90 b/generic3g/specs/FieldSpec.F90 index ee02a3b12492..87cb24d605ba 100644 --- a/generic3g/specs/FieldSpec.F90 +++ b/generic3g/specs/FieldSpec.F90 @@ -872,7 +872,6 @@ subroutine adapt_vertical_grid(this, spec, action, rc) end subroutine adapt_vertical_grid logical function adapter_match_vertical_grid(this, spec, rc) result(match) - class(VerticalGridAdapter), intent(in) :: this class(StateItemSpec), intent(in) :: spec integer, optional, intent(out) :: rc @@ -882,52 +881,10 @@ logical function adapter_match_vertical_grid(this, spec, rc) result(match) match = .false. select type (spec) type is (FieldSpec) - match = same_vertical_grid(spec%vertical_grid, this%vertical_grid, _RC) + match = spec%vertical_grid%is_identical_to(this%vertical_grid) end select _RETURN(_SUCCESS) - - contains - - logical function same_vertical_grid(src_grid, dst_grid, rc) - class(VerticalGrid), intent(in) :: src_grid - class(VerticalGrid), allocatable, intent(in) :: dst_grid - integer, optional, intent(out) :: rc - - same_vertical_grid = .false. - if (.not. allocated(dst_grid)) then - same_vertical_grid = .true. - _RETURN(_SUCCESS) ! mirror grid - end if - - same_vertical_grid = src_grid%same_id(dst_grid) - if (same_vertical_grid) then - _RETURN(_SUCCESS) - end if - - select type(src_grid) - type is(BasicVerticalGrid) - select type(dst_grid) - type is(BasicVerticalGrid) - same_vertical_grid = (src_grid%get_num_levels() == dst_grid%get_num_levels()) - class default - _FAIL("not implemented yet") - end select - type is(FixedLevelsVerticalGrid) - select type(dst_grid) - type is(FixedLevelsVerticalGrid) - same_vertical_grid = (src_grid == dst_grid) - class default - same_vertical_grid = .false. - end select - class default ! ModelVerticalGrid - same_vertical_grid = .false. - ! _FAIL("not implemented yet") - end select - - _RETURN(_SUCCESS) - end function same_vertical_grid - end function adapter_match_vertical_grid function new_TypekindAdapter(typekind) result(typekind_adapter) diff --git a/generic3g/vertical/BasicVerticalGrid.F90 b/generic3g/vertical/BasicVerticalGrid.F90 index 54d2da7bbcf4..a823ec623d2a 100644 --- a/generic3g/vertical/BasicVerticalGrid.F90 +++ b/generic3g/vertical/BasicVerticalGrid.F90 @@ -23,6 +23,7 @@ module mapl3g_BasicVerticalGrid procedure :: get_num_levels procedure :: get_coordinate_field procedure :: can_connect_to + procedure :: is_identical_to procedure :: write_formatted end type BasicVerticalGrid @@ -81,18 +82,17 @@ logical function can_connect_to(this, dst, rc) class(VerticalGrid), intent(in) :: dst integer, optional, intent(out) :: rc - select type(dst) - type is (BasicVerticalGrid) - can_connect_to = (this%get_num_levels() == dst%get_num_levels()) - type is (MirrorVerticalGrid) - can_connect_to = .true. - class default - _FAIL("BasicVerticalGrid can only connect to BasicVerticalGrid, or MirrorVerticalGrid") - end select - - _RETURN(_SUCCESS) + _FAIL("BasicVerticalGrid::can_connect_to - NOT implemented yet") end function can_connect_to + logical function is_identical_to(this, that, rc) + class(BasicVerticalGrid), intent(in) :: this + class(VerticalGrid), allocatable, intent(in) :: that + integer, optional, intent(out) :: rc + + _FAIL("BasicVerticalGrid::is_identical_to - NOT implemented yet") + end function is_identical_to + elemental logical function equal_to(a, b) type(BasicVerticalGrid), intent(in) :: a, b equal_to = a%num_levels == b%num_levels diff --git a/generic3g/vertical/FixedLevelsVerticalGrid.F90 b/generic3g/vertical/FixedLevelsVerticalGrid.F90 index cd59462dae99..67e02577351e 100644 --- a/generic3g/vertical/FixedLevelsVerticalGrid.F90 +++ b/generic3g/vertical/FixedLevelsVerticalGrid.F90 @@ -27,6 +27,7 @@ module mapl3g_FixedLevelsVerticalGrid procedure :: get_num_levels procedure :: get_coordinate_field procedure :: can_connect_to + procedure :: is_identical_to procedure :: write_formatted end type FixedLevelsVerticalGrid @@ -118,6 +119,35 @@ logical function can_connect_to(this, dst, rc) _RETURN(_SUCCESS) end function can_connect_to + logical function is_identical_to(this, that, rc) + class(FixedLevelsVerticalGrid), intent(in) :: this + class(VerticalGrid), allocatable, intent(in) :: that + integer, optional, intent(out) :: rc + + logical :: same_id + + is_identical_to = .false. + + ! Mirror grid + if (.not. allocated(that)) then + is_identical_to = .true. + _RETURN(_SUCCESS) ! mirror grid + end if + + ! Same id + is_identical_to = this%same_id(that) + if (is_identical_to) then + _RETURN(_SUCCESS) + end if + + select type(that) + type is(FixedLevelsVerticalGrid) + is_identical_to = (this == that) + end select + + _RETURN(_SUCCESS) + end function is_identical_to + subroutine write_formatted(this, unit, iotype, v_list, iostat, iomsg) class(FixedLevelsVerticalGrid), intent(in) :: this integer, intent(in) :: unit diff --git a/generic3g/vertical/MirrorVerticalGrid.F90 b/generic3g/vertical/MirrorVerticalGrid.F90 index 98f04424f815..2c6048962a87 100644 --- a/generic3g/vertical/MirrorVerticalGrid.F90 +++ b/generic3g/vertical/MirrorVerticalGrid.F90 @@ -26,6 +26,7 @@ module mapl3g_MirrorVerticalGrid procedure :: get_num_levels procedure :: get_coordinate_field procedure :: can_connect_to + procedure :: is_identical_to procedure :: write_formatted end type MirrorVerticalGrid @@ -81,6 +82,18 @@ logical function can_connect_to(this, dst, rc) _UNUSED_DUMMY(dst) end function can_connect_to + logical function is_identical_to(this, that, rc) + class(MirrorVerticalGrid), intent(in) :: this + class(VerticalGrid), allocatable, intent(in) :: that + integer, optional, intent(out) :: rc + + is_identical_to = .false. + + _RETURN(_SUCCESS) + _UNUSED_DUMMY(this) + _UNUSED_DUMMY(that) + end function is_identical_to + subroutine write_formatted(this, unit, iotype, v_list, iostat, iomsg) class(MirrorVerticalGrid), intent(in) :: this integer, intent(in) :: unit diff --git a/generic3g/vertical/ModelVerticalGrid.F90 b/generic3g/vertical/ModelVerticalGrid.F90 index 33dfe9caeb06..d290b417384c 100644 --- a/generic3g/vertical/ModelVerticalGrid.F90 +++ b/generic3g/vertical/ModelVerticalGrid.F90 @@ -35,6 +35,7 @@ module mapl3g_ModelVerticalGrid procedure :: get_num_levels procedure :: get_coordinate_field procedure :: can_connect_to + procedure :: is_identical_to procedure :: write_formatted ! subclass-specific methods @@ -48,6 +49,14 @@ module mapl3g_ModelVerticalGrid procedure new_ModelVerticalGrid_basic end interface ModelVerticalGrid + interface operator(==) + module procedure equal_ModelVerticalGrid + end interface operator(==) + + interface operator(/=) + module procedure not_equal_ModelVerticalGrid + end interface operator(/=) + ! TODO: ! - Ensure that there really is a vertical dimension @@ -179,8 +188,6 @@ logical function can_connect_to(this, dst, rc) class(VerticalGrid), intent(in) :: dst integer, optional, intent(out) :: rc - integer :: status - if (this%same_id(dst)) then can_connect_to = .true. _RETURN(_SUCCESS) @@ -198,4 +205,51 @@ logical function can_connect_to(this, dst, rc) _RETURN(_SUCCESS) end function can_connect_to + logical function is_identical_to(this, that, rc) + class(ModelVerticalGrid), intent(in) :: this + class(VerticalGrid), allocatable, intent(in) :: that + integer, optional, intent(out) :: rc + + is_identical_to = .false. + + ! Mirror grid + if (.not. allocated(that)) then + is_identical_to = .true. + _RETURN(_SUCCESS) ! mirror grid + end if + + ! Same id + is_identical_to = this%same_id(that) + if (is_identical_to) then + _RETURN(_SUCCESS) + end if + + select type(that) + type is(ModelVerticalGrid) + is_identical_to = (this == that) + end select + + _RETURN(_SUCCESS) + end function is_identical_to + + impure elemental logical function equal_ModelVerticalGrid(a, b) result(equal) + type(ModelVerticalGrid), intent(in) :: a, b + + equal = a%standard_name == b%standard_name + if (.not. equal) return + equal = (a%get_units() == b%get_units()) + if (.not. equal) return + equal = (a%num_levels == b%num_levels) + if (.not. equal) return + equal = (a%short_name_edge == b%short_name_edge) + if (.not. equal) return + equal = (a%short_name_center == b%short_name_center) + end function equal_ModelVerticalGrid + + impure elemental logical function not_equal_ModelVerticalGrid(a, b) result(not_equal) + type(ModelVerticalGrid), intent(in) :: a, b + + not_equal = .not. (a==b) + end function not_equal_ModelVerticalGrid + end module mapl3g_ModelVerticalGrid diff --git a/generic3g/vertical/VerticalGrid.F90 b/generic3g/vertical/VerticalGrid.F90 index 1fdf5c66076c..307814540b6c 100644 --- a/generic3g/vertical/VerticalGrid.F90 +++ b/generic3g/vertical/VerticalGrid.F90 @@ -15,6 +15,7 @@ module mapl3g_VerticalGrid procedure(I_get_num_levels), deferred :: get_num_levels procedure(I_get_coordinate_field), deferred :: get_coordinate_field procedure(I_can_connect_to), deferred :: can_connect_to + procedure(I_is_identical_to), deferred :: is_identical_to procedure(I_write_formatted), deferred :: write_formatted generic :: write(formatted) => write_formatted @@ -59,6 +60,13 @@ logical function I_can_connect_to(this, dst, rc) result(can_connect_to) integer, optional, intent(out) :: rc end function I_can_connect_to + logical function I_is_identical_to(this, that, rc) result(is_identical_to) + import VerticalGrid + class(VerticalGrid), intent(in) :: this + class(VerticalGrid), allocatable, intent(in) :: that + integer, optional, intent(out) :: rc + end function I_is_identical_to + subroutine I_write_formatted(this, unit, iotype, v_list, iostat, iomsg) import VerticalGrid class(VerticalGrid), intent(in) :: this