! $Id: CouplerMod.F90,v 1.13 2011/06/30 05:58:24 theurich Exp $ ! !------------------------------------------------------------------------- !BOP ! ! \subsection{Source for 2-way Coupler Component CouplerMod.F90} ! ! !DESCRIPTION: ! The Coupler Component provides two-way coupling between the Injector ! and FlowSolver Models. During initialization this Component is ! responsible for setting that data "is needed" from the export state ! of each model. In its run routine it calls route to transfer the ! needed data directly from one Component's export state to the other ! Component's import state. ! ! !EOP module XduneXbeachCouplerMod ! ESMF Framework module - defines ESMF data types and procedures use ESMF use libxbeach_module use introspection_module implicit none private type(ESMF_RouteHandle),save :: fromXBeach_rh, fromXDune_rh ! Public entry point public XduneXbeachCoupler_register contains !------------------------------------------------------------------------------ !BOPI ! !IROUTINE: Xdunexbeachcoupler_register - public SetServices entry point ! !INTERFACE: subroutine XduneXbeachCoupler_register(comp, rc) ! ! !ARGUMENTS: type(ESMF_CplComp) :: comp integer, intent(out) :: rc ! ! !DESCRIPTION: ! User-supplied setservices routine. ! ! The arguments are: ! \begin{description} ! \item[comp] ! Component. ! \item[rc] ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors, ! otherwise {\tt ESMF\_FAILURE}. ! \end{description} ! !EOPI ! because none of the arguments to this subroutine will ever be optional, ! go ahead and set rc to an initial return code before using it below. ! (this makes some eager error-checking compilers happy.) rc = ESMF_FAILURE ! Register the callback routines. call ESMF_CplCompSetEntryPoint(comp, ESMF_METHOD_INITIALIZE, userRoutine=xdunexbeachcoupler_init, rc=rc) if(rc /= ESMF_SUCCESS) call ESMF_Finalize(endflag=ESMF_END_ABORT, rc=rc) call ESMF_CplCompSetEntryPoint(comp, ESMF_METHOD_RUN, userRoutine=xdunexbeachcoupler_run, rc=rc) if(rc /= ESMF_SUCCESS) call ESMF_Finalize(endflag=ESMF_END_ABORT, rc=rc) call ESMF_CplCompSetEntryPoint(comp, ESMF_METHOD_FINALIZE, userRoutine=xdunexbeachcoupler_final, rc=rc) if(rc /= ESMF_SUCCESS) call ESMF_Finalize(endflag=ESMF_END_ABORT, rc=rc) print *, "XDune-XBeach coupler module: Registered Initialize, Run, and Finalize routines" end subroutine XduneXbeachCoupler_register !------------------------------------------------------------------------------ !BOPI ! !IROUTINE: xdunexbeachcoupler_init - xdunexbeachcoupler init routine ! !INTERFACE: subroutine xdunexbeachcoupler_init(comp, importState, exportState, clock, rc) ! ! !ARGUMENTS: type(ESMF_CplComp) :: comp type(ESMF_State) :: importState, exportState type(ESMF_Clock) :: clock integer, intent(out) :: rc character(len=100) :: name, tofieldname, fromfieldname, & importstatename, exportstatename type(ESMF_Field) :: tofield, fromfield integer(ESMF_KIND_I4), pointer :: indices(:,:) real(ESMF_KIND_R8), pointer :: weights(:) ! ! !DESCRIPTION: ! User-supplied init routine. ! ! The arguments are: ! \begin{description} ! \item[comp] ! Component. ! \item[importState] ! Nested state object containing import data. ! \item[exportState] ! Nested state object containing export data. ! \item[clock] ! External clock. ! \item[rc] ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors, ! otherwise {\tt ESMF\_FAILURE}. ! \end{description} ! !EOPI call ESMF_StateGet(importstate, name=name, rc=rc) ! The import of this component is the export of XBeach/Dune component and vice versa if (trim(name) .eq. 'XBeach Export') then fromfieldname = "dzbexport" tofieldname = "dhimport" elseif (trim(name) .eq. 'XDune Export') then fromfieldname = "dhexport" tofieldname = "dzbimport" end if ! Get Fields... call ESMF_StateGet(importState, fromfieldname , field=fromfield, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, & file=__FILE__)) & call ESMF_Finalize(endflag=ESMF_END_ABORT) call ESMF_StateGet(exportState, tofieldname, field=tofield, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, & file=__FILE__)) & call ESMF_Finalize(endflag=ESMF_END_ABORT) call ESMF_FieldWrite(tofield, trim(tofieldname) // 'init' // '.nc',append=.false., rc=rc) call ESMF_FieldWrite(fromfield, trim(fromfieldname) // 'init' // '.nc',append=.false., rc=rc) ! call ESMF_FieldRegridRelease(fromXDune_rh, rc=rc) ! if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & ! line=__LINE__, & ! file=__FILE__)) & ! call ESMF_Finalize(endflag=ESMF_END_ABORT) ! call ESMF_FieldRegridRelease(fromXBeach_rh, rc=rc) ! if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & ! line=__LINE__, & ! file=__FILE__)) & ! call ESMF_Finalize(endflag=ESMF_END_ABORT) if (trim(name) .eq. 'XBeach Export') then write(*,*) 'Regridstore ', trim(fromfieldname), '->', trim(tofieldname) call ESMF_FieldRegridStore(fromfield, tofield, & unmappedaction=ESMF_UNMAPPEDACTION_IGNORE, & routehandle=fromXBeach_rh, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, & file=__FILE__)) & call ESMF_Finalize(endflag=ESMF_END_ABORT) ! Indices and weights now contain the mapping between dzbdt to dh_dt else if(trim(name) .eq. 'XDune Export') then write(*,*) 'Regridstore ', trim(fromfieldname), '->', trim(tofieldname) call ESMF_FieldRegridStore(fromfield, tofield, & unmappedaction=ESMF_UNMAPPEDACTION_IGNORE, & routehandle=fromXDune_rh, weights=weights, indices=indices, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, & file=__FILE__)) & call ESMF_Finalize(endflag=ESMF_END_ABORT) end if ! ! Local variables rc = ESMF_SUCCESS print *, "Xdunexbeachcoupler Init returning" end subroutine xdunexbeachcoupler_init !------------------------------------------------------------------------------ !BOPI ! !IROUTINE: xdunexbeachcoupler_run - xdunexbeachcoupler run routine ! !INTERFACE: subroutine xdunexbeachcoupler_run(comp, importState, exportState, clock, rc) ! ! !ARGUMENTS: type(ESMF_CplComp) :: comp type(ESMF_State) :: importState, exportState type(ESMF_Clock) :: clock integer, intent(out) :: rc type(ESMF_Grid) :: grid type(ESMF_Array) :: coordarray integer(ESMF_KIND_I4),pointer :: srcMaskValues(:), dstMaskValues(:) integer(ESMF_KIND_I4), pointer :: indices(:,:) real(ESMF_KIND_R8), pointer :: weights(:) real(ESMF_KIND_R8), allocatable :: zscl(:) ! zs at the coastline character(len=100) :: name, tofieldname, fromfieldname, & tototalfieldname, fromtotalfieldname, importstatename, exportstatename real :: tstop, tstart integer(ESMF_KIND_I8) :: advanceCount, i, j integer(ESMF_KIND_I4) :: timeslice type(ESMF_Field) :: tofield, fromfield, tototalfield, fromtotalfield, & zsexportfield, wetzexportfield, zsmeanimportfield real(ESMF_KIND_r8), pointer :: farraydblptr(:,:) integer(ESMF_KIND_i4), pointer :: farrayintptr(:,:) integer(ESMF_KIND_I4), allocatable :: lastwetz(:) integer :: nx, ny, wetzcell real(ESMF_KIND_R8) :: zscumul, zsmean ! This should be a vector.... (n locations) real(ESMF_KIND_R8), pointer :: zsmeanfdbl_f_ptr(:) ! ! !DESCRIPTION: ! User-supplied run routine. ! ! The arguments are: ! \begin{description} ! \item[comp] ! Component. ! \item[importState] ! Nested state object containing import data. ! \item[exportState] ! Nested state object containing export data. ! \item[clock] ! External clock. ! \item[rc] ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors, ! otherwise {\tt ESMF\_FAILURE}. ! \end{description} ! !EOPI rc = getintparameter('nx', nx) rc = getintparameter('ny', ny) ! Define the field names. call ESMF_StateGet(importstate, name=name, rc=rc) ! The import of this component is the export of XBeach/Dune component and vice versa if (trim(name) .eq. 'XBeach Export') then fromfieldname = "dzbexport" tofieldname = "dhimport" fromtotalfieldname = "zbexport" tototalfieldname = "himport" elseif (trim(name) .eq. 'XDune Export') then fromfieldname = "dhexport" tofieldname = "dzbimport" fromtotalfieldname = "hexport" tototalfieldname = "zbimport" end if ! Get Fields... call ESMF_StateGet(importState, fromfieldname , field=fromfield, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, & file=__FILE__)) & call ESMF_Finalize(endflag=ESMF_END_ABORT) call ESMF_StateGet(exportState, tofieldname, field=tofield, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, & file=__FILE__)) & call ESMF_Finalize(endflag=ESMF_END_ABORT) call ESMF_StateGet(importState, fromtotalfieldname , field=fromtotalfield, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, & file=__FILE__)) & call ESMF_Finalize(endflag=ESMF_END_ABORT) call ESMF_StateGet(exportState, tototalfieldname, field=tototalfield, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, & file=__FILE__)) & call ESMF_Finalize(endflag=ESMF_END_ABORT) ! Get Fields... ! Compute the field bundle route handler (mapping) ! The routehandler could be reused but we're not because it might change on different wind directions call ESMF_ClockGet(clock, advanceCount=advanceCount, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, & file=__FILE__)) & call ESMF_Finalize(endflag=ESMF_END_ABORT) ! This is a cast from a 8byte integer to a 4 byte integer.... timeslice = transfer(advanceCount,timeslice) + 1 call ESMF_FieldWrite(fromfield, trim(fromfieldname) // '.nc',append=.false., timeslice=timeslice, rc=rc) call ESMF_FieldWrite(fromtotalfield, trim(fromtotalfieldname) // '.nc',append=.false., timeslice=timeslice, rc=rc) ! Process wetz and zsexport into -> importfield??? if (trim(name) .eq. 'XBeach Export') then call ESMF_StateGet(importState, 'zsexport', field=zsexportfield, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, & file=__FILE__)) & call ESMF_Finalize(endflag=ESMF_END_ABORT) call ESMF_FieldWrite(zsexportfield, 'zsexport.nc',append=.false., timeslice=timeslice, rc=rc) call ESMF_FieldGet(zsexportfield, localDe=0, farrayPtr=farraydblptr, & rc=rc) !totalLBound=ftlb, totalUBound=ftub, totalCount=ftc, if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, & file=__FILE__)) & call ESMF_Finalize(endflag=ESMF_END_ABORT) call ESMF_StateGet(importState, 'wetzexport', field=wetzexportfield, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, & file=__FILE__)) & call ESMF_Finalize(endflag=ESMF_END_ABORT) call ESMF_FieldWrite(wetzexportfield, 'wetzexport.nc',append=.false., timeslice=timeslice, rc=rc) call ESMF_FieldGet(wetzexportfield, localDe=0, farrayPtr=farrayintptr, & rc=rc) !totalLBound=ftlb, totalUBound=ftub, totalCount=ftc, if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, & file=__FILE__)) & call ESMF_Finalize(endflag=ESMF_END_ABORT) allocate(lastwetz(nx+1)) allocate(zscl(nx+1)) lastwetz = 0 zscl = 0 ! Searching for waterlevel at seaside waterline do i=1,(nx+1) ! find the first dry cell do j=2,(ny+1) wetzcell = farrayintptr(i,j) if(wetzcell .eq. 0) exit enddo lastwetz(i) = j ! save the waterlevel in the last wet cell zscl(i) = farraydblptr(i,j-1) enddo zscumul = 0 do i=1,(nx+1) zscumul = zscumul + zscl(i) enddo ! calculate mean waterlevel at waterline zsmean = zscumul / (nx+1) call ESMF_StatePrint(exportState, rc=rc) ! Store in import field of dune (in the exportstate of this coupler) call ESMF_StateGet(exportState, "zsmeanimport" , field=zsmeanimportfield, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, & file=__FILE__)) & call ESMF_Finalize(endflag=ESMF_END_ABORT) call ESMF_FieldGet(zsmeanimportfield, farrayPtr=zsmeanfdbl_f_ptr, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, & file=__FILE__)) & call ESMF_Finalize(endflag=ESMF_END_ABORT) zsmeanfdbl_f_ptr = zsmean call ESMF_FieldWrite(zsmeanimportfield, 'zsmeanimport.nc',append=.false., timeslice=timeslice, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, & file=__FILE__)) & call ESMF_Finalize(endflag=ESMF_END_ABORT) end if ! if (trim(name) .eq. 'XBeach Export') then ! call ESMF_FieldValidate(tofield, rc=rc) ! if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & ! line=__LINE__, & ! file=__FILE__)) & ! call ESMF_Finalize(endflag=ESMF_END_ABORT) ! call ESMF_FieldValidate(fromfield, rc=rc) ! if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & ! line=__LINE__, & ! file=__FILE__)) & ! call ESMF_Finalize(endflag=ESMF_END_ABORT) ! write(*,*) 'Regridstore ', trim(fromfieldname), '->', trim(tofieldname) ! call ESMF_FieldRegridStore(fromfield, tofield, & ! unmappedaction=ESMF_UNMAPPEDACTION_IGNORE, & ! routehandle=fromXBeach_rh, rc=rc) ! if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & ! line=__LINE__, & ! file=__FILE__)) & ! call ESMF_Finalize(endflag=ESMF_END_ABORT) ! ! Indices and weights now contain the mapping between dzbdt to dh_dt ! ! else if(trim(name) .eq. 'XDune Export') then ! write(*,*) 'Regridstore ', trim(fromfieldname), '->', trim(tofieldname) ! call ESMF_FieldRegridStore(fromfield, tofield, & ! ! indices=indices, weights=weights, & ! unmappedaction=ESMF_UNMAPPEDACTION_IGNORE, & ! routehandle=fromXDune_rh, rc=rc) ! if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & ! line=__LINE__, & ! file=__FILE__)) & ! call ESMF_Finalize(endflag=ESMF_END_ABORT) ! ! Indices and weights now contain the mapping between dzbdt to dh_dt ! ! I hope this uses all the regridding options, like the weights and stuff.... ! call ESMF_FieldRegrid(fromfield, tofield, & ! fromXDune_rh, rc=rc) ! end if if (trim(name) .eq. 'XBeach Export') then ! I hope this uses all the regridding options, like the weights and stuff.... write(*,*) 'Regridding', trim(fromfieldname), '->', trim(tofieldname) call ESMF_FieldRegrid(fromfield, tofield, & fromXBeach_rh, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, & file=__FILE__)) & call ESMF_Finalize(endflag=ESMF_END_ABORT) elseif (trim(name) .eq. 'XDune Export') then ! I hope this uses all the regridding options, like the weights and stuff.... write(*,*) 'Regridding', trim(fromfieldname), '->', trim(tofieldname) call ESMF_FieldRegrid(fromfield, tofield, & fromXDune_rh, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, & file=__FILE__)) & call ESMF_Finalize(endflag=ESMF_END_ABORT) end if call ESMF_FieldWrite(tototalfield, trim(tototalfieldname) // '.nc',append=.false., timeslice=timeslice, rc=rc) call ESMF_FieldWrite(tofield, trim(tofieldname) // '.nc',append=.false., timeslice=timeslice, rc=rc) ! rc has the last error code already rc = ESMF_SUCCESS end subroutine xdunexbeachcoupler_run !------------------------------------------------------------------------------ !BOPI ! !IROUTINE: xdunexbeachcoupler_final - finalization routine ! !INTERFACE: subroutine xdunexbeachcoupler_final(comp, importState, exportState, clock, rc) ! ! !ARGUMENTS: type(ESMF_CplComp) :: comp type(ESMF_State) :: importState, exportState type(ESMF_Clock) :: clock integer, intent(out) :: rc ! ! !DESCRIPTION: ! User-supplied finalize routine. ! ! The arguments are: ! \begin{description} ! \item[comp] ! Component. ! \item[importState] ! Nested state object containing import data. ! \item[exportState] ! Nested state object containing export data. ! \item[clock] ! External clock. ! \item[rc] ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors, ! otherwise {\tt ESMF\_FAILURE}. ! \end{description} ! !EOPI print *, "Xdunexbeachcoupler Final starting" ! none of the arguments to this subroutine will ever be optional, so ! go ahead and set rc to an initial return code before using it below. ! (this makes some eager error-checking compilers happy.) rc = ESMF_FAILURE ! Release some things here? rc = ESMF_SUCCESS print *, "Xdunexbeachcoupler Final returning" end subroutine xdunexbeachcoupler_final end module XdunexbeachcouplerMod