MATLAB 利用YALMIP+Gurobi 求解线性规划 -多无人机扫描覆盖

使用要点

1. 创建决策变量
2. 设置目标函数
3. 添加约束条件
4. 参数配置
5. 求解问题

问题描述

假设M个无人机的任务是尽快覆盖一组由 P 顶点表示的多边形凸区域，假设每架无人机的最大飞行时间是有限的，并且是预先知道的。每架无人机的都配备了一个指向下方的机载摄像头，无人机需要使用它来感知 P 指定的整个区域。具体问题如下：

1. 最大限度减少覆盖区域 P 所需时间的无人机数量 m <= M
2. 以最短完成任务时间为目标规划每个 UAV 路径

变量定义

C i j C_{ij} Cij​ represents the traversing cost of the edge ( i , j ) (i,j) (i,j) between nodes i i i and j j j.

for  i=1:numberOfVertices
	for j=1:numberOfVertices
		if i~=j
			C(i,j)=norm(V(i,:)-V(j,:))/(uavSpeed*1000/60);
		end
	end	
end

X i j k ∈ { 0 , 1 } {X_{ij}^{k}}\in\{0,1\} Xijk​∈{0,1} represents whether or not the k k k-th UAV is going to fly from vertex i i i to vertex j j j.

V i j k ∈ R {V_{ij}^{k}}\in\R Vijk​∈R represent the flight speed of the k k k-th UAV while flying from vertex i i i to vertex j j j.

t s ∈ R {t_{s}}\in\R ts​∈R be the individual setup time.

L k {L_k} Lk​ be the battery duration of UAV number k k k.

m ∈ N {m}\in\N m∈N be the variable that represents the number of UAVs designed for a mission

M ∈ N {M}\in\N M∈N be the total number of UAVs available.

O ∈ N {O}\in\N O∈N be the constant number of UAV operators.

N ∈ N {N}\in\N N∈N be the number of nodes of the graph.

d k d_{k} dk​ represents the extra time necessary to launch UAV k k k.

每个无人机的最短时间可表示为：
T k = ∑ i = 1 N ∑ j = 1 N C i j V i j k X i j k + d k T_{k}=\sum_{i=1}^{N} \sum_{j=1}^{N} \frac{C_{i j}}{V_{i j}^{k}} X_{i j}^{k}+d_{k} Tk​=i=1∑N​j=1∑N​Vijk​Cij​​Xijk​+dk​
线性规划问题定义：
m i n ( V ) min(\mathcal{V}) min(V)
∑ i = 1 N ∑ j = 1 N C i j V i j k X i j k + d k ≤ V , k = 1 , … , M \sum_{i=1}^{N} \sum_{j=1}^{N} \frac{C_{i j}}{V_{i j}^{k}} X_{i j}^{k}+d_{k} \leq \mathcal{V}, k =1, \ldots, M i=1∑N​j=1∑N​Vijk​Cij​​Xijk​+dk​≤V,k=1,…,M
t s ⌈ k O ⌉ ∑ j = 1 N X 1 j k = d k , k = 1 , … , M t_{s} {\lceil\frac{k}{O}\rceil \sum_{j=1}^{N} X_{1 j}^{k}=d_{k}, k=1, \ldots, M } ts​⌈Ok​⌉j=1∑N​X1jk​=dk​,k=1,…,M

01 | 创建决策变量

变量设置常用的三种形式：

• sdqvar() 设置实型变量
• intvar() 设置整型变量
• binvar() 设置0-1变量

举例：P = sdqvar(n,m) 表示用 n 行 m 列定义一个矩阵（或标量）P

v = sdpvar(1,1);
X = binvar(numberOfVertices,numberOfVertices,uavNumber,'full');
u = sdpvar(numberOfVertices,1);
m = intvar(1,1);
for k = 1:uavNumber
	for i = 1:numberOfVertices
		X(i,i,k) = 0;
	end
	if uav(k) == 0
		X(:,:,k) = zeros(numberOfVertices);
	end
end

02 | 设置目标函数

YALMIP 默认求解最小值问题，如果问题的目标函数是求解最大值，则需要在目标函数前加负号
目标函数即求最小值 v

03 | 添加约束条件

先初始化约束条件为空，再逐个将约束条件添加

constraints = [];
% 限制1 - 总成本
for k = 1:uavNumber
    sum1 = 0;
    for i = 1:numberOfVertices
        for j = 1:numberOfVertices
            if i ~= j
                sum1 = sum1 + C(i,j)*sum(X(i,j,k));
            end
        end
    end
    constraints = [constraints, v >= sum1 + sum(X(1,:,k))*ceil(k/O)*uavSetupTime];
    % 限制9 - 飞行时间
    constraints = [constraints, sum1 <= uavFlightTime];
end
% 限制2 - 每个顶点必须被访问一次
for j = 2:numberOfVertices
    constraints = [constraints, sum(sum(X(:,j,:))) == 1];
end
% 限制3 - 到达给定节点的节点与离开该节点的节点相同。
for k = 1:uavNumber
    for p = 1:numberOfVertices
        constraints = [constraints, sum(X(:,p,k)) - sum(X(p,:,k)) == 0];
    end
end
% 限制4&5 - 无人机数量
constraints = [constraints, sum(sum(X(1,:,:))) == m];
constraints = [constraints, m <= sum(uav)];
% 限制7 - 循环
for i = 2:numberOfVertices
    for j = 2:numberOfVertices
        if i ~= j
            constraints = [constraints, u(i) - u(j) + numberOfVertices*sum(X(i,j,:)) <= numberOfVertices - 1];
        end
    end
end
% 限制 8 - 强制每条线路仅由一架无人机沿一个方向穿越
for i = 2:2:numberOfVertices
   constraints = [constraints, 1 == sum(X(i,i+1,:)) + sum(X(i+1,i,:))];
end
% 限制 10 - 避免 UA Vs 沿着不平行于覆盖范围的边缘穿过覆盖区域行
for i = 2:2:numberOfVertices
    constraints = [constraints, sum(X(i,i+1,:)) == sum(sum(X(i,3:2:numberOfVertices,:)))];
    constraints = [constraints, sum(X(i+1,i,:)) == sum(sum(X(i+1,2:2:numberOfVertices,:)))];
end

04 | 参数设置

使用 sdpsetting() 函数进行参数设置

% 使用 gurobi 最小化 v 受限于约束变量中包含的约束
options = sdpsettings('verbose',1,'gurobi.Threads',4);

05 | 求解问题

首先使用 solvesdp() 函数求解
再使用 double() 函数将求解的变量转换为实数

solvesdp(constraints,v,options);
X = round(double(X));