concurrent_esmf_test.F90

! This code is part of the SCHISM-ESMF interface. It is a main
! program for running the SCHISM component concurrent to a
! dummy component (in src/model/atmosphere_cmi_esmf.F90).
!
! @copyright (C) 2018, 2019, 2020-2021 Helmholtz-Zentrum Geesthacht
! @author Richard Hofmeister richard.hofmeister@hereon.de
! @author Carsten Lemmen carsten.lemmen@hereon.de
!
! @license Apache License, Version 2.0 (the "License");
! you may not use this file except in compliance with the License.
! You may obtain a copy of the License at
!
!                 http://www.apache.org/licenses/LICENSE-2.0
!
! Unless required by applicable law or agreed to in writing, software
! distributed under the License is distributed on an "AS IS" BASIS,
! WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
! See the License for the specific language governing permissions and
! limitations under the License.
!
! This version accounts for halo(ghost) zone, because ESMF by default
! partitions among nodes instead of elements

#define ESMF_CONTEXT  line=__LINE__,file=ESMF_FILENAME,method=ESMF_METHOD
#define ESMF_ERR_PASSTHRU msg="SCHISM subroutine call returned error"
#undef ESMF_FILENAME
#define ESMF_FILENAME "concurrent_esmf_test.F90"

#define _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(X) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, ESMF_CONTEXT, rcToReturn=X)) call ESMF_Finalize(rc=localrc, endflag=ESMF_END_ABORT)

