How to solve simultaneously a system of ODEs containing both initial and boundary value problems

I have a system of first order ODEs to solve, most of the ODEs are initial value problems , and about two of these ODEs are boundary value problems. For simplicity I have given an illustration below.

I'm currently using ode15s solver to implement the solution for the process as an initial value problem, however the last ODE is best solved as a boundary value problem, where the upper limit of it's boundary is a function of the solution of ODE, dy(1) i.e. y(1).

dy(1) = f{ y(1) }, y(0) = a (constant)

dy(2) = f{ y(1), y(2), y(3)}, y(0) = b (constant)

dy(3) = f{ y(3), y(2) }, y(0) = c (constant) & y(end) = y(1) (variable)

My implementation process involves the use of a script that contains all necessary constant parameters (p), initial condition values (IC) & plots and another script which contains all the ODEs

Therefore the solution syntax to run the process is as illustrated below:

[time,Y] = ode15s(@(t,y) model03cBatchPhBubMod(t,y,p),[0,t_end],IC)

Answer 1

[Alan Stevens]

Your first two odes could be solved by ode15s, since they don't depend on y(3), then, having solved them, you could investigate function bvp4c to solve the third ode.

[Samuel Emebu]

Thanks for your comment, but y(3) depends on y(2) and I have to solve them together to get a correct solution . I can't solve the third ODE separately.

[Alan Stevens]

Having solved for y(2) as a function of the independent variable, you can use its values whenever they are needed in the y(3) ode.

[Samuel Emebu]

Sorry, I have Edited/corrected my question. dy(2) and dy(3) are both function of each other.

[Alan Stevens]

Ok. In that case only y(1) can be done with ode15s on its own. Have a look at bvp4c for the others.

[Samuel Emebu]

So I have to create two function files (for IVP and BVP)? if so how do I connect the variables y(1), y(2) and y(3) across my scripts without getting undefined variable error?

[Alan Stevens]

All the functions can be written in one script. The functions must come at the end of the script. Call the IVP at the beginning of the script first, then go on to call the BVP (but the calling section will still be written before all the functions).

% Data...
tspan = ...
y1init = ...
[t, y1] = ode15s(@y1fn ...)

% data for bvp functions
% call bvp functions
 
% extract /plot/display... results

function y = y1fn(...)
    ...
    end
function y = bvpinit(...)
    ...
    
end

... etc.

.

[Samuel Emebu]

I'm currently working on it, the equations are much. I will give feedback on my success.

Thanks