AU-BCE-EE / ABM

R model for anaerobic microbial degradation of organic matter with multiple microbial groups
1 stars 2 forks source link

Add Instantanuous rates in output #55

Open fdalby opened 1 month ago

fdalby commented 1 month ago

A problem was found when running simulations with and without the times argument. The CH4_emis_rate output was different at the same times. This occurs because output rates are currently estimated from the output (where rows have already been removed to match "times").

We could add an argument to decide on instantanuous or average rates in the abm output

fdalby commented 1 month ago

@sashahafner this commit fixes the instant rates issue. https://github.com/AU-BCE-EE/ABM/commit/67f2c258c6e548d8f0ed96c0abf485d2c7d0fca7

But we still need to fix some of the unit conversions. e.g. CH4_emis_cum_COD, where COD_load_cum is estimated from COD_load_rate and diff(dat$time)

sashahafner commented 1 month ago

Nice! I can work on the COD bit now. But I do not think we need it to go into rates(). COD loading rate is constant in an interval, so these should be accurate:

  dat$COD_load_rate <- dat$COD_conc_fresh * dat$slurry_prod_rate
  dat$COD_load_cum <- cumsum(dat$COD_load_rate * c(0, diff(dat$time))) + dat$COD_conc[1] * dat$slurry_mass[1]

I want to check that the c(0, is at the correct end (vs. c(diff(...), 0). Also earlier I have thought about whether we should have an option to exclude initial COD. Ultimately we might want to group all these options together.

abm() does not allow variable COD (VSd, etc.) concentration in feed over time (at least the main version doesn't), correct?

sashahafner commented 1 month ago

Oh, maybe I am wrong, if we can have output intervals that include multiple loading rates, because of the use of times. Checking.

sashahafner commented 1 month ago

Oh, yup. I see it in abm_variable():

image

Darn it.

sashahafner commented 1 month ago

So, we could use rates() or move the all the loading rate and cumulative loading calculations from abm to abm_regular and abm_variable in order to apply them to dat when it still has a single slurry loading rate for each interval. I like the second option, because abm is already slow and rates is already huge, and because in the second option, we could make a calc_loading function that simplifies the code in abm and both abm_*. What do yo think @fdalby ?

fdalby commented 1 month ago

Nice! I can work on the COD bit now. But I do not think we need it to go into rates(). COD loading rate is constant in an interval, so these should be accurate:

  dat$COD_load_rate <- dat$COD_conc_fresh * dat$slurry_prod_rate
  dat$COD_load_cum <- cumsum(dat$COD_load_rate * c(0, diff(dat$time))) + dat$COD_conc[1] * dat$slurry_mass[1]

I want to check that the c(0, is at the correct end (vs. c(diff(...), 0). Also earlier I have thought about whether we should have an option to exclude initial COD. Ultimately we might want to group all these options together.

abm() does not allow variable COD (VSd, etc.) concentration in feed over time (at least the main version doesn't), correct?

good if we can get it out of rates(). Yes abm does allow variable fresh conc. I think there is a demo in the vignette about it actually.

sashahafner commented 1 month ago

Oh, if it allows variable fresh concentrations with interpolation, we will need to use rates for cumulative values and loading rates. Darn. I'll try to look later today but am not sure I'll have much time. Will check in first via Whatsapp if I do so we don't duplicate efforts.

fdalby commented 1 month ago

So, we could use rates() or move the all the loading rate and cumulative loading calculations from abm to abm_regular and abm_variable in order to apply them to dat when it still has a single slurry loading rate for each interval. I like the second option, because abm is already slow and rates is already huge, and because in the second option, we could make a calc_loading function that simplifies the code in abm and both abm_*. What do yo think @fdalby ?

yes...I agree moving to abm_variable and regular is better way. I made the rates_calc argument, perhaps we could use this instead of yet another argument? Sine instant rates would alwasy go with instant loading output also?

fdalby commented 1 month ago

Oh, if it allows variable fresh concentrations with interpolation, we will need to use rates for cumulative values and loading rates. Darn. I'll try to look later today but am not sure I'll have much time. Will check in first via Whatsapp if I do so we don't duplicate efforts.

I guess you are right. We could solve it by just implementing another state variable COD_load_cum...

fdalby commented 1 month ago

Oh, if it allows variable fresh concentrations with interpolation, we will need to use rates for cumulative values and loading rates. Darn. I'll try to look later today but am not sure I'll have much time. Will check in first via Whatsapp if I do so we don't duplicate efforts.

I guess you are right. We could solve it by just implementing another state variable COD_load_cum...

sashahafner commented 1 month ago

Above commit seems to work, but I have not yet sorted out all the calculations in abm() after the lsoda() results. For future: adding a state variable is pretty simple but we have to remember to update:

  1. Initial state variable vector in abm.R around L 255: https://github.com/AU-BCE-EE/ABM/blob/ef657f41657445d0a2e473a181526df7c5279bd5/R/abm.R#L255
  2. Integer used to select only state variables in output from lsoda, currently 29, around L 220 in abm_variable.R: https://github.com/AU-BCE-EE/ABM/blob/ef657f41657445d0a2e473a181526df7c5279bd5/R/abm_variable.R#L220
  3. Same as 2 but in abm_regular, around L83: https://github.com/AU-BCE-EE/ABM/blob/ef657f41657445d0a2e473a181526df7c5279bd5/R/abm_regular.R#L83
sashahafner commented 1 month ago

Question for @fdalby:

Is it OK to just have total COD loading or do we need it for all substrates?

https://github.com/AU-BCE-EE/ABM/blob/ef657f41657445d0a2e473a181526df7c5279bd5/R/rates.R#L282

And are all these units COD anyway?

I guess I should add it separately for VS, and maybe degradable COD and degradable VS.

sashahafner commented 1 month ago

Above commit comments out some of the derived variables in question. I think we should eliminate a lot of them, but add the cumulative and rate variables needed to calculate them.

fdalby commented 1 month ago

Question for @fdalby:

Is it OK to just have total COD loading or do we need it for all substrates?

https://github.com/AU-BCE-EE/ABM/blob/ef657f41657445d0a2e473a181526df7c5279bd5/R/rates.R#L282

And are all these units COD anyway?

I guess I should add it separately for VS, and maybe degradable COD and degradable VS.

Hey. We need COD and C loaded also. I don't think we need each species separately, but the C load needs individuals inside rates. I do have the COD/C ratios in pars$COD_conv. So can easily be calculated. All the units are in gCOD except VSd_A and VSnd_A which is in g.

fdalby commented 1 month ago

Question for @fdalby: Is it OK to just have total COD loading or do we need it for all substrates? https://github.com/AU-BCE-EE/ABM/blob/ef657f41657445d0a2e473a181526df7c5279bd5/R/rates.R#L282

And are all these units COD anyway? I guess I should add it separately for VS, and maybe degradable COD and degradable VS.

Hey. We need COD and C loaded also. I don't think we need each species separately, but the C load needs individuals inside rates. I do have the COD/C ratios in pars$COD_conv. So can easily be calculated. All the units are in gCOD except VSd_A and VSnd_A which is in g.

fixed COD_load_cum. did not add C yet.