Skip to content

Commit

Permalink
Merge pull request #1155 from ckoven/new_patch_insertion
Browse files Browse the repository at this point in the history 
New patch insertion method
  • Loading branch information
@glemieux
glemieux authored Oct 28, 2024
2 parents f5bd625 + ec61b06 commit 4ffab32
Show file tree
Hide file tree
Showing 4 changed files with 125 additions and 167 deletions.
216 changes: 82 additions & 134 deletions biogeochem/EDPatchDynamicsMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -3815,156 +3815,104 @@ subroutine InsertPatch(currentSite, newPatch)

! !LOCAL VARIABLES:
type (fates_patch_type), pointer :: currentPatch
integer :: insert_method ! Temporary dev
logical :: found_landuselabel_match
integer, parameter :: unordered_lul_groups= 1
integer, parameter :: primaryland_oldest_group = 2
integer, parameter :: numerical_order_lul_groups = 3
integer, parameter :: age_order_only = 4

! Insert new patch case options:
! Option 1: Group the landuse types together, but the group order doesn't matter
! Option 2: Option 1, but primarylands are forced to be the oldest group
! Option 3: Option 1, but groups are in numerical order according to land use label index integer
! (i.e. primarylands=1, secondarylands=2, ..., croplands=5)
! Option 4: Don't group the patches by land use label. Simply add new patches to the youngest end.

! Hardcode the default insertion method. The options developed during FATES V1 land use are
! currently being held for potential future usage.
insert_method = primaryland_oldest_group

! Start from the youngest patch and work to oldest, regarless of insertion_method
currentPatch => currentSite%youngest_patch
logical :: patch_inserted

! For the three grouped cases, if the land use label of the youngest patch on the site
! is a match to the new patch land use label, simply insert it as the new youngest.
! This is applicable to the non-grouped option 4 method as well.
if (currentPatch%land_use_label .eq. newPatch%land_use_label ) then
newPatch%older => currentPatch
newPatch%younger => null()
currentPatch%younger => newPatch
currentSite%youngest_patch => newPatch
else

! If the current site youngest patch land use label doesn't match the new patch
! land use label then work through the list until you find the matching type.
! Since we've just checked the youngest patch, move to the next patch and
! initialize the match flag to false.
found_landuselabel_match = .false.
currentPatch => currentPatch%older
select case(insert_method)
! The goal here is to have patches ordered in a specific way. That way is to have them
! looped as the following, where LU refers to the land use label, PFT refers to the
! nocomp PFT label, and Y and O refer to continuous patch ages.
!
! LU1 ---- LU2 ---- LU3 -- etc
! / \ / \ / \
! PFT1 --- PFT2 | PFT1 --- PFT2 | PFT1 --- PFT2 -- etc
! / \ / \ / \ / \ / \ / \
! O - Y O - Y O - Y O - Y O - Y O - Y -- etc

! Option 1 - order of land use label groups does not matter
case (unordered_lul_groups)
! I.e. to treat land use as the outermost loop element, then nocomp PFT as next loop element,
! and then age as the innermost loop element. Visualizing the above as a linked list patches:

do while(associated(currentPatch) .and. .not. found_landuselabel_match)
if (currentPatch%land_use_label .eq. newPatch%land_use_label) then
found_landuselabel_match = .true.
else
currentPatch => currentPatch%older
end if
end do
! LU1/PFT1/O <-> LU1/PFT1/Y <-> LU1/PFT2/O <- ... -> LU3/PFT2/O <-> LU3/PFT2/Y

! In the case where we've found a land use label matching the new patch label,
! insert the newPatch will as the youngest patch for that land use type.
if (associated(currentPatch)) then
newPatch%older => currentPatch
newPatch%younger => currentPatch%younger
currentPatch%younger%older => newPatch
currentPatch%younger => newPatch
else
! In the case in which we get to the end of the list and haven't found
! a landuse label match simply add the new patch to the youngest end.
newPatch%older => currentSite%youngest_patch
newPatch%younger => null()
currentSite%youngest_patch%younger => newPatch
currentSite%youngest_patch => newPatch
endif

