Passing parameters to dde23 slows it down

I am having a set of Delay Differential Equations which I am trying to solve using dde23. When I define the function in the form,

function change = eqn(t, X, Z)
%Defining parameter 'k' here
%Inserting the set of equations here.
end

, the code runs.

However if I try to pass a parameter into the function

function change = eqn(t, X, Z, k)
%Inserting the set of equations here.
end

And accordingly modify dde23 line in the main code:

K = 1;
sol = dde23(@(t, x, z) eqn(t, x, z, K), tau, init_vals, tspan, options);

The code takes a lot of time to run. If anyone knows why and how I can improve this, please let me know!

I have to be able to pass the parameter to the function since I am going to use dde23 in a loop with a range of the parameter values.

[Walter Roberson]

Can you post your two different versions of the code ?

Answer 1

[Peter O]

Anonymous functions carry some overhead. I wonder if the way it's parsing the arguments to dde23 is causing it to rebuild the anonymous function on each subcall.

Does defining the anonymous function ahead of sending it as an argument help?

for ix=1:N
    k = rand(1);
    dfcn = @(t,x,y,z) eqn(t,x,y,z,k);
    sol = dde23(dfcn, tau, init_vals, tspan, options)
end

[Sneha Srikanth]

I tried what you have suggested. The code still takes a really long time.

[Peter O]

Okay -- something else to try: Profile the code. Turn on the profiler using the instructions at the link below. Does the report it generates give you any insight on where it's spending the time?

https://www.mathworks.com/help/matlab/matlab_prog/profiling-for-improving-performance.html

[Walter Roberson]

sol = dde23(@(t, x, z) eqn(t, x, z, K), tau, init_vals, tspan, options);

would not rebuild the function every time. It would evaluate each of the arguments in turn, getting out a function handle for the first argument, and that function handle would be passed to dde23 -- the same as if you had constructed the function handle before the function.

[Sneha Srikanth]

I found the mistake. It was a mathematical error in the set of equations that caused this. To be precise, it was a missing cosine term which, if present, would have scaled down the output of the function eqn a lot. I still don't get why the code ran faster when the parameters were defined within the function. Won't the mistake be the same regardless of passing the parameters or defining them inside itself?

I am putting the function's code below.

function change = k_eqn(t, state, Z, K)
xf = 0.3;
k = pi*(1:3)';
uf = cos(xf*k')*Z(1:3); % the mistake was: uf = k'*Z(1:3); 
U = state(1:3);
P = state(4:end);
damp = [0.02; 0.01];  
changeU = -k.*P;
changeP =  -2*damp.*k.*P + k.*U + K.*sin(k*xf).*(sqrt(abs(1/3+uf))-sqrt(1/3));
change = [changeU; changeP];
end

[Peter O]

I'm still a little confused, but I think there's another thread to pull.

When you said there was an error in a formula that led to a change in the size of the derivatives, my first impression is that the wrong parameters alter the stiffness of the problem. When, the derivatives become more sensitive to changes in their values, dde23 has to take smaller steps to keep its error tolerances within bounds, which leads to (a lot) more function calls and longer time to run.

https://en.wikipedia.org/wiki/Stiff_equation

But you're saying that in both of your cases you had the error in computing k and there was a noticeable difference in overhead? Since Z seems to affect the parameters, where does it get calculated and is it dependent on anything? There also seems to be a big K and a lowercase k in your solutions. Is everything supposed to be one k? Finally, assuming state is 6 elements long, in the code you posted there seems to be an element-wise [2x1] x [3x1] x [3x1] operation in changeP, which MATLAB should throw an error on for bad dimensions. Is damp missing a third value there?

[Sneha Srikanth]

There were some more irrelevant calculation steps which I had replaced with constants. I missed out the third element of damp while replacing. My bad, sorry.

k is different from K, Sorry for the lack of context. This is a problem about the velocity and pressure fluctuations in a pipe in which a heater is present at position xf from one end of the pipe.

k is a vector of [pi; 2*pi; 3*pi] while K is the heater power.

uf can be thought of as velocity of air at a position xf in a pipe.

The velocity u at a position x in the pipe is U1*cos(pi*x) + U2*cos(2*pi*x) + U3*cos(3*pi*x). Similar expression for pressure p in the pipe is present with U1, U2, etc. replaced by P1, P2 and cosine replaced by sine.

There is a delay in the heat interaction between the heater and any point in the pipe. That has a time delay tau.

Z is the vector containing the state ([U1; U2; U3; P1; P2; P3]) at time (t-tau). I am finding the velocity of air at the heater location xf at time t-tau since it affects the state (velocity and pressure components) at t. Z is calculated by dde23 itself