【三维装箱】基于matlab求解三维装箱优化问题【含Matlab源码 1194期】

一、获取代码方式

获取代码方式1：
通过订阅紫极神光博客付费专栏，凭支付凭证，私信博主，可获得此代码。

获取代码方式2：
通过紫极神光博客主页开通CSDN会员，凭支付凭证，私信博主，可获得此代码。

获取代码方式3：
完整代码已上传我的资源：【三维装箱】基于matlab求解三维装箱优化问题【含Matlab源码 1194期】

备注：开通CSDN会员，仅只能免费获得1份代码（有效期为开通日起，三天内有效）；
订阅紫极神光博客付费专栏，可免费获得1份代码（有效期为订阅日起，三天内有效）；

二、三维装箱简介

1 前言
三维装箱问题（Three Dimensional Container Loading Poblem）指在满足容积限制、外形几何限制和稳定性限制等条件的情况下，把一定数量体积较小的物品放入具有较大容量的一个或多个箱子，达到所用箱子数量最少、空间利用率最高、稳定性最好、装载价值最高、容重比最高等目的组合优化的问题。三维装箱问题广泛存在于运输和包装等行业的物资装载过程中，因此，如何针对实际应用需求设计具有较高操作性和装载效率的三维装箱求解算法，对于提高社会经济效益具有积极作用。

自20世纪60年代开始，三维装箱问题便因切割库存问题而提出，引起了学者关注，并被证明是具有较高复杂度的NP难问题。经过多年发展，三维装箱问题的求解不仅在枚举法、分支定界法等确定性算法上取得了进步。随着元启发式算法、遗传算法、粒子群算法等现代启发式算法的崛起，三维装箱问题的研究更加深入，逐步出现了带有多种现实约束和求解目标的复杂算法。其中，近年来具有代表性的相关研究如下所述。

Lodi等[9]受到一维装箱问题求解算法的启发，针对三维装箱问题提出了基于“先高后面积”（Height First-Area Second）的分层策略的结构式禁忌搜索算法，实验证明算法在分层堆放的情况下能得到较好的结果。Crainic等考虑了装载盒子时箱子剩余空间的描述和备选定位点的重要作用，为寻找每次装载时的合适空间位置和最优定位点，提出了极点（Extreme Points,EPS）的概念和计算方法，并设计了基于极点的三维装箱启发式算法。Karabulut[12]和Kang[13]等先后以规则长方体箱装载最大量规则长方体盒为研究目标，在设计深度、底部和左部方向优先装载（Deepest Bottom Left with Fill）物资装载策略的基础上，提出了基于遗传算法的三维装箱策略求解算法，2个研究的区别在于后者在前者的基础上，对装载剩余空间进行了更为细致的描述。Ramos等重点考虑了运输行业货物装载时静止时的纵向稳定性和行进时的横向稳定性，并提出了具有较高操作性的稳定性评价指标和基于顺序的货物装载启发式算法。

尽管三维装箱问题研究的深度和广度在不断延伸，相关研究的模型建立、算法设计和实际应用也取得了众多成果，但是一些现实约束仍难以满足，解决方法也存在较大的局限性。其中，现有求解三维装箱问题较为常用的遗传算法，在交叉算子设计、变异算子设计和进化性能等方面表现不足。

1 问题描述
1.1 基本问题描述
除球体、非规则立方体等形状盒子的装箱特例外，经典的三维装箱问题（见图1）可描述为：用一定数量的长宽高分别为Li,Wi和Hi的箱子（Container），按照一定顺序和装载规则根据不同转向逐个装入一定数量的长、宽、高分别为lj,wj和hj的盒子（Box）中，实现箱子容量的利用率最大、盒子堆放的稳定性最好、装载价值最高或容重比最恰当等其一或多种目标，并满足外形约束、容积约束、稳定性约束、盒子属性相关性约束、载重约束等实际限制。

