# PRACTICA 3: OPTIMIZACION CON MAXIMA

## MÉTODOS NUMÉRICOS

MAXIMA tiene implementados una serie de comandos para resolver problemas de optimizacion. 

Entre ellos se incluyen comandos para optimizacion con restricciones (fmin_cobyla, augmented_lagrangian_method), 
sin restricciones (lbfgs) y lineal (simplex). 
    
En casi todos los algoritmos numericos de optimizacion se procede de la misma forma: 

- Se toma un valor o aproximacion inicial: sería el valor que se elige en primera instancia como solución del problema.

- A partir de la primera aproximación, se construye una sucesion de aproximaciones sucesivas de forma que el valor 
  de la funcion objetivo va decreciendo o creciendo, segun el objetivo del problema sea de minimizar o maximizar, respectivamente.

- Se utiliza un criterio de terminacion que indica que la aproximacion esta suficientemente cerca de la solucion exacta del problema.

### Optimizacion sin restricciones.

La funcion empleada por MAXIMA para resolver problemas de optimizacion sin restricciones es:


    Funcion MAXIMA  : lbfgs
    Descripcion           : Minimizacion sin restricciones
    Algoritmo                : 	Cuasi Newton BFGS


La funcion lbfgs implementa el llamado algoritmo L-BFGS (metodo de Broyden-Fletcher-Goldfarb-Shanno 
con limitacion de memoria) para resolver el problema de minimizacion sin restricciones siguiente:

    min 	f(x)
    x real

Antes de poder utilizar este comando, es necesario cargar la libreria correspondiente:

	load(lbfgs);

Hay dos formas de utilizar lbfgs:

    lbfgs(Fobj, X, X0, epsilon, iprint)                 :  	Calculo automatico del gradiente de la función objetivo Fobj.
    lbfgs([Fobj, Grad], X, X0, epsilon, iprint)   :  El usuario debe incluir el gradiente de la función objetivo en la variable Grad.

Todos los parametros son obligatorios, salvo Grad.

    Fobj      : Funcion objetivo del problema.
    X            : Variables de decision del problema.
    X0          : Aproximacion inicial.
    epsilon : Error (criterio de parada).
    iprint     : Parametro que controla los mensajes que emite la funcion durante el proceso.

    Grad: Si este argumento esta presente, debe contener el gradiente de la funcion objetivo Fobj 
          respecto de las variables X. Puede ser una lista o una funcion que devuelva una lista con igual 
          numero de elementos que X. Si esta variable no esta presente, el gradiente se calcula 
          automaticamente mediante derivacion simbolica.

El cometido de estos comandos es encontrar una aproximacion a la solucion del problema de minimizacion 
sin restricciones cuya funcion objetivo viene dada por Fobj y variables de decision dadas en X, 
para ello se toma como aproximacion inicial la contenida en el vector X0, de tal manera que el criterio de parada es que 
se cumpla la siguiente desigualdad

    || grad(f(x)) || < epsilon*max(1, ||x|| )

siendo ||⋅|| la norma cuadratica o euclidea del vector correspondiente.

Los valores que toma el argumento iprint se describen a continuacion:

    iprint[1]      : Frecuencia de aparicion de los mensajes de información.
        Valor <0 : No se envian mensajes. Solo la solucion final.
        Valor =0 : Mensaje unicamente en la primera y ultima iteraciones.
        Valor >0 : Imprimer un mensaje con informacion del algoritmo cada iprint[1] iteraciones.

    
    iprint[2]    : Cantidad de información que se muestra.

        Valor 0 : Imprime contador de iteraciones
                  - Numero de evaluaciones de la funcion objetivo.
                  - Valor de la funcion objetivo.
                  - Modulo del gradiente de la funcion objetivo.
                  - Tamaño de paso.
        Valor 1 : Igual que para el valor 0, pero incluye ademas X0, 
                  y el gradiente de la funcion objetivo evaluado en X0.
        Valor 2 : Igual que para el valor 1, pero incluyendo ademas los valores
                  de X en cada iteracion del proceso.
        Valor 3 : Igual que para el valor 2, incluyendo el gradiente de la funcion
                  objetivo en cada iteracion.

