Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add DLND to moab coupler #6567

Merged
merged 8 commits into from
Sep 22, 2024
Merged
Show file tree
Hide file tree
Changes from 7 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
181 changes: 180 additions & 1 deletion components/data_comps/dlnd/src/dlnd_comp_mod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -29,6 +29,10 @@ module dlnd_comp_mod
use dlnd_shr_mod , only: domain_fracname ! namelist input
use dlnd_shr_mod , only: nullstr

#ifdef HAVE_MOAB
use seq_comm_mct, only : mlnid ! id of moab lnd app
use iso_c_binding
#endif
! !PUBLIC TYPES:
implicit none
private ! except
Expand Down Expand Up @@ -100,6 +104,12 @@ subroutine dlnd_comp_init(Eclock, x2l, l2x, &
scmMode, scmlat, scmlon)

! !DESCRIPTION: initialize dlnd model
#ifdef HAVE_MOAB
use iMOAB, only: iMOAB_DefineTagStorage, iMOAB_GetDoubleTagStorage, &
iMOAB_SetIntTagStorage, iMOAB_SetDoubleTagStorage, &
iMOAB_ResolveSharedEntities, iMOAB_CreateVertices, &
iMOAB_GetMeshInfo, iMOAB_UpdateMeshInfo, iMOAB_WriteMesh
iulian787 marked this conversation as resolved.
Show resolved Hide resolved
#endif
implicit none

