Hi everyone,

I am building a model of a growing cell population. I want to make use of callback to evaluate the ‘division’ event, similar to the example in this tutorial - Event Handling and Callback Functions · DifferentialEquations.jl

However unlike in the tutorial, I want to maintain a spatial relationship between parent and daughter cells - if u represents a 1-D array of cells and their protein concentration, then when the protein concentration of a specific cell idx hits 1, it divides and its daughter cell is placed next to it at u[idx+1], shifting all the other cells one space along (as if we called insert!(u,idx,val)), rather than at the end of the array as in the tutorial.

In addition, I want the possibility for multiple cells to be able to divide - suppose we start with two cells with equal starting concentrations, then when one reaches a protein concentration of 1 the other will have also (assuming identical protein dynamics), meaning that 2 new cells should be created rather than just 1 (as is the behaviour in the tutorial’s model).

I have tried the below code, and it doesn’t produce the expected behaviour.

```
const α = 0.3
function f(du,u,p,t)
for i in 1:length(u)
du[i] = α*u[i]
end
end
function condition(u,t,integrator)
1 - maximum(u) # has any cell hit 1?
end
function affect!(integrator)
u = integrator.u
idx = findall(y-> y == 1,u) # find all cells which have hit concentration == 1
resize!(integrator,length(u)+length(idx)) # resize solution to accomodate all new divisions
for id in reverse(idx) # for each cell hitting 1, divide:
for i in reverse(id:length(u)-1) # this just replicates the behaviour of insert!(u,idx+1,u[idx]/0.5) and u[idx] = 0.5*u[idx]
if i == id + 1
u[i+1] = u[i]
u[i] = 0.5*u[id]
elseif i == id
u[i] = 0.5*u[id]
else
u[i+1] = u[i]
end
end
end
nothing
end
callback = ContinuousCallback(condition,affect!)
u0 = [0.2]
tspan = (0.0,10.0)
prob = ODEProblem(f,u0,tspan)
sol = solve(prob,Tsit5(),callback=callback);
plot(sol.t,map((x)->length(x),sol[:]),lw=3,
ylabel="Number of Cells",xlabel="Time")
```

and

```
ts = range(0, stop=10, length=100)
plot(ts,map((x)->x[1],sol.(ts)),lw=3,
ylabel="Amount of X in Cell 1",xlabel="Time")
```

I might have the wrong approach here - I did think about using VectorContinuousCallback where each event in out was a particular cell dividing. However, this would require the length of out to be dynamic (the same as u) which I don’t think is supported.

I am new to this part of diffeq and so any help would be appreciated! Thanks