Outflow

Roland Haas <roland.haas@physcis.gatech.edu>

August 15 2009

Abstract

Outflow calculates the flow of rest mass density across a SphericalSurface, eg. and apparent horizon or a sphere at “infinity”.

1 Introduction

Hydrodynamic simulations should conserve eg. the total rest mass in the system (outside he event horizon, that is). This thorn allow to measure the flux of rest mass across a given SphericalSurface.

2 Physical System

The Valencia formalism defines its $D = \rho W$ variable such that one can define a total rest mass density as: \begin {displaymath} M = \int \sqrt {\gamma } D d^3x \end {displaymath} and from the EOM of $D$ (see [1]) \begin {displaymath} \frac {\partial }{\partial x^0} (\sqrt {\gamma }D )+ \frac {\partial }{\partial x^i} (\sqrt {\gamma } \alpha D (v^i- \beta ^i/\alpha )) = 0 \end {displaymath} one obtains:

\begin {eqnarray} \dot M &=& - \int _V \frac {\partial }{\partial x^i} (\sqrt \gamma \alpha D(v^i - \beta ^i/\alpha )) d^3x \\ &=& -\int _{\partial V} \sqrt \gamma \alpha D (v^i-\beta ^i/\alpha ) d\sigma _i \end {eqnarray}

where $\sigma _i$ is the ordinary flat space directed surface element of the enclosing surface, eg. \begin {displaymath} d\sigma _i = \hat r_i r^2 \sin \theta d\theta d\phi \end {displaymath} with $\hat r_i = [\cos \phi \sin \theta , \sin \phi \sin \theta , \cos (\theta )]$ for a sphere of radius r.

For a generic SphericalSurface parametrized by $\theta $ and $\phi $ one has:

\begin {eqnarray} x &=& \bar r(\theta , \phi ) \cos \phi \sin \theta \\ y &=& \bar r(\theta , \phi ) \sin \phi \sin \theta \\ z &=& \bar r(\theta , \phi ) \cos \theta \end {eqnarray}

where $\bar r$ is the isotropic radius. Consequently the surface element is

\begin {eqnarray} d \sigma _i &=& \left ( \frac {\partial \vec {\bar r}}{\partial \theta } \times \frac {\partial \vec {\bar r}}{\partial \phi } \right )_i \\ &=& \bar r^2 \sin \theta \hat {\bar r}_i - \frac {\partial \bar r}{\partial \theta } \bar r \sin \theta \hat \theta _i - \frac {\partial \bar r}{\partial \phi } \bar r \hat \phi _i \end {eqnarray}

where $\hat {\bar r}$, $\hat \theta $ and $\hat \phi $ are the flat space standard unit vectors on the sphere [3].

3 Numerical Implementation

We implement the surface integral by interpolating the required quantities ($g_{ij}$, $\rho $, $v^i$, $\beta ^i$, $\alpha $) onto the spherical surface and then integrate using a fourth order convergent Newton-Cotes formula.

For the $\theta $ direction SphericalSurfaces defines grid points such that \begin {displaymath} \theta _i = -(n_\theta - 1/2) \Delta _\theta + i \Delta _\theta \qquad 0 \le i < N_\theta -1 \end {displaymath} where $N_\theta $ is the total number of intervals (sf_ntheta), $n_\theta $ is the number of ghost zones in the $\theta $ direction (nghosttheta) and $\Delta _\theta = \frac {\pi }{N_\theta -2n\theta }$. Since with this \begin {displaymath} \theta _i \in [\Delta _\theta /2, \pi - \Delta _\theta /2] \qquad \mbox {for $i \in [n_\theta , N_\theta -n_\theta -1]$} \end {displaymath} we do not have grid points at the end of the interval $\{0,\pi \}$ we derive an open extended Newton-Cotes formula from Eq. 4.1.14 of [2] and a third order accurate extrapolative rule (see Maple worksheet). We find

\begin {eqnarray} \int _{x_0}^{x_{N-1}} f(x) dx &\approx & h \Bigl \{ \frac {13}{12} f_{1/2} + \frac {7}{8} f_{3/2} + \frac {25}{24} f_{5/2} + f_{7/2} + f_{9/2} + \cdots + f_{N-1-7/2} \nonumber \\ && + f_{N-1-9/2} + \frac {25}{24} f_{N-1-5/2} + \frac {7}{8} f_{N-1-3/2} + \frac {13}{12} f_{N-1-1/2} \Bigr \} + O(1/N^4) \end {eqnarray}