图1 三维装箱问题示意
根据文献对三维装箱问题的分类法，该研究以SBSBPP类（Single Bin-Size Bin Packing Problem，有限数量单类型箱体三维装载问题）为研究对象，研究如何在单个规则长方体箱子中装载尽量多的盒子（这些盒子具有相似度较小且分散较大的特点），以最大限度地提高箱子的空间利用率。此研究旨在解决每个盒子应以何种顺序、何种方向放置在具体空间位置的问题，即提供具体的盒子装载方案。此过程需满足如下现实约束。
1）外形约束。进行装载时，所有盒子的几何顶点不能位于箱子之外，装进箱子的盒子在空间上不能发生重叠。
2）容积约束。装载盒子的总体积不能超过箱子的容积。
3）稳定性约束。装载盒子的堆放需满足重心要求，不能出现悬空或重心偏移的情况。
基本假设：盒子质量均匀分布，其质心位于其几何中心；不考虑盒子总质量超载的问题；盒子装载时应与箱子的边长方向平行（即不斜放）。

1.2 箱体最大剩余空间描述
最大剩余空间（Allowable Packing Area，简写为APA）指进行盒子装载前的最大长方体形闲置空间，每个最大剩余空间可用最靠近原点的长方体顶点三维坐标和最远离原点的长方体顶点三维坐标合并表示，全部的最大剩余空间即构成最大剩余空间集（简写为APAs）。盒子的装载是在箱子的最大剩余空间中按照一定规则进行的，最大剩余空间的刻画关系到盒子的放置和空间容纳能力，对于盒子的装载至关重要。该研究采用文献中的I-DBLF算法，完成最大剩余空间的描述。每装入一个箱子，其最大剩余空间将按照盒子的边界和x,y,z轴进行划分和更新。如图2a所示，空置箱子的最大剩余空间只有1个，表示为(0,0,0),(Li,Wi,Hi)；在临原点位置装入一个长、宽、高分别为lj,wj和hj的盒子后，最大剩余空间更新为3个（见图2b—d），最大剩余空间集为{((Li-lj,0,0),(Li,Wi,Hi)),((0,0,Hi-hj),(Li,Wi,Hi)),((0,Wi-wj,0),(Li,Wi,Hi))}。最大剩余空间的计算过程较为复杂，因篇幅所限，此处不做赘述，具体可见文献。

图2 剩余空间示意

三、部分源代码

close all;
clear;
clc;

%% set data
Dbox = [10, 10, 10];
Dobj = [1, 1.5, 1;
        1, 1.5, 1];

PobjI = [7, 5, 5, 0, 0, 0;
         3, 5, 5, 0, 0, 0];

PobjD = [3, 3, -3, 0, 0, 0;
         8, 8, 13, 0, 0, 0];

objN = size(Dobj, 1);
Rcost = 10;



%% shape & bounding box & start point



[x, f] = genFunction(Rcost, PobjD);

Ar = [1 0 0; -1 0 0; 0 1 0; 0 -1 0; 0 0 1; 0 0 -1];
syms vx vy vz;
Vcoor = [vx; vy; vz];
for i=1:objN
    br = [Dobj(i,1); Dobj(i,1); Dobj(i,2); Dobj(i,2); Dobj(i,3); Dobj(i,3)];
    Gc{i} = Ar * Vcoor - br;
end

G = zeros(0);
for i=1:objN
    G = [G; genGi(i, Dbox, Dobj, x)];
end

%%
tic;
[X0,~] = zoutendijk1(f, x, G, Gc, Dbox, Dobj, PobjI, PobjD);
toc
figure;
plotResult(X0, Dbox, Dobj, PobjI, PobjD);

 












function [ f ] = genCollision( Gc, X, x0, Dobj, o1, o2 )

syms vx vy vz;
x = [vx vy vz];

Gc1 = Gc{o1};
Gc2 = Gc{o2};
a1 = Dobj(o1,1);
b1 = Dobj(o1,2);
c1 = Dobj(o1,3);
a2 = Dobj(o2,1);
b2 = Dobj(o2,2);
c2 = Dobj(o2,3);

