Esta página aún no se ha traducido para esta versión. Puede ver la versión más reciente de esta página en inglés.

Programación cuadrática con restricción de límite, basada en Solver

Abrir script en vivo

Este ejemplo muestra cómo determinar la forma de una carpa de circo resolviendo un problema de optimización cuadrática. La carpa está formada por material pesado y elástico, y se asienta en una forma que tiene una energía potencial mínima sujeta a restricciones. Una discretización del problema conduce a un problema de programación cuadrática limitado.

Para obtener una versión basada en problemas de este ejemplo, vea.Programación cuadrática con restricción de límite, basada en problemas

Definición de problema

Considere construir una carpa de circo para cubrir un lote cuadrado. La carpa tiene cinco postes cubiertos con un material pesado y elástico. El problema es encontrar la forma natural de la carpa. Modele la forma como la altura () de la carpa en posición.xpp

La energía potencial del material pesado elevado a la altura es, para una constante que es proporcional al peso del material.xcxc Para este problema, elija = 1/3000.c

La energía potencial elástica de una pieza del material

<math display="block">
<mrow>
<msub>
<mrow>
<mi>E</mi>
</mrow>
<mrow>
<mstyle mathvariant="normal">
<mrow>
<mi>s</mi>
<mi>t</mi>
<mi>r</mi>
<mi>e</mi>
<mi>t</mi>
<mi>c</mi>
<mi>h</mi>
</mrow>
</mstyle>
</mrow>
</msub>
</mrow>
</math>
es aproximadamente proporcional a la segunda derivada de la altura del material, multiplicado por la altura. Usted puede aproximar la segunda derivada por la aproximación de la diferencia finita de 5 puntos (asuma que los pasos de la diferencia finita son del tamaño 1). Dejar
<math display="block">
<mrow>
<mi>Δ</mi>
<mi>x</mi>
</mrow>
</math>
 representan un desplazamiento de 1 en la primera dirección de coordenadas y
<math display="block">
<mrow>
<mi>Δ</mi>
<mi>y</mi>
</mrow>
</math>
 representan un desplazamiento de 1 en la segunda dirección de coordenadas.

<math display="block">
<mrow>
<msub>
<mrow>
<mi>E</mi>
</mrow>
<mrow>
<mstyle mathvariant="normal">
<mrow>
<mrow>
<mi>s</mi>
<mi>t</mi>
<mi>r</mi>
<mi>e</mi>
<mi>t</mi>
<mi>c</mi>
<mi>h</mi>
</mrow>
</mrow>
</mstyle>
</mrow>
</msub>
<mo stretchy="false">(</mo>
<mi>p</mi>
<mo stretchy="false">)</mo>
<mo>=</mo>
<mo maxsize="2.4" minsize="2.4">(</mo>
<mo>-</mo>
<mn>1</mn>
<mo maxsize="2.4" minsize="2.4">(</mo>
<mi>x</mi>
<mo stretchy="false">(</mo>
<mi>p</mi>
<mo>+</mo>
<msub>
<mrow>
<mi>Δ</mi>
</mrow>
<mrow>
<mi>x</mi>
</mrow>
</msub>
<mo stretchy="false">)</mo>
<mo>+</mo>
<mi>x</mi>
<mo stretchy="false">(</mo>
<mi>p</mi>
<mo>-</mo>
<msub>
<mrow>
<mi>Δ</mi>
</mrow>
<mrow>
<mi>x</mi>
</mrow>
</msub>
<mo stretchy="false">)</mo>
<mo>+</mo>
<mi>x</mi>
<mo stretchy="false">(</mo>
<mi>p</mi>
<mo>+</mo>
<msub>
<mrow>
<mi>Δ</mi>
</mrow>
<mrow>
<mi>y</mi>
</mrow>
</msub>
<mo stretchy="false">)</mo>
<mo>+</mo>
<mi>x</mi>
<mo stretchy="false">(</mo>
<mi>p</mi>
<mo>-</mo>
<msub>
<mrow>
<mi>Δ</mi>
</mrow>
<mrow>
<mi>y</mi>
</mrow>
</msub>
<mo stretchy="false">)</mo>
<mo maxsize="2.4" minsize="2.4">)</mo>
<mo>+</mo>
<mn>4</mn>
<mi>x</mi>
<mo stretchy="false">(</mo>
<mi>p</mi>
<mo stretchy="false">)</mo>
<mo maxsize="2.4" minsize="2.4">)</mo>
<mi>x</mi>
<mo stretchy="false">(</mo>
<mi>p</mi>
<mo stretchy="false">)</mo>
<mo>.</mo>
</mrow>
</math>

La forma natural de la carpa minimiza la energía potencial total. Al discretización del problema, usted encuentra que la energía potencial total para minimizar es la suma sobre todas las posiciones dep

