El ciclo for en Python

El ciclo for es una estructura cíclica cuya finalidad es ejecutar un código una cierta cantidad de veces. De manera general, en Python, la sintaxis de for es:

for k in iterable:
    # hacer algo

Donde iterable es justamente eso, un iterable de Python que puede ser desde una lista hasta un string. Para la i-ésima iteración la variable de ciclo k adoptará el valor en la i-ésima posición de iterable.

A continuación vemos un ejemplo muy básico:

In [5]:
numeros = [10,20,30,40,50,60]
for numero in numeros:
    print(numero)
10
20
30
40
50
60

Observe que lo único que se hace es imprimir cada elemento que compone la lista numeros.

Revisemos ahora el siguiente ejemplo:

In [7]:
numeros = [10,20,30,40,50,60]
for numero in numeros:
    print(numero, numero**2)
10 100
20 400
30 900
40 1600
50 2500
60 3600

Note que la única diferencia radica, en que además del propio número, también imprimos su cuadrado.

Naturalmente, no sólo podemos "recorrer" listas, sino cualquier elemento que sea iterable, como es el caso de los strings.

In [9]:
nombres = ["Ana", "Lucas", "Catalina", "Javier"]
for nombre in nombres:
    print("Hola " + nombre)
Hola Ana
Hola Lucas
Hola Catalina
Hola Javier

En ocasiones tenemos iterables dentro de otro iterable. Es posible acceder a cada elemento del sub-iterable utilizando múltiple variables de ciclo. Por ejemplo:

In [22]:
puntos = [(0,0), (3,5), (-1,4), (10,2)]
for x,y in puntos:
    print("Px = {0},  Py ={1}".format(x,y))
Px = 0,  Py =0
Px = 3,  Py =5
Px = -1,  Py =4
Px = 10,  Py =2

Una utilidad muy común es la de la función enumerate:

In [24]:
notas = [10,8,10,9,10,7]
for k,nota in enumerate(notas):
    print("La nota en posición {0} es: {1}".format(k,nota))
La nota en posición 0 es: 10
La nota en posición 1 es: 8
La nota en posición 2 es: 10
La nota en posición 3 es: 9
La nota en posición 4 es: 10
La nota en posición 5 es: 7

De manera muy breve, la variable de ciclo k devuelve la posición del elemento para cada iteración.

Algunos ejemplos utilizando ciclos for

Una tabla de multiplicar

In [27]:
n = 2
for k in range(1,11):
    print("{0} x {1} = {2}".format(n,k,n*k))
2 x 1 = 2
2 x 2 = 4
2 x 3 = 6
2 x 4 = 8
2 x 5 = 10
2 x 6 = 12
2 x 7 = 14
2 x 8 = 16
2 x 9 = 18
2 x 10 = 20

Un triángulo de caracteres

In [30]:
n = 10
for k in range(1,n+1):
    print("*"*k)
*
**
***
****
*****
******
*******
********
*********
**********

Contando vocales

In [46]:
frase = "Anita lava la tina"
k = 0
for letra in frase.lower():
    if letra in ("a","e","i","o","u"):
        k += 1

print("Hay {0} vocales".format(k))
Hay 8 vocales

La suma de los primeros cien enteros positivos

In [47]:
suma = 0
for num in range(1,101):
    suma += num

print(suma)
5050

Gráficas de pastel con Matplotlib

Las gráficas de barras nos sirven para representar porcentajes y proporciones. En Python podemos utilizar la librería Matplotlib para desarrollar este tipo de gráficas.

Matplotlib dispone de la función pie, cuya sintaxis depende del grado de personalización y control que se requiera sobre la gráfica de pastel a dibujar.

Para ejemplificar el uso de esta función vamos a suponer que se tienen los siguientes datos sobre algunas personas que tienen cierta cantidad de manzanas en su poder:

Nombre Manzanas
Ana 20
Juan 10
Diana 25
Catalina 30

Para representar el porcentaje del total del cual dispone cada uno, podemos trazar una gráfica de pastel. Para ello realizamos lo siguiente:

In [20]:
import matplotlib.pyplot as plt

manzanas = [20,10,25,30]
nombres = ["Ana","Juan","Diana","Catalina"]
plt.pie(manzanas, labels=nombres)
plt.show()

Observe que lo primero que hacemos es importar la librería matplotlib, enseguida, en utilizando listas definimos los nombres y el número de manzanas correspondientes. Luego, la función pie acepta un primer argumento que contiene los valores absolutos de cada ítem, además, de un keyword argument labels que contiene las etiquetas correspondientes.

El aspecto achatado de la gráfica se puede solucionar utilizando la función axis.

In [31]:
import matplotlib.pyplot as plt

manzanas = [20,10,25,30]
nombres = ["Ana","Juan","Diana","Catalina"]
plt.pie(manzanas, labels=nombres)
plt.axis("equal")
plt.show()

Indicando el porcentaje

El porcentaje correspondiente a cada ítem se puede indicar mediante el argumento autopct:

In [32]:
import matplotlib.pyplot as plt

manzanas = [20,10,25,30]
nombres = ["Ana","Juan","Diana","Catalina"]
plt.pie(manzanas, labels=nombres, autopct="%0.1f %%")
plt.axis("equal")
plt.show()

Cambiando los colores

Los combinación de colores se puede especificar de manera manual, pasando una lista de color en formato hexadecimal o RGB.

In [88]:
import matplotlib.pyplot as plt

manzanas = [20,10,25,30]
nombres = ["Ana","Juan","Diana","Catalina"]
colores = ["#EE6055","#60D394","#AAF683","#FFD97D","#FF9B85"]
plt.pie(manzanas, labels=nombres, autopct="%0.1f %%", colors=colores)
plt.axis("equal")
plt.show()

Los colores también se pueden determinar y autocalcular utilizando un mapa de color específico. Enseguida se muestra un ejemplo donde la variación es sobre colores en tonos azules.

In [89]:
import matplotlib.pyplot as plt
from matplotlib import cm
from matplotlib import colors

manzanas = [20,10,25,30]
nombres = ["Ana","Juan","Diana","Catalina"]

normdata = colors.Normalize(min(manzanas), max(manzanas))
colormap = cm.get_cmap("Blues")
colores =colormap(normdata(manzanas))

plt.pie(manzanas, labels=nombres, autopct="%0.1f %%", colors=colores)
plt.axis("equal")
plt.show()

Extrayendo rebanadas

Es posible también segmentar o separar del bloque una o más de las rebanadas de la gráfica de pastel. Para ello se debe pasar una lista o tupla con valores entre 0 y n que indican el desfase respecto al centro, 0 indica ningún desfase y n un desfase equivalente a n*r, donde r es el radio de la gráfica de pastel.

