1 Introduction

This thorn provides routines to implement radiative boundary conditions as described in section VI of [2]: \begin {equation} f = f_0 + u(r - v t)/r\mbox {,} \label {eqn:SpacetimeX_NewRadX_bc} \end {equation} in its differential form \begin {equation} \partial _t f + v \partial _r f - v\,(f - f_0) = 0\mbox {,} \label {eqn:SpacetimeX_NewRadX_diff_bc} \end {equation} where \(v\) is the asymptotic propagation speed of the field \(f\) and \(f_0\) is its asymptotic value.

2 Physical System

Quoting [2]:

At the outer boundary we use a radiation (Sommerfeld) boundary condition. We start from the assumption that near the boundary all fields behave as spherical waves, namely we impose the condition \begin {equation} f = f_0 + u(r-vt)/r . \label {eq:radiation} \end {equation} Where \(f_0\) is the asymptotic value of a given dynamical variable (typically 1 for the lapse and diagonal metric components, and zero for everything else), and \(v\) is some wave speed. If our boundary is sufficiently far away one can safely assume that the speed of light is 1, so \(v=1\) for most fields. However, the gauge variables can easily propagate with a different speed implying a different value of \(v\).

In practice, we do not use the boundary condition (??) as it stands, but rather we use it in differential form: \begin {equation} \partial _t f + v \partial _r f - v \, (f - f_0)/r = 0 . \end {equation} Since our code is written in Cartesian coordinates, we transform the last condition to \begin {equation} \frac {x_i}{r} \, \partial _t f + v \partial _i f + \frac {v x_i}{r^2} \, \left ( f - f_0 \right ) = 0 . \end {equation} We finite difference this condition consistently to second order in both space and time and apply it to all dynamic variables (with possible different values of \(f_0\) and \(v\)) at all boundaries.

There is a final subtlety in our boundary treatment. Wave propagation is not the only reason why fields evolve near a boundary. Simple infall of the coordinate observers will cause some small evolution as well, and such evolution is poorly modeled by a propagating wave. This is particularly important at early times, when the radiative boundary condition introduces a bad transient effect. In order to minimize the error at our boundaries introduced by such non-wavelike evolution, we allow for boundary behavior of the form: \begin {equation} f = f_0 + u(r-vt)/r + h(t)/r^n , \label {eq:radpower} \end {equation} with \(h\) a function of \(t\) alone and \(n\) some unknown power. This leads to the differential equation

