trilinos / Trilinos

Primary repository for the Trilinos Project
https://trilinos.org/
Other
1.19k stars 565 forks source link

intrepid2: slightly inaccurate values for `CubatureDirectTriSymmetric` should be rewritten #12260

Closed romintomasetti closed 1 week ago

romintomasetti commented 11 months ago

@trilinos/intrepid2

@mperego

Motivation

Some of our unit tests failed because HCURL on TRI now uses CubatureDirectTriSymmetric. This is expected because some of the values have been written "inaccurately" (sorry for the rude word), though their exact value is known. This is clear in the code excerpt below wherein one should use 1./3. or 1./6., and most importantly ensure that values that should be the same are actually equal (e.g. 0.3333333333333333 is not equal to 0.3333333333333334).

https://github.com/trilinos/Trilinos/blob/57a0cf441dca123b038d6d046e5abfa9c52d8a81/packages/intrepid2/src/Discretization/Integration/Intrepid2_CubatureDirectTriSymmetricDef.hpp#L115-L132

This was properly done e.g. here:

https://github.com/trilinos/Trilinos/blob/57a0cf441dca123b038d6d046e5abfa9c52d8a81/packages/intrepid2/src/Discretization/Integration/Intrepid2_CubatureDirectTetSymmetricDef.hpp#L117

I can file a PR for that :wink:

romintomasetti commented 11 months ago

@maartenarnst

mperego commented 11 months ago

Wow, you have very sensitive tests. It would be great if you can push a PR for improving the accuracy of the quadrature points/weights. Ideally, it would great to compute those quadratures at higher precision, but we used values provided in a paper rather than computing them from scratch.

maartenarnst commented 11 months ago

It's PR #12261. It uses the rational numbers for the two lowest-order quadrature rules, where it's kind of immediate to know which rational numbers it should be.

mperego commented 11 months ago

B.t.w. are you doing a bfb comparison in your test? I'd be alarmed if a 10^-15 error on the quad points would result in errors well above machine precision

maartenarnst commented 11 months ago

Thanks for approving!

I'm not sure what you mean with a "bfb" comparison.

What happened is that we wrote quite a few unit tests that are a mix between documenting and testing functionalities of Intrepid2. One of those unit tests is checking values taken by HCurl basis functions on the TRI. We use those defined in Intrepid2_HCURL_TRI_In_FEM.hpp. And after we updated our Trilinos version, this unit test detected that these values no longer match the expected values up to the tolerance that we had set. We chose a tolerance that was high.

We suspect that the reason for the changed behaviour is that the code in Intrepid2_HCURL_TRI_In_FEM.hpp now uses the CubatureDirectTriSymmetric quadrature rule, and we suspect that this change has an impact on the computation of the matrices that are used in the basis transformations to ultimately get to the basis functions.

mperego commented 11 months ago

@maartenarnst, with bfb comparison I meant bit-for-bit comparison, that is, if you check that the computed values exactly match the reference expected values. I'm a bit surprised that your test fails for such small differences in the quadrature point values, even if you use tight tolerances.

mperego commented 1 week ago

I think this issue was addressed by PR #12261. Closing.