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

Add InverseJacobianInertialToFluidCompute #5811

Merged

Conversation

ffoucart
Copy link
Contributor

@ffoucart ffoucart commented Feb 27, 2024

Proposed changes

Implement a ComputeTag that calculates the inverse jacobian of the map between the inertial frame and an orthonormal tetrad comoving with the fluid.

Upgrade instructions

Code review checklist

  • The code is documented and the documentation renders correctly. Run
    make doc to generate the documentation locally into BUILD_DIR/docs/html.
    Then open index.html.
  • The code follows the stylistic and code quality guidelines listed in the
    code review guide.
  • The PR lists upgrade instructions and is labeled bugfix or
    new feature if appropriate.

Further comments

I needed to add a new frame (Frame::Fluid) -- should it be in the same file as the other frames, or defined separately in the MC code?

@ffoucart ffoucart force-pushed the InverseJacobianInertialToFluid branch from c2a263b to 9304927 Compare February 28, 2024 18:05
@ffoucart
Copy link
Contributor Author

ffoucart commented Feb 28, 2024

Changed the accuracy required in the test; the default worked in Release mode on wheeler (at least every time I tried it) but consistently failed in Debug mode and for many of the configurations run by Jenkins.

gr::Tags::SpatialMetric<DataVector, 3> >;

static void function(gsl::not_null<return_type*> inv_jacobian,
const tnsr::I<DataVector, 3>& spatial_velocity,
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Fix formatting.

(spatial_velocity.get(d) - shift.get(d) / get(lapse));
}

// Then, the other members of the tetrad are constructed using Gram-Shmidt
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

"Schmidt"

// Temporary memory allocation
DataVector temp_dot_product = get(lapse) * 0.0;
std::array<DataVector, 4> tetrad_vector{temp_dot_product, temp_dot_product,
temp_dot_product, temp_dot_product};
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Instead of using this temporary, write the result directly into the inv_jacobian matrix. That saves four memory allocations, which is good for both speed and memory use.

// Then, the other members of the tetrad are constructed using Gram-Shmidt

// Temporary memory allocation
DataVector temp_dot_product = get(lapse) * 0.0;
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

auto temp_dot_product = make_with_value<DataVector>(lapse, 0.0);

@ffoucart ffoucart force-pushed the InverseJacobianInertialToFluid branch 2 times, most recently from f315183 to 78bb9f6 Compare February 29, 2024 00:01
Copy link
Member

@wthrowe wthrowe left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Other changes look good. Squash, including this last request.


// d = 1,2,3 for tetrad components built from x,y,z
for (size_t d = 1; d < 4; d++) {
// Base vector for Gram-Shmidt
for (size_t i = 0; i < 4; i++) {
gsl::at(tetrad_vector, i) = (i == d ? 1.0 : 0.0);
inv_jacobian->get(i, d) = make_with_value<DataVector>(lapse, 0.0);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do set_number_of_grid_points(&inv_jacobian->get(i, d), lapse) (from Utilities/SetNumberOfGridPoints.hpp) and then = 0.0 (or the ?: expression you had before). This avoids a DataVector allocation. Might need a make_not_null.

@ffoucart ffoucart force-pushed the InverseJacobianInertialToFluid branch from 78bb9f6 to 76e8d92 Compare February 29, 2024 20:57
@ffoucart
Copy link
Contributor Author

Squashed with last requested change.

Copy link
Member

@wthrowe wthrowe left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Test failed in CI:

2024-02-29T23:23:33.9252167Z Unit.Evolution.Particles.MonteCarlo.InverseJacobianInertialToFluid
2024-02-29T23:23:33.9253216Z -------------------------------------------------------------------------------
2024-02-29T23:23:33.9255545Z ../../../../../../tests/Unit/Evolution/Particles/MonteCarlo/Test_InverseJacobianInertialToFluid.cpp:22
2024-02-29T23:23:33.9256689Z ...............................................................................
2024-02-29T23:23:33.9257210Z 
2024-02-29T23:23:33.9257566Z ../../../../../../tests/Unit/Framework/TestingFramework.hpp:116: FAILED:
2024-02-29T23:23:33.9258285Z   CHECK( a == appx(b) )
2024-02-29T23:23:33.9258704Z with expansion:
2024-02-29T23:23:33.9259514Z   -0.00000000000212888 == Approx( 0.0 )
2024-02-29T23:23:33.9260039Z with messages:
2024-02-29T23:23:33.9260650Z   Seed is: 2708572239 from ../../../../../../tests/Unit/Evolution/Particles/
2024-02-29T23:23:33.9261452Z   MonteCarlo/Test_InverseJacobianInertialToFluid.cpp:26
2024-02-29T23:23:33.9262207Z   ../../../../../../tests/Unit/Evolution/Particles/MonteCarlo/
2024-02-29T23:23:33.9263027Z   Test_InverseJacobianInertialToFluid.cpp:88: dot_product ==
2024-02-29T23:23:33.9263727Z   expected_dot_product
2024-02-29T23:23:33.9264480Z   a := (-1.55761e-14,3.35154e-17,1.85832e-16,-1.83013e-16,-2.12888e-12)
2024-02-29T23:23:33.9265165Z   b := (0,0,0,0,0)
2024-02-29T23:23:33.9265402Z 
2024-02-29T23:23:33.9265916Z 0.001 s: Unit.Evolution.Particles.MonteCarlo.InverseJacobianInertialToFluid
2024-02-29T23:23:33.9266854Z ===============================================================================
2024-02-29T23:23:33.9267485Z test cases:   1 |   0 passed | 1 failed
2024-02-29T23:23:33.9268385Z assertions: 127 | 126 passed | 1 failed

@ffoucart
Copy link
Contributor Author

ffoucart commented Mar 1, 2024

Is there a way to force the use of that seed in the test locally?
I could just raise the tolerance a bit (it fails the comparison by 2e-12 when the tolerance is 1e-12), but unless it's somehow using weird random values for the metric/velocity, I'd expect higher accuracy.

@wthrowe
Copy link
Member

wthrowe commented Mar 1, 2024

2024-02-29T23:23:33.9260650Z   Seed is: 2708572239 from ../../../../../../tests/Unit/Evolution/Particles/
2024-02-29T23:23:33.9261452Z   MonteCarlo/Test_InverseJacobianInertialToFluid.cpp:26

That line should have a MAKE_GENERATOR call on it. That macro takes an optional second argument of the seed.

@ffoucart ffoucart force-pushed the InverseJacobianInertialToFluid branch from 76e8d92 to e8242e3 Compare March 1, 2024 19:28
@ffoucart
Copy link
Contributor Author

ffoucart commented Mar 1, 2024

I changed the epsilon from 1e-12 to 5e-12.
The failing point had a Lorentz factor of 20.8. The maximum allowed by the RNG is 21.1, and it makes sense that the calculation would loose accuracy for high Lorentz factors so... in my view it's ok for now.

Copy link
Member

@wthrowe wthrowe left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

OK, sounds good.

@wthrowe wthrowe merged commit 477c3b2 into sxs-collaboration:develop Mar 5, 2024
22 checks passed
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