# Librería de cálculo simbólico
import sympy as sym
# Para imprimir en formato TeX
from sympy import init_printing; init_printing(use_latex='mathjax')

sym.var('x', real = True)
f = x**2
f

df = sym.diff(f, x)
df

x_c = sym.solve(df, x)
x_c[0]

"""
Explanation: Optimización de funciones escalares diferenciables con Sympy


Mediante optimización se obtienen soluciones elegantes tanto en teoría como en ciertas aplicaciones. 
La teoría de optimización usa elementos comenzando con cálculo elemental y álgebra lineal básica, y luego se extiende con análisis funcional y convexo.
Las aplicaciones en optimización involucran ciencia, ingeniería, economía, finanzas e industria.
El amplio y creciente uso de la optimización lo hace escencial para estudiantes y profesionales de cualquier rama de la ciencia y la tecnología.


Referencia
- http://www.math.uwaterloo.ca/~hwolkowi//henry/reports/talks.d/t06talks.d/06msribirs.d/optimportance.shtml
Algunas aplicaciones son:

Ingeniería
Encontrar la composición de equilibrio de una mezcla de diferentes átomos.
Planeación de ruta para un robot (o vehículo aéreo no tripulado).
Distribución óptima de recursos.
Distribución de rutas de vuelo.
Encontrar una dieta óptima.
Optimización financiera
Administración de riesgos.

En esta clase veremos aspectos básicos de optimización. En específico, veremos cómo obtener máximos y mínimos de una función escalar de una variable (como en cálculo diferencial).
Basamos todos los resultados en los siguientes teoremas:
1. Teorema de Fermat (análisis)
Si una función $f(x)$ alcanza un máximo o mínimo local en $x=c$, y si la derivada $f'(c)$ existe en el punto $c$, entonces $f'(c) = 0$.
Ejemplo
Sabemos que la función $f(x)=x^2$ tiene un mínimo global en $x=0$, pues
$$f(x)=x^2\geq0,\qquad\text{y}\qquad f(x)=x^2=0 \qquad\text{si y solo si}\qquad x=0.$$
End of explanation
"""


import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline

f_num = sym.lambdify([x], f, 'numpy')
x_vec = np.linspace(-5, 5, 100)

plt.plot(x_vec, f_num(x_vec))
plt.xlabel('$x$')
plt.ylabel('$x^2$')
plt.show()

"""
Explanation: Veamos la gráfica...
End of explanation
"""


def f(x):
    return x**2

f_sym = f(x)
f_sym

df = sym.diff(f(x), x)
df

x_c = sym.solve(df, x)
x_c[0]

"""
Explanation: Otra manera de hacer lo anterior
End of explanation
"""


f = x**2
#d2f = sym.diff(f, x, x)
d2f = sym.diff(f, x, 2)
d2f

d2f>0

"""
Explanation: El converso del teorema anterior no es cierto.
Actividad
Considere $g(x)=x^3$. 
- Usando sympy, muestre que $g'(0)=0$.
- Sin embargo, descartar que $x=0$ es un extremo de $g(x)$ viendo su gráfica.
2. Criterio de la segunda derivada
Sea $f(x)$ una función tal que $f’(c)=0$ y cuya segunda derivada existe en un intervalo abierto que contiene a $c$.
- Si $f’’(c)>0$, entonces $f(c)$ es un mínimo relativo.
- Si $f’’(c)<0$, entonces $f(c)$ es un máximo relativo.
- Si $f’’(c)=0$, entonces el criterio no decide. 
Ejemplo
Mostrar, usando sympy, que la función $f(x)=x^2$ tiene un mínimo relativo en $x=0$.
Ya vimos que $f'(0)=0$. Notemos que:
End of explanation
"""


f = x**2-6*x
f

df = sym.diff(f, x)
df

x_c = sym.solve(df, x)
x_c