! Option 2 - primaryland group must be on the oldest end
case (primaryland_oldest_group)
! Mapping this setup onto the existing "older/younger" scheme means that lower number
! land use and pft labels are considered "older". Note that this matches the current
! initialization scheme in which patches are created and linked in increasing pft
! numerical order starting from 1. This also aligns with the current set_patchno scheme
! in which patches are given an indexable number for the API iteration loops.

! The way to accomplsh this most simply is to define a pseudo-age that includes all of the
! above info and sort the patches based on the pseudo-age. i.e. take some number larger
! than any patch will ever reach in actual age. Then take the LU, multiply it by the big
! number squared, add it to the pft number multiplied by the big number, and add to the age.
! And lastly to sort using that instead of the actual age.

do while(associated(currentPatch) .and. .not. found_landuselabel_match)
if (currentPatch%land_use_label .eq. newPatch%land_use_label) then
found_landuselabel_match = .true.
else
currentPatch => currentPatch%older
end if
end do
! If land use is turned off or nocomp is turned off, then this should devolve to the prior
! behavior of just age sorting.

! In the case where we've found a land use label matching the new patch label,
! insert the newPatch will as the youngest patch for that land use type.
if (associated(currentPatch)) then
newPatch%older => currentPatch
newPatch%younger => currentPatch%younger
currentPatch%younger%older => newPatch
currentPatch%younger => newPatch
else
! In the case in which we get to the end of the list and haven't found
! a landuse label match.

! If the new patch is primaryland add it to the oldest end of the list
if (newPatch%land_use_label .eq. primaryland) then
newPatch%older => null()
newPatch%younger => currentSite%oldest_patch
currentSite%oldest_patch%older => newPatch
currentSite%oldest_patch => newPatch
else
! If the new patch land use type is not primaryland and we are at the
! oldest end of the list, add it to the youngest end
newPatch%older => currentSite%youngest_patch
newPatch%younger => null()
currentSite%youngest_patch%younger => newPatch
currentSite%youngest_patch => newPatch
endif
endif
patch_inserted = .false.

if (GetPseudoPatchAge(newPatch) .le. GetPseudoPatchAge(currentSite%youngest_patch)) then

! Option 3 - groups are numerically ordered with primaryland group starting at oldest end.
case (numerical_order_lul_groups)
! insert new patch at the head of the linked list
newPatch%older => currentSite%youngest_patch
newPatch%younger => null()
currentSite%youngest_patch%younger => newPatch
currentSite%youngest_patch => newPatch

! If the youngest patch landuse label number is greater than the new
! patch land use label number, the new patch must be inserted somewhere
! in between oldest and youngest
do while(associated(currentPatch) .and. .not. found_landuselabel_match)
if (currentPatch%land_use_label .eq. newPatch%land_use_label .or. &
currentPatch%land_use_label .lt. newPatch%land_use_label) then
found_landuselabel_match = .true.
else
currentPatch => currentPatch%older
endif
end do
patch_inserted = .true.
else if (GetPseudoPatchAge(newPatch) .ge. GetPseudoPatchAge(currentSite%oldest_patch)) then

! In the case where we've found a landuse label matching the new patch label
! insert the newPatch will as the youngest patch for that land use type.
if (associated(currentPatch)) then
! insert new patch at the end of the linked list
newPatch%younger => currentSite%oldest_patch
newPatch%older => null()
currentSite%oldest_patch%older => newPatch
currentSite%oldest_patch => newPatch

newPatch%older => currentPatch
newPatch%younger => currentPatch%younger
patch_inserted = .true.
else
! new patch has a pseudo-age somewhere within the linked list. find the first patch which
! has a pseudo age older than it, and put it ahead of that patch
currentPatch => currentSite%youngest_patch
do while (associated(currentPatch) .and. ( .not. patch_inserted) )
if (GetPseudoPatchAge(newPatch) .lt. GetPseudoPatchAge(currentPatch)) then
newPatch%older => currentPatch
newPatch%younger => currentPatch%younger
currentPatch%younger%older => newPatch
currentPatch%younger => newPatch

else

