return Out, rad -- (center [Point], R [float])
end
+-- compute the list of the circumcircle of a triangulation
+function listCircumCenter(listPoints,triangulation)
+ list = {}
+ for j=1,#triangulation do
+ A = listPoints[triangulation[j][1]]
+ B = listPoints[triangulation[j][2]]
+ C = listPoints[triangulation[j][3]]
+ center, radius = circoncircle(A,B,C)
+ table.insert(list,{x=center.x,y=center.y,r=radius})
+ end
+ return list
+end
+
+-- find the three neighbour triangles of T
+function findNeighbour(T,i,triangulation)
+ -- T : triangle
+ -- i : index of T in triangualation
+ -- triangulation
+
+ list = {}
+ -- define the three edge
+ e1 = {T[1],T[2]}
+ e2 = {T[2],T[3]}
+ e3 = {T[3],T[1]}
+ for j=1,#triangulation do
+ if j~= i then
+ if(edgeInTriangle(e1,triangulation[j])) then
+ table.insert(list,j)
+ end
+ if(edgeInTriangle(e2,triangulation[j])) then
+ table.insert(list,j)
+ end
+ if(edgeInTriangle(e3,triangulation[j])) then
+ table.insert(list,j)
+ end
+ end
+ end
+ return list
+end
+
+-- test if edge are the same (reverse)
+function equalEdge(e1,e2)
+ if(((e1[1] == e2[1]) and (e1[2] == e2[2])) or ((e1[1] == e2[2]) and (e1[2] == e2[1]))) then
+ return true
+ else
+ return false
+ end
+end
+
+-- test if the edge belongs to the list
+function edgeInList(e,listE)
+ output = false
+ for i=1,#listE do
+ if(equalEdge(e,listE[i])) then
+ output = true
+ end
+ end
+ return output
+end
+
+-- build the edges of the Voronoi diagram with a given triangulation
+function buildVoronoi(listPoints, triangulation)
+ listCircumCircle = listCircumCenter(listPoints, triangulation)
+ listVoronoi = {}
+ for i=1,#listCircumCircle do
+ listN = findNeighbour(triangulation[i],i,triangulation)
+ for j=1,#listN do
+ edge = {i,listN[j]}
+ if( not edgeInList(edge, listVoronoi)) then
+ table.insert(listVoronoi, edge)
+ end
+ end
+ end
+ return listVoronoi
+end
-------------------------- TeX
-- build the list of points
return listPoints
end
+
+-- trace Voronoi with MP
+function traceVoronoiMP(listPoints, triangulation,listVoronoi, points, tri,styleD,styleV)
+ if(styleD == "dashed") then
+ sDelaunay = "dashed evenly"
+ else
+ sDelaunay = ""
+ end
+ if(styleV == "dashed") then
+ sVoronoi = "dashed evenly"
+ else
+ sVoronoi = ""
+ end
+ listCircumC = listCircumCenter(listPoints,triangulation)
+ output = "";
+ output = output .. " pair MeshPoints[];"
+ for i=1,#listPoints do
+ output = output .. "MeshPoints[".. i .. "] = (" .. listPoints[i].x .. "," .. listPoints[i].y .. ")*u;"
+ end
+ output = output .. " pair CircumCenters[];"
+ for i=1,#listCircumC do
+ output = output .. "CircumCenters[".. i .. "] = (" .. listCircumC[i].x .. "," .. listCircumC[i].y .. ")*u;"
+ end
+ if(tri=="show") then
+ for i=1,#triangulation do
+ PointI = listPoints[triangulation[i][1]]
+ PointJ = listPoints[triangulation[i][2]]
+ PointK = listPoints[triangulation[i][3]]
+ if(triangulation[i].type == "bbox") then
+ output = output .. "draw (".. PointI.x ..",".. PointI.y ..")*u--("..PointJ.x..",".. PointJ.y ..")*u--("..PointK.x..",".. PointK.y ..")*u--cycle "..sDelaunay.." withcolor \\luameshmpcolorBbox;"
+ else
+ output = output .. "draw (".. PointI.x ..",".. PointI.y ..")*u--("..PointJ.x..",".. PointJ.y ..")*u--("..PointK.x..",".. PointK.y ..")*u--cycle "..sDelaunay.." withcolor \\luameshmpcolor;"
+ end
+ end
+ end
+ for i=1,#listVoronoi do
+ PointI = listCircumC[listVoronoi[i][1]]
+ PointJ = listCircumC[listVoronoi[i][2]]
+ output = output .. "draw (".. PointI.x ..",".. PointI.y ..")*u--("..PointJ.x..",".. PointJ.y ..")*u "..sVoronoi.." withcolor \\luameshmpcolorVoronoi;"
+ end
+ if(points=="points") then
+ j=1
+ for i=1,#listPoints do
+ if(listPoints[i].type == "bbox") then
+ output = output .. "dotlabel.llft (btex $\\MeshPoint^{*}_{"..j.."}$ etex, (" .. listPoints[i].x ..",".. listPoints[i].y .. ")*u ) withcolor \\luameshmpcolorBbox ;"
+ j=j+1
+ else
+ output = output .. "dotlabel.llft (btex $\\MeshPoint_{" .. i .. "}$ etex, (" .. listPoints[i].x ..",".. listPoints[i].y .. ")*u ) withcolor \\luameshmpcolor ;"
+ end
+ end
+ for i=1,#listCircumC do
+ output = output .. "dotlabel.llft (btex $\\CircumPoint_{" .. i .. "}$ etex, (" .. listCircumC[i].x ..",".. listCircumC[i].y .. ")*u ) withcolor \\luameshmpcolorVoronoi ;"
+ end
+ else
+ j=1
+ for i=1,#listPoints do
+ if(listPoints[i].type == "bbox") then
+ output = output .. "drawdot (" .. listPoints[i].x ..",".. listPoints[i].y .. ")*u withcolor \\luameshmpcolorBbox withpen pencircle scaled 3;"
+ j=j+1
+ else
+ output = output .. "drawdot (" .. listPoints[i].x ..",".. listPoints[i].y .. ")*u withcolor \\luameshmpcolor withpen pencircle scaled 3;"
+ end
+ end
+ for i=1,#listCircumC do
+ output = output .. "drawdot (" .. listCircumC[i].x ..",".. listCircumC[i].y .. ")*u withcolor \\luameshmpcolorVoronoi withpen pencircle scaled 3;"
+ end
+ end
+
+ return output
+end
+
+
+-- trace Voronoi with TikZ
+function traceVoronoiTikZ(listPoints, triangulation,listVoronoi, points, tri,color,colorBbox,colorVoronoi,styleD,styleV)
+ if(styleD == "dashed") then
+ sDelaunay = ",dashed"
+ else
+ sDelaunay = ""
+ end
+ if(styleV == "dashed") then
+ sVoronoi = ",dashed"
+ else
+ sVoronoi = ""
+ end
+ listCircumC = listCircumCenter(listPoints,triangulation)
+ output = ""
+ for i=1,#listPoints do
+ output = output .. "\\coordinate (MeshPoints".. i .. ") at (" .. listPoints[i].x .. "," .. listPoints[i].y .. ");"
+ end
+ for i=1,#listCircumC do
+ output = output .. "\\coordinate (CircumPoints".. i .. ") at (" .. listCircumC[i].x .. "," .. listCircumC[i].y .. ");"
+ end
+ if(tri=="show") then
+ for i=1,#triangulation do
+ PointI = listPoints[triangulation[i][1]]
+ PointJ = listPoints[triangulation[i][2]]
+ PointK = listPoints[triangulation[i][3]]
+ if(triangulation[i].type == "bbox") then
+ output = output .. "\\draw[color="..colorBbox..sDelaunay.."] (".. PointI.x ..",".. PointI.y ..")--("..PointJ.x..",".. PointJ.y ..")--("..PointK.x..",".. PointK.y ..")--cycle;"
+ else
+ output = output .. "\\draw[color="..color..sDelaunay.."] (".. PointI.x ..",".. PointI.y ..")--("..PointJ.x..",".. PointJ.y ..")--("..PointK.x..",".. PointK.y ..")--cycle;"
+ end
+ end
+ end
+ for i=1,#listVoronoi do
+ PointI = listCircumC[listVoronoi[i][1]]
+ PointJ = listCircumC[listVoronoi[i][2]]
+ output = output .. "\\draw[color="..colorVoronoi..sVoronoi.."] (".. PointI.x ..",".. PointI.y ..")--("..PointJ.x..",".. PointJ.y ..");"
+ end
+ if(points=="points") then
+ j=1
+ for i=1,#listPoints do
+ if(listPoints[i].type == "bbox") then
+ output = output .. "\\draw[color="..colorBbox.."] (" .. listPoints[i].x ..",".. listPoints[i].y .. ") node {$\\bullet$} node[anchor=north east] {$\\MeshPoint^*_{" .. j .. "}$};"
+ j=j+1
+ else
+ output = output .. "\\draw[color="..color.."] (" .. listPoints[i].x ..",".. listPoints[i].y .. ") node {$\\bullet$} node[anchor=north east] {$\\MeshPoint_{" .. i .. "}$};"
+ end
+ end
+ for i=1,#listCircumC do
+ output = output .. "\\draw[color="..colorVoronoi.."] (" .. listCircumC[i].x ..",".. listCircumC[i].y .. ") node {$\\bullet$} node[anchor=north east] {$\\CircumPoint_{" .. i .. "}$};"
+ end
+ else
+ j=1
+ for i=1,#listPoints do
+ if(listPoints[i].type == "bbox") then
+ output = output .. "\\draw[color="..colorBbox.."] (" .. listPoints[i].x ..",".. listPoints[i].y .. ") node {$\\bullet$} ;"
+ j=j+1
+ else
+ output = output .. "\\draw[color="..color.."] (" .. listPoints[i].x ..",".. listPoints[i].y .. ") node {$\\bullet$} ;"
+ end
+ end
+ for i=1,#listCircumC do
+ output = output .. "\\draw[color="..colorVoronoi.."] (" .. listCircumC[i].x ..",".. listCircumC[i].y .. ") node {$\\bullet$};"
+ end
+ end
+ return output
+end
+
+
+
+-- buildVoronoi with MP
+function buildVoronoiMPBW(chaine,mode,points,bbox,scale,tri,styleD,styleV)
+ listPoints = buildList(chaine, mode)
+ triangulation = BowyerWatson(listPoints,bbox)
+ listVoronoi = buildVoronoi(listPoints, triangulation)
+ output = traceVoronoiMP(listPoints,triangulation,listVoronoi,points,tri,styleD,styleV)
+ output = "\\leavevmode\\begin{mplibcode}beginfig(0);u:="..scale.. ";" .. output .."endfig;\\end{mplibcode}"
+ tex.sprint(output)
+end
+
+
+-- buildVoronoi with TikZ
+function buildVoronoiTikZBW(chaine,mode,points,bbox,scale,tri,color,colorBbox,colorVoronoi,styleD,styleV)
+ listPoints = buildList(chaine, mode)
+ triangulation = BowyerWatson(listPoints,bbox)
+ listVoronoi = buildVoronoi(listPoints, triangulation)
+ output = traceVoronoiTikZ(listPoints,triangulation,listVoronoi,points,tri,color,colorBbox,colorVoronoi,styleD,styleV)
+ output = "\\noindent\\begin{tikzpicture}[x=" .. scale .. ",y=" .. scale .."]" .. output .."\\end{tikzpicture}" tex.sprint(output)
+end
+
+
+-- buildVoronoi with MP
+function buildVoronoiMPBWinc(chaine,beginning, ending,mode,points,bbox,scale,tri,styleD,styleV)
+ listPoints = buildList(chaine, mode)
+ triangulation = BowyerWatson(listPoints,bbox)
+ listVoronoi = buildVoronoi(listPoints, triangulation)
+ output = traceVoronoiMP(listPoints,triangulation,listVoronoi,points,tri,styleD,styleV)
+ output = "\\leavevmode\\begin{mplibcode}u:="..scale..";"..beginning .. output .. ending .. "\\end{mplibcode}"
+ tex.sprint(output)
+end
+
+
+-- buildVoronoi with TikZ
+function buildVoronoiTikZBWinc(chaine,beginning, ending,mode,points,bbox,scale,tri,color,colorBbox,colorVoronoi)
+ listPoints = buildList(chaine, mode,styleD,styleV)
+ triangulation = BowyerWatson(listPoints,bbox)
+ listVoronoi = buildVoronoi(listPoints, triangulation)
+ output = traceVoronoiTikZ(listPoints,triangulation,listVoronoi,points,tri,color,colorBbox,colorVoronoi,styleD,styleV)
+ output = "\\noindent\\begin{tikzpicture}[x=" .. scale .. ",y=" .. scale .."]" ..beginning.. output..ending .."\\end{tikzpicture}"
+ tex.sprint(output)
+end
+
+
+
-- trace a triangulation with TikZ
function traceMeshTikZ(listPoints, triangulation,points,color,colorBbox)
output = ""
+ for i=1,#listPoints do
+ output = output .. "\\coordinate (MeshPoints".. i .. ") at (" .. listPoints[i].x .. "," .. listPoints[i].y .. ");"
+ end
for i=1,#triangulation do
PointI = listPoints[triangulation[i][1]]
PointJ = listPoints[triangulation[i][2]]
for i=1,#listPoints do
output = output .. "MeshPoints[".. i .. "] = (" .. listPoints[i].x .. "," .. listPoints[i].y .. ")*u;"
end
-
for i=1,#triangulation do
PointI = listPoints[triangulation[i][1]]
PointJ = listPoints[triangulation[i][2]]
-- print points of the mesh
function tracePointsTikZ(listPoints,points,color,colorBbox)
output = "";
+ for i=1,#listPoints do
+ output = output .. "\\coordinate (MeshPoints".. i .. ") at (" .. listPoints[i].x .. "," .. listPoints[i].y .. ");"
+ end
if(points=="points") then
j=1
for i=1,#listPoints do
-- build the triangulation
triangulation = BowyerWatson(listPoints,bbox)
badTriangles = buildBadTriangles(P,triangulation)
+ for i=1,#listPoints do
+ output = output .. "\\coordinate (MeshPoints".. i .. ") at (" .. listPoints[i].x .. "," .. listPoints[i].y .. ");"
+ end
if(step == "badT") then
-- draw all triangle
for i=1,#triangulation do
output = "\\begin{mplibcode}u:="..scale..";"..beginning .. output .. ending .. "\\end{mplibcode}"
tex.sprint(output)
end
+
+
+function split(pString, pPattern)
+ local Table = {} -- NOTE: use {n = 0} in Lua-5.0
+ local fpat = "(.-)" .. pPattern
+ local last_end = 1
+ local s, e, cap = pString:find(fpat, 1)
+ while s do
+ if s ~= 1 or cap ~= "" then
+ table.insert(Table,cap)
+ end
+ last_end = e+1
+ s, e, cap = pString:find(fpat, last_end)
+ end
+ if last_end <= #pString then
+ cap = pString:sub(last_end)
+ table.insert(Table, cap)
+ end
+ return Table
+end
+
+function readGmsh(file)
+ io.input(file) -- open the file
+ text=io.read("*all")
+ local lines = split(text,"\n+") -- all the lines
+ listPoints={}
+ triangulation ={}
+ boolNodes = false
+ Jnodes = 0
+ boolElements = false
+ Jelements = 0
+ J=0
+ for i=1,#lines-J do
+ if(lines[i+J] == "$EndNodes") then
+ boolNodes = false
+ -- go to the next line
+ end
+ if(boolNodes) then -- we are in the Nodes environment
+ xy=split(lines[i+J]," +")
+ table.insert(listPoints,{x=tonumber(xy[2]),y=tonumber(xy[3])})
+ end
+ if(lines[i+J] == "$Nodes") then
+ boolNodes = true
+ -- go to the next line
+ J=J+1
+ end
+ if(lines[i+J] == "$EndElements") then
+ boolElements = false
+ -- go to the next line
+ end
+ if(boolElements) then -- we are in the Nodes environment
+ xy=split(lines[i+J]," +")
+ if(tonumber(xy[2]) == 2) then -- if the element is a triangle
+ nbrTags = xy[3]+1
+ table.insert(triangulation,{tonumber(xy[2+nbrTags+1]),tonumber(xy[2+nbrTags+2]),tonumber(xy[2+nbrTags+3])})
+ end
+ end
+ if(lines[i+J] == "$Elements") then
+ boolElements = true
+ -- go to the next line
+ J=J+1
+ end
+ end
+ return listPoints, triangulation
+end
+
+function drawGmshMP(file,points,scale)
+ listPoints,triangulation = readGmsh(file)
+ output = traceMeshMP(listPoints,triangulation,points)
+ output = "\\leavevmode\\begin{mplibcode}beginfig(0);u:="..scale.. ";" .. output .."endfig;\\end{mplibcode}"
+ tex.sprint(output)
+end
+
+function drawGmshMPinc(file,beginning,ending,points,scale)
+ listPoints,triangulation = readGmsh(file)
+ output = traceMeshMP(listPoints,triangulation,points)
+ output = "\\leavevmode\\begin{mplibcode}u:="..scale..";"..beginning .. output .. ending .. "\\end{mplibcode}"
+ tex.sprint(output)
+end
+
+
+
+--
+function drawGmshTikZ(file,points,scale,color)
+ listPoints,triangulation = readGmsh(file)
+ output = traceMeshTikZ(listPoints, triangulation,points,color,colorBbox)
+ output = "\\noindent\\begin{tikzpicture}[x=" .. scale .. ",y=" .. scale .."]" .. output .."\\end{tikzpicture}"
+ tex.sprint(output)
+end
+
+--
+function drawGmshTikZinc(file,beginning, ending,points,scale,color)
+ listPoints,triangulation = readGmsh(file)
+ output = traceMeshTikZ(listPoints, triangulation,points,color,colorBbox)
+ output = "\\noindent\\begin{tikzpicture}[x=" .. scale .. ",y=" .. scale .."]" ..beginning.. output..ending .."\\end{tikzpicture}"
+ tex.sprint(output)
+end
+
+
+-- buildVoronoi with MP
+function gmshVoronoiMP(file,points,scale,tri,styleD,styleV)
+ listPoints,triangulation = readGmsh(file)
+ listVoronoi = buildVoronoi(listPoints, triangulation)
+ output = traceVoronoiMP(listPoints,triangulation,listVoronoi,points,tri,styleD,styleV)
+ output = "\\leavevmode\\begin{mplibcode}beginfig(0);u:="..scale.. ";" .. output .."endfig;\\end{mplibcode}"
+ tex.sprint(output)
+end
+
+
+-- buildVoronoi with TikZ
+function gmshVoronoiTikZ(file,points,scale,tri,color,colorVoronoi,styleD,styleV)
+ listPoints,triangulation = readGmsh(file)
+ listVoronoi = buildVoronoi(listPoints, triangulation)
+ output = traceVoronoiTikZ(listPoints,triangulation,listVoronoi,points,tri,color,colorBbox,colorVoronoi,styleD,styleV)
+ output = "\\noindent\\begin{tikzpicture}[x=" .. scale .. ",y=" .. scale .."]" .. output .."\\end{tikzpicture}" tex.sprint(output)
+end
+
+
+-- buildVoronoi with MP
+function gmshVoronoiMPinc(file,beginning, ending,points,scale,tri,styleD,styleV)
+ listPoints,triangulation = readGmsh(file)
+ listVoronoi = buildVoronoi(listPoints, triangulation)
+ output = traceVoronoiMP(listPoints,triangulation,listVoronoi,points,tri,styleD,styleV)
+ output = "\\leavevmode\\begin{mplibcode}u:="..scale..";"..beginning .. output .. ending .. "\\end{mplibcode}"
+ tex.sprint(output)
+end
+
+
+-- buildVoronoi with TikZ
+function gmshVoronoiTikZinc(file,beginning, ending,points,scale,tri,color,colorVoronoi,styleD,styleV)
+ listPoints,triangulation = readGmsh(file)
+ listVoronoi = buildVoronoi(listPoints, triangulation)
+ output = traceVoronoiTikZ(listPoints,triangulation,listVoronoi,points,tri,color,colorBbox,colorVoronoi,styleD,styleV)
+ output = "\\noindent\\begin{tikzpicture}[x=" .. scale .. ",y=" .. scale .."]" ..beginning.. output..ending .."\\end{tikzpicture}"
+ tex.sprint(output)
+end