@@ -26,6 +26,7 @@ module EDCanopyStructureMod
2626 use FatesInterfaceMod , only : hlm_use_planthydro
2727 use FatesInterfaceMod , only : numpft
2828 use FatesPlantHydraulicsMod, only : UpdateH2OVeg,InitHydrCohort, RecruitWaterStorage
29+ use EDTypesMod , only : maxCohortsPerPatch
2930
3031 use PRTGenericMod, only : leaf_organ
3132 use PRTGenericMod, only : all_carbon_elements
@@ -58,6 +59,9 @@ module EDCanopyStructureMod
5859 real (r8 ), parameter :: area_target_precision = 1.0E-11_r8 ! Area conservation must be within this tolerance
5960 real (r8 ), parameter :: area_check_precision = 1.0E-9_r8 ! Area conservation checks must be within this tolerance
6061
62+ real (r8 ), parameter :: similar_height_tol = 1.0E-3_r8 ! I think trees that differ by 1mm
63+ ! can be roughly considered the same right?
64+
6165
6266 ! 10/30/09: Created by Rosie Fisher
6367 ! 2017/2018: Modifications and updates by Ryan Knox
@@ -327,7 +331,9 @@ subroutine DemoteFromLayer(currentSite,currentPatch,i_lyr)
327331 integer , intent (in ) :: i_lyr ! Current canopy layer of interest
328332
329333 ! !LOCAL VARIABLES:
330- type (ed_cohort_type), pointer :: currentCohort,copyc
334+ type (ed_cohort_type), pointer :: currentCohort
335+ type (ed_cohort_type), pointer :: copyc
336+ type (ed_cohort_type), pointer :: nextc ! The next cohort in line
331337 integer :: i_cwd ! Index for CWD pool
332338 real (r8 ) :: cc_loss ! cohort crown area loss in demotion (m2)
333339 real (r8 ) :: leaf_c ! leaf carbon [kg]
@@ -342,10 +348,14 @@ subroutine DemoteFromLayer(currentSite,currentPatch,i_lyr)
342348 real (r8 ) :: remainder_area_hold
343349 real (r8 ) :: sumweights
344350 real (r8 ) :: sumweights_old
351+ real (r8 ) :: sumequal ! for rank-ordered same-size cohorts
352+ ! this tallies their excluded area
345353 real (r8 ) :: arealayer ! the area of the current canopy layer
346354 integer :: exceedance_counter ! when seeking to rebalance demotion exceedance
347355 ! keep a loop counter to check for hangs
348-
356+ logical :: tied_size_with_neighbors
357+ type (ed_cohort_type), pointer :: cohort_tosearch_relative_to, cohort_tocompare_to
358+ real (r8 ) :: total_crownarea_of_tied_cohorts
349359
350360 ! First, determine how much total canopy area we have in this layer
351361
@@ -359,25 +369,99 @@ subroutine DemoteFromLayer(currentSite,currentPatch,i_lyr)
359369 ! In that case, we need to work out which cohorts to demote.
360370 ! We go in order from shortest to tallest for ranked demotion
361371
362- sumweights = 0.0_r8
372+ sumweights = 0.0_r8
363373 currentCohort = > currentPatch% shortest
364374 do while (associated (currentCohort))
365375
366376 call carea_allom(currentCohort% dbh,currentCohort% n, &
367377 currentSite% spread,currentCohort% pft,currentCohort% c_area)
368378
369379 if ( currentCohort% canopy_layer == i_lyr)then
380+
370381 if (ED_val_comp_excln .ge. 0.0_r8 ) then
382+
383+ ! ----------------------------------------------------------
371384 ! normal (stochastic) case. weight cohort demotion by
372385 ! inverse size to a constant power
386+ ! ----------------------------------------------------------
387+
373388 currentCohort% excl_weight = &
374- currentCohort% n/ (currentCohort% dbh** ED_val_comp_excln)
389+ currentCohort% n/ (currentCohort% hite** ED_val_comp_excln)
390+ sumweights = sumweights + currentCohort% excl_weight
391+
375392 else
393+
394+ ! -----------------------------------------------------------
376395 ! Rank ordered deterministic method
377- currentCohort% excl_weight = &
378- max (min (currentCohort% c_area, demote_area - sumweights ), 0._r8 )
396+ ! -----------------------------------------------------------
397+ ! If there are cohorts that have the exact same height (which is possible, really)
398+ ! we don't want to unilaterally promote/demote one before the others.
399+ ! So we <>mote them as a unit
400+ ! now we need to go through and figure out how many equal-size cohorts there are.
401+ ! then we need to go through, add up the collective crown areas of all equal-sized
402+ ! and equal-canopy-layer cohorts,
403+ ! and then demote from each as if they were a single group
404+
405+ total_crownarea_of_tied_cohorts = currentCohort% c_area
406+
407+ tied_size_with_neighbors = .false.
408+ nextc = > currentCohort% taller
409+ do while (associated (nextc))
410+ if ( abs (nextc% hite - currentCohort% hite) < similar_height_tol ) then
411+ if ( nextc% canopy_layer .eq. currentCohort% canopy_layer ) then
412+ tied_size_with_neighbors = .true.
413+ total_crownarea_of_tied_cohorts = &
414+ total_crownarea_of_tied_cohorts + nextc% c_area
415+ end if
416+ else
417+ exit
418+ endif
419+ nextc = > nextc% taller
420+ end do
421+
422+ if ( tied_size_with_neighbors ) then
423+
424+ currentCohort% excl_weight = &
425+ max (0.0_r8 ,min (currentCohort% c_area, &
426+ (currentCohort% c_area/ total_crownarea_of_tied_cohorts) * &
427+ (demote_area - sumweights) ))
428+
429+ sumequal = currentCohort% excl_weight
430+
431+ nextc = > currentCohort% taller
432+ do while (associated (nextc))
433+ if ( abs (nextc% hite - currentCohort% hite) < similar_height_tol ) then
434+ if (nextc% canopy_layer .eq. currentCohort% canopy_layer ) then
435+ ! now we know the total crown area of all equal-sized,
436+ ! equal-canopy-layer cohorts
437+ nextc% excl_weight = &
438+ max (0.0_r8 ,min (nextc% c_area, &
439+ (nextc% c_area/ total_crownarea_of_tied_cohorts) * &
440+ (demote_area - sumweights) ))
441+ sumequal = sumequal + nextc% excl_weight
442+ end if
443+ else
444+ exit
445+ endif
446+ nextc = > nextc% taller
447+ end do
448+
449+ ! Update the current cohort pointer to the last similar cohort
450+ ! Its ok if this is not in the right layer
451+ if (associated (nextc))then
452+ currentCohort = > nextc% shorter
453+ else
454+ currentCohort = > currentPatch% tallest
455+ end if
456+ sumweights = sumweights + sumequal
457+
458+ else
459+ currentCohort% excl_weight = &
460+ max (min (currentCohort% c_area, demote_area - sumweights ), 0._r8 )
461+ sumweights = sumweights + currentCohort% excl_weight
462+ end if
463+
379464 endif
380- sumweights = sumweights + currentCohort% excl_weight
381465 endif
382466 currentCohort = > currentCohort% taller
383467 enddo
@@ -653,12 +737,16 @@ subroutine PromoteIntoLayer(currentSite,currentPatch,i_lyr)
653737 ! !LOCAL VARIABLES:
654738 type (ed_cohort_type), pointer :: currentCohort
655739 type (ed_cohort_type), pointer :: copyc
740+ type (ed_cohort_type), pointer :: nextc ! the next cohort, or used for looping
741+ ! cohorts against the current
656742
657743
658744 real (r8 ) :: promote_area
659745 real (r8 ) :: newarea
660746 real (r8 ) :: sumweights
661747 real (r8 ) :: sumweights_old
748+ real (r8 ) :: sumequal ! for tied cohorts, the sum of weights in
749+ ! their group
662750 integer :: exceedance_counter
663751 real (r8 ) :: remainder_area
664752 real (r8 ) :: remainder_area_hold
@@ -670,6 +758,10 @@ subroutine PromoteIntoLayer(currentSite,currentPatch,i_lyr)
670758 real (r8 ) :: sapw_c ! sapwood carbon [kg]
671759 real (r8 ) :: store_c ! storage carbon [kg]
672760 real (r8 ) :: struct_c ! structure carbon [kg]
761+
762+ logical :: tied_size_with_neighbors
763+ type (ed_cohort_type), pointer :: cohort_tosearch_relative_to, cohort_tocompare_to
764+ real (r8 ) :: total_crownarea_of_tied_cohorts
673765
674766 call CanopyLayerArea(currentPatch,currentSite% spread,i_lyr,arealayer_current)
675767 call CanopyLayerArea(currentPatch,currentSite% spread,i_lyr+1 ,arealayer_below)
@@ -728,14 +820,83 @@ subroutine PromoteIntoLayer(currentSite,currentPatch,i_lyr)
728820 call carea_allom(currentCohort% dbh,currentCohort% n,currentSite% spread, &
729821 currentCohort% pft,currentCohort% c_area)
730822 if (currentCohort% canopy_layer == i_lyr+1 )then ! look at the cohorts in the canopy layer below...
823+
731824 if (ED_val_comp_excln .ge. 0.0_r8 ) then
732825 ! normal (stochastic) case, as above.
733- currentCohort% prom_weight = currentCohort% n* currentCohort% dbh** ED_val_comp_excln
826+ currentCohort% prom_weight = currentCohort% n* currentCohort% hite** ED_val_comp_excln
827+ sumweights = sumweights + currentCohort% prom_weight
734828 else
735- currentCohort% prom_weight = max (min (currentCohort% c_area, &
736- promote_area - sumweights ), 0._r8 )
829+
830+ ! ------------------------------------------------------------------
831+ ! Rank ordered deterministic method
832+ ! If there are cohorts that have the exact same height (which is possible, really)
833+ ! we don't want to unilaterally promote/demote one before the others.
834+ ! So we <>mote them as a unit
835+ ! now we need to go through and figure out how many equal-size cohorts there are.
836+ ! then we need to go through, add up the collective crown areas of all equal-sized
837+ ! and equal-canopy-layer cohorts,
838+ ! and then demote from each as if they were a single group
839+ ! ------------------------------------------------------------------
840+
841+ total_crownarea_of_tied_cohorts = currentCohort% c_area
842+ tied_size_with_neighbors = .false.
843+ nextc = > currentCohort% shorter
844+ do while (associated (nextc))
845+ if ( abs (nextc% hite - currentCohort% hite) < similar_height_tol ) then
846+ if ( nextc% canopy_layer .eq. currentCohort% canopy_layer ) then
847+ tied_size_with_neighbors = .true.
848+ total_crownarea_of_tied_cohorts = &
849+ total_crownarea_of_tied_cohorts + nextc% c_area
850+ end if
851+ else
852+ exit
853+ endif
854+ nextc = > nextc% shorter
855+ end do
856+
857+ if ( tied_size_with_neighbors ) then
858+
859+ currentCohort% prom_weight = &
860+ max (0.0_r8 ,min (currentCohort% c_area, &
861+ (currentCohort% c_area/ total_crownarea_of_tied_cohorts) * &
862+ (promote_area - sumweights) ))
863+ sumequal = currentCohort% prom_weight
864+
865+ nextc = > currentCohort% shorter
866+ do while (associated (nextc))
867+ if ( abs (nextc% hite - currentCohort% hite) < similar_height_tol ) then
868+ if (nextc% canopy_layer .eq. currentCohort% canopy_layer ) then
869+ ! now we know the total crown area of all equal-sized,
870+ ! equal-canopy-layer cohorts
871+ nextc% prom_weight = &
872+ max (0.0_r8 ,min (nextc% c_area, &
873+ (nextc% c_area/ total_crownarea_of_tied_cohorts) * &
874+ (promote_area - sumweights) ))
875+ sumequal = sumequal + nextc% prom_weight
876+ end if
877+ else
878+ exit
879+ endif
880+ nextc = > nextc% shorter
881+ end do
882+
883+ ! Update the current cohort pointer to the last similar cohort
884+ ! Its ok if this is not in the right layer
885+ if (associated (nextc))then
886+ currentCohort = > nextc% taller
887+ else
888+ currentCohort = > currentPatch% shortest
889+ end if
890+ sumweights = sumweights + sumequal
891+
892+ else
893+ currentCohort% prom_weight = &
894+ max (min (currentCohort% c_area, promote_area - sumweights ), 0._r8 )
895+ sumweights = sumweights + currentCohort% prom_weight
896+
897+ end if
898+
737899 endif
738- sumweights = sumweights + currentCohort% prom_weight
739900 endif
740901 currentCohort = > currentCohort% shorter
741902 enddo ! currentCohort
0 commit comments