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.
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
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
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
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')
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.x
x
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)));
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;
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.
Trace la solución encontrada por el solucionador de optimización.
surfl(sol.x); axis tight; view([-20,30]);