! In the case were we get to the end, the new patch
! must be numerically the smallest, so put it at the oldest position
newPatch%older => null()
newPatch%younger => currentSite%oldest_patch
currentSite%oldest_patch%older => newPatch
currentSite%oldest_patch => newPatch
currentPatch%younger => newPatch

patch_inserted = .true.
endif
currentPatch => currentPatch%older
end do
end if

! Option 4 - always add the new patch as the youngest regardless of land use label
case (age_order_only)
! Set the current patch to the youngest patch
newPatch%older => currentSite%youngest_patch
newPatch%younger => null()
currentSite%youngest_patch%younger => newPatch
currentSite%youngest_patch => newPatch
end select
if ( .not. patch_inserted) then
! something has gone wrong. abort.
write(fates_log(),*) 'something has gone wrong in the patch insertion, because no place to put the new patch was found'
call endrun(msg=errMsg(sourcefile, __LINE__))
end if

end subroutine InsertPatch
end subroutine InsertPatch

! =====================================================================================

function GetPseudoPatchAge(CurrentPatch) result(pseudo_age)

! Purpose: we want to sort the patches in a way that takes into account both their
! continuous and categorical variables. Calculate a pseudo age that does this, by taking
! the integer labels, multiplying these by large numbers, and adding to the continuous age.
! Note to ensure that lower integer land use label and pft label numbers are considered
! "younger" (i.e higher index patchno) in the linked list, they are summed and multiplied by
! negative one. The patch age is still added normally to this negative pseudoage calculation
! as a higher age will result in a less negative number correlating with an "older" patch.

type (fates_patch_type), intent(in), pointer :: CurrentPatch
real(r8) :: pseudo_age
real(r8), parameter :: max_actual_age = 1.e4 ! hard to imagine a patch older than 10,000 years
real(r8), parameter :: max_actual_age_squared = 1.e8

pseudo_age = -1.0_r8 * (real(CurrentPatch%land_use_label,r8) * max_actual_age_squared + &
real(CurrentPatch%nocomp_pft_label,r8) * max_actual_age) + CurrentPatch%age

end function GetPseudoPatchAge

! =====================================================================================

Expand Down
10 changes: 8 additions & 2 deletions main/EDInitMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -706,6 +706,12 @@ subroutine init_patches( nsites, sites, bc_in)
do el=1,num_elements
call SiteMassStock(sites(s),el,sites(s)%mass_balance(el)%old_stock, &
biomass_stock,litter_stock,seed_stock)
! Initialize the integrated flux balance diagnostics
! No need to initialize the instantaneous states, those are re-calculated
sites(s)%iflux_balance(el)%iflux_liveveg = &
(biomass_stock + seed_stock)*area_inv
sites(s)%iflux_balance(el)%iflux_litter = litter_stock * area_inv

end do
call set_patchno(sites(s))
enddo
Expand Down Expand Up @@ -972,7 +978,7 @@ subroutine init_patches( nsites, sites, bc_in)
do el=1,num_elements
call SiteMassStock(sites(s),el,sites(s)%mass_balance(el)%old_stock, &
biomass_stock,litter_stock,seed_stock)

! Initialize the integrated flux balance diagnostics
! No need to initialize the instantaneous states, those are re-calculated
sites(s)%iflux_balance(el)%iflux_liveveg = &
Expand All @@ -983,7 +989,7 @@ subroutine init_patches( nsites, sites, bc_in)

call set_patchno(sites(s))

enddo sites_loop !s
enddo sites_loop
end if

! zero all the patch fire variables for the first timestep
Expand Down
5 changes: 4 additions & 1 deletion main/EDMainMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -834,7 +834,10 @@ subroutine ed_update_site( currentSite, bc_in, bc_out, is_restarting )
enddo

! Check to see if the time integrated fluxes match the state
call CheckIntegratedMassPools(currentSite)
! Dont call this if we are restarting, it will double count the flux
if(.not.is_restarting)then
call CheckIntegratedMassPools(currentSite)
end if

! The HLMs need to know about nutrient demand, and/or
! root mass and affinities
Expand Down
Loading

0 comments on commit 4ffab32

Please sign in to comment.