La funcion lbfgs devuelve una serie de valores que serian los siguientes

    I                           : Numero de iteraciones. Se incrementa tras cada busqueda lineal.
    NFN                    : Numero de evaluaciones de la funcion objetivo.
    FUNC                 : Valor de la funcion objetivo al final de cada iteracion.
    GNORM             : Modulo del gradiente de la funcion objetivo al final de cada iteracion.
    STEPLENGTH : Parametro interno del algoritmo de busqueda. 

EJEMPLO 1: Encuentra el minimo de la funcion f(x) = x sin(x+1) utilizando 
X0=0 como aproximacion inicial y un error de 1/10^4

SOLUCION:

lbfgs(x*sin(x+1),[x],[0],1e-4,[1,0]);

Podemos comprobar que se trata de un minimo local, ya que si realizamos la representacion de la
funcion en el intervalo.

wxplot2d(x*sin(x+1),[x,-10,10]);

podemos observar que el minimo global de la funcion en el intervalo [-10,10], se encuentra cercano 
al punto 10, bastante alejado de la solucion que se ha encontrado con el comando lbfgs. 

NOTA: Debido a la presencia de la funcion periodica sin(*), es posible deducir que la funcion no tendra 
ni maximo, ni minimo global sobre la recta real R.

La solucion obtenida por la funcion lbfgs depende del punto de partida elegido, por ejemplo si 
cambiamos el punto de partida anterior por X0=2, y volvemos a utilizar la funcion lbfgs:

lbfgs(x*sin(x+1),[x],[2],1e-4,[1,0]);

se comprueba, efectivamente que esta respuesta, no coincide con la obtenida antes.
Observemos que en la grafica de la funcion es otro minimo local. Esta claro que el punto de partida
es fundamental a la hora de obtener una u otra solucion. 

El valor devuelto por el comando de MAXIMA, en este caso, lbfgs, es aquel que cumple  unos determinados criterios.

EJEMPLO 2: Podemos emplear la misma instruccion para funciones de varias variables, especificando en cada caso 
un valor inicial para cada variable, por ejemplo, para localizar un minimo local de la funcion 

f(x,y)=cos(x^2-3y)+sen(x^2+y^2) 

cerca del punto (1,1) utilizaremos la expresion:

	lbfgs( cos(x^2-3*y) + sin(x^2+y^2), [x,y], [1,1], 1e-3, [-1,0]);

En este caso tanto el vector X=[x,y], como el vector X0=[1,1] tienen dos componentes.

Si representamos graficamente la funcion:

wxplot3d(cos(x^2-3*y) + sin(x^2+y^2), [x,-2,2], [y,-2,2]);

/* Se puede utilizar el siguiente comando:

plot3d(cos(x^2-3*y) + sin(x^2+y^2), [x,-2,2], [y,-2,2]);

para abrir una ventana independiente y poder rotar la figura 3D */

podemos ver que la funcion tiene infinitos maximos y minimos y por tanto, como en el ejemplo en una variable, 
la respuesta dependera del punto de partida. 

En este caso tambien es interesante representar graficamente los isocontornos de la funcion:

	contour_plot(cos(x^2-3*y)+sin(x^2+y^2), [x,-2,2], [y,-2,2], [legend,false], [gnuplot_preamble,"set cntrparam levels 20"]);

donde se pueden observar las zonas donde se encuentran los extremos locales.

Si no se especifica, el gradiente se calcula automaticamente y de forma simbolica, 
no obstante es posible incluirlo como argumento en la funcion lbfgs, como en el siguiente ejemplo:

EJEMPLO 3: Encuentra mediante lbfgs el minimo de la funcion

F(x,y,z)=(x - 5)^2 + (y - 3)^4 + (z - 2)^6 

Incluye en la resolucion el gradiente de la funcion objetivo.

Solucion: 