or in Cartesian coordinates \begin {equation} \frac {x_i}{r} \, \partial _t f + v \partial _i f + \frac {v x_i}{r^2} \, \left ( f - f_0 \right ) \simeq \frac {x_i h'(t)}{r^{n+1}} . \label {eqn:SpacetimeX_NewRadX_num_bc} \end {equation}

This expression still contains the unknown function \(h'(t)\). Having chosen a value of \(n\), one can evaluate the above expression one point away from the boundary to solve for \(h'(t)\), and then use this value at the boundary itself. Empirically, we have found that taking \(n=3\) almost completely eliminates the bad transient caused by the radiative boundary condition on its own.

3 Numerical Implementation

The thorn implements the boundary condition  (??) for the evolution. The original Carpet-based version of NewRad provides an additional routine ExtrapolateGammas, which can be used to extrapolate a grid function from the interior to the outer boundaries at \(t=0\). Originally, this was used in the BSSN evolution thorn ML_BSSN on the contracted Christoffel symbols \(g^{jk} \Gamma ^i_{jk}\) (Xt1…Xt3 in the code). This extrapolation routine is not yet ported over to the current version of NewRadX.

3.1 Radiative boundary condition

The derivatives \(\partial _i f\) appearing in (??) are approximated as 2^nd order accurate finite difference expressions using centred expressions \begin {equation} f'(x) \approx \frac {f(x-\Delta x) - f(x+\Delta x)}{2\Delta x} \end {equation} and one sided expressions \begin {equation} f'(x) \approx \pm \frac {3 f(x) - 4 f(x\mp \Delta x) + f(x\mp 2\Delta x)}{2\Delta x}\mbox {,} \end {equation} for boundaries around the \(\pm \) coordinate direction.

Finally we try to extrapolate for the part of the boundary that does not behave as a pure wave (i.e. Coulomb type terms caused by infall of the coordinate lines).

This we do by comparing the source term nghostzones grid points away from the boundary (which we already have), to what we would have obtained if we had used the boundary condition there. The difference gives us an idea of the missing part and we extrapolate that to the boundary assuming a power-law decay.

4 Using This Thorn

The thorn provides the routine NewRadX_Apply, which should be scheduled in the ODESolvers_RHS bin to apply radiative boundary conditions on evolved variables. To use it in an evolution thorn, it can either be called inside the RHS evaluation routine, or scheduled in a separate function right after the RHS routine.

The call signature is:

void NewRadX_Apply                               \
        (const cGH *restrict const cctkGH,       \
         const Loop::GF3D2<const CCTK_REAL> var, \
         Loop::GF3D2<CCTK_REAL> rhs,             \
         CCTK_REAL var0,                         \
         CCTK_REAL v0,                           \
         CCTK_REAL radpower)

where var0 corresponds to \(f_0\), v0 is \(v\) and radpower is the assumed decay exponent \(n\).

4.1 Obtaining This Thorn

The thorn is part of the SpacetimeX arrangement available on GitHub.

git clone git@github.com:lwJi/SpacetimeX.git

4.2 Basic Usage

For an overview of the implementation and usage of this thorn, please see [3] for a presentation by Alexandru Dima and Cheng-Hsin Cheng during the weekly CarpetX developer meetings.

4.3 Special Behaviour

In contrast with the Carpet-version of NewRad, which updates the RHS variables on the physical outer boundary of the simulation domain, NewRadX updates them on the so-called “outermost interior points”. These points are defined as the interior points which are located nghostzones or fewer grid points from the outer boundary. Therefore, the boundary conditions are effectively applied on a set of points that are labelled as “interior” in CarpetX. It is necessary to introduce this change because imposing a wave-like time evolution on the CarpetX “proper” boundary points is incompatible with the synchronization conducted by AMReX.

Currently, this thorn only supports single-patch, Cartesian grids provided by CarpetX, but not multi-patch coordinate systems.

4.4 Interaction With Other Thorns

This thorn requires thorn Loop which is part of CarpetX.

4.5 Examples

The simplest example is likely to inspect the source code for TestNewRadX, which is a thorn from SpacetimeX that tests the boundary condition on the evolution of a scalar pulse on flat background. Here the boundary conditions are applied in the RHS calculation routine TestNewRadX_RHS, where NewRadX_Apply is called after the main loop that updates the RHS of the state variables in the interior.

<configuration.ccl>
REQUIRES Loop NewRadX

<interface.ccl>
USES INCLUDE HEADER: newradx.hxx

<schedule.ccl>
SCHEDULE TestNewRadX_RHS IN ODESolvers_RHS
{
  LANG: C
  READS: state(everywhere)
  WRITES: rhs(interior)
  SYNC: rhs
} "Calculate scalar wave RHS"

<test_newradx.cxx>
#include <newradx.hxx>

extern "C" void TestNewRadX_RHS(CCTK_ARGUMENTS) {
  DECLARE_CCTK_ARGUMENTSX_TestNewRadX_RHS;
  DECLARE_CCTK_PARAMETERS;

  /* First loop over the interior to update RHS of the state variables */
  grid.loop_int_device<0, 0, 0>(grid.nghostzones,
                                [=] CCTK_DEVICE(const Loop::PointDesc &p)
                                    CCTK_ATTRIBUTE_ALWAYS_INLINE {
                                      uu_rhs(p.I) = vv(p.I);
                                      vv_rhs(p.I) = Laplacian_4thorder(uu, p);
                                    });

  /* After RHS grid functions are updated, apply boundary conditions */
  if (test_use_newradx) {
    NewRadX::NewRadX_Apply(cctkGH, uu, uu_rhs, 0, 1, n_falloff);
    NewRadX::NewRadX_Apply(cctkGH, vv, vv_rhs, 0, 1, n_falloff + 1);
  }

}

Other examples of CarpetX evolution thorns that use NewRadX include CanudaX_BSSNMoL[4] and Z4c[5].

4.6 Support and Feedback

Please use the Einstein Toolkit mailing list users@einsteintoolkit.org to report issues and ask for help.

5 History

5.1 Thorn Source Code

• Initial version of method likely by Miguel Alcubierre;

• Carpet version of code by Erik Schnetter;

• CarpetX version of code by Cheng-Hsin Cheng, Jake Doherty, Michail Chabanov, Alexandru Dima.

5.2 Thorn Documentation

• Initial version by Roland Haas using text by Miguel, Erik.

• Current version by Alexandru Dima and Cheng-Hsin Cheng.

5.3 Acknowledgements

This documentation is based on the original one for the Carpet version of NewRad, which copies text from comments in the source code as well as the paper [2].

References

6 Parameters

7 Interfaces

General

Implements:

newradx

Inherits:

carpetx

Adds header:

newradx.hxx

Uses header:

loop_device.hxx

8 Schedule

This section lists all the variables which are assigned storage by thorn SpacetimeX/NewRadX. 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

NONE