! !INPUT/OUTPUT PARAMETERS:
Expand Down Expand Up @@ -135,6 +145,18 @@ subroutine dlnd_comp_init(Eclock, x2l, l2x, &
character(nec_len) :: nec_str ! elevation class, as character string
character(*), parameter :: domain_fracname_unset = 'null'

#ifdef HAVE_MOAB
character*400 tagname
real(R8) latv, lonv
integer iv, tagindex, ilat, ilon
real(R8), allocatable, target :: data(:)
integer(IN), pointer :: idata(:) ! temporary
real(R8), dimension(:), allocatable :: moab_vert_coords ! temporary
#ifdef MOABDEBUG
character*100 outfile, wopts
#endif
#endif

!--- formats ---
character(*), parameter :: F00 = "('(dlnd_comp_init) ',8a)"
character(*), parameter :: F0L = "('(dlnd_comp_init) ',a, l2)"
Expand Down Expand Up @@ -256,6 +278,119 @@ subroutine dlnd_comp_init(Eclock, x2l, l2x, &

call t_stopf('dlnd_initmctdom')

#ifdef HAVE_MOAB
ilat = mct_aVect_indexRA(ggrid%data,'lat')
ilon = mct_aVect_indexRA(ggrid%data,'lon')
allocate(moab_vert_coords(lsize*3))
do iv = 1, lsize
lonv = ggrid%data%rAttr(ilon, iv) * SHR_CONST_PI/180.
latv = ggrid%data%rAttr(ilat, iv) * SHR_CONST_PI/180.
rljacob marked this conversation as resolved.
Show resolved Hide resolved
moab_vert_coords(3*iv-2)=COS(latv)*COS(lonv)
moab_vert_coords(3*iv-1)=COS(latv)*SIN(lonv)
moab_vert_coords(3*iv )=SIN(latv)
enddo

! create the vertices with coordinates from MCT domain
ierr = iMOAB_CreateVertices(mlnid, lsize*3, 3, moab_vert_coords)
if (ierr .ne. 0) &
call shr_sys_abort('Error: fail to create MOAB vertices in data lnd model')

tagname='GLOBAL_ID'//C_NULL_CHAR
ierr = iMOAB_DefineTagStorage(mlnid, tagname, &
0, & ! dense, integer
iulian787 marked this conversation as resolved.
Show resolved Hide resolved
1, & ! number of components
tagindex )
if (ierr .ne. 0) &
call shr_sys_abort('Error: fail to retrieve GLOBAL_ID tag ')

! get list of global IDs for Dofs
call mct_gsMap_orderedPoints(gsMap, my_task, idata)

ierr = iMOAB_SetIntTagStorage ( mlnid, tagname, lsize, &
0, & ! vertex type
idata)
if (ierr .ne. 0) &
call shr_sys_abort('Error: fail to set GLOBAL_ID tag ')

ierr = iMOAB_ResolveSharedEntities( mlnid, lsize, idata );
if (ierr .ne. 0) &
call shr_sys_abort('Error: fail to resolve shared entities')

deallocate(moab_vert_coords)
deallocate(idata)
vijaysm marked this conversation as resolved.
Show resolved Hide resolved

ierr = iMOAB_UpdateMeshInfo( mlnid )
if (ierr .ne. 0) &
call shr_sys_abort('Error: fail to update mesh info ')

allocate(data(lsize))
ierr = iMOAB_DefineTagStorage( mlnid, "area:aream:frac:mask"//C_NULL_CHAR, &
1, & ! dense, double
1, & ! number of components
tagindex )
if (ierr > 0 ) &
call shr_sys_abort('Error: fail to create tag: area:aream:frac:mask' )

data(:) = ggrid%data%rAttr(mct_aVect_indexRA(ggrid%data,'area'),:)
tagname='area'//C_NULL_CHAR
ierr = iMOAB_SetDoubleTagStorage ( mlnid, tagname, lsize, &
0, & ! set data on vertices
data)
if (ierr > 0 ) &
call shr_sys_abort('Error: fail to get area tag ')

! set the same data for aream (model area) as area
! data(:) = ggrid%data%rAttr(mct_aVect_indexRA(ggrid%data,'aream'),:)
tagname='aream'//C_NULL_CHAR
ierr = iMOAB_SetDoubleTagStorage ( mlnid, tagname, lsize, &
0, & ! set data on vertices
data)
if (ierr > 0 ) &
call shr_sys_abort('Error: fail to set aream tag ')

data(:) = ggrid%data%rAttr(mct_aVect_indexRA(ggrid%data,'mask'),:)
tagname='mask'//C_NULL_CHAR
ierr = iMOAB_SetDoubleTagStorage ( mlnid, tagname, lsize, &
0, & ! set data on vertices
data)
if (ierr > 0 ) &
call shr_sys_abort('Error: fail to set mask tag ')

data(:) = ggrid%data%rAttr(mct_aVect_indexRA(ggrid%data,'frac'),:)
tagname='frac'//C_NULL_CHAR
ierr = iMOAB_SetDoubleTagStorage ( mlnid, tagname, lsize, &
0, & ! set data on vertices
data)
if (ierr > 0 ) &
call shr_sys_abort('Error: fail to set frac tag ')

deallocate(data)

! define tags
ierr = iMOAB_DefineTagStorage( mlnid, trim(seq_flds_x2l_fields)//C_NULL_CHAR, &
1, & ! dense, double
1, & ! number of components
tagindex )
if (ierr > 0 ) &
call shr_sys_abort('Error: fail to create seq_flds_x2l_fields tags ')

ierr = iMOAB_DefineTagStorage( mlnid, trim(seq_flds_l2x_fields)//C_NULL_CHAR, &
1, & ! dense, double
1, & ! number of components
tagindex )
if (ierr > 0 ) &
call shr_sys_abort('Error: fail to create seq_flds_l2x_fields tags ')
#ifdef MOABDEBUG
! debug test
outfile = 'LndDataMesh.h5m'//C_NULL_CHAR
wopts = ';PARALLEL=WRITE_PART'//C_NULL_CHAR !
! write out the mesh file to disk
ierr = iMOAB_WriteMesh(mlnid, trim(outfile), trim(wopts))
if (ierr .ne. 0) then
call shr_sys_abort(subname//' ERROR in writing data mesh lnd ')
endif
#endif
#endif
!----------------------------------------------------------------------------
! Initialize MCT attribute vectors
!----------------------------------------------------------------------------
Expand Down Expand Up @@ -339,8 +474,15 @@ subroutine dlnd_comp_run(EClock, x2l, l2x, &
inst_suffix, logunit, case_name)

! !DESCRIPTION: run method for dlnd model
implicit none
#ifdef HAVE_MOAB
#ifdef MOABDEBUG
use iMOAB, only: iMOAB_WriteMesh
#endif
use seq_flds_mod , only: seq_flds_l2x_fields
use seq_flds_mod , only: moab_set_tag_from_av
#endif

implicit none
! !INPUT/OUTPUT PARAMETERS:
type(ESMF_Clock) , intent(in) :: EClock
type(mct_aVect) , intent(inout) :: x2l
Expand All @@ -366,6 +508,17 @@ subroutine dlnd_comp_run(EClock, x2l, l2x, &
integer(IN) :: nu ! unit number
logical :: write_restart ! restart now
character(len=18) :: date_str
#ifdef HAVE_MOAB
real(R8), allocatable, target :: datam(:)
type(mct_list) :: temp_list
integer :: size_list, index_list, lsize
type(mct_string) :: mctOStr !
character*400 tagname, mct_field
#ifdef MOABDEBUG
integer :: cur_dlnd_stepno, ierr
character*100 outfile, wopts, lnum
#endif
#endif

character(*), parameter :: F00 = "('(dlnd_comp_run) ',8a)"
character(*), parameter :: F04 = "('(dlnd_comp_run) ',2a,2i8,'s')"
Expand Down Expand Up @@ -464,6 +617,32 @@ subroutine dlnd_comp_run(EClock, x2l, l2x, &

call t_stopf('DLND_RUN')

#ifdef HAVE_MOAB
lsize = mct_avect_lsize(l2x) ! is it the same as mct_avect_lsize(avstrm) ?
allocate(datam(lsize)) !
call mct_list_init(temp_list ,seq_flds_l2x_fields)
size_list=mct_list_nitem (temp_list)
do index_list = 1, size_list
call mct_list_get(mctOStr,index_list,temp_list)
mct_field = mct_string_toChar(mctOStr)
tagname= trim(mct_field)//C_NULL_CHAR
call moab_set_tag_from_av(tagname, l2x, index_list, mlnid, datam, lsize) ! loop over all a2x fields, not just a few
enddo
call mct_list_clean(temp_list)
deallocate(datam) ! maybe we should keep it around, deallocate at the final only?

#ifdef MOABDEBUG
call seq_timemgr_EClockGetData( EClock, stepno=cur_dlnd_stepno )
write(lnum,"(I0.2)")cur_dlnd_stepno
outfile = 'dlnd_comp_run_'//trim(lnum)//'.h5m'//C_NULL_CHAR
wopts = 'PARALLEL=WRITE_PART'//C_NULL_CHAR
ierr = iMOAB_WriteMesh(mlnid, outfile, wopts)
if (ierr > 0 ) then
write(logunit,*) 'Failed to write data lnd component state '
endif
#endif
#endif

end subroutine dlnd_comp_run

!===============================================================================
Expand Down
26 changes: 22 additions & 4 deletions components/data_comps/dlnd/src/lnd_comp_mct.F90
Original file line number Diff line number Diff line change
Expand Up @@ -16,7 +16,11 @@ module lnd_comp_mct
use dlnd_comp_mod , only: dlnd_comp_init, dlnd_comp_run, dlnd_comp_final
use dlnd_shr_mod , only: dlnd_shr_read_namelists
use seq_flds_mod , only: seq_flds_x2l_fields, seq_flds_l2x_fields

#ifdef HAVE_MOAB
use seq_comm_mct, only : mlnid ! iMOAB app id for lnd
use iso_c_binding
use iMOAB , only: iMOAB_RegisterApplication
#endif
! !PUBLIC TYPES:
implicit none
private ! except
Expand Down Expand Up @@ -52,7 +56,9 @@ module lnd_comp_mct

!===============================================================================
subroutine lnd_init_mct( EClock, cdata, x2l, l2x, NLFilename )

#ifdef HAVE_MOAB
use shr_stream_mod, only: shr_stream_getDomainInfo, shr_stream_getFile
#endif
! !DESCRIPTION: initialize dlnd model
implicit none

Expand Down Expand Up @@ -146,13 +152,25 @@ subroutine lnd_init_mct( EClock, cdata, x2l, l2x, NLFilename )
!----------------------------------------------------------------------------
! Initialize dlnd
!----------------------------------------------------------------------------

#ifdef HAVE_MOAB
ierr = iMOAB_RegisterApplication(trim("DLND")//C_NULL_CHAR, mpicom, compid, mlnid)
if (ierr .ne. 0) then
write(logunit,*) subname,' error in registering data lnd comp'
call shr_sys_abort(subname//' ERROR in registering data lnd comp')
endif
#endif
call dlnd_comp_init(Eclock, x2l, l2x, &
seq_flds_x2l_fields, seq_flds_l2x_fields, &
SDLND, gsmap, ggrid, mpicom, compid, my_task, master_task, &
inst_suffix, inst_name, logunit, read_restart, &
scmMode, scmlat, scmlon)

#ifdef HAVE_MOAB
if (my_task == master_task) then
call seq_infodata_PutData( infodata, lnd_domain=SDLND%domainFile) ! we use the same one for regular case
iulian787 marked this conversation as resolved.
Show resolved Hide resolved
! in regular case, it is copied from fatmlndfrc ; so we don't know if it is data land or not
write(logunit,*), ' use this land domain file: ', SDLND%domainFile
endif
#endif
!----------------------------------------------------------------------------
! Fill infodata that needs to be returned from dlnd
!----------------------------------------------------------------------------
Expand Down
1 change: 1 addition & 0 deletions driver-moab/main/component_mod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -748,6 +748,7 @@ subroutine component_init_areacor_moab (comp, mbccid, mbcxid, seq_flds_c2x_fluxe
lsize = comp(1)%mblsize
allocate(areas (lsize, 3)) ! lsize is along grid; read mask too
allocate(factors (lsize, 2))
factors = 1.0 ! initialize with 1.0 all factors; then maybe correct them
iulian787 marked this conversation as resolved.
Show resolved Hide resolved
! get areas
tagname='area:aream:mask'//C_NULL_CHAR
arrsize = 3 * lsize
Expand Down
10 changes: 5 additions & 5 deletions driver-moab/main/prep_lnd_mod.F90
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
module prep_lnd_mod

use shr_kind_mod , only: r8 => SHR_KIND_R8
use shr_kind_mod , only: R8 => SHR_KIND_R8
use shr_kind_mod , only: cs => SHR_KIND_CS
use shr_kind_mod , only: cl => SHR_KIND_CL
use shr_kind_mod , only: cxx => SHR_KIND_CXX
Expand Down Expand Up @@ -162,7 +162,7 @@ subroutine prep_lnd_init(infodata, atm_c2_lnd, rof_c2_lnd, glc_c2_lnd, iac_c2_ln
integer nrflds ! number of rof fields projected on land
integer arrsize ! for setting the r2x fields on land to 0
integer ent_type ! for setting tags
real (kind=r8) , allocatable :: tmparray (:) ! used to set the r2x fields to 0
real (kind=R8) , allocatable :: tmparray (:) ! used to set the r2x fields to 0

#endif
character(*), parameter :: subname = '(prep_lnd_init)'
Expand Down Expand Up @@ -222,7 +222,6 @@ subroutine prep_lnd_init(infodata, atm_c2_lnd, rof_c2_lnd, glc_c2_lnd, iac_c2_ln
call seq_map_init_rcfile(mapper_Fr2l, rof(1), lnd(1), &
'seq_maps.rc','rof2lnd_fmapname:','rof2lnd_fmaptype:',samegrid_lr, &
string='mapper_Fr2l initialization',esmf_map=esmf_map_flag)
end if
! symmetric of l2r, from prep_rof
#ifdef HAVE_MOAB
! Call moab intx only if land and river are init in moab
Expand Down Expand Up @@ -370,7 +369,7 @@ subroutine prep_lnd_init(infodata, atm_c2_lnd, rof_c2_lnd, glc_c2_lnd, iac_c2_ln
arrsize = nrflds*mlsize
allocate (tmparray(arrsize)) ! mlsize is the size of local land
! do we need to zero out others or just river ?
tmparray = 0._r8
tmparray = 0._R8
ierr = iMOAB_SetDoubleTagStorage(mblxid, tagname, arrsize , ent_type, tmparray)
if (ierr .ne. 0) then
write(logunit,*) subname,' cant zero out r2x tags on land'
Expand All @@ -381,6 +380,7 @@ subroutine prep_lnd_init(infodata, atm_c2_lnd, rof_c2_lnd, glc_c2_lnd, iac_c2_ln
end if ! if ((mbrxid .ge. 0) .and. (mblxid .ge. 0))
! endif HAVE_MOAB
#endif
end if
call shr_sys_flush(logunit)

if (atm_c2_lnd) then
Expand Down Expand Up @@ -693,7 +693,7 @@ subroutine prep_lnd_mrg_moab (infodata)
#endif
#ifdef MOABCOMP
character(CXX) :: tagname, mct_field
real(r8) :: difference
real(R8) :: difference
vijaysm marked this conversation as resolved.
Show resolved Hide resolved
type(mct_list) :: temp_list
integer :: size_list, index_list, ent_type
type(mct_string) :: mctOStr !
Expand Down
Loading