#define ESMF_METHOD "main"
program main

  use esmf
  use schism_esmf_cap, only: schismSetServices => SetServices
  use atmosphere_cmi_esmf, only: atmosSetServices => SetServices

  implicit none

  interface
    function clockCreateFrmParam(filename, rc)
      use esmf
      character(len=ESMF_MAXSTR), intent(in) :: filename
      integer(ESMF_KIND_I4), intent(out) :: rc
      type(ESMF_Clock) :: clockCreateFrmParam
    end function clockCreateFrmParam
  end interface

  type(ESMF_GridComp) :: schism_component
  type(ESMF_GridComp) :: atmos_component

  type(ESMF_State) :: schism_import, schism_export
  type(ESMF_State) :: atmos_import, atmos_export

  type(ESMF_TimeInterval) :: timestep
  type(ESMF_Time) :: start_time, stop_time
  type(ESMF_Clock) :: clock

  type(ESMF_Field) :: field, field_in, field_out
  type(ESMF_RouteHandle) :: routehandle_air2sea, routehandle_sea2air
  type(ESMF_Vm) :: vm

  integer(ESMF_KIND_I4) :: rc, petCount, i, inum, localrc
  integer, allocatable :: petlist_schism(:), petlist_atmos(:)
  real(ESMF_KIND_R8), pointer :: ptr1d(:)
  logical :: isPresent
  character(len=ESMF_MAXSTR) :: filename

  call ESMF_Initialize(defaultCalKind=ESMF_CALKIND_GREGORIAN, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  ! Inquire the parallel environment about available
  ! resources, and partition the environment to use
  ! all but one PET for SCHISM, and last core for the
  ! dummy atmosphere component

  call ESMF_VMGetGlobal(vm, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  call ESMF_VMGet(vm, petCount=petCount, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  allocate (petlist_schism(max(1, petcount - 1)), stat=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  do i = 1, max(1, petcount - 1)
    petlist_schism(i) = i - 1 !0 based
  end do
  allocate (petlist_atmos(1))
  petlist_atmos(1) = petcount - 1

  ! Create both components on their respective parallel
  ! environment provided by each petList, then register
  ! the components entry points.

  schism_component = ESMF_GridCompCreate(name='schism_component', &
                                         petList=petlist_schism, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  atmos_component = ESMF_GridCompCreate(name='atmosphere_component', &
                                        petList=petlist_atmos, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  call ESMF_GridCompSetServices(schism_component, &
                                schismSetServices, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  call ESMF_GridCompSetServices(atmos_component, &
                                atmosSetServices, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  ! Create states for exchange of information between
  ! components.
  schism_export = ESMF_StateCreate(name='schism export state', rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  schism_import = ESMF_StateCreate(name='schism import state', rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  atmos_export = ESMF_StateCreate(name='atmos export state', rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  atmos_import = ESMF_StateCreate(name='atmos import state', rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  ! Create a clock, this should ideally be done with a configuration
  ! file but is hardcoded here for now. The timestep defined here
  ! is the coupling timestep.

  ! call ESMF_TimeIntervalSet(timestep, s=3600, rc=localrc)
  ! _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)
  !
  ! call ESMF_TimeSet(start_time, yy=2018, mm=4, dd=10, rc=localrc)
  ! _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)
  !
  ! stop_time = start_time + 2*timestep
  !
  ! clock = ESMF_ClockCreate(timestep, start_time, stopTime=stop_time, &
  !   name='main clock',rc=localrc)
  ! _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  filename = './global.nml'
  clock = clockCreateFrmParam(filename, localrc)

  ! Initialize SCHISM
  call ESMF_GridCompInitialize(schism_component, &
                               importState=schism_import, &
                               exportState=schism_export, clock=clock, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  call ESMF_GridCompInitialize(atmos_component, importState=atmos_import, &
                               exportState=atmos_export, clock=clock, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  ! Make sure that all states are reconciled across the
  ! entire VM
  call ESMF_StateReconcile(schism_import, vm=vm, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  call ESMF_StateReconcile(atmos_export, vm=vm, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  ! Within the schism component, the following fields are
  ! defined for import and export (y-component not implemented yet)
  call ESMF_StateGet(schism_import, 'wind_x-velocity_in_10m_height', &
                     field=field_out, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  call ESMF_StateGet(atmos_export, 'wind_x-velocity', &
                     field=field_in, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  ! Precompute the weights
  call ESMF_FieldRegridStore(field_in, field_out, &
                             regridMethod=ESMF_REGRIDMETHOD_BILINEAR, &
                             routehandle=routehandle_air2sea, unmappedaction=ESMF_UNMAPPEDACTION_IGNORE, &
                             rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  ! Loop over coupling timesteps until stopTime
  do while (.not. (ESMF_ClockIsStopTime(clock)))

    !> Directly manipulate the fields from import and export
    !> states.  In less basic applications, this should be
    !> handled by a mediator component.
    call ESMF_FieldRegrid(field_in, field_out, &
                          routeHandle=routehandle_air2sea, rc=localrc)
    _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

    call ESMF_GridCompRun(atmos_component, importState=atmos_import, &
                          exportState=atmos_export, clock=clock, rc=localrc)
    _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

    call ESMF_GridCompRun(schism_component, importState=schism_import, &
                          exportState=schism_export, clock=clock, rc=localrc)
    _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

    call ESMF_ClockAdvance(clock, rc=localrc)
    _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)
  end do

  !> Clean up
  call ESMF_GridCompFinalize(schism_component, importState=schism_import, &
                             exportState=schism_export, clock=clock, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  call ESMF_ClockDestroy(clock, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  call ESMF_StateDestroy(schism_import, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  call ESMF_StateDestroy(schism_export, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  call ESMF_StateDestroy(atmos_import, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  call ESMF_StateDestroy(atmos_export, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  call ESMF_GridCompDestroy(schism_component, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  call ESMF_GridCompDestroy(atmos_component, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  call ESMF_Finalize(rc=localrc)

end program main

function clockCreateFrmParam(filename, rc) result(clock)

  use esmf
  implicit none

  character(len=ESMF_MAXSTR), intent(in) :: filename
  integer(ESMF_KIND_I4), intent(out) :: rc
  type(ESMF_Clock) :: clock

  logical :: isPresent
  integer(ESMF_KIND_I4) :: unit, localrc
  type(ESMF_Time) :: stopTime, startTime
  type(ESMF_TimeInterval) :: timeStep

  integer(ESMF_KIND_I4) :: start_year = 2000, start_month = 1, start_day = 1
  integer(ESMF_KIND_I4) :: start_hour = 0, rnday = 2 !rnday in hours
  namelist /global/ start_year, start_month, start_day, start_hour, rnday

  inquire (file=filename, exist=isPresent)
  if (isPresent) then

    call ESMF_UtilIOUnitGet(unit, rc=localrc)
    _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

    open (unit, file=filename, iostat=localrc)
    _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

    read (unit, nml=global, iostat=localrc)
    _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

    close (unit)
  end if

  ! Set day as timestep temporarily to count later to stop time
  call ESMF_TimeSet(startTime, yy=start_year, mm=start_month, dd=start_day, &
                    h=start_hour, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  call ESMF_TimeIntervalSet(timeStep, h=rnday, rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

  stopTime = startTime + timeStep

  ! Only now define the coupling timestep as fraction of full timeStep (24
  ! coupling steps in total here)
  timeStep = timeStep / 24

  clock = ESMF_ClockCreate(timeStep, startTime, stopTime=stopTime, &
                           name='main clock', rc=localrc)
  _SCHISM_LOG_AND_FINALIZE_ON_ERROR_(rc)

end function clockCreateFrmParam