Hopf Bifurcation diagram for 3D system

Question

kdv0 2024-7-11

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2136538-hopf-bifurcation-diagram-for-3d-system

编辑： Francisco J. Triveno Vargas 2024-7-15

Hello I want to plot the Hopf bifurcation diagram, I have this code but I'm getting nothing other than straight lines.

% Parameters
r=3;
a=0.2;
b=2;
gamma=0.3;
m=0.92;
c=2.5;
alpha=10;
p=100.001;
q=0.01;
% Define the range of parameter values for bifurcation diagram
eta_values = linspace(0, 3, 1000);  % Range of 'eta' values
% Arrays to store results
x_result = zeros(size(eta_values));
y_result = zeros(size(eta_values));
z_result = zeros(size(eta_values));
% Initial conditions (for detecting the Hopf bifurcation)
x0 = 1000;
y0 = 200;
z0 = 600;
% Iterate over each 'eta' value
for i = 1:length(eta_values)
    eta = eta_values(i);
    
    % Solve the system using ODE45
    tspan = [0 1000];  % Time span for integration
    odefun = @(t, Y) [r*Y(1)*(1-Y(1)/Y(3))-a*(1-m)*Y(1)*Y(2)/(1+gamma*(1-m)*Y(1));
                      b*(1-m)*Y(1)*Y(2)/(1+gamma*(1-m)*Y(1))-c*Y(2);
                      (Y(3)-alpha)*(p-Y(3))/(q+Y(3));];
    [~, Y] = ode45(odefun, tspan, [x0; y0; z0]);
    
    % Check for Hopf bifurcation (change in stability)
    % Look for oscillatory behavior in the solution
    if i > 1
        if abs(Y(end, 1) - x_result(i-1)) > 0.01
            fprintf('Hopf bifurcation detected at eta = %f\n', eta);
        end
    end
    
    % Record the final values (or a characteristic of the solution)
    x_result(i) = Y(end, 1);
    y_result(i) = Y(end, 2);
    z_result(i) = Y(end, 3);
end
% Plot the bifurcation diagram
figure;
plot(eta_values, x_result, '.k', 'MarkerSize', 1);
hold on;
plot(eta_values, y_result, '.b', 'MarkerSize', 1);
plot(eta_values, z_result, '.r', 'MarkerSize', 1);
xlabel('\eta');
ylabel('Steady State Values');
legend('x', 'y', 'z');
title('Hopf Bifurcation Diagram');
grid on;

4 个评论
显示 2更早的评论隐藏 2更早的评论

Umar 2024-7-11

Hi kdvO,

After analyzing your code, the issue causing the plotted lines to appear straight instead of showing the expected Hopf bifurcation behavior lies in how the detection of Hopf bifurcation is implemented within the loop. The current logic for detecting the bifurcation relies on comparing the difference between the last value and the current value of `x_result` to a threshold of 0.01. However, this approach may not accurately capture the oscillatory behavior characteristic of a Hopf bifurcation.

To fix this, you need to modify the condition for detecting the Hopf bifurcation by considering both the real and imaginary parts of the solution. Specifically, you can check for changes in the sign of the imaginary part to identify oscillatory behavior. Here's how you can update the code:

% Check for Hopf bifurcation (change in stability)

% Look for oscillatory behavior in the solution

if i > 1

    if sign(imag(Y(end, 1))) ~= sign(imag(x_result(i-1)))
        fprintf('Hopf bifurcation detected at eta = %f\n', eta);
    end
end

By comparing the signs of the imaginary parts of `Y(end, 1)` and `x_result(i-1)`, you can detect when there is a change in stability, indicative of a Hopf bifurcation. Additionally, to enhance visualization and better observe the bifurcation behavior, you might consider increasing the number of points in `eta_values` or adjusting other parameters related to your system dynamics.

Once you make these changes, rerun your code to generate an accurate Hopf bifurcation diagram that showcases the expected behavior with oscillations and dynamic changes in stability. This adjustment should address the issue of getting straight lines in your plot and provide a more informative representation of the system's behavior. Good luck!

Umar 2024-7-11

Hi kdv0,

I was experimenting with your code and was able to print Hopf Bifurcation diagram for 3D system. One potential issue in the code is the initialization of the x0, y0, and z0 variables for the initial conditions. These values are set too high, which might lead to numerical instability during the integration process. Please see attached plot along with snippet code.