PHYSICAL PHENOMENON

$\Downarrow$

MATHEMATICAL MODEL

$\Downarrow$

DISCRETE ALGEBRAIC APPROXIMATION

$\Downarrow$

NUMERICAL ALGORITHMS

$\Downarrow$

SIMULATION PROGRAM

$\Downarrow$

COMPUTER EXPERIMENT

SIMULACIONES NUMERICAS COSMOLOGICAS

Motivaciones

-La GRAVITACION es la interacción dominate en la formación y evolución de objetos astronómicos. (cúmulos estelares, galaxias , cúmulos de galaxias y la estructura en gran escala del Universo).

-En el contexto cosmológico es una maquinaria de hacer crecer perturbaciones. Evolución altamente no lineal.

-Desafortunadamente, no hay solución analítica general con $N > 2$.

-Simulaciones numéricas cosmológicas: algo asi como el laboratorio de la cosmologia

-Testeo de aproximaciones lineales (o de orden mayor) aproximación lineal de Zeldovich (1972), aproximación de segundo orden de Scoccimarro (1997), formalismo de Press-Schechter.

RECETAS DE COCINA PARA SIMULACIONES N-CUERPOS COSMOLOGICAS

-softening grande para minimizar el two-body relaxation.

-paso de tiempo pequeño para eliminar el two-body scattering de los errores de integración $\Delta t < 0.1 (G \rho)^{-1/2}$

-errores en las fuerzas menores que el 0.5 %.

-integrador de segundo orden y reversible temporalmente.

-volumen grande para modelizar todos los efectos no-lineales. Los filamentos afectan la evolución de los objetos mas peque nos).

- redshift inicial grande para que todas las escalas representadas esten aun en su regimen lineal.

N-CUERPOS

La ley de Newton para una distribución discreta de $N_P$ partículas es:

$${\bf F}_i = - G \sum_{i \ne j=1}^{N_p} \frac{{\bf r}_i- {\bf r}_j}{\vert{\bf r}_i-{\bf r}_j\vert^3} m_j$$ (1)

definiendo

$$\Phi_i = - G \sum_{i \ne j=1}^{N_p} \frac{m_j}{\vert{\bf r}_i-{\bf r}_j\vert}$$ (2)

se tiene que

$${\bf F}_i=-{\bf\nabla} \Phi_i$$ (3)

$${\bf F}_{ij}=-{\bf F}_{ji}$$ (4)

$${\bf F}_{ii}=0$$ (5)

$softening$

$${\bf F}_i = - G \sum_{i \ne j=1}^{N_p} \frac{{\bf r}_i- {\bf r}_j}{(\epsilon ^2+\vert{\bf r}_i-{\bf r}_j\vert^2)^{3/2}} m_j$$ (6)

Para una distribución continua:

$$\Phi({\bf r}) = - G \int \frac{\rho({\bf r}')}{\vert{\bf r}-{\bf r}'\vert} dV'$$ (7)

aplicando la Transformada de Fourier

$$\stackrel {\sim} {f}({\bf k}) = \int f({\bf r}) e^{-i{\bf k}.{\bf r}} dV$$ (8)

a ambos miembros y aplicando el teorema de convolución

$$\stackrel {\sim} {\Phi}({\bf k}) = \stackrel {\sim} {\rho}({\bf k}) \stackrel {\sim} {F}({\bf k})$$ (9)

donde

$$F(\vert{\bf r}-{\bf r}'\vert)=-\frac{G}{\vert{\bf r}-{\bf r}'\vert}$$ (10)

es la función de Green

$$k_{max}=\frac{2\pi}{\lambda_{min}}=\frac{2\pi}{2\lambda_{grid}}=\frac{2\pi}{L} \frac{N}{2}$$ (11)

$$k_{min}=\frac{2\pi}{\lambda_{max}}=\frac{2\pi}{L}$$ (12)

ESQUEMA PARTICLE MESH

1) Asignación de masas de las partículas a los nodos del grid $\rightarrow \rho$ (NGP, CIC, TSC)

$$\rho({\bf r}_p)=\sum _{i=1}^{N_p} m_i W({\bf r}_i-{\bf r}_p)$$ (13)

2) Solución de la ecuación de Poisson $\nabla^2 \phi =4 \pi G \rho$ en los nodos:

