-
Notifications
You must be signed in to change notification settings - Fork 92
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
Wedge Support #581
Wedge Support #581
Conversation
Btw the fail in the CI is unrelated to this PR, but I found it when improving the tests. The way the quintic triangle is evaluated is not really stable (numerically). |
Cool, this is an important case for improved generality. |
CI is broken due to the deprecation. Refinement leads to decrease in the L_infty norm on Poisson problems with manufactured solutions, indicating convergence. Edit: Test coverage for the new elements is incomplete right now and will be completed when the review progresses. Edit 2: Linear wedges are exportable to VTK. To visually inspect the quadratic wedges please use FerriteViz#do/crincle-clip (Ferrite-FEM/FerriteViz.jl#56) together with this code to make wedges functional function FerriteViz.decompose!(coord_offset, coord_matrix, ref_coord_matrix, triangle_offset, triangle_matrix, grid, cell::Ferrite.AbstractCell{3,N,5}) where {N}
faces = Ferrite.faces(cell)
(coord_offset, triangle_offset) = FerriteViz.decompose!(coord_offset, coord_matrix, ref_coord_matrix, triangle_offset, triangle_matrix, grid, Ferrite.Cell{3,3,1}(faces[1]))
(coord_offset, triangle_offset) = FerriteViz.decompose!(coord_offset, coord_matrix, ref_coord_matrix, triangle_offset, triangle_matrix, grid, Ferrite.Cell{3,4,1}(faces[2]))
(coord_offset, triangle_offset) = FerriteViz.decompose!(coord_offset, coord_matrix, ref_coord_matrix, triangle_offset, triangle_matrix, grid, Ferrite.Cell{3,4,1}(faces[3]))
(coord_offset, triangle_offset) = FerriteViz.decompose!(coord_offset, coord_matrix, ref_coord_matrix, triangle_offset, triangle_matrix, grid, Ferrite.Cell{3,4,1}(faces[4]))
(coord_offset, triangle_offset) = FerriteViz.decompose!(coord_offset, coord_matrix, ref_coord_matrix, triangle_offset, triangle_matrix, grid, Ferrite.Cell{3,3,1}(faces[5]))
(coord_offset, triangle_offset)
end
FerriteViz.ntriangles(cell::Ferrite.AbstractCell{3,N,5}) where {N} = 1+4+4+4+1
FerriteViz.getfaceip(ip::Ferrite.Lagrange{3, RefPrism, order}, local_face_idx::Int) where {order} = local_face_idx in 2:4 ? Ferrite.Lagrange{2, RefCube, order}() : Ferrite.Lagrange{2, RefTetrahedron, order}()
FerriteViz.linear_face_cell(cell::Ferrite.Cell{3,N,5}, local_face_idx::Int) where N = local_face_idx in 2:4 ? Ferrite.Quadrilateral3D(Ferrite.faces(cell)[local_face_idx]) : Ferrite.Cell{3,3,1}(Ferrite.faces(cell)[local_face_idx]) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Cool!
I looked through, out of interest, and spotted some minor things.
The faulty assumption
has to be fixed in a subsequent PR. I could not find a simple way to propagate the shell dof elimination through the constraint handler without breaking dof elimination for mixed grids. |
Codecov ReportPatch coverage:
📣 This organization is not using Codecov’s GitHub App Integration. We recommend you install it so Codecov can continue to function properly for your repositories. Learn more Additional details and impacted files@@ Coverage Diff @@
## master #581 +/- ##
==========================================
- Coverage 92.66% 92.65% -0.02%
==========================================
Files 28 29 +1
Lines 4270 4441 +171
==========================================
+ Hits 3957 4115 +158
- Misses 313 326 +13
... and 1 file with indirect coverage changes Help us with your feedback. Take ten seconds to tell us how you rate us. Have a feature suggestion? Share it here. ☔ View full report in Codecov by Sentry. |
I think #477 handles shell interpolations better. Maybe merging that PR will makes life easier in this PR? |
It would definitely make things easier, but the underlying issue that the semantics of edge and face is different in different sections of the codebase stays. |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This should summarize the remaining issues and clarify some stuff @fredrikekre
src/interpolations.jl
Outdated
@@ -154,7 +179,7 @@ indices of [`reference_coordinates(::Interpolation)`](@ref). | |||
value(ip::Interpolation, i::Int, ξ::Vec) | |||
|
|||
""" | |||
reference_coordinates(::Interpolation) | |||
Ferrite.reference_coordinates(::Interpolation) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
When should we include the Ferrite.
prefix in the doc string?
TODO for me: Put answer into devdocs.
…tion is unstable.
2f7d8ee
to
92ba486
Compare
Thanks for working through this ! |
This patch reworks the `AbstractCell` interface (refer to #677 for motivation and discussion for alternatives to this patch). Specifically, this patch - introduces the reference dimension as a type parameter on `AbstractRefShape`, - introduces `RefInterval`, `RefTriangle`, `RefQuadrilateral`, and `RefHexahedron` as new subtypes of `AbstractRefShape{refdim}` (note that interpolations are not updated to this system yet since it can be done in a separate patch), - deprecates the old `Cell{refdim, nnodes, nfaces}` struct, - introduces new structs for all supported cell types, e.g. `struct Triangle <: AbstractCell{RefTriangle}`, - defines the unique vertex/edge/face identifiers for DoF-distribution in terms of the reference shape. The new parametrization makes it possible to dispatch on the reference shape instead of union of concrete implementations, i.e. `::AbstractCell{RefTriangle}` instead of `::Union{Triangle, QuadraticTriangle}`, and similarly for the reference dimension. One drawback is that it is no longer possible to use "anonymous celltypes" after this patch, i.e. using, e.g., `Cell{3, 20, 6}` and defining methods for that instead of explicitly naming it `SerendipityQuadraticHexahedron` or similar. Instead, you now have to define a struct and some methods. Arguably this is a good thing -- nobody understands the meaning of `Cell{3, 20, 6}`. For normal usage this patch is not breaking (see e.g. no required changes for examples), unless one dispatches on the `Cell` struct directly. Specific things that remain to be decided: - `struct Triangle <: AbstractCell{2, RefTriangle}` or `struct Triangle <: AbstractRefCell{2, RefTriangle} <: AbstractCell`? I don't know if it is too restrictive to "force" `{refdim, refshape}` on all subtypes of `AbstractCell`. On the other hand, if one implements a non-refshape based element one can define a dummy reference shape. For a real-world example, looking at e.g. https://github.com/kimauth/FerriteCohesiveZones.jl/blob/main/src/cohesive_cell.jl I think that would simply be implemented as ```julia struct CohesiveQuadrilateral <: AbstractCell{2, RefQuadrilateral} nodes::NTuple{4, Int} end ``` which I don't see a direct problem with, other than the fact that the cell will now have 4 edges (according to the fallback definitions). To solve this, one would have to implement a method for `faces(::CohesiveQuadrilateral)` that only return the upper and lower edge. I don't think this is too much to ask for custom elements like this. Alternatively, with #581 merged it is possible to distribute different number of dofs for different edges, so this can also be solved on the (default) interpolation level by distributing 0 dofs on the left and right edge. - `AbstractCell{refdim, refshape}` or just `AbstractCell{refshape}`? The dimension is implicit from the shape now, and you can dispatch on `AbstractCell{<:AbstractRefShape{refdim}} where refdim`. Personally I am in favor of including the dimension, even if redundant. It does not cost anything, and sometimes you want to dispatch on the dimension, sometimes on the refshape. Fixes #677.
After #581 `vertexdof_indices` should return a tuple of tuples, but since `length(::Int) == length(::Tuple{Int}) == 1` this old code worked anyway. This patch removes the method to rely on the generic `Lagrange` method instead.
After #581 `vertexdof_indices` should return a tuple of tuples, but since `length(::Int) == length(::Tuple{Int}) == 1` this old code worked anyway. This patch removes the method to rely on the generic `Lagrange` method instead.
This PR contains a fully functional code to support linear and quadratic wedges (see #580) for both dof handlers. Addresses #272 as a side-effect.