R1V = [a1 b1 c1; -a1 -b1 c1; -a1 b1 c1; a1 -b1 c1;
       a1 b1 -c1; -a1 -b1 -c1; -a1 b1 -c1; a1 -b1 -c1];
 
R2V = [a2 b2 c2; -a2 -b2 c2; -a2 b2 c2; a2 -b2 c2;
       a2 b2 -c2; -a2 -b2 -c2; -a2 b2 -c2; a2 -b2 -c2];
 
% R1V(1,:) = [-a1, -b1];
% R1V(2,:) = [a1, b1];
% R1V(3,:) = [-a1, b1];
% R1V(4,:) = [a1, -b1];
% 
% R2V(1,:) = [-a2, -b2];
% R2V(2,:) = [a2, b2];
% R2V(3,:) = [-a2, b2];
% R2V(4,:) = [a2, -b2];

newR1V = R1V;
newR2V = R2V;
for i=1:size(R1V,1)
    newR1V(i,:) = rt(R1V(i,:), transpose(x0(6*(o1-1)+1:6*(o1-1)+3)), transpose(x0(6*(o2-1)+1:6*(o2-1)+3)), x0(6*(o1-1)+4:6*(o1-1)+6), x0(6*(o2-1)+4:6*(o2-1)+6));
    newR2V(i,:) = rt(R2V(i,:), transpose(x0(6*(o2-1)+1:6*(o2-1)+3)), transpose(x0(6*(o1-1)+1:6*(o1-1)+3)), x0(6*(o2-1)+4:6*(o2-1)+6), x0(6*(o1-1)+4:6*(o1-1)+6));
end

temp = zeros(size(newR2V,1),3);
minEdged = zeros(0,3);
for i=1:size(Gc1,1)
    for j=1:size(newR2V,1)
        temp(j,:) = [subs(Gc1(i), x, newR2V(j,:)), i, j];
    end
    [~, vec] = sort(temp);
    temp1 = temp(vec(:,1),:);
    for k=1:size(newR2V,1)
        if abs(temp1(k,1)-temp1(k+1,1)) < 0.1
            continue
        else
            break
        end
    end
    minEdged = [minEdged; temp1(1:k,:)];
end
[~, vec] = sort(minEdged, 'descend');
temp1 = minEdged(vec(:,1),:);
for k=1:size(temp1,1)
    if abs(temp1(k,1)-temp1(k+1,1)) < 0.1
        continue
    else
        break
    end
end
maxEdged1 = temp1(1:k,:);

minEdged = zeros(0,3);        % value % edge num % v num
for i=1:size(Gc1,1)
    for j=1:size(newR2V,1)
        temp(j,:) = [subs(Gc2(i), x, newR1V(j,:)), i, j];
    end
    [~, vec] = sort(temp);
    temp1 = temp(vec(:,1),:);
    for k=1:size(temp1,1)
        if abs(temp1(k,1)-temp1(k+1,1)) < 0.1
            continue
        else
            break
        end
    end
    minEdged = [minEdged; temp1(1:k,:)];   
end       
[~, vec] = sort(minEdged, 'descend');
temp1 = minEdged(vec(:,1),:);
for k=1:size(temp1,1)
    if abs(temp1(k,1)-temp1(k+1,1)) < 0.1
        continue
    else
        break
    end
end
maxEdged2 = temp1(1:k,:);
      

f = zeros(0,0);
if abs(maxEdged1(1,1)-maxEdged2(1,1)) < 0.1
    for i=1:size(maxEdged1,1)
        f = [f; -subs(Gc1(maxEdged1(i,2)), x, rt(R2V(maxEdged1(i,3),:), X(6*(o2-1)+1:6*(o2-1)+3), X(6*(o1-1)+1:6*(o1-1)+3), X(6*(o2-1)+4:6*(o2-1)+6), X(6*(o1-1)+4:6*(o1-1)+6)))];
    end
    for i=1:size(maxEdged2,1)
        f = [f; -subs(Gc2(maxEdged2(i,2)), x, rt(R1V(maxEdged2(i,3),:), X(6*(o1-1)+1:6*(o1-1)+3), X(6*(o2-1)+1:6*(o2-1)+3), X(6*(o1-1)+4:6*(o1-1)+6), X(6*(o2-1)+4:6*(o2-1)+6)))];
    end