For the $\phi $ direction SphericalSurfaces defines grid points such that \begin {displaymath} \phi _i = -n_\phi \Delta _\phi + i \Delta _\phi \qquad 0 \le i < N_\phi -1 \end {displaymath} where $N_\phi $ is the total number of intervals (sf_nphi), $n_\phi $ is the number of ghost zones in the $\phi $ direction (nghostphi) and $\Delta _\phi = \frac {\pi }{N_\phi -2n\phi }$. With this \begin {displaymath} \phi _i \in [0, 2 \pi - \Delta _\phi ] \qquad \mbox {for $i \in [n_\phi , N_\phi -n_\phi -1]$} \end {displaymath} we use a simple extended trapezoid rule to achieve spectral convergence due to the periodic nature of $\phi $ (note: $x_N = x_0$)

\begin {eqnarray} \int _{x_0}^{x_N} f(x) dx \approx & h \sum _{i=0}^{N-1} f_i \end {eqnarray}

The derivatives of $\vec {\bar r}$ along $\theta $ and $\phi $ are obtained numerically and require at least two ghost zones in $\theta $ and $\phi $.

4 Using This Thorn

Right now surface can only be prescribed by SphericalSurfaces, the flux through each surface is output in in a file outflow_det_%d.asc

4.1 Interaction With Other Thorns

Takes care to schedule itself after SphericalSurfaces_HasBeenSet.

4.2 Examples

See the parameter file in the test directory. For spherical symmetric infall the flux through all detectors should be equal (since rest mass must be conserved).

References

[1] José A. Font, “Numerical Hydrodynamics and Magnetohydrodynamics in General Relativity”, Living Rev. Relativity 11, (2008), 7. URL (cited on August 15. 2009): http://www.livingreviews.org/lrr-2008-7

[2] Press, W. H. (2002). Numerical recipes in C++: the art of scientific computing. Array Cambridge, UK: Cambridge University Press.

[3] Charles W. Misner, Kip S. Thorne, and John Archibald Wheeler (1973). “Gravitation”. New York: W.H. Freeman and Company.

5 Parameters


compute_every	Scope: private	INT

Description: How frequently to compute the outflow

Range		Default: 1
1:*	number of iterations


compute_every_det	Scope: private	INT

Description: How frequently to compute the outflow

Range		Default: -1
-1	take from compute_every
	do not compute
1:*	number of iterations


coord_system	Scope: private	STRING

Description: What is the coord system?

Range		Default: cart3d
.*	Any smart string will do


extra_variables	Scope: private	STRING

Description: extra (scalar) variables to project onto the spherical surface

Range		Default: (none)
.*	Any Cactus variables


interpolator_name	Scope: private	STRING

Description: Which interpolator should I use

Range		Default: Hermite polynomial interpolation
.+	Any nonempty string


interpolator_pars	Scope: private	STRING

Description: Parameters for the interpolator

Range		Default: order=3
.*	”Any string that Util_TableSetFromStr ing() will take”


num_detectors	Scope: private	INT

Description: how many detectors do we have

Range		Default: (none)
0:*	number of detectors


num_thresholds	Scope: private	INT

Description: how many thresholds to use

Range		Default: (none)
	ignore thresholds, output only one flux
1:20	how many to use


out_format	Scope: private	STRING

Description: Which format for Scalar floating-point number output

Range		Default: .19g
see [1] below	output with given precision in exponential / floating point notation

[1]

\^({\textbackslash}.[1]?[0-9])?[EGefg]\$


output_2d_data	Scope: private	BOOLEAN

Description: output 2d data on spherical surface

		Default: no


output_relative_coordinates	Scope: private	BOOLEAN

Description: output coordinates on spheroid relative to sf_centroid, absolute coordinates otherwise

		Default: no


override_radius	Scope: private	BOOLEAN

Description: do not take radius values from spherical surface

		Default: no


rad_rescale	Scope: private	REAL

Description: Factor for scaling radius for this detector surface (e.g. for buffer around AH surface)

