diff --git a/src/Dofs/DofHandler.jl b/src/Dofs/DofHandler.jl index e5ea543fd9..f5a4755c74 100644 --- a/src/Dofs/DofHandler.jl +++ b/src/Dofs/DofHandler.jl @@ -644,7 +644,7 @@ described therein. # typestable way. n_copies = step(dofs) @assert n_copies > 0 - if adjust_during_distribution && !orientation.regular + if adjust_during_distribution && orientation.flipped # Reverse the dofs for the path dofs = reverse(dofs) end @@ -661,12 +661,12 @@ end Returns the unique representation of an edge and its orientation. Here the unique representation is the sorted node index tuple. The -orientation is `true` if the edge is not flipped, where it is `false` +orientation is `false` if the edge is not flipped, where it is `true` if the edge is flipped. """ function sortedge(edge::Tuple{Int, Int}) a, b = edge - return a < b ? (edge, PathOrientationInfo(true)) : ((b, a), PathOrientationInfo(false)) + return a < b ? (edge, PathOrientationInfo(false)) : ((b, a), PathOrientationInfo(true)) end """ @@ -739,7 +739,7 @@ function sortface(face::Tuple{Int, Int, Int}) b, c = minmax(b, c) a, c = minmax(a, c) a, b = minmax(a, b) - return (a, b, c), SurfaceOrientationInfo() # TODO fill struct + return (a, b, c), SurfaceOrientationInfo(face) end @@ -760,7 +760,7 @@ function sortface(face::Tuple{Int, Int, Int, Int}) b, c = minmax(b, c) a, c = minmax(a, c) a, b = minmax(a, b) - return (a, b, c), SurfaceOrientationInfo() # TODO fill struct + return (a, b, c), SurfaceOrientationInfo(face) end diff --git a/src/FEValues/InterfaceValues.jl b/src/FEValues/InterfaceValues.jl index 64ea4dc3fd..e06c993832 100644 --- a/src/FEValues/InterfaceValues.jl +++ b/src/FEValues/InterfaceValues.jl @@ -419,8 +419,8 @@ Return the relative orientation info for facet B with regards to facet A. Relative orientation is computed using a [`OrientationInfo`](@ref) for each side of the interface. """ function InterfaceOrientationInfo(cell_a::AbstractCell{RefShapeA}, cell_b::AbstractCell{RefShapeB}, facet_a::Int, facet_b::Int) where {RefShapeA <: AbstractRefShape, RefShapeB <: AbstractRefShape} - OI_a = OrientationInfo(facets(cell_a)[facet_a]) - OI_b = OrientationInfo(facets(cell_b)[facet_b]) + OI_a = orientation_info(facets(cell_a)[facet_a]) + OI_b = orientation_info(facets(cell_b)[facet_b]) flipped = OI_a.flipped != OI_b.flipped shift_index = OI_b.shift_index - OI_a.shift_index return InterfaceOrientationInfo{RefShapeA, RefShapeB}(flipped, shift_index, OI_b.shift_index, facet_a, facet_b) diff --git a/src/Grid/grid.jl b/src/Grid/grid.jl index 2d3ca08552..2daa8a823c 100644 --- a/src/Grid/grid.jl +++ b/src/Grid/grid.jl @@ -664,8 +664,28 @@ end ################################# # @TODO merge this code with into the logic in `ConstraintHandler`. +function get_edge_direction(edgenodes::NTuple{2, Int}) + positive = edgenodes[2] > edgenodes[1] + return ifelse(positive, 1, -1) +end + +function get_face_direction(facenodes::NTuple{N, Int}) where {N} + N > 2 || throw(ArgumentError("A face must have at least 3 nodes")) + min_idx = argmin(facenodes) + if min_idx == 1 + positive = facenodes[2] < facenodes[end] + elseif min_idx == length(facenodes) + positive = facenodes[1] < facenodes[end - 1] + else + positive = facenodes[min_idx + 1] < facenodes[min_idx - 1] + end + return ifelse(positive, 1, -1) +end + """ PathOrientationInfo + PathOrientationInfo(flipped:Bool) + PathOrientationInfo(edgenodes::NTuple{2, Int}) Orientation information for 1D entities. @@ -674,117 +694,73 @@ associated to the vertices. To give an example, the oriented path ``` 1 ---> 2 ``` -is called *regular*, indicated by `regular=true`, while the oriented path +is called *regular*, indicated by `flipped=false`, while the oriented path ``` 2 ---> 1 ``` -is called *inverted*, indicated by `regular=false`. +is called *inverted*, indicated by `flipped=true`. """ struct PathOrientationInfo - regular::Bool # Indicator whether the orientation is regular or inverted. + flipped::Bool # Indicator whether the orientation is regular or inverted. + shift_index::Int + + function PathOrientationInfo(flipped::Bool) + return new(flipped, 0) + end +end + +function PathOrientationInfo(edgenodes::NTuple{2, Int}) + return PathOrientationInfo(get_edge_direction(edgenodes) < 0) end """ SurfaceOrientationInfo + SurfaceOrientationInfo(flipped:Bool, shift_index::Int) + SurfaceOrientationInfo(facenodes::NTuple{N, Int}) Orientation information for 2D entities. Such an entity can be possibly flipped (i.e. the defining vertex order is reverse to the spanning vertex order) and the vertices can be rotated against each other. Take for example the faces ``` -1---2 2---3 +4---3 3---2 | A | | B | -4---3 1---4 +1---2 4---1 ``` which are rotated against each other by 90° (shift index is 1) or the faces ``` -1---2 2---1 +4---3 2---3 | A | | B | -4---3 3---4 +1---2 1---4 ``` which are flipped against each other. Any combination of these can happen. The combination to map this local face to the defining face is encoded with this data structure via ``rotate \\circ flip`` where the rotation is indiced by the shift index. - !!!NOTE TODO implement me. """ struct SurfaceOrientationInfo - #flipped::Bool - #shift_index::Int -end - - -@doc raw""" - InterfaceOrientationInfo - -Orientation information for 1D and 2D entities. -The orientation is defined by the indices of the grid nodes -associated to the vertices. To give an example, the oriented path -``` -1 ---> 2 -``` -is called *regular*, indicated by `flipped=false`, while the oriented path -``` -2 ---> 1 -``` -is called *inverted*, indicated by `flipped=true`. - -2D entities can be flipped (i.e. the defining vertex order is reverse to the -spanning vertex order) and the vertices can be rotated against each other. - -The reference entity is a one with it's first node is the lowest index vertex -and its vertices span counter-clock-wise. -Take for example the faces -``` -1 2 -| \ | \ -| \ | \ -| A \ | B \ -| \ | \ -2-----3 3-----1 -``` -which are rotated against each other by 240° after tranfroming to an -equilateral triangle (shift index is 2). Or the faces -``` -3 2 -| \ | \ -| \ | \ -| A \ | B \ -| \ | \ -2-----1 3-----1 -``` -which are flipped against each other. -""" -struct OrientationInfo flipped::Bool shift_index::Int end -function OrientationInfo(edgenodes::NTuple{2, Int}) - return OrientationInfo(get_edge_direction(edgenodes) < 0, 0) -end - -function OrientationInfo(facenodes::NTuple{N, Int}) where {N} +function SurfaceOrientationInfo(facenodes::NTuple{N, Int}) where {N} min_idx = argmin(facenodes) shift_index = min_idx - 1 flipped = get_face_direction(facenodes) < 0 - return OrientationInfo(flipped, shift_index) + return SurfaceOrientationInfo(flipped, shift_index) end -function get_edge_direction(edgenodes::NTuple{2, Int}) - positive = edgenodes[2] > edgenodes[1] - return ifelse(positive, 1, -1) + +""" + orientation_info(edgenodes::NTuple{2, Int}) + orientation_info(facenodes::NTuple{N, Int}) + +Helper function that returns a `PathOrientationInfo` or `SurfaceOrientationInfo` as appropriate for the size of the node tuple. +""" +function orientation_info(edgenodes::NTuple{2, Int}) + return PathOrientationInfo(edgenodes) end -function get_face_direction(facenodes::NTuple{N, Int}) where {N} - N > 2 || throw(ArgumentError("A face must have at least 3 nodes")) - min_idx = argmin(facenodes) - if min_idx == 1 - positive = facenodes[2] < facenodes[end] - elseif min_idx == length(facenodes) - positive = facenodes[1] < facenodes[end - 1] - else - positive = facenodes[min_idx + 1] < facenodes[min_idx - 1] - end - return ifelse(positive, 1, -1) +function orientation_info(facenodes::NTuple{N, Int}) where {N} + return SurfaceOrientationInfo(facenodes) end diff --git a/test/runtests.jl b/test/runtests.jl index 6ba1a352dc..71909bcb69 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -18,6 +18,7 @@ include("test_utils.jl") # Unit tests include("test_collectionsofviews.jl") include("test_refshapes.jl") +include("test_orientations.jl") include("test_interpolations.jl") include("test_cellvalues.jl") include("test_facevalues.jl") diff --git a/test/test_orientations.jl b/test/test_orientations.jl new file mode 100644 index 0000000000..bede46b4c5 --- /dev/null +++ b/test/test_orientations.jl @@ -0,0 +1,75 @@ +@testset "Orientations" begin + @testset "PathOrientationInfo" begin + # Path [a, b], flipped + test_paths = [ + ((1, 2), false), # 1 -> 2, regular + ((2, 1), true), # 2 -> 1, inverted + ((5, 1), true), # 5 -> 1, inverted + ((3, 5), false), # 3 -> 5, regular + ] + + for (path, flipped) in test_paths + orientation = Ferrite.PathOrientationInfo(path) + + @test orientation.flipped == flipped + end + end + + @testset "SurfaceOrientationInfo" begin + # Surface [facenodes], flipped, shift_index + test_surfaces = [ + # 2 + # | \ + # 3 - 1, regular, rotation = 0 deg + ((1, 2, 3), false, 0), + # 1 + # | \ + # 2 - 3, regular, rotation = 120 deg + ((3, 1, 2), false, 1), + # 3 + # | \ + # 1 - 2, regular, rotation = 240 deg + ((2, 3, 1), false, 2), + # 3 + # | \ + # 2 - 1, inverted, rotation = 0 deg + ((1, 3, 2), true, 0), + # 1 + # | \ + # 3 - 2, inverted, rotation = 120 deg + ((2, 1, 3), true, 1), + # 2 + # | \ + # 1 - 3, inverted, rotation = 240 deg + ((3, 2, 1), true, 2), + # 4 - 3 + # | | + # 1 - 2, regular, rotation = 0 deg + ((1, 2, 3, 4), false, 0), + # 2 - 3 + # | | + # 1 - 4, inverted, rotation = 0 deg + ((1, 4, 3, 2), true, 0), + # 3 - 2 + # | | + # 4 - 1, regular, rotation = 90 deg + ((4, 1, 2, 3), false, 1), + # 1 - 2 + # | | + # 4 - 3, inverted, rotation = 270 deg + ((4, 3, 2, 1), true, 3), + ] + + + for (surface, flipped, shift_index) in test_surfaces + orientation = Ferrite.SurfaceOrientationInfo(surface) + + @test orientation.flipped == flipped + @test orientation.shift_index == shift_index + end + end + + @assert typeof(Ferrite.orientation_info((1, 2))) <: Ferrite.PathOrientationInfo + @assert typeof(Ferrite.orientation_info((1, 2, 3))) <: Ferrite.SurfaceOrientationInfo + @assert typeof(Ferrite.orientation_info((1, 2, 3, 4))) <: Ferrite.SurfaceOrientationInfo +end