Undefined behavior in pk1 when KA=CL/V

[toddyoder]

The concentrations in the closed form pk1 model appear to be either all zeros or NaN when KA is chosen to be equal to CL/V.

library(dplyr)
library(mrgsolve)

# KA != CL/V
modlib("pk1") %>% 
  param(CL = 0.391, V = 30, KA = 0.013) %>% 
  ev(amt = 100) %>% 
  mrgsim() %>% 
  as_tibble %>% 
  pull(CENT) %>% 
  summary

# KA = CL/V
modlib("pk1") %>% 
  param(CL = 0.39, V = 30, KA = 0.013) %>% 
  ev(amt = 100) %>% 
  mrgsim() %>% 
  as_tibble %>% 
  pull(CENT) %>% 
  summary

# KA = CL/V
modlib("pk1") %>% 
  param(CL = 1, V = 20, KA = 0.05) %>% 
  ev(amt = 100) %>% 
  mrgsim() %>% 
  as_tibble %>% 
  pull(CENT) %>% 
  summary

Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
0.000   6.371  12.868  12.162  18.312  22.829 

Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
0       0       0       0       0       0 

Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
0       0       0       0       0       0      24 

[kylebaron]

mrgsolve/src/odeproblem.cpp

Line 387 in 3d5d701

a[0] = ka/(ka-alpha[0]);