dde_test.jl

using DifferentialEquations
using Plots

function all_model(dz,z,h,p,t)
    b0,k0,T1,L1,g1l,a1,Al,t3,t4,et1,et2,m1,ka,kel,kpt,Kt,erel,km,vpt,kte,Km,d,g2l,g0,k,x,th1,th2,Rm,R50,a,n,tau1,tau2 = p
    hist1 = h(p,t-tau1)[1]
    hist2 = h(p,t-tau1)[3]
    hist3 = h(p,t-tau2)[1]
    hist4 = h(p,t-tau2)[3]
    hist5 = h(p,t-tau1)[7]
    ke(τ) = k0*τ/(1+τ)
    be(x,y) = b0*y/((1+x^m1)*(1+y))
    dz[1] = -g0*z[1]/(1+z[3]^a)-Rm*z[7]*z[1]/(R50+z[7])-(et1+et2)*ke(z[3])*z[1]-(1-et1-et2)*be(z[1],z[3])*z[1]+
        2*(1-et1-et2)*be(hist1,hist2)*z[4]*hist1 + et1*ke(hist2)*z[4]*hist1;
    dz[2] = -g2l*z[2] + Al*ke(hist4)*hist3;
    dz[3] = -k*z[3]+a1/(1+z[2]^n);
    dz[4] = z[4]*(-g0/(1+z[3]^a)-Rm*z[7]/(R50+z[7]) + g0/((1+hist2)^a)+Rm*hist5/(R50+hist5));
    dz[5] = -ka*z[5]+d;
    dz[6] = ka*z[5]-kel*z[6]-kpt*(1-erel)*z[6]/(Kt+z[6])-km*erel*z[6]/(Km+z[6]);
    dz[7] = vpt*kpt*(1-erel)*z[6]/(Kt+z[6])-kte*z[7];
end

function solver()
    b0=2.4;
    k0=0.18;
    T1=0.0782;
    L1=84.0;
    g1l=0.06;
    a1 = 1;
    Al =  563;
    t3 = 0.8;
    t4 = 3.5;
    et1 = 0.3;
    et2 = 0.3;
    m1 = 2;
    ka = 4.8;
    kel = 5;
    kpt = 29.8;
    Kt = 4.04*10^5;
    erel = 0.5;
    km = 655.8;
    vpt = 1;
    kte = 0.0714;
    Km  = 3.28*10^5;
    d = 60;
    g2l = 10^(-3);
    g0 = 0.01;
    k = 0.6;
    x = 0.3;
    th1 = 0.5;
    th2 = 0.36;
    Rm = 0.0022;
    R50 = 82.2;
    a = 0.83;
    n=7 ;
    tau1 = 2.8;
    tau2 = 6;
    lags=[tau1,tau2];
    p0 = (b0,k0,T1,L1,g1l,a1,Al,t3,t4,et1,et2,m1,ka,kel,kpt,Kt,erel,km,vpt,kte,Km,d,g2l,g0,k,x,th1,th2,Rm,R50,a,n,tau1,tau2);
    zz4(x,y)= exp((-g0/(1+x^a)-Rm*y/(R50+y))*2.8);
    u0 = [10^(-3),0,1.66,zz4(1.66,0.00619705),12.5,11.9975,0.00619705];
    #h0(p,t) = z0;
    h0(p,t) = u0;
    tspan=(0.0,1000.0);
    prob0 = DDEProblem(all_model, u0, h0, tspan, p0, constant_lags=lags)
    alg = MethodOfSteps(BS3())
    sol0 = solve(prob0,alg,progress=true)
    plot(sol0)
end