F(x, y, z):=(x - 5)^2 + (y - 3)^4 + (z - 2)^6;
define(F_grad(x,y,z), map(lambda([t], diff(F(x,y,z),t)), [x,y,z]));
epsilon : 1e-4;
x0 : [0,0,0];
solucion : lbfgs([F,F_grad], [x,y,z], x0, epsilon, [0, 0]);

Teniendo en cuenta que estamos buscando el minimo de una funcion que es 
siempre positiva, por ser suma de cuadrados, la solucion optima del problema 
se encuentra en el unico punto donde se anula la funcion: 

x_exacta = (5,3,2)


	F(x,y,z)≥0=F(5,3,2)

Podemos conocer el error cometido al tomar como solucion el vector anterior:

error_absoluto = || x_aprox - x_exacta || 

x_aprox:[rhs(solucion[1]),rhs(solucion[2]),rhs(solucion[3])];
x_exacta:[5,3,2];
error_absoluto : abs(x_aprox-x_exacta);
/* Aquí calculamos la norma del error absoluto */
error_absoluto_norma : sqrt(dotproduct(transpose(error_absoluto),transpose(error_absoluto)));

Y podemos comprobar que se cumple el criterio de parada establecido en la fucion

||grad(f(x))|| < epsilon*max(1,||x||)

/* ||grad(F(x))|| */
v:F_grad(x_aprox[1],x_aprox[2],x_aprox[3]);
norma_gradiente: sqrt(dotproduct(transpose(v),transpose(v)));

/* ||x|| */
u:x_aprox;
norma_solucion: sqrt(dotproduct(transpose(u),transpose(u)));

/*  ||grad(F(x))|| < epsilon*max(1, ||X|| ) */
is(norma_gradiente < epsilon*max(1,norma_solucion));

Es posible definir las funciones empleadas como variables:

kill(all)$
load(lbfgs); 

F : (x - 5)^2 + (y - 3)^4 + (z - 2)^6;
F_grad: [diff(F,x), diff(F,y), diff(F,z)];
solucion : lbfgs([F,F_grad], [x,y,z], x0, epsilon, [0, 0]);

Hay dos variables globales relacionadas con la funcion lbfgs y que describimos a continuacion:

  ∙  lbfgs_nfeval_max (Valor por defecto: 100)

     La variable lbfgs_nfeval_max almacena el numero maximo de evaluaciones que lbfgs utiliza en la funcion objetivo. 
     Cuando se alcanza el valor lbfgs_nfeval_max, la funcion lbfgs devuelve el resultado de la ultima iteracion. 
     Este limite evita que se pueda entrar en un bucle debido a la precision finita del ordenador.


  ∙  lbfgs_ncorrections (Valor por defecto: 25)

     La variable lbfgs_ncorrections esta relacionada con el algoritmo de optimizacion empleado y tiene que ver 
     con una aproximacion a la inversa de la matriz hessiana en cada punto.

A continuacion se incluyen algunos ejemplos del uso de la funcion lbgfs.

EJEMPLO 4: Resuelve el siguiente problema con MAXIMA

            max exp(-x^2)

con un error  10^-3 y punto de partida x0=[0,5].

load (lbfgs);
fobjetivo(x) := exp(-x^2);
x0: [0.5];
epsilon:1e-3;
output:[1,0];
lbfgs(-fobjetivo(x),[x],x0,epsilon,output);

Gráficamente:

wxplot2d(fobjetivo(x),[x,-10,10]);

Tambien es posible definir la funcion en forma de variable:

	fobjetivo : exp(-x^2);

En este caso la llamada a la funcion lbfgs seria:

lbfgs(-fobjetivo,[x],x0,epsilon,output);

EJEMPLO 5: Resuelve el siguiente problema con MAXIMA

        min x^2 + y^2,

con un error de 10^-3 y punto de partida (0.5,0.5).

fobjetivo(x,y) := x^2 +  y^2;
x0: [0.5, 0.5];
epsilon:1e-3;
output:[1,0];
lbfgs(fobjetivo(x,y),[x,y],x0,epsilon,output);

Como antes tambien es posible definir la funcion como una variable:

	fobjetivo:x^2+y^2;

y en este caso la llamada a la funcion lbfgs seria:

	lbfgs(fobjetivo, [x,y], x0, epsilon, output);

