MAPL interface to MITgcm
From Maplcode.org
Contents |
[edit] Background
The MITgcm (MIT General Circulation Model) is a numerical model designed for study of the atmosphere, ocean, and climate. Its non-hydrostatic formulation enables it to simulate fluid phenomena over a wide range of scales; its adjoint capability enables it to be applied to parameter and state estimation problems. By employing fluid isomorphisms, one hydrodynamical kernel can be used to simulate flow in both the atmosphere and ocean.
The code is implemented in Fortran and its implementation is constrained by the capabilities of adjoint mode capable autpmatic differentiation tools. Maintaining compatibility with automatic differentiation precludes the widespread use of truly dynamic, instantiable programming practices. Slides dicussing the coupling of MITgcm and GEOS through MAPL can be found here. The main elements for fitting the code within MAPL (including enabling limited multiple instantiation capabilities) are described below. The steps are described in some detail below as a way to illustrate how to approach fitting a non-instantiable, legacy numerical kernel into a system like MAPL in a way that (i) can exploit a lot of the powerful features of MAPL and (ii) does not break other constraints on the numerical kernel implementation.
[edit] Plug code MIT_GEOS5PlugMod
The plug code follows fairly standard MAPL practices for a leaf component. It contains SetServices, Initialize, Run and Finalize functions.
[edit] Plug SetServices
The SetServices code registers the imports and exports, associating them with the gridded component passed down into SetServices. In SetServices for MIT_GEOS5PlugMod we use the type shown below.
! Variables for setting MAPL import/export specs
TYPE MAPL_SPEC
CHARACTER(len=ESMF_MAXSTR) :: short_name
CHARACTER(len=ESMF_MAXSTR) :: long_name
CHARACTER(len=ESMF_MAXSTR) :: units
INTEGER :: dims
INTEGER :: vlocation
END TYPE
TYPE( MAPL_SPEC ), POINTER :: imports(:)
TYPE( MAPL_SPEC ), POINTER :: exports(:)
This type is used to create arrays of MAPL array specifications used to register the plug imports and exports. These are then used in the executable portion of SetServices to register the imports and exports as follows.
[edit] Import registration
! Imports and exports specification
! ---------------------------------
nimports = 7
allocate(imports(nimports))
! ---------------------------------------------
imports = (/ &
mstate('TAUX', 'Agrid_eastward_stress_on_skin', 'N m-2', MAPL_DimsHorzOnly,MAPL_VLocationNone), &
mstate('TAUY', 'Agrid_northward_stress_on_skin', 'N m-2', MAPL_DimsHorzOnly,MAPL_VLocationNone), &
mstate('PS', 'Surface Atmospheric Pressure', 'Pa', MAPL_DimsHorzOnly,MAPL_VLocationNone), &
mstate('SWHEAT','solar_heating_rate', 'W m-2', MAPL_DimsHorzVert,MAPL_VLocationCenter), &
mstate('QFLX', 'freshwater_flux_from_skin_to_ocean', 'kg m-2 s-1',MAPL_DimsHorzOnly,MAPL_VLocationNone), &
mstate('HFLX', 'turbulent_heat_flux_from_skin_to_ocean','W m-2', MAPL_DimsHorzOnly,MAPL_VLocationNone), &
mstate('SFLX', 'salt_flux_from_skin_to_ocean', 'N m-2', MAPL_DimsHorzOnly,MAPL_VLocationNone) &
/)
DO I=1,NIMPORTS
CALL MAPL_StateAddImportSpec(GC, &
SHORT_NAME = imports(i)%short_name, &
LONG_NAME = imports(i)%long_name, &
UNITS = imports(i)%units, &
DIMS = imports(i)%dims, &
VLOCATION = imports(i)%vlocation, &
RC =status); VERIFY_(STATUS)
ENDDO
[edit] Export registration
nexports =5
allocate(exports(nexports))
exports = (/ &
mstate('US', 'top_layer_Agrid_eastward_velocity', 'm s-1', MAPL_DimsHorzOnly,MAPL_VLocationNone), &
mstate('VS', 'top_layer_Agrid_northward_velocity', 'm s-1', MAPL_DimsHorzOnly,MAPL_VLocationNone), &
mstate('TS', 'top_layer_temperature', 'K', MAPL_DimsHorzOnly,MAPL_VLocationNone), &
mstate('SS', 'top_layer_salinity', 'psu', MAPL_DimsHorzOnly,MAPL_VLocationNone), &
mstate('MASK', 'ocean mask at t-points', '1', MAPL_DimsHorzVert,MAPL_VLocationCenter) &
/)
DO I=1,5
call MAPL_StateAddExportSpec(GC, &
SHORT_NAME = exports(i)%short_name, &
LONG_NAME = exports(i)%long_name, &
UNITS = exports(i)%units, &
DIMS = exports(i)%dims, &
VLOCATION = exports(i)%vlocation, &
RC=STATUS )
VERIFY_(STATUS)
ENDDO
[edit] Plug Initialize
The Initialize routine for the MITgcm GEOS5 plug includes a collection of MITgcm state types that are used to provide limited internal state support. These types were introduced as part of the process of fitting within MAPL, in order to support instantiation. An example of one of these types is shown below
TYPE MITGCM_DYNVARS
SEQUENCE
_RL , POINTER :: uVel( :,:,:,:,:) => NULL()
_RL , POINTER :: vVel( :,:,:,:,:) => NULL()
_RL , POINTER :: wVel( :,:,:,:,:) => NULL()
_RL , POINTER :: theta( :,:,:,:,:) => NULL()
_RL , POINTER :: salt( :,:,:,:,:) => NULL()
_RL , POINTER :: guNM1( :,:,:,:,:) => NULL()
_RL , POINTER :: gvNM1( :,:,:,:,:) => NULL()
_RL , POINTER :: gtNM1( :,:,:,:,:) => NULL()
_RL , POINTER :: gsNM1( :,:,:,:,:) => NULL()
_RL , POINTER :: gU( :,:,:,:,:) => NULL()
_RL , POINTER :: gv( :,:,:,:,:) => NULL()
_RL , POINTER :: uVelD( :,:,:,:,:) => NULL()
_RL , POINTER :: vVelD( :,:,:,:,:) => NULL()
_RL , POINTER :: uNM1( :,:,:,:,:) => NULL()
_RL , POINTER :: vNM1( :,:,:,:,:) => NULL()
_RL , POINTER :: etaN( :,:,:,:) => NULL()
_RL , POINTER :: etaNM1( :,:,:,:) => NULL()
END TYPE
Its contents mirror the content These types are collected together to provide an internal state object. In Initialize are pointer to this collection of types
TYPE(MITGCM_ISTATE_CONTAINER) :: mitgcmIState(1)
is declared and passed down into an internal MITgcm routine
CALL DRIVER_INIT( mitgcmIState(1)%p, myComm)
