Change initial condition by set_state! DynamicalSystems.jl

Hello!
I’m trying to change the initial condition using the set_state! function and I integrate the system, but when I do this, I see that the initial state has not changed. Why doesn’t the initial condition change? Should I use reinit! ?

Thank for your helps!
Code:

using StaticArrays, DifferentialEquations, DynamicalSystems

@inbounds function TM6_glial_ECM(var, par, t)

    g = log( 1.0 + exp( ( (par[9] + par[6] * var[5]) * var[3] * var[2] * var[1] + par[11]) / par[5] ) );
    U = par[10] + par[12] / ( 1.0 + exp( -50.0 * ( var[4] - par[25] ) ) );
    Hy = 1.0 / ( 1.0 + exp( -20.0 * ( var[2] - par[26] ) ) );
    Hecm = par[17] - (par[17] - par[18]) / (1.0 + exp( -(var[1] - par[20]) / par[19] ) );
    Hp =  par[21] - (par[21] - par[22]) / (1.0 + exp( -(var[1] - par[23]) / par[24]) );

    dE = (-var[1] + par[5] * g) / par[1];
    dx = (1.0 - var[2]) / par[2]  -var[3] * var[2] * var[1];
    du = (U - var[3]) / par[3]  + U * (1.0 - var[3]) * var[1];
    dy = (-var[4]) / par[4] + par[13] * Hy;
    decm = -( par[7] + par[16] * var[6] ) * var[5] + par[14] * Hecm; 
    dp = -par[8] * var[6] + par[15] * Hp;

    return SVector(dE, dx, du, dy, decm, dp);
end

function TM6_glial_ECM_get_params()
    τ = 0.013; τD = 0.15; τF = 1.0; τy = 1.8;   
    α = 1.5; αecm = 0.001; αp = 0.01;
    J = 3.07; U0 = 0.3; ΔU0 = 0.305; 
    β = 0.438; βp = 0.01; βecm = 0.01;
    ecm0 = 0.0; ecm1 = 1.0; kecm = 0.15; θecm = 25.6;
    p0 = 0.0; p1 = 1.0; kp = 0.05; γp = 0.1; θp = 26.0; 
    ythr = 0.5; xthr = 0.9;
    αE = 5.0;
    I0 = -1.741;

    param = [τ, τD, τF, τy, α, αE, αecm, αp, J, U0, I0, ΔU0, β, βecm, βp, γp, ecm0, ecm1, kecm, θecm, p0, p1, θp, kp, ythr, xthr];
    return param;
end

t_tt = 3000;
tol = 1e-15;
integrator_setting = (alg = Vern9(), adaptive = true, abstol = tol, reltol = tol);

u0start = [0.9445509341100914, 0.74116702856987, 0.7361196042973006, 0.0646914552140727, 0.15145764079879162, 0.0009327645775731449];
params = TM6_glial_ECM_get_params();
ds = CoupledODEs(TM6_glial_ECM, u0start, params, diffeq = integrator_setting);

ds

set_parameter!(ds, 6, 4.2);
set_state!(ds, zeros(6));

ds

tr, _ = trajectory(ds, 500);
tr[1]

I don’t understand… The picture you paste clearly shows that the state does change. First and second time you print the system show different states.

?

Yes, but when I calculate the trajectory, the first point from the trajectory is not the new state that I set

ok what’s the versions you are using.

I use this version DynamicalSystems v3.2.1

Oh, of course, there is no bug. Have a look at the documentation of trajectory function: DynamicalSystemsBase.jl · DynamicalSystemsBase.jl

The default state is the initial_state(ds). of course, the initial_state does not change when setting the current state. you need to explicitly give a starting state as a third argument to trajectory.

Thank you for help

1 Like

Can I please ask another question here without creating a new topic?

better start a new topic! This one has been solved!