Representación grafica:

wxplot3d(fobjetivo,[x,-2,2],[y,-2,2]);

contour_plot(fobjetivo, [x,-2,2], [y,-2,2], [legend,false], [gnuplot_preamble,"set cntrparam levels 20"]);

### Optimizacion con restricciones

Las funciones empleadas por MAXIMA para resolver problemas de optimizacion con restricciones son:


fmin_cobyla:
    Descripcion : 	Minimizacion con restricciones.
    Algoritmo   : 	Aproximaciones lineales sucesivas.


augmented_lagrangian_method:
    Descripcion : 	Minimizacion con restricciones h(x)=0	
    Algoritmo   : Funcion de penalizacion.

### Paquete COBYLA

El nombre de COBYLA es la abreviatura en ingles de Constrained Optimization BY Linear Aproximacion, es decir, 
optimizacion con restricciones mediante aproximaciones lineales.

COBYLA se utiliza para resolver problemas del tipo

    min  	f(x)
    Sujeto a	 	g(x)≥0

Las restricciones en forma de igualdades se pueden sustituir por pares de desigualdades:

    h(x) = 0 se substituiria por las dos restricciones h(x)>=0 y -h(x)>=0

No obstante la interfaz MAXIMA para COBYLA admite restricciones de igualdad transformandolas 
luego internamente a pares de desigualdades.

El algoritmo hace uso de aproximaciones lineales, tanto de la funcion objetivo como de las 
restricciones; tales aproximaciones se hacen mediante interpolacion lineal de n+1 puntos 
en el espacio de variables. Los puntos de interpolacion son los vertices de un simplex, que es la 
extension a dimension n de un triangulo en dimension 2. Asi un simplex en dimension 3 seria un tetraedro. 

Hay un parametro (RHO) que controla el tamaño del simplex y se reduce automaticamente desde una valor 
inicial (RHOBEG) a un valor final (RHOEND). Para cada valor del parametro RHO la subrutina intenta 
alcanzar un vector adecuado de variables para el tamaño actual, reduciendose entonces RHO hasta alcanzar 
el valor de RHOEND. Por eso, tanto a RHOBEG como a RHOEND se les deben asignar valores razonables, aunque esto
requeriria cierto trabajo empirico previo.

La rutina trata cada restriccion individualmente, en lugar de considerarlas todas dentro de una funcion conjunta de penalizacion.

Antes del uso de las funciones de este paquete, se debe cargar la libreria correspondiente:

	load(fmin_cobyla);

Podemos utilizar la funcion fmin_cobyla de dos formas:

    fmin_cobyla(FObj, X, X0)
	  
  	  fmin_cobyla(FObj, X, X0, optional_args)

La funcion devuelve una aproximacion del valor minimo de la funcion FObj respecto de las variables X, 
sujeta a un conjunto opcional de restricciones. X0 es una lista que contiene la aproximacion inicial en X.  
FObj debe ser una expresion ordinaria, no valen nombres de funciones ni expresiones tipo lambda 
(ver manual de MAXIMA). X es un vector de variables de n componentes.
    
El parametro optional_args hace referencia a argumentos adicionales, que se especifican en la forma 

    symbol = value

