MATLAB error using ODE15s: Function must return a column vector in

Question

I am solving system of four differential equations using ODE15s solver in MATLAB. I want to store values of N1 N2 N3 N4 for each seven different values of W and PION. I also want to plot N1 N2 N3 N4 vs time for each W and PION. Here is my MATLAB code:

% Code to solve system of four simultaneous differential equation

clc

N0 = [9.66e17;0;0;0];   % Intital conditions of N1, N2, N3, N4 (9.66e17,0,0,0)
[T, N] = ode15s(@spl_eq,[0 126e-15], N0);  % Intital conditions of time from 0 to 126e-15
plot(T,N(:,4),T,N(:,3),T,N(:,2),T,N(:,1))   % Plotting of N1 N2 N3 N4 vs Time

function dN = spl_eq( t,N )

t = 126e-15;

dN = zeros(4,7);    % Blank array to store N1 N2 N3 N4 from solution for each values of seven values of W and PION
Q = 10e10;      % Constant value
A = 10e8;       % Constant value
Ppd = 10e10;    % Constant value

% W and PION are not constant and they have seven values as row vectors given below:

W = [0.2716 0.4243 0.6110 0.8317 1.0862 1.3748 1.6972]*1e9;
PION = [0.7072 0.8874 1.0254 1.2375 1.4143 1.5911 1.7679]*1e11;


for j = 1:length(W)
    dN(1) = (-W(j)*N(1)) + ((Q+A)*N(2)) + ((W(j)+Q+A)*N(3)) + 0;            % first equation 
    dN(2) = 0 + (-(Q+A)*N(2)) + ((Q+A)*N(3)) + 0;                             % second equation 
    dN(3) = (W(j)*N(1)) + 0 - (W(j)+(2*Q)+(2*A)+ Ppd + PION(j))*N(3) + 0;     % third equation 
    dN(4) = 0 + 0 + PION(j)*N(3) + 0;                                           % fourth equation 
end

end

But I am getting the following errors:

Error using odearguments
SPL_EQ must return a column vector.
Error in ode15s (line 148)
    odearguments(odeIsFuncHandle, odeTreatAsMFile, solver_name, ode, tspan, y0, options, varargin);

Error in test_soln_rate (line 6)
[T, N] = ode15s(@spl_eq,[0 126e-15], N0);  % Intital conditions of time from 0 to 126e-15

Answer 1

If your model is somewhat dx = f(x, p), for extra parameter p, please check documentation on this regard: https://www.mathworks.com/help/matlab/ref/ode15s.html .

Your parameters must be either time-constant or a timeseries p(t). In case of timeseries p(t), you must use the closure concept: https://www.mathworks.com/matlabcentral/fileexchange/18223-closures-in-the-matlab-language