基于锥规划的分段线性质量-弹簧系统能量最小化方法gydF4y2Ba

打开实时脚本gydF4y2Ba

这个例子展示了如何找到悬挂在两个锚点上的质量-弹簧系统的平衡位置。弹簧具有分段线性拉伸力。该系统包括gydF4y2Ba $ngydF4y2Ba$ 二维的质量。质量gydF4y2Ba $我gydF4y2Ba$ 连接到弹簧gydF4y2Ba $我gydF4y2Ba$ 而且gydF4y2Ba $我gydF4y2Ba +gydF4y2Ba 1gydF4y2Ba$ ．弹簧gydF4y2Ba $1gydF4y2Ba$ 而且gydF4y2Ba $ngydF4y2Ba +gydF4y2Ba 1gydF4y2Ba$ 也连接到单独的锚点。在这种情况下，弹簧的零力长度gydF4y2Ba $我gydF4y2Ba$ 长度是正的gydF4y2Ba $lgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba$ ，弹簧产生力gydF4y2Ba $kgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba 问gydF4y2Ba$ 当拉伸到一定长度时gydF4y2Ba $问gydF4y2Ba +gydF4y2Ba lgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba$ ．问题是找到质量的最小势能结构，其中势能来自重力和非线性弹簧的拉伸。平衡发生在最小能量结构。gydF4y2Ba

这幅图显示了五个弹簧和四个物体悬挂在两个锚点上。gydF4y2Ba

质量的势能gydF4y2Ba $米gydF4y2Ba$ 在高度gydF4y2Ba $hgydF4y2Ba$ 是gydF4y2Ba $米gydF4y2Ba ggydF4y2Ba hgydF4y2Ba$ ,在那里gydF4y2Ba $ggydF4y2Ba$ 是地球上的引力常数。同样，具有弹簧常数的理想线性弹簧的势能gydF4y2Ba $kgydF4y2Ba$ 拉长gydF4y2Ba $问gydF4y2Ba$ 是gydF4y2Ba $kgydF4y2Ba {问gydF4y2Ba}^{2gydF4y2Ba} /gydF4y2Ba 2gydF4y2Ba$ ．目前的模型是弹簧不理想，但有一个非零的静息长度gydF4y2Ba $lgydF4y2Ba$ ．gydF4y2Ba

这个例子的数学基础来自Lobo, Vandenberg, Boyd和LebretgydF4y2Ba[1]gydF4y2Ba．有关此示例的基于问题的版本，请参见gydF4y2Ba基于问题的圆锥规划最小化分段线性质量-弹簧系统能量gydF4y2Ba．gydF4y2Ba

数学公式gydF4y2Ba

质量的位置gydF4y2Ba $我gydF4y2Ba$ 是gydF4y2Ba $xgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba$ ，为横坐标gydF4y2Ba ${xgydF4y2Ba}_{1gydF4y2Ba} （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba$ 纵坐标gydF4y2Ba ${xgydF4y2Ba}_{2gydF4y2Ba} （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba$ ．质量gydF4y2Ba $我gydF4y2Ba$ 的重力有势能吗gydF4y2Ba $ggydF4y2Ba 米gydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba {xgydF4y2Ba}_{2gydF4y2Ba} （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba$ ．春天的势能gydF4y2Ba $我gydF4y2Ba$ 是gydF4y2Ba $kgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba （gydF4y2Ba dgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba -gydF4y2Ba lgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba {）gydF4y2Ba}^{2gydF4y2Ba} /gydF4y2Ba 2gydF4y2Ba$ ,在那里gydF4y2Ba $dgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba$ 弹簧的长度在质量之间吗gydF4y2Ba $我gydF4y2Ba$ 和质量gydF4y2Ba $我gydF4y2Ba -gydF4y2Ba 1gydF4y2Ba$ ．设锚点1为质量0的位置，锚点2为质量的位置gydF4y2Ba $ngydF4y2Ba +gydF4y2Ba 1gydF4y2Ba$ ．由前面的能量计算可知，弹簧的势能gydF4y2Ba $我gydF4y2Ba$ 是gydF4y2Ba

