[Users] Implementing a simple 1+1 BSSN in Python

Tue May 17 16:53:35 CDT 2022

Chris

That's good to hear. (I assumed, same as you, that the Einstein Toolkit
would evolve the boundary points as well in this case...)

-erik

On Tue, May 17, 2022 at 5:01 AM Chris Stevens <
chris.stevens at canterbury.ac.nz> wrote:

> Hi Ian, Erik,
>
> Thanks for your feedback. I've finally tracked the problem to what I
> *thought* the ET was doing at the boundary vs what was *actually
> happening*. My par file has
>
> Coordinates::outer_boundary_size = 1
>
> and together with the none BC, the initial values of the system variables
> remain as such at these points. This is different to my Python code, which
> is evolving the outer boundary points (in preparation for the SAT method).
> Thus, I set the time derivative of all fields to zero at the outer boundary
> points and switched back to RK4 (Euler is indeed unstable!) and everything
> works now. So, the issue can be summarised as my misunderstanding of what
> is happening in the ET.
>
> Thanks again for your helpful comments!
>
> Chris
>
>
>
>
>
>
>
> *Dr Chris Stevens*
>
> *Lecturer in Applied Mathematics*
>
> Rm 602, Jack Erskine building
>
> School of Mathematics and Statistics
>
> T: +64 3 369 0396 (Internal 90396)
>
> University of Canterbury | Te Whare Wānanga o Waitaha
>
> Private Bag 4800, Christchurch 8140, New Zealand
>
> http://www.chrisdoesmaths.com
>
> *Director*
> SCRI Ltd
> http://www.scri.co.nz
>
> ------------------------------
> *From:* Ian Hinder <ian.hinder at manchester.ac.uk>
> *Sent:* 17 May 2022 19:32
> *To:* Chris Stevens <chris.stevens at canterbury.ac.nz>
> *Cc:* Erik Schnetter <schnetter at gmail.com>; Einstein Toolkit Users <
> users at einsteintoolkit.org>
> *Subject:* Re: [Users] Implementing a simple 1+1 BSSN in Python
>
> You don't often get email from ian.hinder at manchester.ac.uk. Learn why
> this is important <https://aka.ms/LearnAboutSenderIdentification>
> Hi Chris,
>
> - Note that the Euler method you are using is unstable, and will probably
> amplify any roundoff errors you have.
> - Does everything work if you use periodic boundary conditions?
>
> Erik said:
>
> Given your description, I assume that all derivatives would be identically
> zero, that there would not even be a floating-point round-off error. Is
> that so? If not, why not?
>
>
> I think this would be true for something like (u_{j+1} - u_j) / dx, but it
> might not be true for ( u_{j+2} - 2u_{j+1} + u_j ) / dx^2 which would
> evaluate to ( u - 2u + u ) / dx^2, and depending on how the compiler
> implements it, (u - 2 u) for example might not be exactly representable in
> floating point, and might not cancel exactly the (+u). If u were zero, I
> would agree.
>
> On 14 May 2022, at 04:22, Chris Stevens <chris.stevens at canterbury.ac.nz>
> wrote:
>
> Hi Erik,
>
> Thanks for your quick response!
>
> I use the none bc type for now, rather than the more complicated radiative
> bc. I use a stencil that leans over as you approach the boundary. Currently
> this is a 4th order interior 3rd order sbp fd operator of Strand.
>
> By looking at the output on the boundary in both codes, there are round
> off level numbers appearing after the first Euler step in multiple
> variables. There is no use comparing these exactly as they're so small but
> note that the same variables have these in both codes. These come from
> taking spatial derivatives.
>
> Sorry, the lapse is initially one.
>
> I have compared output by saving the data in both codes. It seems like,
> since the first few time steps everything is round off errors, it is hard
> to see exactly where things start to go wrong.
>
> I will try your small grid and see if that yields anything more obvious.
>
> Thanks,
>
> Chris
>
>
>
>
> ------------------------------
> *From:* Erik Schnetter <schnetter at gmail.com>
> *Sent:* Saturday, May 14, 2022 3:09:17 PM
> *To:* Chris Stevens <chris.stevens at canterbury.ac.nz>
> *Cc:* Einstein Toolkit Users <users at einsteintoolkit.org>
> *Subject:* Re: [Users] Implementing a simple 1+1 BSSN in Python
>
> Chris
>
> What do you mean by "Boundary conditions are not applied"? What do you do
> at the boundaries? Are you using one-sided differencing there?
>
> Given your description, I assume that all derivatives would be identically
> zero, that there would not even be a floating-point round-off error. Is
> that so? If not, why not?
>
> Are you really setting the lapse to zero initially? It should be one for
> Minkowski, same as e.g. g_xx.
>
> What happens if you choose a very small grid size (e.g. just 10 points)
> and compare the initial conditions between the Einstein Toolkit and the
> Python code? Are they identical? You would need to compare ASCII output and
> compare all the digits.
>
> What happens if you take a single Euler time step? Are the values still
> all identical? Try to find out and which iteration the first grid point
> differs between the two codes.
>
> -erik
>
>
>
> On Fri, May 13, 2022 at 9:59 PM Chris Stevens <
> chris.stevens at canterbury.ac.nz> wrote:
>
> Hi all,
>
> I've been having trouble implementing the BSSN equations, as given by
> CTGamma, in Python using the package COFFEE, which I am very familiar with:
>
> https://gitlab.com/thebarista/coffee
> https://doi.org/10.1016/j.softx.2019.100283
>
> The main comments are:
>
> -- For simplicity, I am starting with a 1+1 code (set 2/3 spatial
> derivatives to be exactly zero).
> -- The 1+log slicing for the lapse and a zero shift.
> -- The "phi" conformal factor type
> -- The evolution equations are written to Python by a very slight
> modification of the .m file in CTGEvolution. I've confirmed this .m file
> outputs the same file that is in the src directory.
> -- The vanishing of the trace of A is enforced after each full timestep in
> the same way as CTGamma.
> -- Boundary conditions are not applied.
> -- Fourth order SBP operators that lean over toward the boundaries are
> used.
> -- Minkowski initial data with all fields zero except \gamma11/22/33.
> -- Euler step is used for time integration to easier compare between
> Python and the Einstein toolkit.
>
> This leads to a stable evolution in the Einstein toolkit with the attached
> .par file using a simple one patch Cartesian grid. However, in Python, I
> get that fields diverge, starting at the boundaries and then propagating
> inward (including with RK4). This behaviour does not go away with higher
> resolution or reduced CFL. It is however, fixed completely by setting the
> second spatial derivative of \alpha to be exactly zero. Trying many
> different SBP FD operators that have proved successful for other projects
> has not changed the behaviour.
>
> The COFFEE code has been used and tested in a variety of projects and so
> there should not be a problem here. Thus, the only thing I can think of, is
> that the Einstein toolkit is doing something that I am not realizing, that
> is helping it remain stable. Hopefully, somebody on this mailing list might
> have an idea?
>
> Thanks in advance,
>
> Chris
>
> <Outlook-zot52qe4.png>
>
>
>
> <Outlook-sdnm0er0.png>
>
>
> *Dr Chris Stevens*
> *Lecturer in Applied Mathematics*
> Rm 602, Jack Erskine building
> School of Mathematics and Statistics
> T: +64 3 369 0396 (Internal 90396)
> University of Canterbury | Te Whare Wānanga o Waitaha
> Private Bag 4800, Christchurch 8140, New Zealand
> http://www.chrisdoesmaths.com
>
> *Director*
> SCRI Ltd
> http://www.scri.co.nz
>
> _______________________________________________
> Users mailing list
> Users at einsteintoolkit.org
> http://lists.einsteintoolkit.org/mailman/listinfo/users
>
>
>
> --
> Erik Schnetter <schnetter at gmail.com>
> http://www.perimeterinstitute.ca/personal/eschnetter/
>
> _______________________________________________
> Users mailing list
> Users at einsteintoolkit.org
> http://lists.einsteintoolkit.org/mailman/listinfo/users
>
>
> --
> Ian Hinder
> Research Software Engineer
> University of Manchester, UK
>
>

-- 
Erik Schnetter <schnetter at gmail.com>
http://www.perimeterinstitute.ca/personal/eschnetter/
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://lists.einsteintoolkit.org/pipermail/users/attachments/20220517/e6b96b72/attachment-0001.html 
-------------- next part --------------
A non-text attachment was scrubbed...
Name: Outlook-ldmmb54u.png
Type: image/png
Size: 13436 bytes
Desc: not available
Url : http://lists.einsteintoolkit.org/pipermail/users/attachments/20220517/e6b96b72/attachment-0002.png 
-------------- next part --------------
A non-text attachment was scrubbed...
Name: Outlook-kg41ljz3.png
Type: image/png
Size: 17337 bytes
Desc: not available
Url : http://lists.einsteintoolkit.org/pipermail/users/attachments/20220517/e6b96b72/attachment-0003.png 