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

Support JRA55-do #61

Merged
merged 12 commits into from
Aug 31, 2023
Merged
4 changes: 3 additions & 1 deletion CDEPS/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -16,14 +16,14 @@ list(APPEND cdeps_drof_src_files
)

list(APPEND cdeps_datm_src_files
CDEPS/datm/atm_comp_nuopc.F90
CDEPS/datm/datm_datamode_cfsr_mod.F90
CDEPS/datm/datm_datamode_clmncep_mod.F90
CDEPS/datm/datm_datamode_core2_mod.F90
CDEPS/datm/datm_datamode_cplhist_mod.F90
CDEPS/datm/datm_datamode_era5_mod.F90
CDEPS/datm/datm_datamode_gefs_mod.F90
CDEPS/datm/datm_datamode_jra_mod.F90
extra_sources/datm_datamode_jra55do_mod.F90
)


Expand Down Expand Up @@ -56,3 +56,5 @@ foreach(LIB cdeps_docn cdeps_dice cdeps_dwav cdeps_drof cdeps_datm)
target_include_directories(${LIB} PUBLIC $<BUILD_INTERFACE:${CMAKE_CURRENT_BINARY_DIR}/fox/include>)
target_link_libraries(${LIB} PUBLIC share cmeps esmf PIO::PIO_Fortran FoX_dom cdeps_common)
endforeach()

add_patched_source(cdeps_datm CDEPS/datm/atm_comp_nuopc.F90)
321 changes: 321 additions & 0 deletions CDEPS/extra_sources/datm_datamode_jra55do_mod.F90
Original file line number Diff line number Diff line change
@@ -0,0 +1,321 @@
module datm_datamode_jra55do_mod

use ESMF , only : ESMF_State, ESMF_StateGet, ESMF_SUCCESS, ESMF_LogWrite, ESMF_LOGMSG_INFO
use ESMF , only : ESMF_MeshGet
use ESMF , only : ESMF_StateItem_Flag, ESMF_STATEITEM_NOTFOUND, operator(/=)
use NUOPC , only : NUOPC_Advertise
use shr_kind_mod , only : r8=>shr_kind_r8, i8=>shr_kind_i8, cl=>shr_kind_cl, cs=>shr_kind_cs
use shr_sys_mod , only : shr_sys_abort
use shr_cal_mod , only : shr_cal_date2julian
use shr_const_mod , only : shr_const_tkfrz, shr_const_pi, shr_const_rdair
use dshr_strdata_mod , only : shr_strdata_get_stream_pointer, shr_strdata_type
use dshr_methods_mod , only : dshr_state_getfldptr, dshr_fldbun_getfldptr, dshr_fldbun_regrid, chkerr
use dshr_mod , only : dshr_restart_read, dshr_restart_write
use dshr_strdata_mod , only : shr_strdata_type
use dshr_fldlist_mod , only : fldlist_type, dshr_fldlist_add

implicit none
private ! except

public :: datm_datamode_jra55do_advertise
public :: datm_datamode_jra55do_init_pointers
public :: datm_datamode_jra55do_advance
public :: datm_datamode_jra55do_restart_write
public :: datm_datamode_jra55do_restart_read

! export state pointers
real(r8), pointer :: Sa_z(:) => null()
real(r8), pointer :: Sa_tbot(:) => null()
real(r8), pointer :: Sa_ptem(:) => null()
real(r8), pointer :: Sa_shum(:) => null()
real(r8), pointer :: Sa_dens(:) => null()
real(r8), pointer :: Sa_pbot(:) => null()
real(r8), pointer :: Sa_pslv(:) => null()
real(r8), pointer :: Faxa_rainc(:) => null()
real(r8), pointer :: Faxa_rainl(:) => null()
real(r8), pointer :: Faxa_snowc(:) => null()
real(r8), pointer :: Faxa_snowl(:) => null()
real(r8), pointer :: Faxa_swndr(:) => null()
real(r8), pointer :: Faxa_swndf(:) => null()
real(r8), pointer :: Faxa_swvdr(:) => null()
real(r8), pointer :: Faxa_swvdf(:) => null()
real(r8), pointer :: Faxa_swnet(:) => null()
real(r8), pointer :: Faxa_ndep(:,:) => null()

! stream data
real(r8), pointer :: strm_prrn(:) => null() ! Rainfall flux
real(r8), pointer :: strm_prsn(:) => null() ! Snowfall flux
real(r8), pointer :: strm_swdn(:) => null()

