Bounding posterior output of pharmacokinetic model

A small thing, the reference you’ve provided, using your notation has:

\mu(T) = \frac{D K_a (e^{-K_e T} - e^{-K_a T})}{C_l (K_a - K_e)}

(in contrast to what you’ve given, K_e is missing from the numerator), but your formula for T_p matches the correct version (your formula for P does not).

Wolfram alpha checks for the derivative and zero of the derivative.

Anyway, I think the main problem is that E in the form you’ve given can easily be larger than 1. Note that the absolute magnitude of \mu and P depends on D/C_l. But T_p does not depend on either D or C_l so I can make E arbitrarily large by increasing D/C_l while keeping the other parameters intact.