<math display="block">
<mrow>
<msub>
<mrow>
<mi>E</mi>
</mrow>
<mrow>
<mstyle mathvariant="normal">
<mrow>
<mrow>
<mi>s</mi>
<mi>t</mi>
<mi>r</mi>
<mi>e</mi>
<mi>t</mi>
<mi>c</mi>
<mi>h</mi>
</mrow>
</mrow>
</mstyle>
</mrow>
</msub>
<mo stretchy="false">(</mo>
<mi>p</mi>
<mo stretchy="false">)</mo>
</mrow>
</math>
+ ( ).cxp

Esta energía potencial es una expresión cuadrática en la variable.x

Especifique la condición de contorno que la altura de la carpa en los bordes es cero. Los postes de la carpa tienen una sección transversal de 1 por 1 unidad, y la carpa tiene un tamaño total de 33-por-33 unidades. Especifique la altura y la ubicación de cada polo. Trazar la región de lote cuadrado y postes de carpa.

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 tight view([-20,30]); title('Tent Poles and Region to Cover')

Crear condiciones de contorno

La matriz define los límites inferiores de la solución.heightx Para restringir la solución a cero en el contorno, establezca el límite superior en cero en el contorno.ub

boundary = false(size(height)); boundary([1,33],:) = true; boundary(:,[1,33]) = true; ub = inf(size(boundary)); % No upper bound on most of the region ub(boundary) = 0;

Crear matrices de funciones objetivas

La formulación del problema es minimizarquadprog

<math display="block">
<mrow>
<mfrac>
<mrow>
<mn>1</mn>
</mrow>
<mrow>
<mn>2</mn>
</mrow>
</mfrac>
<msup>
<mrow>
<mi>x</mi>
</mrow>
<mrow>
<mi>T</mi>
</mrow>
</msup>
<mi>H</mi>
<mi>x</mi>
<mo>+</mo>
<msup>
<mrow>
<mi>f</mi>
</mrow>
<mrow>
<mi>T</mi>
</mrow>
</msup>
<mi>x</mi>
</mrow>
</math>
.

En este caso, el término lineal

<math display="block">
<mrow>
<msup>
<mrow>
<mi>f</mi>
</mrow>
<mrow>
<mi>T</mi>
</mrow>
</msup>
<mi>x</mi>
</mrow>
</math>
corresponde a la energía potencial de la altura del material. Por lo tanto, especifique = 1/3000 para cada componente de.fx

f = ones(size(height))/3000;

Cree la matriz de diferencias finitas que representa

<math display="block">
<mrow>
<msub>
<mrow>
<mi>E</mi>
</mrow>
<mrow>
<mstyle mathvariant="normal">
<mrow>
<mrow>
<mi>s</mi>
<mi>t</mi>
<mi>r</mi>
<mi>e</mi>
<mi>t</mi>
<mi>c</mi>
<mi>h</mi>
</mrow>
</mrow>
</mstyle>
</mrow>
</msub>
</mrow>
</math>
mediante la función.delsq La función devuelve una matriz dispersa con entradas de 4 y-1 correspondientes a las entradas de 4 y-1 en la fórmula paradelsq
<math display="block">
<mrow>
<msub>
<mrow>
<mi>E</mi>
</mrow>
<mrow>
<mstyle mathvariant="normal">
<mrow>
<mrow>
<mi>s</mi>
<mi>t</mi>
<mi>r</mi>
<mi>e</mi>
<mi>t</mi>
<mi>c</mi>
<mi>h</mi>
</mrow>
</mrow>
</mstyle>
</mrow>
</msub>
<mo stretchy="false">(</mo>
<mi>p</mi>
<mo stretchy="false">)</mo>
</mrow>
</math>
. Multiplique la matriz devuelta por 2 para que resuelva el programa cuadrático con la función de energía dada porquadprog
<math display="block">
<mrow>
<msub>
<mrow>
<mi>E</mi>
</mrow>
<mrow>
<mstyle mathvariant="normal">
<mrow>
<mrow>
<mi>s</mi>
<mi>t</mi>
<mi>r</mi>
<mi>e</mi>
<mi>t</mi>
<mi>c</mi>
<mi>h</mi>
</mrow>
</mrow>
</mstyle>
</mrow>
</msub>
</mrow>
</math>
.

H = delsq(numgrid('S',33+2))*2;

Ver la estructura de la matriz.H La matriz funciona, lo que significa que la matriz se convierte en un vector por indexación lineal.x(:)x

spy(H); title('Sparsity Structure of Hessian Matrix');

Ejecute Optimization Solver

Resuelva el problema llamando.quadprog

x = quadprog(H,f,[],[],[],[],height,ub);
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.

Solución de trazado

Remodele la solución a una matriz.xS A continuación, trace la solución.

S = reshape(x,size(height)); surfl(S); axis tight; view([-20,30]);

Temas relacionados

Documentación de Optimization Toolbox
Soporte
8 MATLAB Cheat Sheets for Data Science

8 MATLAB Cheat Sheets for Data Science

Download now