"""
Explanation: Por tanto, por el criterio de la segunda derivada, $f(0)=0$ es un mínimo relativo (en efecto, el mínimo global).
Actividad
¿Qué pasa con $g(x)=x^3$ al intentar utilizar el criterio de la segunda derivada? (usar sympy).
3. Método para determinar extremos absolutos de una función continua y=f(x) en [a,b]

Evaluar $f$ en los extremos $x=a$ y $x=b$.
Determinar todos los valores críticos $c_1, c_2, c_3, \dots, c_n$ en $(a,b)$.
Evaluar $f$ en todos los valores críticos.
El más grande y el más pequeño de los valores de la lista $f(a), f(b), f(c_1), f(c_2), \dots, f(c_n)$ son el máximo absoluto y el mínimo absoluto, respectivamente, de f en el intervalo [a,b].

Ejemplo
Determinar los extremos absolutos de $f(x)=x^2-6x$ en $\left[0,5\right]$.
Obtenemos los puntos críticos de $f$ en $\left[0,5\right]$:
End of explanation
"""


f.subs(x, 0), f.subs(x, 5), f.subs(x, x_c[0])

"""
Explanation: Evaluamos $f$ en los extremos y en los puntos críticos:
End of explanation
"""


f_num = sym.lambdify([x], f, 'numpy')
x_vec = np.linspace(0, 5, 100)

plt.figure(figsize=(8,6))
plt.plot(x_vec, f_num(x_vec), 'k', label = '$y=f(x)$')
plt.plot([0], [0], '*r', label = '$(0,0=\max_{0\leq x\leq 5} f(x))$')
plt.plot([3], [-9], '*b', label = '$(3,-9=\min_{0\leq x\leq 5} f(x))$')
plt.legend(loc='best')
plt.xlabel('x')
plt.show()

"""
Explanation: Concluimos que el máximo absoluto de $f$ en $\left[0,5\right]$ es $0$ y se alcanza en $x=0$, y que el mínimo absoluto es $-9$ y se alcanza en $x=3$.
End of explanation
"""


sym.var('x y')
x, y

def f(x, y):
    return x**2 + y**2

dfx = sym.diff(f(x,y), x)
dfy = sym.diff(f(x,y), y)
dfx, dfy

xy_c = sym.solve([dfx, dfy], [x, y])
xy_c

x_c, y_c = xy_c[x], xy_c[y]

d2fx = sym.diff(f(x,y), x, 2)
d2fy = sym.diff(f(x,y), y, 2)
dfxy = sym.diff(f(x,y), x, y)

Jf = sym.Matrix([[d2fx, dfxy], [dfxy, d2fy]])
Jf.eigenvals()

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')

x = np.linspace(-2, 2, 100)
y = x
X, Y = np.meshgrid(x, y)

ax.plot_surface(X, Y, f(X, Y))
ax.plot([0], [0], [0], '*r')

"""
Explanation: Actividad
Determinar los valores extremos absolutos de $h(x)=x^3-3x$ en $\left[-2.2,1.8\right]$, usando sympy. Mostrar en una gráfica.
En varias variables...
El procedimiento es análogo. 
Si una función $f:\mathbb{R}^n\to\mathbb{R}$ alcanza un máximo o mínimo local en $\boldsymbol{x}=\boldsymbol{c}\in\mathbb{R}^n$, y $f$ es diferenciable en el punto $\boldsymbol{x}=\boldsymbol{c}$, entonces $\left.\frac{\partial f}{\partial \boldsymbol{x}}\right|_{\boldsymbol{x}=\boldsymbol{c}}=\boldsymbol{0}$ (todas las derivadas parciales en el punto $\boldsymbol{x}=\boldsymbol{c}$ son cero).
Criterio de la segunda derivada: para ver si es máximo o mínimo, se toma la segunda derivada (matriz jacobiana) y se verifica definición negativa o positiva, respectivamente.
Si se restringe a cierta región, hay ciertas técnicas. La más general, pero también la más compleja es la de multiplicadores de Lagrange.
End of explanation
"""


