基于问题的有界约束二次规划

打开实时脚本

这个例子展示了如何通过求解二次优化问题来确定马戏团帐篷的形状。帐篷由沉重的、有弹性的材料构成，并形成一个受限制的最小势能的形状。将问题离散化，得到有界-有约束二次规划问题。

问题定义

考虑建造一个马戏团帐篷来覆盖一个广场。帐篷有五根杆子，上面覆盖着一种沉重的、有弹性的材料。问题是要找到帐篷的自然形状。将形状建模为高度x（p)帐篷的位置p．

重力:重物升到高度时的势能x是残雪，对于常数c这与材料的重量成正比。对于这个问题，请选择c= 1/3000。

一块材料的弹性势能 $E_{年代 t r e t c h}$ 大约正比于材料高度的二阶导数，乘以高度。你可以用五点有限差分近似来近似二阶导数(假设有限差分步长为1) $Δ x$ 表示在第一个坐标方向上移位1，和 $Δ y$ 表示在第二个坐标方向上的位移为1。

$E_{年代 t r e t c h} （ p ）＝（ - 1 （ x （ p + Δ_{x} ） + x （ p - Δ_{x} ） + x （ p + Δ_{y} ） + x （ p - Δ_{y} ）） + 4 x （ p ）） x （ p ）．$

帐篷的自然形状使总势能最小化。通过离散化问题，你会发现要最小化的总势能是所有位置的和p的 $E_{年代 t r e t c h} （ p ）$ +残雪（p）.

势能是变量中的二次表达式x．

指定边界条件，即帐篷边缘的高度为零。帐篷杆的横截面为1 × 1单位，帐篷的总尺寸为33 × 33单位。指定每个极点的高度和位置。划出方形地段和帐篷杆。

高度= 0 (33);高度(6:7,6:7)= 0.3;高度(26:27,26:27)= 0.3;身高(6:7,26:27)= 0.3;身高(26:27,6:7)= 0.3;高度(16:17,16:17)= 0.5;colormap(灰色);surfl(高度)轴紧视图([-20,30]);标题(“帐篷杆和覆盖区域”）

图中包含一个轴对象。标题为Tent Poles和Region to Cover的axis对象包含一个类型为surface的对象。

公式优化问题

创建优化变量x表示材料的高度。

X = optimvar(“x”、大小(高度));

集x在方形域的边界上求零。

border = false(size(height));Boundary ([1,33]，:) = true;Boundary (:，[1,33]) = true;x.LowerBound(边界)= 0;x.UpperBound(boundary) = 0;

计算每个点的弹性势能。首先，计算区域内部的势能，其中有限差分不超过包含解的区域。

L = size(height,1);peStretch = optimexpr(L,L);将peStretch初始化为0 (L,L)interior = 2:(L-1);peStretch(内部，内部)= (-1*(x(内部-1，内部)+ x(内部+ 1，内部)).．.+ x(室内，室内- 1)+ x(室内，室内+ 1)+ 4*x(室内，室内).．.. * x(室内,室内);

因为解在区域的边缘被限制为0，所以不需要包括其余的项。所有的项都有倍数x,x在边缘处为零。作为参考，如果你想使用不同的边界条件，下面是一个注释掉的势能版本。

% peStretch(1，内部)= (-1*(x(1，内部-1)+ x(1，内部+ 1)+ x(2，内部))…% + 4*x(1，内部)).*x(1，内部);% peStretch (L,内政部)= (1 * (x (L,室内- 1)+ x (L,室内+ 1)+ x (L - 1,室内))…% + 4*x(L，内部)).*x(L，内部);% peStretch(内部，1)= (-1*(x(内部- 1,1)+ x(内部+ 1,1)+ x(内部，2))…*x(interior,1)).*x(interior,1);% peStretch(内部,L) = (1 * (x(室内- 1,L) + x(内部+ 1,L) +(内部,L - 1))…% + 4*x(内部，L)).*x(内部，L);% peStretch (1, - 1) = (1 * (x (2, 1) + x(1、2))+ 4 * x(1,1))。* x (1,1);% peStretch (L) = (1 * (x (2 L) + x (1, L - 1)) + 4 * x(1升))。* x(1升);% peStretch (L, 1) = (1 * (x (L, 2) + x (L - 1, - 1)) + 4 * x (L, 1))。* x (L, 1);% peStretch (L, L) = (1 * (x (L - 1, L) + x (L, L - 1)) + 4 * x (L, L)。* x (L, L);

定义由材料高度引起的势能，即x / 3000．

peHeight = x/3000;

创建一个名为tentproblem．包括目标函数的表达式，它是所有位置上两个势能的和。

Tentproblem = optimproblem(“目标”，sum(sum(peStretch + peHeight)));

设置约束

的值之上设置解的约束高度矩阵。该矩阵在大多数位置为零，代表地面，并包括每个帐篷杆在其位置的高度。

Htcons = x >= height;tentproblem.Constraints.htcons = htcons;

运行优化求解器

解决问题。忽略结果语句“Your Hessian is not对称”。解决发出这个声明是因为从问题形式到二次矩阵的内部转换不能保证矩阵是对称的。

Sol =解决(帐篷问题);

用quadprog解决问题。你的黑纱是不对称的。重置H = (H + H) / 2。最小值满足约束条件。优化完成是因为目标函数在可行方向上不递减，在最优性容差值范围内，约束条件满足在约束容差值范围内。

策划解决方案

绘制由优化求解器找到的解。

surfl (sol.x);轴紧；视图([-20,30]);

图中包含一个轴对象。axis对象包含一个surface类型的对象。