[Users] Grid structure inconsistency & errors in surface mask locations

Sat Apr 17 20:32:00 CDT 2021

On Sat, Apr 17, 2021 at 8:32 PM Konrad Topolski
<k.topolski2 at student.uw.edu.pl> wrote:
>
> Changing the shiftout did not help unfortunately.
>
> Here's the link to the output, errors, parameter file and my running script:
> Let me know if that's a preferable format.
>
> https://gist.github.com/konrad-topolski/57af20b18539cd7bcb4ae2d11993144f
>
> If it helps troubleshooting in any way, when I set up fully nonsymmetric version (i.e. zmin=-120) (no reflection symmetry in z), the error persists but now the last coordinate is affected, as in:
>
>    CCTK_InterpLocalUniform():
>         interpolation point is either outside the grid,
>         or inside but too close to the grid boundary!
>         (this may be caused by a global interpolation with
>          driver::ghost_size too small)
>         0-origin interpolation point number pt=2963 of N_interp_points=2964
>         interpolation point (x,y,z)=(2.63796e+243,1.58133e+243,3.77779e+243)
>         grid x_min(delta_x)x_max = -122.5(1.25)122.5
>         grid y_min(delta_y)y_max = -122.5(1.25)122.5
>         grid z_min(delta_z)z_max = -122.5(1.25)3.75

Ah, these coordinates come from the interpolator! I thought this was
the grid structure description output by Carpet. That's fine, then;
the local interpolator sees only part of the grid at a time.

In the output I see this:

INFO (AHFinderDirect): AH 1/3: r=0.105016 at (2.999966,-0.002037,0.000000)
INFO (AHFinderDirect): AH 1/3: area=3.141592763 m_irreducible=0.2500000044

which means that the first AH is found fine; it has the right mass.
The coordinate radius is 0.1, which is a pointer as to what grid
resolution you need there. For starters, I'd say the resolution there
needs to be 0.1/15 (or 0.1/20); you can then fine-tune this later
depending on how things converge.

Information about the second horizon is not visible since it's
calculated on the second MPI process. You can either use the "-roe"
command-line option to Cactus to see that output, or look for
"BH_diagnostics" files. If you have 2 such files, then both horizons
are found, and the horizon properties are described there.

When the interpolator reports an error, it is called from the function
"QuasiLocalMeasures::qlm_calculate". This thorn examines horizons
after they are found, and e.g. calculates their spin. (It also
calculates many other quantities, but most people are only interested
in the spin.) It shouldn't be called when a horizon is not found.
Maybe that's the error – a miscommunication between AHFinderDirect,
SphericalSurface, and QuasiLocalMeasures. I'll check the parameter
file, brb...

... am back.

There are 3 surfaces (good); one for each horizon.
There are 3 horizons; each is stored into the correct surface.
(Horizons are numbered starting from 1, surfaces starting from 0.).
Also correct.
QuasiLocalMeasures examines 3 surfaces. Also correct.

This looks wrong:
INFO (SphericalSurface): SphericalSurface[0]: finest refinement-level
that fully contains this surface:  sf_maxreflevel = 1920103026
I don't know where this comes from. Maybe someone else can explain?

This line
AHF BH_diagnostics::save[2] found_flag=0
indicates that the second horizon wasn't found. I would start
debugging there. You set the AH finder to be quite verbose (good for
debugging!), and output such as
INFO (AHFinderDirect):    proc 1/horizon 2:it 19 r_grid=0.634 ||Theta||=2.1e+00
indicates that the AH finder is currently examining a surface at
r=0.634 (this looks very large; the horizon should be much smaller),
and the expansion of that surface is 2.1 (that's very large; a horizon
has expansion 0). You need to find a better initial guess.

I looked at the AH finding iterations:
$ grep 'proc 1/horizon 2' output.out
INFO (AHFinderDirect):    proc 1/horizon 2:it 1 r_grid=0.125 ||Theta||=1.1e+01
INFO (AHFinderDirect):    proc 1/horizon 2:it 2 r_grid=0.137 ||Theta||=9.9e+00
INFO (AHFinderDirect):    proc 1/horizon 2:it 3 r_grid=0.151 ||Theta||=9.0e+00
INFO (AHFinderDirect):    proc 1/horizon 2:it 4 r_grid=0.166 ||Theta||=8.2e+00
INFO (AHFinderDirect):    proc 1/horizon 2:it 5 r_grid=0.182 ||Theta||=7.5e+00
INFO (AHFinderDirect):    proc 1/horizon 2:it 6 r_grid=0.199 ||Theta||=6.8e+00
INFO (AHFinderDirect):    proc 1/horizon 2:it 7 r_grid=0.219 ||Theta||=6.2e+00
INFO (AHFinderDirect):    proc 1/horizon 2:it 8 r_grid=0.240 ||Theta||=5.7e+00
INFO (AHFinderDirect):    proc 1/horizon 2:it 9 r_grid=0.263 ||Theta||=5.2e+00
INFO (AHFinderDirect):    proc 1/horizon 2:it 10 r_grid=0.288 ||Theta||=4.7e+00
INFO (AHFinderDirect):    proc 1/horizon 2:it 11 r_grid=0.315 ||Theta||=4.3e+00
INFO (AHFinderDirect):    proc 1/horizon 2:it 12 r_grid=0.344 ||Theta||=3.9e+00
INFO (AHFinderDirect):    proc 1/horizon 2:it 13 r_grid=0.377 ||Theta||=3.6e+00
INFO (AHFinderDirect):    proc 1/horizon 2:it 14 r_grid=0.412 ||Theta||=3.3e+00
INFO (AHFinderDirect):    proc 1/horizon 2:it 15 r_grid=0.450 ||Theta||=3.0e+00
INFO (AHFinderDirect):    proc 1/horizon 2:it 16 r_grid=0.491 ||Theta||=2.8e+00
INFO (AHFinderDirect):    proc 1/horizon 2:it 17 r_grid=0.535 ||Theta||=2.5e+00
INFO (AHFinderDirect):    proc 1/horizon 2:it 18 r_grid=0.583 ||Theta||=2.3e+00
INFO (AHFinderDirect):    proc 1/horizon 2:it 19 r_grid=0.634 ||Theta||=2.1e+00
INFO (AHFinderDirect):    proc 1/horizon 2:it 20 r_grid=0.689 ||Theta||=2.0e+00

The initial guess is probably too small (maybe r=0.3 is better?), and
the AH finder never latches onto a horizon surface. Mind that a
surface can have an arbitrary shape, and the output shows only the
average radius. If the surface becomes distorted, but the actual
horizon is spherical, then the AH finder won't find it.

If you find it difficult to find the horizon initially (this can be a
tedious process), then you can set up a run and look for 50 horizons,
all with different initial radii, or with different initial shapes
(see "coordinate ellipsoid"). Look also at
<https://einsteintoolkit.org/thornguide/EinsteinAnalysis/AHFinderDirect/documentation.html>
for the gory details.

Good luck hunting for the horizon!

-erik

-- 
Erik Schnetter <schnetter at cct.lsu.edu>
http://www.perimeterinstitute.ca/personal/eschnetter/