In [102]:
import matplotlib.pyplot as plt

manzanas = [20,10,25,30]
nombres = ["Ana","Juan","Diana","Catalina"]
colores = ["#EE6055","#60D394","#AAF683","#FFD97D","#FF9B85"]
desfase = (0, 0, 0, 0.1)
plt.pie(manzanas, labels=nombres, autopct="%0.1f %%", colors=colores, explode=desfase)
plt.axis("equal")
plt.show()

Analizando datos del CONACYT: De los investigadores

En lo siguiente se extraen algunos datos correspondientes a los investigadores que pertenecen al SNI.

In [16]:
import warnings
warnings.filterwarnings("ignore", 'This pattern has match groups')
import pandas as pd
import matplotlib.pyplot as plt
import matplotlib as mpl
plt.style.use("ggplot")
mpl.rcParams['xtick.labelsize'] = 6 
mpl.rcParams['ytick.labelsize'] = 6 
%matplotlib inline

# Datos tomados de: http://datosabiertos.conacyt.gob.mx/publico/default.aspx
filename = "data/inv_sni.csv"
data = pd.read_csv(filename, encoding="latin1") # leyendo datos
data.drop(data.columns[range(8,17)],axis=1,inplace=True) # Limpiando datos

Generalidades

El número de total de investigadores

In [17]:
NUM_INV = len(data)

Investigadores eméritos

In [18]:
emeritos = data['EMERITO'].value_counts().E
emeritos
Out[18]:
108

Por niveles

In [19]:
niveles = data['NIVEL'].value_counts().to_frame()
niveles
Out[19]:
NIVEL
1 12775
C 4575
2 3964
3 2002

Por género

In [20]:
genero = data["GENERO"].value_counts().to_frame("INVESTIGADORES")
genero["PORCENTAJE"] = 100*genero["INVESTIGADORES"]/(NUM_INV)
genero
Out[20]:
INVESTIGADORES PORCENTAJE
H 14969 64.200549
M 8347 35.799451

Por estados

In [21]:
estados = data['ESTADO'].value_counts()
x = range(len(estados))
plt.bar(x, estados.get_values(), 0.5, align="center")
plt.xticks(x, estados.keys(), rotation="vertical")
plt.xlim(min(x)-0.5, max(x)+0.5);
In [22]:
ninv = sum(estados.values) # No. total de investigadores
prc_df = estados.to_frame("INVESTIGADORES")
prc_df["Porcentaje"] = 100*prc_df['INVESTIGADORES']/ninv
prc_df
Out[22]:
INVESTIGADORES Porcentaje
DISTRITO FEDERAL 7831 33.587819
M?XICO, EDO. DE 1361 5.837444
JALISCO 1197 5.134034
NO DISPONIBLE 1013 4.344842
MORELOS 1008 4.323397
NUEVO LEON 962 4.126099
PUEBLA 881 3.778683
GUANAJUATO 801 3.435557
BAJA CALIFORNIA 715 3.066695
VERACRUZ 701 3.006648
MICHOACAN 687 2.946601
QUERETARO 617 2.646365
SAN LUIS POTOSI 573 2.457645
YUCATAN 552 2.367575
SONORA 475 2.037315
SINALOA 401 1.719923
CHIHUAHUA 381 1.634141
COAHUILA 340 1.458289
HIDALGO 323 1.385374
CHIAPAS 278 1.192365
OAXACA 270 1.158053
BAJA CALIFORNIA SUR 250 1.072271
ZACATECAS 201 0.862106
TAMAULIPAS 196 0.840661
COLIMA 192 0.823504
AGUASCALIENTES 173 0.742012
TABASCO 158 0.677675
DURANGO 155 0.664808
TLAXCALA 145 0.621917
CAMPECHE 131 0.561870
QUINTANA ROO 129 0.553292
NAYARIT 116 0.497534
GUERRERO 102 0.437487

Por instituciones

Las 20 instituciones con mayor número de investigadores

In [34]:
inst = data['INSTITUCION'].value_counts()
n = 20
x = range(len(inst))[:n]
plt.bar(x, inst.get_values()[:n], 0.5, align="center")
plt.xticks(x, inst.keys()[:n], rotation=90)
plt.xlim(min(x)-0.5, max(x)+0.5);

Sobre las instituciones del Tecnológico Nacional de México (TECNM)

Lo subsiguiente presenta el porcentaje de representación de las instituciones correspondientes al Tecnológico Nacional de México: es decir, todos los institutos tecnológicos centralizados y descentralizados, más el CENIDET y CIIDET.

In [24]:
TECNOLOGICOS = sum(inst[inst.index.str.contains("INSTITUTO TECNOLOGICO (SUP|DE)")])
CENIDET = sum(inst[inst.index.str.contains("CENTRO NACIONAL DE INV*")])
CIIDET = sum(inst[inst.index.str.contains("INTERDISCIPLINARIO [\w*\s*\.]* EDUCACION")])
porcentaje_tecnm = 100*(TECNOLOGICOS + CENIDET + CIIDET)/NUM_INV
porcentaje_tecnm
Out[24]:
2.5347400926402472

Los 20 tecnológicos con mayor cantidad de investigadores son:

In [29]:
(inst[inst.index.str.contains("INSTITUTO TECNOLOGICO (SUP|DE)")][:20]).to_frame("No. de Investigadores")
Out[29]:
No. de Investigadores
INSTITUTO TECNOLOGICO DE CELAYA 45
INSTITUTO TECNOLOGICO DE SONORA 43
INSTITUTO TECNOLOGICO DE TIJUANA 37
INSTITUTO TECNOLOGICO DE CIUDAD MADERO 25
INSTITUTO TECNOLOGICO DE MORELIA 23
INSTITUTO TECNOLOGICO DE VERACRUZ 20
INSTITUTO TECNOLOGICO DE TOLUCA 20
INSTITUTO TECNOLOGICO DE AGUASCALIENTES 19
INSTITUTO TECNOLOGICO DE CONKAL 16
INSTITUTO TECNOLOGICO DE ORIZABA 13
INSTITUTO TECNOLOGICO DE DURANGO 13
INSTITUTO TECNOLOGICO DE MERIDA 11
INSTITUTO TECNOLOGICO SUPERIOR DE IRAPUATO 11
INSTITUTO TECNOLOGICO DE TEPIC 11
INSTITUTO TECNOLOGICO DE LA LAGUNA 11
INSTITUTO TECNOLOGICO DE OAXACA 10
INSTITUTO TECNOLOGICO DE ROQUE GUANAJUATO 9
INSTITUTO TECNOLOGICO DE TUXTLA GUTIERREZ 9
INSTITUTO TECNOLOGICO DE CULIACAN 9
INSTITUTO TECNOLOGICO DE SALTILLO 9
In [ ]:
 

