geco-bern / rsofun

Implements the Simulating Optimal FUNctioning framework for site-scale simulations of ecosystem processes, including model calibration. It contains Fortran 90 modules for the P-model, SPLASH, and BiomeE models.
https://geco-bern.github.io/rsofun/
GNU General Public License v3.0
25 stars 28 forks source link

Temporal instability of p-model #26

Closed khufkens closed 3 years ago

khufkens commented 3 years ago

When calling the p-model with the same input parameters several times the returned results are not the same.

@stineb are there any non-deterministic parts in the p-model (optimized), relying on internal optimization? I'm not sure which parameter set to use with the additions of tau_acclim_tempstress etc. However, even if parameters are nonsensical the output should be stable (if deterministic).

Worked example returning jmax25

library(rsofun)

set.seed(1)

# load parameters
params_modl <- list(
  kphio           = 0.09423773,
  soilm_par_a     = 0.33349283,
  soilm_par_b     = 1.45602286,
  tau_acclim_tempstress = 10,
  par_shape_tempstress  = 0.0
)

# read in demo data of package
df <- p_model_drivers

# loop over 100 runs return jmax25
values <- lapply(1:100,function(x){
  # run the SOFUN Fortran P-model
  mod <- run_pmodel_f_bysite( 
    df$sitename[1],
    df$params_siml[[1]],
    df$siteinfo[[1]],
    df$forcing[[1]], 
    df$df_soiltexture[[1]],
    params_modl = params_modl,
    makecheck = FALSE
  )

  mod$jmax25
})
stineb commented 3 years ago

As long as par_shape_tempstress is zero, this takes no effect and values of tau_acclim_tempstress. I cannot reproduce this. My R session crashes when running the above code after merging the pull request from your fork.

Instability of Jmax. This should not happen. Is it the only variable that is unstable? Jmax is not tested separately in the benchmarking.

khufkens commented 3 years ago

Just checked the unit tests, seems that segfault in mac and windows. I only changed the CFLAG and _Bool to Int to avoid errors as commented in the code.

stineb commented 3 years ago

This may be a problem: In sofun_r.f90:

logical(kind=c_bool), intent(in) :: soilmstress

and in wrappersc.c:

    int    *soilmstress,

Should probably be received in Fortran as a integer (integer(kind=c_int)) if passed on in wrappers as integer. And what about R? Now we have

as.logical(params_siml$soilmstress)

Shouldn't it be converted to integer (0 or 1) already here?

khufkens commented 3 years ago

Not sure how if() statements are dealt with in Fortran and if 0 / 1 i accepted as identical to .TRUE. / .FALSE.

as.logical(TRUE)

is TRUE

as.integer(TRUE)

1

khufkens commented 3 years ago

See the note on Boleans in this post

https://www.r-bloggers.com/2018/12/the-need-for-speed-part-1-building-an-r-package-with-fortran-or-c/

Basically it is asking for trouble, depending on the compiler and platform used.

khufkens commented 3 years ago

Although I doubt that this is the issue here.

The runs on my old setup worked on both MacOS and Windows.

https://github.com/bluegreen-labs/rsofun_reference/runs/3132618063?check_suite_focus=true

khufkens commented 3 years ago

Missing tmin and tmax here:

https://github.com/stineb/rsofun/blob/463af7b984d4541e03a5994c8b0069d37aadf674/R/run_pmodel_f_bysite.R#L58

Need to update the driver file as well it seems.

khufkens commented 3 years ago

Fixed as per https://github.com/stineb/rsofun/pull/27