Error in using ode23t to solve a differential equation on a graph/network

Hi, I am trying to solve a differential equation about a substance spreading on a graph/network.

I am using the ode23t solver, but I am getting the following errors:

Error using odearguments (line 113)
Inputs must be floats, namely single or double.

Error in ode23t (line 143)
    odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);

Error in fr_net (line 48)
[T,H] =  ode23t(@(t,h) MYODE(t,h,G,N,e,n,z,initialnode,w,Q,lambda),tspan,h,opts);

Do you have any idea how to fix it?

Here following, I tried to summarise the relevant parts of the code:

% network
s = [1 1 2 2 3 3 4 5 5 6  6 7  7  8  9  9 10 10 11 11 12 13 14 15];
t = [2 5 3 6 4 7 8 6 9 7 10 8 11 12 10 13 11 14 12 15 16 14 15 16];
G = graph(s,t);
N = numnodes(G); % 16 nodes
e = table2array(G.Edges);
n = table2array(G.Nodes);

% networks coordinates (x,y,z)
x = ...
y = ...
z = ...
G.Nodes.X = x'; G.Nodes.Y = y'; G.Nodes.Z = z';

% inputs/initial conditions/parameters
initialnode = 1;
h = repelem(0.01,N); 
h(initialnode) = 0.5; 
w = repelem(1,N); 
Q(initialnode) = 0; 
lambda = 0.3;

% ode23t solver
tspan = [0 3600]; % t = time
opts = odeset('MaxStep',3600);
[T,H] =  ode23t(@(t,h) MYODE(t,h,G,N,e,n,z,initialnode,w,Q,lambda),tspan,h,opts);

% MYODE 
function dhdt = MYODE(t,h,G,N,e,n,z,initialnode,w,Q,lambda),tspan,h,opts);
    for i = 1 : N
        j = randi([1 N]);
        % some operation/equation (I use cells for Q, since it could "store" several numbers)
        Qcell{i} = ... % some operation/equation
        dhdt_tmp{i} = cell2mat(Qcell{i}) / function_ij; % differential equation
    end
    dhdt = cell2mat(dhdt_tmp);
end

Answer 1

[Walter Roberson]

function dhdt = MYODE(t,h,G,N,e,n,z,initialnode,w,Q,lambda),tspan,h,opts);

That is an illegal function declaration. I suspect you made a mistake in your copying.

You appear to be building dhdt_tmp as a row cell array. cell2mat() would convert that to a row vector, but the output needs to be a column vector.

Inputs must be floats, namely single or double.

At that polnt in the code, that error would result if tspan or h are not numeric (which they certainly look to be) or if the ode called on the first time and first boundary condition, returns something that is not numeric.

[Sim]

Thanks a lot Walter for your kind reply!

• Sorry, but I cannot understand which part of MYODE function would be an "illegal function declaration". Which part should be modified?
• About the dhdt_tmp, I just commented both the function MYODE and the ode23t call for checking if dhdt_tmp, and consequently dhdt, were row or column arrays. I did a small modification to dhdt_tmp to change it from a row array to a column array:

% row-to-col array
dhdt_tmp = dhdt_tmp';

% since the dhdt_tmp 's cells could contain more than one number, 
% I just take the first ones and I store them in dhdt
idx = ~cellfun('isempty',dhdt_tmp);
dhdt = zeros(size(dhdt_tmp));
dhdt(idx) = cellfun(@(v)v(1),dhdt_tmp(idx));

And I got this (however, the previous error is still present):

>> dhdt_tmp

dhdt_tmp =

  16×1 cell array

    {2×1 double }
    {2×1 double }
    {2×1 double }
    {[ 0.014678]}
    {2×1 double }
    {2×1 double }
    {2×1 double }
    {[-0.039071]}
    {2×1 double }
    {2×1 double }
    {2×1 double }
    {[-0.026954]}
    {[-0.018813]}
    {[-0.037621]}
    {[-0.050246]}
    {0×0 double }

>> dhdt

