Improvements udim2

src/CompScienceMeshes.jl

@@ -38,6 +38,7 @@ export connectivity, cellpairs # marked for deprecation
 export translate, translate!, rotate, rotate!
 export fliporientation!, fliporientation
 export weld, union
+export permutate
 
 # mesh refinement
 export barycentric_refinement, bisecting_refinement
@@ -142,6 +143,7 @@ include("primitives/primitives.jl")
 include("baryref.jl")
 include("subdivision.jl")
 include("weld.jl")
+include("permutatie_mesh.jl")
 
 include("mapper.jl")
 include("restrict.jl")

src/charts.jl

@@ -217,7 +217,7 @@ function _normals(tangents::SVector{2,SVector{2,T}}, ::Type{Val{0}}) where {T}
 
     t = tangents[1]
     s = tangents[2]
-    v = (t[1]*s[2] - t[2]*s[1])/2
+    v = abs(t[1]*s[2] - t[2]*s[1])/2
     # n[3] = tangents[1] × tangents[2]
     # l = norm(n)
 

src/mesh.jl

@@ -331,7 +331,7 @@ Base.getindex(m::AbstractMesh, I::Vector{Int}) = Mesh(vertices(m), cells(m)[I])
 
 Returns the boundary of `mesh` as a mesh of lower dimension.
 """
-function boundary(mesh)
+function boundary(mesh, include_interior_points=true)
 
     D = dimension(mesh)
 
@@ -366,7 +366,22 @@ function boundary(mesh)
     end
 
     resize!(bnd_edges, i-1)
-    bnd = Mesh(vertices(mesh), bnd_edges)
+    if include_interior_points
+        return Mesh(vertices(mesh), bnd_edges)
+    else
+        originaltonew = Dict{Int,Int}()
+        sizehint!(originaltonew, i-1)
+        bnd_edges_new = Vector{I}(undef, length(bnd_edges))
+        for (e,edge) in enumerate(bnd_edges)
+            get!(originaltonew, edge.indices[1], e)
+        end
+        for (e,edge) in enumerate(bnd_edges)
+            bnd_edges_new[e] =SimplexGraph(
+                get(originaltonew, edge.indices[1], -1), 
+                get(originaltonew, edge.indices[2], -1)) 
+        end
+        return Mesh(vertices(mesh)[collect(keys(originaltonew))], bnd_edges_new)
+    end
 end
 
 

src/overlap.jl

@@ -106,3 +106,34 @@ function overlap(p::Simplex{3,3,0,4,T}, q::Simplex{3,3,0,4,T}) where T
     all(0+tol .<= u .<= 1-tol) && return true
     return false
 end
+
+
+function overlap(p::Simplex{2,2,0,3,T}, q::Simplex{2,2,0,3,T}) where T
+
+    #   tol = sqrt(eps(T))
+    tol = 1e3 * eps(T)
+
+  # Are the patches in the same plane?
+  u1 = q.tangents[1]
+  u2 = q.tangents[2]
+  v = p.vertices[1] - q.vertices[2]
+
+  for i in 1:3
+    a = p.vertices[mod1(i+1,3)]
+    b = p.vertices[mod1(i+2,3)]
+    c = p.vertices[i]
+    t = b - a
+    m = StaticArrays.SVector{2,T}(t[2],-t[1])
+
+    sp = zeros(T,3); sp[i] = dot(c-a, m)
+    sq = T[dot(q.vertices[j]-a, m) for j in 1:3]
+
+    minp, maxp = extrema(sp)
+    minq, maxq = extrema(sq)
+
+    maxq <= minp + tol && return false
+    maxp <= minq + tol && return false
+  end
+
+  return true
+end

src/permutatie_mesh.jl

@@ -0,0 +1,113 @@
+"""
+permutate_mesh permutate the vertices of a mesh, while keeping the same cells.
+
+Permutation is represented as a Vector v which sets the v[i]-th vertex at the i-th place.
+
+mesh can be permutated from:
+- Vector or SVector
+- Permutations.Permutation
+- another mesh with same universedimension
+- another mesh with different universedimension, only works for combination of udims 2,3.
+
+In the last two, control that the other mesh is a submesh of the starting mesh. In either direct
+
+
+"""
+
+function permutate(mesh::Mesh, σ::Permutations.Permutation)
+    mesh.vertices = mesh.vertices[σ.data]
+    for (i, a) in enumerate(mesh.faces)
+        ind = σ'.data[a.indices]
+        mesh.faces[i] = SimplexGraph(ind...)
+    end
+end
+
+function permutate(mesh::Mesh, σ::Union{Vector{<:Integer}, StaticArrays.SVector{Int32, <:Integer}})
+    @assert numvertices(mesh) == length(σ)
+    permutate(mesh, Permutations.Permutation(σ))
+end
+
+function permutate(X::Mesh{U}, Y::Mesh{U}) where {U}
+    tol = sqrt(eps(coordtype(X)))
+    permut = Vector{Int32}()
+    temp = collect(1:numvertices(X))
+    for p in Y.vertices
+        index = findfirst(isapprox(p;atol = tol), X.vertices)
+        @assert !isnothing(index)
+
+        temp[index] = 0
+        append!(permut, index)
+    end
+
+    for i in temp
+        if i !=0
+            push!(permut, i)
+    end end
+    permutate(X, Permutations.Permutation(permut))
+    X
+end
+
+function permutate(X::Mesh{2}, Y::Mesh{3})
+    tol = sqrt(eps(coordtype(X)))
+
+    # assert that the z-coordinate for the vertices in X is constant.
+
+    permut = Vector{Int32}()  
+    temp = collect(1:numvertices(X))
+
+    for p in Y.vertices
+        q = zeros(eltype(X.vertices[begin]), 2)
+        q[begin] = p[begin]
+        q[2] = p[2]
+
+
+        index = findfirst(isapprox(q;atol = tol), X.vertices)
+        @assert !isnothing(index)   # vertex from Y exist in X
+        @assert temp[index] != 0    # vertex from Y unique in X
+
+        temp[index] = 0
+        append!(permut, index)
+    end
+
+    for i in temp
+        if i !=0
+            push!(permut, i)
+    end end
+
+    permutate(X, Permutations.Permutation(permut))
+    X
+end
+
+function permutate(X::Mesh{3}, Y::Mesh{2})
+    tol = sqrt(eps(coordtype(X)))
+
+    # assert that the z-coordinate for the vertices in X is constant.
+    for v in X.vertices
+        @assert v[3] == X.vertices[1][3]
+    end
+
+    permut = Vector{Int32}()  
+    temp = collect(1:numvertices(X))
+
+    for p in Y.vertices
+        q = zeros(eltype(X.vertices[begin]), 3)
+        q[begin] = p[begin]
+        q[2] = p[2]
+        q[3] = X.vertices[1][3]
+
+        index = findfirst(isapprox(q;atol = tol), X.vertices)
+        @assert !isnothing(index)   # vertex from Y exist in X
+        @assert temp[index] != 0    # vertex from Y unique in X
+
+        temp[index] = 0
+        append!(permut, index)
+    end
+
+    for i in temp
+        if i !=0
+            push!(permut, i)
+    end end
+
+    permutate(X, Permutations.Permutation(permut))
+    X
+end

src/plotlyjs_glue.jl

@@ -44,6 +44,44 @@ function __init__()
             return s
         end
 
+        @eval function patch(Γ::CompScienceMeshes.AbstractMesh{2}, fcr=nothing;
+            caxis=nothing, showscale=true, color="red", kwargs...)
+
+            v = vertexarray(Γ)
+            c = cellarray(Γ)
+
+            x = v[:,1];    y = v[:,2];    z = zeros(length(y))
+            i = c[:,1].-1; j = c[:,2].-1; k = c[:,3].-1
+
+
+            if fcr != nothing
+                m, M = extrema(fcr)
+                if caxis != nothing
+                    m, M = caxis
+                end
+
+                s = PlotlyBase.mesh3d(;
+                    x=x, y=y, z=z,
+                    i=i, j=j, k=k,    
+                    intensitymode="cell",
+                    intensity=fcr,
+                    colorscale="Viridis",
+                    showscale=showscale,
+                    cmin=m,
+                    cmax=M,
+                    kwargs...
+                )
+            else
+                s = PlotlyBase.mesh3d(;
+                    x=x, y=y, z=z,
+                    i=i, j=j, k=k,
+                    color=color,
+                    kwargs...
+                )
+            end
+            return s
+        end
+
         @eval function patch(a::Vector{<:Simplex}; kwargs...)
             vertices = reduce(vcat, [v.vertices for v in a])
             faces = collect(SVector(3*(i-1)+1, 3*(i-1)+2, 3*(i-1)+3) for i in 1:length(a))

test/runtests.jl

@@ -32,12 +32,14 @@ include("test_geometry.jl")
 include("test_patches.jl")
 include("test_submesh.jl")
 include("test_overlap.jl")
+include("test_overlap2D.jl")
 include("test_intersect.jl")
 include("test_sh_intersection.jl")
 include("test_isinside.jl")
 include("test_jctweld.jl")
 include("test_isinclosure.jl")
 include("test_permute_vertices.jl")
+include("test_permutate_mesh.jl")
 
 include("test_trgauss.jl")
 include("test_trdunavant.jl")

test/test_mesh.jl

@@ -82,4 +82,23 @@ T = Float64
     Σᵀ = connectivity(edges, faces, identity)
     Σ = connectivity(faces, edges, identity)
     @test norm(Σᵀ - Σ', Inf) == 0
+
+    rectangle = meshrectangle(T(1.0), T(1.0), T(0.5));
+    bnd = boundary(rectangle, true)
+
+    @test dimension(bnd) == 1
+    @test numvertices(bnd) == 9
+    @test numcells(bnd) == 8
+
+    bnd = boundary(rectangle, false)
+
+    @test dimension(bnd) == 1
+    @test numvertices(bnd) == 8
+    @test numcells(bnd) == 8
+
+
+    m = meshrectangle(T(1.0),T(1.0),T(0.5))
+    CompScienceMeshes.rotate!(m, [1,0,0]*T(π/2))
+
+    @test CompScienceMeshes.isoriented(boundary(m, false))
 # end

test/test_overlap2D.jl

@@ -0,0 +1,32 @@
+using Test
+using CompScienceMeshes
+using StaticArrays
+
+p = simplex(
+    point(0,0),
+    point(1,0),
+    point(0,1))
+
+q1 = simplex(
+    point(0.6, 0.6),
+    point(1.6, 0.6),
+    point(0.6, 1.6),
+)
+
+q2 = simplex(
+    point(0.4, 0.4),
+    point(1.4, 0.4),
+    point(0.4, 1.4),
+)
+
+@test overlap(p, q1) == false
+@test overlap(p, q2) == true
+
+
+# test a case where the segments are:
+#   not of unit length
+#   colinear and opposite
+#   meet in a common point
+ch1 = simplex(point(1/3,0), point(1/3,1/3))
+ch2 = simplex(point(1/3,1/3), point(1/3,2/3))
+@test !overlap(ch1, ch2)

test/test_patches.jl

@@ -1,6 +1,7 @@
 using Test
 using CompScienceMeshes
 
+# universe dimension 3
 for T in [Float32, Float64]
     local mesh = meshrectangle(T(1.0), T(1.0), T(1.0))
     # local faces = skeleton(mesh, 2)
@@ -19,3 +20,22 @@ for T in [Float32, Float64]
         point(T, 0.0, 0.0, -1.0)]
     @test p.volume == T(0.5)
 end
+
+# universe dimension 2
+for T in [Float32, Float64]
+    local mesh = meshrectangle(T(1.0), T(1.0), T(1.0),2)
+    # local faces = skeleton(mesh, 2)
+    # local verts = cellvertices(mesh, 1)
+    p = chart(mesh, cells(mesh)[1])
+    #p = simplex(verts)
+
+    @test p.vertices == [
+        point(T, 0.0, 0.0),
+        point(T, 0.0, 1.0),
+        point(T, 1.0, 0.0)]
+    @test p.tangents == [
+        point(T, -1.0, 0.0),
+        point(T, -1.0, 1.0)]
+    @test p.normals == []
+    @test p.volume == T(0.5)
+end