Analizando datos del CONACYT: De los becarios nacionales

En este pequeño post vamos a analizar (en realidad no, sólo a trazar grafiquitas con algunos datos de interés) datos del Consejo Nacional de Ciencia y Tecnología (CONACYT), que hasta hace unos días me listaba entre sus becarios nacionales. Pero bueno, ¿qué es el CONACYT?, según wikipedia:

El Consejo Nacional de Ciencia y Tecnología (Conacyt) es un organismo público descentralizado del gobierno federal mexicano dedicado a promover y estimular el desarrollo de la ciencia y la tecnología en ese país. Tiene la responsabilidad oficial para elaborar las políticas de ciencia y tecnología nacionales.

Por medio del Conacyt es posible para los estudiantes conseguir apoyo económico a fin de realizar estudios de posgrado (maestría o doctorado) en universidades con reconocida excelencia académica dentro y fuera del país.

En resumen el CONACYT posibilita que muchos de nosotros (jóvenes mexicanos) ingresemos a un programa de posgrado, mismo que debe pertenecer a un padrón de posgrados de calidad.

En lo subsiguiente vamos a tomar algunos datos públicos del CONACYT referentes a los becarios nacionales y trataremos de segmentar y/o clasificar la información obtenida, básicamente en lo que concierne

In [1]:
import pandas as pd
import matplotlib.pyplot as plt
import matplotlib as mpl
plt.style.use("ggplot")
mpl.rcParams['xtick.labelsize'] = 6 
mpl.rcParams['ytick.labelsize'] = 6 
%matplotlib inline

# Datos tomados de: http://datosabiertos.conacyt.gob.mx/publico/default.aspx
filename = "data/becas_nacionales.csv"
data = pd.read_csv(filename)

Por género

In [33]:
genero = data['Genero'].value_counts()
print("Mujeres: {0}".format(genero.get_value("F")))
print("Hombres: {0}".format(genero.get_value("M")))
Mujeres: 10954
Hombres: 12646

Por estados

In [34]:
estados = data['Entidad Federativa'].value_counts()
x = range(len(estados))
plt.bar(x, estados.get_values(), 0.5, align="center")
plt.xticks(x, estados.keys(), rotation="vertical")
plt.xlim(min(x)-0.5, max(x)+0.5);

Áreas de conocimientos

In [35]:
acon = data['Area de Conocimiento'].value_counts()
x = range(len(acon))
plt.bar(x, acon.get_values(), 0.3, align="center")
plt.xticks(x, acon.keys(), rotation=90)
plt.xlim(min(x)-0.5, max(x)+0.5);

Por instituciones

In [36]:
inst = data['Institucion'].value_counts()
n = 10
x = range(len(inst))[:n]
plt.bar(x, inst.get_values()[:n], 0.5, align="center")
plt.xticks(x, inst.keys()[:n], rotation=90)
plt.xlim(min(x)-0.5, max(x)+0.5);
In [31]:
data['Nivel'].value_counts()
Out[31]:
2. MAE    17034
1. DOC     5270
3. ESP     1296
Name: Nivel, dtype: int64

Máximos y mínimos de una función

Supongamos que en $x=x_1$ la derivada de la función $y=f(x)$ se reduce a cero, es decir, $f'(x)=0$. Admitamos, además que existe la segunda derivada, $f''(x)$, y es continua sobre cierta vecindad del punto $x_1$. Para este caso es válido el siguiente teorema:

Teorema 1. Si $f'(x_1)=0$, entonces en $x=x_1$ la función tiene máximo cuando $f''(x_1)<0$, y, un mínimo cuando $f''(x_1)>0$.

De acuerdo al teorema anterior, para caracterizar los puntos críticos de una función $f(x)$ es necesario utilizar el criterio de la segunda derivada. Entendiendo que los puntos críticos se obtiene de resolver la ecuación $f'(x)=0$.

En Python es posible realizar cálculo simbólico mediante la librería SymPy.


Ejemplo 1. Calcular y caracterizar los puntos críticos de la función $f(x) = x^3 - x$

Primero, y como siempre, importamos los módulos necesarios.

In [48]:
from sympy import *
from sympy.plotting import plot
from sympy.abc import x
init_printing()

Definimos la función con la cual operaremos:

In [22]:
f = x**3 - x

Acto seguido, calculamos la primera y segunda derivada de $f(x)$ mediante la función diff:

In [49]:
df = diff(f, x) # 1era. derivada
d2f = diff(f,x,2) # Segunda derivada
df, d2f
Out[49]:
$$\left ( 3 x^{2} - 1, \quad 6 x\right )$$

Los puntos críticos los calculamos resolviendo la ecuación $f'(x)=0$

In [50]:
pc = solve(Eq(df, 0))
pc
Out[50]:
$$\left [ - \frac{\sqrt{3}}{3}, \quad \frac{\sqrt{3}}{3}\right ]$$

Para determinar si se trata de un mínimo o máximo, utilizamos el criterio de la segunda derivada, sustituyendo los puntos críticos en $f''(x)$, es decir:

In [51]:
d2f.subs(x,pc[0]) # Primer punto crítico
Out[51]:
$$- 2 \sqrt{3}$$
In [52]:
d2f.subs(x,pc[1]) # Segundo punto crítico
Out[52]:
$$2 \sqrt{3}$$

Con esto determinamos, acorde al teorema 1, que el primer punto crítico $\left(\frac{-\sqrt{3}}{3}\right)$ es un mínimo y el segundo $\left(\frac{\sqrt{3}}{3}\right)$ un máximo. Podemos comprobarlo trazando la gráfica correspondiente:

In [54]:
plot(f, (x, -1, 1))
Out[54]:
<sympy.plotting.plot.Plot at 0x94ca890>

Podemos automatizar un poco este proceso definiendo una función que realice el procedimiento descrito.

In [84]:
def maxminf(f):
    """ Calcula los máximos y mínimos de una función f(x) """
    df = diff(f, x) # 1era. derivada
    d2f = diff(f, x, 2) # 2da. derivada
    pcs = solve(Eq(df,0)) # puntos críticos
    for p in pcs:
        if d2f.subs(x,p)>0: 
            tipo="Min"
        elif d2f.subs(x,p)<0: 
            tipo="Max"
        else: 
            tipo="Indefinido"
        print("x = %f (%s)"%(p,tipo))
In [85]:
maxminf(x**3 - x)
x = -0.577350 (Max)
x = 0.577350 (Min)

