kinnala / scikit-fem

Simple finite element assemblers
https://scikit-fem.readthedocs.io
BSD 3-Clause "New" or "Revised" License
470 stars 76 forks source link

Lowest order integration rules for tets #1110

Closed kinnala closed 1 month ago

kinnala commented 3 months ago

Right now the lowest order integration rule for tets has 4 points. The following rule is suitable for P1 Laplacian though:

quadrature=(np.array([[.25],[.25],[.25]]), np.array([1/6]))

It should be added to get_quadrature.

adtzlr commented 1 month ago

Great, this saves some computation time and memory for linear elasticity.

The only downside (if it is really considered to be one) of it is that projections do not work with this lowest-order quadrature. But it is so obvious that this should not be a problem. I think a typical scikit-fem user has enough knowledge.

kinnala commented 1 month ago

I think it won't be a problem since usually it would not enabled by default, user has to write Basis(..., intorder=1). Perhaps it will be enabled by default if initializing a basis for Basis(..., ElementTetP0())). I will need to check this.

kinnala commented 1 month ago

So basically, before this change, if user initializes Basis(mesh, ElementTetP0()) it chose automatically the rule with 4 integration points while, after the change, it would automatically choose the rule with 1 integration point. So there might be some surprises involved if a user relied on this kind of behaviour.

kinnala commented 1 month ago

I am considering of circumventing the aforementioned issue by increasing ElementTetP0.maxdeg from 0 to 1.

kinnala commented 1 month ago

Nah, I will let the speed benefit show there.