Los argumentos opcionales que se reconocen son:

    ∙constraints

        Lista de restricciones en forma de desigualdades e igualdades que debe satisfacer X. 
        Las restricciones de desigualdad deben ser de la forma g(X) >= h(X) o g(X) <= h(X). 
        Las restricciones de igualdad deben ser de la forma g(X) = h(X).

    ∙  rhobeg

        Valor inicial de la variable interna RHO, que controla el tamaño del simplex. Su valor por defecto es 1.0.

    ∙  rhoend

        El valor final deseado para el parametro RHO. Es aproximadamente la precision de las variables. 
        Su valor por defecto es 1e-6.

    ∙  iprint

        Nivel de informacion durante el proceso de optimizacion. Los posibles valores son: 
            Valor 0: Sin informacion de salida, es el valor por defecto.
            Valor 1: Resumen al final de todos los calculos.
            Valor 2: Se van mostrando los nuevos valores de RHO y SIGA, incluyendo el vector de variables.
            Valor 3: Como en 2, pero la informacion se muestra cuandos se calcula F(X).

    ∙  maxfun

        Numero maximo de evaluaciones de la funcion, el proceso termina si se supera este valor que por defecto es 1000.


    El resultado que devuelve fmin_cobyla es un vector cuyas componentes estan descritas en la siguiente tabla:
    
    . 	Los valores de las variables en las que se alcanza el valor minimo. 	Es una lista de elementos de la forma 
      var = value para cada una de las variables listadas en X

    . 	El valor minimo de la funcion objetivo.

    . 	El numero de evaluaciones de la funcion.

    . 	Codigo de respuesta con los siguientes significados:
        Valor 0: No ha habido errores durante el proceso.
        Valor 1: Se ha alcanzado el numero maximo de evaluaciones de la funcion.
        Valor 2: Los errores de redondeo han impedido el avance del proceso.

EJEMPLO 6: Resuelve el siguiente problema con MAXIMA

    Minimizar                  	xy
    Sujeto a                     1-x^2-y^2 >= 0

load (fmin_cobyla);
fobjetivo(x,y) := x*y;
x0: [1, 1];
fmin_cobyla(fobjetivo(x,y),[x,y], x0, constraints=[x^2+y^2<=1], iprint=1);

Tambien es posible definir la funcion como una variable de esta forma

	fobjetivo : x*y;

en este caso la llamada a la funcion fmin_cobyla seria:

	fmin_cobyla(fobjetivo,[x,y], x0, constraints=[x^2+y^2<=1], iprint=1);

### Funcion augmented_lagrangian_method

La funcion augmented_lagrangian_method minimiza el llamado "lagrangiano aumentado" (o penalizado) 
haciendo uso del algoritmo LBFGS, se utiliza para problemas del tipo:

    Minimizar     	f(x)
    Sujeto a	        h(x) = 0
    
    que se transforma en un problema sin restricciones de la forma

    min f(x) + l*h(x)^2

Como el sumatorio esta formado por terminos positivos, para cualquier valor de x que no sea factible
(es decir h(x)≠0) contribuira con un valor positivo y sera penalizado frente a puntos que si son factibles (h(x)=0).

Antes de su uso, se debe cargar la libreria correspondiente:

	load(augmented_lagrangian);

El uso de augmented_lagrangian_method viene dado en la siguiente lista:


∙ augmented_lagrangian_method(FObj, X, R, X0)

∙ augmented_lagrangian_method(FObj, X, R, X0, optional_args)

∙ augmented_lagrangian_method([FObj, Grad], X, R, X0)

∙ augmented_lagrangian_method([FObj, Grad], X, R, X0, optional_args)


La funcion devuelve una aproximacion del valor minimo de la funcion objetivo FObj respecto de las variables 
dadas en el vector X y manteniendo las restricciones definidas en R igual a cero. La lista X0 contiene 
las aproximaciones iniciales para X.

Si Grad esta presente en la llamada a la funcion, se interpreta como el gradiente de la funcion FObj 
respecto de las variables X y estara representado como una lista de tantas expresiones como variables tenga X. 
Si el argumento Grad no esta, se calcula de forma automatica.

El argumento optional_args hace referencia a argumentos adicionales, los cuales se especifican de la forma 

symbol = value

Estos argumentos son:

    niter            : 	Numero de iteraciones del algoritmo	
    lbfgs_tolerance  : 	Tolerancia del algoritmo	
    iprint           : Lista de dos enteros que controlan la frecuencia de mensajes durante el proceso
                       asi como la cantidad de informacion en cada mensaje (ver funcion lbfgs).
    lambda           : Valor inicial de lambda que sera utilizado para calcular el lagrangiano aumentado.

EJEMPLO 7: Resuelve el siguiente problema con MAXIMA y la funcion augmented_lagrangian_method

    Minimizar    x^2 + y^2
    Sujeto a      	 x + y = 1

Solucion:

