Use a new _ParametricGeometry for specializations

[mikeingold]

Changes

• Implements an internal _ParametricGeometry <: Meshes.Geometry type to serve as a simple wrapper for custom parametric functions.
• Uses this new geometry type to replace specialization methods that used hard-coded domain transformations and analytic methods. The new methods simply generate a _ParametricGeometry with domain-transformed parametric function, and then use the standard _integral code paths. This new process was applied to specialization methods for:
  • BezierCurve
  • Line
  • Plane
  • Ray
  • Tetrahedron
  • Triangle
• Deprecates jacobian(::BezierCurve, ts, ::Analytical) method - performance was inferior to the generic FiniteDifference method.
• Changes to internal utility functions:
  • Removes now-obsolete _guarantee_analytical
  • Defines a _parametric(::Geometry, ts...) used for specializations
  • _zeros([T,] N) and _ones([T,] N) return an NTuple{N,T} of zeros/ones (faster than allocating a vector with Base.zeros/ones)
• Added a new Specializations section to benchmarks
  • Moved out BezierCurve from normal section, reducing redundancy.
  • Added benchmark for Triangle.

Results

• Normalizes support for the six specialization types listed above, i.e. they can now be expected to behave like any other curve/surface/volume.
• Closes Integral on Tetrahedron not fully supported #40 by adding full support for integrating Tetrahedron types.
• Mostly completes Use new Parametric Geometry types for internal domain transformations #128: will update to indicate that we’re waiting on a proper implementation upstream to replace _ParametricGeometry in the future.
• Significant performance increases
  • Massive improvement (~75x faster) for integrals over a BezierCurve.
  • Significant improvement (up to 4x) for scalar integrands with HAdaptiveCubature rule.
• Significantly reduces code and math complexity.
• Incurs a performance regression (~4x slower) for integrals over a Triangle due to new requirement to call differential inside integrand (vs baked-in analytical), but simultaneously discovered and fixed a bug where Unitful functions couldn't be integrated with the old version.

TODO

• Update Documenter specializations page

[mikeingold]

I'm trying out a slightly (moderately (strongly)) hacky solution. I don't love it, but it gets the analytical jacobian for BezierCurve working, so let's see how performance compares.

[mikeingold]

Update: The hacky solution was able to get the jacobian(bezier, ts, ::Analytical) back in use through a _ParametricGeometry wrapper. Surprisingly, performance was not significantly different vs using the FiniteDifference method: ~80x faster than main with Analytical but ~85x with FiniteDifference. It seems like something else entirely about the old integral methods for BezierCurve was just plain slower. I'm just going to call this one dead and deprecate this analytical method (and un-roll the latest hacky test changes).

[mikeingold]

cc: @JoshuaLampert, @kylebeggs

Just need to update the docs and I think this one is otherwise just about done.

[JoshuaLampert]

So this means has_analytical always returns false, which means we don't need the Analytical differentiation method anymore? Would it be the easiest to simply remove the option entirely? Then we would also not need default_method anymore, but simply make FiniteDifference always the default. What do you think?
The new method for integrating Triangles is ~5 times slower now. What do you think about that?
Also: Now we are using the finite difference approximation also for the geometries using the _ParametricGeometry, right? So why don't we lose accuracy if we use an approximation to the Jacobian instead of the exact Jacobian as before? So why are the tests not failing?

[kylebeggs]

Sorry, I have not been able to keep up with this PR. However, if this is true:

The new method for integrating Triangles is ~5 times slower now. What do you think about that?

I would urge to fix this before merging. Just my opinion but performance should be near the top in priority in a package like this.

[mikeingold]

So this means has_analytical always returns false, which means we don't need the Analytical differentiation method anymore? Would it be the easiest to simply remove the option entirely? Then we would also not need default_method anymore, but simply make FiniteDifference always the default. What do you think?

I think this is definitely worth considering. I had visions of defining analytical methods for geometries like Box, Segment, etc where they're fairly simple to derive. However, I'm not sure they'd actually be significantly faster than FiniteDifference, or any more accurate than Auto(AD), so it might not be worth the extra complexity/maintenance/testing.

The new method for integrating Triangles is ~5 times slower now. What do you think about that?

I'm a little disappointed that performance decreased for Triangle, granted eval times are still well into the microsecond regime, which is still reasonably fast. On balance with the significant reduction in complexity (no more _integral-alikes outside of src/integrals.jl), fixed handling of Unitful integrands, etc I think it might be worth the trade.

Also: Now we are using the finite difference approximation also for the geometries using the _ParametricGeometry, right? So why don't we lose accuracy if we use an approximation to the Jacobian instead of the exact Jacobian as before? So why are the tests not failing?

This is a good question. I haven't done any real numerical sensitivity analysis but I suspect that most practical geometries have parametric functions that evolve fairly slowly, so the sensitivity of the result to changes in step size is likely reduced, i.e. 1e-6 step size -> <1e-6 loss of accuracy. This would almost definitely fail to hold in the extreme where some less-well-behaved geometry has been defined, e.g. a fractal.

[mikeingold]

@kylebeggs, I didn't see your comment until after posting mine. I'm going to need some time to dial in the implementation with this one, so (if you've got the time) go ahead and press forward with #129 (`AutoEnzyme`) and we can plan to merge that one first.

…

On Thu, Nov 14, 2024 at 8:25 AM Kyle Beggs ***@***.***> wrote: Sorry, I have not been able to keep up with this PR. However, if this is true: The new method for integrating Triangles is ~5 times slower now. What do you think about that? I would urge to fix this before merging. Just my opinion but performance should be near the top in priority in a package like this. — Reply to this email directly, view it on GitHub <#131 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/AA657OUMNXB3Z6AIAPY3SEL2ASQGDAVCNFSM6AAAAABRQ5Y7ZSVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDINZWGM2TGOJXHA> . You are receiving this because you were assigned.Message ID: ***@***.***>

[mikeingold]

Performance with this change is oddly kind of a mixed bag. There's actually a couple of ideas contained in this PR, so I think I'm going to partition it into multiple smaller changes to make them easier to review and reason about.