La función anterior puede arrojar errores cuando se tienen raíces complejas, lo cual podría considerarse mediante una estructura try-except o bien verificando si el punto crítico es un valor real.

wxdx, una GUI para cálculo elemental

En este blog ya hemos tratado algunas veces cómo desarrollar interfaces gráficas utilizando wxPython (véase aquí). En este post vamos a presentar el desarrollo de una mini-aplicación que calcula límites, derivadas e integrales utilizando SymPy como motor de cálculo simbólico.

En la siguiente imagen se muestra la interfaz gráfica resultante:

Ahora vamos a describir un poco el diseño y funcionamiento de la aplicación. Como puede notarse la aplicación está basada en la clase wxNotebook, teniéndose cuatro páginas: Límites, Derivadas, Integrales indefinidas y Ayuda. Cada página del Notebook contiene un instancia heredada de wxPanel, misma que contiene todos los controles necesarios en cada uno de los casos.

El proceso de cálculo es básicamente como sigue: se lee una función $f(x)$ introducida en los wxTextCtrl correspondientes y se hace un preprocesamiento mínimo, enseguida se transforma el string leído en una expresión de SymPy, para luego llevar a cabo la operación requerida. El resultado de la operación dado por SymPy se muestra en un canvas de Matplotlib, en el cual se renderiza la expresión LaTeX del resultado obtenido.

El código de la aplicación completa se adjunta a continuación.

Resolviendo problemas de optimización con Python

Ejemplo. Se requiere diseñar un contenedor de desechos, cuya forma será la de un prisma rectangular de dimensiones x, y, z como se muestra en la figura. El volumen del contenedor deberá ser de 0.08 $m^3$. Por requerimientos, la altura debe estar en el rango 0.6-0.8 m, y los lados de la base deben medir al menos 0.1 m. Luego, la idea es diseñar el contenedor de tal forma que se utilice la menor cantidad de material posible. Nota: El contenedor está abierto en la parte superior.

La cantidad de material utilizada para fabricar el contenedor está directamente relacionada con el área del prisma rectangular, donde esta corresponde a la suma de las áreas para cada una de las caras, luego, basándonos en el esquema de la figura anterior, se puede expresar el área en función de sus dimensiones como sigue:

$$ A(x_1,x_2,x_3) = 2x_1x_3 + 2x_2x_3 + x_1x_2 $$

De la restricción de volumen se tiene que:

$$ V = x_1 x_2 x_3 = 0.08\,\, \text{m}^3 \,\,\,\,\,\, \rightarrow \,\,\,\,\,\, x_3 = \frac{V}{x_1 x_2} $$

Entonces, sustituyendo $x_3$ en la primera ecuación podemos reescribirla en términos de $x_1$ y $x_2$ , como sigue:

$$ A(x_1,x_2) = 2\frac{V}{x_2} + 2\frac{V}{x_1} + x_1x_2 $$

Se pide que $0.6 \leq x_3 \leq 0.8$, lo cual se puede expresar como dos ecuaciones de restricción:

$$ g_1(x_1,x_2) = \frac{V}{x_1x_2} - 0.6 \geq 0 $$$$ g_2(x_1,x_2) = \frac{V}{x_1x_2} - 0.8 \leq 0 $$

Además, para los requirimientos dimensionales de la base se tiene:

$$ g_3(x_1,x_2) = x_1 - 0.1 \geq 0 $$$$ g_4(x_1,x_2) = x_2 - 0.1 \geq 0 $$

Óptimo no restringido

Un valor óptimo no restringido puede calcularse mediante herramientas elementales de cálculo variacional. Para obtener un punto crítico se resuelve el sistema de ecuaciones resultante de igualar a cero el gradiente de la función, es decir:

$$ \left[ \vec{\nabla} f \right] = [0] $$

Esto lo podemos calcular utilizando SymPy. Importamos la librería y las funciones hessian y zeros del módulo matrices:

In [1]:
from __future__ import division
import sympy as sym
from sympy.matrices import hessian, zeros
sym.init_printing(use_latex="mathjax")

def gradient(f, varls):
    """ Calcula el gradiente de una función f(x1,x2,...)"""
    n = len(varls)
    G = zeros(n,1)
    for i in range(n):
        G[i] = f.diff(varls[i])
    return G

Definimos las variables simbólicas a utilizar y la función $A(x_1,x_2)$

In [2]:
x1,x2,V = sym.symbols("x_1,x_2,V")
A = 2*(V/x2) + 2*(V/x1) + x1*x2

Calculamos el gradiente de la función y establecemos el sistema de ecuaciones:

In [3]:
G = gradient(A.subs(V,0.08), (x1,x2))
sym.Eq(G,zeros(2,1))
Out[3]:
$$\left[\begin{matrix}x_{2} - \frac{0.16}{x_{1}^{2}}\\x_{1} - \frac{0.16}{x_{2}^{2}}\end{matrix}\right] = \left[\begin{matrix}0\\0\end{matrix}\right]$$

Resolviendo el sistema de ecuaciones, se tiene:

In [4]:
sol = sym.solve(sym.Eq(G,zeros(2,1)), (x1,x2))
sol
Out[4]:
$$\left [ \left ( 0.542883523318981, \quad 0.542883523318981\right ), \quad \left ( -0.271441761659491 - 0.470150922490239 i, \quad -0.271441761659491 - 0.470150922490239 i\right ), \quad \left ( -0.271441761659491 + 0.470150922490239 i, \quad -0.271441761659491 + 0.470150922490239 i\right )\right ]$$

Se obtienen dos soluciones complejas y una real, además, la solución real presenta dos valores iguales.

Hasta ahora sabemos que la solución obtenida es un punto crítico, pero no podemos asegurar si es un máximo, un mínimo o un punto de silla. Luego, para comprobar si el punto calculado es un mínimo se puede utilizar el siguiente teorema:

Teorema 1. Condición necesaria: si $f(x)$ tiene un mínimo local en $x^*$ entonces:

$$ \frac{\partial f(x^*)}{\partial x_i} = 0 ; \,\,\,\, i=1,2,...,n $$

Condición necesaria de segundo orden: si $f(x)$ tiene un mínimo local en $x^*$, entonces la matriz Hessiana

$$ H(x^*) = \left[ \frac{\partial^2 f}{\partial x_i \partial x_j} \right]_{(nxn)} $$

es positiva definida o semidefinida en el punto $x^*$.

Condición suficiente de segundo orden: si la matriz $H(x^*)$ es positiva definida en el punto estacionario $x^*$, entonces $x^*$ es un mínimo local para la función $f(x)$.

Una matriz es positiva definida si sus eigenvalores son estrictamente positivos, es decir: $\lambda_i > 0$