$EgydF4y2Ba ngydF4y2Ba egydF4y2Ba rgydF4y2Ba ggydF4y2Ba ygydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba ＝gydF4y2Ba \frac{kgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba （gydF4y2Ba 为gydF4y2Ba xgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba -gydF4y2Ba xgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba -gydF4y2Ba 1gydF4y2Ba ）gydF4y2Ba 为gydF4y2Ba -gydF4y2Ba lgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba {）gydF4y2Ba}^{2gydF4y2Ba}}{2gydF4y2Ba}$ ．gydF4y2Ba

将势能问题重新表述为二阶锥问题需要引入一些新的变量，如Lobo中所述gydF4y2Ba[1]gydF4y2Ba．创建变量gydF4y2Ba $tgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba$ 等于这项的平方根gydF4y2Ba $EgydF4y2Ba ngydF4y2Ba egydF4y2Ba rgydF4y2Ba ggydF4y2Ba ygydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba$ ．gydF4y2Ba

$tgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba ＝gydF4y2Ba \sqrt{\frac{kgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba （gydF4y2Ba 为gydF4y2Ba xgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba -gydF4y2Ba xgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba -gydF4y2Ba 1gydF4y2Ba ）gydF4y2Ba 为gydF4y2Ba -gydF4y2Ba lgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba {）gydF4y2Ba}^{2gydF4y2Ba}}{2gydF4y2Ba}} ．gydF4y2Ba$

让gydF4y2Ba $egydF4y2Ba$ 是单位列向量gydF4y2Ba $［gydF4y2Ba \binom{0gydF4y2Ba}{1gydF4y2Ba} ］gydF4y2Ba$ ．然后gydF4y2Ba ${xgydF4y2Ba}_{2gydF4y2Ba} （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba ＝gydF4y2Ba {egydF4y2Ba}^{TgydF4y2Ba} xgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba$ ．问题变成了gydF4y2Ba

