User: knarf Date: 2011/10/03 12:06 AM

Modified: / ET.tex /examples/tov/ ham_conv.pdf, ham_conv.py

Log: add convergence plot to paper

File Changes:

Directory: / ============

File [modified]: ET.tex Delta lines: +49 -26 =================================================================== --- ET.tex 2011-10-03 05:05:49 UTC (rev 140) +++ ET.tex 2011-10-03 05:06:51 UTC (rev 141) @@ -2557,8 +2557,8 @@ The examples in the previous subsections did not include the evolution of matter within a relativistic spacetime. One interesting test of a coupled matter-spacetime evolution is to measure the eigenfrequencies of a stable TOV -star. These eigenfrequencies can be compared to values known from linear -perturbation theory. +star~\todo{cite papers also doing this}. These eigenfrequencies can be compared +to values known from linear perturbation theory.

This test uses an equillibrium configuration of a self-gravitating fluid sphere, described by a polytropic equation of state. This one-dimensional @@ -2568,21 +2568,11 @@ \codename{McLachlan} and the hydrodynamics evolution system implemented in \codename{GRHydro}.

-Despite the system being in equilibrium, the solution experiences oscillations -mainly due to perturbations from two sources, the interpolation of the initial -data onto the evolution grid, and due to the interaction of the star and the -atmosphere. The former is damped quickly, and can be reduced by increasing the -resolution of the evolution grid. It can be seen as initial spike in the -evolution of the central density of the star. The remaining perturbations, -mainly from the atmosphere treatment, are present during the whole evolution, -continuously exciting oscillations. Given enough evolution time, the -frequencies of these oscillations can be measured with satisfactory accuracy. - Specifically for this test, we set up a stable TOV star described by a polytropic equation of state $p=K\varrho^\Gamma$ with $K=100$ and $\Gamma=2$, and an initial central density of $\varrho_0=1.28\times10^{-3}$. This model can be taken to represent a non-rotating neutron star with a mass of -$M=1.4\mathrm{M}_\odot$. The computational domain is a cube of length +$M=1.4\mathrm{M}_\odot$. The computational domain is a cube of length $640\mathrm{M}$ and a base resolution of $2\mathrm{M}$ ($4\mathrm{M}$, $8\mathrm{M}$) in each dimension. Four additional grids refine the region around the star located at the origin, each doubling the resolution, with sizes @@ -2593,10 +2583,10 @@ In figure~\ref{fig:tov_rho_max} we show the evolution of the central density of the star, over an evolution time of $1300\mathrm{M}$. The initial spike is due to the perturbation of the solution resulting from the interpolation onto the -evolution grid. These oscillations are quickly damped. Most of the remaining -oscillations are due to perturbations resulting from the treatment of the -atmosphere. The mean of the central density shows a drift, which converges -to zero with increasing resolution.\todo{show} +evolution grid. The remaining oscillations are due to the interaction of the +star and the artificial atmosphere and are present during the whole evolution. +Given enough evolution time, the frequencies of these oscillations can be +measured with satisfactory accuracy.

\begin{figure} \label{fig:tov_rho_max} @@ -2609,15 +2599,18 @@ are three different resolutions.} \end{figure}

-In figure~\ref{fig:tov_mode_spectrum} we show the power spectral density (PSD) of the central -density oscillations computed from a full 3D relativistic hydrodynamics -simulation, compared to the coresponding frequencies as obtained with -perturbative techniques, kindly provided to use by Kentaro Takami and computed -using the method described in~\cite{Yoshida01}. The PSD was computed using the entire -time series, removing the linear trend, averaging over Hanning windows overlapping half the -signal length after padding the signal to five time it's length. -The agreement of the frequencies of the fundamental mode and the first three overtones -is clearly visible. +In figure~\ref{fig:tov_mode_spectrum} we show the power spectral density (PSD) +of the central density oscillations computed from a full 3D relativistic +hydrodynamics simulation, compared to the corresponding frequencies as obtained +with perturbative techniques, kindly provided to use by Kentaro Takami and +computed using the method described in~\cite{Yoshida01}. The PSD was computed +using the entire time series of the high-resolution run, removing the linear +trend, averaging over Hanning windows overlapping half the signal length after +padding the signal to five time it's length. The agreement of the frequencies +of the fundamental mode and the first three overtones is clearly visible, but +we are not able to show the same for the fourth overtone, at this resolution. +We expect to be able to resolve higher overtones when using even higher +resolution, but did not try because of the involved computational cost.

\begin{figure} \label{fig:tov_mode_spectrum} @@ -2627,7 +2620,37 @@ to the values obtained by perturbation theory.\todo{improve caption}} \end{figure}