unidimensional

$$\frac{\phi_{p-1}-2\phi_{p}+\phi_{p+1}}{H^2} = \rho _{p}$$ (14)

o

FFT

3) Cómputo de la fuerza en los nodos a partir del potencial por diferencias finitas

$${\bf f}_p =\frac{\phi_{p-1}-\phi_{p+1}}{H}$$ (15)

4) Interpolación de la fuerza en los nodos a las posiciones de las partículas

$${\bf f}_i=\sum _{i=1}^{N_g} {\bf f}_p W({\bf r}_i-{\bf r}_p)$$ (16)

5) Integrar la ecuaciones de movimiento

$${\bf v}^{n+1/2}={\bf v}^{n-1/2}+{\bf f}^{n} \Delta t$$ (17)

$${\bf r}^{n+1}={\bf r}^{n}+{\bf v}^{n+1/2} \Delta t$$ (18)

TECNICAS DE N-CUERPOS

-particle-particle:

suma directa particulas.

el CPU-time escala como N(N-1)

Fuerzas de rango corto y largo

pasos de tiempo individuales

softening

No hay hipótesis sobre la geometría del sistema

-tree code:

promedia partículas lejans como un solo objeto en el centro de masa.

Más dificil programación

N log(N)

no hay condiciones de contorno

softening

-particle-mesh:

mas eficiente para N grandes

el tiempo de CPU

-particle-particle-particle-mesh: P$^3$M

es un PP+PM

La parte PP es una corrección a las fuerzas de escalas grandes computadas por un PM.

La parte PP se vuelve muy costosa cuando se desarrolla el clustering

-adaptative particle-particle-particle-mesh: AP$^3$M refinamiento de los grids en las zonas de alta densidad

-expansión de potencial:

pura, armónicos esféricos, centros dobles

expansión truncada de la distribución de densidad en algunas funciones base

falta de flexibilidad.

LIMITES DE RESOLUCION

fuerza Depende del softening y determina la la escala espacial más pequeña donde los resultados son confiables.

masa Depende de la masa por partícula y determina los tamaños de los objetos que podrán ser analizados.

tiempo Depende del paso de integración y debe ir de acuerdo con la resolución de la fuerza. Mayor resoluciín en la fuerza mayor resolución temporal.

TEST DE SIMULACIONES

-comparación con soluciones analíticas de la teoría lineal de perturbaciones y de orden mayor

-colapso esférico

-resultados de simulaciones free-scale

-comparaciones de resultados de diferentes códigos

-Estabilidad de los códigos individualmente, convergencia de los resultados con el paso de tiempo, número de partículas, tamaño del sistema, resolución de la fuerza.

COORDENADAS COMOVILES

$$\stackrel {\sim} {\bf r} = {\bf r} / a$$ (19)

$a(t) \equiv$ factor de expansión $\sim t^{2/3}$

$${\bf v} = \dot {a} \stackrel {\sim} {\bf r} + a \stackrel {\sim} {\bf v}={\bf v}_{H} + {\bf v}_{P}$$ (20)

$$\dot {\bf v} = \stackrel{..}{a} {\bf r}+2 \dot {a} \stackrel... ...m}{\bf v} + a \dot {\stackrel {\sim}{\bf v}}=-{\bf\nabla} \Phi$$ (21)

