From 7150f84ea1fdf88779af8384734f3b38bf20253a Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Mon, 18 Sep 2023 16:32:31 -0600 Subject: [PATCH 01/52] Specify stream_fldFileName_swindow_(start,end) instead of stream_fldFileName_sdate. * Replicate stream_fldFileName_sdate behavior by using same file for stream_fldFileName_swindow_start and stream_fldFileName_swindow_end. * RXCROPMATURITY test looks fine. --- bld/CLMBuildNamelist.pm | 15 +- .../namelist_definition_ctsm.xml | 9 +- cime_config/SystemTests/rxcropmaturity.py | 3 +- src/biogeochem/CNPhenologyMod.F90 | 90 ++++++----- src/biogeochem/CropType.F90 | 15 +- src/cpl/share_esmf/cropcalStreamMod.F90 | 141 +++++++++++++----- src/main/clm_varctl.F90 | 2 +- 7 files changed, 186 insertions(+), 89 deletions(-) diff --git a/bld/CLMBuildNamelist.pm b/bld/CLMBuildNamelist.pm index 2136cfc05c..ce486d9e1c 100755 --- a/bld/CLMBuildNamelist.pm +++ b/bld/CLMBuildNamelist.pm @@ -3962,21 +3962,28 @@ sub setup_logic_cropcal_streams { # Option checks my $generate_crop_gdds = $nl->get_value('generate_crop_gdds') ; my $use_mxmat = $nl->get_value('use_mxmat') ; - my $sdate_file = $nl->get_value('stream_fldFileName_sdate') ; + my $swindow_start_file = $nl->get_value('stream_fldFileName_swindow_start') ; + my $swindow_end_file = $nl->get_value('stream_fldFileName_swindow_end') ; my $gdd_file = $nl->get_value('stream_fldFileName_cultivar_gdds') ; my $mesh_file = $nl->get_value('stream_meshfile_cropcal') ; + if ( ($swindow_start_file eq '' and $swindow_start_file ne '') or ($swindow_start_file ne '' and $swindow_start_file eq '') ) { + $log->fatal_error("When specifying sowing window dates, you must provide both swindow_start_file and swindow_end_file. To specify exact sowing dates, use the same file." ); + } if ( $generate_crop_gdds eq '.true.' ) { if ( $use_mxmat eq '.true.' ) { $log->fatal_error("If generate_crop_gdds is true, you must also set use_mxmat to false" ); } - if ( $sdate_file eq '' ) { - $log->fatal_error("If generate_crop_gdds is true, you must specify stream_fldFileName_sdate") + if ( $swindow_start_file eq '' or $swindow_end_file eq '' ) { + $log->fatal_error("If generate_crop_gdds is true, you must specify stream_fldFileName_swindow_start and stream_fldFileName_swindow_end") + } + if ( $swindow_start_file ne $swindow_end_file ) { + $log->fatal_error("If generate_crop_gdds is true, you must specify exact sowing dates by setting stream_fldFileName_swindow_start and stream_fldFileName_swindow_end to the same file") } if ( $gdd_file ne '' ) { $log->fatal_error("If generate_crop_gdds is true, do not specify stream_fldFileName_cultivar_gdds") } } - if ( $mesh_file eq '' and ( $sdate_file ne '' or $gdd_file ne '' ) ) { + if ( $mesh_file eq '' and ( $swindow_start_file ne '' or $gdd_file ne '' ) ) { $log->fatal_error("If prescribing crop sowing dates and/or maturity requirements, you must specify stream_meshfile_cropcal") } } diff --git a/bld/namelist_files/namelist_definition_ctsm.xml b/bld/namelist_files/namelist_definition_ctsm.xml index 9a0663c34f..caeab9eb48 100644 --- a/bld/namelist_files/namelist_definition_ctsm.xml +++ b/bld/namelist_files/namelist_definition_ctsm.xml @@ -1835,9 +1835,14 @@ Last year to loop over for crop calendar data Simulation year that aligns with stream_year_first_cropcal value - -Filename of input stream data for sowing dates +Filename of input stream data for date (day of year) of start of sowing window + + + +Filename of input stream data for date (day of year) of end of sowing window crop_inst%fertnitro_patch , & ! Input: [real(r8) (:) ] fertilizer nitrogen hui => crop_inst%hui_patch , & ! Input: [real(r8) (:) ] crop patch heat unit index (growing degree-days); set to 0 at sowing and accumulated until harvest - leafout => crop_inst%gddtsoi_patch , & ! Input: [real(r8) (:) ] gdd from top soil layer temperature + leafout => crop_inst%gddtsoi_patch , & ! Input: [real(r8) (:) ] gdd from top soil layer temperature harvdate => crop_inst%harvdate_patch , & ! Output: [integer (:) ] harvest date - croplive => crop_inst%croplive_patch , & ! Output: [logical (:) ] Flag, true if planted, not harvested + croplive => crop_inst%croplive_patch , & ! Output: [logical (:) ] Flag, true if planted, not harvested vf => crop_inst%vf_patch , & ! Output: [real(r8) (:) ] vernalization factor - next_rx_sdate => crop_inst%next_rx_sdate_patch , & ! Inout: [integer (:) ] prescribed sowing date of next growing season this year + next_rx_swindow_start => crop_inst%next_rx_swindow_start_patch , & ! Inout: [integer (:) ] start of prescribed sowing window of next growing season this year + next_rx_swindow_end => crop_inst%next_rx_swindow_end_patch , & ! Inout: [integer (:) ] end of prescribed sowing window of next growing season this year sowing_count => crop_inst%sowing_count , & ! Inout: [integer (:) ] number of sowing events this year for this patch harvest_count => crop_inst%harvest_count , & ! Inout: [integer (:) ] number of harvest events this year for this patch peaklai => cnveg_state_inst%peaklai_patch , & ! Output: [integer (:) ] 1: max allowed lai; 0: not at max @@ -1907,17 +1910,21 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & do k = repr_grain_min, repr_grain_max cnveg_carbonflux_inst%repr_grainc_to_food_thisyr_patch(p,k) = 0._r8 end do - next_rx_sdate(p) = crop_inst%rx_sdates_thisyr_patch(p,1) + next_rx_swindow_start(p) = crop_inst%rx_swindow_starts_thisyr_patch(p,1) + next_rx_swindow_end (p) = crop_inst%rx_swindow_ends_thisyr_patch (p,1) end if - ! Get next sowing date + ! Get dates of next prescribed sowing window (if any) if (sowing_count(p) < mxsowings) then - next_rx_sdate(p) = crop_inst%rx_sdates_thisyr_patch(p,sowing_count(p)+1) + next_rx_swindow_start(p) = crop_inst%rx_swindow_starts_thisyr_patch(p,sowing_count(p)+1) + next_rx_swindow_end (p) = crop_inst%rx_swindow_ends_thisyr_patch (p,sowing_count(p)+1) end if - ! CropType->InitAllocate() initializes next_rx_sdate to -1. It's only changed from that, by cropCalStreamMod->cropcal_interp(), when use_cropcal_rx_sdates is true. So if not using prescribed sowing dates, this boolean will always be false, because jday can never be -1. - do_plant_prescribed = next_rx_sdate(p) == jday .and. sowing_count(p) < mxsowings - + ! CropType->InitAllocate() initializes next_rx_swindow_start to -1. It's only changed from that, by cropCalStreamMod->cropcal_interp(), when use_cropcal_rx_swindows is true. So if not using prescribed sowing windows, this boolean will always be false, because jday can never be -1. + do_plant_prescribed = next_rx_swindow_start(p) == jday .and. sowing_count(p) < mxsowings + ! We only want to plant on a specific day if the prescribed sowing window starts AND ends on the same day. + do_plant_prescribed = do_plant_prescribed .and. next_rx_swindow_start(p) == next_rx_swindow_end(p) + ! BACKWARDS_COMPATIBILITY(wjs/ssr, 2022-02-18) ! When resuming from a run with old code, may need to manually set these. ! Will be needed until we can rely on all restart files have been generated @@ -1931,13 +1938,26 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & end if ! This is outside the croplive check so that the "harvest if planting conditions were met today" conditional works. - is_in_sowing_window = is_doy_in_interval(minplantjday(ivt(p),h), maxplantjday(ivt(p),h), jday) - is_end_sowing_window = jday == maxplantjday(ivt(p),h) + if (use_cropcal_rx_swindows) then + if (sowing_count(p) < mxsowings) then + sowing_window_startdate = crop_inst%rx_swindow_starts_thisyr_patch(p,sowing_count(p)+1) + sowing_window_enddate = crop_inst%rx_swindow_ends_thisyr_patch (p,sowing_count(p)+1) + else ! Do not allow more sowings after reaching max number of sowings (mxsowings) + sowing_window_startdate = 999 + sowing_window_enddate = 999 + end if + else + sowing_window_startdate = minplantjday(ivt(p),h) + sowing_window_enddate = maxplantjday(ivt(p),h) + end if + is_in_sowing_window = is_doy_in_interval(sowing_window_startdate, sowing_window_enddate, jday) + is_end_sowing_window = jday == sowing_window_enddate ! ! Only allow sowing according to normal "window" rules if not using prescribed - ! sowing dates at all, or if this cell had no values in the prescribed sowing - ! date file. - allow_unprescribed_planting = (.not. use_cropcal_rx_sdates) .or. crop_inst%rx_sdates_thisyr_patch(p,1)<0 + ! sowing windows at all, or if this cell had no values in either of the + ! prescribed sowing window files. (Note that if either sowing window start was + ! or end was missing a value, they are both set to negative values.) + allow_unprescribed_planting = (.not. use_cropcal_rx_swindows) .or. crop_inst%rx_swindow_starts_thisyr_patch(p,1)<0 if (sowing_count(p) == mxsowings) then do_plant_normal = .false. do_plant_lastchance = .false. @@ -2140,7 +2160,7 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & do_harvest = .false. fake_harvest = .false. did_plant_prescribed_today = .false. - if (use_cropcal_rx_sdates .and. sowing_count(p) > 0) then + if (use_cropcal_rx_swindows .and. sowing_count(p) > 0) then did_plant_prescribed_today = crop_inst%sdates_thisyr_patch(p,sowing_count(p)) == real(jday, r8) end if @@ -2161,10 +2181,10 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & harvest_reason = HARVEST_REASON_SOWTODAY ! If generate_crop_gdds and this patch has prescribed sowing inputs - else if (generate_crop_gdds .and. crop_inst%rx_sdates_thisyr_patch(p,1) .gt. 0) then - if (next_rx_sdate(p) >= 0) then - ! Harvest the day before the next sowing date this year. - do_harvest = jday == next_rx_sdate(p) - 1 + else if (generate_crop_gdds .and. crop_inst%rx_swindow_starts_thisyr_patch(p,1) .gt. 0) then + if (next_rx_swindow_start(p) >= 0) then + ! Harvest the day before the start of the next sowing window this year. + do_harvest = jday == next_rx_swindow_start(p) - 1 ! ... unless that will lead to growing season length 365 (or 366, ! if last year was a leap year). This would result in idop==jday, @@ -2188,10 +2208,10 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & do_harvest = .false. ! ... unless first sowing next year happens Jan. 1. - ! WARNING: This implementation assumes that sowing dates don't change over time! + ! WARNING: This implementation assumes that sowing windows don't change over time! ! In order to avoid this, you'd have to read this year's AND next year's prescribed - ! sowing dates. - if (crop_inst%rx_sdates_thisyr_patch(p,1) == 1) then + ! sowing window start dates. + if (crop_inst%rx_swindow_starts_thisyr_patch(p,1) == 1) then do_harvest = jday == dayspyr end if @@ -2212,13 +2232,13 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & ! Original harvest rule do_harvest = hui(p) >= gddmaturity(p) .or. idpp >= mxmat - ! Always harvest the day before the next prescribed sowing, if still alive. - ! WARNING: This implementation assumes that sowing dates don't change over time! + ! Always harvest the day before the next prescribed sowing window starts, if still alive. + ! WARNING: This implementation assumes that sowing windows don't change over time! ! In order to avoid this, you'd have to read this year's AND next year's prescribed - ! sowing dates. - ! WARNING: This implementation assumes that all patches use prescribed sowing dates. - if (use_cropcal_rx_sdates) then - will_plant_prescribed_tomorrow = (jday == next_rx_sdate(p) - 1) .or. & + ! sowing windows. + ! WARNING: This implementation assumes that all patches use prescribed sowing windows. + if (use_cropcal_rx_swindows) then + will_plant_prescribed_tomorrow = (jday == next_rx_swindow_start(p) - 1) .or. & (crop_inst%sdates_thisyr_patch(p,1) == 1 .and. & jday == dayspyr) else @@ -2497,7 +2517,7 @@ subroutine PlantCrop(p, leafcn_in, jday, kyr, do_plant_normal, & ! !USES: use clm_varctl , only : use_c13, use_c14 - use clm_varctl , only : use_cropcal_rx_sdates, use_cropcal_rx_cultivar_gdds, use_cropcal_streams + use clm_varctl , only : use_cropcal_rx_swindows, use_cropcal_rx_cultivar_gdds, use_cropcal_streams use clm_varcon , only : c13ratio, c14ratio use clm_varpar , only : mxsowings use pftconMod , only : ntmp_corn, nswheat, nwwheat, ntmp_soybean @@ -2537,7 +2557,7 @@ subroutine PlantCrop(p, leafcn_in, jday, kyr, do_plant_normal, & ivt => patch%itype , & ! Input: [integer (:) ] patch vegetation type croplive => crop_inst%croplive_patch , & ! Output: [logical (:) ] Flag, true if planted, not harvested harvdate => crop_inst%harvdate_patch , & ! Output: [integer (:) ] harvest date - next_rx_sdate => crop_inst%next_rx_sdate_patch , & ! Inout: [integer (:) ] prescribed sowing date of next growing season this year + next_rx_swindow_start => crop_inst%next_rx_swindow_start_patch , & ! Inout: [integer (:) ] start of prescribed sowing window of next growing season this year sowing_count => crop_inst%sowing_count , & ! Inout: [integer (:) ] number of sowing events this year for this patch sowing_reason => crop_inst%sowing_reason_thisyr_patch , & ! Output: [real(r8) (:) ] reason for each sowing this year for this patch gddmaturity => cnveg_state_inst%gddmaturity_patch , & ! Output: [real(r8) (:) ] gdd needed to harvest @@ -2568,10 +2588,10 @@ subroutine PlantCrop(p, leafcn_in, jday, kyr, do_plant_normal, & harvdate(p) = NOT_Harvested sowing_count(p) = sowing_count(p) + 1 - if (use_cropcal_rx_sdates .and. sowing_count(p) < mxsowings) then - next_rx_sdate(p) = crop_inst%rx_sdates_thisyr_patch(p, sowing_count(p)+1) + if (use_cropcal_rx_swindows .and. sowing_count(p) < mxsowings) then + next_rx_swindow_start(p) = crop_inst%rx_swindow_starts_thisyr_patch(p, sowing_count(p)+1) else - next_rx_sdate(p) = -1 + next_rx_swindow_start(p) = -1 endif crop_inst%sdates_thisyr_patch(p,sowing_count(p)) = real(jday, r8) diff --git a/src/biogeochem/CropType.F90 b/src/biogeochem/CropType.F90 index f84c4c1821..564f4babc8 100644 --- a/src/biogeochem/CropType.F90 +++ b/src/biogeochem/CropType.F90 @@ -47,8 +47,10 @@ module CropType character(len=20) :: baset_mapping real(r8) :: baset_latvary_intercept real(r8) :: baset_latvary_slope - integer , pointer :: next_rx_sdate_patch (:) ! prescribed sowing date for the next growing season this year - integer , pointer :: rx_sdates_thisyr_patch (:,:) ! all prescribed sowing dates for this patch this year (day of year) [patch, mxsowings] + integer , pointer :: next_rx_swindow_start_patch (:) ! start of prescribed sowing window for the next growing season this year + integer , pointer :: next_rx_swindow_end_patch (:) ! end of prescribed sowing window for the next growing season this year + integer , pointer :: rx_swindow_starts_thisyr_patch(:,:) ! all prescribed sowing window start dates for this patch this year (day of year) [patch, mxsowings] + integer , pointer :: rx_swindow_ends_thisyr_patch (:,:) ! all prescribed sowing window end dates for this patch this year (day of year) [patch, mxsowings] real(r8), pointer :: rx_cultivar_gdds_thisyr_patch (:,:) ! all cultivar GDD targets for this patch this year (ddays) [patch, mxsowings] real(r8), pointer :: sdates_thisyr_patch (:,:) ! all actual sowing dates for this patch this year (day of year) [patch, mxsowings] real(r8), pointer :: sdates_perharv_patch (:,:) ! all actual sowing dates for crops *harvested* this year (day of year) [patch, mxharvests] @@ -227,8 +229,10 @@ subroutine InitAllocate(this, bounds) allocate(this%cphase_patch (begp:endp)) ; this%cphase_patch (:) = cphase_not_planted allocate(this%sowing_reason_patch (begp:endp)) ; this%sowing_reason_patch (:) = -1 allocate(this%latbaset_patch (begp:endp)) ; this%latbaset_patch (:) = spval - allocate(this%next_rx_sdate_patch(begp:endp)) ; this%next_rx_sdate_patch(:) = -1 - allocate(this%rx_sdates_thisyr_patch(begp:endp,1:mxsowings)) ; this%rx_sdates_thisyr_patch(:,:) = -1 + allocate(this%next_rx_swindow_start_patch(begp:endp)) ; this%next_rx_swindow_start_patch(:) = -1 + allocate(this%next_rx_swindow_end_patch (begp:endp)) ; this%next_rx_swindow_end_patch (:) = -1 + allocate(this%rx_swindow_starts_thisyr_patch(begp:endp,1:mxsowings)); this%rx_swindow_starts_thisyr_patch(:,:) = -1 + allocate(this%rx_swindow_ends_thisyr_patch(begp:endp,1:mxsowings)) ; this%rx_swindow_ends_thisyr_patch (:,:) = -1 allocate(this%rx_cultivar_gdds_thisyr_patch(begp:endp,1:mxsowings)) ; this%rx_cultivar_gdds_thisyr_patch(:,:) = spval allocate(this%sdates_thisyr_patch(begp:endp,1:mxsowings)) ; this%sdates_thisyr_patch(:,:) = spval allocate(this%sdates_perharv_patch(begp:endp,1:mxharvests)) ; this%sdates_perharv_patch(:,:) = spval @@ -242,7 +246,8 @@ subroutine InitAllocate(this, bounds) allocate(this%sowing_count(begp:endp)) ; this%sowing_count(:) = 0 allocate(this%harvest_count(begp:endp)) ; this%harvest_count(:) = 0 - this%rx_sdates_thisyr_patch(:,:) = -1 + this%rx_swindow_starts_thisyr_patch(:,:) = -1 + this%rx_swindow_ends_thisyr_patch (:,:) = -1 this%rx_cultivar_gdds_thisyr_patch(:,:) = spval end subroutine InitAllocate diff --git a/src/cpl/share_esmf/cropcalStreamMod.F90 b/src/cpl/share_esmf/cropcalStreamMod.F90 index 8c2168d01d..91f565771d 100644 --- a/src/cpl/share_esmf/cropcalStreamMod.F90 +++ b/src/cpl/share_esmf/cropcalStreamMod.F90 @@ -13,7 +13,7 @@ module cropcalStreamMod use decompMod , only : bounds_type use abortutils , only : endrun use clm_varctl , only : iulog - use clm_varctl , only : use_cropcal_rx_sdates, use_cropcal_rx_cultivar_gdds, use_cropcal_streams + use clm_varctl , only : use_cropcal_rx_swindows, use_cropcal_rx_cultivar_gdds, use_cropcal_streams use clm_varpar , only : mxpft use perf_mod , only : t_startf, t_stopf use spmdMod , only : masterproc, mpicom, iam @@ -31,12 +31,14 @@ module cropcalStreamMod ! !PRIVATE MEMBER DATA: integer, allocatable :: g_to_ig(:) ! Array matching gridcell index to data index - type(shr_strdata_type) :: sdat_cropcal_sdate ! sdate input data stream + type(shr_strdata_type) :: sdat_cropcal_swindow_start ! sowing window start input data stream + type(shr_strdata_type) :: sdat_cropcal_swindow_end ! sowing window end input data stream type(shr_strdata_type) :: sdat_cropcal_cultivar_gdds ! sdate input data stream - character(len=CS), allocatable :: stream_varnames_sdate(:) + character(len=CS), allocatable :: stream_varnames_sdate(:) ! used for both start and end dates character(len=CS), allocatable :: stream_varnames_cultivar_gdds(:) integer :: ncft ! Number of crop functional types (excl. generic crops) - character(len=CL) :: stream_fldFileName_sdate ! sdate stream filename to read + character(len=CL) :: stream_fldFileName_swindow_start ! sowing window start stream filename to read + character(len=CL) :: stream_fldFileName_swindow_end ! sowing window end stream filename to read character(len=CL) :: stream_fldFileName_cultivar_gdds ! cultivar growing degree-days stream filename to read character(len=*), parameter :: sourcefile = & @@ -81,7 +83,8 @@ subroutine cropcal_init(bounds) stream_year_first_cropcal, & stream_year_last_cropcal, & model_year_align_cropcal, & - stream_fldFileName_sdate, & + stream_fldFileName_swindow_start, & + stream_fldFileName_swindow_end, & stream_fldFileName_cultivar_gdds, & stream_meshfile_cropcal @@ -90,7 +93,8 @@ subroutine cropcal_init(bounds) stream_year_last_cropcal = 1 ! last year in stream to use model_year_align_cropcal = 1 ! align stream_year_first_cropcal with this model year stream_meshfile_cropcal = '' - stream_fldFileName_sdate = '' + stream_fldFileName_swindow_start = '' + stream_fldFileName_swindow_end = '' stream_fldFileName_cultivar_gdds = '' ! Will need modification to work with mxsowings > 1 ncft = mxpft - npcropmin + 1 ! Ignores generic crops @@ -119,7 +123,8 @@ subroutine cropcal_init(bounds) call shr_mpi_bcast(stream_year_first_cropcal , mpicom) call shr_mpi_bcast(stream_year_last_cropcal , mpicom) call shr_mpi_bcast(model_year_align_cropcal , mpicom) - call shr_mpi_bcast(stream_fldFileName_sdate , mpicom) + call shr_mpi_bcast(stream_fldFileName_swindow_start, mpicom) + call shr_mpi_bcast(stream_fldFileName_swindow_end , mpicom) call shr_mpi_bcast(stream_fldFileName_cultivar_gdds, mpicom) call shr_mpi_bcast(stream_meshfile_cropcal , mpicom) @@ -129,7 +134,8 @@ subroutine cropcal_init(bounds) write(iulog,'(a,i8)') ' stream_year_first_cropcal = ',stream_year_first_cropcal write(iulog,'(a,i8)') ' stream_year_last_cropcal = ',stream_year_last_cropcal write(iulog,'(a,i8)') ' model_year_align_cropcal = ',model_year_align_cropcal - write(iulog,'(a,a)' ) ' stream_fldFileName_sdate = ',trim(stream_fldFileName_sdate) + write(iulog,'(a,a)' ) ' stream_fldFileName_swindow_start = ',trim(stream_fldFileName_swindow_start) + write(iulog,'(a,a)' ) ' stream_fldFileName_swindow_end = ',trim(stream_fldFileName_swindow_end) write(iulog,'(a,a)' ) ' stream_fldFileName_cultivar_gdds = ',trim(stream_fldFileName_cultivar_gdds) write(iulog,'(a,a)' ) ' stream_meshfile_cropcal = ',trim(stream_meshfile_cropcal) do n = 1,ncft @@ -139,14 +145,15 @@ subroutine cropcal_init(bounds) write(iulog,*) endif - use_cropcal_rx_sdates = stream_fldFileName_sdate /= '' + ! CLMBuildNamelist checks that both start and end files are provided if either is + use_cropcal_rx_swindows = stream_fldFileName_swindow_start /= '' use_cropcal_rx_cultivar_gdds = stream_fldFileName_cultivar_gdds /= '' - use_cropcal_streams = use_cropcal_rx_sdates .or. use_cropcal_rx_cultivar_gdds + use_cropcal_streams = use_cropcal_rx_swindows .or. use_cropcal_rx_cultivar_gdds - ! Initialize the cdeps data type sdat_cropcal_sdate + ! Initialize the cdeps data type sdat_cropcal_swindow_start ! NOTE: stream_dtlimit 1.5 didn't work for some reason - if (use_cropcal_rx_sdates) then - call shr_strdata_init_from_inline(sdat_cropcal_sdate, & + if (use_cropcal_rx_swindows) then + call shr_strdata_init_from_inline(sdat_cropcal_swindow_start, & my_task = iam, & logunit = iulog, & compname = 'LND', & @@ -155,7 +162,7 @@ subroutine cropcal_init(bounds) stream_meshfile = trim(stream_meshfile_cropcal), & stream_lev_dimname = 'null', & stream_mapalgo = trim(cropcal_mapalgo), & - stream_filenames = (/trim(stream_fldFileName_sdate)/), & + stream_filenames = (/trim(stream_fldFileName_swindow_start)/), & stream_fldlistFile = stream_varnames_sdate, & stream_fldListModel = stream_varnames_sdate, & stream_yearFirst = stream_year_first_cropcal, & @@ -165,7 +172,36 @@ subroutine cropcal_init(bounds) stream_taxmode = 'extend', & stream_dtlimit = 1.0e30_r8, & stream_tintalgo = cropcal_tintalgo, & - stream_name = 'sowing date data', & + stream_name = 'sowing window start data', & + rc = rc) + if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) then + call ESMF_Finalize(endflag=ESMF_END_ABORT) + end if + end if + + ! Initialize the cdeps data type sdat_cropcal_swindow_end + ! NOTE: stream_dtlimit 1.5 didn't work for some reason + if (use_cropcal_rx_swindows) then + call shr_strdata_init_from_inline(sdat_cropcal_swindow_end, & + my_task = iam, & + logunit = iulog, & + compname = 'LND', & + model_clock = model_clock, & + model_mesh = mesh, & + stream_meshfile = trim(stream_meshfile_cropcal), & + stream_lev_dimname = 'null', & + stream_mapalgo = trim(cropcal_mapalgo), & + stream_filenames = (/trim(stream_fldFileName_swindow_end)/), & + stream_fldlistFile = stream_varnames_sdate, & + stream_fldListModel = stream_varnames_sdate, & + stream_yearFirst = stream_year_first_cropcal, & + stream_yearLast = stream_year_last_cropcal, & + stream_yearAlign = model_year_align_cropcal, & + stream_offset = cropcal_offset, & + stream_taxmode = 'extend', & + stream_dtlimit = 1.0e30_r8, & + stream_tintalgo = cropcal_tintalgo, & + stream_name = 'sowing window end data', & rc = rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) then call ESMF_Finalize(endflag=ESMF_END_ABORT) @@ -227,8 +263,12 @@ subroutine cropcal_advance( bounds ) call get_curr_date(year, mon, day, sec) mcdate = year*10000 + mon*100 + day - if (use_cropcal_rx_sdates) then - call shr_strdata_advance(sdat_cropcal_sdate, ymd=mcdate, tod=sec, logunit=iulog, istr='cropcaldyn', rc=rc) + if (use_cropcal_rx_swindows) then + call shr_strdata_advance(sdat_cropcal_swindow_start, ymd=mcdate, tod=sec, logunit=iulog, istr='cropcaldyn', rc=rc) + if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) then + call ESMF_Finalize(endflag=ESMF_END_ABORT) + end if + call shr_strdata_advance(sdat_cropcal_swindow_end, ymd=mcdate, tod=sec, logunit=iulog, istr='cropcaldyn', rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) then call ESMF_Finalize(endflag=ESMF_END_ABORT) end if @@ -275,9 +315,11 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) integer :: n, g integer :: lsize integer :: rc - real(r8), pointer :: dataptr1d_sdate(:) - real(r8), pointer :: dataptr2d_sdate(:,:) + real(r8), pointer :: dataptr1d_swindow_start(:) + real(r8), pointer :: dataptr1d_swindow_end (:) real(r8), pointer :: dataptr1d_cultivar_gdds(:) + real(r8), pointer :: dataptr2d_swindow_start(:,:) + real(r8), pointer :: dataptr2d_swindow_end (:,:) real(r8), pointer :: dataptr2d_cultivar_gdds(:,:) !----------------------------------------------------------------------- @@ -288,30 +330,39 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) ! Place all data from each type into a temporary 2d array lsize = bounds%endg - bounds%begg + 1 - ! Read prescribed sowing dates from input files - allocate(dataptr2d_sdate(lsize, ncft)) - if (use_cropcal_rx_sdates) then + ! Read prescribed sowing window start dates from input files + allocate(dataptr2d_swindow_start(lsize, ncft)) + allocate(dataptr2d_swindow_end (lsize, ncft)) + if (use_cropcal_rx_swindows) then ! Starting with npcropmin will skip generic crops do n = 1, ncft - call dshr_fldbun_getFldPtr(sdat_cropcal_sdate%pstrm(1)%fldbun_model, trim(stream_varnames_sdate(n)), & - fldptr1=dataptr1d_sdate, rc=rc) + call dshr_fldbun_getFldPtr(sdat_cropcal_swindow_start%pstrm(1)%fldbun_model, trim(stream_varnames_sdate(n)), & + fldptr1=dataptr1d_swindow_start, rc=rc) + if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) then + call ESMF_Finalize(endflag=ESMF_END_ABORT) + end if + call dshr_fldbun_getFldPtr(sdat_cropcal_swindow_end%pstrm(1)%fldbun_model, trim(stream_varnames_sdate(n)), & + fldptr1=dataptr1d_swindow_end, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) then call ESMF_Finalize(endflag=ESMF_END_ABORT) end if ! Note that the size of dataptr1d includes ocean points so it will be around 3x larger than lsize ! So an explicit loop is required here do g = 1,lsize - + ! If read-in value is invalid, allow_unprescribed_planting in CropPhenology() - if (dataptr1d_sdate(g) <= 0 .or. dataptr1d_sdate(g) > 365) then - dataptr1d_sdate(g) = -1 + if (dataptr1d_swindow_start(g) <= 0 .or. dataptr1d_swindow_start(g) > 365 & + .or. dataptr1d_swindow_end(g) <= 0 .or. dataptr1d_swindow_end(g) > 365) then + dataptr1d_swindow_start(g) = -1 + dataptr1d_swindow_end (g) = -1 end if - - dataptr2d_sdate(g,n) = dataptr1d_sdate(g) + + dataptr2d_swindow_start(g,n) = dataptr1d_swindow_start(g) + dataptr2d_swindow_end (g,n) = dataptr1d_swindow_end (g) end do end do - - ! Set rx_sdate for each gridcell/patch combination + + ! Set sowing window for each gridcell/patch combination do fp = 1, num_pcropp p = filter_pcropp(fp) ivt = patch%itype(p) @@ -320,27 +371,35 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) n = ivt - npcropmin + 1 ! vegetated pft ig = g_to_ig(patch%gridcell(p)) - crop_inst%rx_sdates_thisyr_patch(p,1) = dataptr2d_sdate(ig,n) + crop_inst%rx_swindow_starts_thisyr_patch(p,1) = dataptr2d_swindow_start(ig,n) + crop_inst%rx_swindow_ends_thisyr_patch (p,1) = dataptr2d_swindow_end (ig,n) ! Sanity check: Should only read in valid values - if (crop_inst%rx_sdates_thisyr_patch(p,1) > 365) then - write(iulog,'(a,i0,a,i0)') 'cropcal_interp(): Crop patch (ivt ',ivt,') has dataptr2d prescribed sowing date ',& - crop_inst%rx_sdates_thisyr_patch(p,1) + if (crop_inst%rx_swindow_starts_thisyr_patch(p,1) > 365) then + write(iulog,'(a,i0,a,i0)') 'cropcal_interp(): Crop patch (ivt ',ivt,') has dataptr2d prescribed sowing window start date ',& + crop_inst%rx_swindow_starts_thisyr_patch(p,1) call ESMF_Finalize(endflag=ESMF_END_ABORT) end if - - ! Only for first sowing date of the year + if (crop_inst%rx_swindow_ends_thisyr_patch(p,1) > 365) then + write(iulog,'(a,i0,a,i0)') 'cropcal_interp(): Crop patch (ivt ',ivt,') has dataptr2d prescribed sowing window end date ',& + crop_inst%rx_swindow_ends_thisyr_patch(p,1) + call ESMF_Finalize(endflag=ESMF_END_ABORT) + end if + + ! Only for first sowing window of the year ! The conditional here is to ensure nothing weird happens if it's called incorrectly on day 365 if (crop_inst%sdates_thisyr_patch(p,1) <= 0) then - crop_inst%next_rx_sdate_patch(p) = crop_inst%rx_sdates_thisyr_patch(p,1) + crop_inst%next_rx_swindow_start_patch(p) = crop_inst%rx_swindow_starts_thisyr_patch(p,1) + crop_inst%next_rx_swindow_end_patch (p) = crop_inst%rx_swindow_ends_thisyr_patch (p,1) end if else - write(iulog,'(a,i0)') 'cropcal_interp(), rx_sdates: Crop patch has ivt ',ivt + write(iulog,'(a,i0)') 'cropcal_interp(), prescribed sowing windows: Crop patch has ivt ',ivt call ESMF_Finalize(endflag=ESMF_END_ABORT) endif end do - end if ! use_cropcal_rx_sdates - deallocate(dataptr2d_sdate) + end if ! use_cropcal_rx_swindows + deallocate(dataptr2d_swindow_start) + deallocate(dataptr2d_swindow_end) allocate(dataptr2d_cultivar_gdds(lsize, ncft)) if (use_cropcal_rx_cultivar_gdds) then diff --git a/src/main/clm_varctl.F90 b/src/main/clm_varctl.F90 index 7fe5d24bba..1ec698b68a 100644 --- a/src/main/clm_varctl.F90 +++ b/src/main/clm_varctl.F90 @@ -320,7 +320,7 @@ module clm_varctl !---------------------------------------------------------- logical, public :: use_cropcal_streams = .false. - logical, public :: use_cropcal_rx_sdates = .false. + logical, public :: use_cropcal_rx_swindows = .false. logical, public :: use_cropcal_rx_cultivar_gdds = .false. !---------------------------------------------------------- From 87c1d21f710d66dff48632bb2a790cc5abd1b1cf Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Mon, 18 Sep 2023 16:44:04 -0600 Subject: [PATCH 02/52] Added SWINDOW_STARTS and SWINDOW_ENDS outputs. --- src/biogeochem/CNPhenologyMod.F90 | 14 +++++++++++++- src/biogeochem/CropType.F90 | 14 ++++++++++++++ 2 files changed, 27 insertions(+), 1 deletion(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index 15863774a7..369f304c19 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -1892,6 +1892,8 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & harvest_count(p) = 0 do s = 1, mxsowings crop_inst%sdates_thisyr_patch(p,s) = -1._r8 + crop_inst%swindow_starts_thisyr_patch(p,s) = -1._r8 + crop_inst%swindow_ends_thisyr_patch (p,s) = -1._r8 crop_inst%sowing_reason_thisyr_patch(p,s) = -1._r8 end do do s = 1, mxharvests @@ -1914,7 +1916,9 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & next_rx_swindow_end (p) = crop_inst%rx_swindow_ends_thisyr_patch (p,1) end if - ! Get dates of next prescribed sowing window (if any) + ! Get dates of next prescribed sowing window (if any). + ! This can probably be converted to a recursive subroutine. + ! TODO: Rework this to not rely on sowing count, because if you missed a sowing window it's not "sowing count" if (sowing_count(p) < mxsowings) then next_rx_swindow_start(p) = crop_inst%rx_swindow_starts_thisyr_patch(p,sowing_count(p)+1) next_rx_swindow_end (p) = crop_inst%rx_swindow_ends_thisyr_patch (p,sowing_count(p)+1) @@ -1938,6 +1942,7 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & end if ! This is outside the croplive check so that the "harvest if planting conditions were met today" conditional works. + ! TODO: Rework this to not rely on sowing count, because if you missed a sowing window it's not "sowing count" if (use_cropcal_rx_swindows) then if (sowing_count(p) < mxsowings) then sowing_window_startdate = crop_inst%rx_swindow_starts_thisyr_patch(p,sowing_count(p)+1) @@ -1953,6 +1958,13 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & is_in_sowing_window = is_doy_in_interval(sowing_window_startdate, sowing_window_enddate, jday) is_end_sowing_window = jday == sowing_window_enddate ! + ! Save these diagnostic variables only on the first day of the window to ensure that windows spanning the new year aren't double-counted. + ! TODO: Rework this to not rely on sowing count, because if you missed a sowing window it's not "sowing count" + if (jday == sowing_window_startdate .and. sowing_count(p) < mxsowings) then + crop_inst%swindow_starts_thisyr_patch(p,sowing_count(p)+1) = sowing_window_startdate + crop_inst%swindow_ends_thisyr_patch (p,sowing_count(p)+1) = sowing_window_enddate + end if + ! ! Only allow sowing according to normal "window" rules if not using prescribed ! sowing windows at all, or if this cell had no values in either of the ! prescribed sowing window files. (Note that if either sowing window start was diff --git a/src/biogeochem/CropType.F90 b/src/biogeochem/CropType.F90 index 564f4babc8..2239d3e9bd 100644 --- a/src/biogeochem/CropType.F90 +++ b/src/biogeochem/CropType.F90 @@ -53,6 +53,8 @@ module CropType integer , pointer :: rx_swindow_ends_thisyr_patch (:,:) ! all prescribed sowing window end dates for this patch this year (day of year) [patch, mxsowings] real(r8), pointer :: rx_cultivar_gdds_thisyr_patch (:,:) ! all cultivar GDD targets for this patch this year (ddays) [patch, mxsowings] real(r8), pointer :: sdates_thisyr_patch (:,:) ! all actual sowing dates for this patch this year (day of year) [patch, mxsowings] + real(r8), pointer :: swindow_starts_thisyr_patch(:,:) ! all sowing window start dates for this patch this year (day of year) [patch, mxsowings] + real(r8), pointer :: swindow_ends_thisyr_patch (:,:) ! all sowing window end dates for this patch this year (day of year) [patch, mxsowings] real(r8), pointer :: sdates_perharv_patch (:,:) ! all actual sowing dates for crops *harvested* this year (day of year) [patch, mxharvests] real(r8), pointer :: syears_perharv_patch (:,:) ! all actual sowing years for crops *harvested* this year (day of year) [patch, mxharvests] real(r8), pointer :: hdates_thisyr_patch (:,:) ! all actual harvest dates for this patch this year (day of year) [patch, mxharvests] @@ -235,6 +237,8 @@ subroutine InitAllocate(this, bounds) allocate(this%rx_swindow_ends_thisyr_patch(begp:endp,1:mxsowings)) ; this%rx_swindow_ends_thisyr_patch (:,:) = -1 allocate(this%rx_cultivar_gdds_thisyr_patch(begp:endp,1:mxsowings)) ; this%rx_cultivar_gdds_thisyr_patch(:,:) = spval allocate(this%sdates_thisyr_patch(begp:endp,1:mxsowings)) ; this%sdates_thisyr_patch(:,:) = spval + allocate(this%swindow_starts_thisyr_patch(begp:endp,1:mxsowings)) ; this%swindow_starts_thisyr_patch(:,:) = spval + allocate(this%swindow_ends_thisyr_patch (begp:endp,1:mxsowings)) ; this%swindow_ends_thisyr_patch (:,:) = spval allocate(this%sdates_perharv_patch(begp:endp,1:mxharvests)) ; this%sdates_perharv_patch(:,:) = spval allocate(this%syears_perharv_patch(begp:endp,1:mxharvests)) ; this%syears_perharv_patch(:,:) = spval allocate(this%hdates_thisyr_patch(begp:endp,1:mxharvests)) ; this%hdates_thisyr_patch(:,:) = spval @@ -307,6 +311,16 @@ subroutine InitHistory(this, bounds) avgflag='I', long_name='actual crop sowing dates; should only be output annually', & ptr_patch=this%sdates_thisyr_patch, default='inactive') + this%swindow_starts_thisyr_patch(begp:endp,:) = spval + call hist_addfld2d (fname='SWINDOW_STARTS', units='day of year', type2d='mxsowings', & + avgflag='I', long_name='crop sowing window start dates; should only be output annually', & + ptr_patch=this%swindow_starts_thisyr_patch, default='inactive') + + this%swindow_ends_thisyr_patch(begp:endp,:) = spval + call hist_addfld2d (fname='SWINDOW_ENDS', units='day of year', type2d='mxsowings', & + avgflag='I', long_name='crop sowing window end dates; should only be output annually', & + ptr_patch=this%swindow_ends_thisyr_patch, default='inactive') + this%sdates_perharv_patch(begp:endp,:) = spval call hist_addfld2d (fname='SDATES_PERHARV', units='day of year', type2d='mxharvests', & avgflag='I', long_name='actual sowing dates for crops harvested this year; should only be output annually', & From d9b7cdf22a74bb26b328e5137309b3cf88755d87 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Mon, 18 Sep 2023 16:54:21 -0600 Subject: [PATCH 03/52] CropPhenology() can now handle missed sowing windows. # Conflicts: # src/biogeochem/CNPhenologyMod.F90 --- src/biogeochem/CNPhenologyMod.F90 | 8 ++++++++ 1 file changed, 8 insertions(+) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index 369f304c19..68e3c446ef 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -1922,6 +1922,14 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & if (sowing_count(p) < mxsowings) then next_rx_swindow_start(p) = crop_inst%rx_swindow_starts_thisyr_patch(p,sowing_count(p)+1) next_rx_swindow_end (p) = crop_inst%rx_swindow_ends_thisyr_patch (p,sowing_count(p)+1) + ! Handle a missed sowing window + if (next_rx_swindow_start(p) < next_rx_swindow_end(p) .and. jday > next_rx_swindow_end(p)) then + sowing_count(p) = sowing_count(p) + 1 + if (sowing_count(p) < mxsowings) then + next_rx_swindow_start(p) = crop_inst%rx_swindow_starts_thisyr_patch(p,sowing_count(p)+1) + next_rx_swindow_end (p) = crop_inst%rx_swindow_ends_thisyr_patch (p,sowing_count(p)+1) + end if + end if end if ! CropType->InitAllocate() initializes next_rx_swindow_start to -1. It's only changed from that, by cropCalStreamMod->cropcal_interp(), when use_cropcal_rx_swindows is true. So if not using prescribed sowing windows, this boolean will always be false, because jday can never be -1. From c3fd7f3a6eada961aa625399f03b94b888cc1085 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Tue, 19 Sep 2023 13:53:05 -0600 Subject: [PATCH 04/52] Prevent crops from being sown twice in one sowing window. --- src/biogeochem/CNPhenologyMod.F90 | 12 +++++++++++- src/biogeochem/CropType.F90 | 27 +++++++++++++++++++++++++++ 2 files changed, 38 insertions(+), 1 deletion(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index 68e3c446ef..6336c926fd 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -1726,7 +1726,7 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & use clm_time_manager , only : get_prev_calday, get_curr_days_per_year, is_beg_curr_year use clm_time_manager , only : get_average_days_per_year use clm_time_manager , only : get_prev_date - use clm_time_manager , only : is_doy_in_interval + use clm_time_manager , only : is_doy_in_interval, is_end_curr_day use pftconMod , only : ntmp_corn, nswheat, nwwheat, ntmp_soybean use pftconMod , only : nirrig_tmp_corn, nirrig_swheat, nirrig_wwheat, nirrig_tmp_soybean use pftconMod , only : ntrp_corn, nsugarcane, ntrp_soybean, ncotton, nrice @@ -1964,6 +1964,10 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & sowing_window_enddate = maxplantjday(ivt(p),h) end if is_in_sowing_window = is_doy_in_interval(sowing_window_startdate, sowing_window_enddate, jday) + if (crop_inst%sown_in_this_window(p) .and. .not. is_in_sowing_window) then + ! Probably unnecessary since it's set to false in last timestep of sowing window at the end of CropPhenology() + crop_inst%sown_in_this_window(p) = .false. + end if is_end_sowing_window = jday == sowing_window_enddate ! ! Save these diagnostic variables only on the first day of the window to ensure that windows spanning the new year aren't double-counted. @@ -2388,6 +2392,11 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & endif end if ! croplive + ! At the end of the sowing window, AFTER we've done everything crop-related, set this to false + if (is_end_sowing_window .and. is_end_curr_day()) then + crop_inst%sown_in_this_window(p) = .false. + end if + end do ! prognostic crops loop end associate @@ -2603,6 +2612,7 @@ subroutine PlantCrop(p, leafcn_in, jday, kyr, do_plant_normal, & ! impose limit on growing season length needed ! for crop maturity - for cold weather constraints croplive(p) = .true. + crop_inst%sown_in_this_window(p) = .true. idop(p) = jday iyop(p) = kyr harvdate(p) = NOT_Harvested diff --git a/src/biogeochem/CropType.F90 b/src/biogeochem/CropType.F90 index 2239d3e9bd..f1e0bba7cf 100644 --- a/src/biogeochem/CropType.F90 +++ b/src/biogeochem/CropType.F90 @@ -47,6 +47,7 @@ module CropType character(len=20) :: baset_mapping real(r8) :: baset_latvary_intercept real(r8) :: baset_latvary_slope + logical , pointer :: sown_in_this_window (:) ! patch flag. True if the crop has already been sown during the current sowing window. False otherwise or if not in a sowing window. integer , pointer :: next_rx_swindow_start_patch (:) ! start of prescribed sowing window for the next growing season this year integer , pointer :: next_rx_swindow_end_patch (:) ! end of prescribed sowing window for the next growing season this year integer , pointer :: rx_swindow_starts_thisyr_patch(:,:) ! all prescribed sowing window start dates for this patch this year (day of year) [patch, mxsowings] @@ -231,6 +232,7 @@ subroutine InitAllocate(this, bounds) allocate(this%cphase_patch (begp:endp)) ; this%cphase_patch (:) = cphase_not_planted allocate(this%sowing_reason_patch (begp:endp)) ; this%sowing_reason_patch (:) = -1 allocate(this%latbaset_patch (begp:endp)) ; this%latbaset_patch (:) = spval + allocate(this%sown_in_this_window(begp:endp)) ; this%sown_in_this_window(:) = .false. allocate(this%next_rx_swindow_start_patch(begp:endp)) ; this%next_rx_swindow_start_patch(:) = -1 allocate(this%next_rx_swindow_end_patch (begp:endp)) ; this%next_rx_swindow_end_patch (:) = -1 allocate(this%rx_swindow_starts_thisyr_patch(begp:endp,1:mxsowings)); this%rx_swindow_starts_thisyr_patch(:,:) = -1 @@ -605,6 +607,31 @@ subroutine Restart(this, bounds, ncid, cnveg_state_inst, flag) end if deallocate(temp1d) + allocate(temp1d(bounds%begp:bounds%endp)) + if (flag == 'write') then + do p= bounds%begp,bounds%endp + if (this%sown_in_this_window(p)) then + temp1d(p) = 1 + else + temp1d(p) = 0 + end if + end do + end if + call restartvar(ncid=ncid, flag=flag, varname='sown_in_this_window', xtype=ncd_log, & + dim1name='pft', & + long_name='Flag that patch was sown already during the current sowing window', & + interpinic_flag='interp', readvar=readvar, data=temp1d) + if (flag == 'read') then + do p= bounds%begp,bounds%endp + if (temp1d(p) == 1) then + this%sown_in_this_window(p) = .true. + else + this%sown_in_this_window(p) = .false. + end if + end do + end if + deallocate(temp1d) + call restartvar(ncid=ncid, flag=flag, varname='harvdate', xtype=ncd_int, & dim1name='pft', long_name='harvest date', units='jday', nvalid_range=(/1,366/), & interpinic_flag='interp', readvar=readvar, data=this%harvdate_patch) From 1b55382b1d4ab52114cb5e4b4397058d1884201a Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 20 Sep 2023 16:19:14 -0600 Subject: [PATCH 05/52] Actually USE sown_in_this_window. --- src/biogeochem/CNPhenologyMod.F90 | 1 + 1 file changed, 1 insertion(+) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index 6336c926fd..d85ac7e015 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -2018,6 +2018,7 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & gdd820(p) /= spval end if do_plant = do_plant_prescribed .or. do_plant_normal .or. do_plant_lastchance + do_plant = do_plant .and. .not. crop_inst%sown_in_this_window(p) did_plant = .false. ! Once outputs can handle >1 planting per year, remove 2nd condition. From dcc55d604ab3a31c9813d31653f572d892b15aa4 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Mon, 25 Sep 2023 22:11:35 -0600 Subject: [PATCH 06/52] Use paramfile if either rx sowing window value is missing. --- src/biogeochem/CNPhenologyMod.F90 | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index d85ac7e015..076fd7e27e 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -1959,7 +1959,10 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & sowing_window_startdate = 999 sowing_window_enddate = 999 end if - else + end if + ! Fall back on parameter file in case either prescribed sowing window value is missing, or if not using prescribed sowing windows at all. + ! TODO: Test whether ".or. .not. use_cropcal_rx_swindows" is necessary. I think not: in that case, rx_swindow_starts_thisyr_patch and rx_swindow_ends_thisyr_patch are initialized to and stay -1. + if (min(sowing_window_startdate, sowing_window_enddate) <= 0 .or. .not. use_cropcal_rx_swindows) then sowing_window_startdate = minplantjday(ivt(p),h) sowing_window_enddate = maxplantjday(ivt(p),h) end if From 6e110401181331b207ed4e6ff27f48f8d64abf5e Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 27 Sep 2023 08:56:12 -0600 Subject: [PATCH 07/52] Add subroutines get_doy_tomorrow() and get_swindow(), plus unit tests. Currently unused. --- src/biogeochem/CNPhenologyMod.F90 | 99 +++++++++++++++++++ .../test/CNPhenology_test/test_CNPhenology.pf | 94 ++++++++++++++++++ src/cpl/share_esmf/cropcalStreamMod.F90 | 6 ++ src/utils/clm_time_manager.F90 | 18 ++++ .../test_clm_time_manager.pf | 13 +++ 5 files changed, 230 insertions(+) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index 076fd7e27e..a3c65e6fb9 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -58,6 +58,7 @@ module CNPhenologyMod public :: CNPhenologySetParams ! Set the parameters explicitly for unit tests public :: SeasonalDecidOnset ! Logical function to determine is seasonal decidious onset should be triggered public :: SeasonalCriticalDaylength ! Critical day length needed for Seasonal decidious offset + public :: get_swindow ! !PRIVITE MEMBER FIUNCTIONS: private :: CNPhenologyClimate ! Get climatological everages to figure out triggers for Phenology @@ -1712,6 +1713,104 @@ subroutine CNStressDecidPhenology (num_soilp, filter_soilp , & end subroutine CNStressDecidPhenology + + !----------------------------------------------------------------------- + subroutine get_swindow(jday, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + ! !DESCRIPTION: + ! Determine when the "next" sowing window is. This is either the sowing window we are + ! currently in or, if not in a sowing window, the next one that will occur. + + ! !USES: + use clm_time_manager , only : get_curr_days_per_year, is_doy_in_interval, get_doy_tomorrow + ! !ARGUMENTS: + integer, intent(in) :: jday ! Day of year + integer, dimension(:), intent(in) :: rx_starts, rx_ends ! All prescribed sowing window start and end dates for this patch + integer, intent(in) :: param_start, param_end ! Sowing window start and end dates from parameter file + integer, intent(out) :: w ! Index of "next" sowing window + integer, intent(out) :: start_w, end_w ! Start and end dates of "next" sowing window + logical, intent(out) :: sowing_window_starts_tomorrow + ! + ! !LOCAL VARIABLES + integer :: next_swindow_start + integer :: i, x + integer :: jday_tomorrow + integer :: mxsowings_in ! Due to unit testing, we can't assume the length of the rx sowing window arrays is mxsowings as set in clm_varpar + + ! Initialize + w = -1 + start_w = -1 + end_w = -1 + + ! Get info + jday_tomorrow = get_doy_tomorrow(jday) + mxsowings_in = size(rx_starts) + + ! If no sowing windows are prescribed, use the values from the parameter file. + if (maxval(rx_starts) < 1) then + w = 1 + start_w = param_start + end_w = param_end + + ! Otherwise, if today is after the latest sowing window end date, use the first sowing window. This works only if sowing windows that span the new year are located at index w = 1. + else if (jday > maxval(rx_ends)) then + w = 1 + start_w = rx_starts(w) + end_w = rx_ends(w) + end if + + ! If we already got sowing window dates, do this and exit. + if (w > 0) then + sowing_window_starts_tomorrow = start_w == jday_tomorrow + return + end if + + ! Otherwise, use the first prescribed sowing window we find whose end is >= today. + do w = 1, mxsowings_in + ! If nothing prescribed at this w, stop looking and exit loop. + if (min(rx_starts(w), rx_ends(w)) < 0) then + exit + end if + + if (jday <= rx_ends(w)) then + start_w = rx_starts(w) + end_w = rx_ends(w) + exit + end if + end do + + ! Ensure that a window was found + if (start_w < 1 .or. end_w < 1) then + call endrun(msg="get_swindow(): No sowing window found") + end if + + ! Get the start date of the NEXT sowing window (not including the sowing window we're currently in, if any) + if (is_doy_in_interval(start_w, end_w, jday)) then + next_swindow_start = -1 + x = w + i = 0 ! For checking infinite loop + do while (next_swindow_start < 1) + ! Check for infinite loop + i = i + 1 + if (i > mxsowings_in + 1) then + call endrun("Infinite loop in get_swindow()") + end if + + x = x + 1 + if (x > mxsowings_in) then + x = 1 + end if + next_swindow_start = rx_starts(x) + end do + else + next_swindow_start = start_w + end if + + ! Does the NEXT sowing window start tomorrow? + sowing_window_starts_tomorrow = next_swindow_start == jday_tomorrow + + end subroutine get_swindow + + !----------------------------------------------------------------------- subroutine CropPhenology(num_pcropp, filter_pcropp , & waterdiagnosticbulk_inst, temperature_inst, crop_inst, canopystate_inst, cnveg_state_inst , & diff --git a/src/biogeochem/test/CNPhenology_test/test_CNPhenology.pf b/src/biogeochem/test/CNPhenology_test/test_CNPhenology.pf index 135364e474..ebab50ad13 100644 --- a/src/biogeochem/test/CNPhenology_test/test_CNPhenology.pf +++ b/src/biogeochem/test/CNPhenology_test/test_CNPhenology.pf @@ -293,6 +293,100 @@ contains end do end subroutine check_crit_dayl_dependsonlatnveg + @Test + subroutine test_get_swindow_startend(this) + use clm_time_manager, only : get_curr_days_per_year + class(TestCNPhenology), intent(inout) :: this + integer, dimension(3), parameter :: rx_starts = (/1, 150, -1/) + integer, dimension(3), parameter :: rx_ends = (/45, 180, -1/) + integer, parameter :: param_start = 200 + integer, parameter :: param_end = 250 + integer :: w, start_w, end_w + logical :: sowing_window_starts_tomorrow + + call get_swindow(1, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + @assertEqual(1, w) + @assertEqual(1, start_w) + @assertEqual(45, end_w) + @assertFalse(sowing_window_starts_tomorrow) + + call get_swindow(45, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + @assertEqual(1, w) + @assertEqual(1, start_w) + @assertEqual(45, end_w) + @assertFalse(sowing_window_starts_tomorrow) + + call get_swindow(149, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + @assertEqual(2, w) + @assertEqual(150, start_w) + @assertEqual(180, end_w) + @assertTrue(sowing_window_starts_tomorrow) + + call get_swindow(175, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + @assertEqual(2, w) + @assertEqual(150, start_w) + @assertEqual(180, end_w) + @assertFalse(sowing_window_starts_tomorrow) + + call get_swindow(229, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + @assertEqual(1, w) + @assertEqual(1, start_w) + @assertEqual(45, end_w) + @assertFalse(sowing_window_starts_tomorrow) + + call get_swindow(get_curr_days_per_year(), rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + @assertEqual(1, w) + @assertEqual(1, start_w) + @assertEqual(45, end_w) + @assertTrue(sowing_window_starts_tomorrow) + end subroutine test_get_swindow_startend + @Test + subroutine test_get_swindow_endstart(this) + use clm_time_manager, only : get_curr_days_per_year + class(TestCNPhenology), intent(inout) :: this + integer, dimension(3), parameter :: rx_starts = (/360, 150, -1/) + integer, dimension(3), parameter :: rx_ends = (/45, 180, -1/) + integer, parameter :: param_start = 200 + integer, parameter :: param_end = 250 + integer :: w, start_w, end_w + logical :: sowing_window_starts_tomorrow + + call get_swindow(1, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + @assertEqual(1, w) + @assertEqual(360, start_w) + @assertEqual(45, end_w) + @assertFalse(sowing_window_starts_tomorrow) + + call get_swindow(45, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + @assertEqual(1, w) + @assertEqual(start_w, 360) + @assertEqual(end_w, 45) + @assertFalse(sowing_window_starts_tomorrow) + + call get_swindow(149, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + @assertEqual(2, w) + @assertEqual(150, start_w) + @assertEqual(180, end_w) + @assertTrue(sowing_window_starts_tomorrow) + + call get_swindow(175, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + @assertEqual(2, w) + @assertEqual(150, start_w) + @assertEqual(180, end_w) + @assertFalse(sowing_window_starts_tomorrow) + + call get_swindow(229, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + @assertEqual(1, w) + @assertEqual(360, start_w) + @assertEqual(45, end_w) + @assertFalse(sowing_window_starts_tomorrow) + + call get_swindow(get_curr_days_per_year(), rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + @assertEqual(1, w) + @assertEqual(360, start_w) + @assertEqual(45, end_w) + @assertFalse(sowing_window_starts_tomorrow) + end subroutine test_get_swindow_endstart end module test_CNPhenology diff --git a/src/cpl/share_esmf/cropcalStreamMod.F90 b/src/cpl/share_esmf/cropcalStreamMod.F90 index 91f565771d..5f91cc391e 100644 --- a/src/cpl/share_esmf/cropcalStreamMod.F90 +++ b/src/cpl/share_esmf/cropcalStreamMod.F90 @@ -389,6 +389,7 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) ! Only for first sowing window of the year ! The conditional here is to ensure nothing weird happens if it's called incorrectly on day 365 if (crop_inst%sdates_thisyr_patch(p,1) <= 0) then + ! TODO: Add handling of mid-year restarts crop_inst%next_rx_swindow_start_patch(p) = crop_inst%rx_swindow_starts_thisyr_patch(p,1) crop_inst%next_rx_swindow_end_patch (p) = crop_inst%rx_swindow_ends_thisyr_patch (p,1) end if @@ -397,6 +398,11 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) call ESMF_Finalize(endflag=ESMF_END_ABORT) endif end do + + ! TODO: Ensure that, if mxsowings > 1, sowing windows are ordered such that ENDS are monotonically increasing. This is necessary because of how get_swindow() works. + + ! TODO: Fail if a sowing window start date is prescribed without an end date (or vice versa) + end if ! use_cropcal_rx_swindows deallocate(dataptr2d_swindow_start) deallocate(dataptr2d_swindow_end) diff --git a/src/utils/clm_time_manager.F90 b/src/utils/clm_time_manager.F90 index 5c65f5decd..14572935c3 100644 --- a/src/utils/clm_time_manager.F90 +++ b/src/utils/clm_time_manager.F90 @@ -42,6 +42,7 @@ module clm_time_manager get_prev_calday, &! return calendar day at beginning of current timestep get_calday, &! return calendar day from input date get_calendar, &! return calendar + get_doy_tomorrow, &! return next day of year get_average_days_per_year,&! return the average number of days per year for the given calendar get_curr_days_per_year, &! return the days per year for year as of the end of the current time step get_prev_days_per_year, &! return the days per year for year as of the beginning of the current time step @@ -1274,6 +1275,23 @@ end function get_calendar !========================================================================================= + function get_doy_tomorrow(doy_today) result(doy_tomorrow) + + !--------------------------------------------------------------------------------- + ! Given a day of the year (doy_today), return the next day of the year + + integer, intent(in) :: doy_today + integer :: doy_tomorrow + + if (doy_today == get_curr_days_per_year()) then + doy_tomorrow = 1 + else + doy_tomorrow = doy_today + 1 + end if + end function get_doy_tomorrow + + !========================================================================================= + real(r8) function get_average_days_per_year() !--------------------------------------------------------------------------------- diff --git a/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf b/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf index 78565fd54d..0665e5f6d0 100644 --- a/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf +++ b/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf @@ -658,4 +658,17 @@ contains end subroutine check_is_today_in_doy_interval + @Test + subroutine check_get_doy_tomorrow(this) + class(TestTimeManager), intent(inout) :: this + + ! We don't care about the actual date here; we just want to enable get_curr_days_per_year() + call unittest_timemgr_setup(dtime=dtime, use_gregorian_calendar=.true.) + call set_date(yr=2009, mon=10, day=28, tod=dtime) + + @assertEqual(get_doy_tomorrow(1), 2) + @assertEqual(get_doy_tomorrow(150), 151) + @assertEqual(get_doy_tomorrow(get_curr_days_per_year()), 1) + end subroutine check_get_doy_tomorrow + end module test_clm_time_manager From 9104f88cd2c4e8691ae7f704e9447a5e882bf2af Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 27 Sep 2023 10:07:11 -0600 Subject: [PATCH 08/52] Use get_swindow() in CropPhenology(). * Allows significant cleanup there and hopefully fixes a bug. * Removed from CropType: next_rx_swindow_start, next_rx_swindow_end. --- src/biogeochem/CNPhenologyMod.F90 | 84 ++++++------------------- src/biogeochem/CropType.F90 | 4 -- src/cpl/share_esmf/cropcalStreamMod.F90 | 8 --- 3 files changed, 20 insertions(+), 76 deletions(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index a3c65e6fb9..ff01a3a2f3 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -1862,6 +1862,7 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & integer h ! hemisphere indices integer s ! growing season indices integer k ! grain pool indices + integer w ! sowing window index integer idpp ! number of days past planting integer mxmat ! maximum growing season length integer kyr ! current year @@ -1870,6 +1871,7 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & integer mcsec ! seconds of day (0, ..., seconds/day) integer sowing_window_startdate ! date (day of year) of first day of sowing window integer sowing_window_enddate ! date (day of year) of last day of sowing window + logical sowing_window_starts_tomorrow ! Is tomorrow the start of a sowing window? real(r8) harvest_reason real(r8) dayspyr ! days per year in this year real(r8) avg_dayspyr ! average number of days per year @@ -1887,7 +1889,6 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & logical do_harvest ! Are harvest conditions satisfied? logical fake_harvest ! Dealing with incorrect Dec. 31 planting logical did_plant_prescribed_today ! Was the crop sown today? - logical will_plant_prescribed_tomorrow ! Is tomorrow a prescribed sowing day? logical vernalization_forces_harvest ! Was the crop killed by freezing during vernalization? !------------------------------------------------------------------------ @@ -1916,8 +1917,6 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & harvdate => crop_inst%harvdate_patch , & ! Output: [integer (:) ] harvest date croplive => crop_inst%croplive_patch , & ! Output: [logical (:) ] Flag, true if planted, not harvested vf => crop_inst%vf_patch , & ! Output: [real(r8) (:) ] vernalization factor - next_rx_swindow_start => crop_inst%next_rx_swindow_start_patch , & ! Inout: [integer (:) ] start of prescribed sowing window of next growing season this year - next_rx_swindow_end => crop_inst%next_rx_swindow_end_patch , & ! Inout: [integer (:) ] end of prescribed sowing window of next growing season this year sowing_count => crop_inst%sowing_count , & ! Inout: [integer (:) ] number of sowing events this year for this patch harvest_count => crop_inst%harvest_count , & ! Inout: [integer (:) ] number of harvest events this year for this patch peaklai => cnveg_state_inst%peaklai_patch , & ! Output: [integer (:) ] 1: max allowed lai; 0: not at max @@ -2011,30 +2010,16 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & do k = repr_grain_min, repr_grain_max cnveg_carbonflux_inst%repr_grainc_to_food_thisyr_patch(p,k) = 0._r8 end do - next_rx_swindow_start(p) = crop_inst%rx_swindow_starts_thisyr_patch(p,1) - next_rx_swindow_end (p) = crop_inst%rx_swindow_ends_thisyr_patch (p,1) end if - ! Get dates of next prescribed sowing window (if any). - ! This can probably be converted to a recursive subroutine. - ! TODO: Rework this to not rely on sowing count, because if you missed a sowing window it's not "sowing count" - if (sowing_count(p) < mxsowings) then - next_rx_swindow_start(p) = crop_inst%rx_swindow_starts_thisyr_patch(p,sowing_count(p)+1) - next_rx_swindow_end (p) = crop_inst%rx_swindow_ends_thisyr_patch (p,sowing_count(p)+1) - ! Handle a missed sowing window - if (next_rx_swindow_start(p) < next_rx_swindow_end(p) .and. jday > next_rx_swindow_end(p)) then - sowing_count(p) = sowing_count(p) + 1 - if (sowing_count(p) < mxsowings) then - next_rx_swindow_start(p) = crop_inst%rx_swindow_starts_thisyr_patch(p,sowing_count(p)+1) - next_rx_swindow_end (p) = crop_inst%rx_swindow_ends_thisyr_patch (p,sowing_count(p)+1) - end if - end if - end if + ! Get dates of current or next sowing window. + call get_swindow(jday, crop_inst%rx_swindow_starts_thisyr_patch(p,:), crop_inst%rx_swindow_ends_thisyr_patch(p,:), minplantjday(ivt(p),h), maxplantjday(ivt(p),h), w, sowing_window_startdate, sowing_window_enddate, sowing_window_starts_tomorrow) - ! CropType->InitAllocate() initializes next_rx_swindow_start to -1. It's only changed from that, by cropCalStreamMod->cropcal_interp(), when use_cropcal_rx_swindows is true. So if not using prescribed sowing windows, this boolean will always be false, because jday can never be -1. - do_plant_prescribed = next_rx_swindow_start(p) == jday .and. sowing_count(p) < mxsowings - ! We only want to plant on a specific day if the prescribed sowing window starts AND ends on the same day. - do_plant_prescribed = do_plant_prescribed .and. next_rx_swindow_start(p) == next_rx_swindow_end(p) + ! We only want to plant on a specific day if the prescribed sowing window starts AND ends on the same day. Also make sure we haven't planted yet today. + ! TODO: Change last condition to `maxval(crop_inst%sdates_perharv_patch(p,:)) /= jday` + do_plant_prescribed = sowing_window_startdate == jday .and. & + sowing_window_enddate == jday .and. & + sowing_count (p) < mxsowings ! BACKWARDS_COMPATIBILITY(wjs/ssr, 2022-02-18) ! When resuming from a run with old code, may need to manually set these. @@ -2048,35 +2033,19 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & crop_inst%sdates_thisyr_patch(p,1) = real(idop(p), r8) end if + ! Are we currently in a sowing window? ! This is outside the croplive check so that the "harvest if planting conditions were met today" conditional works. - ! TODO: Rework this to not rely on sowing count, because if you missed a sowing window it's not "sowing count" - if (use_cropcal_rx_swindows) then - if (sowing_count(p) < mxsowings) then - sowing_window_startdate = crop_inst%rx_swindow_starts_thisyr_patch(p,sowing_count(p)+1) - sowing_window_enddate = crop_inst%rx_swindow_ends_thisyr_patch (p,sowing_count(p)+1) - else ! Do not allow more sowings after reaching max number of sowings (mxsowings) - sowing_window_startdate = 999 - sowing_window_enddate = 999 - end if - end if - ! Fall back on parameter file in case either prescribed sowing window value is missing, or if not using prescribed sowing windows at all. - ! TODO: Test whether ".or. .not. use_cropcal_rx_swindows" is necessary. I think not: in that case, rx_swindow_starts_thisyr_patch and rx_swindow_ends_thisyr_patch are initialized to and stay -1. - if (min(sowing_window_startdate, sowing_window_enddate) <= 0 .or. .not. use_cropcal_rx_swindows) then - sowing_window_startdate = minplantjday(ivt(p),h) - sowing_window_enddate = maxplantjday(ivt(p),h) - end if is_in_sowing_window = is_doy_in_interval(sowing_window_startdate, sowing_window_enddate, jday) if (crop_inst%sown_in_this_window(p) .and. .not. is_in_sowing_window) then - ! Probably unnecessary since it's set to false in last timestep of sowing window at the end of CropPhenology() + ! TODO: Test whether this is necessary, since it's set to false in last timestep of sowing window at the end of CropPhenology(). crop_inst%sown_in_this_window(p) = .false. end if is_end_sowing_window = jday == sowing_window_enddate ! ! Save these diagnostic variables only on the first day of the window to ensure that windows spanning the new year aren't double-counted. - ! TODO: Rework this to not rely on sowing count, because if you missed a sowing window it's not "sowing count" - if (jday == sowing_window_startdate .and. sowing_count(p) < mxsowings) then - crop_inst%swindow_starts_thisyr_patch(p,sowing_count(p)+1) = sowing_window_startdate - crop_inst%swindow_ends_thisyr_patch (p,sowing_count(p)+1) = sowing_window_enddate + if (jday == sowing_window_startdate) then + crop_inst%swindow_starts_thisyr_patch(p,w) = sowing_window_startdate + crop_inst%swindow_ends_thisyr_patch (p,w) = sowing_window_enddate end if ! ! Only allow sowing according to normal "window" rules if not using prescribed @@ -2308,10 +2277,11 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & harvest_reason = HARVEST_REASON_SOWTODAY ! If generate_crop_gdds and this patch has prescribed sowing inputs + ! TODO: Shouldn't this also apply to paramfile sowing windows? else if (generate_crop_gdds .and. crop_inst%rx_swindow_starts_thisyr_patch(p,1) .gt. 0) then - if (next_rx_swindow_start(p) >= 0) then - ! Harvest the day before the start of the next sowing window this year. - do_harvest = jday == next_rx_swindow_start(p) - 1 + if (sowing_window_startdate >= 0) then + ! Harvest the day before the start of the next sowing window. + do_harvest = sowing_window_starts_tomorrow ! ... unless that will lead to growing season length 365 (or 366, ! if last year was a leap year). This would result in idop==jday, @@ -2363,21 +2333,13 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & ! WARNING: This implementation assumes that sowing windows don't change over time! ! In order to avoid this, you'd have to read this year's AND next year's prescribed ! sowing windows. - ! WARNING: This implementation assumes that all patches use prescribed sowing windows. - if (use_cropcal_rx_swindows) then - will_plant_prescribed_tomorrow = (jday == next_rx_swindow_start(p) - 1) .or. & - (crop_inst%sdates_thisyr_patch(p,1) == 1 .and. & - jday == dayspyr) - else - will_plant_prescribed_tomorrow = .false. - end if - do_harvest = do_harvest .or. will_plant_prescribed_tomorrow + do_harvest = do_harvest .or. sowing_window_starts_tomorrow if (hui(p) >= gddmaturity(p)) then harvest_reason = HARVEST_REASON_MATURE else if (idpp >= mxmat) then harvest_reason = HARVEST_REASON_MAXSEASLENGTH - else if (will_plant_prescribed_tomorrow) then + else if (sowing_window_starts_tomorrow) then harvest_reason = HARVEST_REASON_SOWTOMORROW end if endif @@ -2689,7 +2651,6 @@ subroutine PlantCrop(p, leafcn_in, jday, kyr, do_plant_normal, & ivt => patch%itype , & ! Input: [integer (:) ] patch vegetation type croplive => crop_inst%croplive_patch , & ! Output: [logical (:) ] Flag, true if planted, not harvested harvdate => crop_inst%harvdate_patch , & ! Output: [integer (:) ] harvest date - next_rx_swindow_start => crop_inst%next_rx_swindow_start_patch , & ! Inout: [integer (:) ] start of prescribed sowing window of next growing season this year sowing_count => crop_inst%sowing_count , & ! Inout: [integer (:) ] number of sowing events this year for this patch sowing_reason => crop_inst%sowing_reason_thisyr_patch , & ! Output: [real(r8) (:) ] reason for each sowing this year for this patch gddmaturity => cnveg_state_inst%gddmaturity_patch , & ! Output: [real(r8) (:) ] gdd needed to harvest @@ -2721,11 +2682,6 @@ subroutine PlantCrop(p, leafcn_in, jday, kyr, do_plant_normal, & harvdate(p) = NOT_Harvested sowing_count(p) = sowing_count(p) + 1 - if (use_cropcal_rx_swindows .and. sowing_count(p) < mxsowings) then - next_rx_swindow_start(p) = crop_inst%rx_swindow_starts_thisyr_patch(p, sowing_count(p)+1) - else - next_rx_swindow_start(p) = -1 - endif crop_inst%sdates_thisyr_patch(p,sowing_count(p)) = real(jday, r8) this_sowing_reason = 0._r8 diff --git a/src/biogeochem/CropType.F90 b/src/biogeochem/CropType.F90 index f1e0bba7cf..6d101df456 100644 --- a/src/biogeochem/CropType.F90 +++ b/src/biogeochem/CropType.F90 @@ -48,8 +48,6 @@ module CropType real(r8) :: baset_latvary_intercept real(r8) :: baset_latvary_slope logical , pointer :: sown_in_this_window (:) ! patch flag. True if the crop has already been sown during the current sowing window. False otherwise or if not in a sowing window. - integer , pointer :: next_rx_swindow_start_patch (:) ! start of prescribed sowing window for the next growing season this year - integer , pointer :: next_rx_swindow_end_patch (:) ! end of prescribed sowing window for the next growing season this year integer , pointer :: rx_swindow_starts_thisyr_patch(:,:) ! all prescribed sowing window start dates for this patch this year (day of year) [patch, mxsowings] integer , pointer :: rx_swindow_ends_thisyr_patch (:,:) ! all prescribed sowing window end dates for this patch this year (day of year) [patch, mxsowings] real(r8), pointer :: rx_cultivar_gdds_thisyr_patch (:,:) ! all cultivar GDD targets for this patch this year (ddays) [patch, mxsowings] @@ -233,8 +231,6 @@ subroutine InitAllocate(this, bounds) allocate(this%sowing_reason_patch (begp:endp)) ; this%sowing_reason_patch (:) = -1 allocate(this%latbaset_patch (begp:endp)) ; this%latbaset_patch (:) = spval allocate(this%sown_in_this_window(begp:endp)) ; this%sown_in_this_window(:) = .false. - allocate(this%next_rx_swindow_start_patch(begp:endp)) ; this%next_rx_swindow_start_patch(:) = -1 - allocate(this%next_rx_swindow_end_patch (begp:endp)) ; this%next_rx_swindow_end_patch (:) = -1 allocate(this%rx_swindow_starts_thisyr_patch(begp:endp,1:mxsowings)); this%rx_swindow_starts_thisyr_patch(:,:) = -1 allocate(this%rx_swindow_ends_thisyr_patch(begp:endp,1:mxsowings)) ; this%rx_swindow_ends_thisyr_patch (:,:) = -1 allocate(this%rx_cultivar_gdds_thisyr_patch(begp:endp,1:mxsowings)) ; this%rx_cultivar_gdds_thisyr_patch(:,:) = spval diff --git a/src/cpl/share_esmf/cropcalStreamMod.F90 b/src/cpl/share_esmf/cropcalStreamMod.F90 index 5f91cc391e..8ad162783c 100644 --- a/src/cpl/share_esmf/cropcalStreamMod.F90 +++ b/src/cpl/share_esmf/cropcalStreamMod.F90 @@ -385,14 +385,6 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) crop_inst%rx_swindow_ends_thisyr_patch(p,1) call ESMF_Finalize(endflag=ESMF_END_ABORT) end if - - ! Only for first sowing window of the year - ! The conditional here is to ensure nothing weird happens if it's called incorrectly on day 365 - if (crop_inst%sdates_thisyr_patch(p,1) <= 0) then - ! TODO: Add handling of mid-year restarts - crop_inst%next_rx_swindow_start_patch(p) = crop_inst%rx_swindow_starts_thisyr_patch(p,1) - crop_inst%next_rx_swindow_end_patch (p) = crop_inst%rx_swindow_ends_thisyr_patch (p,1) - end if else write(iulog,'(a,i0)') 'cropcal_interp(), prescribed sowing windows: Crop patch has ivt ',ivt call ESMF_Finalize(endflag=ESMF_END_ABORT) From e5a3dc0c1732c2fab8793ef66e61136964ea01f2 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 27 Sep 2023 10:41:26 -0600 Subject: [PATCH 09/52] Added a TODO. --- src/biogeochem/CNPhenologyMod.F90 | 1 + 1 file changed, 1 insertion(+) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index ff01a3a2f3..fc680d74c9 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -2043,6 +2043,7 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & is_end_sowing_window = jday == sowing_window_enddate ! ! Save these diagnostic variables only on the first day of the window to ensure that windows spanning the new year aren't double-counted. + ! TODO: Change this to the LAST day of the window, to ensure the same ordering as required for rx sowing window input files. if (jday == sowing_window_startdate) then crop_inst%swindow_starts_thisyr_patch(p,w) = sowing_window_startdate crop_inst%swindow_ends_thisyr_patch (p,w) = sowing_window_enddate From b0107fda54367d82dfa6c6c46435b755da3953a1 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Thu, 28 Sep 2023 21:01:13 -0600 Subject: [PATCH 10/52] Bugfix to allow_unprescribed_planting. --- src/biogeochem/CNPhenologyMod.F90 | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index fc680d74c9..3ab8ee1eae 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -2050,10 +2050,8 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & end if ! ! Only allow sowing according to normal "window" rules if not using prescribed - ! sowing windows at all, or if this cell had no values in either of the - ! prescribed sowing window files. (Note that if either sowing window start was - ! or end was missing a value, they are both set to negative values.) - allow_unprescribed_planting = (.not. use_cropcal_rx_swindows) .or. crop_inst%rx_swindow_starts_thisyr_patch(p,1)<0 + ! sowing dates. + allow_unprescribed_planting = sowing_window_startdate /= sowing_window_enddate if (sowing_count(p) == mxsowings) then do_plant_normal = .false. do_plant_lastchance = .false. From 1693c0d829428ed59f92a7cce1a8adb6d2345fa3 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Fri, 29 Sep 2023 09:05:10 -0600 Subject: [PATCH 11/52] Reverted behavior change re: HARVEST_REASON_SOWTOMORROW. Now only happens if patch actually has a prescribed sowing DAY, instead of happening at the start of the window. This had caused answer changes for sugarcane. --- src/biogeochem/CNPhenologyMod.F90 | 21 +++++++++++++++------ 1 file changed, 15 insertions(+), 6 deletions(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index 3ab8ee1eae..7b9471c67a 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -1826,6 +1826,7 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & use clm_time_manager , only : get_average_days_per_year use clm_time_manager , only : get_prev_date use clm_time_manager , only : is_doy_in_interval, is_end_curr_day + use clm_time_manager , only : get_doy_tomorrow use pftconMod , only : ntmp_corn, nswheat, nwwheat, ntmp_soybean use pftconMod , only : nirrig_tmp_corn, nirrig_swheat, nirrig_wwheat, nirrig_tmp_soybean use pftconMod , only : ntrp_corn, nsugarcane, ntrp_soybean, ncotton, nrice @@ -1877,12 +1878,14 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & real(r8) avg_dayspyr ! average number of days per year real(r8) crmcorn ! comparitive relative maturity for corn real(r8) ndays_on ! number of days to fertilize + logical has_rx_sowing_date ! does the crop have a single sowing date instead of a window? logical is_in_sowing_window ! is the crop in its sowing window? logical is_end_sowing_window ! is it the last day of the crop's sowing window? logical sowing_gdd_requirement_met ! has the gridcell historically been warm enough to support the crop? logical do_plant_normal ! are the normal planting rules defined and satisfied? logical do_plant_lastchance ! if not the above, what about relaxed rules for the last day of the planting window? logical do_plant_prescribed ! is today the prescribed sowing date? + logical do_plant_prescribed_tomorrow ! is tomorrow the prescribed sowing date? logical do_plant ! are we planting in this time step for any reason? logical did_plant ! did we plant the crop in this time step? logical allow_unprescribed_planting ! should crop be allowed to be planted according to sowing window rules? @@ -2017,9 +2020,15 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & ! We only want to plant on a specific day if the prescribed sowing window starts AND ends on the same day. Also make sure we haven't planted yet today. ! TODO: Change last condition to `maxval(crop_inst%sdates_perharv_patch(p,:)) /= jday` - do_plant_prescribed = sowing_window_startdate == jday .and. & - sowing_window_enddate == jday .and. & + ! TODO: ¿Allow use of NON-prescribed sowing with one-day-long windows? + has_rx_sowing_date = sowing_window_startdate == sowing_window_enddate + do_plant_prescribed = has_rx_sowing_date .and. & + sowing_window_startdate == jday .and. & sowing_count (p) < mxsowings + do_plant_prescribed_tomorrow = \ + has_rx_sowing_date .and. & + sowing_window_startdate == get_doy_tomorrow(jday) .and. & + sowing_count (p) < mxsowings ! BACKWARDS_COMPATIBILITY(wjs/ssr, 2022-02-18) ! When resuming from a run with old code, may need to manually set these. @@ -2051,7 +2060,7 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & ! ! Only allow sowing according to normal "window" rules if not using prescribed ! sowing dates. - allow_unprescribed_planting = sowing_window_startdate /= sowing_window_enddate + allow_unprescribed_planting = .not. has_rx_sowing_date if (sowing_count(p) == mxsowings) then do_plant_normal = .false. do_plant_lastchance = .false. @@ -2280,7 +2289,7 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & else if (generate_crop_gdds .and. crop_inst%rx_swindow_starts_thisyr_patch(p,1) .gt. 0) then if (sowing_window_startdate >= 0) then ! Harvest the day before the start of the next sowing window. - do_harvest = sowing_window_starts_tomorrow + do_harvest = do_plant_prescribed_tomorrow ! ... unless that will lead to growing season length 365 (or 366, ! if last year was a leap year). This would result in idop==jday, @@ -2332,13 +2341,13 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & ! WARNING: This implementation assumes that sowing windows don't change over time! ! In order to avoid this, you'd have to read this year's AND next year's prescribed ! sowing windows. - do_harvest = do_harvest .or. sowing_window_starts_tomorrow + do_harvest = do_harvest .or. do_plant_prescribed_tomorrow if (hui(p) >= gddmaturity(p)) then harvest_reason = HARVEST_REASON_MATURE else if (idpp >= mxmat) then harvest_reason = HARVEST_REASON_MAXSEASLENGTH - else if (sowing_window_starts_tomorrow) then + else if (do_plant_prescribed_tomorrow) then harvest_reason = HARVEST_REASON_SOWTOMORROW end if endif From bb577533b0486553ea6a439ffdac2037a7ad67d6 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Fri, 29 Sep 2023 12:30:23 -0600 Subject: [PATCH 12/52] Removed unused logical sowing_window_starts_tomorrow. --- src/biogeochem/CNPhenologyMod.F90 | 15 ++------ .../test/CNPhenology_test/test_CNPhenology.pf | 38 ++++++------------- 2 files changed, 15 insertions(+), 38 deletions(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index 7b9471c67a..bc93810d92 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -1715,7 +1715,7 @@ end subroutine CNStressDecidPhenology !----------------------------------------------------------------------- - subroutine get_swindow(jday, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + subroutine get_swindow(jday, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w) ! !DESCRIPTION: ! Determine when the "next" sowing window is. This is either the sowing window we are ! currently in or, if not in a sowing window, the next one that will occur. @@ -1728,7 +1728,6 @@ subroutine get_swindow(jday, rx_starts, rx_ends, param_start, param_end, w, star integer, intent(in) :: param_start, param_end ! Sowing window start and end dates from parameter file integer, intent(out) :: w ! Index of "next" sowing window integer, intent(out) :: start_w, end_w ! Start and end dates of "next" sowing window - logical, intent(out) :: sowing_window_starts_tomorrow ! ! !LOCAL VARIABLES integer :: next_swindow_start @@ -1750,17 +1749,13 @@ subroutine get_swindow(jday, rx_starts, rx_ends, param_start, param_end, w, star w = 1 start_w = param_start end_w = param_end + return ! Otherwise, if today is after the latest sowing window end date, use the first sowing window. This works only if sowing windows that span the new year are located at index w = 1. else if (jday > maxval(rx_ends)) then w = 1 start_w = rx_starts(w) end_w = rx_ends(w) - end if - - ! If we already got sowing window dates, do this and exit. - if (w > 0) then - sowing_window_starts_tomorrow = start_w == jday_tomorrow return end if @@ -1805,9 +1800,6 @@ subroutine get_swindow(jday, rx_starts, rx_ends, param_start, param_end, w, star next_swindow_start = start_w end if - ! Does the NEXT sowing window start tomorrow? - sowing_window_starts_tomorrow = next_swindow_start == jday_tomorrow - end subroutine get_swindow @@ -1872,7 +1864,6 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & integer mcsec ! seconds of day (0, ..., seconds/day) integer sowing_window_startdate ! date (day of year) of first day of sowing window integer sowing_window_enddate ! date (day of year) of last day of sowing window - logical sowing_window_starts_tomorrow ! Is tomorrow the start of a sowing window? real(r8) harvest_reason real(r8) dayspyr ! days per year in this year real(r8) avg_dayspyr ! average number of days per year @@ -2016,7 +2007,7 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & end if ! Get dates of current or next sowing window. - call get_swindow(jday, crop_inst%rx_swindow_starts_thisyr_patch(p,:), crop_inst%rx_swindow_ends_thisyr_patch(p,:), minplantjday(ivt(p),h), maxplantjday(ivt(p),h), w, sowing_window_startdate, sowing_window_enddate, sowing_window_starts_tomorrow) + call get_swindow(jday, crop_inst%rx_swindow_starts_thisyr_patch(p,:), crop_inst%rx_swindow_ends_thisyr_patch(p,:), minplantjday(ivt(p),h), maxplantjday(ivt(p),h), w, sowing_window_startdate, sowing_window_enddate) ! We only want to plant on a specific day if the prescribed sowing window starts AND ends on the same day. Also make sure we haven't planted yet today. ! TODO: Change last condition to `maxval(crop_inst%sdates_perharv_patch(p,:)) /= jday` diff --git a/src/biogeochem/test/CNPhenology_test/test_CNPhenology.pf b/src/biogeochem/test/CNPhenology_test/test_CNPhenology.pf index ebab50ad13..30fba3505b 100644 --- a/src/biogeochem/test/CNPhenology_test/test_CNPhenology.pf +++ b/src/biogeochem/test/CNPhenology_test/test_CNPhenology.pf @@ -302,43 +302,36 @@ contains integer, parameter :: param_start = 200 integer, parameter :: param_end = 250 integer :: w, start_w, end_w - logical :: sowing_window_starts_tomorrow - call get_swindow(1, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + call get_swindow(1, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w) @assertEqual(1, w) @assertEqual(1, start_w) @assertEqual(45, end_w) - @assertFalse(sowing_window_starts_tomorrow) - call get_swindow(45, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + call get_swindow(45, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w) @assertEqual(1, w) @assertEqual(1, start_w) @assertEqual(45, end_w) - @assertFalse(sowing_window_starts_tomorrow) - call get_swindow(149, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + call get_swindow(149, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w) @assertEqual(2, w) @assertEqual(150, start_w) @assertEqual(180, end_w) - @assertTrue(sowing_window_starts_tomorrow) - call get_swindow(175, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + call get_swindow(175, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w) @assertEqual(2, w) @assertEqual(150, start_w) @assertEqual(180, end_w) - @assertFalse(sowing_window_starts_tomorrow) - call get_swindow(229, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + call get_swindow(229, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w) @assertEqual(1, w) @assertEqual(1, start_w) @assertEqual(45, end_w) - @assertFalse(sowing_window_starts_tomorrow) - call get_swindow(get_curr_days_per_year(), rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + call get_swindow(get_curr_days_per_year(), rx_starts, rx_ends, param_start, param_end, w, start_w, end_w) @assertEqual(1, w) @assertEqual(1, start_w) @assertEqual(45, end_w) - @assertTrue(sowing_window_starts_tomorrow) end subroutine test_get_swindow_startend @Test @@ -350,43 +343,36 @@ contains integer, parameter :: param_start = 200 integer, parameter :: param_end = 250 integer :: w, start_w, end_w - logical :: sowing_window_starts_tomorrow - call get_swindow(1, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + call get_swindow(1, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w) @assertEqual(1, w) @assertEqual(360, start_w) @assertEqual(45, end_w) - @assertFalse(sowing_window_starts_tomorrow) - call get_swindow(45, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + call get_swindow(45, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w) @assertEqual(1, w) @assertEqual(start_w, 360) @assertEqual(end_w, 45) - @assertFalse(sowing_window_starts_tomorrow) - call get_swindow(149, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + call get_swindow(149, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w) @assertEqual(2, w) @assertEqual(150, start_w) @assertEqual(180, end_w) - @assertTrue(sowing_window_starts_tomorrow) - call get_swindow(175, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + call get_swindow(175, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w) @assertEqual(2, w) @assertEqual(150, start_w) @assertEqual(180, end_w) - @assertFalse(sowing_window_starts_tomorrow) - call get_swindow(229, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + call get_swindow(229, rx_starts, rx_ends, param_start, param_end, w, start_w, end_w) @assertEqual(1, w) @assertEqual(360, start_w) @assertEqual(45, end_w) - @assertFalse(sowing_window_starts_tomorrow) - call get_swindow(get_curr_days_per_year(), rx_starts, rx_ends, param_start, param_end, w, start_w, end_w, sowing_window_starts_tomorrow) + call get_swindow(get_curr_days_per_year(), rx_starts, rx_ends, param_start, param_end, w, start_w, end_w) @assertEqual(1, w) @assertEqual(360, start_w) @assertEqual(45, end_w) - @assertFalse(sowing_window_starts_tomorrow) end subroutine test_get_swindow_endstart end module test_CNPhenology From 0585ef24e0cf0200f00d599f9221b75d947d6b98 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Sat, 30 Sep 2023 12:55:40 -0600 Subject: [PATCH 13/52] Resolved a TODO by thinking about it. --- src/biogeochem/CNPhenologyMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index bc93810d92..3f5515d1df 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -2037,7 +2037,7 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & ! This is outside the croplive check so that the "harvest if planting conditions were met today" conditional works. is_in_sowing_window = is_doy_in_interval(sowing_window_startdate, sowing_window_enddate, jday) if (crop_inst%sown_in_this_window(p) .and. .not. is_in_sowing_window) then - ! TODO: Test whether this is necessary, since it's set to false in last timestep of sowing window at the end of CropPhenology(). + ! Although sown_in_this_window is set to false in last timestep of sowing window at the end of CropPhenology(), this extra check may be necessary if sowing windows change. crop_inst%sown_in_this_window(p) = .false. end if is_end_sowing_window = jday == sowing_window_enddate From 0658cab229140adca9ac912caf51443ed1eb80b5 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Sat, 30 Sep 2023 12:56:38 -0600 Subject: [PATCH 14/52] Sowing window outputs now saved on end date. This ensures ordering is identical to how inputs should be. --- src/biogeochem/CNPhenologyMod.F90 | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index 3f5515d1df..12152de52e 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -2043,8 +2043,7 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & is_end_sowing_window = jday == sowing_window_enddate ! ! Save these diagnostic variables only on the first day of the window to ensure that windows spanning the new year aren't double-counted. - ! TODO: Change this to the LAST day of the window, to ensure the same ordering as required for rx sowing window input files. - if (jday == sowing_window_startdate) then + if (jday == sowing_window_enddate) then crop_inst%swindow_starts_thisyr_patch(p,w) = sowing_window_startdate crop_inst%swindow_ends_thisyr_patch (p,w) = sowing_window_enddate end if From ba014b05b860fb4b1ebf9f5f1b6383ecb1e30d3d Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Sat, 30 Sep 2023 13:04:09 -0600 Subject: [PATCH 15/52] Fail if a sowing window start date is prescribed without an end date (or vice versa). --- src/cpl/share_esmf/cropcalStreamMod.F90 | 25 ++++++++++++++++++------- 1 file changed, 18 insertions(+), 7 deletions(-) diff --git a/src/cpl/share_esmf/cropcalStreamMod.F90 b/src/cpl/share_esmf/cropcalStreamMod.F90 index 8ad162783c..113627a10f 100644 --- a/src/cpl/share_esmf/cropcalStreamMod.F90 +++ b/src/cpl/share_esmf/cropcalStreamMod.F90 @@ -323,6 +323,11 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) real(r8), pointer :: dataptr2d_cultivar_gdds(:,:) !----------------------------------------------------------------------- + associate( & + starts => crop_inst%rx_swindow_starts_thisyr_patch, & + ends => crop_inst%rx_swindow_ends_thisyr_patch & + ) + SHR_ASSERT_FL( (lbound(g_to_ig,1) <= bounds%begg ), sourcefile, __LINE__) SHR_ASSERT_FL( (ubound(g_to_ig,1) >= bounds%endg ), sourcefile, __LINE__) @@ -371,18 +376,18 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) n = ivt - npcropmin + 1 ! vegetated pft ig = g_to_ig(patch%gridcell(p)) - crop_inst%rx_swindow_starts_thisyr_patch(p,1) = dataptr2d_swindow_start(ig,n) - crop_inst%rx_swindow_ends_thisyr_patch (p,1) = dataptr2d_swindow_end (ig,n) + starts(p,1) = dataptr2d_swindow_start(ig,n) + ends(p,1) = dataptr2d_swindow_end (ig,n) ! Sanity check: Should only read in valid values - if (crop_inst%rx_swindow_starts_thisyr_patch(p,1) > 365) then + if (starts(p,1) > 365) then write(iulog,'(a,i0,a,i0)') 'cropcal_interp(): Crop patch (ivt ',ivt,') has dataptr2d prescribed sowing window start date ',& - crop_inst%rx_swindow_starts_thisyr_patch(p,1) + starts(p,1) call ESMF_Finalize(endflag=ESMF_END_ABORT) end if - if (crop_inst%rx_swindow_ends_thisyr_patch(p,1) > 365) then + if (ends(p,1) > 365) then write(iulog,'(a,i0,a,i0)') 'cropcal_interp(): Crop patch (ivt ',ivt,') has dataptr2d prescribed sowing window end date ',& - crop_inst%rx_swindow_ends_thisyr_patch(p,1) + ends(p,1) call ESMF_Finalize(endflag=ESMF_END_ABORT) end if else @@ -393,7 +398,11 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) ! TODO: Ensure that, if mxsowings > 1, sowing windows are ordered such that ENDS are monotonically increasing. This is necessary because of how get_swindow() works. - ! TODO: Fail if a sowing window start date is prescribed without an end date (or vice versa) + ! Fail if a sowing window start date is prescribed without an end date (or vice versa) + if (any((starts >= 1 .and. ends < 1) .or. (starts < 1 .and. ends >= 1))) then + write(iulog, *) 'Every prescribed sowing window start date must have a corresponding end date.' + call ESMF_Finalize(endflag=ESMF_END_ABORT) + end if end if ! use_cropcal_rx_swindows deallocate(dataptr2d_swindow_start) @@ -461,6 +470,8 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) deallocate(dataptr2d_cultivar_gdds) + end associate + end subroutine cropcal_interp end module cropcalStreamMod From d0f77fc2c283a66709317208137da62dee7822ee Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Sat, 30 Sep 2023 13:19:54 -0600 Subject: [PATCH 16/52] Ensure that, if mxsowings > 1, sowing windows are ordered such that ENDS are monotonically increasing. --- src/cpl/share_esmf/cropcalStreamMod.F90 | 10 +++++++++- 1 file changed, 9 insertions(+), 1 deletion(-) diff --git a/src/cpl/share_esmf/cropcalStreamMod.F90 b/src/cpl/share_esmf/cropcalStreamMod.F90 index 113627a10f..00f439e265 100644 --- a/src/cpl/share_esmf/cropcalStreamMod.F90 +++ b/src/cpl/share_esmf/cropcalStreamMod.F90 @@ -15,6 +15,7 @@ module cropcalStreamMod use clm_varctl , only : iulog use clm_varctl , only : use_cropcal_rx_swindows, use_cropcal_rx_cultivar_gdds, use_cropcal_streams use clm_varpar , only : mxpft + use clm_varpar , only : mxsowings use perf_mod , only : t_startf, t_stopf use spmdMod , only : masterproc, mpicom, iam use pftconMod , only : npcropmin @@ -396,7 +397,14 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) endif end do - ! TODO: Ensure that, if mxsowings > 1, sowing windows are ordered such that ENDS are monotonically increasing. This is necessary because of how get_swindow() works. + ! Ensure that, if mxsowings > 1, sowing windows are ordered such that ENDS are monotonically increasing. This is necessary because of how get_swindow() works. + if (mxsowings > 1) then + if (any(ends(:,2:mxsowings) <= ends(:,1:mxsowings-1) .and. & + ends(:,2:mxsowings) >= 1)) then + write(iulog, *) 'Sowing window inputs must be ordered such that end dates are monotonically increasing.' + call ESMF_Finalize(endflag=ESMF_END_ABORT) + end if + end if ! Fail if a sowing window start date is prescribed without an end date (or vice versa) if (any((starts >= 1 .and. ends < 1) .or. (starts < 1 .and. ends >= 1))) then From 3422d517e4e39999fea772666039b4c5d14e9149 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Sat, 30 Sep 2023 13:58:54 -0600 Subject: [PATCH 17/52] Simplified some harvest logic. --- src/biogeochem/CNPhenologyMod.F90 | 55 +++++++++---------------------- 1 file changed, 16 insertions(+), 39 deletions(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index 12152de52e..18111994ac 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -2010,7 +2010,6 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & call get_swindow(jday, crop_inst%rx_swindow_starts_thisyr_patch(p,:), crop_inst%rx_swindow_ends_thisyr_patch(p,:), minplantjday(ivt(p),h), maxplantjday(ivt(p),h), w, sowing_window_startdate, sowing_window_enddate) ! We only want to plant on a specific day if the prescribed sowing window starts AND ends on the same day. Also make sure we haven't planted yet today. - ! TODO: Change last condition to `maxval(crop_inst%sdates_perharv_patch(p,:)) /= jday` ! TODO: ¿Allow use of NON-prescribed sowing with one-day-long windows? has_rx_sowing_date = sowing_window_startdate == sowing_window_enddate do_plant_prescribed = has_rx_sowing_date .and. & @@ -2277,44 +2276,22 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & ! If generate_crop_gdds and this patch has prescribed sowing inputs ! TODO: Shouldn't this also apply to paramfile sowing windows? else if (generate_crop_gdds .and. crop_inst%rx_swindow_starts_thisyr_patch(p,1) .gt. 0) then - if (sowing_window_startdate >= 0) then - ! Harvest the day before the start of the next sowing window. - do_harvest = do_plant_prescribed_tomorrow - - ! ... unless that will lead to growing season length 365 (or 366, - ! if last year was a leap year). This would result in idop==jday, - ! which would invoke the "manually setting sowing_count and - ! sdates_thisyr" code. This would lead to crops never getting - ! harvested. Instead, always harvest the day before idop. - if ((.not. do_harvest) .and. & - (idop(p) > 1 .and. jday == idop(p) - 1) .or. & - (idop(p) == 1 .and. jday == dayspyr)) then - do_harvest = .true. - if (do_harvest) then - harvest_reason = HARVEST_REASON_IDOPTOMORROW - end if - else if (do_harvest) then - harvest_reason = HARVEST_REASON_SOWTOMORROW - end if - - else - ! If this patch has already had all its plantings for the year, don't harvest - ! until some time next year. - do_harvest = .false. - - ! ... unless first sowing next year happens Jan. 1. - ! WARNING: This implementation assumes that sowing windows don't change over time! - ! In order to avoid this, you'd have to read this year's AND next year's prescribed - ! sowing window start dates. - if (crop_inst%rx_swindow_starts_thisyr_patch(p,1) == 1) then - do_harvest = jday == dayspyr - end if - - if (do_harvest) then - harvest_reason = HARVEST_REASON_SOWTOMORROW - end if - - endif + ! Harvest the day before the next prescribed sowing. + do_harvest = do_plant_prescribed_tomorrow + + ! ... unless that will lead to growing season length 365 (or 366, + ! if last year was a leap year). This would result in idop==jday, + ! which would invoke the "manually setting sowing_count and + ! sdates_thisyr" code. This would lead to crops never getting + ! harvested. Instead, always harvest the day before idop. + if ((.not. do_harvest) .and. & + (idop(p) > 1 .and. jday == idop(p) - 1) .or. & + (idop(p) == 1 .and. jday == dayspyr)) then + do_harvest = .true. + harvest_reason = HARVEST_REASON_IDOPTOMORROW + else if (do_harvest) then + harvest_reason = HARVEST_REASON_SOWTOMORROW + end if else if (did_plant_prescribed_today) then ! Do not harvest on the day this growing season began; From baed0a736cd1e8a900ba3b0d69838f721c7f2599 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Sun, 1 Oct 2023 16:01:57 -0600 Subject: [PATCH 18/52] Removed a TODO. --- src/biogeochem/CNPhenologyMod.F90 | 1 - 1 file changed, 1 deletion(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index 18111994ac..b625416c54 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -2274,7 +2274,6 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & harvest_reason = HARVEST_REASON_SOWTODAY ! If generate_crop_gdds and this patch has prescribed sowing inputs - ! TODO: Shouldn't this also apply to paramfile sowing windows? else if (generate_crop_gdds .and. crop_inst%rx_swindow_starts_thisyr_patch(p,1) .gt. 0) then ! Harvest the day before the next prescribed sowing. do_harvest = do_plant_prescribed_tomorrow From 5ea05d5532aeba08d425e44ec9f3eea83c759d35 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Mon, 2 Oct 2023 09:18:34 -0600 Subject: [PATCH 19/52] Replaced an EOL backslash with &. --- src/biogeochem/CNPhenologyMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index b625416c54..a4b50735d5 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -2015,7 +2015,7 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & do_plant_prescribed = has_rx_sowing_date .and. & sowing_window_startdate == jday .and. & sowing_count (p) < mxsowings - do_plant_prescribed_tomorrow = \ + do_plant_prescribed_tomorrow = & has_rx_sowing_date .and. & sowing_window_startdate == get_doy_tomorrow(jday) .and. & sowing_count (p) < mxsowings From fd7140e342938f6d17ad4d5201382d75c1871c3b Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Fri, 29 Sep 2023 10:59:27 -0600 Subject: [PATCH 20/52] Added a comment. --- src/biogeochem/CNPhenologyMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index a4b50735d5..b3e82445a5 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -2363,7 +2363,7 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & crop_inst%gddaccum_thisyr_patch(p, harvest_count(p)) = crop_inst%gddaccum_patch(p) crop_inst%hui_thisyr_patch(p, harvest_count(p)) = hui(p) crop_inst%sowing_reason_perharv_patch(p, harvest_count(p)) = real(crop_inst%sowing_reason_patch(p), r8) - crop_inst%sowing_reason_patch(p) = -1 + crop_inst%sowing_reason_patch(p) = -1 ! "Reason for most recent sowing of this patch." So in the line above we save, and here we reset. crop_inst%harvest_reason_thisyr_patch(p, harvest_count(p)) = harvest_reason endif From 1cade9901635da260642b02bccd80b8c8a275127 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Mon, 2 Oct 2023 15:23:10 -0600 Subject: [PATCH 21/52] Add function was_sown_in_this_window() to handle edge cases. Includes unit testing. --- src/biogeochem/CNPhenologyMod.F90 | 51 ++++++++++++- .../test/CNPhenology_test/test_CNPhenology.pf | 75 +++++++++++++++++++ 2 files changed, 122 insertions(+), 4 deletions(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index b3e82445a5..923e2ca053 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -59,6 +59,7 @@ module CNPhenologyMod public :: SeasonalDecidOnset ! Logical function to determine is seasonal decidious onset should be triggered public :: SeasonalCriticalDaylength ! Critical day length needed for Seasonal decidious offset public :: get_swindow + public :: was_sown_in_this_window ! !PRIVITE MEMBER FIUNCTIONS: private :: CNPhenologyClimate ! Get climatological everages to figure out triggers for Phenology @@ -1803,6 +1804,51 @@ subroutine get_swindow(jday, rx_starts, rx_ends, param_start, param_end, w, star end subroutine get_swindow + !----------------------------------------------------------------------- + function was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, sown_in_this_window) + ! !DESCRIPTION: + ! Determine whether the crop was sown in the current sowing window. Although sown_in_this_window is set to false in last timestep of sowing window at the end of CropPhenology(), these extra checks may be necessary if sowing windows change. + ! + ! !USES: + use clm_time_manager , only : is_doy_in_interval + ! !ARGUMENTS: + integer, intent(in) :: sowing_window_startdate, sowing_window_enddate, jday, idop + logical, intent(in) :: sown_in_this_window + ! !LOCAL VARIABLES + logical :: is_in_sowing_window, idop_in_sowing_window + ! !RESULT + logical :: was_sown_in_this_window + + was_sown_in_this_window = sown_in_this_window + + ! If not in a sowing window, sown_in_this_window must be false. + is_in_sowing_window = is_doy_in_interval(sowing_window_startdate, sowing_window_enddate, jday) + if (.not. is_in_sowing_window) then + was_sown_in_this_window = .false. + return + end if + + ! If we're in a sowing window but the day of planting isn't in the active sowing window, we must be in a different sowing window. + idop_in_sowing_window = is_doy_in_interval(sowing_window_startdate, sowing_window_enddate, idop) + if (is_in_sowing_window .and. .not. idop_in_sowing_window) then + was_sown_in_this_window = .false. + return + end if + + ! Sometimes we're in an active sowing window, and the patch was sown between the start and end dates of the window, but not *the currently active* window. + if (sowing_window_startdate <= sowing_window_enddate .and. idop > jday) then + was_sown_in_this_window = .false. + else if (sowing_window_startdate >= sowing_window_enddate) then + if (jday <= sowing_window_enddate .and. idop <= sowing_window_enddate .and. idop > jday) then + was_sown_in_this_window = .false. + else if (jday >= sowing_window_startdate .and. (idop > jday .or. idop <= sowing_window_enddate)) then + was_sown_in_this_window = .false. + end if + end if + + end function was_sown_in_this_window + + !----------------------------------------------------------------------- subroutine CropPhenology(num_pcropp, filter_pcropp , & waterdiagnosticbulk_inst, temperature_inst, crop_inst, canopystate_inst, cnveg_state_inst , & @@ -2035,10 +2081,7 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & ! Are we currently in a sowing window? ! This is outside the croplive check so that the "harvest if planting conditions were met today" conditional works. is_in_sowing_window = is_doy_in_interval(sowing_window_startdate, sowing_window_enddate, jday) - if (crop_inst%sown_in_this_window(p) .and. .not. is_in_sowing_window) then - ! Although sown_in_this_window is set to false in last timestep of sowing window at the end of CropPhenology(), this extra check may be necessary if sowing windows change. - crop_inst%sown_in_this_window(p) = .false. - end if + crop_inst%sown_in_this_window(p) = was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop(p), crop_inst%sown_in_this_window(p)) is_end_sowing_window = jday == sowing_window_enddate ! ! Save these diagnostic variables only on the first day of the window to ensure that windows spanning the new year aren't double-counted. diff --git a/src/biogeochem/test/CNPhenology_test/test_CNPhenology.pf b/src/biogeochem/test/CNPhenology_test/test_CNPhenology.pf index 30fba3505b..b4b65b5b1c 100644 --- a/src/biogeochem/test/CNPhenology_test/test_CNPhenology.pf +++ b/src/biogeochem/test/CNPhenology_test/test_CNPhenology.pf @@ -375,4 +375,79 @@ contains @assertEqual(45, end_w) end subroutine test_get_swindow_endstart + @Test + subroutine test_was_sown_in_this_window_startend(this) + use clm_time_manager , only : is_doy_in_interval + class(TestCNPhenology), intent(inout) :: this + integer, parameter :: sowing_window_startdate = 84 + integer, parameter :: sowing_window_enddate = 205 + integer :: jday, idop + + jday = 100 + + idop = 90 + @assertTrue(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .false.)) + + idop = 100 + @assertTrue(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .false.)) + + idop = 101 + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + + idop = 120 + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + + idop = 360 + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + + jday = 300 + + idop = 90 + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + + idop = 360 + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + end subroutine test_was_sown_in_this_window_startend + + @Test + subroutine test_was_sown_in_this_window_endstart(this) + use clm_time_manager , only : is_doy_in_interval + class(TestCNPhenology), intent(inout) :: this + integer, parameter :: sowing_window_startdate = 205 + integer, parameter :: sowing_window_enddate = 84 + integer :: jday, idop + + jday = 300 + + idop = 60 + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + + idop = 205 + @assertTrue(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .false.)) + + idop = 300 + @assertTrue(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .false.)) + + idop = 301 + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + + jday = 70 + + idop = 60 + @assertTrue(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .false.)) + + idop = 75 + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + + idop = 301 + @assertTrue(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .false.)) + + end subroutine test_was_sown_in_this_window_endstart + end module test_CNPhenology From bab784806a13cd75ad97a0b74860f8d5227ea482 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Mon, 9 Oct 2023 18:21:46 -0600 Subject: [PATCH 22/52] Correct a comment. --- src/biogeochem/CNPhenologyMod.F90 | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index 923e2ca053..c030f24da5 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -2346,10 +2346,10 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & ! Original harvest rule do_harvest = hui(p) >= gddmaturity(p) .or. idpp >= mxmat - ! Always harvest the day before the next prescribed sowing window starts, if still alive. - ! WARNING: This implementation assumes that sowing windows don't change over time! + ! Always harvest the day before the next prescribed sowing date, if still alive. + ! WARNING: This implementation assumes that prescribed sowing dates don't change over time! ! In order to avoid this, you'd have to read this year's AND next year's prescribed - ! sowing windows. + ! sowing dates. do_harvest = do_harvest .or. do_plant_prescribed_tomorrow if (hui(p) >= gddmaturity(p)) then From 16698054724868806949c69eaca38350a03874bf Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Mon, 9 Oct 2023 18:25:13 -0600 Subject: [PATCH 23/52] Remove an unused 'use' variable. --- src/biogeochem/CNPhenologyMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index c030f24da5..98aefc54d3 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -2628,7 +2628,7 @@ subroutine PlantCrop(p, leafcn_in, jday, kyr, do_plant_normal, & ! !USES: use clm_varctl , only : use_c13, use_c14 - use clm_varctl , only : use_cropcal_rx_swindows, use_cropcal_rx_cultivar_gdds, use_cropcal_streams + use clm_varctl , only : use_cropcal_rx_cultivar_gdds, use_cropcal_streams use clm_varcon , only : c13ratio, c14ratio use clm_varpar , only : mxsowings use pftconMod , only : ntmp_corn, nswheat, nwwheat, ntmp_soybean From bc77099aef313a005e7b5b17b689e8c161ce6153 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Mon, 9 Oct 2023 18:29:08 -0600 Subject: [PATCH 24/52] Combined two if-blocks. --- src/cpl/share_esmf/cropcalStreamMod.F90 | 10 ++++------ 1 file changed, 4 insertions(+), 6 deletions(-) diff --git a/src/cpl/share_esmf/cropcalStreamMod.F90 b/src/cpl/share_esmf/cropcalStreamMod.F90 index 00f439e265..120310dd49 100644 --- a/src/cpl/share_esmf/cropcalStreamMod.F90 +++ b/src/cpl/share_esmf/cropcalStreamMod.F90 @@ -151,9 +151,9 @@ subroutine cropcal_init(bounds) use_cropcal_rx_cultivar_gdds = stream_fldFileName_cultivar_gdds /= '' use_cropcal_streams = use_cropcal_rx_swindows .or. use_cropcal_rx_cultivar_gdds - ! Initialize the cdeps data type sdat_cropcal_swindow_start - ! NOTE: stream_dtlimit 1.5 didn't work for some reason if (use_cropcal_rx_swindows) then + ! Initialize the cdeps data type sdat_cropcal_swindow_start + ! NOTE: stream_dtlimit 1.5 didn't work for some reason call shr_strdata_init_from_inline(sdat_cropcal_swindow_start, & my_task = iam, & logunit = iulog, & @@ -178,11 +178,9 @@ subroutine cropcal_init(bounds) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) then call ESMF_Finalize(endflag=ESMF_END_ABORT) end if - end if - ! Initialize the cdeps data type sdat_cropcal_swindow_end - ! NOTE: stream_dtlimit 1.5 didn't work for some reason - if (use_cropcal_rx_swindows) then + ! Initialize the cdeps data type sdat_cropcal_swindow_end + ! NOTE: stream_dtlimit 1.5 didn't work for some reason call shr_strdata_init_from_inline(sdat_cropcal_swindow_end, & my_task = iam, & logunit = iulog, & From e161bf5003c96aad5edb0e6aacd18b31600b26de Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Tue, 10 Oct 2023 16:30:07 -0600 Subject: [PATCH 25/52] Add SWINDOW start/end dates to crop testmod. --- cime_config/testdefs/testmods_dirs/clm/crop/user_nl_clm | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/cime_config/testdefs/testmods_dirs/clm/crop/user_nl_clm b/cime_config/testdefs/testmods_dirs/clm/crop/user_nl_clm index dc22a973b2..4679194e9f 100644 --- a/cime_config/testdefs/testmods_dirs/clm/crop/user_nl_clm +++ b/cime_config/testdefs/testmods_dirs/clm/crop/user_nl_clm @@ -10,7 +10,7 @@ hist_fincl2 += 'DYN_COL_SOIL_ADJUSTMENTS_C' ! Annual crop variables on per-sowing/per-harvest axes, per PFT. -hist_fincl3 = 'SDATES', 'SDATES_PERHARV', 'SYEARS_PERHARV', 'HDATES', 'GRAINC_TO_FOOD_PERHARV', 'GRAINC_TO_FOOD_ANN', 'HDATES', 'GDDHARV_PERHARV', 'GDDACCUM_PERHARV', 'HUI_PERHARV', 'SOWING_REASON_PERHARV', 'HARVEST_REASON_PERHARV' +hist_fincl3 = 'SDATES', 'SDATES_PERHARV', 'SYEARS_PERHARV', 'HDATES', 'GRAINC_TO_FOOD_PERHARV', 'GRAINC_TO_FOOD_ANN', 'HDATES', 'GDDHARV_PERHARV', 'GDDACCUM_PERHARV', 'HUI_PERHARV', 'SOWING_REASON_PERHARV', 'HARVEST_REASON_PERHARV', 'SWINDOW_STARTS', 'SWINDOW_ENDS' hist_nhtfrq(3) = 17520 hist_mfilt(3) = 1 hist_type1d_pertape(3) = 'PFTS' From 703de51ceed80464a32f2a9f11705de2569f3a7c Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Tue, 10 Oct 2023 16:31:15 -0600 Subject: [PATCH 26/52] Add sowingWindows testmod. --- .../testmods_dirs/clm/sowingWindows/include_user_mods | 1 + .../testdefs/testmods_dirs/clm/sowingWindows/user_nl_clm | 5 +++++ 2 files changed, 6 insertions(+) create mode 100644 cime_config/testdefs/testmods_dirs/clm/sowingWindows/include_user_mods create mode 100644 cime_config/testdefs/testmods_dirs/clm/sowingWindows/user_nl_clm diff --git a/cime_config/testdefs/testmods_dirs/clm/sowingWindows/include_user_mods b/cime_config/testdefs/testmods_dirs/clm/sowingWindows/include_user_mods new file mode 100644 index 0000000000..fe0e18cf88 --- /dev/null +++ b/cime_config/testdefs/testmods_dirs/clm/sowingWindows/include_user_mods @@ -0,0 +1 @@ +../default diff --git a/cime_config/testdefs/testmods_dirs/clm/sowingWindows/user_nl_clm b/cime_config/testdefs/testmods_dirs/clm/sowingWindows/user_nl_clm new file mode 100644 index 0000000000..61bdf9aa93 --- /dev/null +++ b/cime_config/testdefs/testmods_dirs/clm/sowingWindows/user_nl_clm @@ -0,0 +1,5 @@ +stream_fldFileName_swindow_start = '/glade/u/home/samrabin/inputdata/lnd/clm2/cropdata/calendars/processed/swindow_starts_ggcmi_crop_calendar_phase3_v1.01.2000-2000.20231005_145103.nc' +stream_fldFileName_swindow_end = '/glade/u/home/samrabin/inputdata/lnd/clm2/cropdata/calendars/processed/swindow_ends_ggcmi_crop_calendar_phase3_v1.01.2000-2000.20231005_145103.nc' +stream_meshfile_cropcal = '/glade/u/home/samrabin/inputdata/share/meshes/360x720_120830_ESMFmesh_c20210507_cdf5.nc' +stream_year_first_cropcal = 2000 +stream_year_last_cropcal = 2000 From 30f47a5cbd99d354cac81eadd62ddc7e6194f0a8 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 11 Oct 2023 12:06:57 -0600 Subject: [PATCH 27/52] Remove unused next_swindow_start from get_swindow(). --- src/biogeochem/CNPhenologyMod.F90 | 24 ------------------------ 1 file changed, 24 deletions(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index 98aefc54d3..8d60b16c33 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -1731,8 +1731,6 @@ subroutine get_swindow(jday, rx_starts, rx_ends, param_start, param_end, w, star integer, intent(out) :: start_w, end_w ! Start and end dates of "next" sowing window ! ! !LOCAL VARIABLES - integer :: next_swindow_start - integer :: i, x integer :: jday_tomorrow integer :: mxsowings_in ! Due to unit testing, we can't assume the length of the rx sowing window arrays is mxsowings as set in clm_varpar @@ -1779,28 +1777,6 @@ subroutine get_swindow(jday, rx_starts, rx_ends, param_start, param_end, w, star call endrun(msg="get_swindow(): No sowing window found") end if - ! Get the start date of the NEXT sowing window (not including the sowing window we're currently in, if any) - if (is_doy_in_interval(start_w, end_w, jday)) then - next_swindow_start = -1 - x = w - i = 0 ! For checking infinite loop - do while (next_swindow_start < 1) - ! Check for infinite loop - i = i + 1 - if (i > mxsowings_in + 1) then - call endrun("Infinite loop in get_swindow()") - end if - - x = x + 1 - if (x > mxsowings_in) then - x = 1 - end if - next_swindow_start = rx_starts(x) - end do - else - next_swindow_start = start_w - end if - end subroutine get_swindow From df30d0479afce773d311ad9ad20472d388ebd6d8 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 11 Oct 2023 12:13:07 -0600 Subject: [PATCH 28/52] Change check for invalid dates in get_swindow(). --- src/biogeochem/CNPhenologyMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index 8d60b16c33..c9e6457be3 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -1761,7 +1761,7 @@ subroutine get_swindow(jday, rx_starts, rx_ends, param_start, param_end, w, star ! Otherwise, use the first prescribed sowing window we find whose end is >= today. do w = 1, mxsowings_in ! If nothing prescribed at this w, stop looking and exit loop. - if (min(rx_starts(w), rx_ends(w)) < 0) then + if (min(rx_starts(w), rx_ends(w)) < 1) then exit end if From 18c9eeb8c219baece4e8a3a88074afa7d392edb9 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 11 Oct 2023 12:57:54 -0600 Subject: [PATCH 29/52] Improve comments in get_swindow(). --- src/biogeochem/CNPhenologyMod.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index c9e6457be3..24839bc0cf 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -1758,9 +1758,9 @@ subroutine get_swindow(jday, rx_starts, rx_ends, param_start, param_end, w, star return end if - ! Otherwise, use the first prescribed sowing window we find whose end is >= today. + ! Otherwise, use the first prescribed sowing window we find whose end is >= today. This works only if sowing windows that span the new year are located at index w = 1. do w = 1, mxsowings_in - ! If nothing prescribed at this w, stop looking and exit loop. + ! If nothing prescribed at this w, stop looking and exit loop. Will trigger "No sowing window found" error, which we do not move here because it's possible that no start or end date is < 1. if (min(rx_starts(w), rx_ends(w)) < 1) then exit end if From f058160e7fc3592daec30fd347e4cccf2da97e11 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Tue, 17 Oct 2023 13:58:25 -0600 Subject: [PATCH 30/52] Comment fixes and improvements. --- bld/namelist_files/namelist_definition_ctsm.xml | 4 ++-- src/biogeochem/CNPhenologyMod.F90 | 5 +++-- src/biogeochem/test/CNPhenology_test/test_CNPhenology.pf | 2 ++ 3 files changed, 7 insertions(+), 4 deletions(-) diff --git a/bld/namelist_files/namelist_definition_ctsm.xml b/bld/namelist_files/namelist_definition_ctsm.xml index caeab9eb48..3f494e0711 100644 --- a/bld/namelist_files/namelist_definition_ctsm.xml +++ b/bld/namelist_files/namelist_definition_ctsm.xml @@ -1837,12 +1837,12 @@ Simulation year that aligns with stream_year_first_cropcal value -Filename of input stream data for date (day of year) of start of sowing window +Filename of input stream data for date (day of year) of start of sowing window. Cells with the same sowing window start and end date are always planted on that date, regardless of climatic conditions/history. -Filename of input stream data for date (day of year) of end of sowing window +Filename of input stream data for date (day of year) of end of sowing window. Cells with the same sowing window start and end date are always planted on that date, regardless of climatic conditions/history. Date: Tue, 17 Oct 2023 14:02:25 -0600 Subject: [PATCH 31/52] CropType: Don't allocate things twice. --- src/biogeochem/CropType.F90 | 4 ---- 1 file changed, 4 deletions(-) diff --git a/src/biogeochem/CropType.F90 b/src/biogeochem/CropType.F90 index 6d101df456..29c4717ab3 100644 --- a/src/biogeochem/CropType.F90 +++ b/src/biogeochem/CropType.F90 @@ -248,10 +248,6 @@ subroutine InitAllocate(this, bounds) allocate(this%sowing_count(begp:endp)) ; this%sowing_count(:) = 0 allocate(this%harvest_count(begp:endp)) ; this%harvest_count(:) = 0 - this%rx_swindow_starts_thisyr_patch(:,:) = -1 - this%rx_swindow_ends_thisyr_patch (:,:) = -1 - this%rx_cultivar_gdds_thisyr_patch(:,:) = spval - end subroutine InitAllocate !----------------------------------------------------------------------- From 90bbb16279d7110f3b50f046f721d41f1f3973b2 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Tue, 17 Oct 2023 14:18:26 -0600 Subject: [PATCH 32/52] Add input validity check to get_doy_tomorrow(). --- src/utils/clm_time_manager.F90 | 6 ++++++ .../test/clm_time_manager_test/test_clm_time_manager.pf | 7 +++++++ 2 files changed, 13 insertions(+) diff --git a/src/utils/clm_time_manager.F90 b/src/utils/clm_time_manager.F90 index 14572935c3..35abbfc6bf 100644 --- a/src/utils/clm_time_manager.F90 +++ b/src/utils/clm_time_manager.F90 @@ -1282,6 +1282,12 @@ function get_doy_tomorrow(doy_today) result(doy_tomorrow) integer, intent(in) :: doy_today integer :: doy_tomorrow + character(len=*), parameter :: sub = 'clm::get_doy_tomorrow' + + if ( doy_today < 1 .or. doy_today > 367 )then + write(iulog,*) sub, ' = ', doy_today + call shr_sys_abort( sub//': error doy_today out of range' ) + end if if (doy_today == get_curr_days_per_year()) then doy_tomorrow = 1 diff --git a/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf b/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf index 0665e5f6d0..736d745e3a 100644 --- a/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf +++ b/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf @@ -661,6 +661,8 @@ contains @Test subroutine check_get_doy_tomorrow(this) class(TestTimeManager), intent(inout) :: this + character(len=256) :: expected_msg + integer :: doy_tomorrow ! We don't care about the actual date here; we just want to enable get_curr_days_per_year() call unittest_timemgr_setup(dtime=dtime, use_gregorian_calendar=.true.) @@ -669,6 +671,11 @@ contains @assertEqual(get_doy_tomorrow(1), 2) @assertEqual(get_doy_tomorrow(150), 151) @assertEqual(get_doy_tomorrow(get_curr_days_per_year()), 1) + + doy_tomorrow = get_doy_tomorrow(400) + expected_msg = endrun_msg("clm::get_doy_tomorrow: error doy_today out of range" ) + @assertExceptionRaised(expected_msg) + end subroutine check_get_doy_tomorrow end module test_clm_time_manager From 92b54358f2cad8fb777cf942e0f0285a1d88b6e4 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Tue, 17 Oct 2023 14:48:15 -0600 Subject: [PATCH 33/52] Don't use hardcoded 365 in cropcalStreamMod. --- src/cpl/share_esmf/cropcalStreamMod.F90 | 12 ++++++++---- 1 file changed, 8 insertions(+), 4 deletions(-) diff --git a/src/cpl/share_esmf/cropcalStreamMod.F90 b/src/cpl/share_esmf/cropcalStreamMod.F90 index 120310dd49..e14979a708 100644 --- a/src/cpl/share_esmf/cropcalStreamMod.F90 +++ b/src/cpl/share_esmf/cropcalStreamMod.F90 @@ -299,6 +299,7 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) ! !USES: use CropType , only : crop_type use PatchType , only : patch + use clm_time_manager, only : get_curr_days_per_year use dshr_methods_mod , only : dshr_fldbun_getfldptr ! ! !ARGUMENTS: @@ -311,6 +312,7 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) ! !LOCAL VARIABLES: integer :: ivt, p, ip, ig integer :: nc, fp + integer :: dayspyr integer :: n, g integer :: lsize integer :: rc @@ -334,6 +336,8 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) ! Place all data from each type into a temporary 2d array lsize = bounds%endg - bounds%begg + 1 + dayspyr = get_curr_days_per_year() + ! Read prescribed sowing window start dates from input files allocate(dataptr2d_swindow_start(lsize, ncft)) allocate(dataptr2d_swindow_end (lsize, ncft)) @@ -355,8 +359,8 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) do g = 1,lsize ! If read-in value is invalid, allow_unprescribed_planting in CropPhenology() - if (dataptr1d_swindow_start(g) <= 0 .or. dataptr1d_swindow_start(g) > 365 & - .or. dataptr1d_swindow_end(g) <= 0 .or. dataptr1d_swindow_end(g) > 365) then + if (dataptr1d_swindow_start(g) <= 0 .or. dataptr1d_swindow_start(g) > dayspyr & + .or. dataptr1d_swindow_end(g) <= 0 .or. dataptr1d_swindow_end(g) > dayspyr) then dataptr1d_swindow_start(g) = -1 dataptr1d_swindow_end (g) = -1 end if @@ -379,12 +383,12 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) ends(p,1) = dataptr2d_swindow_end (ig,n) ! Sanity check: Should only read in valid values - if (starts(p,1) > 365) then + if (starts(p,1) > dayspyr) then write(iulog,'(a,i0,a,i0)') 'cropcal_interp(): Crop patch (ivt ',ivt,') has dataptr2d prescribed sowing window start date ',& starts(p,1) call ESMF_Finalize(endflag=ESMF_END_ABORT) end if - if (ends(p,1) > 365) then + if (ends(p,1) > dayspyr) then write(iulog,'(a,i0,a,i0)') 'cropcal_interp(): Crop patch (ivt ',ivt,') has dataptr2d prescribed sowing window end date ',& ends(p,1) call ESMF_Finalize(endflag=ESMF_END_ABORT) From 8d1537c69018e797a492ea64712a479f07d1165c Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Tue, 17 Oct 2023 15:27:34 -0600 Subject: [PATCH 34/52] do_plant_prescribed: Check sown_in_this_window, not sowing_count. --- src/biogeochem/CNPhenologyMod.F90 | 17 ++++++++--------- 1 file changed, 8 insertions(+), 9 deletions(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index bca6200607..60931fb13b 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -2032,16 +2032,21 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & ! Get dates of current or next sowing window. call get_swindow(jday, crop_inst%rx_swindow_starts_thisyr_patch(p,:), crop_inst%rx_swindow_ends_thisyr_patch(p,:), minplantjday(ivt(p),h), maxplantjday(ivt(p),h), w, sowing_window_startdate, sowing_window_enddate) + ! Are we currently in a sowing window? + ! This is outside the croplive check so that the "harvest if planting conditions were met today" conditional works. + is_in_sowing_window = is_doy_in_interval(sowing_window_startdate, sowing_window_enddate, jday) + crop_inst%sown_in_this_window(p) = was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop(p), crop_inst%sown_in_this_window(p)) + is_end_sowing_window = jday == sowing_window_enddate + ! We only want to plant on a specific day if the prescribed sowing window starts AND ends on the same day. Also make sure we haven't planted yet today. ! TODO: ¿Allow use of NON-prescribed sowing with one-day-long windows? has_rx_sowing_date = sowing_window_startdate == sowing_window_enddate do_plant_prescribed = has_rx_sowing_date .and. & sowing_window_startdate == jday .and. & - sowing_count (p) < mxsowings + .not. crop_inst%sown_in_this_window(p) do_plant_prescribed_tomorrow = & has_rx_sowing_date .and. & - sowing_window_startdate == get_doy_tomorrow(jday) .and. & - sowing_count (p) < mxsowings + sowing_window_startdate == get_doy_tomorrow(jday) ! BACKWARDS_COMPATIBILITY(wjs/ssr, 2022-02-18) ! When resuming from a run with old code, may need to manually set these. @@ -2055,12 +2060,6 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & crop_inst%sdates_thisyr_patch(p,1) = real(idop(p), r8) end if - ! Are we currently in a sowing window? - ! This is outside the croplive check so that the "harvest if planting conditions were met today" conditional works. - is_in_sowing_window = is_doy_in_interval(sowing_window_startdate, sowing_window_enddate, jday) - crop_inst%sown_in_this_window(p) = was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop(p), crop_inst%sown_in_this_window(p)) - is_end_sowing_window = jday == sowing_window_enddate - ! ! Save these diagnostic variables only on the last day of the window to ensure that windows spanning the new year aren't double-counted. Doing this on the last day ensures that outputs are ordered as inputs should be. if (jday == sowing_window_enddate) then crop_inst%swindow_starts_thisyr_patch(p,w) = sowing_window_startdate From bc148ebfaf8a67be8914131f50fa6eff302029f3 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Tue, 17 Oct 2023 16:10:32 -0600 Subject: [PATCH 35/52] Remove unneeded check in cropcal_interp(). --- src/cpl/share_esmf/cropcalStreamMod.F90 | 14 ++------------ 1 file changed, 2 insertions(+), 12 deletions(-) diff --git a/src/cpl/share_esmf/cropcalStreamMod.F90 b/src/cpl/share_esmf/cropcalStreamMod.F90 index e14979a708..4976ef74d0 100644 --- a/src/cpl/share_esmf/cropcalStreamMod.F90 +++ b/src/cpl/share_esmf/cropcalStreamMod.F90 @@ -340,7 +340,9 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) ! Read prescribed sowing window start dates from input files allocate(dataptr2d_swindow_start(lsize, ncft)) + dataptr2d_swindow_start(:,:) = -1._r8 allocate(dataptr2d_swindow_end (lsize, ncft)) + dataptr2d_swindow_end(:,:) = -1._r8 if (use_cropcal_rx_swindows) then ! Starting with npcropmin will skip generic crops do n = 1, ncft @@ -381,18 +383,6 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) ig = g_to_ig(patch%gridcell(p)) starts(p,1) = dataptr2d_swindow_start(ig,n) ends(p,1) = dataptr2d_swindow_end (ig,n) - - ! Sanity check: Should only read in valid values - if (starts(p,1) > dayspyr) then - write(iulog,'(a,i0,a,i0)') 'cropcal_interp(): Crop patch (ivt ',ivt,') has dataptr2d prescribed sowing window start date ',& - starts(p,1) - call ESMF_Finalize(endflag=ESMF_END_ABORT) - end if - if (ends(p,1) > dayspyr) then - write(iulog,'(a,i0,a,i0)') 'cropcal_interp(): Crop patch (ivt ',ivt,') has dataptr2d prescribed sowing window end date ',& - ends(p,1) - call ESMF_Finalize(endflag=ESMF_END_ABORT) - end if else write(iulog,'(a,i0)') 'cropcal_interp(), prescribed sowing windows: Crop patch has ivt ',ivt call ESMF_Finalize(endflag=ESMF_END_ABORT) From 96a47bd262cc773e3490607d0e561564043fd03f Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Tue, 17 Oct 2023 16:40:10 -0600 Subject: [PATCH 36/52] Remove some commented-out lines from cropcal_interp(). --- src/cpl/share_esmf/cropcalStreamMod.F90 | 4 ---- 1 file changed, 4 deletions(-) diff --git a/src/cpl/share_esmf/cropcalStreamMod.F90 b/src/cpl/share_esmf/cropcalStreamMod.F90 index 4976ef74d0..1ca4d331c6 100644 --- a/src/cpl/share_esmf/cropcalStreamMod.F90 +++ b/src/cpl/share_esmf/cropcalStreamMod.F90 @@ -436,10 +436,6 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) do fp = 1, num_pcropp p = filter_pcropp(fp) -! if (.not. patch%active(p)) then -! continue -! end if - ivt = patch%itype(p) ! Will skip generic crops if (ivt >= npcropmin) then From 6081633aa85623d38296e66ee7469709a90a2b62 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 18 Oct 2023 09:11:40 -0600 Subject: [PATCH 37/52] Updated paths in sowingWindows/user_nl_clm. --- .../testdefs/testmods_dirs/clm/sowingWindows/user_nl_clm | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/cime_config/testdefs/testmods_dirs/clm/sowingWindows/user_nl_clm b/cime_config/testdefs/testmods_dirs/clm/sowingWindows/user_nl_clm index 61bdf9aa93..d3d922f721 100644 --- a/cime_config/testdefs/testmods_dirs/clm/sowingWindows/user_nl_clm +++ b/cime_config/testdefs/testmods_dirs/clm/sowingWindows/user_nl_clm @@ -1,5 +1,5 @@ -stream_fldFileName_swindow_start = '/glade/u/home/samrabin/inputdata/lnd/clm2/cropdata/calendars/processed/swindow_starts_ggcmi_crop_calendar_phase3_v1.01.2000-2000.20231005_145103.nc' -stream_fldFileName_swindow_end = '/glade/u/home/samrabin/inputdata/lnd/clm2/cropdata/calendars/processed/swindow_ends_ggcmi_crop_calendar_phase3_v1.01.2000-2000.20231005_145103.nc' -stream_meshfile_cropcal = '/glade/u/home/samrabin/inputdata/share/meshes/360x720_120830_ESMFmesh_c20210507_cdf5.nc' +stream_fldFileName_swindow_start = '/glade/p/cesmdata/cseg/inputdata/lnd/clm2/cropdata/calendars/processed/swindow_starts_ggcmi_crop_calendar_phase3_v1.01.2000-2000.20231005_145103.nc' +stream_fldFileName_swindow_end = '/glade/p/cesmdata/cseg/inputdata/lnd/clm2/cropdata/calendars/processed/swindow_ends_ggcmi_crop_calendar_phase3_v1.01.2000-2000.20231005_145103.nc' +stream_meshfile_cropcal = '/glade/p/cesmdata/cseg/inputdata/share/meshes/360x720_120830_ESMFmesh_c20210507_cdf5.nc' stream_year_first_cropcal = 2000 stream_year_last_cropcal = 2000 From 16b73f15caabae283499a4a54a6cb5ca16bf13a5 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 18 Oct 2023 10:12:08 -0600 Subject: [PATCH 38/52] Require allow_invalid_swindow_inputs to avoid error on invalid sowing window read. --- .../namelist_definition_ctsm.xml | 5 ++++ src/cpl/share_esmf/cropcalStreamMod.F90 | 23 +++++++++++++++---- 2 files changed, 23 insertions(+), 5 deletions(-) diff --git a/bld/namelist_files/namelist_definition_ctsm.xml b/bld/namelist_files/namelist_definition_ctsm.xml index 3f494e0711..5c0cc0fa12 100644 --- a/bld/namelist_files/namelist_definition_ctsm.xml +++ b/bld/namelist_files/namelist_definition_ctsm.xml @@ -1835,6 +1835,11 @@ Last year to loop over for crop calendar data Simulation year that aligns with stream_year_first_cropcal value + +By default, a value in stream_fldFileName_swindow_start or _end outside the range [1, 365] (or 366 in leap years) will cause the run to fail. Set this to .true. to instead fall back on the paramfile sowing windows. + + Filename of input stream data for date (day of year) of start of sowing window. Cells with the same sowing window start and end date are always planted on that date, regardless of climatic conditions/history. diff --git a/src/cpl/share_esmf/cropcalStreamMod.F90 b/src/cpl/share_esmf/cropcalStreamMod.F90 index 1ca4d331c6..aa523bea82 100644 --- a/src/cpl/share_esmf/cropcalStreamMod.F90 +++ b/src/cpl/share_esmf/cropcalStreamMod.F90 @@ -38,6 +38,7 @@ module cropcalStreamMod character(len=CS), allocatable :: stream_varnames_sdate(:) ! used for both start and end dates character(len=CS), allocatable :: stream_varnames_cultivar_gdds(:) integer :: ncft ! Number of crop functional types (excl. generic crops) + logical :: allow_invalid_swindow_inputs ! Fall back on paramfile sowing windows in cases of invalid values in stream_fldFileName_swindow_start and _end? character(len=CL) :: stream_fldFileName_swindow_start ! sowing window start stream filename to read character(len=CL) :: stream_fldFileName_swindow_end ! sowing window end stream filename to read character(len=CL) :: stream_fldFileName_cultivar_gdds ! cultivar growing degree-days stream filename to read @@ -84,6 +85,7 @@ subroutine cropcal_init(bounds) stream_year_first_cropcal, & stream_year_last_cropcal, & model_year_align_cropcal, & + allow_invalid_swindow_inputs, & stream_fldFileName_swindow_start, & stream_fldFileName_swindow_end, & stream_fldFileName_cultivar_gdds, & @@ -93,6 +95,7 @@ subroutine cropcal_init(bounds) stream_year_first_cropcal = 1 ! first year in stream to use stream_year_last_cropcal = 1 ! last year in stream to use model_year_align_cropcal = 1 ! align stream_year_first_cropcal with this model year + allow_invalid_swindow_inputs = .false. stream_meshfile_cropcal = '' stream_fldFileName_swindow_start = '' stream_fldFileName_swindow_end = '' @@ -124,6 +127,7 @@ subroutine cropcal_init(bounds) call shr_mpi_bcast(stream_year_first_cropcal , mpicom) call shr_mpi_bcast(stream_year_last_cropcal , mpicom) call shr_mpi_bcast(model_year_align_cropcal , mpicom) + call shr_mpi_bcast(allow_invalid_swindow_inputs, mpicom) call shr_mpi_bcast(stream_fldFileName_swindow_start, mpicom) call shr_mpi_bcast(stream_fldFileName_swindow_end , mpicom) call shr_mpi_bcast(stream_fldFileName_cultivar_gdds, mpicom) @@ -135,6 +139,7 @@ subroutine cropcal_init(bounds) write(iulog,'(a,i8)') ' stream_year_first_cropcal = ',stream_year_first_cropcal write(iulog,'(a,i8)') ' stream_year_last_cropcal = ',stream_year_last_cropcal write(iulog,'(a,i8)') ' model_year_align_cropcal = ',model_year_align_cropcal + write(iulog,'(a,i8)') ' allow_invalid_swindow_inputs = ',allow_invalid_swindow_inputs write(iulog,'(a,a)' ) ' stream_fldFileName_swindow_start = ',trim(stream_fldFileName_swindow_start) write(iulog,'(a,a)' ) ' stream_fldFileName_swindow_end = ',trim(stream_fldFileName_swindow_end) write(iulog,'(a,a)' ) ' stream_fldFileName_cultivar_gdds = ',trim(stream_fldFileName_cultivar_gdds) @@ -360,7 +365,7 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) ! So an explicit loop is required here do g = 1,lsize - ! If read-in value is invalid, allow_unprescribed_planting in CropPhenology() + ! If read-in value is invalid, set to -1. Will be handled later in this subroutine. if (dataptr1d_swindow_start(g) <= 0 .or. dataptr1d_swindow_start(g) > dayspyr & .or. dataptr1d_swindow_end(g) <= 0 .or. dataptr1d_swindow_end(g) > dayspyr) then dataptr1d_swindow_start(g) = -1 @@ -398,10 +403,18 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) end if end if - ! Fail if a sowing window start date is prescribed without an end date (or vice versa) - if (any((starts >= 1 .and. ends < 1) .or. (starts < 1 .and. ends >= 1))) then - write(iulog, *) 'Every prescribed sowing window start date must have a corresponding end date.' - call ESMF_Finalize(endflag=ESMF_END_ABORT) + ! Handle invalid sowing window values + if (any(starts < 1 .or. ends < 1)) then + ! Fail if not allowing fallback to paramfile sowing windows + if ((.not. allow_invalid_swindow_inputs) .and. any(all(starts < 1, dim=2))) then + write(iulog, *) 'At least one crop in one gridcell has invalid prescribed sowing window start date(s). To ignore and fall back to paramfile sowing windows, set allow_invalid_swindow_inputs to .true.' + call ESMF_Finalize(endflag=ESMF_END_ABORT) + + ! Fail if a sowing window start date is prescribed without an end date (or vice versa) + else if (any((starts >= 1 .and. ends < 1) .or. (starts < 1 .and. ends >= 1))) then + write(iulog, *) 'Every prescribed sowing window start date must have a corresponding end date.' + call ESMF_Finalize(endflag=ESMF_END_ABORT) + end if end if end if ! use_cropcal_rx_swindows From e71e8e7deec850dc348b05a395226cd74140416b Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 18 Oct 2023 10:15:08 -0600 Subject: [PATCH 39/52] Only error on invalid swindow if wtgcell > 0. --- src/cpl/share_esmf/cropcalStreamMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/cpl/share_esmf/cropcalStreamMod.F90 b/src/cpl/share_esmf/cropcalStreamMod.F90 index aa523bea82..82dd8d4f82 100644 --- a/src/cpl/share_esmf/cropcalStreamMod.F90 +++ b/src/cpl/share_esmf/cropcalStreamMod.F90 @@ -406,7 +406,7 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) ! Handle invalid sowing window values if (any(starts < 1 .or. ends < 1)) then ! Fail if not allowing fallback to paramfile sowing windows - if ((.not. allow_invalid_swindow_inputs) .and. any(all(starts < 1, dim=2))) then + if ((.not. allow_invalid_swindow_inputs) .and. any(all(starts < 1, dim=2) .and. patch%wtgcell > 0._r8)) then write(iulog, *) 'At least one crop in one gridcell has invalid prescribed sowing window start date(s). To ignore and fall back to paramfile sowing windows, set allow_invalid_swindow_inputs to .true.' call ESMF_Finalize(endflag=ESMF_END_ABORT) From 52d665fd21772f271226601befa6895d9724d51c Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 18 Oct 2023 10:59:15 -0600 Subject: [PATCH 40/52] Only error on invalid swindow if patch is prognostic crop. --- src/cpl/share_esmf/cropcalStreamMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/cpl/share_esmf/cropcalStreamMod.F90 b/src/cpl/share_esmf/cropcalStreamMod.F90 index 82dd8d4f82..70afd263ff 100644 --- a/src/cpl/share_esmf/cropcalStreamMod.F90 +++ b/src/cpl/share_esmf/cropcalStreamMod.F90 @@ -406,7 +406,7 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) ! Handle invalid sowing window values if (any(starts < 1 .or. ends < 1)) then ! Fail if not allowing fallback to paramfile sowing windows - if ((.not. allow_invalid_swindow_inputs) .and. any(all(starts < 1, dim=2) .and. patch%wtgcell > 0._r8)) then + if ((.not. allow_invalid_swindow_inputs) .and. any(all(starts < 1, dim=2) .and. patch%wtgcell > 0._r8 .and. patch%itype >= npcropmin)) then write(iulog, *) 'At least one crop in one gridcell has invalid prescribed sowing window start date(s). To ignore and fall back to paramfile sowing windows, set allow_invalid_swindow_inputs to .true.' call ESMF_Finalize(endflag=ESMF_END_ABORT) From ce7e9ecc86779ecca00390ad67d08cdfe129db7a Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 18 Oct 2023 11:00:21 -0600 Subject: [PATCH 41/52] List crop types with invalid swindows. --- src/cpl/share_esmf/cropcalStreamMod.F90 | 11 +++++++++++ 1 file changed, 11 insertions(+) diff --git a/src/cpl/share_esmf/cropcalStreamMod.F90 b/src/cpl/share_esmf/cropcalStreamMod.F90 index 70afd263ff..d07359ce81 100644 --- a/src/cpl/share_esmf/cropcalStreamMod.F90 +++ b/src/cpl/share_esmf/cropcalStreamMod.F90 @@ -305,6 +305,7 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) use CropType , only : crop_type use PatchType , only : patch use clm_time_manager, only : get_curr_days_per_year + use pftconMod , only : pftname use dshr_methods_mod , only : dshr_fldbun_getfldptr ! ! !ARGUMENTS: @@ -408,6 +409,16 @@ subroutine cropcal_interp(bounds, num_pcropp, filter_pcropp, crop_inst) ! Fail if not allowing fallback to paramfile sowing windows if ((.not. allow_invalid_swindow_inputs) .and. any(all(starts < 1, dim=2) .and. patch%wtgcell > 0._r8 .and. patch%itype >= npcropmin)) then write(iulog, *) 'At least one crop in one gridcell has invalid prescribed sowing window start date(s). To ignore and fall back to paramfile sowing windows, set allow_invalid_swindow_inputs to .true.' + write(iulog, *) 'Affected crops:' + do ivt = npcropmin, mxpft + do fp = 1, num_pcropp + p = filter_pcropp(fp) + if (ivt == patch%itype(p) .and. patch%wtgcell(p) > 0._r8 .and. all(starts(p,:) < 1)) then + write(iulog, *) ' ',pftname(ivt),' (',ivt,')' + exit ! Stop looking for patches of this type + end if + end do + end do call ESMF_Finalize(endflag=ESMF_END_ABORT) ! Fail if a sowing window start date is prescribed without an end date (or vice versa) From 5cb70af593d6732d09e89bbd52c947f9275e63a0 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Thu, 26 Oct 2023 12:50:45 -0600 Subject: [PATCH 42/52] Removed a TODO. --- src/biogeochem/CNPhenologyMod.F90 | 1 - 1 file changed, 1 deletion(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index 60931fb13b..60ab4f6dba 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -2039,7 +2039,6 @@ subroutine CropPhenology(num_pcropp, filter_pcropp , & is_end_sowing_window = jday == sowing_window_enddate ! We only want to plant on a specific day if the prescribed sowing window starts AND ends on the same day. Also make sure we haven't planted yet today. - ! TODO: ¿Allow use of NON-prescribed sowing with one-day-long windows? has_rx_sowing_date = sowing_window_startdate == sowing_window_enddate do_plant_prescribed = has_rx_sowing_date .and. & sowing_window_startdate == jday .and. & From d783396e659d90ddb30ec2ff2dc09e2724e891ac Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Thu, 26 Oct 2023 16:48:18 -0600 Subject: [PATCH 43/52] Fix bug in was_sown_in_this_window() for windows w/ start==end. --- src/biogeochem/CNPhenologyMod.F90 | 6 +-- .../test/CNPhenology_test/test_CNPhenology.pf | 49 +++++++++++++++++++ 2 files changed, 52 insertions(+), 3 deletions(-) diff --git a/src/biogeochem/CNPhenologyMod.F90 b/src/biogeochem/CNPhenologyMod.F90 index 60ab4f6dba..5a3f02a1ee 100644 --- a/src/biogeochem/CNPhenologyMod.F90 +++ b/src/biogeochem/CNPhenologyMod.F90 @@ -1812,10 +1812,10 @@ function was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, return end if - ! Sometimes we're in an active sowing window, and the patch was sown between the start and end dates of the window, but not *the currently active* window. - if (sowing_window_startdate <= sowing_window_enddate .and. idop > jday) then + ! Sometimes we're in an active sowing window, and the patch was sown between the start and end dates of the window, but not *the currently active* window. Note that windows with start==end are not checked here; we always trust the input value of sown_in_this_window in such cases. + if (sowing_window_startdate < sowing_window_enddate .and. idop > jday) then was_sown_in_this_window = .false. - else if (sowing_window_startdate >= sowing_window_enddate) then + else if (sowing_window_startdate > sowing_window_enddate) then if (jday <= sowing_window_enddate .and. idop <= sowing_window_enddate .and. idop > jday) then was_sown_in_this_window = .false. else if (jday >= sowing_window_startdate .and. (idop > jday .or. idop <= sowing_window_enddate)) then diff --git a/src/biogeochem/test/CNPhenology_test/test_CNPhenology.pf b/src/biogeochem/test/CNPhenology_test/test_CNPhenology.pf index 9a3701d68c..9e06bc74e2 100644 --- a/src/biogeochem/test/CNPhenology_test/test_CNPhenology.pf +++ b/src/biogeochem/test/CNPhenology_test/test_CNPhenology.pf @@ -452,4 +452,53 @@ contains end subroutine test_was_sown_in_this_window_endstart + @Test + subroutine test_was_sown_in_this_window_sameday(this) + use clm_time_manager , only : is_doy_in_interval + class(TestCNPhenology), intent(inout) :: this + integer, parameter :: sowing_window_startdate = 205 + integer, parameter :: sowing_window_enddate = 205 + integer :: jday, idop + + ! If today == start == end, we trust whatever the current value of sown_in_this_window is. + jday = 205 + idop = 205 + ! If it's false, then even if idop == jday, it that idop value must be left over from planting in a PREVIOUS year. + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .false.)) + ! The ONLY way was_sown_in_this_window() should return true is if today == start == end == idop AND the current value is true. + @assertTrue(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + + ! There is no other situation where was_sown_in_this_window() should return true. + + jday = 300 + idop = 60 + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .false.)) + + idop = 205 + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .false.)) + + idop = 300 + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .false.)) + + idop = 301 + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + + jday = 70 + + idop = 60 + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .false.)) + + idop = 75 + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + + idop = 301 + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .true.)) + @assertFalse(was_sown_in_this_window(sowing_window_startdate, sowing_window_enddate, jday, idop, .false.)) + + end subroutine test_was_sown_in_this_window_sameday + end module test_CNPhenology From 54e5f247031210f5e082bdfb978d692c83483b2a Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Fri, 27 Oct 2023 13:45:36 -0600 Subject: [PATCH 44/52] More precise error check in get_doy_tomorrow(). --- src/utils/clm_time_manager.F90 | 2 +- src/utils/test/clm_time_manager_test/test_clm_time_manager.pf | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/src/utils/clm_time_manager.F90 b/src/utils/clm_time_manager.F90 index 35abbfc6bf..11b05a171a 100644 --- a/src/utils/clm_time_manager.F90 +++ b/src/utils/clm_time_manager.F90 @@ -1284,7 +1284,7 @@ function get_doy_tomorrow(doy_today) result(doy_tomorrow) integer :: doy_tomorrow character(len=*), parameter :: sub = 'clm::get_doy_tomorrow' - if ( doy_today < 1 .or. doy_today > 367 )then + if ( doy_today < 1 .or. doy_today > get_curr_days_per_year() )then write(iulog,*) sub, ' = ', doy_today call shr_sys_abort( sub//': error doy_today out of range' ) end if diff --git a/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf b/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf index 736d745e3a..c28b87fb2c 100644 --- a/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf +++ b/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf @@ -672,7 +672,7 @@ contains @assertEqual(get_doy_tomorrow(150), 151) @assertEqual(get_doy_tomorrow(get_curr_days_per_year()), 1) - doy_tomorrow = get_doy_tomorrow(400) + doy_tomorrow = get_doy_tomorrow(get_curr_days_per_year() + 1) expected_msg = endrun_msg("clm::get_doy_tomorrow: error doy_today out of range" ) @assertExceptionRaised(expected_msg) From 99e0b916651ccb4ffc72540bbfc043d5ff51ce6d Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Fri, 27 Oct 2023 15:43:46 -0600 Subject: [PATCH 45/52] Added entries in ChangeLog and ChangeSum. --- doc/ChangeLog | 94 +++++++++++++++++++++++++++++++++++++++++++++++++++ doc/ChangeSum | 1 + 2 files changed, 95 insertions(+) diff --git a/doc/ChangeLog b/doc/ChangeLog index c9192e4ffa..657f78cca1 100644 --- a/doc/ChangeLog +++ b/doc/ChangeLog @@ -1,4 +1,98 @@ =============================================================== +Tag name: ctsm5.1.dev147 +Originator(s): samrabin (Sam Rabin, UCAR/TSS, samrabin@ucar.edu) +Date: Fri Oct 27 15:32:54 MDT 2023 +One-line Summary: Add sowing window input files + +Purpose and description of changes +---------------------------------- + +Previously, one could run crops with either (a) sowing windows defined by the hemisphere-specific start and end dates on the paramfile or (b) prescribed sowing dates specified by input file stream_fldFileName_sdate. This PR replaces the latter with two new input files, stream_fldFileName_swindow_start and stream_fldFileName_swindow_end. + + +Significant changes to scientifically-supported configurations +-------------------------------------------------------------- + +Does this tag change answers significantly for any of the following physics configurations? +(Details of any changes will be given in the "Answer changes" section below.) + + [Put an [X] in the box for any configuration with significant answer changes.] + +[ ] clm5_1 + +[ ] clm5_0 + +[ ] ctsm5_0-nwp + +[ ] clm4_5 + + +Notes of particular relevance for users +--------------------------------------- + +Changes to CTSM's user interface (e.g., new/renamed XML or namelist variables): +* Replaces input file stream_fldfilename_sdate (prescribed sowing date) with stream_fldFileName_swindow_start (start of sowing window) and stream_fldFileName_swindow_end (end of sowing window). +* Any gridcell with sowing window start == end will experience prescribed sowing, matching previous behavior with stream_fldfilename_sdate. +* Setting new parameter allow_invalid_swindow_inputs to .true. makes it so that gridcell-crops without values in provided sowing window files will fall back to paramfile sowing windows. Otherwise, such cells will cause an error. + + +Testing summary: +---------------- + + [PASS means all tests PASS; OK means tests PASS other than expected fails.] + + build-namelist tests (if CLMBuildNamelist.pm has changed): + + cheyenne - + + regular tests (aux_clm: https://github.com/ESCOMP/CTSM/wiki/System-Testing-Guide#pre-merge-system-testing): + + cheyenne ---- + izumi ------- + + any other testing (give details below): + * RXCROPMATURITY test passes. + +If the tag used for baseline comparisons was NOT the previous tag, note that here: + + +Answer changes +-------------- + +Changes answers relative to baseline: + + [ If a tag changes answers relative to baseline comparison the + following should be filled in (otherwise remove this section). + And always remove these three lines and parts that don't apply. ] + + Summarize any changes to answers, i.e., + - what code configurations: + - what platforms/compilers: + - nature of change (roundoff; larger than roundoff/same climate; new climate): + + If bitwise differences were observed, how did you show they were no worse + than roundoff? Roundoff differences means one or more lines of code change results + only by roundoff level (because order of operation changes for example). Roundoff + changes to state fields usually grow to greater than roundoff as the simulation progresses. + + If this tag changes climate describe the run(s) done to evaluate the new + climate (put details of the simulations in the experiment database) + - casename: + + URL for LMWG diagnostics output used to validate new climate: + + +Other details +------------- +[Remove any lines that don't apply. Remove entire section if nothing applies.] + +List any externals directories updated (cime, rtm, mosart, cism, fates, etc.): + +Pull Requests that document the changes (include PR ids): +* #2193 (https://github.com/ESCOMP/CTSM/pull/2193) + +=============================================================== +=============================================================== Tag name: ctsm5.1.dev146 Originator(s): glemieux (Gregory Lemieux, LBNL, glemieux@lbl.gov) Date: Tue Oct 24 20:13:17 MDT 2023 diff --git a/doc/ChangeSum b/doc/ChangeSum index f530ef9f7a..b6fc5f0214 100644 --- a/doc/ChangeSum +++ b/doc/ChangeSum @@ -1,5 +1,6 @@ Tag Who Date Summary ============================================================================================================================ + ctsm5.1.dev147 samrabin 10/27/2023 Add sowing window input files ctsm5.1.dev146 glemieux 10/24/2023 FATES cross-grid seed dispersal ctsm5.1.dev145 slevis 10/19/2023 SNICAR snow albedo scheme updates ctsm5.1.dev144 samrabin 10/19/2023 Remove a deprecated shr_mpi_bcast call From e212bc9b1750884e0d05b76b8775da1a05342c3b Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Fri, 27 Oct 2023 16:01:19 -0600 Subject: [PATCH 46/52] Fix formatting in a write() call. --- src/cpl/share_esmf/cropcalStreamMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/cpl/share_esmf/cropcalStreamMod.F90 b/src/cpl/share_esmf/cropcalStreamMod.F90 index d07359ce81..46696eeba9 100644 --- a/src/cpl/share_esmf/cropcalStreamMod.F90 +++ b/src/cpl/share_esmf/cropcalStreamMod.F90 @@ -139,7 +139,7 @@ subroutine cropcal_init(bounds) write(iulog,'(a,i8)') ' stream_year_first_cropcal = ',stream_year_first_cropcal write(iulog,'(a,i8)') ' stream_year_last_cropcal = ',stream_year_last_cropcal write(iulog,'(a,i8)') ' model_year_align_cropcal = ',model_year_align_cropcal - write(iulog,'(a,i8)') ' allow_invalid_swindow_inputs = ',allow_invalid_swindow_inputs + write(iulog,'(a,l1)') ' allow_invalid_swindow_inputs = ',allow_invalid_swindow_inputs write(iulog,'(a,a)' ) ' stream_fldFileName_swindow_start = ',trim(stream_fldFileName_swindow_start) write(iulog,'(a,a)' ) ' stream_fldFileName_swindow_end = ',trim(stream_fldFileName_swindow_end) write(iulog,'(a,a)' ) ' stream_fldFileName_cultivar_gdds = ',trim(stream_fldFileName_cultivar_gdds) From 80a01b07c6a1ae57d19c0cc89f6fb6952a1684fa Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Mon, 30 Oct 2023 10:53:12 -0600 Subject: [PATCH 47/52] Centralize defs of secs_in_day & dtime in test_clm_time_manager. --- .../test_clm_time_manager.pf | 111 +++++++++--------- 1 file changed, 55 insertions(+), 56 deletions(-) diff --git a/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf b/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf index c28b87fb2c..4ac3613164 100644 --- a/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf +++ b/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf @@ -15,7 +15,10 @@ module test_clm_time_manager save real(r8), parameter :: tol = 1.e-13_r8 - integer, parameter :: dtime = 1800 + integer, parameter :: secs_in_day_int = 86400 + real(r8), parameter :: secs_in_day_r8 = 86400._r8 + integer, parameter :: dtime_int = 1800 + real(r8), parameter :: dtime_r8 = 1800._r8 @TestCase type, extends(TestCase) :: TestTimeManager @@ -42,11 +45,11 @@ contains class(TestTimeManager), intent(inout) :: this integer :: step_size - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) step_size = get_step_size() - @assertEqual(dtime, step_size) + @assertEqual(dtime_int, step_size) end subroutine getStepSize_returnsCorrectValue @Test @@ -54,7 +57,7 @@ contains class(TestTimeManager), intent(inout) :: this real(r8) :: calday - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) calday = get_calday(101, 0) @@ -66,7 +69,7 @@ contains class(TestTimeManager), intent(inout) :: this real(r8) :: calday - call unittest_timemgr_setup(dtime=dtime, use_gregorian_calendar=.true.) + call unittest_timemgr_setup(dtime=dtime_int, use_gregorian_calendar=.true.) calday = get_calday(41231, 43200) @@ -78,7 +81,7 @@ contains class(TestTimeManager), intent(inout) :: this real(r8) :: calday - call unittest_timemgr_setup(dtime=dtime, use_gregorian_calendar=.true.) + call unittest_timemgr_setup(dtime=dtime_int, use_gregorian_calendar=.true.) calday = get_calday(41231, 43200, reuse_day_365_for_day_366=.true.) @@ -90,7 +93,7 @@ contains class(TestTimeManager), intent(inout) :: this real(r8) :: calday - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) call set_date(yr=2000, mon=1, day=1, tod=0) @@ -104,7 +107,7 @@ contains class(TestTimeManager), intent(inout) :: this real(r8) :: calday - call unittest_timemgr_setup(dtime=dtime, use_gregorian_calendar=.true.) + call unittest_timemgr_setup(dtime=dtime_int, use_gregorian_calendar=.true.) call set_date(yr=2000, mon=12, day=31, tod=43200) @@ -118,7 +121,7 @@ contains class(TestTimeManager), intent(inout) :: this real(r8) :: calday - call unittest_timemgr_setup(dtime=dtime, use_gregorian_calendar=.true.) + call unittest_timemgr_setup(dtime=dtime_int, use_gregorian_calendar=.true.) call set_date(yr=2000, mon=12, day=31, tod=43200) @@ -133,12 +136,12 @@ contains real(r8) :: calday real(r8) :: calday_expected - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) call set_date(yr=2000, mon=1, day=1, tod=0) calday = get_prev_calday() - calday_expected = 366._r8 - dtime/86400._r8 + calday_expected = 366._r8 - dtime_int/secs_in_day_r8 @assertEqual(calday_expected, calday, tolerance=tol) end subroutine getPrevCalday_jan1Time0_returnsCorrectValue @@ -149,12 +152,12 @@ contains real(r8) :: calday real(r8) :: calday_expected - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) call set_date(yr=2000, mon=1, day=2, tod=0) calday = get_prev_calday() - calday_expected = 2._r8 - dtime/86400._r8 + calday_expected = 2._r8 - dtime_int/secs_in_day_r8 @assertEqual(calday_expected, calday, tolerance=tol) end subroutine getPrevCalday_jan2Time0_returnsCorrectValue @@ -164,7 +167,7 @@ contains class(TestTimeManager), intent(inout) :: this real(r8) :: yearfrac - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) call set_date(yr=2, mon=1, day=1, tod=0) @@ -179,7 +182,7 @@ contains real(r8) :: yearfrac real(r8) :: yearfrac_expected - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) call set_date(yr=2, mon=3, day=1, tod=43200) @@ -195,7 +198,7 @@ contains real(r8) :: yearfrac real(r8) :: yearfrac_expected - call unittest_timemgr_setup(dtime=dtime, use_gregorian_calendar=.true.) + call unittest_timemgr_setup(dtime=dtime_int, use_gregorian_calendar=.true.) call set_date(yr=2000, mon=12, day=31, tod=43200) @@ -209,16 +212,15 @@ contains subroutine getPrevYearfrac_atYearBoundary_returnsLargeValue(this) class(TestTimeManager), intent(inout) :: this real(r8) :: yearfrac - integer, parameter :: secs_in_day = 86400 real(r8) :: yearfrac_expected - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) call set_date(yr=2, mon=1, day=1, tod=0) yearfrac = get_prev_yearfrac() - yearfrac_expected = (365._r8 - real(dtime, r8) / real(secs_in_day, r8)) / 365._r8 + yearfrac_expected = (365._r8 - dtime_r8 / secs_in_day_r8) / 365._r8 @assertEqual(yearfrac_expected, yearfrac) end subroutine getPrevYearfrac_atYearBoundary_returnsLargeValue @@ -226,16 +228,15 @@ contains subroutine getPrevYearfrac_inMiddleOfYear_returnsCorrectValue(this) class(TestTimeManager), intent(inout) :: this real(r8) :: yearfrac - integer, parameter :: secs_in_day = 86400 real(r8) :: yearfrac_expected - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) call set_date(yr=2, mon=3, day=1, tod=43200) yearfrac = get_prev_yearfrac() - yearfrac_expected = (59.5_r8 - real(dtime, r8) / real(secs_in_day, r8)) / 365._r8 + yearfrac_expected = (59.5_r8 - dtime_r8 / secs_in_day_r8) / 365._r8 @assertEqual(yearfrac_expected, yearfrac) end subroutine getPrevYearfrac_inMiddleOfYear_returnsCorrectValue @@ -243,16 +244,15 @@ contains subroutine getPrevYearfrac_leapYearInMiddleOfYear_returnsCorrectValue(this) class(TestTimeManager), intent(inout) :: this real(r8) :: yearfrac - integer, parameter :: secs_in_day = 86400 real(r8) :: yearfrac_expected - call unittest_timemgr_setup(dtime=dtime, use_gregorian_calendar = .true.) + call unittest_timemgr_setup(dtime=dtime_int, use_gregorian_calendar = .true.) call set_date(yr=2000, mon=3, day=1, tod=43200) yearfrac = get_prev_yearfrac() - yearfrac_expected = (60.5_r8 - real(dtime, r8) / real(secs_in_day, r8)) / 366._r8 + yearfrac_expected = (60.5_r8 - dtime_r8 / secs_in_day_r8) / 366._r8 @assertEqual(yearfrac_expected, yearfrac) end subroutine getPrevYearfrac_leapYearInMiddleOfYear_returnsCorrectValue @@ -261,10 +261,9 @@ contains ! This ensures that the correct year is used in determining the number of days in the year class(TestTimeManager), intent(inout) :: this real(r8) :: yearfrac - integer, parameter :: secs_in_day = 86400 real(r8) :: yearfrac_expected - call unittest_timemgr_setup(dtime=dtime, use_gregorian_calendar = .true.) + call unittest_timemgr_setup(dtime=dtime_int, use_gregorian_calendar = .true.) call set_date(yr=2000, mon=1, day=1, tod=0) @@ -272,7 +271,7 @@ contains ! In the following, note that we have 365 and not 366, because the prev_yearfrac uses ! year 1999, which is not a leap year: - yearfrac_expected = (365._r8 - real(dtime, r8) / real(secs_in_day, r8)) / 365._r8 + yearfrac_expected = (365._r8 - dtime_r8 / secs_in_day_r8) / 365._r8 @assertEqual(yearfrac_expected, yearfrac) end subroutine getPrevYearfrac_leapYearAtYearBoundary_returnsCorrectValue @@ -281,7 +280,7 @@ contains class(TestTimeManager), intent(inout) :: this integer :: nstep - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) nstep = get_nstep() @@ -294,7 +293,7 @@ contains integer, parameter :: expected_nstep = 3 integer :: nstep - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) call set_nstep(expected_nstep) @@ -308,9 +307,9 @@ contains class(TestTimeManager), intent(inout) :: this logical :: is_beg - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) - call set_date(yr=2, mon=1, day=1, tod=dtime) + call set_date(yr=2, mon=1, day=1, tod=dtime_int) is_beg = is_beg_curr_year() @@ -322,9 +321,9 @@ contains class(TestTimeManager), intent(inout) :: this logical :: is_beg - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) - call set_date(yr=2, mon=1, day=2, tod=dtime) + call set_date(yr=2, mon=1, day=2, tod=dtime_int) is_beg = is_beg_curr_year() @@ -336,7 +335,7 @@ contains class(TestTimeManager), intent(inout) :: this logical :: is_end - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) call set_date(yr=2, mon=1, day=1, tod=0) @@ -350,7 +349,7 @@ contains class(TestTimeManager), intent(inout) :: this logical :: is_end - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) call set_date(yr=2, mon=1, day=2, tod=0) @@ -364,7 +363,7 @@ contains class(TestTimeManager), intent(inout) :: this logical :: is_first - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) is_first = is_first_step() @@ -376,7 +375,7 @@ contains class(TestTimeManager), intent(inout) :: this logical :: is_first - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) call set_nstep(1) @@ -390,7 +389,7 @@ contains class(TestTimeManager), intent(inout) :: this integer :: nstep - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) nstep = 100 call set_nstep(nstep) @@ -419,7 +418,7 @@ contains real(r8) :: londeg integer :: expected - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) ! Check for local noon at Greenich londeg = 0.0_r8 @@ -446,7 +445,7 @@ contains integer :: secs real(r8) :: londeg - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) londeg = 0.0_r8 do while ( londeg <= 360.0_r8 ) @@ -470,16 +469,16 @@ contains real(r8) :: londeg integer :: expected - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) ! Check for local noon at Greenich for 1 time step ahead londeg = 0.0_r8 - secs = 3600*12 + dtime + secs = 3600*12 + dtime_int call set_date(yr=2018, mon=9, day=3, tod=secs) - expected = secs - dtime - @assertEqual( expected, get_local_time( londeg, offset=-dtime ) ) + expected = secs - dtime_int + @assertEqual( expected, get_local_time( londeg, offset=-dtime_int ) ) londeg = 360.0_r8 - @assertEqual( expected, get_local_time( londeg, offset=-dtime ) ) + @assertEqual( expected, get_local_time( londeg, offset=-dtime_int ) ) end subroutine check_local_time_woffset @@ -490,7 +489,7 @@ contains integer :: secs real(r8) :: londeg - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) ! Check for local noon at Greenich will be true from 11 to 1pm londeg = 0.0_r8 @@ -528,7 +527,7 @@ contains integer :: secs logical :: check - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) londeg = 0.0_r8 secs = get_local_time( londeg ) @@ -549,7 +548,7 @@ contains real(r8) :: londeg integer :: secs - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) londeg = -200.0_r8 secs = get_local_time( londeg ) @@ -567,7 +566,7 @@ contains real(r8) :: londeg integer :: secs - call unittest_timemgr_setup(dtime=dtime) + call unittest_timemgr_setup(dtime=dtime_int) londeg = 400.0_r8 secs = get_local_time( londeg ) @@ -627,7 +626,7 @@ contains integer :: start, end - call unittest_timemgr_setup(dtime=dtime, use_gregorian_calendar=.true.) + call unittest_timemgr_setup(dtime=dtime_int, use_gregorian_calendar=.true.) start = 100 ! April 10 end = 300 ! October 27 @@ -641,7 +640,7 @@ contains @assertFalse(is_today_in_doy_interval(start, end)) ! Test first timestep of interval - call set_date(yr=2009, mon=4, day=10, tod=dtime) + call set_date(yr=2009, mon=4, day=10, tod=dtime_int) @assertTrue(is_today_in_doy_interval(start, end)) ! Test well within interval @@ -653,7 +652,7 @@ contains @assertTrue(is_today_in_doy_interval(start, end)) ! Test first timestep after interval - call set_date(yr=2009, mon=10, day=28, tod=dtime) + call set_date(yr=2009, mon=10, day=28, tod=dtime_int) @assertFalse(is_today_in_doy_interval(start, end)) end subroutine check_is_today_in_doy_interval @@ -662,11 +661,11 @@ contains subroutine check_get_doy_tomorrow(this) class(TestTimeManager), intent(inout) :: this character(len=256) :: expected_msg - integer :: doy_tomorrow + integer :: doy_tomorrow ! Dummy needed for exception tests ! We don't care about the actual date here; we just want to enable get_curr_days_per_year() - call unittest_timemgr_setup(dtime=dtime, use_gregorian_calendar=.true.) - call set_date(yr=2009, mon=10, day=28, tod=dtime) + call unittest_timemgr_setup(dtime=dtime_int, use_gregorian_calendar=.true.) + call set_date(yr=2009, mon=10, day=28, tod=dtime_int) @assertEqual(get_doy_tomorrow(1), 2) @assertEqual(get_doy_tomorrow(150), 151) From 31d84925d80bcbb4818f53f4b0f38a24e4b9bc0d Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Mon, 30 Oct 2023 11:34:24 -0600 Subject: [PATCH 48/52] Better handling of days in year. --- src/utils/clm_time_manager.F90 | 10 +++++++--- .../clm_time_manager_test/test_clm_time_manager.pf | 7 +++++-- 2 files changed, 12 insertions(+), 5 deletions(-) diff --git a/src/utils/clm_time_manager.F90 b/src/utils/clm_time_manager.F90 index 11b05a171a..783493d528 100644 --- a/src/utils/clm_time_manager.F90 +++ b/src/utils/clm_time_manager.F90 @@ -1282,14 +1282,18 @@ function get_doy_tomorrow(doy_today) result(doy_tomorrow) integer, intent(in) :: doy_today integer :: doy_tomorrow + integer :: days_in_year character(len=*), parameter :: sub = 'clm::get_doy_tomorrow' - if ( doy_today < 1 .or. doy_today > get_curr_days_per_year() )then - write(iulog,*) sub, ' = ', doy_today + days_in_year = get_curr_days_per_year() + + if ( doy_today < 1 .or. doy_today > days_in_year )then + write(iulog,*) 'doy_today = ', doy_today + write(iulog,*) 'days_in_year = ', days_in_year call shr_sys_abort( sub//': error doy_today out of range' ) end if - if (doy_today == get_curr_days_per_year()) then + if (doy_today == days_in_year) then doy_tomorrow = 1 else doy_tomorrow = doy_today + 1 diff --git a/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf b/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf index 4ac3613164..fa05b904ca 100644 --- a/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf +++ b/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf @@ -662,16 +662,19 @@ contains class(TestTimeManager), intent(inout) :: this character(len=256) :: expected_msg integer :: doy_tomorrow ! Dummy needed for exception tests + integer :: days_in_year ! We don't care about the actual date here; we just want to enable get_curr_days_per_year() call unittest_timemgr_setup(dtime=dtime_int, use_gregorian_calendar=.true.) call set_date(yr=2009, mon=10, day=28, tod=dtime_int) + days_in_year = get_curr_days_per_year() + @assertEqual(get_doy_tomorrow(1), 2) @assertEqual(get_doy_tomorrow(150), 151) - @assertEqual(get_doy_tomorrow(get_curr_days_per_year()), 1) + @assertEqual(get_doy_tomorrow(days_in_year), 1) - doy_tomorrow = get_doy_tomorrow(get_curr_days_per_year() + 1) + doy_tomorrow = get_doy_tomorrow(days_in_year + 1) expected_msg = endrun_msg("clm::get_doy_tomorrow: error doy_today out of range" ) @assertExceptionRaised(expected_msg) From 818f3ce6aff9f13eee9508cccad17a073ae68316 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Mon, 30 Oct 2023 11:40:28 -0600 Subject: [PATCH 49/52] Add leap-year test of get_doy_tomorrow(). --- .../test_clm_time_manager.pf | 27 +++++++++++++++++-- 1 file changed, 25 insertions(+), 2 deletions(-) diff --git a/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf b/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf index fa05b904ca..fafe7c0f09 100644 --- a/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf +++ b/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf @@ -658,7 +658,7 @@ contains end subroutine check_is_today_in_doy_interval @Test - subroutine check_get_doy_tomorrow(this) + subroutine check_get_doy_tomorrow_noleap(this) class(TestTimeManager), intent(inout) :: this character(len=256) :: expected_msg integer :: doy_tomorrow ! Dummy needed for exception tests @@ -678,6 +678,29 @@ contains expected_msg = endrun_msg("clm::get_doy_tomorrow: error doy_today out of range" ) @assertExceptionRaised(expected_msg) - end subroutine check_get_doy_tomorrow + end subroutine check_get_doy_tomorrow_noleap + + @Test + subroutine check_get_doy_tomorrow_leap(this) + class(TestTimeManager), intent(inout) :: this + character(len=256) :: expected_msg + integer :: doy_tomorrow ! Dummy needed for exception tests + integer :: days_in_year + + call unittest_timemgr_setup(dtime=dtime_int, use_gregorian_calendar=.true.) + + ! Ensure that things work even in the last timestep of a leap year... + call set_date(yr=2008, mon=12, day=31, tod=secs_in_day_int - dtime_int) + days_in_year = get_curr_days_per_year() + @assertEqual(get_doy_tomorrow(days_in_year), 1) + doy_tomorrow = get_doy_tomorrow(days_in_year + 1) + @assertExceptionRaised(expected_msg) + + ! ... as well as in the first timestep after a leap year + call set_date(yr=2009, mon=1, day=1, tod=0) + @assertEqual(get_doy_tomorrow(1), 2) + + end subroutine check_get_doy_tomorrow_leap + end module test_clm_time_manager From c803feac091dd3ef3d1138d8207780fe1869892d Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Mon, 30 Oct 2023 11:57:54 -0600 Subject: [PATCH 50/52] get_doy_tomorrow(): Use get_prev_days_per_year() instead of _curr_. --- src/utils/clm_time_manager.F90 | 2 +- .../test/clm_time_manager_test/test_clm_time_manager.pf | 7 ++++--- 2 files changed, 5 insertions(+), 4 deletions(-) diff --git a/src/utils/clm_time_manager.F90 b/src/utils/clm_time_manager.F90 index 783493d528..5f535b7126 100644 --- a/src/utils/clm_time_manager.F90 +++ b/src/utils/clm_time_manager.F90 @@ -1285,7 +1285,7 @@ function get_doy_tomorrow(doy_today) result(doy_tomorrow) integer :: days_in_year character(len=*), parameter :: sub = 'clm::get_doy_tomorrow' - days_in_year = get_curr_days_per_year() + days_in_year = get_prev_days_per_year() if ( doy_today < 1 .or. doy_today > days_in_year )then write(iulog,*) 'doy_today = ', doy_today diff --git a/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf b/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf index fafe7c0f09..14d244f65c 100644 --- a/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf +++ b/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf @@ -664,11 +664,11 @@ contains integer :: doy_tomorrow ! Dummy needed for exception tests integer :: days_in_year - ! We don't care about the actual date here; we just want to enable get_curr_days_per_year() + ! We don't care about the actual date here; we just want to enable get_prev_days_per_year() call unittest_timemgr_setup(dtime=dtime_int, use_gregorian_calendar=.true.) call set_date(yr=2009, mon=10, day=28, tod=dtime_int) - days_in_year = get_curr_days_per_year() + days_in_year = get_prev_days_per_year() @assertEqual(get_doy_tomorrow(1), 2) @assertEqual(get_doy_tomorrow(150), 151) @@ -691,9 +691,10 @@ contains ! Ensure that things work even in the last timestep of a leap year... call set_date(yr=2008, mon=12, day=31, tod=secs_in_day_int - dtime_int) - days_in_year = get_curr_days_per_year() + days_in_year = get_prev_days_per_year() @assertEqual(get_doy_tomorrow(days_in_year), 1) doy_tomorrow = get_doy_tomorrow(days_in_year + 1) + expected_msg = endrun_msg("clm::get_doy_tomorrow: error doy_today out of range") @assertExceptionRaised(expected_msg) ! ... as well as in the first timestep after a leap year From 0bcec665702da6ba44c527b74ec9fac9fc580d18 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Mon, 30 Oct 2023 13:17:23 -0600 Subject: [PATCH 51/52] Explain use of get_prev_days_per_year(). --- src/utils/clm_time_manager.F90 | 1 + src/utils/test/clm_time_manager_test/test_clm_time_manager.pf | 2 ++ 2 files changed, 3 insertions(+) diff --git a/src/utils/clm_time_manager.F90 b/src/utils/clm_time_manager.F90 index 5f535b7126..955d98057a 100644 --- a/src/utils/clm_time_manager.F90 +++ b/src/utils/clm_time_manager.F90 @@ -1285,6 +1285,7 @@ function get_doy_tomorrow(doy_today) result(doy_tomorrow) integer :: days_in_year character(len=*), parameter :: sub = 'clm::get_doy_tomorrow' + ! Use get_prev_days_per_year() instead of get_curr_days_per_year() because the latter, in the last timestep of a year, actually returns the number of days in the NEXT year. days_in_year = get_prev_days_per_year() if ( doy_today < 1 .or. doy_today > days_in_year )then diff --git a/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf b/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf index 14d244f65c..df8a59de4b 100644 --- a/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf +++ b/src/utils/test/clm_time_manager_test/test_clm_time_manager.pf @@ -668,6 +668,7 @@ contains call unittest_timemgr_setup(dtime=dtime_int, use_gregorian_calendar=.true.) call set_date(yr=2009, mon=10, day=28, tod=dtime_int) + ! Use get_prev_days_per_year() instead of get_curr_days_per_year() because the latter, in the last timestep of a year, actually returns the number of days in the NEXT year. days_in_year = get_prev_days_per_year() @assertEqual(get_doy_tomorrow(1), 2) @@ -691,6 +692,7 @@ contains ! Ensure that things work even in the last timestep of a leap year... call set_date(yr=2008, mon=12, day=31, tod=secs_in_day_int - dtime_int) + ! Use get_prev_days_per_year() instead of get_curr_days_per_year() because the latter, in the last timestep of a year, actually returns the number of days in the NEXT year. days_in_year = get_prev_days_per_year() @assertEqual(get_doy_tomorrow(days_in_year), 1) doy_tomorrow = get_doy_tomorrow(days_in_year + 1) From 806f12414bc15ade6b204dfe6254e89e9a7d8d18 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Mon, 30 Oct 2023 16:53:29 -0600 Subject: [PATCH 52/52] Update ChangeLog and ChangeSum. --- doc/ChangeLog | 39 ++++----------------------------------- doc/ChangeSum | 2 +- 2 files changed, 5 insertions(+), 36 deletions(-) diff --git a/doc/ChangeLog b/doc/ChangeLog index 657f78cca1..76286ce770 100644 --- a/doc/ChangeLog +++ b/doc/ChangeLog @@ -1,7 +1,7 @@ =============================================================== Tag name: ctsm5.1.dev147 Originator(s): samrabin (Sam Rabin, UCAR/TSS, samrabin@ucar.edu) -Date: Fri Oct 27 15:32:54 MDT 2023 +Date: Mon Oct 30 16:53:20 MDT 2023 One-line Summary: Add sowing window input files Purpose and description of changes @@ -43,50 +43,19 @@ Testing summary: build-namelist tests (if CLMBuildNamelist.pm has changed): - cheyenne - + cheyenne - PASS regular tests (aux_clm: https://github.com/ESCOMP/CTSM/wiki/System-Testing-Guide#pre-merge-system-testing): - cheyenne ---- - izumi ------- + cheyenne ---- OK (some diffs in field lists) + izumi ------- OK (some diffs in field lists) any other testing (give details below): * RXCROPMATURITY test passes. -If the tag used for baseline comparisons was NOT the previous tag, note that here: - - -Answer changes --------------- - -Changes answers relative to baseline: - - [ If a tag changes answers relative to baseline comparison the - following should be filled in (otherwise remove this section). - And always remove these three lines and parts that don't apply. ] - - Summarize any changes to answers, i.e., - - what code configurations: - - what platforms/compilers: - - nature of change (roundoff; larger than roundoff/same climate; new climate): - - If bitwise differences were observed, how did you show they were no worse - than roundoff? Roundoff differences means one or more lines of code change results - only by roundoff level (because order of operation changes for example). Roundoff - changes to state fields usually grow to greater than roundoff as the simulation progresses. - - If this tag changes climate describe the run(s) done to evaluate the new - climate (put details of the simulations in the experiment database) - - casename: - - URL for LMWG diagnostics output used to validate new climate: - Other details ------------- -[Remove any lines that don't apply. Remove entire section if nothing applies.] - -List any externals directories updated (cime, rtm, mosart, cism, fates, etc.): Pull Requests that document the changes (include PR ids): * #2193 (https://github.com/ESCOMP/CTSM/pull/2193) diff --git a/doc/ChangeSum b/doc/ChangeSum index b6fc5f0214..5e4f816f64 100644 --- a/doc/ChangeSum +++ b/doc/ChangeSum @@ -1,6 +1,6 @@ Tag Who Date Summary ============================================================================================================================ - ctsm5.1.dev147 samrabin 10/27/2023 Add sowing window input files + ctsm5.1.dev147 samrabin 10/30/2023 Add sowing window input files ctsm5.1.dev146 glemieux 10/24/2023 FATES cross-grid seed dispersal ctsm5.1.dev145 slevis 10/19/2023 SNICAR snow albedo scheme updates ctsm5.1.dev144 samrabin 10/19/2023 Remove a deprecated shr_mpi_bcast call