dhdt =

     0.057469
       -2.816
  -0.00027738
     0.014678
      -3.9819
    -0.028159
     0.014802
    -0.039071
    -0.013186
    -0.033483
   0.00055477
    -0.026954
    -0.018813
    -0.037621
    -0.050246
            0

• Yes, tspan and h are numeric. Therefore, I guess I should check "if the ode called on the first time and first boundary condition, returns something that is not numeric", but I dont know how to debug the ode23t solver when it runs.......

[Walter Roberson]

function dhdt =   MYODE(t,h,G,N,e,n,z,initialnode,w,Q,lambda),tspan,h,opts);
                       1                                    0             -1

The numbers are the active bracket nesting count "after" the character they are below. You have more ) than you have ( on the line.

The code on that line is equivalent to the lines of code

function dhdt =   MYODE(t,h,G,N,e,n,z,initialnode,w,Q,lambda)
disp(tspan)
disp(h)
disp(opts)
);           %this is a syntax error

and displaying tspan and opts would each fail because they are not defined within the scope of the function.

I think it very likely that the line you wanted was just

function dhdt =   MYODE(t,h,G,N,e,n,z,initialnode,w,Q,lambda)

[Walter Roberson]

Put a breakpoint at the beginning of the line

[T,H] =  ode23t(@(t,h) MYODE(t,h,G,N,e,n,z,initialnode,w,Q,lambda),tspan,h,opts);

Before executing the line, execute

result = MYODE(tspan(1), h(:), G, N, e, n, z, initialnode, w, Q, lambda)
class(result)    %needs to be double or single
size(result)     %needs to be a column vector with same number of elements as h has

[Sim]

Thanks a lot!

Oh, yes, I got it... it was a very stupid mistake in copying...Indeed, I copied the ode23t argument, instead of MYODE function. Ok, I need to drink some more coffee...!

About the second reply, I performed the debug you indicated

result = MYODE(tspan(1), h(:), G, N, e, n, z, initialnode, w, Q, lambda)
class(result)    %needs to be double or single
size(result)     %needs to be a column vector with same number of elements as h has
[T,H] =  ode23t(@(t,h) MYODE(t,h,G,N,e,n,z,initialnode,w,Q,lambda),tspan,h,opts);

And I got this:

Undefined function 'neighbors' for input arguments of type 'double'.

Error in fr_net>MYODE (line 64)
    neighb = neighbors(G,i);

Error in fr_net (line 48)
result = MYODE(tspan(1), h(:), G, N, e, n, z, initialnode, w, Q, lambda)

Indeed, inside MYODE I called the function neighbors(G,i) to find the neighbours of the i-node (this is one of the code lines that I omitted just to simplify my question):

function dhdt = MYODE(t,h,G,N,e,n,z,initialnode,w,Q,lambda)
    % the following loop is calculating some equation between i-node and j-node
    for i = 1 : N
        neighb = neighbors(G,i);
        % other stuff
        % other stuff
        % other stuff
    end
end

...So, now I have the doubt that I cannot use graph G as input - and its relative functions - inside MYODE for solving the diff.equation...

[Walter Roberson]

You can use graph objects as extra inputs, no problem. However, at that particular breakpoint for debugging, first ask

class(G)

to make sure it is still a graph object at the time of the call.

[Sim]

ok... doing this:

class(G)
result = MYODE(tspan(1), h(:), G, N, e, n, z, initialnode, w, Q, lambda)

I obtain:

ans =

    'graph'  

...but I still get this error:

Undefined function 'neighbors' for input arguments of type 'double'.

Error in fr_net>MYODE (line 66)
    neighb = neighbors(G,i);

Error in fr_net (line 49)
result = MYODE(tspan(1), h(:), G, N, e, n, z, initialnode, w, Q, lambda)

....Probably I should put G as global variable (or open a new question about this last error....) ?!

[Sim]

Problem solved! There was a typo in the function!

Thank you very much Walter! You gave me the right advices and diagnostic tools to understand my mistakes!

All the best,

Sim