问题定义

考虑建造马戏团帐篷覆盖平方批次。帐篷有五极覆盖着沉重的弹性材料。问题是找到帐篷的自然形状。将形状塑造为高度x(p）帐篷在位置p。

The potential energy of heavy material lifted to heightx是cx，常数c这与材料的重量成比例。对于这个问题，选择c= 1/3000。

材料的弹性势能 $$E_{stretch}$$是approximately proportional to the second derivative of the material height, times the height. You can approximate the second derivative by the five-point finite difference approximation (assume that the finite difference steps are of size 1). Let $$\mathrm{\delta .}x$$表示在第一坐标方向上1的偏移， $$\mathrm{\delta .}y$$represent a shift of 1 in the second coordinate direction.

帐篷的自然形状最小化总势能。通过离散问题，您发现最小化的总潜在能量是所有位置的总和p的 $$E_{stretch}(p)$$+cx(p）。

这种势能是变量中的二次表达x。

Specify the boundary condition that the height of the tent at the edges is zero. The tent poles have a cross section of 1-by-1 unit, and the tent has a total size of 33-by-33 units. Specify the height and location of each pole. Plot the square lot region and tent poles.

height = zeros(33); height(6:7,6:7) = 0.3; height(26:27,26:27) = 0.3; height(6:7,26:27) = 0.3; height(26:27,6:7) = 0.3; height(16:17,16:17) = 0.5; colormap(gray); surfl(height) axis紧的查看（[ -  20,30]）;标题（'Tent Poles and Region to Cover')

制定优化问题

Create an optimization variablexrepresenting the height of the material.

x = optimvar('X'，尺寸（高度））;

Setx在方框的边界上为零。

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

Calculate the elastic potential energy at each point. First, calculate the potential energy in the interior of the region, where the finite differences do not overstep the region containing the solution.

L = size(height,1); peStretch = optimexpr(L,L);% This initializes peStretch to zeros(L,L)interior = 2:(L-1); peStretch(interior,interior) = (-1*(x(interior - 1,interior) + x(interior + 1,interior)。。。+ x(interior,interior - 1) + x(interior,interior + 1)) + 4*x(interior,interior))。。。。*x(interior, interior);

Because the solution is constrained to be 0 at the edges of the region, you do not need to include the remainder of the terms. All terms have a multiple ofx， 和x在边缘为零。出于参考，在您想要使用不同的边界条件的情况下，以下是潜在能量的注释版本。

％Perthetch（1，内部）=（-1 *（x（1，内部 -  1）+ x（1，内部+ 1）+ x（2，内部））...% + 4*x(1,interior)).*x(1,interior);％Pertrach（l，内部）=（-1 *（x（l，内部 -  1）+ x（l，内部+ 1）+ x（l-1，内部））...% + 4*x(L,interior)).*x(L,interior);% peStretch(interior,1) = (-1*(x(interior - 1,1) + x(interior + 1,1) + x(interior,2))...% + 4*x(interior,1)).*x(interior,1);％Pertrch（内部，L）=（-1 *（x（内部 -  1，l）+ x（内部+ 1，l）+ x（内部，l-1））...% + 4*x(interior,L)).*x(interior,L);% peStretch(1,1) = (-1*(x(2,1) + x(1,2)) + 4*x(1,1)).*x(1,1);% peStretch(1,L) = (-1*(x(2,L) + x(1,L-1)) + 4*x(1,L)).*x(1,L);％perthetch（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);

Define the potential energy due to material height, which isX / 3000.。

peheight = x / 3000;

Create an optimization problem named帐篷问题。包括目标函数的表达式，这是所有位置的两个潜在能量的总和。

TentProblus = OptimProblem（'Objective',sum(sum(peStretch + peHeight)));

设置约束

Set the constraint that the solution must lie above the values of theheight矩阵。该矩阵在大多数位置处为零，表示地面，并且包括每个帐篷杆的位置处的位置。

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

Run Optimization Solver

Solve the problem. Ignore the resulting statement "Your Hessian is not symmetric."solve发布此语句，因为问题表单到二次矩阵的内部转换不会确保矩阵是对称的。

sol =解决（帐篷问题）;

Solving problem using quadprog. Your Hessian is not symmetric. Resetting H=(H+H')/2. Minimum found that satisfies the constraints. Optimization completed because the objective function is non-decreasing in feasible directions, to within the value of the optimality tolerance, and constraints are satisfied to within the value of the constraint tolerance.

Plot Solution

绘制优化求解器发现的解决方案。

surfl(sol.x); axis紧的;查看（[ -  20,30]）;