Main Content

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

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 el solucionador de este ejemplo, vea.Programación cuadrática con restricción de límite, basada en Solver

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. Puede aproximar la segunda derivada por la aproximación de la diferencia finita de cinco puntos (supongamos que los pasos de la diferencia finita son de 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')

Formular problema de optimización

Cree una variable de optimización que represente la altura del material.x

x = optimvar('x',size(height));

Se establece en cero en los límites del dominio cuadrado.x

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

Calcule la energía potencial elástica en cada punto. En primer lugar, calcule la energía potencial en el interior de la región, donde las diferencias finitas no exageran la región que contiene la solución.

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);

Dado que la solución está restringida a ser 0 en los bordes de la región, no es necesario incluir el resto de los términos. Todos los términos tienen un múltiplo de, y en el borde es cero.xx Como referencia en caso de que desee utilizar una condición de contorno diferente, la siguiente es una versión comentada de la energía potencial.

% peStretch(1,interior) = (-1*(x(1,interior - 1) + x(1,interior + 1) + x(2,interior))... %     + 4*x(1,interior)).*x(1,interior); % peStretch(L,interior) = (-1*(x(L,interior - 1) + x(L,interior + 1) + x(L-1,interior))... %     + 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); % peStretch(interior,L) = (-1*(x(interior - 1,L) + x(interior + 1,L) + x(interior,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); % 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);

Definir la energía potencial debido a la altura del material, que es.x/3000

peHeight = x/3000;

Cree un problema de optimización denominado.tentproblem Incluya la expresión para la función objetiva, que es la suma de las dos energías potenciales en todas las ubicaciones.

tentproblem = optimproblem('Objective',sum(sum(peStretch + peHeight)));

Establecer restricción

Establezca la restricción de que la solución debe estar por encima de los valores de la matriz.height Esta matriz es cero en la mayoría de los lugares, representando el suelo, e incluye la altura de cada poste de la carpa en su ubicación.

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

Ejecute Optimization Solver

Resuelve el problema. Ignore la declaración resultante "su hessian no es simétrico." emite esta instrucción porque la conversión interna de un problema a una matriz cuadrática no garantiza que la matriz sea simétrica.solve

sol = solve(tentproblem);
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.

Solución de trazado

Trace la solución encontrada por el solucionador de optimización.

surfl(sol.x); 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