Here's the code checking the eigenvalues. As I could see, eigenvalues for t are near zero but eigenvalues for k are quite large.
Should I take different perturbation for k and t?
% Define the equations as anonymous functions eq1 = @(vars, m) (sqrt(vars(1))/sqrt(1-vars(1))) + ((vars(2)*(2*(m^3) + vars(2)^3 - 3*vars(2)*((m+1)^2))^2)/(6*vars(2) - 3*((m-vars(2))^2)) + ((m-vars(2))^2)*(m+2*vars(2))*(m/(3*vars(2)) + 2/3)) / ((m^3 - (vars(2)^2)*(3*m+3) + 2*(vars(2)^3))*(m/(3*vars(2)) - (2*m+vars(2))/(3*(2*vars(2)-((m-vars(2))^2))) + 2/3) - vars(2)*((m-vars(2))^2)*(2*m+vars(2))*((2*m/3) + vars(2)/3)); eq2 = @(vars, m) vars(2) - ((sqrt(vars(1))/sqrt(1-vars(1)))*(((m/vars(2))+2)/3) - (1/3)*(2 + (m/vars(2)) + ((2*m+vars(2))/(((m-vars(2))^2)-2*vars(2))))) / ((sqrt(vars(1))/sqrt(1-vars(1)))*(2*(m^3) - 3*((1+m)^2)*vars(2) + vars(2)^3)/(3*(((m-vars(2))^2)-2*vars(2))) - ((2*m+vars(2))/3));
% Define the range of m values m_values = linspace(0, 1, 100); % 100 values of m between 0 and 1
% Preallocate arrays to store solutions k_solutions = zeros(size(m_values)); t_solutions = zeros(size(m_values)); stabilities = cell(size(m_values));
% Loop over m values and solve the system of equations for each m
for i = 1:length(m_values)
m = m_values(i); % System of equations for current m
system_of_equations = @(vars) [eq1(vars, m); eq2(vars, m)]; % Initial guess for [k, t]
initial_guess = [0.75, 1.5]; % Solve the system of equations using fsolve
options = optimoptions('fsolve', 'Display', 'off');
[solution, ~, exitflag] = fsolve(system_of_equations, initial_guess, options); % Check if fsolve converged
if exitflag > 0
% Store solutions if they meet the criteria
k_solution = solution(1);
t_solution = solution(2); if k_solution > 0.5 && k_solution < 1 && t_solution > 0
k_solutions(i) = k_solution;
t_solutions(i) = t_solution; % Compute the Jacobian matrix using finite differences
delta = 1e-6;
J = zeros(2); % Initialize Jacobian matrix % Compute f_original once outside the loop
f_original = system_of_equations(solution); % Calculate the partial derivatives for the Jacobian matrix
for j = 1:2
perturbed_solution = solution;
perturbed_solution(j) = perturbed_solution(j) + delta;
f_perturbed = system_of_equations(perturbed_solution);
J(:, j) = (f_perturbed - f_original) / delta;
end % Calculate eigenvalues of the Jacobian matrix
eigenvalues = eig(J); % Debug: print eigenvalues for inspection
fprintf('m = %.2f, k = %.4f, t = %.4f, Eigenvalues: %.6f, %.6f\n', m, k_solution, t_solution, eigenvalues(1), eigenvalues(2)); % Analyze stability
if all(real(eigenvalues) < 0)
stabilities{i} = 'Stable';
elseif all(real(eigenvalues) == 0)
% Use center manifold method for stability analysis
% Modify this part with your center manifold method implementation
stabilities{i} = 'Center manifold method';
elseif any(real(eigenvalues) > 0)
stabilities{i} = 'Unstable';
end
else
k_solutions(i) = NaN;
t_solutions(i) = NaN;
stabilities{i} = 'NaN';
end
else
k_solutions(i) = NaN;
t_solutions(i) = NaN;
stabilities{i} = 'NaN';
end
end % Display solutions and stabilities
disp('Solutions for each m:');
for i = 1:length(m_values)
fprintf('m = %.2f, k = %.4f, t = %.4f, Stability = %s\n', m_values(i), k_solutions(i), t_solutions(i), stabilities{i});
end % Plot solutions
figure;
plot(m_values, k_solutions, '-o', 'DisplayName', 'k');
hold on;
plot(m_values, t_solutions, '-o', 'DisplayName', 't');
xlabel('m');
ylabel('Values');
title('Solutions as a Function of m');
legend('k', 't');
hold off; % Plot stability
figure;
plot(m_values, strcmp(stabilities, 'Stable'), '-o', 'DisplayName', 'Stable');
hold on;
plot(m_values, strcmp(stabilities, 'Unstable'), '-o', 'DisplayName', 'Unstable');
xlabel('m');
ylabel('Stability');
title('Stability as a Function of m');
legend('Stable', 'Unstable');
hold off;
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!发生错误
由于页面发生更改，无法完成操作。请重新加载页面以查看其更新后的状态。
Translated by 
您也可以从以下列表中选择网站：
如何获得最佳网站性能
选择中国网站（中文或英文）以获得最佳网站性能。其他 MathWorks 国家/地区网站并未针对您所在位置的访问进行优化。
美洲
欧洲
亚太