load (augmented_lagrangian);
FObj(x,y) := x^2 + 2*y^2;
R : [x+y-1];
x0: [1, 1];
augmented_lagrangian_method(FObj(x,y), [x,y], R, x0, iprint=[-1,0]);

Como en las funciones anteriores, tambien es posible definir la funcion como una variable:

	FObj : x^2+2*y ^2;

en este caso la llamada a la funcion augmented_lagrangian_method seria:

	augmented_lagrangian_method(FObj, [x,y], R, x0, iprint=[-1,0]);

EJEMPLO 8: Resuelve el siguiente problema con MAXIMA y la funcion augmented_lagrangian_method

    Minimizar      	y
    Sujeto a         x^2 + z^2 = 10
                            x + y + z = 3

Solucion:

load (augmented_lagrangian);
FObj(x,y,z) := y;
R : [x ^2+z ^2-10, x+y+z-3];
x0: [1, 1, 1];
augmented_lagrangian_method(FObj(x,y,z), [x,y,z], R, x0, iprint=[-1,0]);

### Optimizacion lineal

El paquete simplex utiliza el algoritmo del mismo nombre para resolver problemas de programacion lineal. 

Para usar las funciones de esta seccion, hay que cargar en primer lugar el paquete correspondiente:

load(simplex);

### Funcion linear_program

    La funcion linear_program es una implementacion del algoritmo simplex y se utiliza para resolver problemas 
de la forma

    Minimizar c_{1} x_{1} + c_{2} x_{2} + ..... + c_{n} x_{n}
    Sujeto a:

        A x = b
        x_{k} >=0
 
donde A es una matriz  de dimensión m x n,  b es un vector de R^m y c es un vector de R^n son vectores.

La instruccion que resuelve ese tipo de problemas es

    	linear_program(A, b, c)

El argumento A es una matriz y los argumentos b y c son vectores, definidos según la sintaxis de MAXIMA.

Si el problema tiene solucion, entonces el valor que devuelve linear_program 
es una lista que contiene un vector x, con la solucion optima y el valor optimo. 

Por el contrario, si el problema no esta acotado devuelve el mensaje 

"Problem not bounded!" 

y si el problema no es factible devuelve el mensaje 

"Problem not feasible!".

EJEMPLO 9: Resuelve el siguiente problema con MAXIMA y la funcion linear_program

    min x_{1} - 2 x_{2}
    Sujeto a
        x_{1} + x_{2} - x_{3} = 1
        2x_{1} - 3x_{2} - x_{4} = 1
        4x_{1} - 5x_{2} = 6

        	x_{j}≥0


Solucion:

load (simplex);
A : matrix([1,1,-1,0], [2,-3,0,-1], [4,-5,0,0]);
b : [1, 1,6];
c : [1, -2, 0, 0];
linear_program(A,b,c);

### Funciones minimize_lp y maximize_lp

La funcion minimize_lp se utiliza para resolver problemas lineales de la forma:


    min    c_{1} x_{1} + c_{2} x_{2} + ..... + c_{n} x_{n}
    Sujeto a
            A_{e} x = b_{e}
   
            A_{g} x >= b_{g}

            A_{l} x <= b_{l} 
            
            x_{k} >= 0


Donde A_{e}, A_{g}, A_{l} son matrices y b_{e}, b_{g} y b_{l} son vectores.

El uso de esta funcion seria de la forma:

    	minimize_lp(obj, con, [pos])

El objetivo es minimizar la funcion lineal (c_{1} x_{1} + c_{2} x_{2} + ..... + c_{n} x_{n}), sujeta a las restricciones 
tambien lineales indicadas en la variable con, estas restricciones se representan como una lista de ecuaciones 
o inecuaciones lineales.

A la hora de resolver el problema en las inecuaciones estrictas se reemplaza > por ≥ y < por ≤. 
El argumento pos es opcional y contiene una lista con las variables de decision que deben ser positivas, 
en otro caso serian libres. 

Si el minimo existe, minimize_lp devuelve una lista que contiene 
el valor minimo de la funcion objetivo y una lista de valores para las variables de decision 
con los que se alcanza el minimo. Si el problema no esta acotado, devuelve el mensaje 