De acuerdo al teorema 1, para que el punto calculado sea un mínimo, los eigenvalores de la matriz Hessiana de $A(x_1, x_2)$ evaluada en el punto crítico deben ser positivos.

Calculando la matriz Hessiana:

In [5]:
H = hessian(A, (x1,x2))
# evaluando en el punto calculado
_x1, _x2 = sol[0][0], sol[0][1]
H = H.subs({x1:_x1, x2:_x2, V:0.08})
H
Out[5]:
$$\left[\begin{matrix}2.0 & 1\\1 & 2.0\end{matrix}\right]$$

Calculando los eigenvalores:

In [6]:
egvals = H.eigenvals()
egvals
Out[6]:
$$\left \{ 1 : 1, \quad 3 : 1\right \}$$

Dado que tanto $\lambda_1 = 1$ como $\lambda_2 = 3$ son positivos, entonces el punto calculado es un mínimo global no restringido.

Ahora, si calculamos el valor de la altura $x_3$:

In [59]:
_x3 = 0.08/(_x1*_x2)
_x3
Out[59]:
$$0.271441761659491$$

Vemos que no cumple con la restricción de altura dada. En lo subsiguiente se abordará el cómo utilizar otros métodos para resolver problemas de optimización con restricciones.

Óptimización con restricciones

La optimización con restricciones implica tomar en cuenta ecuaciones que definen regiones factibles. Uno de los métodos más simples para resolver problemas con restricciones es el método gráfico.

Usando el método gráfico

Resolver un problema de optimización por el método gráfico implica trazar las gráficas de contorno de la función objetivo, así como las regiones definidas por las restricciones, y luego mediante inspección ubicar el punto que optimice la función objetivo. Claro está que es un método que puede resultar poco preciso, pero normalmente sirve para dar una aproximación aceptable y que pueda utilizarse como punto de entrada para algoritmos numéricos como los que se verán en la siguiente sección.

Lo primero que haremos es importar las librerías a utilizar: NumPy y Matplotlib.

In [50]:
import matplotlib.pyplot as plt
import numpy as np
%matplotlib inline

Posteriormente definimos la función objetivo y las restricciones:

In [54]:
V = 0.08
xmin,xmax,ymin,ymax = 0,1,0,1
xx1 = np.linspace(xmin,xmax)
xx2 = np.linspace(ymin,ymax)
x1, x2 = np.meshgrid(xx1, xx2)
f = 2*(V/x2) + 2*(V/x1) + x1*x2
g1 = (V/(x1*x2)) - 0.6
g2 = (V/(x1*x2)) - 0.8
g3 = x1 - 0.1
g4 = x2 - 0.1

Para trazar la gráfica de contorno de $A(x_1,x_2)$ utilizamos la función contour:

In [57]:
plt.contour(x1, x2, f, cmap="viridis")
Out[57]:
<matplotlib.contour.QuadContourSet at 0xbb675f0>

Vemos que la gráfica obtenida no parece muy adecuada, esto debido a que los niveles utilizados por defecto no se ajustan a la función en cuestión. Si evaluamos $A(x_1,x_2)$ en el punto óptimo no restringido de la sección anterior se tiene:

In [63]:
A.subs({"x_1":_x1,"x_2":_x2,"V":0.08})
Out[63]:
$$0.884167559673693$$

Luego, se pueden utilizar niveles personalizados para la función contour, partiendo desde el valor mínimo anterior:

In [76]:
lvl = [0.885, 0.9, 0.95, 1.0, 1.1, 1.25, 1.5, 2.0]
cl = plt.contour(x1, x2, f, levels=lvl, cmap="viridis")
plt.clabel(cl);

Ahora trazamos las gráficas de las restricciones como se muestra enseguida:

In [77]:
plt.contour(x1, x2, g1, levels=[0], linestyles='--', linewidths=2, colors="b")
plt.contour(x1, x2, g2, levels=[0], linestyles='--', linewidths=2, colors='r')
plt.contour(x1, x2, g3, levels=[0], linestyles='--', linewidths=2, colors='k')
plt.contour(x1, x2, g4, levels=[0], linestyles='--', linewidths=2, colors='k');

De la gráfica anterior, la región factible corresponde a la región ubicada entre las líneas rojas y azules correspondientes a las ecuaciones de restricción para la altura $g_1$ y $g_2$. Esto se puede hacer más evidente se utilizamos la función contourf.

In [119]:
plt.contourf(x1, x2, g1, levels=[0,5], colors="#808080", alpha=0.4)
plt.contourf(x1, x2, g2, levels=[-5,0], colors="#808080", alpha=0.4)
plt.contourf(x1, x2, g3, levels=[0,1], colors="#808080", alpha=0.2)
plt.contourf(x1, x2, g4, levels=[0,1], colors="#808080", alpha=0.2);

Donde la región más oscura corresponde a la región factible.

Trazando todo lo anterior en una misma gráfica:

In [122]:
lvl = [0.885, 0.9, 0.95, 1.0, 1.1, 1.25, 1.5, 2.0]
cl = plt.contour(x1, x2, f, levels=lvl, cmap="viridis")
plt.clabel(cl)
plt.contour(x1, x2, g1, levels=[0], linestyles='--', linewidths=2, colors="b")
plt.contour(x1, x2, g2, levels=[0], linestyles='--', linewidths=2, colors='r')
plt.contour(x1, x2, g3, levels=[0], linestyles='--', linewidths=2, colors='k')
plt.contour(x1, x2, g4, levels=[0], linestyles='--', linewidths=2, colors='k');

#ax.plot([0.365],[0.365], 'wo')
#ax.set_xlim(xmin,xmax)
#ax.set_ylim(ymin,ymax)

Por inspección podemos observar que el punto (0.365, 0.365) podría ser un punto óptimo aproximado:

In [137]:
lvl = [0.885, 0.9, 0.95, 1.0, 1.1, 1.25, 1.5, 2.0]
cl = plt.contour(x1, x2, f, levels=lvl, cmap="viridis")
plt.clabel(cl)
plt.contour(x1, x2, g1, levels=[0], linestyles='--', linewidths=2, colors="b")
plt.contour(x1, x2, g2, levels=[0], linestyles='--', linewidths=2, colors='r')
plt.contour(x1, x2, g3, levels=[0], linestyles='--', linewidths=2, colors='k')
plt.contour(x1, x2, g4, levels=[0], linestyles='--', linewidths=2, colors='k');
plt.plot([0.365],[0.365], 'mo', markersize=8)
plt.text(0.365 + 0.01, 0.365 + 0.01, "Q");

Utilizando SciPy

In [153]:
from scipy.optimize import minimize

def A(x):
    v = 0.08
    return 2*(v/x[1]) + 2*(v/x[0]) + x[0]*x[1]
    
