search

krcools / CompScienceMeshes.jl

Computational Geometry Foundations for Finite and Boundary Element Methods
Other
28 stars 22 forks source link

Inconsistent quadrature rules #41

Open surzks opened 3 months ago

surzks commented 3 months ago

Currently, the quadrature rules for triangles are a mixture of Dunavant's rules and some rules based on a tensor-product approach.

To be consistent, I replaced the latter ones with the remaining Dunavant rules in PR #40

Moreover, they are now collected in SVectors, which yields a speed-up.

TODO: These days, there are better quadrature rules available, as can be seen in the QuadPy package. However, it might be better to write a dedicated Julia package to collect these.

sbadrian commented 3 months ago

We just noticed that some of the unit tests of BEAST fail when using the Dunavant quadratures. Even when we increased from order=13 (the maximum with the old quadrature) to 20. Naturally, the old quadratures, being based on tensor products, can integrate functions exactly that exceed order = 13 at the price of many more quadrature points.

Still, I am hesitating to propose reducing the accuracy in the BEAST unit tests. Probably, the best would be, if several quadratures are concurrently available and the user can select them.

Maybe this is possible by some wrapping type? Instead of passing 13 one would pass TensorQuad(13). Of course, integers could be interpreted to default to Dunavant, as this should be what in most cases would be needed.

In future, it would make it easier to integrate other strategies which are potentially better than Dunavant.

krcools commented 2 months ago

Fully agree. This is what I had in mind originally. I suggest to create wrapping types both for tensor and Duvanant and let rules of type Int dispatch to the original inconsistent series for backward compatibility.

sbadrian commented 2 months ago

In PR #43, the Dunavant rules are added as extra rules and we default to the old quadrature. What do you think of the layout of the types?