+Within this test it is also interesting to study the convergence behavior of +the coupled curvature and matter evolution code. One of the variables often +used for this test is the Hamiltonian constraint violation. This violation +vanishes for the continuum problem, but is non-zero and resolution-dependend in +discrete simulations. The expected rate of convergence if the hydrodynamics +code is between $1$ and $2$. It cannot be higher than $2$ because of the used +directional flux-split algorithm which is of second order, and in some cases +the hydrodynamics code is only of first order, in particular at extrema (like +the center of the star), or at the star surface.

+Figure~\ref{fig:tov_ham_conv} shows the order of convergence of the Hamiltonian +constraint violation, using the two highest-resolution runs, at the stellar +center and a coordinate radius of $r=5$ which is about half-way between the +center and the surface. The observed convergence rate is most of the time +between $1.4$ and $1.5$ at the center, and between $1.6$ and $2$ inside the +star. This difference is expected because of the mentioned data-dependend +convergence orders of the underlying hydrodynamics evolution schemes. + +\begin{figure} + \label{fig:tov_ham_conv} + \includegraphics[width=0.9\textwidth]{examples/tov/ham_conv} + \caption{Convergence factor of Hamiltonian constraint violation at $r=0$ and + $r=5$. The observed convergence order of about $1.5$ at the center of + the star is lower then the general second order of the hydrodynamics + evolution scheme. This is expected because the scheme's convergence rate drops to first + order at extremas or shocks, like the stellar center or the star surface. + Consequently, the observed convergence order about half-way between the + stellar center and surface is higher than $1.5$, but most of the time below + $2$.} +\end{figure} + \subsection{Collapse\pages{2 Christian}} Show TOV collapse and BH formation \todo{RH: there are three runs for this for convergence testing. All three will be

Directory: /examples/tov/ =========================

File [modified]: ham_conv.pdf Delta lines: +0 -0 =================================================================== (Binary files differ)

File [modified]: ham_conv.py Delta lines: +16 -9 =================================================================== --- examples/tov/ham_conv.py 2011-10-03 05:05:49 UTC (rev 140) +++ examples/tov/ham_conv.py 2011-10-03 05:06:51 UTC (rev 141) @@ -6,14 +6,17 @@

# stuff xlim = (0,1400) -ylim = (1.1,1.6) +ylim = (1,2.2)

# load data -Ft, F4 = np.loadtxt("tov_4/ml_bssn::ml_ham.x.asc.x0.bz2", comments="#", usecols=(8,12), unpack=True) -Ft, F2 = np.loadtxt("tov_2/ml_bssn::ml_ham.x.asc.x0.bz2", comments="#", usecols=(8,12), unpack=True) +Ft, F04 = np.loadtxt("tov_4/ml_bssn::ml_ham.x.asc.x0.bz2", comments="#", usecols=(8,12), unpack=True) +Ft, F02 = np.loadtxt("tov_2/ml_bssn::ml_ham.x.asc.x0.bz2", comments="#", usecols=(8,12), unpack=True) +Ft, F54 = np.loadtxt("tov_4/ml_bssn::ml_ham.x.asc.x5.bz2", comments="#", usecols=(8,12), unpack=True) +Ft, F52 = np.loadtxt("tov_2/ml_bssn::ml_ham.x.asc.x5.bz2", comments="#", usecols=(8,12), unpack=True)

-for i in range(len(F2)): - F2[i] = math.log(abs(F4[i]/F2[i]),2) +for i in range(len(F02)): + F02[i] = math.log(abs(F04[i]/F02[i]),2) + F52[i] = math.log(abs(F54[i]/F52[i]),2)

# Plot basics fig = plt.figure(figsize=(6, 3)) @@ -21,8 +24,12 @@ ax = fig.add_subplot(1,1,1)

# plot -ax.plot(Ft, F2, linestyle='-', color='black') +ax.plot(Ft, F52, linestyle='-' , color='blue' , label=r'$r=5$') +ax.plot(Ft, F02, linestyle='--', color='black', label=r'$r=0$')

+handles, labels = ax.get_legend_handles_labels() +ax.legend(handles, labels, loc='best') + # plot properties ax.set_xlim(xlim) ax.set_ylim(ylim) @@ -38,9 +45,9 @@ ax2.set_xlim((ax.get_xlim()[0]/M_to_ms, ax.get_xlim()[1]/M_to_ms)) ax2.xaxis.set_major_locator(mticker.MaxNLocator(7)) ax2.xaxis.set_minor_locator(mticker.MaxNLocator(14)) -ax.set_ylabel(r'$\left|H\left(x=0\right)\right|$') -ax.yaxis.set_major_locator(mticker.MaxNLocator(3)) -ax.yaxis.set_minor_locator(mticker.MaxNLocator(6)) +ax.set_ylabel(r'$\left|H\right|$') +ax.yaxis.set_major_locator(mticker.MaxNLocator(7)) +ax.yaxis.set_minor_locator(mticker.MaxNLocator(14)) ax.yaxis.grid(False) set_tick_sizes(ax, 8, 4)