def dxA(x):
    v = 0.08
    dfdx0 = -2*(v/x[0]**2) + x[1]
    dfdx1 = -2*(v/x[1]**2) + x[0]
    return np.array([dfdx0, dfdx1])

cons = ({'type': 'ineq',
        'fun' : lambda x: np.array([x[0] - 0.1]),
        'jac' : lambda x: np.array([1,0])},
        {'type': 'ineq',
        'fun' : lambda x: np.array([x[1] - 0.1]),
        'jac' : lambda x: np.array([0,1])},
        {'type': 'ineq',
        'fun' : lambda x: np.array([0.08/(x[0]*x[1]) - 0.6]),
        'jac' : lambda x: np.array([-0.08/(x[0]**2*x[1]), -0.08/(x[0]*x[1]**2)])},
        {'type': 'ineq',
        'fun' : lambda x: -np.array([0.08/(x[0]*x[1]) - 0.8]),
        'jac' : lambda x: -np.array([-0.08/(x[0]**2*x[1]), -0.08/(x[0]*x[1]**2)])},
       )
        
#~ cons = ({'type': 'ineq',
        #~ 'fun' : lambda x: np.array([x[0] - 0.1])},
        #~ {'type': 'ineq',
        #~ 'fun' : lambda x: np.array([x[1] - 0.1])},
        #~ {'type': 'ineq',
        #~ 'fun' : lambda x: np.array([0.08/(x[0]*x[1]) - 0.6])}
        #~ )
    
res = minimize(A, [0.3,0.3], method='SLSQP', jac= dxA, 
               constraints=cons)
res
Out[153]:
     fun: 1.0096894255035318
     jac: array([-0.83485163, -0.83485163,  0.        ])
 message: 'Optimization terminated successfully.'
    nfev: 5
     nit: 5
    njev: 5
  status: 0
 success: True
       x: array([ 0.36514837,  0.36514837])
In [ ]:
 

Análisis estructural con Python I. Una introducción

El análisis estructural es uno de los aspectos elementales para aquellos que nos dedicamos al diseño mecánico o cuestiones similares. En los cursos universitarios de resistencia de materiales se enseñan algunos métodos analíticos que permiten obtener soluciones rápidas para componentes mecánicos simples. No obstante, cuando las geometrías se complican se hace necesario utilizar el enfoque numérico e implementar una metodología de solución utilizando el método de los elementos finitos, el cual proporciona una solución aproximada que será adecuada en medida de ciertos criterios, tales como el tamaño y tipo de elementos, la física del problema, entre otros.

El propósito del presente minicurso es introducir al lector en el uso de las herramientas numéricas que proporciona Python para la solución de problemas de análisis estructural utilizando el método de los elementos finitos.

El método de los elementos finitos consiste en la discretización de un continuo en pequeños elementos, con la finalidad de establecer un sistema de ecuaciones que describa de manera aproximada el comportamiento individual y global del sistema, pasando por la inclusión de las condiciones de frontera y todas las consideraciones físicas que deriven en la simplificación del problema.

En análisis numérico estructural el método de los desplazamientos o de la rigidez, asume que los desplazamientos nodales son las variables desconocidas y comúnmente se debe resolver una ecuación algebraica del tipo:

$$ K U = F $$

Donde $F$ es el vector de fuerzas nodales, $K$ la matriz global de rigidez y $U$ el vector de desplazamientos nodales.

Dado que normalmente los desplazamientos son las incógnitas, se tiene que:

$$ U = K^{-1} F $$

En el caso de un análisis estático lineal esta ecuación se resuelve como se muestra: calculando la inversa de la matriz rigidez y multiplicando por el vector de fuerzas nodales, para el caso de un análisis no lineal se utilizan métodos iterativos.

La matriz global $K$ se obtiene de ensamblar todas las matrices de rigidez por elemento acorde a la numeración o posición de sus nodos.

Un caso muy simple: elementos resorte

El elemento resorte es el elemento finito más simple, tiene un grado de libertad (por cada nodo): traslación en dirección axial.

La matriz de rigidez para un elemento resorte viene dada por:

$$ k^{(e)} = \begin{bmatrix} k_e & - k_e \\ - k_e & k_e \\ \end{bmatrix} $$

La obtención de la matriz de rigidez puede consultarla en la mayoría de los libros introductorios de elementos finitos, por ejemplo en [1]. En lo anterior $k_e$ es la rigidez característica del resorte.

Para ejemplificar cómo funciona el método de los elementos finitos en elementos de este tipo vamos a resolver el siguiente ejemplo.

Primero, las matrices de rigidez por elemento vienen dadas por:

$$ k^{(1)} = \begin{bmatrix} 1000 & -1000 \\ -1000 & 1000 \\ \end{bmatrix} \,\,\,\,\,\,\,\, ; k^{(2)} = \begin{bmatrix} 2000 & -2000 \\ -2000 & 2000 \\ \end{bmatrix} \,\,\,\,\,\,\,\, ; k^{(3)} = \begin{bmatrix} 3000 & -3000 \\ -3000 & 3000 \\ \end{bmatrix} $$

La matriz global se obtiene ensamblando las matrices elementales utilizando el principio de superposición, es decir, asumiendo que los efectos individuales de cada elemento a la matriz global pueden adicionarse de manera independiente. Para esto se debe expandir la matriz de rigidez elemental y rellenar sólo las posiciones correspondientes a los nodos que componen el elemento.

Por ejemplo, el elemento 1 está conformado por los nodos 1 y 3, entonces:

$$ K^{(1)} = \begin{bmatrix} 1000 & 0 & -1000 & 0\\ 0 & 0 & 0 & 0 \\ -1000 & 0 & 1000 & 0 \\ 0 & 0 & 0 & 0 \\ \end{bmatrix} $$

Mismo caso para los otros elementos:

$$ K^{(2)} = \begin{bmatrix} 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ 0 & 0 & 2000 & -2000 \\ 0 & 0 & -2000 & 2000 \\ \end{bmatrix} $$$$ K^{(3)} = \begin{bmatrix} 0 & 0 & 0 & 0 \\ 0 & 3000 & 0 & -3000 \\ 0 & 0 & 0 & 0 \\ 0 & -3000 & 0 & 3000 \\ \end{bmatrix} $$

Luego, la matriz global de rigidez se obtiene sumando las matrices elementales expandidas:

$$ K = \begin{bmatrix} 1000 & 0 & -1000 & 0 \\ 0 & 3000 & 0 & -3000 \\ -1000 & 0 & 3000 & -2000 \\ 0 & -3000 & -2000 & 5000 \\ \end{bmatrix} $$