elseif (maxEdged1(1,1) < maxEdged2(1,1))
    for i=1:size(maxEdged2,1)
        f = [f; -subs(Gc2(maxEdged2(i,2)), x, rt(R1V(maxEdged2(i,3),:), X(6*(o1-1)+1:6*(o1-1)+3), X(6*(o2-1)+1:6*(o2-1)+3), X(6*(o1-1)+4:6*(o1-1)+6), X(6*(o2-1)+4:6*(o2-1)+6)))];
    end
else
    for i=1:size(maxEdged1,1)
        f = [f; -subs(Gc1(maxEdged1(i,2)), x, rt(R2V(maxEdged1(i,3),:), X(6*(o2-1)+1:6*(o2-1)+3), X(6*(o1-1)+1:6*(o1-1)+3), X(6*(o2-1)+4:6*(o2-1)+6), X(6*(o1-1)+4:6*(o1-1)+6)))];
    end
end

return

end


  

 

  
1
  
2
  
3
  
4
  
5
  
6
  
7
  
8
  
9
  
10
  
11
  
12
  
13
  
14
  
15
  
16
  
17
  
18
  
19
  
20
  
21
  
22
  
23
  
24
  
25
  
26
  
27
  
28
  
29
  
30
  
31
  
32
  
33
  
34
  
35
  
36
  
37
  
38
  
39
  
40
  
41
  
42
  
43
  
44
  
45
  
46
  
47
  
48
  
49
  
50
  
51
  
52
  
53
  
54
  
55
  
56
  
57
  
58
  
59
  
60
  
61
  
62
  
63
  
64
  
65
  
66
  
67
  
68
  
69
  
70
  
71
  
72
  
73
  
74
  
75
  
76
  
77
  
78
  
79
  
80
  
81
  
82
  
83
  
84
  
85
  
86
  
87
  
88
  
89
  
90
  
91
  
92
  
93
  
94
  
95
  
96
  
97
  
98
  
99
  
100
  
101
  
102
  
103
  
104
  
105
  
106
  
107
  
108
  
109
  
110
  
111
  
112
  
113
  
114
  
115
  
116
  
117
  
118
  
119
  
120
  
121
  
122
  
123
  
124
  
125
  
126
  
127
  
128
  
129
  
130
  
131
  
132
  
133
  
134
  
135
  
136
  
137
  
138
  
139
  
140
  
141
  
142
  
143
  
144
  
145
  
146
  
147
  
148
  
149
  
150
  
151
  
152
  
153
  
154
  
155
  
156
  
157
  
158
  
159
  
160
  
161
  
162
  
163
  
164
  
165
  
166
  
167
  
168
  
169
  
170
  
171
  
172
  
173
  
174
 

四、运行结果

五、matlab版本及参考文献

1 matlab版本
2014a

2 参考文献
[1] 包子阳,余继周,杨杉.智能优化算法及其MATLAB实例（第2版）[M].电子工业出版社，2016.
[2]张岩,吴水根.MATLAB优化算法源代码[M].清华大学出版社，2017.
[3]周品.MATLAB 神经网络设计与应用[M].清华大学出版社，2013.
[4]陈明.MATLAB神经网络原理与实例精解[M].清华大学出版社，2013.
[5]方清城.MATLAB R2016a神经网络设计与应用28个案例分析[M].清华大学出版社，2018.

文章来源: qq912100926.blog.csdn.net，作者：海神之光，版权归原作者所有，如需转载，请联系作者。

原文链接：qq912100926.blog.csdn.net/article/details/119493181