Range		Default: 1
(0:*	Anything positive


radius	Scope: private	REAL

Description: spherical radius for this detector surface

Range		Default: -1
-1	Invalid
(0:*	any positive value


surface_index	Scope: private	INT

Description: the indices of the sphericalsurfaces from which we take the surfaces to compute the outflow on

Range		Default: -1
-1	Invalid
0:*	valid surface


surface_name	Scope: private	STRING

Description: the indices of the sphericalsurfaces from which we take the surfaces to compute the outflow on

Range		Default: (none)
	use surface_index
.*	any string


threshold	Scope: private	REAL

Description: compute fluxes with Lorentz factors / specific energy at infinity above these levels only

Range		Default: (none)
:	Any real number


threshold_on_var	Scope: private	KEYWORD

Description: Which variable to use as threshold

Range		Default: w_lorentz
eninf	(Approximate) specific energy at infinity, - u_t - 1
w_lorentz	Lorentz factor


verbose	Scope: private	INT

Description: How verbose do you want it?

Range		Default: 1
0:*	0 is nothing, 1 is interesting, >1 is crazy


maxnphi	Scope: shared from SPHERICALSURFACE	INT


maxntheta	Scope: shared from SPHERICALSURFACE	INT


nghostsphi	Scope: shared from SPHERICALSURFACE	INT


nghoststheta	Scope: shared from SPHERICALSURFACE	INT


sphericalsurfaces_nsurfaces	Scope: shared from SPHERICALSURFACE	INT

6 Interfaces

General

Implements:

outflow

Inherits:

admbase

hydrobase

sphericalsurface

Grid Variables

6.0.1 PRIVATE GROUPS


Group Names	Variable Names	Details

outflow_flux	outflow_flux	compact	0
		description	flux of mass through the detectors
		dimensions	0
		distribution	CONSTANT
		group type	SCALAR
		tags	checkpoint=”no”
		timelevels	1
		vararray_size	num_detectors
		variable type	REAL

fluxdens_projected	fluxdens_projected	compact	0
		description	2D (theta
		description	phi) grid arrays for flux density
		dimensions	2
		distribution	CONSTANT
		group type	ARRAY
		size	SPHERICALSURFACE::MAXNTHETA
		size	SPHERICALSURFACE::MAXNPHI
		tags	checkpoint=”no”
		timelevels	1
		vararray_size	num_detectors
		variable type	REAL

w_lorentz_projected	w_lorentz_projected	compact	0
		description	2D (theta
		description	phi) grid arrays for Lorentz factor
		dimensions	2
		distribution	CONSTANT
		group type	ARRAY
		size	SPHERICALSURFACE::MAXNTHETA
		size	SPHERICALSURFACE::MAXNPHI
		tags	checkpoint=”no”
		timelevels	1
		vararray_size	num_detectors
		variable type	REAL

eninf_projected	eninf_projected	compact	0
		description	2D (theta
		description	phi) grid arrays for the specific energy at infinity
		dimensions	2
		distribution	CONSTANT
		group type	ARRAY
		size	SPHERICALSURFACE::MAXNTHETA
		size	SPHERICALSURFACE::MAXNPHI
		tags	checkpoint=”no”
		timelevels	1
		vararray_size	num_detectors
		variable type	REAL

surface_projections		compact	0
	surface_projection_0	description	2D (theta
		description	phi) grid arrays for points on the spherical surfaces
	surface_projection_1	dimensions	2
	surface_projection_2	distribution	CONSTANT
	surface_projection_3	group type	ARRAY
	surface_projection_4	size	SPHERICALSURFACE::MAXNTHETA
		size	SPHERICALSURFACE::MAXNPHI
	surface_projection_5	tags	checkpoint=”no”
	surface_projection_6	timelevels	1
	surface_projection_7	vararray_size	num_detectors
	surface_projection_8	variable type	REAL

7 Schedule

This section lists all the variables which are assigned storage by thorn EinsteinAnalysis/Outflow. Storage can either last for the duration of the run (Always means that if this thorn is activated storage will be assigned, Conditional means that if this thorn is activated storage will be assigned for the duration of the run if some condition is met), or can be turned on for the duration of a schedule function.

Storage

Always:
outflow_flux
fluxdens_projected
w_lorentz_projected
eninf_projected
surface_projections

Scheduled Functions

CCTK_ANALYSIS

outflow

compute outflow

	Language:	c
	Options:	global
	Type:	function

CCTK_STARTUP

outflow_setup

set up global ompute data structures

	Language:	c
	Options:	global
	Type:	function