Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Rework handling of roundoff domain extent #3199

Merged
merged 9 commits into from
Apr 6, 2023

Conversation

dpgrote
Copy link
Contributor

@dpgrote dpgrote commented Mar 14, 2023

Summary

This reworks previous changes that created the round off domain which was done to avoid round off errors when locating particles in the domain, specifically when handling periodic boundaries. This PR changes how this is done so that the roundoff domain makes as small a change as possible.

Additional background

Previous PRs #2839 and #2950 created a round off domain that has a slightly smaller extent that the actual domain. Any particles that ended up outside the round off domain but inside the actual domain were shifted to be on the boundary of the round off domain. The main purpose was to handle problems with periodic boundary with mixed precision (single precision particles and double precision grid). In that case, it is possible for a particle to have z < lo and z+(hi-lo) > hi. It is also possible for a particle to have z < hi but (z-lo)*dzinv > ihi putting it in a grid cell beyond the domain. Unfortunately, the roundoff domains was being set further inside that was needed and would affect valid particles that would not have those issues. While the problem is small, the effects can be noticeable when examining simulation properties, for example charge conservation in WarpX.

This PR changes the roundoff domain so that only problematic particles are affected. The code now uses std::nextafter to find the exact place where the problems occurs. At the lower end, it ensures that the casted position will be greater than ProbLo. At the upper end, it ensures that (z-lo)*dzinv <= ihi.

One comment is that this assumes that particle grid coordinates will always be calculated at the higher precision, with expressions like iz = (int)(z - lo), without lo being cast to ParticleReal.

Checklist

The proposed changes:

  • fix a bug or incorrect behavior in AMReX
  • add new capabilities to AMReX
  • changes answers in the test suite to more than roundoff level
  • are likely to significantly affect the results of downstream AMReX users
  • include documentation in the code and/or rst files, if appropriate

@atmyers atmyers self-requested a review March 14, 2023 22:07
@atmyers atmyers self-assigned this Mar 14, 2023
@dpgrote dpgrote changed the title Rework handling of round off errors with periodic boundaries [WIP]Rework handling of round off errors with periodic boundaries Mar 16, 2023
@dpgrote dpgrote changed the title [WIP]Rework handling of round off errors with periodic boundaries Rework handling of round off errors with periodic boundaries Mar 21, 2023
@dpgrote dpgrote changed the title Rework handling of round off errors with periodic boundaries Rework handling of roundoff domain extent Mar 21, 2023
@atmyers atmyers enabled auto-merge (squash) April 6, 2023 00:41
Src/Base/AMReX_Geometry.cpp Outdated Show resolved Hide resolved
@atmyers
Copy link
Member

atmyers commented Apr 6, 2023

Thank you for this!

@atmyers atmyers merged commit ce3401f into AMReX-Codes:development Apr 6, 2023
@atmyers atmyers mentioned this pull request Apr 7, 2023
5 tasks
atmyers added a commit that referenced this pull request Apr 7, 2023
PR #3199 broke the Nyx regression tests. The issue is we require the
particle positions to be in [plo, phi). This won't happen if rhi is
exactly equal to phi. The fix is to make sure rhi is always slightly
inside the domain.

EDIT - fixed this in a different way. Now, when applying periodic
boundaries, the particle positions will be set to be strictly less than
rhi.
WeiqunZhang added a commit that referenced this pull request Apr 10, 2023
There was a flaw in amrex::Geometry::computeRoundoffDomain. If probhi is
close to zero, std::nextafter towards negative infinity will be a very
small number (e.g., 1.e-300). So the function will be able to find the
lowest value of roundoff_hi that's outside the domain within a
reasonable number of iterations.

In the new implementation, we use bisection to find the lowest value of
roundoff_lo that's inside the domain and the highest value fo
roundoff_hi that's inside the domain, up to a tolerance.

X-ref:
- #3243
- #3199
guj pushed a commit to guj/amrex that referenced this pull request Jul 13, 2023
## Summary
This reworks previous changes that created the round off domain which
was done to avoid round off errors when locating particles in the
domain, specifically when handling periodic boundaries. This PR changes
how this is done so that the roundoff domain makes as small a change as
possible.

## Additional background
Previous PRs AMReX-Codes#2839 and AMReX-Codes#2950 created a round off domain that has a
slightly smaller extent that the actual domain. Any particles that ended
up outside the round off domain but inside the actual domain were
shifted to be on the boundary of the round off domain. The main purpose
was to handle problems with periodic boundary with mixed precision
(single precision particles and double precision grid). In that case, it
is possible for a particle to have `z < lo` and `z+(hi-lo) > hi`. It is
also possible for a particle to have `z < hi` but `(z-lo)*dzinv > ihi`
putting it in a grid cell beyond the domain. Unfortunately, the roundoff
domains was being set further inside that was needed and would affect
valid particles that would not have those issues. While the problem is
small, the effects can be noticeable when examining simulation
properties, for example charge conservation in WarpX.

This PR changes the roundoff domain so that only problematic particles
are affected. The code now uses `std::nextafter` to find the exact place
where the problems occurs. At the lower end, it ensures that the casted
position will be greater than ProbLo. At the upper end, it ensures that
`(z-lo)*dzinv <= ihi`.

One comment is that this assumes that particle grid coordinates will
always be calculated at the higher precision, with expressions like `iz
= (int)(z - lo)`, without `lo` being cast to `ParticleReal`.

Co-authored-by: atmyers <[email protected]>
guj pushed a commit to guj/amrex that referenced this pull request Jul 13, 2023
PR AMReX-Codes#3199 broke the Nyx regression tests. The issue is we require the
particle positions to be in [plo, phi). This won't happen if rhi is
exactly equal to phi. The fix is to make sure rhi is always slightly
inside the domain.

EDIT - fixed this in a different way. Now, when applying periodic
boundaries, the particle positions will be set to be strictly less than
rhi.
guj pushed a commit to guj/amrex that referenced this pull request Jul 13, 2023
There was a flaw in amrex::Geometry::computeRoundoffDomain. If probhi is
close to zero, std::nextafter towards negative infinity will be a very
small number (e.g., 1.e-300). So the function will be able to find the
lowest value of roundoff_hi that's outside the domain within a
reasonable number of iterations.

In the new implementation, we use bisection to find the lowest value of
roundoff_lo that's inside the domain and the highest value fo
roundoff_hi that's inside the domain, up to a tolerance.

X-ref:
- AMReX-Codes#3243
- AMReX-Codes#3199
WeiqunZhang pushed a commit to WeiqunZhang/amrex that referenced this pull request Sep 28, 2023
  * Rework handling of roundoff domain extent (AMReX-Codes#3199)

  * Fix periodic boundary bug in AMReX-Codes#3199 (AMReX-Codes#3243)

  * Fix Roundoff Domain (AMReX-Codes#3247)

  * Add Tests/RoundoffDomain (AMReX-Codes#3248)

  * Roundoff errors in computeRoundoffDomain (AMReX-Codes#3255)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants