octubre 2012

You are currently browsing the monthly archive for octubre 2012.

Toda aplicación es habitual que tenga un fichero de configuración. Lo primero que hace el programa cuando empieza a funcionar es leer este fichero en el que se especifican una serie de parámetros que determinan una cierta configuración de funcionamiento.

Cuando somos nosotros los que estamos desarrollando una aplicación, la pregunta que se nos pasa por la cabeza es: ¿Qué formato debe tener un fichero de configuración? Ciertamente, dado que nosotros somos los desarrolladores, tenemos total libertad para determinar cual sera éste. Sin embargo, el hecho de seleccionar uno u otro nos permitirá aproximar nuestra aplicación a la manera estándar de funcionar el software, teniendo también la posibilidad de utilizar herramientas desarrolladas para estos estándares de facto.

Volviendo a nuestra libertad de elección, la primera manera que se nos ocurre es tener un fichero de texto con un formato determinado. Lo único que tenemos que hacer es desarrollar un código capaz de leer y escribir ficheros con este formato.

De unos años a esta parte, los formatos mas utilizados para los ficheros de configuración son los ficheros INI, XML, YAML y JSON.

INI: es el mas sencillo. Solo admite dos niveles y no esta preparado para binarios.

;config valSPH
...
numParticles=10000
kernel=cubic
odeSolver=rungeKutta23
...

XML: muy extendido. Soporta múltiples niveles, binarios… Hay muchas herramientas para parsearlo, por ejemplo Xerces de Apache, o también es fácil programarse uno basado en SAX.


<!-- config valSPH -->
...
<numParticles>10000</numParticles>
<kernel>cubic</cubic>
<odeSolver>rungeKutta23<odeSolver>
...

JSON: Dispone de MIME type, “application/json” con codificación y decodificación nativa en navegadores. Muy utilizado en projectos Ajax.

{

...
"numParticles": 10000,
"kernel": "cubic",
"odeSolver": "rungeKutta23",
...

}

JAML: nace como una generalización de JSON. Existe soporte en C++ con yaml-cpp

---# config valSPH
...
numParticles:  10000
kernel:        cubic
odeSolver:     rungeKutta23
...

Tags: , , , , , , ,

Actualizamos la minianimación que teniamos incluyendo los algoritmos y las estructuras de datos necesarias para la NNPS.

Con el fin de poder comprobar la determinación de vecinos, el color de cada péndulo depende, en un instante dado, del número de vecinos que tiene a una distancia menor o igual a $latex h$.

En la siguiente animación se utiliza una LinkedList:

Tags: , , ,

En el método SPH necesitamos conocer cual es su conjunto de vecinos de una partícula $latex a$. Este viene determinado por la smoothing length $latex h$, uno de los parámetros de las funciones kernel.

Supongamos que tenemos $latex N$ partículas. Si $latex N$ es lo suficientemente pequeño, lo único que tenemos que hacer es, para cada partícula, recorrer el resto de partículas calculando su distancia a nuestra partícula y tener en consideración solo aquellas que esten a distancia menor que $latex h$ (all-pair search). La complejidad del algorítmo es $latex O(N^2)$.

Sin embargo, cuando el número de partículas crece, y que lo ha de hacer pues en eso se basa la potencia del método, esta manera de trabajar es mejorable.

A continuación veremos algunas opciones, pero la idea general es que necesitamos tener una descomposición espacial de nuestro escenario que nos permita encontrar de manera eficiente el vecindario de cada una de nuestras partículas y que tendremos que mantener actualizada a medida que evoluciona el sistema.

Otra opción es la linked-list. Se trata de superponer una malla con un tamaño de celda $latex 2h$, de manera que, dada una partícula, solo tendremos que buscar las partículas que se encuentran en las celdas adyacentes. En C++ tenemos las clases list y vector. La primera por nombre y la segunda porque la clase vector puede crecer dinámicamente. Solo podremos utilizarlas si la $latex h$ es constante para todas las partículas. La complejidad, con un número suficientemente pequeño de partículas por celda, es $latex O(N)$.

Finalmente, podemos trabajar con árboles. Los árboles permiten complejidades genéricas en búsquedas de $latex O(N log N)$ y permiten trabajar con $latex h$ diferentes para cada partícula. En astrofísica, por ejemplo, podríamos tener un arbol binario de búsqueda, pues tenemos que simular objetos autogravitantes y necesitamos resolver un problema de n-cuerpos, que parece ser que es lo que mas tiempo de CPU consume. Para hacer esto de manera eficiente necesitamos un BST. Se trata de utilizar este mismo para realizar las búsquedas. La clase map de C++ se implementa con un arbol binario balanceado (AVL).

Existen otras estructuras de datos que permiten hacer descomposiones espaciales que podria ser interesante analizar (UB-tree, Octrees, BST, kd-tree, R-tree).

Tags: , , , , , , ,

En la siguiente simulación tenemos $latex 1000$ partículas. Cada partícula $latex a$ representa un péndulo de masa $latex m_a$ y posición $latex (x_a(t),y_a(t),z_a(t))$, con $latex z_a(t)$ constante, sometido a las fuerzas:

  1. $latex F_1 = -m_ag$, la fuerza de la gravedad, con $latex g$ constante.
  2. $latex F_2$, la fuerza elática de la suspensión, de sentido opuesta a la posición del péndulo y de magnitud $latex k_1d$, con $latex k_1 > 0$ y $latex d$ la variación de la longitud del péndulo respecto de su longitud en reposo $latex l$.
  3. $latex F_3$, la fuerza de fricción, con sentido opuesta a la velocidad del péndulo y magnitud $latex k_2v$, con $latex v = sqrt{(x_a’)^2+(y_a’)^2}$ y $latex k_2 > 0$.

La ley de Newton $latex F = ma$ queda, para cada uno de los péndulos, es decir, para cada partícula, como el sistema de dos ecuaciones diferenciales de segundo orden:

$latex m_a begin{bmatrix} x_a”(t) \ y_a”(t) end{bmatrix} = F_1 + F_2 + F_3$

que es:

$latex m_ax_a” = -k_1 (1-frac{l}{sqrt{x_a^2+y_a^2}})x_a – k_2frac{x’_a}{sqrt{x_a’^2+y_a’^2}}$

$latex m_ay_a” = -k_1 (1-frac{l}{sqrt{x_a^2+y_a^2}})y_a – k_2frac{y’_a}{sqrt{x_a’^2+y_a’^2}}$

Podemos expresar este sistema como un sistemad de ecuaciones de primer orden introduciendo incongnitas adicionales: $latex u_a = x_a’$ y $latex v_a = y_a’$, con lo que nos queda:

$latex left {
begin{array}{rcl}
x_a’ &=& u_a \
y_a’ &=& v_a \
m_au_a’ & = & k_1 (frac{l}{sqrt{x_a^2+y_a^2}}-1)x_a – k_2frac{u_a}{sqrt{u_a^2+v_a^2}} \
m_av_a’ & = & k_1 (frac{l}{sqrt{x_a^2+y_a^2}}-1)y_a – k_2frac{v_a}{sqrt{u_a^2+v_a^2}} – m_ag \
end{array}
right .
$

En la siguiente simulación tomamos $latex m_a$ aleatoria para cada partícula $latex a$, por lo que tendremos que resolver el sistema de EDOs numéricamente para cada una de ellas. Además, tendremos $latex l=1$, $latex k_1 = 10000$ y $latex k_2=1$ en todos los péndulos. Cada uno de éstos se distribuye a lo largo de una espilar cilíndrica y las velocidades iniciales son, para todas ellas, $latex u_{a_0}=0$, $latex v_{a_0} = -0.7$ y $latex v_{z_0} = 0$. El intervalo de tiempo es desde $latex t_i = 0$ hasta $latex t_f=10$ que se divide en $latex 400$ subintervalos.

El resultado es el siguiente:

[youtube http://www.youtube.com/watch?v=emirYBIds6o?rel=0&w=420&h=315]

Tags: , , ,

Para poder generar los call graph y los caller graph desde doxygen necesitamos tener instalado Graphviz.

Aquí tenemos la direccion desde donde nos lo podemos bajar. En el caso de mac, es un pkg que instala el fichero binario dot en el directorio /src/local/bin. Tendremos que añadir este path en la variable DOT_PATH en la pestaña Expert de doxygen.

A continuación, unos ejemplos del tipo de diagramas generados:

Además, en la segunda captura podemos observar que las clases pueden mostrarse con estilo UML.

Tags: , , , ,

A medida que codificamos y probamos nuestras clases las iremos documentando mediante doxygen.

La instalación y utilización básica de doxygen es extremadamente sencilla. El programa está disponible aquí.

A continuación tenemos un ejemplo del formato de los comentarios añadidos y el tipo de salida generado.

Tags: ,

A continuación nuestra primera animación creada con VisIt a partir de ficheros silo generados a partir de llamadas a métodos de nuestras primeras clases.

Temporalmente, aunque es mejor utilizar nuestras partículas simplemente como nodos de interpolación a la hora de visualizar, las mostraremos como puntos.

A partir de este momento, podremos tener una referencia visual de nuestras partículas que nos ayudará mucho a la hora de evaluar como de bien estamos haciendo las cosas.

Sin mas dilación, nuestra minianimación:

[youtube http://www.youtube.com/watch?v=-bgdirYegrc?rel=0&w=420&h=315]

No es mas que un conjunto de mil partículas desplazandose sobre un cilindro. La animación consta de 100 ficheros .silo y hemos generado el mpeg con el wizard de VisIt.

Tags: , ,

Ya dedicamos un post a VisIt. En el creamos una animación a partir de un conjunto de ficheros existentes.

¿Como generamos ficheros .silo desde nuestro código? En primer lugar, necesitamos tener las librerias correspondientes, que las podemos conseguir ejecutando:

./build_visit --console --no-visit --no-thirdparty  --thirdparty-path /usr/local --silo --hdf5 --szip

desde el terminal, y donde build_visit es un fichero que podemos conseguir aqui.

A continuación, necesitamos incluir la libreria <silo.h> en nuestro fichero y añadir:

PLATFORM=i386-apple-darwin10_gcc-4.2

SZIP_DIR=/usr/local/szip/2.1/$(PLATFORM)
SZIP_CPPFLAGS=-I$(SZIP_DIR)/include
SZIP_LDFLAGS=-L$(SZIP_DIR)/lib
SZIP_LIBS=-lsz

HDF5_DIR=/usr/local/hdf5/1.8.4/$(PLATFORM)
HDF5_CPPFLAGS=-I$(HDF5_DIR)/include $(SZIP_CPPFLAGS)
HDF5_LDFLAGS=-L$(HDF5_DIR)/lib $(SZIP_LDFLAGS)
HDF5_LIBS=-lhdf5 $(SZIP_LIBS) -lz

SILO_DIR=/usr/local/silo/4.6.2/$(PLATFORM)
SILO_CPPFLAGS=-I$(SILO_DIR)/include $(HDF5_CPPFLAGS)
SILO_LDFLAGS=-L$(SILO_DIR)/lib $(HDF5_LDFLAGS)
SILO_LIBS=-lsiloh5 $(HDF5_LIBS) -lm

LDFLAGS=$(LDFLAGS) $(SILO_LDFLAGS)
LIBS=$(SILO_LIBS)
CPPFLAGS=$(CPPFLAGS) $(SILO_CPPFLAGS)

a nuestro Makefile

Tags: ,

Cuando eramos niños, nos pasabamos el curso esperando la llegada de las vacaciones. Sin embargo, a medida que pasaban, cada vez tenías mas ganas de volver a las clases.

En esto de programar pasa un poco lo mismo. Cuando estas metido en un proyecto en el que te pasas dias y dias escribiendo, probando código y peleandote con los ordenadores, lo que mas desea uno en el mundo es acabar con aquello. Sin embargo, una vez terminado y pasado un tiempo, después de disfrutar de un merecido descanso y del placer y de la potencia de ver aquello funcionando, nos entra el gusanillo de volver a la carga.

Existen muchos modelos para el ciclo de vida del software, es decir, diferentes maneras de enfrentarte al reto de como desarrollar sofware. Sin embargo, con el tiempo y despues de acumular muchas “horas de vuelo”, uno ya tiene mas o menos definida su pequeña estrategia. Son cosas de sentido común, pero que no va mal tenerlas presente.

Para empezar, es muy importante tener una visión global del sistema, de manera que, nuestro primer diseño, es interesante que lo abarque todo. Para ello van bien, por ejemplo, realizar el diagrama de casos de uso, el diagrama de clases, y algunos diagramas de comunicación y de secuencia. Además, existe software que sera capaz de generarnos, a partir de ellos, la arquitectura de la aplicación, es decir, el conjunto de todos los ficheros de clases que necesitaremos y la interrelación entre ellos.

Sin embargo, en el momento de empezar a implementar el software, a escribir el programa, es mejor seguir una estrategia en espiral.

cuya idea de fondo es sencilla: no pretender programarlo todo de golpe sino hacerlo de manera incremental. Para ello, empezamos seleccionando un pequeño conjunto del global de los objetivos de nuestro software, seleccionaremos las clases que nos permitiran asumirlos, las implementaremos y las probaremos. Llegados a este punto habremos dado una vuelta en la espiral ya que estaremos en condiciones de volver a empezar la misma rutina con otro pequeño conjunto de objetivos pero con la diferencia de que nuestro programa habrá crecido con respecto a la vuelta anterior y lo que allí nos propusimos estará funcionando.

Pues nada, teniendo ésto presente y sin mas preámbulos, nos ponemos el mono de programador y, por fin, ¡a programar!

Tags: ,

En el artículo A modified SPH approach for fluids with large density differences de F. Ott y E. Schnetter se explica una nueva modificación del método SPH que permite la interacción entre fluidos con densidades muy diferentes.

En la aproximación del SPH estandar tenemos:

$latex rho_i = sum_j m_j W_{ij}$

con $latex W_{ij}:=W(boldsymbol{x}_i – boldsymbol{x}_j)$, o también:

$latex frac{d}{dt}rho_i = sum_j m_j boldsymbol{v}_{ij} cdot nabla W_{ij}$

donde $latex boldsymbol{v}_{ij}:=(boldsymbol{v}_i – boldsymbol{v}_j)$ y $latex nabla W_{ij}:=(nabla W)(boldsymbol{x}_i – boldsymbol{x}_j)$, que tiene la ventaja de que los valores iniciales para $latex rho_i$ pueden escogerse libremente permitiendo la modelización de discontinuidades.

Ellos proponen utilizar,

$latex frac{d}{dt}rho_i = m_i sum_j boldsymbol{v}_{ij} cdot nabla W_{ij}$

deducida suavizando la densidad de partículas $latex n_i = frac{1}{V_i}$ en lugar de la densidad de masa $latex rho_i$, donde $latex V_i = frac{m_i}{rho_i}$ y $latex frac{d}{dt}m_i = 0$.

justificandolo con tests positivos en tres problemas: ecuación de advección, onda sonora y tubo de choque.

Tags: , ,

En el review que Rosswog hace sobre el método SPH, en especial en sus aplicaciones a la astrofísica, hay un apartado dedicado al tratamiento de los shocks. El tratamiento de los strong shocks es uno de los puntos débiles del método SPH (penetración de partículas parcialmente resuelto mediante XSPH). En este apartado comenta que, básicamente, exiten dos maneras de tratarlos:

  1. Utilizar soluciones (analíticas) de un problema de Riemann entre dos partículas adyacentes.
  2. Ampliar la discontinuidad hasta una longitud resoluble numéricamente, para lo que se añaden terminos disipativos artificiales extra al flujo.

Comenta que la mayoria de las implementaciones de SPH utilizan esta segunda aproximación, consiguiendo esta amplicación o mediante el método de discretización utilizado, utilizando métodos que modifican las ecuaciones originales de una manera que nos puede venir bien, tipo el esquema Lax; o mediante la adición de terminos ad-hoc, algo bastante habitual, de tipo presión.

Pero también existen métodos que implementan la primera, y cita dos artículos:

  1. El artículo Reformulation of Smoothed Particle Hydrodynamics with Riemann Solver de Inutsuka y Shu-Ichiro.
  2. El artículo Implementations and tests of Godunov-type particle hydrodynamics de Cha y Whitworth.

En el primero comenta como consiguen, tras algunos intentos fallidos,  introducir Riemann Solver en las ecuaciones de manera similar a como se consigue en los métodos basados en malla, manteniendo la formulación conservativa de las mismas. Para empezar, comentan que mientras el error en la mayoria de los métodos basados en malla proviene de la truncación de la serie de de Taylor, en SPH, si consideramos:

$latex langle f rangle(boldsymbol{x}) := int f(boldsymbol{x’})W(boldsymbol{x’}-boldsymbol{x},h)dboldsymbol{x’}$

como la convolución de la función $latex f(x)$ con el kernel. Mediante Taylor, de $latex f$ alrededor de $latex x=x_i$, podemos ver que:

$latex langle frangle(boldsymbol{x}_i) – f(boldsymbol{x}_i) = frac{h^2_{ef}}{4} nabla^2 f(boldsymbol{x}) + O(h^4_{ef})$.

donde

$latex h^2_{ef}:=2 int x^2 W(boldsymbol{x},h)dboldsymbol{x}$,

por lo que el error proviene de la convolución con la función simétrica $latex W$.

También se tiene:

$latex langle frac{partial f}{partial x} rangle (boldsymbol{x}) = int frac{partial f(boldsymbol{x’})}{partial x’} W(boldsymbol{x’}-boldsymbol{x},h)dboldsymbol{x’}$

$latex = – int f(boldsymbol{x’}) frac{partial}{partial x’} W(boldsymbol{x’}-boldsymbol{x},h)dboldsymbol{x’}$

$latex = frac{partial}{partial x} int f(x’)W(x’-x,h)dx’$

De manera que $latex nabla f = nabla langle f rangle$. Si definimos:

$latex rho(boldsymbol{x}):=sum_j m_j W(boldsymbol{x}-boldsymbol{x}_j,h)$,

entonces:

$latex sum_j frac{m_j}{rho(boldsymbol{x})}W(boldsymbol{x}-boldsymbol{x}_j,h) = 1$,

$latex sum_j m_j nabla big [ frac{W(boldsymbol{x}-boldsymbol{x}_j,h)}{rho(boldsymbol{x})} big ]= 0$

que son la clave del método, pues ahora se puede hacer:

$latex f_i:= langle f rangle (boldsymbol{x}_i) = int f(boldsymbol{x’})W(boldsymbol{x’}-boldsymbol{x}_i,h)dboldsymbol{x’}$

$latex = int sum_j m_j frac{f(boldsymbol{x}’)}{rho(boldsymbol{x}’)}W(boldsymbol{x}’-boldsymbol{x}_i,h)W(boldsymbol{x}’-boldsymbol{x}_j,h)dboldsymbol{x}’$

$latex sum_j Delta f_{i,j}$

donde:

$latex Delta f_{i,j}:=int m_j frac{f(boldsymbol{x}’)}{rho(boldsymbol{x}’)}W(boldsymbol{x}’-boldsymbol{x}_i,h)W(boldsymbol{x}’-boldsymbol{x}_j,h)dboldsymbol{x}’$.

En el método SPH estandar nos encontramos, por una parte, con:

$latex sum_j frac{m_j}{rho(boldsymbol{x}_j)}W(boldsymbol{x}-boldsymbol{x}_j,h) approx 1$

y por otra, con:

$latex Delta f_{i,j} approx m_j frac{f(boldsymbol{x}_j)}{rho(boldsymbol{x}_j)}W(boldsymbol{x}_i-boldsymbol{x}_j,h)$

y que, evidentemente, son aproximaciones mas pobres a las anteriormente obtenidas.

A continuación derivan las ecuaciones para la posición:

$latex ddot{boldsymbol{x}}:=int frac{dboldsymbol{v}(boldsymbol{x})}{dt}W(boldsymbol{x}-boldsymbol{x}_i,h)dboldsymbol{x}$

$latex = – int frac{1}{rho(boldsymbol{x})} nabla P(x) W(boldsymbol{x}-boldsymbol{x}_i,h)dboldsymbol{x}$

y para la energía

$latex dot{u}_i := int frac{du(boldsymbol{x})}{dt} W(boldsymbol{x}-boldsymbol{x}_i,h)dboldsymbol{x}$

$latex = sum_j m_j int frac{P}{rho^2}[boldsymbol{v}-dot{boldsymbol{x}}_i] cdot [nabla W(boldsymbol{x}-boldsymbol{x}_i,h)W(boldsymbol{x}-boldsymbol{x}_j,h) – W(boldsymbol{x}-boldsymbol{x}_i,h) nabla W(boldsymbol{x}-boldsymbol{x}_j,h) ] dboldsymbol{x}$,

que se modifica a

$latex$

para tener en cuenta la disipación

En el segundo artículo define e implementa el método Godunov-type Particle Hydrodynamics (GPH). Básicamente GPH = SPH + Godunov upwind schema.

Tags: , , , ,

Hablemos de generalizaciones, cosa esencial en matemáticas:

  • Desde el punto de vista del análisis funcional, ¿qué pasa si nuestro espacio tiene asociada una métrica no definida positiva?
  • En el caso finito, ¿qué pasa cuando en lugar de trabajar con $latex mathbb{R}^n$ o $latex mathbb{C}^n$ tengo variedades diferenciables mas generales?
  • En analisis funcional tengo espacios de funciones donde estas cumplen una propiedad asociada con una medida y especificada normalmente mediante una integral. ¿Puedo tener espacios de variedades cumpliendo propiedades de este tipo? (pues tiene sentido hablar de integración en variedades orientadas)
  • Al hablar de variedades de Riemann, asociamos una métrica a una variedad para poder hablar de distancias, areas, angulos, etc. En realidad, asociamos un producto escalar, un tensor 2 veces covariante, del que deriva una norma a partir de la que podemos especificar una distancia.

Tags: , ,

Teorema de Laurent.

Definición (serie de Laurent): Sea $latex z_0 in mathbb{C}$. Una serie de Laurent en $latex z_0$ es formalmente una expresión de la forma:

$latex sum_{n in mathbb{Z}} a_n (z-z_0)^n$ con $latex a_n in mathbb{C}$.

Llamaremos parte analítica a $latex sum_{n=0}^infty a_n (z-z_0)^n$ y parte principal a $latex sum_{n=1}^infty a_n (z-z_0)^{-n}$. Diremos que la serie es convegente en $latex z$ si lo son su parte analítica y su parte principal.

Definición (anillo): Llamaremos anillo centrado en $latex z_0$ de radio menor $latex r$ y radio mayor $latex R$ a:

$latex A(z_0;r,R):= { z in mathbb{C}: r < |z-z_0| < R}$.

En el caso que $latex r=0$ entonces $latex A(z_0;0,R) = D(z_0,R)-{z_0} = D'(z_o,R)$ tenemos un disco perforado. Si $latex r=0 y R=infty$ entonces $latex A(z_o;0,infty) = mathbb{C}-{z_0}$. Finalmente, si $latex r neq 0$ y $latex R = infty$ entonces $latex A(z_0,r,infty = mathbb{C}-overline{D(z_o,r)}$.

Teorema (de Laurent): Si $latex f in mathcal{H}(A(z_0;r,R))$ entonces existe $latex a_n in mathbb{C}$ tal que:

$latex f(z)= sum_{n in mathbb{Z}} a_n(z-z_0)^n$ para todo $z in A(z_0;r,R)$.

demostración:

Clasificación de singularidades aisladas

Teorema de los residuos.

Definición (función meromorfa): Sea $latex f:U-A longrightarrow mathbb{C}$ una función analítica con $latex U subset mathbb{C}$ un abierto y $latex A$ el conjunto de singularidades aisladas de $latex f$. Diremos que $latex f$ es una función meromorfa en $latex U$.

Teorema (de los residuos)

Integrales tipo I, II, III, IV, V i VI

Tags: , , , , , , ,

Postulado 1: La representación del estado de un sistema físico.

P1. La máxima información posible sobre un sistema físico en un instante dado $latex t$ es su estado cuántico, que se representa por un vector (ket) de norma $latex 1$ y de fase arbitraria en un espacio de Hilbert separable.

Superposición en un producto tensorial de espacios de Hilbert

Postulado 2: La representación de las magnitudes medibles.

P2. Cualquier magnitud del sistema que en pricipio se pueda medir tiene asociada un operador lineal autoadjunto definido sobre el espacio vectorial de los estados. La totalidad de los autovalores recibe el nombre de espectro y los autovectores definen una base del espacio de Hilbert.

Conjunto completo de observables compatibles

Postulado 3: El resultado de la medida.

P3. El resultado de una medida de un observable $latex A$ sobre un estado $latex |psi rangle$ solo puede ser uno de sus autovalores. La probabilidad de que al medir $latex A$, de descomposicón espectral $latex A=sum_i lambda_i Pi_i,, (lambda_i neq lambda_j)$, sobre un estado $latex |psi rangle$ obtengamos el resultado $latex lambda_i$ es:

$latex mathcal{P}_{| psi rangle} (A:lambda_i) = ||Pi_i | psi rangle ||^2 = langle psi | Pi_i | psi rangle$

donde $latex A:lambda_i$ indica que la medida del observable $latex A$ nos da el valor $latex lambda_i$

Dualidad onda-particula

Incertidumbre de Heisenberg

Postulado 4: La proyección, reducción o colapso.

P4. Cualquier estado $latex | phi rangle$ sobre el que se realiza una medida del observable $latex A$ que es filtrante respecto del autovalor $latex lambda_i$, pasa a estar inmediatamente despues de la medida en el estado

$latex frac{Pi_i | phi rangle}{||Pi_i | phi rangle||}$

con probabilidad $latex ||Pi_i | phi rangle||^2$, o se destruye durante el proceso de medida.

La no localidad de la mecánica cuántica

Postulado 5: La evolución temporal.

P5. Entre medidas, el sistema evoluciona según

$latex ihbar frac{partial}{partial t} | varphi(t) rangle = H | varphi (t) rangle$

donde $latex H$ es el operador hamiltoniano del sistema.

Postulado 6: los operadores de posición y de momento.

P6. Los operadores de posición, en coordenadas cartesianas, y de momento satisfacen las reglas de conmutación:

$latex [X_i,X_j]=0, , [P_i,P_j]=0, , [X_i,P_j]=i hbar delta_{ij} I$.

La función de onda en el espacio de posiciones

La función de onda en el espacio de momentos

Dimensiones y teorema de Ehrenfest

Los valores de los operadores $latex X_i$ y $latex P_i$ evolucionan según las leyes clásicas pero con fuerzas medias.

Tags: , , , , , , , , , , ,

Ayer asistí al primer Friday’s miniWorkshop del IVICFA (L’institut Valencià d’Investigació Coorporativa de Física Avançada) que trataba sobre supercomputación y computación GRID. Forma parte de un ciclo de seminarios sobre “Fronteras de la Física” y contaba ayer, entre otros atractivos, con la apertura de los mismo por parte de Ian Bird, Project Leader of the Worldwide LHC Computing GRID (WLCG), con “Petabyte scale computing for the LHC”, o el seminario “Modeling complex solar magnetodynamical phenomena using supercomputing and visualization techniques” de Fernando Moreno Insertis, del Instituto de Astrofísica de Canarias (IAC) y PI of the European Solaire Network.

Aunque ya hicimos algunos comentarios en un post anterior, aprovechamos la ocasión para volver a hablar sobre estos temas y así aclararnos las ideas.

Un computador computa, por lo que un supercomputador supercomputa ;-). Como el elemento fundamental para computar en la arquitectura Von Neumann son las CPUs, en un supercomputador dispondremos de multiples de éstas. El otro elemento fundamental en la arquitectura Von Neumann, que es la que las convierte en máquinas de proposito general, es la memoria. En función de como las CPUs comparten esta memoria, nos encontramos con dos modelos de supercomputación: memoria compartida y memoria distribuida. Las dos interfaces que nos permiten programar en estos ambientes son OpenMP y MPI respectivamente.

Tenemos que aclarar que, cuando hablamos de supercomputación, nos estamos refiriendo a que estas entidades, múltiples procesadores y memorias, existen físicamente, pues podemos encontrarnos todos estos elementos de manera virtual en ordenadores mas sencillos con sistemas operativos multiproceso.

La última tendencia en arquitecturas paralelas son las arquitecturas vectoriales, arquitecturas SIMD: una instrucción múltiples datos, debido a su enorme desarrollo en la evolución de las tarjetas gráficas, que es donde aparecen de manera natural al tener que aplicar la misma operación a múltiples píxeles. CUDA es un estandar de facto en la extensión de estas operaciones a cualquier tipo de datos no necesariamente gráficos.

Por tanto, en el mundo de los supercomputadores actuales es fácil encontrarnos con arquitecutras paralelas heterogeneas en donde conviven simultanemente la memoria compartida con la memoria distribuida de múltiples unidades de proceso que soportan un juego de instrucciones vectorial. La paralelización de estos códigos se realiza ad-hoc, pero, al igual que sucedió con el triunfo de las arquitecturas RISC sobre las CISC, la última palabra es posible que la tengan los compiladores, aunque es un hecho objetivo que la paralelización automática es instrínsecamente muy compleja.

Un compilador, aunque se suele asociar a los traductores de lenguajes de alto nivel a lenguaje máquina, en realidad son traductores entre dos lenguajes, sean del tipo que sean, por lo que es totalmente lógico pensar en compiladores de programas secuenciales en un lenguaje determinado a programas paralelos heterogeneos en el mismo lenguaje. Mientras se programó en lenguaje máquina, las arquitecturas CISC dominaron el mercado. Con la aparición de los lenguajes de alto nivel y de compiladores potentes para los mismos que solo utilizaban un subconjunto muy reducido del amplísimo conjunto de instrucciones disponibles, entraron en escena las arquitecturas RISC y desplazaron a las primeros.

Por otra parte, ¿cúal es la propuesta del GRID computing? Pués su aparición es, al igual que por ejemplo los protocolos TCP/IP de Internet, una solución de facto a un problema existente y es, por una parte, el hecho de que existan en un momento determinado multitud de recursos repartidos a lo largo del planeta, y por otra, a la posibilidad de compartirlos de manera transparente por multitud de usuarios igualmente dispersos.

Por ejemplo, la realidad es que un grupo de la Universidad de las Palmas de Gran Canaria tiene un supercomputador fruto de sus necesidades en un determinado momento y otro grupo de la University of Tasmania tiene otro con los mismos recursos. Si en un momento dado se dan cuenta de que podrían compartir sus recursos, lo cual a priori siempre es positivo desde el punto de vista de la teoría de juegos, la tecnología grid ofrece una capa de abstracción, la middleware grid, gracias a la cual pasamos a tener un único sistema con el total de los recursos.

Pensando trivialmente, solo por aclarar ideas, y sabiendo que es un caso totalmente irreal, si ambas máquinas solo se utilizasen durante las ocho horas laborables locales, supongamos de 8h a 18h, obviamente los dos grupos salen muy beneficiados, pués con una diferencia horaria de 10 horas no habría conflictos de acceso y los dos pasarían a disponer del doble de recursos de los que tenían.

Por tanto, los supercomputadores existen y existirán, pues son la mejor solución a nivel local, pero la tecnología GRID, sin entrar en consideraciones sobre Cloud, es la mejor solución a nivel global, ya que permite la interconexión transparente de estos óptimos locales que pueden llegar a ser muy heterogéneos.

Tags: , , , , , , , , , , , , , , , , ,

El teorema de los residuos es una parte fundamental de la variable compleja y es una generalización del teorema integral de Cauchy y la fórmula integral de Cauchy.

Fórmula integral de Cauchy

Sean $latex Omega subset mathbb{C}$ un abierto simplemente conexo ($latex I(gamma,z) = 0$ si $latex z notin Omega$), $latex f in mathcal{H}(Omega)$ y $latex gamma$ un ciclo tal que $latex gamma^* subset Omega$. Entonces, para $latex z in Omega-gamma^*$:

$latex frac{1}{2 pi i} int_gamma frac{f(u)}{u-z}du = f(z) I(gamma,z)$.

demostración:

Teorema integral de Cauchy.

Sean $latex Omega subset mathbb{C}$ un abierto simplemente conexo ($latex I(gamma,z) = 0$ si $latex z notin Omega$), $latex f in mathcal{H}(Omega)$ y $latex gamma$ un ciclo tal que $latex gamma^* subset Omega$. Entonces, para $latex z in Omega-gamma^*$:

$latex int_gamma f(u)du = 0$

demostración:

Fórmula integral de Cauchy para las derivadas.

Sean $latex Omega subset mathbb{C}$ un abierto simplemente conexo ($latex I(gamma,z) = 0$ si $latex z notin Omega$), $latex f in mathcal{H}(Omega)$ y $latex gamma$ un ciclo tal que $latex gamma^* subset Omega$. Entonces, para $latex z in Omega-gamma^*$:

$latex frac{n!}{2 pi i} int_gamma frac{f(u)}{(u-z)^{n+1}}du = f^{n)}(z) I(gamma,z)$

demostración:

Tags: , , ,

Para formular matemáticamente la mecánica cuántica necesitamos hablar de espacios de Hilbert y de operadores lineales. La rama de las matemáticas que trata estos temas es el análisis funcional.

Procederemos a dar cada una de las definiciones en el momento que las necesitemos, de manera que, como lo que nos interesan son los espacios de Hilbert y los operadores lineales, definiremos estos en primer lugar. Sin embargo, como estas definiciones se basan en otras definiciones, que por abstracción adquieren identidad propia, iremos introduciendo las mismas a medida que nos vayan apareciendo.

Espacios de Hilbert

Definición (Espacio de Hilbert): Un espacio de Hilbert es un espacio prehilbertiano completo.

A esto nos referiamos, pués de momento, nos hemos quedado tal y como estabamos, pues no sabemos que significa ser prehilbertiano y tampoco completo. Estos nuevos conceptos aparecen porque al estudiar los espacios de Hilbert y abstraer ciertas partes, estas adquieren interes per se.

Definición (Espacio prehilbertiano): Un espacio prehilbertiano es un espacio vectorial $latex H$ sobre un cuerpo $latex mathbb{K} = mathbb{R}, mathbb{C}$ en el que tenemos definido un producto escalar.

Definición (Producto escalar): Un producte escalar en un $latex mathbb{K}$-espacio vectorial $latex H$ es una aplicación:

$latex langle cdot , cdot rangle: H times H longrightarrow mathbb{K}$

tal que es:

  1. lineal en la primera componente: $latex langle x + y , z rangle = langle x, z rangle + langle y, z rangle$ y $$
  2. simétrica en $latex mathbb{R}$,  $latex langle x,y rangle = langle y,x rangle$, o hermítica en $latex mathbb{C}$, $latex langle x,y rangle = overline{ langle y,x rangle}$.
  3. definida positiva: $latex langle x,x rangle geq 0$ y $latex langle x,x rangle = 0 Leftrightarrow x=0$

Ejemplos

  • $latex (mathbb{R}^n,langle cdot, cdot rangle)$ con $latex langle x,y rangle = sum_{i=1}^n x_i y_i$ es un espacio prehilbertiano.
  • $latex (mathbb{C}^n, langle cdot,cdot rangle)$ con $latex langle z,w rangle = sum_{i=1}^n z_i overline{w_i}$ es un espacio prehilbertiano.
  • Sea $latex A$ un conjunto igual a $latex { 1, ldots , n}, mathbb{N}, mathbb{Z}$. El espacio $latex ell^2(A)$ es un espacio prehilbertiano con el producto escalar $latex langle x,y rangle := sum_{alpha in A} x_alpha overline{y_alpha}$ con $latex x = { x_alpha }_{alpha in A}$ y $latex y = { y_alpha }_{alpha in A}$. En el caso que $latex A = mathbb{N}$ o $latex mathbb{Z}$ escribiremos $latex ell^2(A) = ell^2$. Si $latex A = {1,ldots,n}$ entonces $latex ell^2(A) equiv mathbb{K}^n$ que tambien se denota por $latex ell^2(n)$.

Vamos a demostrar este último ejemplo. Probaremos, en primer lugar, que $latex ell^2$ es un espacio vectorial, y posteriormente, que $latex langle cdot,cdot rangle$ es un producto interior.

Definición (Notación de Dirac): Una notación alternativa para el producto escalar introducida por Dirac y ampliamente utilizada en mecánica cuántica por su versatilidad es la notación bra-ket. Podemos escribir el producto escalara en notación matricial:

$latex langle cdot | cdot rangle : $

Espacios de Banach

Para empezar, recordaremos que un espacio normado es un par $latex (E,|| cdot ||_E)$ formado por un espacio vectorial $latex E$ sobre un cuerpo $latex K$, que puede ser $latex mathbb{R}$ o $latex mathbb{C}$, y una norma $latex ||cdot||_E$ sobre $latex E$, que es una función:

$latex ||cdot||_E : E longrightarrow [0,+infty[$

que satisface las propiedades de:

  1. separación: $latex ||x||_E = 0 Leftrightarrow x = 0$.
  2. homogeneidad: $latex ||lambda x|| = |lambda|||x||_E$.
  3. desigualdad triangular:

Tags: , ,

Uno de los grandes logros de la variable compleja es su efectividad a la hora de permitir demostrar resultados muy alejadas de su aparente alcance. Un claro ejemplo de esto es su aplicación a la demostración del teorema fundamental del algebra. A continuación presentamos una demostración del mismo junto con unos resultados previos que necesitaremos para la misma.

Desigualdad de Cauchy

Sea $latex f$ una función analítica, es decir, $latex f = sum_{n=0}^infty a_n(z-z_0)^n$ con $latex z in D(z_0,R)$. Sea $latex 0 < r < R$ y

$latex M(r) := max (|f(z)|:|z-z_0|=r)$.

Entonces

$latex |a_n| leq frac{M(r)}{r^n}$

demostración:

$latex a_n = frac{f^{n)}(z_0)}{n!}$

ya que $latex f'(z) = sum_{n=1}^infty n,a_n (z-z_0)^{n-1}$, por lo que $latex f'(z_0) = 1 cdot a_1$. De la misma manera, $latex f”(z) = sum_{n=2}^infty (n-1)n,a_n (z-z_0)^{n-2}$, por lo que $latex f”(z_0) = 2 cdot 1 cdot a_2$, etc.

Si aplicamos ahora la fórmula integral de Cauchy para las derivadas:

$latex f^{n)}(z) = frac{n!}{2 pi i} int_{C(z_0,R)} frac{f(u)}{(u-z)^{n+1}} du,, forall z in D(z_0,R)$

tenemos:

$latex a_n = frac{1}{2 pi i} int_{C(z_0,r)} frac{f(z)}{(z-z_0)^{n+1}} dz$

por lo que si aplicamos la propiedad de que $latex |int_gamma f(z) dz| leq L(gamma) , max { |f(z)| : z in gamma([a,b]) }$, nos queda:

$latex |a_n| leq big | frac{1}{2 pi i} big | (2 pi r) max_{|z-z_0|=r} big | frac{f(z)}{(z-z_0)^{n+1}} big |$.

Como $latex |frac{1}{2 pi i}| = frac{1}{2 pi}$ y

$latex max_{|z-z_0|=r} big | frac{f(z)}{(z-z_0)^{n+1}} big | = max_{|z-z_0|=r} frac{|f(z)|}{r^{n+1}}$

entonces:

$latex |a_n| leq frac{r}{r^{n+1}} max_{|z-z_0|=r} |f(z)| = frac{M(r)}{r^n} ,Box$

Teorema de Liouville

Sea $latex f in mathcal{H}(mathbb{C})$, $latex f$ no constante. Entonces $latex f$ no está acotada.

demostración:

Vamos a demostrar el reciproco, es decir, que si $latex fin mathcal{H}(mathbb{C})$ y $latex f$ está acotada, entonces $latex f$ es necesariamente constante.

Para empezar, como $latex f$ es holomorfa, entonces es analítica y

$latex f(z) = sum_{n=0}^infty a_n z^n, , forall z in mathbb{C}$.

Fijamos $latex 0 < r < +infty$. Si le aplicamos ahora la desigualdad de Cauchy, nos queda

$latex |a_n| leq frac{M(r)}{r^n}$.

que, como $latex f$ está acotada, es decir, $latex exists C > 0: |f(z)| leq C, , forall z in mathbb{C}$, que no depende de $latex r$, nos queda:

$latex |a_n| leq frac{C}{r^n}$.

Finalmente, si $latex n geq 1$ entonces $latex |a_n| leq lim_{r rightarrow infty} frac{C}{r^n} = 0$, por lo que $latex f(z) = a_0$ que implica que $latex f$ es constante $latex forall z in mathbb{C} , Box$

Lo que nos está diciendo el teorema de Liouville, por ejemplo, es que la función $latex sin(z)$ no está acotada (a diferencia de lo que ocurre en variable real), pues solo lo estan las funciones constantes.

Teorema fundamental del Álgebra

Sea $latex p(z) = a_0 + a_1 z + a_2 z^2 + ldots + a_n z^n$ un polinomio no constante. Entonces la ecuación $latex p(z)=0$ tiene solución.

demostración:

Vamos a proceder por reducción al absurdo. Supondremos que $latex p(z) neq 0, , forall z in mathbb{C}$ y consideraremos la función $latex f(z) = frac{1}{p(z)}$. Es fácil ver que $latex f(z) in mathcal{H}(mathbb{C})$ y que no es constante.

Veremos que $latex f(z)$ está acotada, lo que supondra una contradicción con Liouville. Efectivamente, para $latex p(z) = a_0 + a_1 z + a_2 z^2 + ldots + a_n z^n$, como $latex n geq 1$ y $latex a_n neq 0$, tiene sentido escribir $latex p(z)$ como:

$latex p(z) = z^n (a_n + frac{a_{n-1}}{z} + ldots + frac{a_0}{z^n})$

con lo que, si tomamos valores absolutos y utilizando que $latex |a+b| geq |a|-|b|$, nos queda:

$latex |p(z)| geq |z|^n (|a_n| – frac{|a_{n-1}|}{|z|} – ldots – frac{|a_0|}{|z|^n})$

$latex ,Box$

Proposición:

Sea $latex p(z)$ un polinomio de grado $latex n geq 1$. Entonces:

$latex p(z) = a(z-z_1)^{k_1}(z-z_2)^{k_2}ldots (z-z_m)^{k_m}$ con $latex sum_{i=1}^m k_i = n$.

demostración: (por inducción sobre $latex n$)

a) Si $latex n=1$ entonces $latex p(z)=az+b$ con $latex a neq 0$. Sea $latex z_1$ tal que $latex p(z_1) =0$. Entonces:

$latex p(z) = az -az_1 + az_1 + b = a(z-z_1)$, pues $latex az_1 + b = 0$.

b) Suponemos que para grado $latex leq n$ es cierto y sea $latex p(z)$ un polinomio de grado $latex n+1$.

Tags: , ,

Fórmula integral de Cauchy para las derivadas.

Sea $latex f(z) in mathcal{H}(A)$ y $latex overline{D(z_0,R)} subset A$. Entonces:

$latex f^{n)}(z) = frac{n!}{2 pi i} int_{C(z_0,R)} frac{u}{u-z}du, , forall z in D(z_0,R)$.

demostración:

Sabemos que $latex f(xi) = frac{1}{2 pi i} int_{C(z_0,R)} frac{f(u)}{u-xi}du$ siempre que $latex |xi – z_0|<R$.

Fijamos ahora $latex z in D(z_0,R)$ y escogemos $latex 0 < r < R: D(z,r) subset D(z_0,R)$. Tenemos un teorema que nos asegura que si $latex |xi – z| < r$ entonces:

$latex f(xi) = frac{1}{2 pi i} sum_{n=0}^infty (int_{C(z_0,R)} frac{f(u)}{(u-z)^{n+1}} du) (xi – z)^n$,

por lo que:

$latex frac{f^{n)}(z)}{n!} = frac{1}{2 pi i} int_{C(z_0,R)} frac{u}{u-z}du$

Aplicaciones:

1. Calcular

$latex int_{C(0,frac{1}{2})} frac{e^z}{z(1-z)^2}dz$.

En este caso, la única singularidad que queda dentro de la circunferencia es el $latex 0$. Tenemos que aislarla, por lo que escribimos la integral de la siguiente forma:

$latex int_{C(0,frac{1}{2})} frac{frac{e^z}{(1-z)^2}}{z-0} dz = int_{C(0,frac{1}{2})} frac{f(z)}{(z-0)}$

con $latex f(z):=frac{e^z}{(1-z)^2}$. Ahora tenemos que $latex f(z)$ es holomorfa en todos los puntos excepto en el $latex 1$, que no está dentro de la circunferencia que consideramos en la integral i, por tanto, podemos aplicar el teorema integral de Cauchy para las derivadas y obtenemos

$latex int_{C(0,frac{1}{2})} frac{f(z)}{z} = 2 pi i f(0) = 2 pi i$

ya que $latex f(0) = frac{e^0}{(1-0)^2} = 1$.

2. Calcular

$latex int_{C(1,frac{1}{2})} frac{e^z}{z(1-z)^2}dz$.

Ahora, la singularidad que queda dentro de la circunferencia es el $latex 1$ y tenemos que aislarlo para poder aplicar el teorema integral para derivadas. Escribimos

$latex int_{C(1,frac{1}{2})} frac{frac{e^z}{z}}{(z-1)^2} = 2 pi i f'(1)$

aplicando el teorema. De esta manera, como $latex f'(z) = frac{ze^z – e^z}{z^2}$ y $latex f'(1)=0$, donde $latex f(z):=frac{e^z}{z}$, nos queda que la integral vale $latex 0$.

3. Calcular

$latex I:=int_{C(0,2)} frac{e^z}{z(1-z)^2}dz$.

Ahora las dos singularidades están dentro de la circunferencia considerada. Lo que haremos es decomponer $latex frac{1}{z(z-1)^2}$ en fracciones simples:

$latex frac{1}{z(z-1)^2} = frac{A}{z} + frac{B}{z-1} + frac{C}{(z-1)^2}$.

Esto nos lleva a un sistema cuya solución es $latex A=1=C, B=-1$, por lo que podemos escribir:

$latex I = int_{C(0,2)} frac{e^z}{z} dz – int_{C(0,2)} frac{e^z}{z-1} dz + int_{C(0,2)} frac{e^z}{(z-1)^2}$

de manera que podemos aplicar la fórmula integral de Cauchy para las derivadas, con $latex f(z):=e^z$ que es holomorfa en todo el plano, a cada integral:

$latex I = 2 pi i (f(0) – f(1) + f'(1)) = 2 pi i (1 – e + e) = 2 pi i$.

Tags: ,

FireStats icon Powered by FireStats