Make disturbance rates happen in parallel

[ckoven]

A thing that's been on our to-do list for a while is to improve the disturbance logic so that it doesn't have to decide which is the dominant type of disturbance in a given timestep. E.g. see here: #266 (comment) and subsequent comments and also point 5 in #450 (comment) for discussion on proposed fix. Marking this as an issue here to hopefully come back to this and try to fix it soon. It seems to me like the simplest solution might be to handle the overlap and possible disturbed area > 100% problems in the subroutine disturbance_rates and then loop over disturbance types instead of picking the largest in subroutine spawn_patches?

[rosiealice]

Just thinking about this a bit, and adding my voice to the idea from the perspective of fire, where we can in principle have low level fires every day that wouldn't result in any mortaliy or disturbance if they were smaller than the physiological mortality, abnd vice-versa... Then if we add low level logging in too, we might end up substantially under estimating disturbance and turnover rates, i think the idea of dealing with overlap in disturbance rates should work fine.

Do we assume that disturbances are non-overlapping until they are > 100% in total, or do we impose some kind of conditional probability thing? (assuming, for example, that logger might avoid forests that are currently on fire...)

[rgknox]

If we assume that the disturbances are uniformly distributed in space, then the fraction of ground Fj experienced jointly by both of two disturbances with area fractions F1 and F2 would simply be their multiplier, right? Fj = F1 * F2, Ftotal = Fj + (F1-Fj) + (F2-Fj) And for three different entities.. I feel like we could figure this out be drawing out a cube with overlapping fractions in the 3 dimensions, but there would be three regions where 2 types overlap, and 1 region where 3 types overlap.

[rosiealice]

I guess we could write some algorithm for that type of logic. Would perhaps not be trivially extensible to n dimensional disturbance types, but there's probably code that does that on the internet somewhere though.

…

On Thu, 2 May 2019, 19:26 Ryan Knox, ***@***.***> wrote: If we assume that the disturbances are uniformly distributed in space, then the fraction of ground Fj experienced jointly by both of two disturbances with area fractions F1 and F2 would simply be their multiplier, right? Fj = F1 * F2, Ftotal = Fj + (F1-Fj) + (F2-Fj) And for three different entities.. I feel like we could figure this out be drawing out a cube with overlapping fractions in the 3 dimensions, but there would be three regions where 2 types overlap, and 1 region where 3 types overlap. — You are receiving this because you commented. Reply to this email directly, view it on GitHub <#490 (comment)>, or mute the thread <https://github.com/notifications/unsubscribe-auth/ADC2YQZQ7O5SCOZGAXJCFSTPTMP27ANCNFSM4G2DOZFA> .

[ckoven]

Just wanted to note here that I have some code that I think will fix this issue, in https://github.com/ckoven/fates/tree/parallel_disturbance2, basically by looping over disturbance types rather than evaluating which is the largest and executing only that. It doesn't try to deal with the joint disturbance problem though. I will try to test this over the next few days, but if anyone wants to take a look to see if it makes sense, and/or try it once I confirm that it is basically working, let me know. E.g., this could be relevant to your harvest code @sshu88.

[sshu88]

@ckoven I will look and try the code after your test.

[ckoven]

@sshu88 ok I ran the code on that branch for a couple years in a global 4x5 run as a test, and weirdly it seems to just work correctly. the key thing is that now the history variable fates_disturbance_rate_potential equals the sum of fates_disturbance_rate_fire and fates_disturbance_rate_treefall (which it didn't before). But that is all without logging, because I don't have the forcing file handy for running with harvest globally. If you could try incorporating into your code and running with the harvest active to see if it makes things behave better, I'd appreciate it. In the meantime I will also shortly submit a PR from that branch to get it in the queue.

[ckoven]

@sshu88 ok I ran the code on that branch for a couple years in a global 4x5 run as a test, and weirdly it seems to just work correctly. the key thing is that now the history variable fates_disturbance_rate_potential equals the sum of fates_disturbance_rate_fire and fates_disturbance_rate_treefall (which it didn't before). But that is all without logging, because I don't have the forcing file handy for running with harvest globally. If you could try incorporating into your code and running with the harvest active to see if it makes things behave better, I'd appreciate it. In the meantime I will also shortly submit a PR from that branch to get it in the queue.

[sshu88]

@ckoven Thank you Charlie! Will merge your work and test with logging.

[sshu88]

@ckoven I tested the code and passed the same global 4x5 test (fates_disturbance_rate_potential = fates_disturbance_rate_logging + fates_disturbance_rate_treefall).
One comment on the code:

       ! if the sum of all disturbance rates is such that they will exceed total patch area on this day, then reduce them all proportionally.
       if ( sum(currentPatch%disturbance_rates(:)) .gt. 1.0_r8 ) then
          tempsum = sum(currentPatch%disturbance_rates(:))
          do i_dist = 1,N_DIST_TYPES
             currentPatch%disturbance_rates(i_dist) = currentPatch%disturbance_rates(i_dist) / tempsum
          end do
       endif

Here the disturbance rates are modified if their sum exceeds one. But the diagnostic output fates_disturbance_rate_potential is calculated before the change. Maybe we can change the calculation sequence?

[ckoven]

Hi @sshu88 thanks for testing, and glad that it works with logging! I guess that, assuming everything goes well with #875 and it gets integrated, I actually might suggest that we delete the diagnostic variable fates_disturbance_rate_potential entirely, since it will become redundant once we are no longer losing part of the disturbance rate.