"Problem not bounded!" 

y si el problema no es factible, devuelve el mensaje 

"Problem not feasible!".
    

Si no se especifica lo contrario, las variables de decision serian libres, para indicar que 
todas las variables de decision son no negativas, se le puede asignar el valor true a la variable 
del sistema nonegative_lp. 

Si solo algunas de las variables de decision fueran positivas entonces 
se pueden incluir restricciones de la forma x≥0, aunque es mas eficiente utilizar el argumento pos.

La funcion maximize_lp es identica en funcionamiento a minimize_lp pero en este caso el objetivo del 
problema seria maximizar.

Variables opcionales

  ∙  nonegative_lp:

        Si nonegative_lp es true, entonces todas las variables de decision pasadas a minimize_lp 
        y a maximize_lp se suponen positivas. El valor por defecto es false.    

EJEMPLO 10: Resuelve mediante minimize_lp o maximize_lp, cuando corresponda:

A)  Minimizar x + y

    sujeto a

        3x +  y = 0
         x + 2y >= 2

minimize_lp(x+y, [3*x+y=0, x+2*y>2]);

B) Minimizar x + y

   Sujeto a

        3x +  y > 0
         x + 2y > 2
         x,y >= 0

minimize_lp(x+y,[3*x+y>0,x+2*y>2]),nonegative_lp=true;

/* O bien */

minimize_lp(x+y,[3*x+y>0,x+2*y>2],[x,y]);

C)  Minimizar 2x + y

    Sujeto a
         x +  y <= 2 
        -x +  y <= 3
        3x + 2y <= 10
        x libre
        y >= 0

maximize_lp(2*x+y,[x+y<=2,-x+y<=3,3*x+2*y<=10],[y]);┊

D) Minimizar -10x + 4y

   Sujeto a
       
        x -  y <= 2
       5x - 2y <= 16
       x,y >= 0

minimize_lp(-10*x+4*y,[x-y<=2, 5*x-2*y<=16]),nonegative_lp=true;

/* O tambien */

minimize_lp(-10*x+4*y,[x-y<=2, 5*x-2*y<=16],[x,y]);

E) Maximizar 2x + 5y

   Sujeto a
       
       x + y >= 4
       x     >= 2
       x, y >= 0

maximize_lp(2*x+5*y,[x+y>=4, x>=2]),nonegative_lp=true;

/* O tambien */

maximize_lp(2*x+5*y,[x+y>=4, x>=2],[x,y]);

F) Minimizar x+y

    Sujeto a
    
         x +  y <= 1
        4x + 2y >= 6
        x, y >= 0

minimize_lp(x+y, [x+y<=1,4*x+2*y>=6]), nonegative_lp=true;

/* O tambien */

minimize_lp(x+y, [x+y<=1,4*x+2*y>=6],[x,y]);

## MÉTODOS ANALÍTICOS

Es posible utilizar los métodos analíticos para resolver problemas de optimización utilizando
los puntos que cumplen las condiciones de KKT.

Definimos la función objetivo como una variable MAXIMA:

fobjetivo : x+y;

Definimos las restricciones del problema en forma estándar (h(x)=0 o g(x)<=0), 
tambien como variables:

r1 : x^2+y^2-1;
r2 : x^2/4 + 2*y^2 - 1;

Definimos la función lagrangiana

lagrangiano : fobjetivo + m1*r1 + m2*r2;

Planteamos el conjunto de todas las ecuaciones

Condición Estacionaria

ec1: diff(lagrangiano,x)=0;
ec2: diff(lagrangiano,y)=0;

Condiciones de holgura

ec3: m1*r1 = 0;
ec4: m2*r2 = 0;

Resolvemos el sistema en todas las variables (decisión y multiplicadores):

solucion1:solve([ec1,ec2,ec3,ec4],[x,y,m1,m2]);

El problema es que MAXIMA no nos da todas las soluciones y el resto hay que hacerlo "a mano"

solucion2:solve([r1,r2],[x,y]);
solucion3:solve([ec1,ec2],[m1,m2]);

