Compute steady states as a funciton of k, ODE

I will write a function called "compute_states" that takes as input a value for k and then returns as output ALL real-valued steady states for the ODE. The function will be assessed by running "Ceq = compute_states(k)" for many different values of k. I also want to be sure to filter out any complex-valued steady states. I know that "roots" and "imag" might be helpful.

Input: a scalar for k

Output: a vector of steady states, where the lenght of the vector is the number of all real-valued steady states.

ODE: dC/dt = 0.1*(C-20)*(23-C)*(C-26) +k

What I have so far:

Function:

function vectCeq = compute_states(k)
    % note that vectCeq is a vector of all real-valued steady states
    
    
    
    
    
    
end

Code to call my function:

k=0; % try for different values of k
vectCeq = compute_states(k)

Answer 1

[Anurag]

Hi Erik,

I understand that you need to compute all real valued steady states of given ordinary differential equation. For that you need to set the right-hand side of the equation to zero and solve for “C”.

Here is an example code snippet that you can refer to:

function vectCeq = compute_states(k)
    % Define the polynomial coefficients (using coeffs of the equation you
    % have written) 
    coefficients = [0.1,−6.9,145−k,−520];
    % Find the roots of the polynomial equation
    roots_eq = roots(coefficients);
    % Filter out complex roots and keep only real roots
    real_roots = roots_eq(imag(roots_eq) == 0);
    % Remove duplicate roots (if any)
    vectCeq = unique(real_roots);
    % Display the results
    disp(['For k = ', num2str(k), ', real-valued steady states:']);
    disp(vectCeq);
  
end
k = 0; % or any other value
vectCeq = compute_states(k); 

This function calculates the real-valued steady states for the given value of k and displays the results.

Documentation link for the function used:

https://in.mathworks.com/help/matlab/ref/roots.html

Hope this helps!

Regards,

Anurag

[Tommy]

May I ask how you got those coefficients (Line 4)? I'm not getting the same

[Sam Chak]

Hi @Erik Eriksson

If @Anurag's solution was what you have been seeking in the last few days, please consider clicking 'Accept' ✔ on the answer to close the issue.

Hi @Tommy,

Those coefficients are derived from the manually solved "expanded form" of the cubic polynomial.