Update of upstream patch version of TetGen and API improvements (#25)

Update of upstream patch version of TetGen and API improvements

* Update of upstream patch version of TetGen to "August 18, 2018". Formally this still appears to be 1.5.1 , and it is manifest in TetGen_jll 1.5.1+1
* Catch errors thrown by TetGen in a rc variable, throw julia error if rc!=0
* Handle unsuitable callback to allow local refinement.

Co-authored-by: Juergen Fuhrmann <[email protected]>

Author: j-fu

Project.toml

@@ -1,20 +1,21 @@
 name = "TetGen"
 uuid = "c5d3f3f7-f850-59f6-8a2e-ffc6dc1317ea"
 authors = ["SimonDanisch <[email protected]>", "Juergen Fuhrmann <[email protected]"]
-version = "1.0.0"
+version = "1.1.0"
 
 [deps]
 DocStringExtensions = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae"
 GeometryBasics = "5c1252a2-5f33-56bf-86c9-59e7332b4326"
-Libdl = "8f399da3-3557-5675-b5ff-fb832c97cbdb"
+LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
+Printf = "de0858da-6303-5e67-8744-51eddeeeb8d7"
 StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
 TetGen_jll = "b47fdcd6-d2c1-58e9-bbba-c1cee8d8c179"
 
 [compat]
 DocStringExtensions = "^0.8"
 GeometryBasics = "0.2, 0.3"
 StaticArrays = "0.11, 0.12"
-TetGen_jll = "^1.5.1"
+TetGen_jll = "1.5.1"
 julia = "1.3"
 
 [extras]

README.md

@@ -9,7 +9,7 @@
 
 The `TetGen.jl` package is a Julia wrapper for the C++ project [TetGen](https://wias-berlin.de/software/index.jsp?id=TetGen&lang=1). This wrapper enables TetGen based tetrahedral meshing, and (constrained) 3D Delaunay and Voronoi tesselation.
 
-## Example
+## Example using GeometryBasics datatypes
 
 ```julia
 using TetGen
@@ -51,6 +51,55 @@ Plotted with Makie:
 ![image](https://user-images.githubusercontent.com/1010467/82307971-69252000-99c1-11ea-8b82-e3a206381bd3.png)
 
 
+## Example using plain Julia arrays
+
+```julia
+using TetGen
+let
+    tetunsuitable() do pa,pb,pc,pd
+        vol=det(hcat(pb-pa,pc-pa,pd-pa))/6
+        center=0.25*(pa+pb+pc+pd)-[0.5,0.5,0.5]
+        vol> 0.05*norm(center)^2.5
+    end
+
+    input=TetGen.RawTetGenIO{Cdouble}()
+    input.pointlist=[0 0 0;  
+                     1 0 0;
+                     1 1 0;
+                     0 1 0;
+                     0 0 1;  
+                     1 0 1;
+                     1 1 1;
+                     0 1 1]'
+
+    TetGen.facetlist!(input,[1 2 3 4;
+                             5 6 7 8;
+                             1 2 6 5;
+                             2 3 7 6;
+                             3 4 8 7;
+                             4 1 5 8]')
+    tetrahedralize(input, "pQa")
+end
+```
+
+Output:
+
+```julia
+RawTetGenIO(
+numberofpoints=169,
+numberofedges=27,
+numberoftrifaces=112,
+numberoftetrahedra=809,
+pointlist'=[0.0 1.0 … 0.500059730245037 0.4996534431688176; 0.0 0.0 … 0.5074057466787957 0.49707528530503103; 0.0 0.0 … 0.5033015055704277 0.4953177845338027],
+tetrahedronlist'=Int32[34 47 … 15 143; 6 24 … 143 15; 58 52 … 154 150; 70 73 … 168 168],
+trifacelist'=Int32[3 58 … 99 22; 19 6 … 22 8; 78 70 … 158 158],
+trifacemarkerlist'=Int32[-1, -1, -1, -1, -1, -1, -1, -1, -1, -1  …  -1, -1, -1, -1, -1, -1, -1, -1, -1, -1],
+edgelist'=Int32[3 5 … 70 157; 18 24 … 6 32],
+edgemarkerlist'=Int32[-1, -1, -1, -1, -1, -1, -1, -1, -1, -1  …  -1, -1, -1, -1, -1, -1, -1, -1, -1, -1],
+)
+```
+
+
 ## [Contributing](https://github.com/JuliaGeometry/TetGen.jl/blob/master/CONTRIBUTING.md)   
 
 

examples/examples.jl

@@ -1,4 +1,5 @@
 using TetGen
+using LinearAlgebra
 
 """
      random_delaunay(;npoints=20)
@@ -38,6 +39,44 @@ function cube(;vol=1)
     tetrahedralize(input, "pQa$(vol)")
 end
 
+
+
+"""
+   cube_localref()
+
+   Tetrahedralization of cube with local refinement callback
+"""
+
+function cube_localref()
+
+    tetunsuitable() do pa,pb,pc,pd
+        vol=det(hcat(pb-pa,pc-pa,pd-pa))/6
+        center=0.25*(pa+pb+pc+pd)-[0.5,0.5,0.5]
+        vol> 0.05*norm(center)^2.5
+    end
+
+    input=TetGen.RawTetGenIO{Cdouble}()
+    input.pointlist=[0 0 0;  
+                     1 0 0;
+                     1 1 0;
+                     0 1 0;
+                     0 0 1;  
+                     1 0 1;
+                     1 1 1;
+                     0 1 1]'
+
+    TetGen.facetlist!(input,[1 2 3 4;
+                             5 6 7 8;
+                             1 2 6 5;
+                             2 3 7 6;
+                             3 4 8 7;
+                             4 1 5 8]')
+    tetrahedralize(input, "pQa")
+end
+
+
+
+
 """
    prism(;vol=1)
 

src/TetGen.jl

@@ -1,12 +1,16 @@
 module TetGen
 using DocStringExtensions
 using TetGen_jll
+
 using GeometryBasics
 using GeometryBasics: Polygon, MultiPolygon, Point, LineFace, Polytope, Line,
     Simplex, connect, Triangle, NSimplex, Tetrahedron,
     TupleView, TriangleFace, SimplexFace, LineString, Mesh, TetrahedronP, TriangleP,
     NgonFace, Ngon, faces, coordinates, metafree, meta, faces
 
+
+using Printf
+
 using StaticArrays
 
 
@@ -18,6 +22,8 @@ include("api.jl")
 
 
 export tetrahedralize
+export tetunsuitable
+export TetGenError
 export RawTetGenIO, facetlist!, RawFacet
 export numberofpoints,numberoftetrahedra,numberoftrifaces,numberofedges
 export volumemesh,surfacemesh

src/cpptetgenio.jl

@@ -57,5 +57,80 @@ struct CPPTetGenIO{T}
     numberofedges::Cint
 end
 
-tetrahedralize(input::CPPTetGenIO{Float64}, command::String)=ccall((:tetrahedralizef64, libtet), CPPTetGenIO{Float64}, (CPPTetGenIO{Float64}, Cstring), input, command)
+"""
+   Error struct for TetGen
+"""
+struct TetGenError <: Exception
+    rc::Cint
+end
+
+"""
+   Show TetGen error, messages have been lifted
+   from TetGen 
+"""
+function Base.show(io::IO, e::TetGenError)
+    if e.rc==1
+        println(io,"TetGen error $(e.rc): out of memory."); 
+    elseif e.rc==2
+        println(io,"TetGen error $(e.rc): internal error.")
+    elseif e.rc==3
+        println(io,"TetGen error $(e.rc): a self-intersection was detected. Hint: use -d option to detect all self-intersections."); 
+    elseif e.rc==4
+        println(io,"TetGen error $(e.rc): a very small input feature size was detected. Hint: use -T option to set a smaller tolerance.")
+    elseif e.rc==5
+        println(io,"TetGen error $(e.rc): two very close input facets were detected. Hint: use -Y option to avoid adding Steiner points in boundary.\n");
+    elseif e.rc==10
+        println(io,"TetGen error $(e.rc): an input error was detected.\n"); 
+    else
+        println(io,"TetGen error $(e.rc): unknown error.\n"); 
+    end    
+end
+
+
+"""
+$(SIGNATURES)
+
+Tetrahedralization with error handling
+"""
+function tetrahedralize(input::CPPTetGenIO{Float64}, command::String)
+    rc=Cint[0]
+    output = ccall((:tetrahedralize2_f64, libtet), CPPTetGenIO{Float64}, (CPPTetGenIO{Float64}, Cstring,  Ptr{Cint}), input, command, rc)
+    if rc[1]!=0
+        throw(TetGenError(rc[1]))
+    end
+    output
+end
 
+
+"""
+   Trivial Julia tetunsuitable function
+"""
+my_jl_tetunsuitable=(pa,pb,pc,pd)->0
+
+
+"""
+   Tetunsuitable function called from C wrapper
+"""
+function jl_wrap_tetunsuitable(pa::Ptr{Float64}, pb::Ptr{Float64}, pc::Ptr{Float64}, pd::Ptr{Float64})
+    pax=Base.unsafe_wrap(Array,pa,(3,),own=false)
+    pbx=Base.unsafe_wrap(Array,pb,(3,),own=false)
+    pcx=Base.unsafe_wrap(Array,pc,(3,),own=false)
+    pdx=Base.unsafe_wrap(Array,pd,(3,),own=false)
+    Cint(my_jl_tetunsuitable(pax,pbx,pcx,pdx))
+end
+
+
+"""
+   Set tetunsuitable function called from C wrapper.
+   Setting this function is valid only for one subsequent
+   call to tetrahedralize     
+"""
+function tetunsuitable(unsuitable::Function;check_signature=true)
+    if check_signature
+        unsuitable(rand(3),rand(3),rand(3),rand(3))
+    end
+    global my_jl_tetunsuitable
+    my_jl_tetunsuitable=unsuitable
+    c_wrap_tetunsuitable=@cfunction(jl_wrap_tetunsuitable, Cint, (Ptr{Cdouble},Ptr{Cdouble},Ptr{Cdouble},Ptr{Cdouble}))
+    ccall((:tetunsuitable_callback,libtet),Cvoid,(Ptr{Cvoid},),c_wrap_tetunsuitable)
+end

src/jltetgenio.jl

@@ -329,8 +329,11 @@ function Base.unsafe_convert(::Type{CPPTetGenIO{T}}, x::Tuple{CPPTetGenIO{T}, Ve
 end
 
 
-
 function tetrahedralize(input::JLTetGenIO{Float64}, command::String)
-    cres = ccall((:tetrahedralizef64, libtet), CPPTetGenIO{Float64}, (CPPTetGenIO{Float64}, Cstring), input, command)
+    rc=Cint[0]
+    cres = ccall((:tetrahedralize2_f64, libtet), CPPTetGenIO{Float64}, (CPPTetGenIO{Float64}, Cstring, Ptr{Cint}), input, command,rc)
+    if rc[1]!=0
+        throw(TetGenError(rc[1]))
+    end
     return convert(JLTetGenIO, cres)
 end

src/rawtetgenio.jl

@@ -419,8 +419,7 @@ function CPPTetGenIO(tio::RawTetGenIO{T}) where T
             edgemarkerlist=pointer(tio.edgemarkerlist)
         end            
     end
-    
-    
+
     # Create struct
     CPPTetGenIO{T}(firstnumber,
                    mesh_dim,
@@ -469,10 +468,10 @@ function RawTetGenIO(ctio::CPPTetGenIO{T}) where T
         tio.pointlist = convert(Array{T,2}, Base.unsafe_wrap(Array, ctio.pointlist, (3,Int(ctio.numberofpoints)), own=true))
     end
     if ctio.numberofpointattributes>0  && ctio.pointattributelist!=C_NULL
-        tio.pointattributelist=convert(Array{T,2}, Base.unsafe_wrap(Array, ctio.pointattributelistlist, (Int(ctio.numberofpointattributes),Int(ctio.numberofpoints)), own=true))
+        tio.pointattributelist=convert(Array{T,2}, Base.unsafe_wrap(Array, ctio.pointattributelist, (Int(ctio.numberofpointattributes),Int(ctio.numberofpoints)), own=true))
     end
     if ctio.numberofpointmtrs>0  && ctio.pointmtrlist!=C_NULL
-        tio.pointmtrlist=convert(Array{T,2}, Base.unsafe_wrap(Array, ctio.pointmtrlistlist, (Int(ctio.numberofpointmtrs),Int(ctio.numberofpoints)), own=true))
+        tio.pointmtrlist=convert(Array{T,2}, Base.unsafe_wrap(Array, ctio.pointmtrlist, (Int(ctio.numberofpointmtrs),Int(ctio.numberofpoints)), own=true))
     end
     if ctio.pointmarkerlist!=C_NULL
         tio.pointmarkerlist=convert(Array{Cint,1}, Base.unsafe_wrap(Array, ctio.pointmarkerlist, (Int(ctio.numberofpoints)), own=true))
@@ -626,7 +625,12 @@ Tetrahedralize input.
 """
 function tetrahedralize(input::RawTetGenIO{Float64}, flags::String)
     cinput,flist,plist=CPPTetGenIO(input)
-    coutput = ccall((:tetrahedralizef64, libtet), CPPTetGenIO{Float64}, (CPPTetGenIO{Float64}, Cstring), cinput, flags)
+    rc=Cint[0]
+    coutput = ccall((:tetrahedralize2_f64, libtet), CPPTetGenIO{Float64}, (CPPTetGenIO{Float64}, Cstring, Ptr{Cint}),
+                    cinput, flags, rc)
+    if rc[1]!=0
+        throw(TetGenError(rc[1]))
+    end
     RawTetGenIO(coutput)
 end