Quedando el sistema de ecuaciones como:

$$\begin{bmatrix} 1000 & 0 & -1000 & 0 \\ 0 & 3000 & 0 & -3000 \\ -1000 & 0 & 3000 & -2000 \\ 0 & -3000 & -2000 & 5000 \end{bmatrix} \begin{bmatrix} u_1 \\ u_2 \\ u_3 \\ u_4 \end{bmatrix} = \begin{bmatrix} 0 \\ 0 \\ 0 \\ 5000 \end{bmatrix} $$
In [1]:
import numpy as np
import numpy.linalg as la

# Datos iniciales
k1 = 1000.0
k2 = 2000.0
k3 = 3000.0
P = 5000.0

# Matrices por elemento
K1 = np.array([[k1,-k1],[-k1,k1]])
K2 = np.array([[k2,-k2],[-k2,k2]])
K3 = np.array([[k3,-k3],[-k3,k3]])

# Matriz global 
K = np.array([[  K1[0,0],        0,           K1[0,1],                0],
               [      0 ,  K3[0,0],                 0,          K3[0,1]],
               [ K1[1,0],        0,   K1[1,1]+K2[0,0],          K2[0,1]],
               [       0,  K3[1,0],           K2[1,0],  K2[1,1]+K3[1,1]]])

F = np.array([0, 0, 0, P])

# Condiciones de frontera
# Nodos 1 y 2 conocidos -> UX = 0
KS = K[2:,2:]
FS = F[2:]

# Resolviendo
USOL = la.solve(KS, FS)

# Vector de desplazamientos
USOL = np.concatenate(([0,0],USOL))

# Obteniendo las fuerzas nodales
NF = np.dot(K,USOL)

# Presentando los resultados
for nodo in range(4):
    print("%g  UX = %-8.4f    FX = %-8.4f"%(nodo+1, USOL[nodo], NF[nodo]))
1  UX = 0.0000      FX = -909.0909
2  UX = 0.0000      FX = -4090.9091
3  UX = 0.9091      FX = 0.0000  
4  UX = 1.3636      FX = 5000.0000

Utilizando NuSA: una librería para análisis estructural

NuSA es una librería Python para el análisis estructural, que facilita el planteamiento y la solución de este tipo de análisis mediante una colección de clases que reciben como dato de entrada las características elementales de un modelo: coordenadas modales, propiedades del material, condiciones de frontera, etc., y retorna valores de salida básicos como desplazamientos y fuerzas.

Para testear las capacidades de NuSA vamos a resolver el ejemplo del elemento resorte. Lo primero es importar las clases que usaremos: Node, Spring y SpringModel.

In [25]:
from nusa.core import Node
from nusa.element import Spring
from nusa.model import SpringModel

Luego creamos un modelo de tipo Spring utilizando la clase SpringModel:

In [26]:
m1 = SpringModel("Modelo 01")
print(m1)
Model: Modelo 01
Nodes: 0
Elements: 0

Ahora creamos los nodos que conformarán el elemento resorte:

In [27]:
n1 = Node((0,0))
n2 = Node((0,0))
n3 = Node((0,0))
n4 = Node((0,0))

En este caso no es necesario especificar las coordenadas nodales dado que un elemento resorte sólo necesita la rigidez para la formulación, así que se pueden dejar ambos nodos con coordenadas (0,0).

Enseguida se define un elemento de tipo SpringElement, cuyos datos de entrada son una tupla con los nodos que le conforman y la rigidez del elemento.

In [28]:
e1 = Spring((n1,n3),1000)
e2 = Spring((n3,n4),2000)
e3 = Spring((n4,n2),3000)

Una vez se han definido los nodos y elementos se procede a agregar estos al modelo creado.

In [40]:
m1.add_node(n1)
m1.add_node(n2)
m1.add_node(n3)
m1.add_node(n4)

m1.add_element(e1)
m1.add_element(e2)
m1.add_element(e3)

Luego, establecemos las condiciones de frontera y la carga externa aplicada. Finalmente utilizamos el método solve e imprimimos los resultados obtenidos.

In [41]:
m1.add_constraint(n1, ux=0)
m1.add_constraint(n2, ux=0)
m1.add_force(n4, (5000,))

m1.solve()

for nodo in m1.get_nodes():
    print("%s  UX = %-8.4f    FX = %-8.4f"%(nodo.label,nodo.ux,nodo.fx))
0  UX = 0.0000      FX = -909.0909
1  UX = 0.0000      FX = -4090.9091
2  UX = 0.9091      FX = 0.0000  
3  UX = 1.3636      FX = 5000.0000

Referencias

  • [1] Logan, D. L. (1986). A first course in the finite element method. Boston: PWS Engineering.
  • [2] Zienkiewicz, O. C., Taylor, R. L., Zhu, J. Z., Zienkiewicz, O. C., & Zienkiewicz, O. C. (2005). The finite element method: Its basis and fundamentals. Oxford: Elsevier Butterworth-Heinemann.
In [ ]:
 

Calculando integrales con Python y SymPy

Las integrales son unos de los conceptos básicos en la formación matemática de un ingeniero, es en términos básicos la operación inversa de la derivación. Pero, además del concepto puramente matemático, las integrales tienen múltiples interpretaciones geométricas y físicas.

En un curso ordinario de cálculo se nos enseñan métodos para resolver de forma analítica funciones que sean integrables. Por ejemplo todos sabemos que la integral de una función constante será:

$$ \int a\,dx = ax + C $$

Y lo sabemos porque nos hemos aprendido reglas básicas de integración y por supuesto a indentificar el tipo de función. Actualmente existen paquetes de álgebra simbólica que son capaces de realizar esta tarea: identificar el caso que se tiene y aplicar métodos computacionales, hasta cierto grado complejos, para determinar la solución.

Y claro, SymPy es uno de esos sistemas de álgebra computacional (CAS), en el que solo necesitamos escribir nuestra función a integrar, utilizar por ahí alguna rutina y obtener un resultado rápidamente. Pero claro, para ello debemos aprender mínimamente la sintaxis y eso es justo lo que veremos enseguida.

Integrales simples

Vamos a ver cómo resolver integrales simples indefinidas, si, de esas que vemos en un primer curso. Para resolverlas tendremos que utilizar la función integrate. Por ejemplo se tiene la siguiente función $f(x)=x^2-3x+2 $.

Como primer paso debemos importar lo que necesitaremos del paquete SymPy:

In [1]:
from sympy import integrate, init_printing
from sympy.abc import x
init_printing(use_latex="mathjax")

Del módulo abc importamos la variable simbólica x e integrate para resolver nuestra integral. Ahora, podemos guardar la función a integrar en una variable o bien pasarla directamente como argumento:

