module CNDVDriverMod !----------------------------------------------------------------------- ! !DESCRIPTION: ! Note that this module was created simply to contain the subroutine dv ! which cannot cannot be in CNDVMod due to circular dependencies ! ! !USES: use shr_kind_mod , only : r8 => shr_kind_r8 use shr_log_mod , only : errMsg => shr_log_errMsg use abortutils , only : endrun use decompMod , only : bounds_type use atm2lndType , only : atm2lnd_type use CNDVType , only : dgvs_type use CNVegCarbonStateType , only : cnveg_carbonstate_type use CNVegCarbonFluxType , only : cnveg_carbonflux_type use clm_varcon , only : grlnd use LandunitType , only : lun use PatchType , only : patch ! ! !PUBLIC TYPES: implicit none private ! ! !PUBLIC MEMBER FUNCTIONS: public :: CNDVDriver public :: CNDVHist ! ! !PRIVATE MEMBER FUNCTIONS: private set_dgvm_filename character(len=*), parameter, private :: sourcefile = & __FILE__ !----------------------------------------------------------------------- contains !----------------------------------------------------------------------- subroutine CNDVDriver(bounds, & num_natvegp, filter_natvegp, kyr, & atm2lnd_inst, cnveg_carbonflux_inst, cnveg_carbonstate_inst, dgvs_inst) ! ! !DESCRIPTION: ! Drives the annual dynamic vegetation that works with CN ! ! !USES: use CNDVLightMod , only : Light use CNDVEstablishmentMod , only : Establishment ! ! !ARGUMENTS: type(bounds_type) , intent(in) :: bounds integer , intent(inout) :: num_natvegp ! number of naturally-vegetated patches in filter integer , intent(inout) :: filter_natvegp(:) ! filter for naturally-vegetated patches integer , intent(in) :: kyr ! used in routine climate20 below type(atm2lnd_type) , intent(inout) :: atm2lnd_inst type(cnveg_carbonflux_type) , intent(in) :: cnveg_carbonflux_inst type(cnveg_carbonstate_type) , intent(inout) :: cnveg_carbonstate_inst type(dgvs_type) , intent(inout) :: dgvs_inst ! ! !LOCAL VARIABLES: integer :: p ! patch index !----------------------------------------------------------------------- associate( & fpcgrid => dgvs_inst%fpcgrid_patch , & ! Input: [real(r8) (:) ] foliar projective cover on gridcell (fraction) agdd20 => dgvs_inst%agdd20_patch , & ! Output: [real(r8) (:) ] 20-yr running mean of agdd tmomin20 => dgvs_inst%tmomin20_patch , & ! Output: [real(r8) (:) ] 20-yr running mean of tmomin agdd => dgvs_inst%agdd_patch , & ! Input: [real(r8) (:) ] accumulated growing degree days above 5 t_mo_min => atm2lnd_inst%t_mo_min_patch , & ! Output: [real(r8) (:) ] annual min of t_mo (Kelvin) leafcmax => cnveg_carbonstate_inst%leafcmax_patch & ! Output: [real(r8) (:) ] (gC/m2) ann max leaf C ) ! ************************************************************************* ! S. Levis version of LPJ's routine climate20: 'Returns' tmomin20 & agdd20 ! for use in routine bioclim, which I have placed in routine Establishment ! Instead of 20-yr running mean of coldest monthly temperature, ! use 20-yr running mean of minimum 10-day running mean ! ************************************************************************* do p = bounds%begp, bounds%endp if (kyr == 2) then ! slevis: add ".and. start_type==arb_ic" here? tmomin20(p) = t_mo_min(p) ! NO, b/c want to be able to start dgvm agdd20(p) = agdd(p) ! w/ clmi file from non-dgvm simulation end if tmomin20(p) = (19._r8 * tmomin20(p) + t_mo_min(p)) / 20._r8 agdd20(p) = (19._r8 * agdd20(p) + agdd(p) ) / 20._r8 end do ! Rebuild filter of present natually-vegetated patches after Kill() num_natvegp = 0 do p = bounds%begp,bounds%endp if (dgvs_inst%present_patch(p)) then num_natvegp = num_natvegp + 1 filter_natvegp(num_natvegp) = p end if end do ! Returns fpcgrid and nind call Light(bounds, num_natvegp, filter_natvegp, & cnveg_carbonstate_inst, dgvs_inst) ! Returns updated fpcgrid, nind, crownarea, and present. Due to updated ! present, we do not use the natveg filter in this subroutine. call Establishment(bounds, & atm2lnd_inst, cnveg_carbonflux_inst, cnveg_carbonstate_inst, dgvs_inst) ! Reset dgvm variables needed in next yr (too few to keep subr. dvreset) do p = bounds%begp,bounds%endp leafcmax(p) = 0._r8 t_mo_min(p) = 1.0e+36_r8 end do end associate end subroutine CNDVDriver !----------------------------------------------------------------------- subroutine CNDVHist(bounds, dgvs_inst) ! ! !DESCRIPTION: ! Write CNDV history file ! ! !USES: use shr_const_mod , only : SHR_CONST_CDAY use shr_sys_mod , only : shr_sys_getenv use clm_varpar , only : maxpatch_pft use clm_varctl , only : caseid, ctitle, finidat, fsurdat, paramfile, iulog use clm_varcon , only : spval use clm_time_manager, only : get_ref_date, get_nstep, get_curr_date, get_curr_time use domainMod , only : ldomain use fileutils , only : get_filename use spmdMod , only : masterproc use ncdio_pio ! ! !ARGUMENTS: type(bounds_type), intent(in) :: bounds type(dgvs_type) , intent(in) :: dgvs_inst ! ! !LOCAL VARIABLES: character(len=256) :: dgvm_fn ! dgvm history filename type(file_desc_t) :: ncid ! netcdf file id integer :: ncprec ! output precision integer :: g,p,l ! indices integer :: ier ! error status integer :: mdcur, mscur, mcdate ! outputs from get_curr_time integer :: yr,mon,day,mcsec ! outputs from get_curr_date integer :: hours,minutes,secs ! hours,minutes,seconds of hh:mm:ss integer :: nstep ! time step integer :: nbsec ! seconds components of a date integer :: dimid ! dimension, variable id real(r8):: time ! current time character(len=256) :: str ! temporary string character(len= 8) :: curdate ! current date character(len= 8) :: curtime ! current time character(len= 10) :: basedate ! base date (yyyymmdd) character(len= 8) :: basesec ! base seconds real(r8) , pointer :: rbuf2dg (:,:) ! Input: [real(r8) (:,:)] temporary !----------------------------------------------------------------------- associate(& fpcgrid => dgvs_inst%fpcgrid_patch , & ! Input: [real(r8) (:)] foliar projective cover on gridcell (fraction) nind => dgvs_inst%nind_patch & ! Input: [real(r8) (:)] number of individuals (#/m**2) ) allocate(rbuf2dg(bounds%begg:bounds%endg,maxpatch_pft), stat=ier) if (ier /= 0) call endrun(msg='histCNDV: allocation error for rbuf2dg'//& errMsg(sourcefile, __LINE__)) ! Set output precision ncprec = ncd_double ! ----------------------------------------------------------------------- ! Create new netCDF file. File will be in define mode ! ----------------------------------------------------------------------- dgvm_fn = set_dgvm_filename() call ncd_pio_createfile(ncid, trim(dgvm_fn)) ! ----------------------------------------------------------------------- ! Create global attributes. ! ----------------------------------------------------------------------- str = 'CF1.0' call ncd_putatt (ncid, ncd_global, 'conventions', trim(str)) call getdatetime(curdate, curtime) str = 'created on ' // curdate // ' ' // curtime call ncd_putatt(ncid, ncd_global,'history', trim(str)) call shr_sys_getenv('LOGNAME', str, ier) if (ier /= 0) call endrun(msg='error: LOGNAME environment variable not defined'//& errMsg(sourcefile, __LINE__)) call ncd_putatt (ncid, ncd_global, 'logname', trim(str)) call shr_sys_getenv('HOST', str, ier) call ncd_putatt (ncid, ncd_global, 'host', trim(str)) str = 'Community Land Model: CLM3' call ncd_putatt (ncid, ncd_global, 'source', trim(str)) str = '$Name$' call ncd_putatt (ncid, ncd_global, 'version', trim(str)) str = '$Id$' call ncd_putatt (ncid, ncd_global, 'revision_id', trim(str)) str = ctitle call ncd_putatt (ncid, ncd_global, 'case_title', trim(str)) str = caseid call ncd_putatt (ncid, ncd_global, 'case_id', trim(str)) str = get_filename(fsurdat) call ncd_putatt(ncid, ncd_global, 'Surface_dataset', trim(str)) str = 'arbitrary initialization' if (finidat /= ' ') str = get_filename(finidat) call ncd_putatt(ncid, ncd_global, 'Initial_conditions_dataset', trim(str)) str = get_filename(paramfile) call ncd_putatt(ncid, ncd_global, 'PFT_physiological_constants_dataset', trim(str)) ! ----------------------------------------------------------------------- ! Define dimensions. ! ----------------------------------------------------------------------- if (ldomain%isgrid2d) then call ncd_defdim (ncid, 'lon' ,ldomain%ni, dimid) call ncd_defdim (ncid, 'lat' ,ldomain%nj, dimid) else call ncd_defdim (ncid, 'gridcell', ldomain%ns, dimid) end if call ncd_defdim (ncid, 'pft' , maxpatch_pft , dimid) call ncd_defdim (ncid, 'time', ncd_unlimited, dimid) call ncd_defdim (ncid, 'string_length', 80 , dimid) ! ----------------------------------------------------------------------- ! Define variables ! ----------------------------------------------------------------------- ! Define coordinate variables (including time) if (ldomain%isgrid2d) then call ncd_defvar(ncid=ncid, varname='lon', xtype=ncprec, dim1name='lon', & long_name='coordinate longitude', units='degrees_east') call ncd_defvar(ncid=ncid, varname='lat', xtype=ncprec, dim1name='lat', & long_name='coordinate latitude', units='degrees_north') end if call get_curr_time(mdcur, mscur) call get_ref_date(yr, mon, day, nbsec) hours = nbsec / 3600 minutes = (nbsec - hours*3600) / 60 secs = (nbsec - hours*3600 - minutes*60) write(basedate,80) yr,mon,day 80 format(i4.4,'-',i2.2,'-',i2.2) write(basesec ,90) hours, minutes, secs 90 format(i2.2,':',i2.2,':',i2.2) str = 'days since ' // basedate // " " // basesec time = mdcur + mscur/SHR_CONST_CDAY call ncd_defvar(ncid=ncid, varname='time', xtype=ncd_double, dim1name='time', & long_name='time', units=str) ! Define surface grid (coordinate variables, latitude, longitude, surface type). if (ldomain%isgrid2d) then call ncd_defvar(ncid=ncid, varname='longxy', xtype=ncprec, & dim1name='lon', dim2name='lat', & long_name='longitude', units='degrees_east') call ncd_defvar(ncid=ncid, varname='latixy', xtype=ncprec, & dim1name='lon', dim2name='lat', & long_name='latitude', units='degrees_north') call ncd_defvar(ncid=ncid, varname='landmask', xtype=ncd_int, & dim1name='lon', dim2name='lat', & long_name='land/ocean mask (0.=ocean and 1.=land)') else call ncd_defvar(ncid=ncid, varname='longxy', xtype=ncprec, & dim1name='gridcell',& long_name='longitude', units='degrees_east') call ncd_defvar(ncid=ncid, varname='latixy', xtype=ncprec, & dim1name='gridcell',& long_name='latitude', units='degrees_north') call ncd_defvar(ncid=ncid, varname='landmask', xtype=ncd_int, & dim1name='gridcell', & long_name='land/ocean mask (0.=ocean and 1.=land)') end if ! Define time information call ncd_defvar(ncid=ncid, varname='mcdate', xtype=ncd_int, dim1name='time',& long_name='current date (YYYYMMDD)') call ncd_defvar(ncid=ncid, varname='mcsec', xtype=ncd_int, dim1name='time',& long_name='current seconds of current date', units='s') call ncd_defvar(ncid=ncid, varname='mdcur', xtype=ncd_int, dim1name='time',& long_name='current day (from base day)') call ncd_defvar(ncid=ncid, varname='mscur', xtype=ncd_int, dim1name='time',& long_name='current seconds of current day', units='s') call ncd_defvar(ncid=ncid, varname='nstep', xtype=ncd_int, dim1name='time',& long_name='time step', units='s') ! Define time dependent variables if (ldomain%isgrid2d) then call ncd_defvar(ncid=ncid, varname='FPCGRID', xtype=ncprec, & dim1name='lon', dim2name='lat', dim3name='pft', dim4name='time', & long_name='plant functional type cover', units='fraction of vegetated area', & missing_value=spval, fill_value=spval) call ncd_defvar(ncid=ncid, varname='NIND', xtype=ncprec, & dim1name='lon', dim2name='lat', dim3name='pft', dim4name='time', & long_name='number of individuals', units='individuals/m2 vegetated land', & missing_value=spval, fill_value=spval) else call ncd_defvar(ncid=ncid, varname='FPCGRID', xtype=ncprec, & dim1name='gridcell', dim2name='pft', dim3name='time', & long_name='plant functional type cover', units='fraction of vegetated area', & missing_value=spval, fill_value=spval) call ncd_defvar(ncid=ncid, varname='NIND', xtype=ncprec, & dim1name='gridcell', dim2name='pft', dim3name='time', & long_name='number of individuals', units='individuals/m2 vegetated land', & missing_value=spval, fill_value=spval) end if call ncd_enddef(ncid) ! ----------------------------------------------------------------------- ! Write variables ! ----------------------------------------------------------------------- ! Write surface grid (coordinate variables, latitude, longitude, surface type). call ncd_io(ncid=ncid, varname='longxy' , data=ldomain%lonc, flag='write', & dim1name=grlnd) call ncd_io(ncid=ncid, varname='latixy' , data=ldomain%latc, flag='write', & dim1name=grlnd) call ncd_io(ncid=ncid, varname='landmask', data=ldomain%mask, flag='write', & dim1name=grlnd) ! Write current date, current seconds, current day, current nstep call get_curr_date(yr, mon, day, mcsec) mcdate = yr*10000 + mon*100 + day nstep = get_nstep() call ncd_io(ncid=ncid, varname='mcdate', data=mcdate, nt=1, flag='write') call ncd_io(ncid=ncid, varname='mcsec' , data=mcsec , nt=1, flag='write') call ncd_io(ncid=ncid, varname='mdcur' , data=mdcur , nt=1, flag='write') call ncd_io(ncid=ncid, varname='mscur' , data=mcsec , nt=1, flag='write') call ncd_io(ncid=ncid, varname='nstep' , data=nstep , nt=1, flag='write') call ncd_io(ncid=ncid, varname='time' , data=time , nt=1, flag='write') ! Write time dependent variables to CNDV history file ! The if .not. ifspecial statment below guarantees that the m index will ! always lie between 1 and maxpatch_pft rbuf2dg(bounds%begg : bounds%endg, :) = 0._r8 do p = bounds%begp,bounds%endp g = patch%gridcell(p) l = patch%landunit(p) if (.not. lun%ifspecial(l)) rbuf2dg(g,patch%mxy(p)) = fpcgrid(p)*100._r8 end do call ncd_io(ncid=ncid, varname='FPCGRID', dim1name=grlnd, data=rbuf2dg, & nt=1, flag='write') rbuf2dg(bounds%begg : bounds%endg, :) = 0._r8 do p = bounds%begp,bounds%endp g = patch%gridcell(p) l = patch%landunit(p) if (.not. lun%ifspecial(l)) rbuf2dg(g,patch%mxy(p)) = nind(p) end do call ncd_io(ncid=ncid, varname='NIND', dim1name=grlnd, data=rbuf2dg, & nt=1, flag='write') ! Deallocate dynamic memory deallocate(rbuf2dg) !------------------------------------------------------------------ ! Close and archive netcdf CNDV history file !------------------------------------------------------------------ call ncd_pio_closefile(ncid) if (masterproc) then write(iulog,*)'(histCNDV): Finished writing CNDV history dataset ',& trim(dgvm_fn), 'at nstep = ',get_nstep() end if end associate end subroutine CNDVHist !----------------------------------------------------------------------- character(len=256) function set_dgvm_filename () ! ! !DESCRIPTION: ! Determine initial dataset filenames ! ! !USES: use clm_varctl , only : caseid, inst_suffix use clm_time_manager , only : get_curr_date ! ! !ARGUMENTS: implicit none ! ! !LOCAL VARIABLES: character(len=256) :: cdate !date char string integer :: day !day (1 -> 31) integer :: mon !month (1 -> 12) integer :: yr !year (0 -> ...) integer :: sec !seconds into current day !----------------------------------------------------------------------- call get_curr_date (yr, mon, day, sec) write(cdate,'(i4.4,"-",i2.2,"-",i2.2,"-",i5.5)') yr,mon,day,sec set_dgvm_filename = "./"//trim(caseid)//".clm2"//trim(inst_suffix)//& ".hv."//trim(cdate)//".nc" end function set_dgvm_filename !----------------------------------------------------------------------- subroutine BuildNatVegFilter(bounds, num_natvegp, filter_natvegp, dgvs_inst) ! ! !DESCRIPTION: ! Reconstruct a filter of naturally-vegetated Patches for use in DGVM ! ! !ARGUMENTS: type(bounds_type) , intent(in) :: bounds integer , intent(out) :: num_natvegp ! number of patches in naturally-vegetated filter integer , intent(out) :: filter_natvegp(:) ! patch filter for naturally-vegetated points type(dgvs_type) , intent(in) :: dgvs_inst ! ! !LOCAL VARIABLES: integer :: p !----------------------------------------------------------------------- num_natvegp = 0 do p = bounds%begp,bounds%endp if (dgvs_inst%present_patch(p)) then num_natvegp = num_natvegp + 1 filter_natvegp(num_natvegp) = p end if end do end subroutine BuildNatVegFilter end module CNDVDriverMod