PecanProject / pecan

The Predictive Ecosystem Analyzer (PEcAn) is an integrated ecological bioinformatics toolbox.
www.pecanproject.org
Other
202 stars 234 forks source link

High Rd0 values break ED2 runs #1861

Closed istfer closed 6 years ago

istfer commented 6 years ago

Context

I've had successful ED runs (with reasonable outputs) with a particular drivers/ED2IN/model combination in the last months.

With the latest pull, those runs started to give me stepsize underflow errors (leaf temperature is off-track / wood temperature is off track etc.). But some of those runs did finish (although a really small percent finishes without error), and when I look at the outputs, for example, NEE values are all positive..

(I was testing the new ED2 pull when this happened, so I compiled back the older version and recreated the same errors, so that's not due to recent ED2 changes)

Description

After some troubleshooting, found out that Rd0 values written to the configs.xml (10<Rd0<70) of my unsuccessful runs were orders of magnitude higher than the ones in my successful runs (0.06<Rd0<0.3). And the culprit was this line where we write Rd0 values as leaf_respiration_rate_m2.

Carrying our correspondence with @serbinsh here:

Are there incorrect Rd0 in BETY for your PFT?

None of my PFTs have associated priors on Rd0

As far as I can see if that line isn't in write.configs then you wont get Rd0 values written out to the ED2 xml param files?

Here in write.configs.ed we read defaults from history files:

> history$Rd0
 [1] 0.4375000 0.2718750 0.1812500 0.0906250 0.2653500 0.1645750 0.1645750 0.0658300
 [9] 0.2956126 0.2530929 0.1012372 0.2653500 0.2653500 0.4375000 0.4375000 0.2718750
[17] 0.2265625

and then use the defaults if no other value is provided. That's how I get Rd0 values otherwise.

Did something happen to the Arrhenius scaling function to cause the erroneous Rd0? I guess the fact that ED2 crashes with the large Rd0 means that ED2 is looking for that param. I guess not if Vcmax (Vm0) values are reasonable. Do your Vcmax values come out of write configs with reasonable numbers?

Arrhenius scaling function was modified but not changing its outputs.

My Vcmax values are a bit high, but nothing unreasonable I guess. I am using the prior on Bety

I am also using the generic prior that was already on Bety for leaf_respiration_rate_m2 with this description "Distribution and parameters derived from the NASA Forest Functional Types (FFT) project. Dataset consists of leaves sampled across temperate and boreal". But ED2 doesn't read it anyway, see #1628 .

Problem only starts when it's written as Rd0 which ED does read.

Possible Implementation

I think that particular line is not the problem itself, it is consistent with how Rd0 would have been calculated in ED2 if we didn't provide anything.

I think the prior on leaf_respiration_rate_m2 is not good or it's not the same thing as Rd0 in ED2. Though, as it is being constrained by MA, maybe it's the latter?

What do you think @serbinsh?

istfer commented 6 years ago

There is another leaf respiration rate variable in ED

leaf_resp ! Leaf respiration rate   [µmol/m²/s]

It's being calculated from rd_nocorr

         !------ Correct rd. -----------------------------------------------------------!
         aparms(ib)%leaf_resp = rd_nocorr / (tlow_fun * thigh_fun)
         !------------------------------------------------------------------------------!

which is calculated from rd0:

         !----- Find Rd using the Arrhenius equation, with no correction. -----------------!
         rd_nocorr = greeness * arrhenius(met(ib)%leaf_temp,thispft(ib)%rd0                &
                                         ,thispft(ib)%rd_hor)
         !---------------------------------------------------------------------------------!

Not sure if this is the right track, but these make me think that leaf_respiration_rate_m2 is not the same thing as Rd0

serbinsh commented 6 years ago

@istfer are you using the prior that is unif(0,100)?? That wont work as it has a median of 50!! If you don't have data to update that prior then the values will be far too high. Rd0 or leaf_respiration_rate is generally less than 3 umols/m2/s. I would suggest using another prior.

That said, I am still not certain I understand the issue. It seems that the values of Rd0 or leaf_respiration_rate in ed2 history are correct adn ED2's XML tag is for leaf respiration, is it not? That is, when it reads parameters out of the XML at execution it is looking for a tag unless that has changed.

Here are all the leaf respiration related tags I see in my history files for an arbitrary PFT

      <dark_respiration_factor>
                 0.0281375311
      </dark_respiration_factor>
      <Rd_low_temp>
                 4.7136998177
      </Rd_low_temp>
      <Rd_high_temp>
                45.0000000000
      </Rd_high_temp>
      <Rd_decay_e>
                 0.4000000060
      </Rd_decay_e>
      <Rd_hor>
              3000.0000000000
      </Rd_hor>
      <Rd_q10>
                 2.4000000954
      </Rd_q10>
      <Rd0>
                 1.2661544085
      </Rd0>
istfer commented 6 years ago

Hmm, yes I was using that prior, and using random effects = TRUE in MA, so it wasn't updated with data because there wasn't enough data.

The issue was that, I was reading good Rd0 values from the history file, but then bad leaf_respiration_rate_m2 values (due to bad prior, with no data to update it) were overwriting it after you added the line. So values in my configs.xml become like this:

   ...
  <Rd_hor>3000</Rd_hor>
  <Rd_q10>2.4000000954</Rd_q10>
  <Rd0>47.3481320964846</Rd0>
  <D0>0.0160000008</D0> 
   ...

Thanks @serbinsh, which prior are you using for leaf_respiration_rate_m2?

Then, should we get rid of this ridiculously wide prior and associate ED PFTs with a better one? @mdietze?

istfer commented 6 years ago

I now changed the prior on leaf_respiration_rate_m2 for default temperate.Early_Hardwood from ID:189 to ID:197 and closing this issue.