In [2]:
f = x**2 - 3*x + 2
integrate(f)
Out[2]:
$$\frac{x^{3}}{3} - \frac{3 x^{2}}{2} + 2 x$$

En este caso no hemos tenido incovenientes, porque en la expresión a integrar sólo existe una variable simbólica, pero si la expresión tuviese más de una, habría que especificar de manera explícita la variable respecto a la cual se integra, de lo contrario Python nos mostrará un error:

In [3]:
from sympy.abc import a,b,c
f = a*x**2+b*x+c
integrate(f)
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-3-476839d3c49d> in <module>()
      1 from sympy.abc import a,b,c
      2 f = a*x**2+b*x+c
----> 3 integrate(f)

~\Anaconda3\lib\site-packages\sympy\integrals\integrals.py in integrate(*args, **kwargs)
   1289     risch = kwargs.pop('risch', None)
   1290     manual = kwargs.pop('manual', None)
-> 1291     integral = Integral(*args, **kwargs)
   1292 
   1293     if isinstance(integral, Integral):

~\Anaconda3\lib\site-packages\sympy\integrals\integrals.py in __new__(cls, function, *symbols, **assumptions)
     73             return function._eval_Integral(*symbols, **assumptions)
     74 
---> 75         obj = AddWithLimits.__new__(cls, function, *symbols, **assumptions)
     76         return obj
     77 

~\Anaconda3\lib\site-packages\sympy\concrete\expr_with_limits.py in __new__(cls, function, *symbols, **assumptions)
    375                     " more than one free symbol, an integration variable should"
    376                     " be supplied explicitly e.g., integrate(f(x, y), x)"
--> 377                     % function)
    378             limits, orientation = [Tuple(s) for s in free], 1
    379 

ValueError:  specify dummy variables for a*x**2 + b*x + c. If the integrand contains more than one free symbol, an integration variable should be supplied explicitly e.g., integrate(f(x, y), x)

Pues eso, si intentamos integrar la función $f(x)=ax^2+bx+c$ sin especificar la variable de integración, Python nos mandará un error que es bastante sugerente al respecto. Así, lo correcto sería:

In [4]:
integrate(f, x)
Out[4]:
$$\frac{a x^{3}}{3} + \frac{b x^{2}}{2} + c x$$

Integrales definidas

Una integral definida usualmente se utiliza para calcular el área bajo la curva de una función en un intervalo finito. En SymPy, para calcular una integral definida se utiliza la función integrate, considerando el hecho que deben adicionarse los límites de evaluación mediante la sintaxis:

integrate(fun, (var, a, b))

Donde fun es la función, var la variable respecto a la cual se integra, a el límite inferior y b el límite superior.

Para ejemplificar vamos a resolver la siguiente integral definida:

$$ \int_0^{\frac{\pi}{2}} \cos x \,\, dx $$
In [5]:
from sympy import cos,pi

integrate(cos(x), (x,0,pi/2.0))
Out[5]:
$$1$$

Otro ejemplo:

$$ \int_0^5 x \, dx $$
In [6]:
integrate(x, (x,0,5))
Out[6]:
$$\frac{25}{2}$$

Integrales múltiples

Ahora vamos a resolver integrales dobles (la sintaxis/metodología de resolución que se revisará aplica para cualquier integral múltiple). Por ejemplo vamos a resolver la integral dada por:

$$ \int_a^b \int_c^d \, dy \, dx $$

Recuerde que este tipo de integrales múltiples se resuelven de forma iterada, yendo de dentro hacia afuera, es decir, para la integral anterior se procedería:

$$ I_1 = \int_c^d \, dy \qquad \rightarrow \qquad I = \int_a^b \, I_1 \, dx $$

En Python/SymPy hacemos exactamente lo mismo:

In [9]:
from sympy.abc import x,y,z,a,b,c,d
from sympy import simplify
In [10]:
I1 = integrate(1, (y,c,d))
simplify( integrate(I1, (x,a,b) ) )
Out[10]:
$$\left(a - b\right) \left(c - d\right)$$

Almacenando imágenes como código utilizando img2py

En esta breve entrada vamos a hablar de cómo utilizar la función ìmg2py de wxPython, para generar códigos a partir de imágenes, y posteriormente poder utilizarlas en conjunto conPyEmbeddedImage.



Lo primero que debemos hacer es importar la función del módulo wx.tools.img2py:

from wx.tools.img2py import img2py

Y posteriormente sólo llamaremos a la función pasando como argumentos el archivo de la imagen y el archivo Python de salida, por ejemplo:

img2py("icono.png", "icono.py")

Con lo anterior se obtendrá un archivo Python icono.py que contiene algo parecido a lo siguiente:

#----------------------------------------------------------------------
# This file was generated by C:\Users\User\Desktop\LABPro\_blogs_\Posts\test_app.py
#
from wx.lib.embeddedimage import PyEmbeddedImage

icono = PyEmbeddedImage(
"iVBORw0KGgoAAAANSUhEUgAAAA0AAAANCAYAAABy6+R8AAAABHNCSVQICAgIfAhkiAAAAAlw"
"SFlzAAAN1wAADdcBQiibeAAAABl0RVh0U29mdHdhcmUAd3d3Lmlua3NjYXBlLm9yZ5vuPBoA"
"AABlSURBVCiRY2TAChryGBgYJmKXY2BgwiWBD7D819JqQxfMfHHXZMY7ZdyaGBgYKtEFNdk+"
"4bWJLOeR5ydsgtrsnxiSBe5hiG/9IrX2xR+OD4z/tbT+E20FI6M249Wr1+jnp0GuCQD2OxXd"
"g2khLwAAAABJRU5ErkJggg==")

Luego, ese código generado lo podremos utilizar posteriormente en nuestras aplicaciones, por ejemplo para un icono:

# -*- coding:utf-8 -*-
import wx
from icono import icono

if __name__=='__main__':
app = wx.App()
fr = wx.Frame(None, -1, u"wxPython")
fr.SetIcon(icono.GetIcon())
fr.Show()
app.MainLoop()

¿Y si necesito embeber varias imágenes en un mismo archivo Python?

Si, es bastante probable que necesite esto, hacer una colección de imagenes embebidas en un mismo archivo Python, la cuestión es muy parecida, sólo tenemos que agregar el keyword argumentappend a la función img2py para evitar que borre lo que se ha colocado anteriormente, un pequeño ejemplo:

import glob
from wx.tools.img2py import img2py

if __name__=='__main__':
for img in glob.glob("img/*.png"):
img2py(img,"iconos.py", append=True)

Lo anterior incrusta el código correspondiente a todas las imágenes PNG contenidas en la carpeta img/.