如何在Matlab中结合微分求解非线性方程组（fsolve应用）

原代码核心问题是符号运算与数值求解不兼容：fsolve要求目标函数接收数值输入并输出数值，但代码中混用了符号函数diff与数值变量，且主脚本错误传递符号变量给fsolve，导致类型冲突。同时手动求导无法适配100个方程的规模，必须实现自动求导。

方法一：符号工具箱自动生成数值函数（推荐，精度高、效率优） #

先通过符号运算推导原函数加导数的表达式，再转换为可直接用于fsolve的数值匿名函数，彻底规避符号与数值的混用问题。

% 主脚本
syms x(1:2) % 扩展到100个变量时，改为syms x(1:100)
x_sym = x;

% 定义原符号函数F（扩展时需通用化表达式）
denom_sym = 1 + sum(exp(-x_sym + [2, 1])); % [2,1]为对应每个方程的常数项，扩展时替换为长度100的常数向量
F_sym = exp(-x_sym + [2, 1])./denom_sym;

% 计算每个F(i)对x(i)的导数（仅取雅可比矩阵对角元）
dF_dx = arrayfun(@(i) diff(F_sym(i), x_sym(i)), 1:length(x_sym));

% 构造新方程组：原函数+对应导数
F_new_sym = F_sym + dF_dx;

% 将符号表达式转换为数值匿名函数
fun = matlabFunction(F_new_sym, 'Vars', {x_sym});

% 调用fsolve求解
x0 = [0, 0]; % 扩展时对应100维初始向量
x_sol = fsolve(fun, x0);
disp('求解结果：');
disp(x_sol);
disp('验证新方程组的值：');
disp(fun(x_sol));

方法二：数值差分计算导数（无需符号工具箱） #

通过对每个变量施加微小扰动，用数值方法近似导数，适合无法使用符号工具箱的场景。

% 主脚本
fun = @root2d_with_num_deriv;
x0 = [0, 0];
x_sol = fsolve(fun, x0);
disp('求解结果：');
disp(x_sol);
disp('验证新方程组的值：');
disp(fun(x_sol));

% 带数值导数的目标函数
function F = root2d_with_num_deriv(x)
    n = length(x);
    c = [2, 1]; % 对应每个方程的常数项，扩展时替换为长度100的常数向量
    h = 1e-8; % 差分步长，平衡精度与数值稳定性
    
    % 计算原函数值（向量运算，扩展时无需修改）
    denom = 1 + sum(exp(-x + c));
    F = exp(-x + c)./denom;
    
    % 计算每个F(i)对x(i)的数值导数
    for i = 1:n
        x_pert = x;
        x_pert(i) = x_pert(i) + h;
        denom_pert = 1 + sum(exp(-x_pert + c));
        F_pert_i = exp(-x_pert(i) + c(i))/denom_pert;
        dF_dxi = (F_pert_i - F(i))/h;
        F(i) = F(i) + dF_dxi;
    end
end

扩展到100个方程的关键细节 #

1. 方法一扩展：只需将syms x(1:2)改为syms x(1:100)，并将常数向量[2,1]替换为长度100的对应常数数组，其余代码无需修改。
2. 方法二扩展：仅需更新常数向量c为100维数组，向量运算部分天然支持大规模变量，无需额外修改循环逻辑。

内容的提问来源于stack exchange，提问作者Wengritzky