matlab：涉及复杂函数图像的交点求解

在MATLAB中求解两个图像的交点是一个常见的需求。本文将通过一个示例，展示如何求解两个图像的交点，并提供相应的MATLAB代码。

画出图像

首先，我们需要绘制两个图像，以便直观地看到它们的交点。以下是绘制图像的MATLAB代码：

% 定义符号变量
syms x1 x2;% 上边界方程
eq1 = 10 + 110 * (0.8 + 0.05 + 0.4 * sin(4 * atan2(x2, x1))^16)^2 - (x1 + x2) == 0;% 绘制图形
figure;% 使用 fimplicit 绘制上边界
fimplicit(@(x1, x2) 10 + 110 * (b + 0.05 + 0.4 * sin(4 * atan2(x2, x1))^16)^2 - (x1 + x2), [0, 100, 0, 100]);
hold on;% 使用 fimplicit 绘制 y = 100 - x1
fimplicit(@(x1, x2) x1 + x2 - 100, [0, 100, 0, 100]);% 设置图例和标题
legend('Upper Boundary', 'y = 100 - x1');
title('Plot of Equations');
xlabel('x1');
ylabel('x2');
hold off;

绘制的图像如下所示：
在这里插入图片描述
从图中可以看出，存在四个实数域交点。

求解交点方法1：solve不加求解范围（失败）

首先尝试使用solve函数来求解交点：

% 定义符号变量
syms x1 x2;% 定义参数
b = 0.8;
l = atan2(x2, x1);% 定义方程组
eq1 = 10 + 110 * (b + 0.05 + 0.4 * sin(4 * l)^16)^2 - (x1 + x2);
eq2 = x1 + x2 - 100;% 解方程组
[sol_x1, sol_x2] = solve([eq1 == 0, eq2 == 0], [x1, x2]);% 转换为数值解
sol_x1 = double(sol_x1);
sol_x2 = double(sol_x2);% 筛选实数解
real_solutions = [sol_x1, sol_x2];
real_solutions = real_solutions(imag(real_solutions(:, 1)) == 0 & imag(real_solutions(:, 2)) == 0, :);% 输出实数解
disp('Real solutions (x1, x2):');
disp(real_solutions);

输出为空值，这表明在解的过程中遇到了问题。

求解交点方法2：solve加求解范围（成功）

% 定义符号变量
syms x1 x2;% 求解范围
assume(x1>=0&x1<=100)
assume(x2>=0&x2<=100)% 定义参数
b = 0.8;
l = atan2(x2, x1);% 定义方程组
eq1 = 10 + 110 * (b + 0.05 + 0.4 * sin(4 * l)^16)^2 - (x1 + x2);
eq2 = x1 + x2 - 100;% 解方程组
[sol_x1, sol_x2] = solve([eq1 == 0, eq2 == 0], [x1, x2]);% 转换为数值解
sol_x1 = double(sol_x1);
sol_x2 = double(sol_x2);% 筛选实数解
real_solutions = [sol_x1, sol_x2];
real_solutions = real_solutions(imag(real_solutions(:, 1)) == 0 & imag(real_solutions(:, 2)) == 0, :);% 输出实数解
disp('Real solutions (x1, x2):');
disp(real_solutions);

输出：
在这里插入图片描述

求解交点方法3：fsolve（成功）

下面通过数值方法而不是符号方法来找到解，使用 fsolve（数值求解函数）：

% 定义匿名函数
func = @(x) [10 + 110 * (b + 0.05 + 0.4 * sin(4 * atan2(x(2), x(1)))^16)^2 - (x(1) + x(2)), x(1) + x(2) - 100];% 设置选项以使用较大的初始搜索范围
options = optimoptions('fsolve', 'Display', 'off', 'MaxFunctionEvaluations', 6000, 'MaxIterations', 4000);% 存储解
solutions = [];% 尝试多个随机初始猜测
for i = 1:100initial_guess = rand(1, 2) * 100; % 生成0到100之间的随机初始猜测[sol, fval, exitflag, output] = fsolve(func, initial_guess, options);% 只有当fsolve成功收敛时才记录解if exitflag > 0 && all(abs(fval) < 1e-6)solutions = [solutions; sol];end
end% 去除重复的解，考虑数值误差
solutions = round(solutions, 3); % 四舍五入到三位小数
solutions = unique(solutions, 'rows', 'stable');% 过滤掉不在感兴趣区域的解
solutions = solutions(all(solutions >= 0 & solutions <= 100, 2), :);% 输出数值解
disp('Numerical solutions (x1, x2):');
disp(solutions);