how to write a variable parameter inside recall function?

In the following code, I want to enter rho as a linear function, and rho at the same time is a parameter inside the kinetics function.

I tried that, but I got error regarding this line.

>> [t, nc] = ode45(@(t,nc) kinetics(t,nc,rho,beta,betasum), tspan, nc0)

>>beta = [0.00021; 0.00141; 0.00127; 0.00255; 0.00074; 0.00027];
>> betasum = sum(beta);
>> [rho] = (@(rho) InReactivity(t));
>> c0 = zeros(6,1);
>> tspan = [0, 0.1];
>> nc0 = [n0; c0];
>> [t, nc] = ode45(@(t,nc) kinetics(t,nc,rho,beta,betasum), tspan, nc0);

rho function:

function [rho] = InReactivity(t)
if 0<t<1
    rho = (0.002/betasum)+t
else
    rho = 0.002/betasum
end
end

Kinetics function

function dncdt = kinetics(~,nc,rho,beta,betasum)
        
        L = 0.0001;
        lam = [0.0126; 0.0337; 0.111; 0.301; 1.14; 3.01];
        n = nc(1);
        c = nc(2:7);
        dndt = (rho - betasum)/L + sum(lam.*c);
        dcdt = beta*n/L - lam.*c;
        
        dncdt = [dndt; dcdt];        
end

Answer 1

[Simran]

Hi @Hadeer Abdullah,

I see that you are encountering the error in

[t, nc] = ode45(@(t,nc) kinetics(t,nc,rho,beta,betasum), tspan, nc0)  

line of your code.

The problem arises from how you are defining and using “rho”. You can take following steps to fix it:

1.) Instead of defining it as a variable, you can define it as a function that takes time “t” as input. This will allow it to be evaluated at each time step during the ODE solution process.

2.) Pass the “InReactivity” function as an argument, instead of depending on “betasum” as a global variable.

This is the corrected code script:

rho = @(t) InReactivity(t, betasum); 

3.) Lastly, use a "function handle" to pass “rho” as a time-dependent function to “ode45”.

Here is the plot between “neutron density” vs “time”, I got by following the above workflow:

You can refer to the following documentation more information:

Anonymous Function: https://www.mathworks.com/help/releases/R2021a/matlab/matlab_prog/anonymous-functions.html

“ode45”:

https://www.mathworks.com/help/releases/R2021a/matlab/ref/ode45.html

"function handle":

https://www.mathworks.com/help/releases/R2021a/matlab/matlab_prog/creating-a-function-handle.html

Hope this helps!