! othe module arrays
real(R8), pointer :: yc(:) ! array of model latitudes

! constants
real(R8) , parameter :: tKFrz = SHR_CONST_TKFRZ
real(R8) , parameter :: rdair = SHR_CONST_RDAIR ! dry air gas constant ~ J/K/kg
real(R8) , parameter :: degtorad = SHR_CONST_PI/180.0_R8
real(R8) , parameter :: phs_c0 = 0.298_R8
real(R8) , parameter :: dLWarc = -5.000_R8

character(*), parameter :: nullstr = 'null'
character(*), parameter :: rpfile = 'rpointer.atm'
character(*), parameter :: u_FILE_u = &
__FILE__

!===============================================================================
contains
!===============================================================================

subroutine datm_datamode_jra55do_advertise(exportState, fldsexport, flds_scalar_name, &
flds_co2, flds_wiso, flds_presaero, flds_presndep, rc)

! input/output variables
type(esmf_State) , intent(inout) :: exportState
type(fldlist_type) , pointer :: fldsexport
character(len=*) , intent(in) :: flds_scalar_name
logical , intent(in) :: flds_co2
logical , intent(in) :: flds_wiso
logical , intent(in) :: flds_presaero
logical , intent(in) :: flds_presndep
integer , intent(out) :: rc

! local variables
type(fldlist_type), pointer :: fldList
!-------------------------------------------------------------------------------

rc = ESMF_SUCCESS

call dshr_fldList_add(fldsExport, trim(flds_scalar_name))
call dshr_fldList_add(fldsExport, 'Sa_z' )
call dshr_fldList_add(fldsExport, 'Sa_u' )
call dshr_fldList_add(fldsExport, 'Sa_v' )
call dshr_fldList_add(fldsExport, 'Sa_ptem' )
call dshr_fldList_add(fldsExport, 'Sa_dens' )
call dshr_fldList_add(fldsExport, 'Sa_pslv' )
call dshr_fldList_add(fldsExport, 'Sa_tbot' )
call dshr_fldList_add(fldsExport, 'Sa_pbot' )
call dshr_fldList_add(fldsExport, 'Sa_shum' )
call dshr_fldList_add(fldsExport, 'Faxa_rainc' )
call dshr_fldList_add(fldsExport, 'Faxa_rainl' )
call dshr_fldList_add(fldsExport, 'Faxa_snowc' )
call dshr_fldList_add(fldsExport, 'Faxa_snowl' )
call dshr_fldList_add(fldsExport, 'Faxa_swndr' )
call dshr_fldList_add(fldsExport, 'Faxa_swvdr' )
call dshr_fldList_add(fldsExport, 'Faxa_swndf' )
call dshr_fldList_add(fldsExport, 'Faxa_swvdf' )
call dshr_fldList_add(fldsExport, 'Faxa_swnet' )
call dshr_fldList_add(fldsExport, 'Faxa_lwdn' )
call dshr_fldList_add(fldsExport, 'Faxa_swdn' )

if (flds_co2) then
call dshr_fldList_add(fldsExport, 'Sa_co2prog')
call dshr_fldList_add(fldsExport, 'Sa_co2diag')
end if
if (flds_presaero) then
call dshr_fldList_add(fldsExport, 'Faxa_bcph' , ungridded_lbound=1, ungridded_ubound=3)
call dshr_fldList_add(fldsExport, 'Faxa_ocph' , ungridded_lbound=1, ungridded_ubound=3)
call dshr_fldList_add(fldsExport, 'Faxa_dstwet' , ungridded_lbound=1, ungridded_ubound=4)
call dshr_fldList_add(fldsExport, 'Faxa_dstdry' , ungridded_lbound=1, ungridded_ubound=4)
end if
if (flds_presndep) then
call dshr_fldList_add(fldsExport, 'Faxa_ndep', ungridded_lbound=1, ungridded_ubound=2)
end if
if (flds_wiso) then
call dshr_fldList_add(fldsExport, 'Faxa_rainc_wiso', ungridded_lbound=1, ungridded_ubound=3)
call dshr_fldList_add(fldsExport, 'Faxa_rainl_wiso', ungridded_lbound=1, ungridded_ubound=3)
call dshr_fldList_add(fldsExport, 'Faxa_snowc_wiso', ungridded_lbound=1, ungridded_ubound=3)
call dshr_fldList_add(fldsExport, 'Faxa_snowl_wiso', ungridded_lbound=1, ungridded_ubound=3)
call dshr_fldList_add(fldsExport, 'Faxa_shum_wiso' , ungridded_lbound=1, ungridded_ubound=3)
end if

