using ParameterizedFunctions, OrdinaryDiffEq, Test

f = @ode_def begin
    dGut    = -Ka*Gut
    dCent   =  Ka*Gut - CL*(Cent/Vc) -Q*(Cent/Vc) + Q*(Periph/Vp)
    dPeriph =  Q*(Cent/Vc)  - Q*(Periph/Vp)
    dResp   =  Kin*(1-(IMAX*(Cent/Vc)^γ/(IC50^γ+(Cent/Vc)^γ)))  - Kout*Resp
end Ka CL Vc Q Vp Kin Kout IC50 IMAX γ

p = (Ka = 2.769674683505292, CL = 0.8315548845614537, Vc = 45.87339297906569,
      Q = 1.4544520084387496, Vp = 9.460961128900388, Kin = 6.360418843010778,
      Kout = 3.524234094210465, IC50 = 0.4483840673123028, IMAX = 1.0, γ = 2.0,)

u0 = [0.35793622925128044,
      0.7123008259214687,
      0.10685666303615937,
      0.00]
prob = ODEProblem(f, u0, (0,26.), p)
sol = solve(prob, AutoTsit5(Rodas5(), dtfac=2))

@test sol([21.2])[1][4] ≈ 1.803 atol = 1e-2

#=
using Plots
plot(sol, vars=4, lab="dtfac=1")
sol = solve(prob, AutoTsit5(Rodas5(), dtfac=2))
plt1 = plot!(sol, vars=4, lab="dtfac=2")
plt2 = plot(sol.t, sol.alg_choice, lab="alg_choice")
plot(plt1, plt2, layout=(2,1))
plot!(;title="AutoTsit5(Rodas5())",dpi=200)
=#
