Решение уравнений с ограничениями неравенства

fsolve решает систему нелинейных уравнений. Однако она не позволяет включать какие-либо ограничения, даже связанные ограничения. Как решить систему нелинейных уравнений при наличии ограничений?

Решение, удовлетворяющее вашим ограничениям, не гарантировано. На самом деле, проблема может не иметь никакого решения, даже такого, которое не удовлетворяет вашим ограничениям. Однако существуют методы поиска решений, удовлетворяющих вашим ограничениям.

Чтобы проиллюстрировать методы, рассмотрим, как решить уравнения

где компоненты $$x$$ должны быть неотрицательными. Уравнения имеют четыре решения:

Только одно решение удовлетворяет ограничениям, а именно $$x = (\mathrm{}\mathrm{}\mathrm{}\mathrm{10,20}$$).

fbnd вспомогательная функция в конце этого примера вычисляет $$F(x$$) численно.

Использовать разные начальные точки

Как правило, система $$N$$ уравнений в $$N$$ переменных имеет изолированные решения, что означает, что каждое решение не имеет соседних соседей, которые также являются решениями. Таким образом, одним из способов поиска решения, удовлетворяющего некоторым ограничениям, является создание нескольких начальных точек. x0, а затем выполнить fsolve начиная с каждого x0.

Для этого примера, чтобы найти решение для системы уравнений $$F(x)\mathrm{}$$= 0, возьмем 10 случайных точек, которые обычно распределены со средним значением 0 и стандартным отклонением 100.

rng default % For reproducibility
N = 10; % Try 10 random start points
pts = 100*randn(N,2); % Initial points are rows in pts
soln = zeros(N,2); % Allocate solution
opts = optimoptions('fsolve','Display','off');
for k = 1:N
    soln(k,:) = fsolve(@fbnd,pts(k,:),opts); % Find solutions
end

Перечислите решения, удовлетворяющие ограничениям.

idx = soln(:,1) >= 0 & soln(:,2) >= 0;
disp(soln(idx,:))

   10.0000   20.0000
   10.0000   20.0000
   10.0000   20.0000
   10.0000   20.0000
   10.0000   20.0000

Использовать различные алгоритмы

fsolve имеет три алгоритма. Каждый из них может привести к различным решениям.

Для этого примера возьмем x0 = [1,9] и проверьте решение, которое возвращает каждый алгоритм.

x0 = [1,9];
opts = optimoptions(@fsolve,'Display','off',...
    'Algorithm','trust-region-dogleg');
x1 = fsolve(@fbnd,x0,opts)

x1 = 1×2

   -1.0000   -2.0000

opts.Algorithm = 'trust-region';
x2 = fsolve(@fbnd,x0,opts)

x2 = 1×2

   -1.0000   20.0000

opts.Algorithm = 'levenberg-marquardt';
x3 = fsolve(@fbnd,x0,opts)

x3 = 1×2

    0.9523    8.9941

Здесь все три алгоритма находят разные решения для одной и той же начальной точки. Ни один из них не удовлетворяет ограничениям. Сообщаемое «решение» x3 это даже не решение, а просто локальная стационарная точка.

Использовать lsqnonlin с границами

lsqnonlin пытается минимизировать сумму квадратов компонентов в векторной функции $$F(x$$). Поэтому он пытается решить уравнение $$F( x) \mathrm{=}$$0. Также,lsqnonlin принимает связанные ограничения.

Сформулировать пример проблемы для lsqnonlin и решить его.

lb = [0,0];
rng default
x0 = 100*randn(2,1);
[x,res] = lsqnonlin(@fbnd,x0,lb)

Local minimum found.

Optimization completed because the size of the gradient is less than
the value of the optimality tolerance.

x = 2×1

   10.0000
   20.0000

res = 2.4783e-25

В этом случае lsqnonlin сходится к решению, удовлетворяющему ограничениям. Вы можете использовать lsqnonlin с помощью инструментария глобальной оптимизации MultiStart для автоматического поиска по многим начальным точкам. См. раздел MultiStart Using lsqcurvefit или lsqnonlin (Global Optimization Toolbox).

Задать уравнения и неравенства как fmincon Ограничения

Можно переформулировать проблему и использовать fmincon следующим образом:

fminconstr вспомогательная функция в конце этого примера реализует нелинейные ограничения. Решите проблему с ограничениями.

lb = [0,0]; % Lower bound constraint
rng default % Reproducible initial point
x0 = 100*randn(2,1);
opts = optimoptions(@fmincon,'Algorithm','interior-point','Display','off');
x = fmincon(@(x)0,x0,[],[],[],[],lb,[],@fminconstr,opts)

x = 2×1

   10.0000
   20.0000

В этом случае fmincon решает проблему из начальной точки.

Вспомогательные функции

Этот код создает fbnd функция помощника.

function F = fbnd(x)

F(1) = (x(1)+1)*(10-x(1))*(1+x(2)^2)/(1+x(2)^2+x(2));
F(2) = (x(2)+2)*(20-x(2))*(1+x(1)^2)/(1+x(1)^2+x(1));
end

Этот код создает fminconstr функция помощника.

function [c,ceq] = fminconstr(x)

c = []; % No nonlinear inequality
ceq = fbnd(x); % fsolve objective is fmincon nonlinear equality constraints
end