fldlist => fldsExport ! the head of the linked list
do while (associated(fldlist))
call NUOPC_Advertise(exportState, standardName=fldlist%stdname, rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call ESMF_LogWrite('(datm_comp_advertise): Fr_atm'//trim(fldList%stdname), ESMF_LOGMSG_INFO)
fldList => fldList%next
enddo

end subroutine datm_datamode_jra55do_advertise

!===============================================================================
subroutine datm_datamode_jra55do_init_pointers(exportState, sdat, rc)

! input/output variables
type(ESMF_State) , intent(inout) :: exportState
type(shr_strdata_type) , intent(in) :: sdat
integer , intent(out) :: rc

! local variables
integer :: n
integer :: lsize
integer :: spatialDim ! number of dimension in mesh
integer :: numOwnedElements ! size of mesh
real(r8), pointer :: ownedElemCoords(:) ! mesh lat and lons
type(ESMF_StateItem_Flag) :: itemFlag
character(len=*), parameter :: subname='(datm_datamode_jra55do_init_pointers): '
!-------------------------------------------------------------------------------

rc = ESMF_SUCCESS

lsize = sdat%model_lsize

call ESMF_MeshGet(sdat%model_mesh, spatialDim=spatialDim, numOwnedElements=numOwnedElements, rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
allocate(ownedElemCoords(spatialDim*numOwnedElements))
allocate(yc(numOwnedElements))
call ESMF_MeshGet(sdat%model_mesh, ownedElemCoords=ownedElemCoords)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
do n = 1,numOwnedElements
yc(n) = ownedElemCoords(2*n)
end do

call shr_strdata_get_stream_pointer( sdat, 'Faxa_prrn' , strm_prrn , rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call shr_strdata_get_stream_pointer( sdat, 'Faxa_prsn' , strm_prsn , rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call shr_strdata_get_stream_pointer( sdat, 'Faxa_swdn' , strm_swdn , rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return

call dshr_state_getfldptr(exportState, 'Sa_z' , fldptr1=Sa_z , rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call dshr_state_getfldptr(exportState, 'Sa_tbot' , fldptr1=Sa_tbot , rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call dshr_state_getfldptr(exportState, 'Sa_pbot' , fldptr1=Sa_pbot , rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call dshr_state_getfldptr(exportState, 'Sa_pslv' , fldptr1=Sa_pslv , rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call dshr_state_getfldptr(exportState, 'Sa_ptem' , fldptr1=Sa_ptem , rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call dshr_state_getfldptr(exportState, 'Sa_shum' , fldptr1=Sa_shum , rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call dshr_state_getfldptr(exportState, 'Sa_dens' , fldptr1=Sa_dens , rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call dshr_state_getfldptr(exportState, 'Faxa_rainc' , fldptr1=Faxa_rainc , rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call dshr_state_getfldptr(exportState, 'Faxa_rainl' , fldptr1=Faxa_rainl , rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call dshr_state_getfldptr(exportState, 'Faxa_snowc' , fldptr1=Faxa_snowc , rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call dshr_state_getfldptr(exportState, 'Faxa_snowl' , fldptr1=Faxa_snowl , rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call dshr_state_getfldptr(exportState, 'Faxa_swvdr' , fldptr1=Faxa_swvdr , rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call dshr_state_getfldptr(exportState, 'Faxa_swvdf' , fldptr1=Faxa_swvdf , rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call dshr_state_getfldptr(exportState, 'Faxa_swndr' , fldptr1=Faxa_swndr , rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call dshr_state_getfldptr(exportState, 'Faxa_swndf' , fldptr1=Faxa_swndf , rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call dshr_state_getfldptr(exportState, 'Faxa_swnet' , fldptr1=Faxa_swnet , rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return

call ESMF_StateGet(exportState, 'Faxa_ndep', itemFlag, rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
if (itemflag /= ESMF_STATEITEM_NOTFOUND) then
call dshr_state_getfldptr(exportState, 'Faxa_ndep', fldptr2=Faxa_ndep, rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
end if

! erro check
if (.not. associated(strm_prrn) .or. .not. associated(strm_prsn) .or. .not. associated(strm_swdn)) then
call shr_sys_abort(trim(subname)//'ERROR: prrn, prsn and swdn must be in streams for JRA55-do')
endif

end subroutine datm_datamode_jra55do_init_pointers

!===============================================================================
subroutine datm_datamode_jra55do_advance(exportstate, target_ymd, target_tod, model_calendar, rc)

! input/output variables
type(ESMF_State) , intent(inout) :: exportState
integer , intent(in) :: target_ymd
integer , intent(in) :: target_tod
character(len=*) , intent(in) :: model_calendar
integer , intent(out) :: rc

! local variables
integer :: n
integer :: lsize
real(R8) :: avg_alb ! average albedo
real(R8) :: rday ! elapsed day
real(R8) :: cosFactor ! cosine factor
character(len=*), parameter :: subname='(datm_datamode_jra55do_advance): '
!-------------------------------------------------------------------------------

rc = ESMF_SUCCESS

lsize = size(Sa_z)

call shr_cal_date2julian(target_ymd, target_tod, rday, model_calendar)
rday = mod((rday - 1.0_R8),365.0_R8)
cosfactor = cos((2.0_R8*SHR_CONST_PI*rday)/365 - phs_c0)

do n = 1,lsize
Sa_z(n) = 10.0_R8
Sa_pbot(n) = Sa_pslv(n)
Sa_ptem(n) = Sa_tbot(n)

! density computation for JRA55-do forcing
Sa_dens(n) = Sa_pbot(n)/(rdair*Sa_tbot(n)*(1 + 0.608*Sa_shum(n)))

! precipitation data
Faxa_rainc(n) = 0.0_R8 ! default zero
Faxa_snowc(n) = 0.0_R8
Faxa_snowl(n) = strm_prsn(n) ! Snowfall flux
Faxa_rainl(n) = strm_prrn(n) ! Rainfall flux

! radiation data - fabricate required swdn components from net swdn
Faxa_swvdr(n) = strm_swdn(n)*(0.28_R8)
Faxa_swndr(n) = strm_swdn(n)*(0.31_R8)
Faxa_swvdf(n) = strm_swdn(n)*(0.24_R8)
Faxa_swndf(n) = strm_swdn(n)*(0.17_R8)

! radiation data - compute net short-wave based on LY08 latitudinally-varying albedo
avg_alb = ( 0.069 - 0.011*cos(2.0_R8*yc(n)*degtorad ) )
Faxa_swnet(n) = strm_swdn(n)*(1.0_R8 - avg_alb)
enddo ! lsize

if (associated(Faxa_ndep)) then
! convert ndep flux to units of kgN/m2/s (input is in gN/m2/s)
Faxa_ndep(:,:) = Faxa_ndep(:,:) / 1000._r8
end if

end subroutine datm_datamode_jra55do_advance

!===============================================================================
subroutine datm_datamode_jra55do_restart_write(case_name, inst_suffix, ymd, tod, &
logunit, my_task, sdat)

! input/output variables
character(len=*) , intent(in) :: case_name
character(len=*) , intent(in) :: inst_suffix
integer , intent(in) :: ymd ! model date
integer , intent(in) :: tod ! model sec into model date
integer , intent(in) :: logunit
integer , intent(in) :: my_task
type(shr_strdata_type) , intent(inout) :: sdat
!-------------------------------------------------------------------------------

call dshr_restart_write(rpfile, case_name, 'datm', inst_suffix, ymd, tod, &
logunit, my_task, sdat)

end subroutine datm_datamode_jra55do_restart_write

!===============================================================================
subroutine datm_datamode_jra55do_restart_read(rest_filem, inst_suffix, logunit, my_task, mpicom, sdat)

! input/output arguments
character(len=*) , intent(inout) :: rest_filem
character(len=*) , intent(in) :: inst_suffix
integer , intent(in) :: logunit
integer , intent(in) :: my_task
integer , intent(in) :: mpicom
type(shr_strdata_type) , intent(inout) :: sdat
!-------------------------------------------------------------------------------

call dshr_restart_read(rest_filem, rpfile, inst_suffix, nullstr, logunit, my_task, mpicom, sdat)

end subroutine datm_datamode_jra55do_restart_read

end module datm_datamode_jra55do_mod
Loading