Using vector of vectors for initial guess in TwoPointBVProblem

The following works:

using BoundaryValueDiffEq
using OrdinaryDiffEq

#System Constants
ss = 1 #excitatory parameter
sj = 0 #inhibitory parameter
glb = 0.05
el = -70
gnab = 3
ena = 50
gkb = 5
ek = -90
gtb = 2 
et = 90 
gex = 0.1 
vex = 0
gsyn = 0.13
vsyn = -85
iext = 0.41
eps = 1
qht = 2.5

#System functions
function f(v,h,r)
    -(glb*(v-el)+gnab*(1/(1+exp(-(v+37)/7)))^3*h*(v-ena)+gkb*(0.75*(1-h))^4*(v-ek)+gtb*(1/(1+exp(-(v+60)/6.2)))^2*r*(v-et))-gex*ss*(v-vex)-gsyn*sj*(v-vsyn)+iext
end

function g(v,h)
    eps*((1/(1+exp((v+41)/4)))-h)/(1/((0.128*exp(-(v+46)/18))+(4/(1+exp(-(v+23)/5)))))
end

function k(v,r)
    qht*((1/(1+exp((v+84)/4)))-r)/((28+exp(-(v+25)/10.5)))
end

#Dynamical System
function TC!(du,u,p,t)
    v,h,r = u
    
    du[1] = dv = f(v,h,r)
    du[2] = dh = g(v,h)
    du[3] = dr = k(v,r)
end

#Finding initial guesses by forward integration
T = 7.588145762136627 #orbit length
u0 = [-40.296570996984855, 0.7298857398191566, 0.0011680534089275774]
tspan = (0.0,T)
prob = ODEProblem(TC!,u0,tspan,dt=0.01)
sol = solve(prob, Rodas4P(), reltol = 1e-12, abstol = 1e-12, saveat = 0.5)

#The BVP set up
function bc_po!(residual, u, p, t)
    residual[1] = u[1][1] - u[end][1]
    residual[2] = u[1][2] - u[end][2]
    residual[3] = u[1][3] - u[end][3]
end

#This is the part of the code that has problems
bvp1 = TwoPointBVProblem(TC!, bc_po!, sol.u, tspan)
sol6 = solve(bvp1, GeneralMIRK4(), dt = 0.5) 

with the patch Fix vector of vector initial conditions with GeneralMIRK4 by ChrisRackauckas · Pull Request #82 · SciML/BoundaryValueDiffEq.jl · GitHub

2 Likes