$$\frac{\partial \Phi}{\partial x}=\frac {\partial \Phi}{\part... ...\frac{\partial \Phi}{\partial \stackrel{\sim}{x}} \frac{1}{a}$$ (22)

luego

$${\bf\nabla} \Phi = \stackrel {\sim} {\bf\nabla } \Phi /a$$ (23)

definiendo

$$\stackrel {\sim}{\Phi} = a \Phi + \frac{1}{2} a^2 \stackrel{..}{a} \stackrel{\sim}{r}^2$$ (24)

$$2 \frac{\dot a}{a} \stackrel {\sim}{\bf v} + \dot {\stackrel... ...{\bf v}} =- \stackrel{\sim}{\bf\nabla}{\stackrel {\sim}{\Phi}}$$ (25)

ECUACION DE POISSON

EN COORDENADAS COMOVILES

$${\bf\nabla} ^2 \Phi = 4 \pi G \rho$$ (26)

$${\bf\nabla} ^2 = \stackrel {\sim} {\bf\nabla} ^2 /a^2$$ (27)

$$\stackrel {\sim} {\bf\nabla} ^2 \stackrel{\sim}{\Phi}= a^3 {\bf\nabla}^2 {\Phi} +3 a^2 \stackrel {..}{a}$$ (28)

y como

$$\stackrel {..}{a} = - G M /a^2 = - \frac{4}{3} \pi G a \rho_b$$ (29)

$$\stackrel {\sim} {\bf\nabla}^2 \stackrel{\sim}{\Phi}= 4 \pi G (\stackrel{\sim}{\rho} -\stackrel{\sim}{\rho}_b)$$ (30)

GENERACION DE CONDICIONES INICIALES

Aproximación de Zeldovich

Es un método lineal analítico

$${\bf r} = a(t)({\bf q } - b(t) {\bf\Psi}({\bf q}))$$ (31)

$${\bf v} = \frac {\dot a }{a} {\bf r} -a \dot b {\bf\Psi}$$ (32)

en la aproximación lineal

$$\delta = {\bf\nabla} . {\bf\Psi}$$ (33)

donde

$$\delta({\bf q})=\frac {\rho({\bf q})-\rho_b}{\rho_b}= \sum \delta_{\bf k} e^{i {\bf k}.{\bf q}}$$ (34)

con

$$P(k)=<\vert\delta_{\bf k}\vert^2>$$ (35)

SIMULACIONES COSMOLOGICAS HIDRODINAMICAS

Motivaciones

-La gravitación describe correctamente solamente el comportamiento de la materia oscura.

-Simulaciones $\rightarrow$ materia oscura

-Observaciones $\rightarrow$ bariones

-La gravitación NO describe cuando y donde se forman las galaxias.

-En pequeña escala ( $<$ 5Mpc) los efectos del gas son dominantes.

-Las simulaciones deben incluir una descripción "adecuada" del gas.

-Testear si los modelos producen galaxias y sistemas de galaxias realísticas.

-Correlación entre la distribución de materia oscura y las galaxias (parámetro de $biasing$ como función de la escala)

-Formación y evolución de galaxias y sistemas de galaxias

-X-Rays y efecto Sunayev-Zeldovich

-Propiedades de absorción del IGM (Ly $\alpha$ forest)

DIFERENCIAS ENTRE GAS Y MATERIA OSCURA

-El gas NO se interpenetra como la materia oscura en una colisión

-Las ondas de choque transforman su energía cinética en energía térmica

-La ecuación de movimiento incluye gradiente de presión

-Pierde energia por enfriamiento radiativo

-Puede convertirse en estrellas

ECUACIONES PARA EL GAS

$$\frac{d {\bf v}}{ d t} = - {\bf\nabla}{\Phi} - \frac{{\bf\nabla }{\Phi}}{\rho}$$ (36)

donde

$$P=P+P_{vis}$$ (37)

$$P=(\gamma-1) \rho \epsilon \; \; (\gamma=5/3)$$ (38)

$$\frac{\partial \rho}{\partial t} = - {\bf\nabla} . (\rho {\bf v})$$ (39)

$$\frac{\partial \epsilon}{\partial t} = - \frac {P}{\rho} {\bf\nabla} . {\bf v}+ \frac{\Lambda(\epsilon,\rho)}{\rho}$$ (40)

METODOS COMPUTACIONALES

-Método Euleriano:

Resuelve las ecuaciones del gas en un grid fijo por diferencias finitas.

Puede incluir diferente número de celdas en las regiones de alta densidad.

-Método Lagrangiano:

Utiliza generalmente el algoritmo SPH (Smoothd Particle Hydrodynamics)

Se representa un elemento de fluido por partículas que llevan información dinámica ( ${\bf r} , {\bf v}$) y térmica ($\rho, T, P$). Elimina la necesidad de un grid.

Se computan la densidad y gradientes de presión por sumas suavizadas sobre las partiículas vecinas.

Usa viscosidad artificial para evitar la interpenetración de partículas.

La longitud de suvizado escala con la densidad, mayor resolución en las zonas de alta densidad. .

$$\rho = \sum m_i W(r_i/h)$$ (41)

$${\bf\nabla} P = \sum P_i {\bf\nabla} W(r_i/h)$$ (42)

$W(r/h) \equiv$ kernel interpolatorio esférico.

función derivable con un máximo en $r=0$ y un decaimiento rápido. No es necesario tomar la contribución de todos los vecinos (N $>$ 30 ).

Shapiro et al. (1996) propone un kernel asimétrico.

GRAVITACION + HIDRODINAMICA

PP+SPH

(implementación paralela)

TREE+SPH

(Hernquist & Katz, Steinmetz 1996, Navarro & White 1993)

HPM+SPH

( Shapiro et al. 1996)

P$^3$M+SPH

(Evrard 1988)

AP$^3$M+SPH

(Couchman et al. 1995, Tissera, Lambas & Abadi 1997)

CONSIDERACIONES NUMERICAS

-El paso de tiempo debe satisfacer la condición de Courant $\Delta t < 0.3 h/c$ donde $h$ es la longitud se suavizado y $c^2 \equiv \frac{\partial P}{\partial \rho} = (\gamma -1) \epsilon$

-La resolución espacial y temporal están estan ligadas. En general se usa un mínimo $\Delta t$ entonces hay que agrandar $h$

-El $cooling$ es muy problemático, porque el tiempo de cooling puede ser muy pequeño

$$t_{cool}=\frac{\epsilon}{\Lambda} << t_{dyn}= (G\rho)^{-1/2}$$ (43)

TESTS DE SIMULACIONES HIDRODINAMICAS

-Colapso esférico con soluciones de autosimilaridad , Bertschinger (1989). Habilidad de reproducir correctamente los $shocks$

-Idem + $cooling$ autosimilar (Abadi, Bower & Navarro, 1998)

-Tests de convergencia de los códigos individuales y entre ellos (Frenk et al. 1998 comparación de la formación de un cúmulo de galaxias con diferentes códigos )

FORMACION ESTELAR

-Parte fundamental y crítica ya que determina la estructura de las "galaxias" que se forman en una simulación.

-Una región forma estrellas si:

es un sistema colapsante ${\bf\nabla}.{\bf v} < 0$

la presión no puede soportar a la gravitación $(T_{sound} > T_{dyn})$ Inestabilidad de Jeans

se enfria rápidamente y se mantiene fría a medida que colapsa $T_{cool}<<T_{dyn}$

-Por la forma de la función de $cooling$ las condiciones anteriores se satisfacen con un criterio de densidad simple $\rho > \rho_{crit} \simeq 7 \times 10^{-26} g \; cm^{-3}$

-Tasa de formación estelar $-\dot \rho_{gas} = \dot \rho_{star} = \rho_{gas}/t_{\star}$ con $t_{\star}$ escala de tiempo característica

-La masa de las estrellas estan distribuidas por la IMF.

-$Feedback$: como las explosiones de supernova afectan al medio. La energia se inyecta al gas que la rodea como energia cinética y/o térmica. (Extremadamente incierto). Puede suprimir la formación estelar de halos pequeños $V_c < 100 km/s$.

Bibliografia

-Hockney & Eastwood, Computer Simulations Using Particles 1981

-Efstathiou G., Davis M., Frenk C.S. & White S.D.M., Numerical Techniques for Large Cosmological N-Body Simulations, ApJS 57, 241 1985

-Katz N., Weinberg D. & Hernquist L., Cosmological Simulations with TreeSPH, ApJ Supp, 105, 19. astroph/9509107 1996

-Merchán, M., Trabajo Especial, Facultad de Matemática Astronomía y Física, Universidad Nacional de Córdoba, 1995

-Steinmetz M., Simulating Galaxy Formation, astro-ph/9512013, 1996

-Yepes G., From Quantum Fluctuations to Cosmological Structures, ASP Conference Series, 126, 279, 1997

- Weinberg D.,Katz N., & Hernquist L., Simulating Cosmic Structure Formation astro-ph/9708213, 1997

UW HPCC/ESS Project

http://www-hpcc.astro.washington.edu/

Grand Challenge Cosmology Consortium

http://zeus.ncsa.uiuc.edu:8080/GC3_Home_Page.html

The Virgo Consortium

http://star-www.dur.ac.uk/ frazerp/virgo/virgo.html