Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Incorrect (?) calculation of per-age-class outputs #1205

Open
samsrabin opened this issue May 21, 2024 · 10 comments · May be fixed by #1252
Open

Incorrect (?) calculation of per-age-class outputs #1205

samsrabin opened this issue May 21, 2024 · 10 comments · May be fixed by #1252

Comments

@samsrabin
Copy link
Contributor

samsrabin commented May 21, 2024

There are (at least) two different ways that per-age-class (_AP) outputs are be calculated in update_history_dyn().

When I look at such a variable, I expect the value in each age class to be the mean across all patches in that age class, weighted by patch area. This is how it's done for LAI (FATES_LAI_AP), requiring two steps. First, in patchloop:

hio_lai_si_age(io_si,cpatch%age_class) = hio_lai_si_age(io_si,cpatch%age_class) &
   + sum(cpatch%tlai_profile(:,:,:)) * cpatch%area

Then you need to divide by the total age class area. Since that's accumulated during patchloop, this happens after that loop:

! divide so-far-just-summed but to-be-averaged patch-age-class
! variables by patch-age-class area to get mean values
do ipa2 = 1, nlevage
   if (hio_area_si_age(io_si, ipa2) .gt. nearzero) then
         hio_lai_si_age(io_si, ipa2) = hio_lai_si_age(io_si, ipa2) / (hio_area_si_age(io_si, ipa2)*AREA)
         ...
   else
      hio_lai_si_age(io_si, ipa2) = 0._r8
      ...
   endif
end do

However, for FATES_FUEL_AMOUNT_AP, it's all done in one step within patchloop:

 ! Fuel sum [kg/m2]
 hio_fire_sum_fuel_si_age(io_si, cpatch%age_class) = hio_fire_sum_fuel_si_age(io_si, cpatch%age_class) +  &
    cpatch%sum_fuel * cpatch%area * AREA_INV

There's never a division by the total age class area, so this doesn't give the average fuel load in each patch age class (kg/m2) as the units of FATES_FUEL_AMOUNT_AP state. Instead, it gives each age class's contribution to the total (not per-m2) fuel load. So if you do an unweighted sum across all age classes, you'll get the total fuel load. And it turns out, FATES_FUEL_AMOUNT also gives the total fuel load, in contrast to its units of kg/m2. :
screenshot_0862

As expected, this isn't the case for LAI:
screenshot_0864

I don't know if this is intended behavior for FATES_FUEL_AMOUNT_AP, but as I mentioned it's certainly not what I would expect as a user of the output. Of outputs with fates_levage dimension (I haven't checked outputs with age class crossed by something else, like PFT), the following outputs (in addition to FATES_FUEL_AMOUNT_AP) appear to follow the same behavior:

  • FATES_ZSTAR_AP
  • FATES_SECONDAREA_ANTHRODIST_AP
  • FATES_SECONDAREA_DIST_AP
  • FATES_BURNFRAC_AP
  • FATES_CANOPYAREA_AP
  • FATES_VEGC_AP
@samsrabin
Copy link
Contributor Author

I'm realizing now that my code doesn't have the fix to FATES_LAI_AP that Ryan made as part of #1174. Will add that, redo my run, and see if the LAI map looks any different.

@samsrabin
Copy link
Contributor Author

Progress on this is blocked until #1207 is resolved.

@samsrabin
Copy link
Contributor Author

Now that I know I can work around #1207, I can proceed on this.

@samsrabin
Copy link
Contributor Author

samsrabin commented Jul 8, 2024

This is still the case as of tag sci.1.76.4_api.35.1.0 (June 11).
screenshot_1098

@samsrabin
Copy link
Contributor Author

samsrabin commented Jul 8, 2024

Actually... I think #1174 didn't fully fix the LAI averaging problem. It looks fine after 12 months, which I think is why I told @rgknox that #1174 fixed things, but it turns out there's only one age class at that point. After 24 months (3 age classes):

screenshot_1099

That second map—the cross-ageclass mean of LAI with weights from FATES_PATCHAREA_AP—should be all zero. It looks similar with weights from FATES_CANOPYAREA_AP.

@samsrabin
Copy link
Contributor Author

It turns out that FATES_FUEL_AMOUNT(_AP) is a bad example of the bad behavior, because it's always zero for the test runs I've been doing. FATES_VEGC(_AP) is a good example, though—the unweighted cross-ageclass sum of FATES_VEGC_AP is equal within 1e-6 to FATES_VEGC:
screenshot_1101

Note that this actually represents a third way (second wrong way) of doing things.

FATES_VEGC:

this%hvars(ih_totvegc_si)%r81d(io_si) =                      &
     this%hvars(ih_totvegc_si)%r81d(io_si)+ ccohort%n *        &
     total_m / m2_per_ha

FATES_VEGC_AP:

hio_biomass_si_age(io_si,cpatch%age_class) = hio_biomass_si_age(io_si,cpatch%age_class) &
     + total_m * ccohort%n * AREA_INV

/ m2_per_ha == * AREA_INV == * 1e-4, so it turns out that these are doing the exact same thing—no weights applied. I think here that means that FATES_VEGC is wrong but FATES_VEGC_AP is right.

@ckoven
Copy link
Contributor

ckoven commented Jul 8, 2024

@samsrabin my hunch here is that this gets to the issue with averaging ratios rather than taking the ratio of averages? So properties that are conserved, like biomass, behave well when averaging over time, in the sense that the sum of the age-resolved variable equals the total-ecosystem variable. This isn't the case for variables like LAI, that have a ratio in them where both the numerator (leaf area) and denominator (patch area) change in time, and so the ratio isn't conserved. So perhaps a solution is to get rid of a patch-age-resolved LAI variable entirely, and replace it with a patch-age-resolved leaf area variable, which a user could use to calculate LAI by dividing it by FATES_PATCHAREA_AP?

@rgknox
Copy link
Contributor

rgknox commented Jul 11, 2024

@samsrabin , sorry if this is a dumb question, but are you weighting you sum of the _AP variables by patch area? It seems to do this comparison, you would need that, looks like it would be FATES_PATCHAREA_AP.

@rgknox
Copy link
Contributor

rgknox commented Jul 11, 2024

Actually... I think #1174 didn't fully fix the LAI averaging problem. It looks fine after 12 months, which I think is why I told @rgknox that #1174 fixed things, but it turns out there's only one age class at that point. After 24 months (3 age classes):

screenshot_1099

That second map—the cross-ageclass mean of LAI with weights from FATES_PATCHAREA_AP—should be all zero. It looks similar with weights from FATES_CANOPYAREA_AP.

oh, nevermind, I see you do use it here

@samsrabin
Copy link
Contributor Author

samsrabin commented Jul 17, 2024

@samsrabin my hunch here is that this gets to the issue with averaging ratios rather than taking the ratio of averages? So properties that are conserved, like biomass, behave well when averaging over time, in the sense that the sum of the age-resolved variable equals the total-ecosystem variable. This isn't the case for variables like LAI, that have a ratio in them where both the numerator (leaf area) and denominator (patch area) change in time, and so the ratio isn't conserved. So perhaps a solution is to get rid of a patch-age-resolved LAI variable entirely, and replace it with a patch-age-resolved leaf area variable, which a user could use to calculate LAI by dividing it by FATES_PATCHAREA_AP?

@ckoven I think you're right, thanks! And your solution sounds smart to me.

@samsrabin samsrabin linked a pull request Sep 16, 2024 that will close this issue
4 tasks
@samsrabin samsrabin linked a pull request Sep 16, 2024 that will close this issue
4 tasks
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
Status: ❕Todo
3 participants