Comprobamos la factibilidad de cada solución

for sol in solucion1 do  (ldisp(sol), factible1:subst(sol,r1), factible2:subst(sol,r2), factible:is(factible1<=0 and factible2 <=0), ldisp(factible) )$
for sol in solucion2 do  (ldisp(sol), factible1:subst(sol,r1), factible2:subst(sol,r2), factible:is(factible1<=0 and factible2 <=0), ldisp(factible) )$

Vemos que los 4 primeros son factibles, mientras que los 4 últimos son infactibles.

Veamos el signo de los multiplicadores de estos 4 puntos:

/* for sol in solucion2 do (dummy:subst(sol,solucion3),ldisp(dummy),  signo:rhs(dummy[1][1])*rhs(dummy[1][2]), compatible:is(signo >=0),ldisp(compatible));*/
for sol in solucion2 do (dummy:subst(sol,solucion3),ldisp(dummy), ldisp(is(rhs(dummy[1][1])<=0)), ldisp(is(rhs(dummy[1][2])<=0)));

Se ha comparado cada valor con el 0. 

Para la primera pareja de valores el resultado en ambos valores es false y por tanto ambos valores son 
negativos y el punto correspondiente puede ser de máximo.

Las dos siguientes parejas son de distinto signo.

La última pareja de valores nos arroja un resultado de true en ambos casos, luego la pareja de multiplicadores es positiva y el punto puede ser de mínimo.

El estudio del Hessiano y el espacio tangente para estos puntos es un poco largo.

Hessiano : hessian(fobjetivo,[x,y])+ m1 * hessian(r1,[x,y])+m2 * hessian(r2,[x,y]);

h1:Hessiano,solucion1[1];
eigenvalues(h1);

h2:Hessiano,solucion2[2];
eigenvalues(h2);

## Ejercicios para practicar.

PROBLEMA 1. Halla los extremos locales y globales de f(x)=x³-12x+3 en el intervalo [-4,4].

PROBLEMA 2. Estudia los extremos relativos y absolutos de las siguientes funciones:
    a) f(x,y) = x^2 + 3xy - y^2	
    b) f(x,y) = x^3 + 3xy^2 - 15x - 12y	
    c) f(x,y)= senx + seny + cos(x+y)       0 < x,y < 2π

PROBLEMA 3. Determina los puntos de la elipse que se obtienen al cortar el cilindro de 
ecuacion x^2+y^2=1 con el plano x+y+z=1, cuyas distancias al origen sean respectivamente 
maxima y minima.

PROBLEMA 4. Resuelve el problema

Optimizar x^2 + y^2 + z^2
Sujeto a
    x^2 + y^2 + z^2 + xy + yz + zx - 1 = 0
    	x+2y-3z=0

PROBLEMA 5. Halla la minima distancia entre la recta x+y=4 y la circunferencia x^2 + y^2 = 1.

PROBLEMA 6 Resuelve el siguiente problema:

Optimizar -x^2 + 2y + z
Sujeto a
    x + 2y - z=0
    	2y + z = 0

PROBLEMA 7: Resuelve el siguiente problema

Optimizar  	-x + 2y +z
Sujeto a
        	x + 2y - z = 0
        	2y + z = 0

PROBLEMA 8: Resuelve el siguiente problema

Minimizar 	 x + y
Sujeto a
    3x + 2y > 2
    	 x + 4y > 3

PROBLEMA 9: Resuelve el siguiente problema

Minimizar 	 x + 2y		
Sujeto a
    	-5x 	+ y 	 =	7
     x + y >= 26	
    x >= 3
    y >= 4

PROBLEMA 10: Resuelve los siguientes problemas

Minimizar 	 x + y
Sujeto a
    	x + 2y >=3
    	x,  y >= 0 

Minimizar 	x + y
Sujeto a
    	x + 2y = 3
    x,y >= 0

Minimizar 	 x + y
Sujeto a
    x + 2y <= 3
    x, y >= 0

Minimizar 	x + y
Sujeto a
    x + 2y >= 3
    -1 <= x <= 2
    -1 <= y <= 1