在 matlab 中进行优化：找出约束方程

Question

我有 N 行由 y 截距和角度 q 定义。约束是所有 N 条线必须在一个点相交。我可以想出最终得到约束的方程式如下：

Y = tan(q(1))X + y(1)
Y = tan(q(2))X + y(2)
...

如果 N = 3 或 4，我可以手动获得约束，但如果 N 大于 4，我只能获得一个约束。如果 N = 3 或 4，那么当我求解上述 X 方程时，我得到 2 个方程，然后可以将它们设置为彼此相等。如果 N > 4，我得到超过 2 个等于 X 的方程，我不知道如何将它们压缩为一个约束。如果我不能将它们压缩为一个约束并且能够使用动态创建的多个约束（取决于传入的 N）来解决优化问题，那也很好。

为了更好地理解我在做什么，我将展示如何获得 N = 3 的约束。我从以下三个等式开始：

Y = tan(q(1))X + y(1)
Y = tan(q(2))X + y(2)
Y = tan(q(3))X + y(3)

然后我将它们设置为彼此相等并得到这些方程：

tan(q(1))X + y(1) = tan(q(2))X + y(2)
tan(q(2))X + y(2) = tan(q(3))X + y(3)

然后我求解 X 并得到这个约束：

(y(2) - y(1)) / (tan(q(1)) - tan(q(2))) = (y(3) - y(2)) / (tan(q(2)) - tan(q(3)))

注意我有 2 个方程来求解 X。当 N > 4 时，我最终得到的结果超过 2。如果我能够动态创建约束，然后在 MATLAB 中调用一个优化函数来处理多个限制，但到目前为止还没有找到。

Answer 1

您说优化算法需要调整q，以使“真实”问题最小化，同时上述等式也成立。

请注意，fifth Euclid axoim 确保所有线始终与所有其他线相交，除非两个 qs 相等但相应的 y0s 不相等。最后一种情况非常罕见（在浮点上下文中），我将在此处跳过它，但为了增加稳健性，您最终应该包含它。

现在，首先，从矩阵的角度来思考。您的约束可以通过矩阵方程表示：

y = tan(q)*x + y0

其中q、y 和y0 是[Nx1] 矩阵，x 是未知标量。请注意y = c*ones(N,1)，例如，仅包含相同常数的矩阵。这实际上是一个非线性约束——也就是说，它不能表示为

A*q <= b   or   A*q == b

带有A 一些设计矩阵和b 一些解决方案向量。因此，您必须编写一个定义此非线性约束的函数，您可以将其传递给像fmincon 这样的优化器。来自the documentation：

x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon) 最小化非线性不等式 c(x) 或等式 ceq(x) 在 nonlcon 中定义。 fmincon 优化使得 c(x) ≤ 0 和 ceq(x) = 0. 如果不存在边界，则设置 lb = [] 和/或 ub = []。

请注意，您实际上是在朝着正确的方向前进。您可以使用以下等式求解任意一对线 q(n),y0(n) 和 q(m),y0(m) 的交点的 x 位置：

x(n,m) = (y0(n)-y0(m)) / (q(m)-q(n))

您的nonlcon 函数应该为所有 个可能的对n,m 找到x，并检查它们是否所有 相等。你可以像这样方便地做到这一点：

function [c, ceq] = nonlcon(q, y0)
    % not using inequalities
    c = -1; % NOTE: setting it like this will always satisfy this constraint

    % compute tangents 
    tanq = tan(q);

    % compute solutions to x for all pairs 
    x = bsxfun(@minus, y0, y0.') ./ -bsxfun(@minus, tanq, tanq.');

    % equality contraints: they all need to be equal 
    ceq = diff(x(~isnan(x))); % NOTE: if all(ceq==0), converged.

end

请注意，您实际上并没有明确地求解 q（或者根本不需要交叉点的 y 坐标）——这就是 fmincon 的全部工作。

你需要做一些实验，因为有时定义就足够了

x = x(~isnan(x));
ceq = norm(x-x(1)); % e.g., only 1 equality constraint

这会更快（计算的导数更少），但确实需要其他问题

x = x(~isnan(x));
ceq = x-x(1); % e.g., N constraints

或类似的技巧。这实际上取决于问题的其余部分，优化器将找到每种情况的难度。