$\underset{xgydF4y2Ba ，gydF4y2Ba tgydF4y2Ba}{最小值gydF4y2Ba} （gydF4y2Ba \underset{我gydF4y2Ba}{\sumgydF4y2Ba} ggydF4y2Ba 米gydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba {egydF4y2Ba}^{TgydF4y2Ba} xgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba +gydF4y2Ba 为gydF4y2Ba tgydF4y2Ba {为gydF4y2Ba}^{2gydF4y2Ba} ）gydF4y2Ba ．gydF4y2Ba$ (1）gydF4y2Ba

现在考虑gydF4y2Ba $tgydF4y2Ba$ 作为一个自由向量变量，不是由前面的方程给出的gydF4y2Ba $tgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba$ ．结合之间的关系gydF4y2Ba $xgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba$ 而且gydF4y2Ba $tgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba$ 在新的锥约束集合中gydF4y2Ba

$为gydF4y2Ba xgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba -gydF4y2Ba xgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba -gydF4y2Ba 1gydF4y2Ba ）gydF4y2Ba 为gydF4y2Ba -gydF4y2Ba lgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba \leqgydF4y2Ba \sqrt{\frac{2gydF4y2Ba}{kgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba}} tgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba ．gydF4y2Ba$ （2）gydF4y2Ba

目标函数的变量还不是线性的，如所要求的gydF4y2BaconeproggydF4y2Ba．引入一个新的标量变量gydF4y2Ba $ygydF4y2Ba$ ．注意这个不等式gydF4y2Ba $为gydF4y2Ba tgydF4y2Ba {为gydF4y2Ba}^{2gydF4y2Ba} \leqgydF4y2Ba ygydF4y2Ba$ 等于这个不等式吗gydF4y2Ba

$为gydF4y2Ba ［gydF4y2Ba \binom{2gydF4y2Ba tgydF4y2Ba}{1gydF4y2Ba -gydF4y2Ba ygydF4y2Ba} ］gydF4y2Ba 为gydF4y2Ba \leqgydF4y2Ba 1gydF4y2Ba +gydF4y2Ba ygydF4y2Ba$ ．（3）gydF4y2Ba

现在的问题是最小化gydF4y2Ba

$\underset{xgydF4y2Ba ，gydF4y2Ba tgydF4y2Ba ，gydF4y2Ba ygydF4y2Ba}{最小值gydF4y2Ba} （gydF4y2Ba \underset{我gydF4y2Ba}{\sumgydF4y2Ba} ggydF4y2Ba 米gydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba {egydF4y2Ba}^{TgydF4y2Ba} xgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba +gydF4y2Ba ygydF4y2Ba ）gydF4y2Ba$ （4）gydF4y2Ba

受锥上的约束gydF4y2Ba $xgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba$ 而且gydF4y2Ba $tgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba$ (2)和附加锥约束(3)。锥约束(3)确保gydF4y2Ba $为gydF4y2Ba tgydF4y2Ba {为gydF4y2Ba}^{2gydF4y2Ba} \leqgydF4y2Ba ygydF4y2Ba$ ．因此，问题(4)等价于问题(1)。gydF4y2Ba

问题(4)中的目标函数和锥约束适用于求解gydF4y2BaconeproggydF4y2Ba．gydF4y2Ba

MATLAB®配方gydF4y2Ba

定义六个弹簧常数gydF4y2Ba $kgydF4y2Ba$ ，六个长度常数gydF4y2Ba $lgydF4y2Ba$ ，和五个质量gydF4y2Ba $米gydF4y2Ba$ ．gydF4y2Ba

K = 40*(1:6);L = [1 1/2 1 2 1 1/2];M = [2 1 3 2 1];gydF4y2Ba

定义地球上的近似引力常数gydF4y2Ba $ggydF4y2Ba$ ．gydF4y2Ba

G = 9.807;gydF4y2Ba

优化的变量是10个组成部分gydF4y2Ba $xgydF4y2Ba$ 向量的六个分量gydF4y2Ba $tgydF4y2Ba$ 向量，而gydF4y2Ba $ygydF4y2Ba$ 变量。让gydF4y2BavgydF4y2Ba是包含所有这些变量的向量。gydF4y2Ba

(v (1), (2))gydF4y2Ba对应于二维变量gydF4y2Ba $xgydF4y2Ba （gydF4y2Ba 1gydF4y2Ba ）gydF4y2Ba$ ．gydF4y2Ba
(v (3), (4))gydF4y2Ba对应于二维变量gydF4y2Ba $xgydF4y2Ba （gydF4y2Ba 2gydF4y2Ba ）gydF4y2Ba$ ．gydF4y2Ba
(v (5), (6)]gydF4y2Ba对应于二维变量gydF4y2Ba $xgydF4y2Ba （gydF4y2Ba 3.gydF4y2Ba ）gydF4y2Ba$ ．gydF4y2Ba
(v(7)、(8)]gydF4y2Ba对应于二维变量gydF4y2Ba $xgydF4y2Ba （gydF4y2Ba 4gydF4y2Ba ）gydF4y2Ba$ ．gydF4y2Ba
(v(9)、(10))gydF4y2Ba对应于二维变量gydF4y2Ba $xgydF4y2Ba （gydF4y2Ba 5gydF4y2Ba ）gydF4y2Ba$ ．gydF4y2Ba
[v (11): v (16))gydF4y2Ba对应于6维向量gydF4y2Ba $tgydF4y2Ba$ ．gydF4y2Ba
v (17)gydF4y2Ba对应于标量变量gydF4y2Ba $ygydF4y2Ba$ ．gydF4y2Ba

使用这些变量，创建相应的目标函数向量gydF4y2BafgydF4y2Ba．gydF4y2Ba

F = 0 (size(m));F = [F;g*m];F = F (:);F = [F; 0(长度(k)+1,1)];F (end) = 1;gydF4y2Ba

创建与质量(2)之间的弹簧相对应的锥约束gydF4y2Ba

的gydF4y2BaconeproggydF4y2Ba求解器对变量向量使用锥约束gydF4y2Ba $vgydF4y2Ba$ 在表格中gydF4y2Ba

$为gydF4y2Ba {一个gydF4y2Ba}_{年代gydF4y2Ba cgydF4y2Ba} \cdotgydF4y2Ba vgydF4y2Ba -gydF4y2Ba {bgydF4y2Ba}_{年代gydF4y2Ba cgydF4y2Ba} 为gydF4y2Ba \leqgydF4y2Ba {dgydF4y2Ba}_{年代gydF4y2Ba cgydF4y2Ba}^{TgydF4y2Ba} vgydF4y2Ba -gydF4y2Ba γgydF4y2Ba$ ．gydF4y2Ba

在下面的代码中，gydF4y2BaAscgydF4y2Ba矩阵表示项gydF4y2Ba $为gydF4y2Ba xgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba -gydF4y2Ba xgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba -gydF4y2Ba 1gydF4y2Ba ）gydF4y2Ba 为gydF4y2Ba$ ,gydF4y2Ba二元同步通信gydF4y2Ba＝gydF4y2Ba(0, 0)gydF4y2Ba．锥变量gydF4y2BadscgydF4y2Ba＝gydF4y2Ba $\sqrt{2gydF4y2Ba /gydF4y2Ba kgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba}$ 和相应的gydF4y2BaγgydF4y2Ba＝gydF4y2Ba $-gydF4y2Ba lgydF4y2Ba （gydF4y2Ba 我gydF4y2Ba ）gydF4y2Ba ．gydF4y2Ba$

D = 0(1，长度(f));Asc = d;Asc([1 3]) = [1 -1];A2 = circshift(Asc,1);Asc = [Asc;A2];Ml =长度(m);Dbase = 2*ml;BSC = [0;0];gydF4y2Ba为gydF4y2BaI = 2:ml gamma = -l(I);DSC = d;Dsc (dbase + i) =√(2/k(i));conecons(i) = secondordercone(Asc,bsc,dsc,gamma);Asc = circshift(Asc,2,2);gydF4y2Ba结束gydF4y2Ba

通过使用末端质量位置的锚点，创建与末端质量和锚点之间的弹簧相对应的锥约束，如前面的代码中所示。gydF4y2Ba

X0 = [0;5];Xn = [5;4];Asc = 0 (size(Asc));Asc(1，(dbase-1)) = 1;Asc(2,dbase) = 1;BSC = xn;Gamma = -l(ml);DSC = d;Dsc (dbase + ml) =√(2/k(ml));conecons(ml + 1) = secondordercone(Asc,bsc,dsc,gamma); Asc = zeros(size(Asc)); Asc(1,1) = 1; Asc(2,2) = 1; bsc = x0; gamma = -l(1); dsc = d; dsc(dbase + 1) = sqrt(2/k(1)); conecons(1) = secondordercone(Asc,bsc,dsc,gamma);

创建圆锥约束(3)对应于gydF4y2Ba $ygydF4y2Ba$ 变量gydF4y2Ba

$为gydF4y2Ba ［gydF4y2Ba \binom{2gydF4y2Ba tgydF4y2Ba}{1gydF4y2Ba -gydF4y2Ba ygydF4y2Ba} ］gydF4y2Ba 为gydF4y2Ba \leqgydF4y2Ba 1gydF4y2Ba +gydF4y2Ba ygydF4y2Ba$

通过创建矩阵gydF4y2BaAscgydF4y2Ba乘以什么gydF4y2BavgydF4y2Ba向量，给出了向量gydF4y2Ba $［gydF4y2Ba \binom{2gydF4y2Ba tgydF4y2Ba}{-gydF4y2Ba ygydF4y2Ba} ］gydF4y2Ba$ ．的gydF4y2Ba二元同步通信gydF4y2Ba向量对应于项中的常数1gydF4y2Ba $1gydF4y2Ba -gydF4y2Ba ygydF4y2Ba$ ．的gydF4y2BadscgydF4y2Ba向量乘以gydF4y2BavgydF4y2Ba,返回gydF4y2Ba $ygydF4y2Ba$ ．和gydF4y2BaγgydF4y2Ba＝gydF4y2Ba $-gydF4y2Ba 1gydF4y2Ba$ ．gydF4y2Ba

Asc = 2*眼(长(f));Asc(1:dbase，:) = [];Asc(end,end) = -1;bsc = 0 (size(Asc,1)，1);Bsc (end) = -1;DSC = d;Dsc (end) = 1;Gamma = -1;conecons(ml+2) = secondordercone(Asc,bsc,dsc,gamma);gydF4y2Ba

最后，创建对应于gydF4y2Ba $tgydF4y2Ba$ 而且gydF4y2Ba $ygydF4y2Ba$ 变量。gydF4y2Ba

Lb = -inf(size(f));Lb (dbase+1:end) = 0;gydF4y2Ba

解决问题并设计解决方案gydF4y2Ba

问题的公式是完整的。打电话解决问题gydF4y2BaconeproggydF4y2Ba．gydF4y2Ba

[v, fval、exitflag、输出]= coneprog (f, conecons ,[],[],[],[], 磅);gydF4y2Ba

找到最优解。gydF4y2Ba

画出解决点和锚点。gydF4y2Ba

Pp = v(1:2*ml);Pp =重塑(Pp,2，[]);Pp = Pp ';情节(pp(: 1)、pp (:, 2),gydF4y2Ba“罗”gydF4y2Ba)举行gydF4y2Ba在gydF4y2BaXx = [x0,xn]';情节(xx (: 1), xx (:, 2),gydF4y2Ba“ks”gydF4y2Ba) xlim ((x0 (1) -0.5, xn (1) + 0.5]) ylim ([min (pp(:, 2)) -0.5,马克斯(x0 (2), xn (2)) + 0.5]) xxx = [x0; pp; xn);情节(xxx (: 1) xxx (:, 2),gydF4y2Ba“b——”gydF4y2Ba)传说(gydF4y2Ba“计算点”gydF4y2Ba，gydF4y2Ba“锚点”gydF4y2Ba，gydF4y2Ba“泉”gydF4y2Ba，gydF4y2Ba“位置”gydF4y2Ba，gydF4y2Ba“最佳”gydF4y2Ba)举行gydF4y2Ba从gydF4y2Ba

图中包含一个轴对象。axis对象包含3个line类型的对象。这些对象表示计算点、锚点、弹簧。gydF4y2Ba

您可以修改参数的值gydF4y2Ba米gydF4y2Ba，gydF4y2BalgydF4y2Ba,gydF4y2BakgydF4y2Ba看看它们是如何影响解的。你也可以改变质量的数量;该代码从您提供的数据中获取质量数。gydF4y2Ba

参考文献gydF4y2Ba

[1] Lobo, Miguel Sousa, Lieven Vandenberghe, Stephen Boyd和Hervé Lebret。二阶锥规划的应用gydF4y2Ba线性代数及其应用gydF4y2Ba284年,没有。1-3(1998年11月):193-228。gydF4y2Bahttps://doi.org/10.1016/s0024 - 3795 (98) 10032 - 0gydF4y2Ba．gydF4y2Ba

另请参阅gydF4y2Ba

coneproggydF4y2Ba|gydF4y2BasecondorderconegydF4y2Ba