Point attributes working again,

Added RCM node renumbering, 
Some changes to mesh properties

git-svn-id: https://triangle.svn.codeplex.com/svn@68848 0e2699bc-83d4-4a8f-98e7-55e24ab8c7a5

Triangle.NET/TestApp/Controls/AngleHistogram.cs

@@ -95,6 +95,7 @@ namespace MeshExplorer.Controls
             if (maxAngleCount == 0)
             {
                 this.maxAngles = null;
+                return;
             }
 
             this.Invalidate();

Triangle.NET/TestApp/Controls/DarkCheckBox.cs

@@ -82,6 +82,37 @@ namespace MeshExplorer.Controls
             InitializeComponent();
         }
 
+        private void DrawText(Graphics g, Color forecolor, Point location)
+        {
+            if (this.UseCompatibleTextRendering)
+            {
+                //using (StringFormat stringFormat = this.CreateStringFormat())
+                {
+                    if (this.Enabled)
+                    {
+                        g.DrawString(this.Text, base.Font, new SolidBrush(forecolor), location.X, location.Y);
+                    }
+                    else
+                    {
+                        g.DrawString(this.Text, base.Font, new SolidBrush(forecolor), location.X, location.Y);
+                    }
+                }
+            }
+            else
+            {
+                //TextFormatFlags textFormatFlags = this.CreateTextFormatFlags();
+                if (this.Enabled)
+                {
+                    TextRenderer.DrawText((IDeviceContext)g, this.Text, this.Font, location, forecolor);
+                }
+                else
+                {
+                    //forecolor = TextRenderer.DisabledTextColor(this.BackColor);
+                    TextRenderer.DrawText((IDeviceContext)g, this.Text, this.Font, location, forecolor);
+                }
+            }
+        }
+
         #region Control overrides
 
         protected override void OnPaint(PaintEventArgs e)
@@ -153,7 +184,7 @@ namespace MeshExplorer.Controls
             e.Graphics.TextRenderingHint = TextRenderingHint.ClearTypeGridFit;
 
             SizeF szL = e.Graphics.MeasureString(this.Text, base.Font, this.Width);
-            e.Graphics.DrawString(this.Text, base.Font, new SolidBrush(text_color), boxSize + 4, (this.Height - szL.Height) / 2);
+            DrawText(e.Graphics, text_color, new Point(boxSize + 4, (int)((this.Height - szL.Height) / 2) + 1));
 
             if (this.isChecked)
             {

Triangle.NET/TestApp/Controls/DarkTabControl.cs

@@ -16,7 +16,7 @@ namespace MeshExplorer.Controls
     /// <summary>
     /// Summary description for FlatTabControl.
     /// </summary>
-	public class DarkTabControl : System.Windows.Forms.TabControl
+    public class DarkTabControl : System.Windows.Forms.TabControl
     {
         #region Designer
 

Triangle.NET/TestApp/FormMain.Designer.cs

@@ -101,6 +101,8 @@
             this.menuTools = new System.Windows.Forms.ToolStripMenuItem();
             this.menuToolsGen = new System.Windows.Forms.ToolStripMenuItem();
             this.menuToolsCheck = new System.Windows.Forms.ToolStripMenuItem();
+            this.toolStripSeparator4 = new System.Windows.Forms.ToolStripSeparator();
+            this.menuToolsRcm = new System.Windows.Forms.ToolStripMenuItem();
             this.btnMesh = new MeshExplorer.Controls.DarkButton();
             this.renderControl1 = new MeshExplorer.Controls.RendererControl();
             ((System.ComponentModel.ISupportInitialize)(this.splitContainer1)).BeginInit();
@@ -755,7 +757,7 @@
             this.label19.Name = "label19";
             this.label19.Size = new System.Drawing.Size(134, 40);
             this.label19.TabIndex = 0;
-            this.label19.Text = "Beta 2 (2012-06-20)\r\nChristian Woltering\r\nMIT";
+            this.label19.Text = "Beta 2 (2012-07-30)\r\nChristian Woltering\r\nMIT";
             // 
             // label18
             // 
@@ -894,7 +896,9 @@
             // 
             this.menuTools.DropDownItems.AddRange(new System.Windows.Forms.ToolStripItem[] {
             this.menuToolsGen,
-            this.menuToolsCheck});
+            this.menuToolsCheck,
+            this.toolStripSeparator4,
+            this.menuToolsRcm});
             this.menuTools.Name = "menuTools";
             this.menuTools.Size = new System.Drawing.Size(46, 24);
             this.menuTools.Text = "Tools";
@@ -902,7 +906,7 @@
             // menuToolsGen
             // 
             this.menuToolsGen.Name = "menuToolsGen";
-            this.menuToolsGen.Size = new System.Drawing.Size(157, 22);
+            this.menuToolsGen.Size = new System.Drawing.Size(195, 22);
             this.menuToolsGen.Text = "Input Generator";
             this.menuToolsGen.Click += new System.EventHandler(this.menuToolsGenerator_Click);
             // 
@@ -910,10 +914,22 @@
             // 
             this.menuToolsCheck.Enabled = false;
             this.menuToolsCheck.Name = "menuToolsCheck";
-            this.menuToolsCheck.Size = new System.Drawing.Size(157, 22);
+            this.menuToolsCheck.Size = new System.Drawing.Size(195, 22);
             this.menuToolsCheck.Text = "Check Mesh";
             this.menuToolsCheck.Click += new System.EventHandler(this.menuToolsCheck_Click);
             // 
+            // toolStripSeparator4
+            // 
+            this.toolStripSeparator4.Name = "toolStripSeparator4";
+            this.toolStripSeparator4.Size = new System.Drawing.Size(192, 6);
+            // 
+            // menuToolsRcm
+            // 
+            this.menuToolsRcm.Name = "menuToolsRcm";
+            this.menuToolsRcm.Size = new System.Drawing.Size(195, 22);
+            this.menuToolsRcm.Text = "Renumber nodes (RCM)";
+            this.menuToolsRcm.Click += new System.EventHandler(this.menuToolsRcm_Click);
+            // 
             // btnMesh
             // 
             this.btnMesh.Enabled = false;
@@ -1050,6 +1066,8 @@
         private System.Windows.Forms.Label lbNumVert2;
         private System.Windows.Forms.Label label23;
         private Controls.DarkCheckBox cbConformDel;
+        private System.Windows.Forms.ToolStripSeparator toolStripSeparator4;
+        private System.Windows.Forms.ToolStripMenuItem menuToolsRcm;
 
     }
 }

Triangle.NET/TestApp/FormMain.cs

@@ -387,8 +387,6 @@ namespace MeshExplorer
 
             if (sfd.ShowDialog() == DialogResult.OK)
             {
-                mesh.Renumber();
-
                 FileProcessor.Save(sfd.FileName, mesh);
             }
         }
@@ -523,6 +521,19 @@ namespace MeshExplorer
             UpdateLog();
         }
 
+        private void Renumber()
+        {
+            if (mesh == null || settings.ExceptionThrown) return;
+
+            bool tmp = Behavior.Verbose;
+            Behavior.Verbose = true;
+
+            mesh.Renumber(NodeNumbering.CuthillMcKee);
+            ShowLog();
+
+            Behavior.Verbose = tmp;
+        }
+
         private void Smooth()
         {
             if (mesh == null || settings.ExceptionThrown) return;
@@ -674,5 +685,10 @@ namespace MeshExplorer
         {
             this.Close();
         }
+
+        private void menuToolsRcm_Click(object sender, EventArgs e)
+        {
+            Renumber();
+        }
     }
 }

Triangle.NET/TestApp/IO/EpsImage.cs

@@ -261,7 +261,7 @@ namespace MeshExplorer.IO
 
         private void DrawPoints(StreamWriter eps, Mesh mesh, bool label)
         {
-            int n = mesh.NumberOfVertices;
+            int n = mesh.Vertices.Count;
 
             int circle_size = 1;
 

Triangle.NET/TestApp/IO/FileProcessor.cs

@@ -97,10 +97,6 @@ namespace MeshExplorer.IO
             {
                 provider = new JsonFile();
             }
-            else if (ext == ".dat")
-            {
-                provider = new DatFile();
-            }
 
             if (provider == null)
             {

Triangle.NET/TestApp/IO/Formats/DatFile.cs

@@ -1,80 +0,0 @@
-// -----------------------------------------------------------------------
-// <copyright file="DatFile.cs" company="">
-// Christian Woltering, Triangle.NET, http://triangle.codeplex.com/
-// </copyright>
-// -----------------------------------------------------------------------
-
-namespace MeshExplorer.IO.Formats
-{
-    using System;
-    using System.Collections.Generic;
-    using System.Linq;
-    using System.Text;
-    using TriangleNet.IO;
-    using System.IO;
-    using MeshExplorer.Rendering;
-    using TriangleNet.Geometry;
-    using TriangleNet;
-
-    /// <summary>
-    /// Read a polygon from DAT file (used by poly2tri).
-    /// </summary>
-    public class DatFile : IMeshFile
-    {
-        /// <summary>
-        /// Gets the supported file extensions.
-        /// </summary>
-        public string[] Extensions
-        {
-            get { return new string[] { ".dat" }; }
-        }
-
-        public bool ContainsMeshData(string filename)
-        {
-            return false;
-        }
-
-        public InputGeometry Read(string filename)
-        {
-            InputGeometry data = new InputGeometry();
-
-            string line;
-            string[] split;
-
-            using (TextReader reader = new StreamReader(filename))
-            {
-                while ((line = reader.ReadLine()) != null)
-                {
-                    split = line.Split(new char[] { ' ' }, StringSplitOptions.RemoveEmptyEntries);
-
-                    if (split.Length == 2)
-                    {
-                        data.AddPoint(
-                            double.Parse(split[0], Util.Nfi),
-                            double.Parse(split[1], Util.Nfi));
-                    }
-                }
-            }
-
-            int n = data.Count;
-
-            for (int i = 0; i < n; i++)
-            {
-                data.AddSegment(i, (i + 1) % n);
-
-            }
-
-            return data;
-        }
-
-        public Mesh Import(string filename)
-        {
-            throw new NotImplementedException();
-        }
-
-        public void Write(Mesh mesh, string filename)
-        {
-            throw new NotImplementedException();
-        }
-    }
-}

Triangle.NET/TestApp/IO/Formats/JsonFile.cs

@@ -16,6 +16,7 @@ namespace MeshExplorer.IO.Formats
     using TriangleNet.Geometry;
     using TriangleNet;
     using System.Collections;
+    using TriangleNet.Data;
 
     /// <summary>
     /// Read and write JSON files.
@@ -101,10 +102,10 @@ namespace MeshExplorer.IO.Formats
         {
             using (StreamWriter writer = new StreamWriter(filename))
             {
-                int nv = mesh.NumberOfVertices;
+                int nv = mesh.Vertices.Count;
-                int ns = mesh.NumberOfSegments;
+                int ns = mesh.Segments.Count;
                 int nh = mesh.Holes.Count;
-                int ne = mesh.NumberOfTriangles;
+                int ne = mesh.Triangles.Count;
 
                 writer.Write("{");
 
@@ -115,6 +116,11 @@ namespace MeshExplorer.IO.Formats
                 writer.Write("\"dim\":2");
                 writer.Write("}");
 
+                if (mesh.CurrentNumbering == NodeNumbering.None)
+                {
+                    mesh.Renumber(NodeNumbering.Linear);
+                }
+
                 // Write the coordinates
                 if (nv > 0)
                 {
@@ -397,17 +403,47 @@ namespace MeshExplorer.IO.Formats
 
         #region Write helpers
 
-        private void WritePoints(Mesh data, StreamWriter writer, int nv)
+        private void WritePoints(Mesh mesh, StreamWriter writer, int nv)
         {
-            int i = 0;
-
-            StringBuilder markers = new StringBuilder();
             bool useMarkers = false;
 
-            string seperator;
+            StringBuilder markers;
 
             writer.Write("\"points\":{\"data\":[");
-            foreach (var item in data.Vertices)
+
+            if (mesh.CurrentNumbering == NodeNumbering.Linear)
+            {
+                markers = WritePoints(writer, mesh.Vertices, nv, useMarkers);
+            }
+            else
+            {
+                Vertex[] nodes = new Vertex[mesh.Vertices.Count];
+
+                foreach (var node in mesh.Vertices)
+                {
+                    nodes[node.ID] = node;
+                }
+
+                markers = WritePoints(writer, nodes, nv, useMarkers);
+            }
+
+            writer.Write("]");
+            if (useMarkers)
+            {
+                writer.Write(",\"markers\":[" + markers.ToString() + "]");
+            }
+
+            // TODO: writer.Write(",\"attributes\":[]");
+            writer.Write("}");
+        }
+
+        private static StringBuilder WritePoints(StreamWriter writer, IEnumerable<Vertex> nodes, int nv, bool useMarkers)
+        {
+            StringBuilder markers = new StringBuilder();
+
+            int i = 0;
+            string seperator;
+            foreach (var item in nodes)
             {
                 seperator = (i == nv - 1) ? String.Empty : ", ";
 
@@ -424,13 +460,8 @@ namespace MeshExplorer.IO.Formats
 
                 i++;
             }
-            writer.Write("]");
+
-            if (useMarkers)
+            return markers;
-            {
-                writer.Write(",\"markers\":[" + markers.ToString() + "]");
-            }
-            //writer.Write(",\"attributes\":[]");
-            writer.Write("}");
         }
 
         private void WriteHoles(Mesh data, StreamWriter writer, int nh)

Triangle.NET/TestApp/IO/Formats/TriangleFile.cs

@@ -57,7 +57,7 @@ namespace MeshExplorer.IO.Formats
 
         public void Write(Mesh mesh, string filename)
         {
-            if (mesh.NumberOfVertices > 0)
+            if (mesh.Vertices.Count > 0)
             {
                 format.Write(mesh, filename);
             }

Triangle.NET/TestApp/IO/RasterImage.cs

@@ -51,13 +51,13 @@ namespace MeshExplorer.IO
             Bitmap bitmap;
 
             // Check if the specified width is reasonable
-            if (width < 2 * Math.Sqrt(mesh.NumberOfVertices))
+            if (width < 2 * Math.Sqrt(mesh.Vertices.Count))
             {
                 bitmap = new Bitmap(400, 200);
                 Graphics g = Graphics.FromImage(bitmap);
                 g.Clear(colors.Background);
 
-                string message = String.Format("Sorry, I won't render {0} points on such a small image!", mesh.NumberOfVertices);
+                string message = String.Format("Sorry, I won't render {0} points on such a small image!", mesh.Vertices.Count);
 
                 SizeF sz = g.MeasureString(message, SystemFonts.DefaultFont);
 

Triangle.NET/TestApp/IO/SvgImage.cs

@@ -153,7 +153,7 @@ namespace MeshExplorer.IO
 
         private void DrawPoints(StreamWriter svg, Mesh mesh, bool label)
         {
-            int n = mesh.NumberOfVertices;
+            int n = mesh.Vertices.Count;
 
             int circle_size = 1;
 

Triangle.NET/TestApp/Mesh Explorer.csproj

@@ -107,7 +107,6 @@
     <Compile Include="Generators\StarInBox.cs" />
     <Compile Include="IO\EpsImage.cs" />
     <Compile Include="IO\FileProcessor.cs" />
-    <Compile Include="IO\Formats\DatFile.cs" />
     <Compile Include="IO\Formats\JsonFile.cs" />
     <Compile Include="IO\Formats\TriangleFile.cs" />
     <Compile Include="IO\IMeshFile.cs" />

Triangle.NET/TestApp/Rendering/RenderData.cs

@@ -73,7 +73,7 @@ namespace MeshExplorer.Rendering
 
             this.Segments = null;
 
-            int n = mesh.NumberOfVertices;
+            int n = mesh.Vertices.Count;
 
             // Convert points to float
             this.Points = new PointF[n];
@@ -85,7 +85,7 @@ namespace MeshExplorer.Rendering
             {
                 var segs = mesh.Segments;
 
-                List<Edge> segList = new List<Edge>(mesh.NumberOfSegments);
+                List<Edge> segList = new List<Edge>(mesh.Segments.Count);
 
                 foreach (var seg in segs)
                 {

Triangle.NET/TestApp/Settings.cs

@@ -59,7 +59,7 @@ namespace MeshExplorer
             SfdFilterIndex = 1;
 
             OfdFilter = SfdFilter;
-            OfdFilter += "|Polygon data (*.dat)|*.dat";
+            //OfdFilter += "|Polygon data (*.dat)|*.dat";
             //OfdFilter += "|COMSOL mesh (*.mphtxt)|*.mphtxt";
             //OfdFilter += "|AVS UCD data (*.ucd)|*.ucd";
             //OfdFilter += "|VTK data (*.vtk)|*.vtk";

Triangle.NET/Triangle/BadTriQueue.cs

@@ -162,10 +162,6 @@ namespace TriangleNet
             newbad.triangorg = enqorg;
             newbad.triangdest = enqdest;
 
-            Vertex org = enqtri.Org();
-            Vertex dest = enqtri.Dest();
-            Vertex apex = enqtri.Apex();
-
             Enqueue(newbad);
         }
 

Triangle.NET/Triangle/Behavior.cs

@@ -12,7 +12,7 @@ namespace TriangleNet
     /// <summary>
     /// Controls the behavior of the meshing software.
     /// </summary>
-    class Behavior
+    public class Behavior
     {
         /// <summary>
         /// Load behavior defaults.

Triangle.NET/Triangle/Data/Vertex.cs

@@ -27,18 +27,51 @@ namespace TriangleNet.Data
         internal VertexType type;
         internal Otri tri;
 
+        /// <summary>
+        /// Initializes a new instance of the <see cref="Vertex" /> class.
+        /// </summary>
         public Vertex()
-            : this(0, 0, 0)
+            : this(0, 0, 0, 0)
-        { }
+        {
+        }
 
+        /// <summary>
+        /// Initializes a new instance of the <see cref="Vertex" /> class.
+        /// </summary>
+        /// <param name="x">The x coordinate of the vertex.</param>
+        /// <param name="y">The y coordinate of the vertex.</param>
         public Vertex(double x, double y)
-            : this(x, y, 0)
+            : this(x, y, 0, 0)
-        { }
+        {
+        }
 
+        /// <summary>
+        /// Initializes a new instance of the <see cref="Vertex" /> class.
+        /// </summary>
+        /// <param name="x">The x coordinate of the vertex.</param>
+        /// <param name="y">The y coordinate of the vertex.</param>
+        /// <param name="mark">The boundary mark.</param>
         public Vertex(double x, double y, int mark)
+            : this(x, y, mark, 0)
+        {
+        }
+
+        /// <summary>
+        /// Initializes a new instance of the <see cref="Vertex" /> class.
+        /// </summary>
+        /// <param name="x">The x coordinate of the vertex.</param>
+        /// <param name="y">The y coordinate of the vertex.</param>
+        /// <param name="mark">The boundary mark.</param>
+        /// <param name="attribs">The number of point attributes.</param>
+        public Vertex(double x, double y, int mark, int attribs)
             : base(x, y, mark)
         {
             this.type = VertexType.InputVertex;
+
+            if (attribs > 0)
+            {
+                this.attributes = new double[attribs];
+            }
         }
 
         #region Public properties

Triangle.NET/Triangle/Enums.cs

@@ -95,5 +95,5 @@ namespace TriangleNet
     /// <summary>
     /// Node renumbering algorithms.
     /// </summary>
-    public enum NodeNumbering { Linear, CuthillMcKee };
+    public enum NodeNumbering { None, Linear, CuthillMcKee };
 }

Triangle.NET/Triangle/Geometry/InputGeometry.cs

@@ -23,6 +23,9 @@ namespace TriangleNet.Geometry
 
         BoundingBox bounds;
 
+        // Used to check consitent use of point attributes.
+        private int pointAttributes = -1;
+
         /// <summary>
         /// Initializes a new instance of the <see cref="InputGeometry" /> class.
         /// </summary>
@@ -44,6 +47,8 @@ namespace TriangleNet.Geometry
             regions = new List<RegionPointer>();
 
             bounds = new BoundingBox();
+
+            pointAttributes = -1;
         }
 
         /// <summary>
@@ -55,7 +60,7 @@ namespace TriangleNet.Geometry
         }
 
         /// <summary>
-        /// Gets a value indicating whether the geometry should be treated as a PLSG.
+        /// Indicates, whether the geometry should be treated as a PSLG.
         /// </summary>
         public bool HasSegments
         {
@@ -95,7 +100,7 @@ namespace TriangleNet.Geometry
         }
 
         /// <summary>
-        /// Gets the list of input holes.
+        /// Gets the list of regions.
         /// </summary>
         public IEnumerable<RegionPointer> Regions
         {
@@ -111,6 +116,8 @@ namespace TriangleNet.Geometry
             segments.Clear();
             holes.Clear();
             regions.Clear();
+
+            pointAttributes = -1;
         }
 
         /// <summary>
@@ -136,6 +143,45 @@ namespace TriangleNet.Geometry
             bounds.Update(x, y);
         }
 
+        /// <summary>
+        /// Adds a point to the geometry.
+        /// </summary>
+        /// <param name="x">X coordinate.</param>
+        /// <param name="y">Y coordinate.</param>
+        /// <param name="boundary">Boundary marker.</param>
+        /// <param name="attribute">Point attribute.</param>
+        public void AddPoint(double x, double y, int boundary, double attribute)
+        {
+            AddPoint(x, y, 0, new double[] { attribute });
+        }
+
+        /// <summary>
+        /// Adds a point to the geometry.
+        /// </summary>
+        /// <param name="x">X coordinate.</param>
+        /// <param name="y">Y coordinate.</param>
+        /// <param name="boundary">Boundary marker.</param>
+        /// <param name="attribs">Point attributes.</param>
+        public void AddPoint(double x, double y, int boundary, double[] attribs)
+        {
+            if (pointAttributes < 0)
+            {
+                pointAttributes = attribs == null ? 0 : attribs.Length;
+            }
+            else if (attribs == null && pointAttributes > 0)
+            {
+                throw new ArgumentException("Inconsitent use of point attributes.");
+            }
+            else if (attribs != null && pointAttributes != attribs.Length)
+            {
+                throw new ArgumentException("Inconsitent use of point attributes.");
+            }
+
+            points.Add(new Vertex(x, y, boundary) { attributes = attribs });
+
+            bounds.Update(x, y);
+        }
+
         /// <summary>
         /// Adds a hole location to the geometry.
         /// </summary>

Triangle.NET/Triangle/Geometry/Point.cs

@@ -1,6 +1,5 @@
 // -----------------------------------------------------------------------
 // <copyright file="Point.cs" company="">
-// Original Triangle code by Jonathan Richard Shewchuk, http://www.cs.cmu.edu/~quake/triangle.html
 // Triangle.NET code by Christian Woltering, http://triangle.codeplex.com/
 // </copyright>
 // -----------------------------------------------------------------------
@@ -144,6 +143,11 @@ namespace TriangleNet.Geometry
             return (x < other.x || (x == other.x && y < other.y)) ? -1 : 1;
         }
 
+        public override int GetHashCode()
+        {
+            return x.GetHashCode() ^ y.GetHashCode();
+        }
+
         public override string ToString()
         {
             return String.Format("[{0},{1}]", x, y);

Triangle.NET/Triangle/Geometry/RegionPointer.cs

@@ -1,6 +1,5 @@
 // -----------------------------------------------------------------------
 // <copyright file="RegionPointer.cs" company="">
-// Original Triangle code by Jonathan Richard Shewchuk, http://www.cs.cmu.edu/~quake/triangle.html
 // Triangle.NET code by Christian Woltering, http://triangle.codeplex.com/
 // </copyright>
 // -----------------------------------------------------------------------

Triangle.NET/Triangle/IO/DataReader.cs

@@ -92,7 +92,7 @@ namespace TriangleNet.IO
                 // Create the subsegments.
                 for (i = 0; i < mesh.insegments; i++)
                 {
-                    mesh.MakeSubseg(ref subseg);
+                    mesh.MakeSegment(ref subseg);
                     // Mark the subsegment as living.
                     //subseg.ss.subsegs[0].ss = subseg.ss;
                 }

Triangle.NET/Triangle/IO/FileReader.cs

@@ -57,30 +57,31 @@ namespace TriangleNet.IO
         /// <param name="data">The input geometry.</param>
         /// <param name="index">The current vertex index.</param>
         /// <param name="line">The current line.</param>
-        /// <param name="n">Number of point attributes</param>
+        /// <param name="attributes">Number of point attributes</param>
-        static void ReadVertex(InputGeometry data, int index, string[] line, int n)
+        /// <param name="marks">Number of point markers (0 or 1)</param>
+        static void ReadVertex(InputGeometry data, int index, string[] line, int attributes, int marks)
         {
             double x = double.Parse(line[1], nfi);
             double y = double.Parse(line[2], nfi);
             int mark = 0;
+            double[] attribs = attributes == 0 ? null : new double[attributes];
 
             // Read the vertex attributes.
-            for (int j = 0; j < n; j++)
+            for (int j = 0; j < attributes; j++)
             {
                 if (line.Length > 3 + j)
                 {
-                    // TODO:
+                    attribs[j] = double.Parse(line[3 + j]);
-                    //vertex.attributes[j] = double.Parse(line[3 + j]);
                 }
             }
 
             // Read a vertex marker.
-            if (line.Length > 3 + n)
+            if (marks > 0 && line.Length > 3 + attributes)
             {
-                mark = int.Parse(line[3 + n]);
+                mark = int.Parse(line[3 + attributes]);
             }
 
-            data.AddPoint(x, y, mark);
+            data.AddPoint(x, y, mark, attribs);
         }
 
         #endregion
@@ -218,7 +219,7 @@ namespace TriangleNet.IO
                             startIndex = int.Parse(line[0], nfi);
                         }
 
-                        ReadVertex(data, i, line, attributes);
+                        ReadVertex(data, i, line, attributes, nodemarkers);
                     }
                 }
             }
@@ -325,7 +326,7 @@ namespace TriangleNet.IO
                             startIndex = int.Parse(line[0], nfi);
                         }
 
-                        ReadVertex(data, i, line, attributes);
+                        ReadVertex(data, i, line, attributes, nodemarkers);
                     }
                 }
                 else

Triangle.NET/Triangle/IO/FileWriter.cs

@@ -12,6 +12,7 @@ namespace TriangleNet.IO
     using System.Globalization;
     using TriangleNet.Data;
     using TriangleNet.Geometry;
+    using System.Collections.Generic;
 
     /// <summary>
     /// Helper methods for writing Triangle file formats.
@@ -25,12 +26,10 @@ namespace TriangleNet.IO
         /// </summary>
         /// <param name="mesh"></param>
         /// <param name="filename"></param>
-        public static void WriteNodes(Mesh mesh, string filename)
+        public static void Write(Mesh mesh, string filename)
         {
-            using (StreamWriter writer = new StreamWriter(filename))
+            FileWriter.WritePoly(mesh, Path.ChangeExtension(filename, ".poly"));
-            {
+            FileWriter.WriteElements(mesh, Path.ChangeExtension(filename, ".ele"));
-                FileWriter.WriteNodes(mesh, writer);
-            }
         }
 
         /// <summary>
@@ -38,10 +37,20 @@ namespace TriangleNet.IO
         /// </summary>
         /// <param name="mesh"></param>
         /// <param name="filename"></param>
-        private static void WriteNodes(Mesh mesh, StreamWriter writer)
+        public static void WriteNodes(Mesh mesh, string filename)
         {
-            Vertex vertex;
+            using (StreamWriter writer = new StreamWriter(filename))
-            long outvertices = mesh.vertices.Count;
+            {
+                FileWriter.WriteNodes(writer, mesh);
+            }
+        }
+
+        /// <summary>
+        /// Number the vertices and write them to a .node file.
+        /// </summary>
+        private static void WriteNodes(StreamWriter writer, Mesh mesh)
+        {
+            int outvertices = mesh.vertices.Count;
 
             Behavior behavior = mesh.behavior;
 
@@ -50,8 +59,6 @@ namespace TriangleNet.IO
                 outvertices = mesh.vertices.Count - mesh.undeads;
             }
 
-            int index = 0;
-
             if (writer != null)
             {
                 // Number of vertices, number of dimensions, number of vertex attributes,
@@ -59,32 +66,69 @@ namespace TriangleNet.IO
                 writer.WriteLine("{0} {1} {2} {3}", outvertices, mesh.mesh_dim, mesh.nextras,
                     behavior.UseBoundaryMarkers ? "1" : "0");
 
-                foreach (var item in mesh.vertices.Values)
+                if (mesh.numbering == NodeNumbering.None)
                 {
-                    vertex = item;
+                    // If the mesh isn't numbered yet, use linear node numbering.
+                    mesh.Renumber();
+                }
 
-                    if (!behavior.Jettison || vertex.type != VertexType.UndeadVertex)
+                if (mesh.numbering == NodeNumbering.Linear)
+                {
+                    // If numbering is linear, just use the dictionary values.
+                    WriteNodes(writer, mesh.vertices.Values, behavior.UseBoundaryMarkers,
+                        mesh.nextras, behavior.Jettison);
+                }
+                else
+                {
+                    // If numbering is not linear, a simple 'foreach' traversal of the dictionary
+                    // values doesn't reflect the actual numbering. Use an array instead.
+
+                    // TODO: Could use a custom sorting function on dictionary values instead.
+                    Vertex[] nodes = new Vertex[mesh.vertices.Count];
+
+                    foreach (var node in mesh.vertices.Values)
                     {
-                        // Vertex number, x and y coordinates.
+                        nodes[node.id] = node;
-                        writer.Write("{0} {1} {2}", index, vertex.x.ToString(nfi), vertex.y.ToString(nfi));
-
-                        // Write attributes.
-                        for (int j = 0; j < mesh.nextras; j++)
-                        {
-                            writer.Write(" {0}", vertex.attributes[j].ToString(nfi));
-                        }
-
-                        if (behavior.UseBoundaryMarkers)
-                        {
-                            // Write the boundary marker.
-                            writer.Write(" {0}", vertex.mark);
-                        }
-
-                        writer.WriteLine();
-
-                        // Assign array index to vertex ID for later use.
-                        vertex.id = index++;
                     }
+
+                    WriteNodes(writer, nodes, behavior.UseBoundaryMarkers,
+                        mesh.nextras, behavior.Jettison);
+                }
+            }
+        }
+
+        /// <summary>
+        /// Write the vertices to a stream.
+        /// </summary>
+        /// <param name="nodes"></param>
+        /// <param name="writer"></param>
+        private static void WriteNodes(StreamWriter writer, IEnumerable<Vertex> nodes, bool markers,
+            int attribs, bool jettison)
+        {
+            int index = 0;
+
+            foreach (var vertex in nodes)
+            {
+                if (!jettison || vertex.type != VertexType.UndeadVertex)
+                {
+                    // Vertex number, x and y coordinates.
+                    writer.Write("{0} {1} {2}", index, vertex.x.ToString(nfi), vertex.y.ToString(nfi));
+
+                    // Write attributes.
+                    for (int j = 0; j < attribs; j++)
+                    {
+                        writer.Write(" {0}", vertex.attributes[j].ToString(nfi));
+                    }
+
+                    if (markers)
+                    {
+                        // Write the boundary marker.
+                        writer.Write(" {0}", vertex.mark);
+                    }
+
+                    writer.WriteLine();
+
+                    index++;
                 }
             }
         }
@@ -163,7 +207,7 @@ namespace TriangleNet.IO
                 if (writeNodes)
                 {
                     // Write nodes to this file.
-                    FileWriter.WriteNodes(mesh, writer);
+                    FileWriter.WriteNodes(writer, mesh);
                 }
                 else
                 {
@@ -202,10 +246,11 @@ namespace TriangleNet.IO
                 }
 
                 // Holes
+                j = 0;
                 writer.WriteLine("{0}", mesh.holes.Count);
                 foreach (var hole in mesh.holes)
                 {
-                    writer.WriteLine("{0} {1}", hole.X.ToString(nfi), hole.Y.ToString(nfi));
+                    writer.WriteLine("{0} {1} {2}", j++, hole.X.ToString(nfi), hole.Y.ToString(nfi));
                 }
 
                 // Regions

Triangle.NET/Triangle/Mesh.cs

@@ -32,7 +32,6 @@ namespace TriangleNet
 
         // Stack that maintains a list of recently flipped triangles.
         Stack<Otri> flipstack;
-        //FlipStacker lastflip;
 
         // TODO: Check if custom hashmap implementation could be faster.
 
@@ -69,12 +68,11 @@ namespace TriangleNet
         // Triangular bounding box vertices.
         internal Vertex infvertex1, infvertex2, infvertex3;
 
-        // The 'triangle' that occupies all of "outer space."
+        // The 'triangle' that occupies all of 'outer space'.
         internal static Triangle dummytri;
 
-        // The omnipresent subsegment. Referenced by any triangle or
+        // The omnipresent subsegment. Referenced by any triangle or subsegment
-        // subsegment that isn't really connected to a subsegment at
+        // that isn't really connected to a subsegment at that location.
-        // that location.
         internal static Segment dummysub;
 
         // Pointer to a recently visited triangle. Improves point location if
@@ -84,6 +82,9 @@ namespace TriangleNet
         // Controls the behavior of the mesh instance.
         internal Behavior behavior;
 
+        // The current node numbering
+        internal NodeNumbering numbering;
+
         #endregion
 
         #region Public properties
@@ -99,7 +100,7 @@ namespace TriangleNet
         /// <summary>
         /// Gets the mesh vertices.
         /// </summary>
-        public IEnumerable<Vertex> Vertices
+        public ICollection<Vertex> Vertices
         {
             get { return this.vertices.Values; }
         }
@@ -115,7 +116,7 @@ namespace TriangleNet
         /// <summary>
         /// Gets the mesh triangles.
         /// </summary>
-        public IEnumerable<Triangle> Triangles
+        public ICollection<Triangle> Triangles
         {
             get { return this.triangles.Values; }
         }
@@ -123,7 +124,7 @@ namespace TriangleNet
         /// <summary>
         /// Gets the mesh segments.
         /// </summary>
-        public IEnumerable<Segment> Segments
+        public ICollection<Segment> Segments
         {
             get { return this.subsegs.Values; }
         }
@@ -133,21 +134,6 @@ namespace TriangleNet
         /// </summary>
         public int NumberOfInputPoints { get { return invertices; } }
 
-        /// <summary>
-        /// Gets the number of mesh vertices.
-        /// </summary>
-        public int NumberOfVertices { get { return this.vertices.Count; } }
-
-        /// <summary>
-        /// Gets the number of mesh triangles.
-        /// </summary>
-        public int NumberOfTriangles { get { return this.triangles.Count; } }
-
-        /// <summary>
-        /// Gets the number of mesh segments.
-        /// </summary>
-        public int NumberOfSegments { get { return this.subsegs.Count; } }
-
         /// <summary>
         /// Gets the number of mesh edges.
         /// </summary>
@@ -158,6 +144,14 @@ namespace TriangleNet
         /// </summary>
         public bool IsPolygon { get { return this.insegments > 0; } }
 
+        /// <summary>
+        /// Gets the current node numbering.
+        /// </summary>
+        public NodeNumbering CurrentNumbering
+        {
+            get { return numbering; }
+        }
+
         #endregion
 
         /// <summary>
@@ -461,6 +455,8 @@ namespace TriangleNet
         /// </summary>
         public void Smooth()
         {
+            numbering = NodeNumbering.None;
+
             //ISmoother smoother = new CvdSmoother(this);
             //smoother.Smooth();
         }
@@ -478,29 +474,38 @@ namespace TriangleNet
         /// </summary>
         public void Renumber(NodeNumbering num)
         {
+            // Don't need to do anything if the nodes are already numbered.
+            if (num == this.numbering)
+            {
+                return;
+            }
+
             int id;
 
             if (num == NodeNumbering.Linear)
             {
                 id = 0;
-                foreach (var item in this.vertices.Values)
+                foreach (var node in this.vertices.Values)
                 {
-                    item.id = id++;
+                    node.id = id++;
                 }
             }
-            else
+            else if (num == NodeNumbering.CuthillMcKee)
             {
-                //CuthillMcKee rcm = new CuthillMcKee();
+                CuthillMcKee rcm = new CuthillMcKee();
-                //int[] perm_inv = rcm.Renumber(this);
+                int[] perm_inv = rcm.Renumber(this);
 
-                //// Permute the node indices.
+                // Permute the node indices.
-                //foreach (var node in this.vertices.Values)
+                foreach (var node in this.vertices.Values)
-                //{
+                {
-                //    node.id = perm_inv[node.id];
+                    node.id = perm_inv[node.id];
-                //}
+                }
             }
 
-            // Triangles will always be numbered from 0..n
+            // Remember the current numbering.
+            numbering = num;
+
+            // Triangles will always be numbered from 0 to n-1
             id = 0;
             foreach (var item in this.triangles.Values)
             {
@@ -723,6 +728,8 @@ namespace TriangleNet
         /// </summary>
         private void Reset()
         {
+            numbering = NodeNumbering.None;
+
             recenttri.triangle = null; // No triangle has been visited yet.
             undeads = 0;               // No eliminated input vertices yet.
             checksegments = false;     // There are no segments in the triangulation yet.
@@ -816,13 +823,10 @@ namespace TriangleNet
                 throw new Exception("Input must have at least three input vertices.");
             }
 
-            this.nextras = 0; // TODO: points[0].Attributes == null ? 0 : points[0].Attributes.Length;
+            this.nextras = points[0].attributes == null ? 0 : points[0].attributes.Length;
 
             foreach (Vertex vertex in points)
             {
-                // TODO: Set vertex attributes.
-                //vertex.attribs = points[i].Attributes;
-
                 vertex.hash = this.hash_vtx++;
                 vertex.id = vertex.hash;
 
@@ -856,7 +860,7 @@ namespace TriangleNet
         /// Create a new subsegment with orientation zero.
         /// </summary>
         /// <param name="newsubseg">Reference to the new subseg.</param>
-        internal void MakeSubseg(ref Osub newsubseg)
+        internal void MakeSegment(ref Osub newsubseg)
         {
             Segment seg = new Segment();
             seg.hash = this.hash_seg++;
@@ -866,6 +870,7 @@ namespace TriangleNet
 
             subsegs.Add(seg.hash, seg);
         }
+
         #endregion
 
         #region Manipulation
@@ -1485,7 +1490,7 @@ namespace TriangleNet
             if (newsubseg.seg == dummysub)
             {
                 // Make new subsegment and initialize its vertices.
-                MakeSubseg(ref newsubseg);
+                MakeSegment(ref newsubseg);
                 newsubseg.SetOrg(tridest);
                 newsubseg.SetDest(triorg);
                 newsubseg.SetSegOrg(tridest);
@@ -2510,7 +2515,7 @@ namespace TriangleNet
             Vertex leftvertex, rightvertex;
             Vertex newvertex;
             InsertVertexResult success;
-            FindDirectionResult collinear;
+
             double ex, ey;
             double tx, ty;
             double etx, ety;
@@ -2535,13 +2540,23 @@ namespace TriangleNet
                 throw new Exception("Attempt to find intersection of parallel segments.");
             }
             split = (ey * etx - ex * ety) / denom;
+
             // Create the new vertex.
-            newvertex = new Vertex(torg.x + split * (tdest.x - torg.x),
+            newvertex = new Vertex(
-                torg.y + split * (tdest.y - torg.y), splitsubseg.seg.boundary);
+                torg.x + split * (tdest.x - torg.x),
+                torg.y + split * (tdest.y - torg.y),
+                splitsubseg.seg.boundary,
+                this.nextras);
+
             newvertex.hash = this.hash_vtx++;
             newvertex.id = newvertex.hash;
-            // TODO: nextras //(vertex) poolalloc(&m.vertices);
+
             // Interpolate its attributes.
+            for (int i = 0; i < nextras; i++)
+            {
+                newvertex.attributes[i] = torg.attributes[i] + split * (tdest.attributes[i] - torg.attributes[i]);
+            }
+
             vertices.Add(newvertex.hash, newvertex);
 
             // Insert the intersection vertex.  This should always succeed.
@@ -2576,7 +2591,8 @@ namespace TriangleNet
 
             // Inserting the vertex may have caused edge flips.  We wish to rediscover
             // the edge connecting endpoint1 to the new intersection vertex.
-            collinear = FindDirection(ref splittri, endpoint1);
+            FindDirection(ref splittri, endpoint1);
+
             rightvertex = splittri.Dest();
             leftvertex = splittri.Apex();
             if ((leftvertex.x == endpoint1.x) && (leftvertex.y == endpoint1.y))

Triangle.NET/Triangle/NewLocation.cs

@@ -34,28 +34,25 @@ namespace TriangleNet
         /// <summary>
         /// Find a new location for a Steiner point.
         /// </summary>
-        /// <param name="mesh"></param>
         /// <param name="torg"></param>
         /// <param name="tdest"></param>
         /// <param name="tapex"></param>
-        /// <param name="circumcenter"></param>
         /// <param name="xi"></param>
         /// <param name="eta"></param>
         /// <param name="offcenter"></param>
         /// <param name="badotri"></param>
-        public void FindLocation(Vertex torg, Vertex tdest, Vertex tapex,
+        /// <returns></returns>
-                              Vertex circumcenter, ref double xi, ref double eta, bool offcenter, Otri badotri) // TODO: ref circumcenter???
+        public Point FindLocation(Vertex torg, Vertex tdest, Vertex tapex,
+            ref double xi, ref double eta, bool offcenter, Otri badotri)
         {
             // Based on using -U switch, call the corresponding function
             if (behavior.MaxAngle == 0.0)
             {
-                FindNewLocationWithoutMaxAngle(torg, tdest, tapex, circumcenter, ref xi, ref eta, true, badotri);
+                return FindNewLocationWithoutMaxAngle(torg, tdest, tapex, ref xi, ref eta, true, badotri);
-            }
-            else
-            {
-                FindNewLocation(torg, tdest, tapex, circumcenter, ref xi, ref eta, true, badotri);
             }
 
+            // With max angle
+            return FindNewLocation(torg, tdest, tapex, ref xi, ref eta, true, badotri);
         }
 
         /// <summary>
@@ -69,8 +66,8 @@ namespace TriangleNet
         /// <param name="eta"></param>
         /// <param name="offcenter"></param>
         /// <param name="badotri"></param>
-        private void FindNewLocationWithoutMaxAngle(Vertex torg, Vertex tdest, Vertex tapex,
+        private Point FindNewLocationWithoutMaxAngle(Vertex torg, Vertex tdest, Vertex tapex,
-                              Vertex circumcenter, ref double xi, ref double eta, bool offcenter, Otri badotri)
+            ref double xi, ref double eta, bool offcenter, Otri badotri)
         {
             double offconstant = behavior.Offconstant;
 
@@ -740,6 +737,8 @@ namespace TriangleNet
                 }// end of relocation				 
             }// end of almostGood	
 
+            Point circumcenter = new Point();
+
             if (relocated <= 0)
             {
                 circumcenter.x = torg.x + dx;
@@ -754,6 +753,7 @@ namespace TriangleNet
             xi = (yao * dx - xao * dy) * (2.0 * denominator);
             eta = (xdo * dy - ydo * dx) * (2.0 * denominator);
 
+            return circumcenter;
         }
 
         /// <summary>
@@ -767,8 +767,8 @@ namespace TriangleNet
         /// <param name="eta"></param>
         /// <param name="offcenter"></param>
         /// <param name="badotri"></param>
-        private void FindNewLocation(Vertex torg, Vertex tdest, Vertex tapex,
+        private Point FindNewLocation(Vertex torg, Vertex tdest, Vertex tapex,
-                              Vertex circumcenter, ref double xi, ref double eta, bool offcenter, Otri badotri)
+            ref double xi, ref double eta, bool offcenter, Otri badotri)
         {
             double offconstant = behavior.Offconstant;
 
@@ -1300,9 +1300,7 @@ namespace TriangleNet
                                         && (IsBadTriangleAngle(middleAngleCorner.x, middleAngleCorner.y, largestAngleCorner.x, largestAngleCorner.y, petal_slab_inter_x_first, petal_slab_inter_y_first))
                                         && MinDistanceToNeighbor(petal_slab_inter_x_first, petal_slab_inter_y_first, ref neighborotri) > MinDistanceToNeighbor(line_inter_x, line_inter_y, ref neighborotri))
                                     {
-                                        // 							
+                                        // check the neighbor's vertices also, which one if better
-                                        /// check the neighbor's vertices also, which one if better
-
                                         //slab and petal intersection is advised
                                         dxFirstSuggestion = petal_slab_inter_x_first - torg.x;
                                         dyFirstSuggestion = petal_slab_inter_y_first - torg.y;
@@ -1329,26 +1327,19 @@ namespace TriangleNet
                                                 // go back to circumcenter
                                                 dxFirstSuggestion = dx;
                                                 dyFirstSuggestion = dy;
-
-
                                             }
                                             else
                                             {
                                                 // intersection point is suggested
                                                 dxFirstSuggestion = line_inter_x - torg.x;
                                                 dyFirstSuggestion = line_inter_y - torg.y;
-
-
                                             }
-
-
                                         }
                                         else
                                         {// we are not creating a bad triangle
                                             // slab intersection is advised
                                             dxFirstSuggestion = line_result[2] - torg.x;
                                             dyFirstSuggestion = line_result[3] - torg.y;
-
                                         }
                                     }
                                     //------------------------------------------------------//
@@ -1361,8 +1352,6 @@ namespace TriangleNet
                                         // go back to circumcenter
                                         dxFirstSuggestion = dx;
                                         dyFirstSuggestion = dy;
-
-
                                     }
                                     else
                                     {
@@ -1370,11 +1359,8 @@ namespace TriangleNet
                                         // neighbor's circumcenter is suggested
                                         dxFirstSuggestion = voronoiOrInter[2] - torg.x;
                                         dyFirstSuggestion = voronoiOrInter[3] - torg.y;
-
-
                                     }
                                 }
-
                             }
                             else
                             { // there is no voronoi vertex between intersection point and circumcenter
@@ -1397,7 +1383,6 @@ namespace TriangleNet
                                         //slab and petal intersection is advised
                                         dxFirstSuggestion = petal_slab_inter_x_first - torg.x;
                                         dyFirstSuggestion = petal_slab_inter_y_first - torg.y;
-
                                     }
                                     else
                                     { // slab intersection point is further away							
@@ -1421,32 +1406,25 @@ namespace TriangleNet
                                                 // go back to circumcenter
                                                 dxFirstSuggestion = dx;
                                                 dyFirstSuggestion = dy;
-
                                             }
                                             else
                                             {
                                                 // intersection point is suggested
                                                 dxFirstSuggestion = line_inter_x - torg.x;
                                                 dyFirstSuggestion = line_inter_y - torg.y;
-
                                             }
-
-
                                         }
                                         else
                                         {// we are not creating a bad triangle
                                             // slab intersection is advised
                                             dxFirstSuggestion = line_result[2] - torg.x;
                                             dyFirstSuggestion = line_result[3] - torg.y;
-
                                         }
                                     }
                                     //------------------------------------------------------//
-
                                 }
                                 else
                                 {
-
                                     if (IsBadTriangleAngle(largestAngleCorner.x, largestAngleCorner.y, middleAngleCorner.x, middleAngleCorner.y, inter_x, inter_y))
                                     {
                                         //printf("testtriangle returned false! bad triangle\n");	
@@ -1469,16 +1447,12 @@ namespace TriangleNet
                                             // go back to circumcenter
                                             dxFirstSuggestion = dx;
                                             dyFirstSuggestion = dy;
-
-
-
                                         }
                                         else
                                         {
                                             // intersection point is suggested
                                             dxFirstSuggestion = inter_x - torg.x;
                                             dyFirstSuggestion = inter_y - torg.y;
-
                                         }
                                     }
                                     else
@@ -1486,7 +1460,6 @@ namespace TriangleNet
                                         // intersection point is suggested
                                         dxFirstSuggestion = inter_x - torg.x;
                                         dyFirstSuggestion = inter_y - torg.y;
-
                                     }
                                 }
                             }
@@ -1614,8 +1587,6 @@ namespace TriangleNet
                                         // slab and petal intersection is advised
                                         dxSecondSuggestion = petal_slab_inter_x_second - torg.x;
                                         dySecondSuggestion = petal_slab_inter_y_second - torg.y;
-
-
                                     }
                                     else
                                     { // slab intersection point is further away	
@@ -1639,8 +1610,6 @@ namespace TriangleNet
                                                 // go back to circumcenter
                                                 dxSecondSuggestion = dx;
                                                 dySecondSuggestion = dy;
-
-
                                             }
                                             else
                                             {
@@ -1649,15 +1618,12 @@ namespace TriangleNet
                                                 dySecondSuggestion = line_inter_y - torg.y;
 
                                             }
-
-
                                         }
                                         else
                                         {// we are not creating a bad triangle
                                             // slab intersection is advised
                                             dxSecondSuggestion = line_result[2] - torg.x;
                                             dySecondSuggestion = line_result[3] - torg.y;
-
                                         }
                                     }
                                     //------------------------------------------------------//
@@ -1669,16 +1635,12 @@ namespace TriangleNet
                                         // go back to circumcenter
                                         dxSecondSuggestion = dx;
                                         dySecondSuggestion = dy;
-
-
                                     }
                                     else
                                     { // we are not creating a bad triangle
                                         // neighbor's circumcenter is suggested
                                         dxSecondSuggestion = voronoiOrInter[2] - torg.x;
                                         dySecondSuggestion = voronoiOrInter[3] - torg.y;
-
-
                                     }
                                 }
                             }
@@ -1703,7 +1665,6 @@ namespace TriangleNet
                                         // slab and petal intersection is advised
                                         dxSecondSuggestion = petal_slab_inter_x_second - torg.x;
                                         dySecondSuggestion = petal_slab_inter_y_second - torg.y;
-
                                     }
                                     else
                                     { // slab intersection point is further away							;
@@ -1727,28 +1688,23 @@ namespace TriangleNet
                                                 // go back to circumcenter
                                                 dxSecondSuggestion = dx;
                                                 dySecondSuggestion = dy;
-
                                             }
                                             else
                                             {
                                                 // intersection point is suggested
                                                 dxSecondSuggestion = line_inter_x - torg.x;
                                                 dySecondSuggestion = line_inter_y - torg.y;
-
                                             }
-
-
                                         }
                                         else
-                                        {// we are not creating a bad triangle
+                                        {
+                                            // we are not creating a bad triangle
                                             // slab intersection is advised
                                             dxSecondSuggestion = line_result[2] - torg.x;
                                             dySecondSuggestion = line_result[3] - torg.y;
-
                                         }
                                     }
                                     //------------------------------------------------------//
-
                                 }
                                 else
                                 {
@@ -1773,24 +1729,19 @@ namespace TriangleNet
                                             // go back to circumcenter
                                             dxSecondSuggestion = dx;
                                             dySecondSuggestion = dy;
-
-
                                         }
                                         else
                                         {
                                             // intersection point is suggested
                                             dxSecondSuggestion = inter_x - torg.x;
                                             dySecondSuggestion = inter_y - torg.y;
-
                                         }
                                     }
                                     else
                                     {
-
                                         // intersection point is suggested
                                         dxSecondSuggestion = inter_x - torg.x;
                                         dySecondSuggestion = inter_y - torg.y;
-
                                     }
                                 }
                             }
@@ -1827,13 +1778,11 @@ namespace TriangleNet
                             {
                                 dx = dxSecondSuggestion;
                                 dy = dySecondSuggestion;
-
                             }
                             else
                             {
                                 dx = dxFirstSuggestion;
                                 dy = dyFirstSuggestion;
-
                             }
                         }
                         else if (neighborNotFound_first)
@@ -1850,13 +1799,11 @@ namespace TriangleNet
                             {
                                 dx = dxSecondSuggestion;
                                 dy = dySecondSuggestion;
-
                             }
                             else
                             {
                                 dx = dxFirstSuggestion;
                                 dy = dyFirstSuggestion;
-
                             }
                         }
                         else if (neighborNotFound_second)
@@ -1873,13 +1820,11 @@ namespace TriangleNet
                             {
                                 dx = dxSecondSuggestion;
                                 dy = dySecondSuggestion;
-
                             }
                             else
                             {
                                 dx = dxFirstSuggestion;
                                 dy = dyFirstSuggestion;
-
                             }
                         }
                         else
@@ -1896,16 +1841,13 @@ namespace TriangleNet
                             {
                                 dx = dxSecondSuggestion;
                                 dy = dySecondSuggestion;
-
                             }
                             else
                             {
                                 dx = dxFirstSuggestion;
                                 dy = dyFirstSuggestion;
-
                             }
                         }
-
                     }
                     else
                     { // acute : consider other direction
@@ -1923,13 +1865,11 @@ namespace TriangleNet
                             {
                                 dx = dxSecondSuggestion;
                                 dy = dySecondSuggestion;
-
                             }
                             else
                             {
                                 dx = dxFirstSuggestion;
                                 dy = dyFirstSuggestion;
-
                             }
                         }
                         else if (neighborNotFound_first)
@@ -1946,13 +1886,11 @@ namespace TriangleNet
                             {
                                 dx = dxSecondSuggestion;
                                 dy = dySecondSuggestion;
-
                             }
                             else
                             {
                                 dx = dxFirstSuggestion;
                                 dy = dyFirstSuggestion;
-
                             }
                         }
                         else if (neighborNotFound_second)
@@ -1969,13 +1907,11 @@ namespace TriangleNet
                             {
                                 dx = dxSecondSuggestion;
                                 dy = dySecondSuggestion;
-
                             }
                             else
                             {
                                 dx = dxFirstSuggestion;
                                 dy = dyFirstSuggestion;
-
                             }
                         }
                         else
@@ -1992,13 +1928,11 @@ namespace TriangleNet
                             {
                                 dx = dxSecondSuggestion;
                                 dy = dySecondSuggestion;
-
                             }
                             else
                             {
                                 dx = dxFirstSuggestion;
                                 dy = dyFirstSuggestion;
-
                             }
                         }
 
@@ -2006,6 +1940,8 @@ namespace TriangleNet
                 }// end of relocation				 
             }// end of almostGood	
 
+            Point circumcenter = new Point();
+
             if (relocated <= 0)
             {
                 circumcenter.x = torg.x + dx;
@@ -2019,6 +1955,7 @@ namespace TriangleNet
             xi = (yao * dx - xao * dy) * (2.0 * denominator);
             eta = (xdo * dy - ydo * dx) * (2.0 * denominator);
 
+            return circumcenter;
         }
 
         /// <summary>

Triangle.NET/Triangle/Primitives.cs

@@ -19,10 +19,9 @@ namespace TriangleNet
     {
         static double splitter;       // Used to split double factors for exact multiplication.
         static double epsilon;        // Floating-point machine epsilon.
-        static double resulterrbound;
+        //static double resulterrbound;
-        static double ccwerrboundA, ccwerrboundB, ccwerrboundC;
+        static double ccwerrboundA; // ccwerrboundB, ccwerrboundC;
-        static double iccerrboundA, iccerrboundB, iccerrboundC;
+        static double iccerrboundA; // iccerrboundB, iccerrboundC;
-        static double o3derrboundA, o3derrboundB, o3derrboundC;
 
         /// <summary>
         /// Initialize the variables used for exact arithmetic.
@@ -69,16 +68,13 @@ namespace TriangleNet
             } while ((check != 1.0) && (check != lastcheck));
             splitter += 1.0;
             // Error bounds for orientation and incircle tests.
-            resulterrbound = (3.0 + 8.0 * epsilon) * epsilon;
+            //resulterrbound = (3.0 + 8.0 * epsilon) * epsilon;
             ccwerrboundA = (3.0 + 16.0 * epsilon) * epsilon;
-            ccwerrboundB = (2.0 + 12.0 * epsilon) * epsilon;
+            //ccwerrboundB = (2.0 + 12.0 * epsilon) * epsilon;
-            ccwerrboundC = (9.0 + 64.0 * epsilon) * epsilon * epsilon;
+            //ccwerrboundC = (9.0 + 64.0 * epsilon) * epsilon * epsilon;
             iccerrboundA = (10.0 + 96.0 * epsilon) * epsilon;
-            iccerrboundB = (4.0 + 48.0 * epsilon) * epsilon;
+            //iccerrboundB = (4.0 + 48.0 * epsilon) * epsilon;
-            iccerrboundC = (44.0 + 576.0 * epsilon) * epsilon * epsilon;
+            //iccerrboundC = (44.0 + 576.0 * epsilon) * epsilon * epsilon;
-            o3derrboundA = (7.0 + 56.0 * epsilon) * epsilon;
-            o3derrboundB = (3.0 + 28.0 * epsilon) * epsilon;
-            o3derrboundC = (26.0 + 288.0 * epsilon) * epsilon * epsilon;
         }
 
         /// <summary>
@@ -386,7 +382,7 @@ namespace TriangleNet
                               ref double xi, ref double eta)
         {
             double xdo, ydo, xao, yao;
-            double dodist, aodist, dadist;
+            double dodist, aodist;
             double denominator;
             double dx, dy;
 
@@ -399,8 +395,6 @@ namespace TriangleNet
             yao = tapex.y - torg.y;
             dodist = xdo * xdo + ydo * ydo;
             aodist = xao * xao + yao * yao;
-            dadist = (tdest.x - tapex.x) * (tdest.x - tapex.x) +
-                     (tdest.y - tapex.y) * (tdest.y - tapex.y);
 
             if (Behavior.NoExact)
             {

Triangle.NET/Triangle/Quality.cs

@@ -373,7 +373,7 @@ namespace TriangleNet
             double area;
             double dist1, dist2;
 
-            double maxedge, maxangle;
+            double maxangle;
 
             torg = testtri.Org();
             tdest = testtri.Dest();
@@ -463,21 +463,21 @@ namespace TriangleNet
             if ((apexlen > orglen) && (apexlen > destlen))
             {
                 // The edge opposite the apex is longest.
-                maxedge = apexlen;
+                // maxedge = apexlen;
                 // Find the cosine of the angle at the apex.
                 maxangle = (orglen + destlen - apexlen) / (2 * Math.Sqrt(orglen) * Math.Sqrt(destlen));
             }
             else if (orglen > destlen)
             {
                 // The edge opposite the origin is longest.
-                maxedge = orglen;
+                // maxedge = orglen;
                 // Find the cosine of the angle at the origin.
                 maxangle = (apexlen + destlen - orglen) / (2 * Math.Sqrt(apexlen) * Math.Sqrt(destlen));
             }
             else
             {
                 // The edge opposite the destination is longest.
-                maxedge = destlen;
+                // maxedge = destlen;
                 // Find the cosine of the angle at the destination.
                 maxangle = (apexlen + orglen - destlen) / (2 * Math.Sqrt(apexlen) * Math.Sqrt(orglen));
             }
@@ -600,7 +600,6 @@ namespace TriangleNet
             double split;
             double multiplier, divisor;
             bool acuteorg, acuteorg2, acutedest, acutedest2;
-            int dummy;
 
             // Note that steinerleft == -1 if an unlimited number
             // of Steiner points is allowed.
@@ -727,7 +726,13 @@ namespace TriangleNet
                     }
 
                     // Create the new vertex.
-                    newvertex = new Vertex(); // TODO: mesh.nextras
+                    newvertex = new Vertex(
+                        eorg.x + split * (edest.x - eorg.x),
+                        eorg.y + split * (edest.y - eorg.y),
+                        currentenc.Mark(),
+                        mesh.nextras);
+
+                    newvertex.type = VertexType.SegmentVertex;
 
                     newvertex.hash = mesh.hash_vtx++;
                     newvertex.id = newvertex.hash;
@@ -741,9 +746,6 @@ namespace TriangleNet
                             + split * (edest.attributes[i] - eorg.attributes[i]);
                     }
 
-                    newvertex.x = eorg.x + split * (edest.x - eorg.x);
-                    newvertex.y = eorg.y + split * (edest.y - eorg.y);
-
                     if (!Behavior.NoExact)
                     {
                         // Roundoff in the above calculation may yield a 'newvertex'
@@ -764,9 +766,6 @@ namespace TriangleNet
                         }
                     }
 
-                    newvertex.mark = currentenc.Mark();
-                    newvertex.type = VertexType.SegmentVertex;
-
                     // Check whether the new vertex lies on an endpoint.
                     if (((newvertex.x == eorg.x) && (newvertex.y == eorg.y)) ||
                         ((newvertex.x == edest.x) && (newvertex.y == edest.y)))
@@ -791,9 +790,9 @@ namespace TriangleNet
                         mesh.steinerleft--;
                     }
                     // Check the two new subsegments to see if they're encroached.
-                    dummy = CheckSeg4Encroach(ref currentenc);
+                    CheckSeg4Encroach(ref currentenc);
                     currentenc.NextSelf();
-                    dummy = CheckSeg4Encroach(ref currentenc);
+                    CheckSeg4Encroach(ref currentenc);
                 }
 
                 // Set subsegment's origin to NULL. This makes it possible to detect dead 
@@ -828,7 +827,7 @@ namespace TriangleNet
         {
             Otri badotri = default(Otri);
             Vertex borg, bdest, bapex;
-            Vertex newvertex;
+            Point newloc; // Location of the new vertex
             double xi = 0, eta = 0;
             InsertVertexResult success;
             bool errorflag;
@@ -846,7 +845,6 @@ namespace TriangleNet
             {
                 errorflag = false;
                 // Create a new vertex at the triangle's circumcenter.
-                newvertex = new Vertex(); // TODO: mesh.nextras
 
                 // Using the original (simpler) Steiner point location method
                 // for mesh refinement.
@@ -854,19 +852,17 @@ namespace TriangleNet
                 // reset VertexType?
                 if (behavior.FixedArea || behavior.VarArea)
                 {
-                    Point tmp = Primitives.FindCircumcenter(borg, bdest, bapex, ref xi, ref eta, behavior.Offconstant);
+                    newloc = Primitives.FindCircumcenter(borg, bdest, bapex, ref xi, ref eta, behavior.Offconstant);
-                    newvertex.x = tmp.x;
-                    newvertex.y = tmp.y;
                 }
                 else
                 {
-                    newLocation.FindLocation(borg, bdest, bapex, newvertex, ref xi, ref eta, true, badotri);
+                    newloc = newLocation.FindLocation(borg, bdest, bapex, ref xi, ref eta, true, badotri);
                 }
 
                 // Check whether the new vertex lies on a triangle vertex.
-                if (((newvertex.x == borg.x) && (newvertex.y == borg.y)) ||
+                if (((newloc.x == borg.x) && (newloc.y == borg.y)) ||
-                    ((newvertex.x == bdest.x) && (newvertex.y == bdest.y)) ||
+                    ((newloc.x == bdest.x) && (newloc.y == bdest.y)) ||
-                    ((newvertex.x == bapex.x) && (newvertex.y == bapex.y)))
+                    ((newloc.x == bapex.x) && (newloc.y == bapex.y)))
                 {
                     if (Behavior.Verbose)
                     {
@@ -876,6 +872,11 @@ namespace TriangleNet
                 }
                 else
                 {
+                    // The new vertex must be in the interior, and therefore is a
+                    // free vertex with a marker of zero.
+                    Vertex newvertex = new Vertex(newloc.x, newloc.y, 0, mesh.nextras);
+                    newvertex.type = VertexType.FreeVertex;
+
                     for (int i = 0; i < mesh.nextras; i++)
                     {
                         // Interpolate the vertex attributes at the circumcenter.
@@ -883,10 +884,6 @@ namespace TriangleNet
                             + xi * (bdest.attributes[i] - borg.attributes[i])
                             + eta * (bapex.attributes[i] - borg.attributes[i]);
                     }
-                    // The new vertex must be in the interior, and therefore is a
-                    // free vertex with a marker of zero.
-                    newvertex.type = VertexType.FreeVertex;
-                    //newvertex.mark = 0;
 
                     // Ensure that the handle 'badotri' does not represent the longest
                     // edge of the triangle.  This ensures that the circumcenter must

Triangle.NET/Triangle/Tools/QualityMeasure.cs

@@ -1,5 +1,6 @@
 // -----------------------------------------------------------------------
 // <copyright file="QualityMeasure.cs" company="">
+// Original Matlab code by John Burkardt, Florida State University
 // Triangle.NET code by Christian Woltering, http://triangle.codeplex.com/
 // </copyright>
 // -----------------------------------------------------------------------
@@ -239,6 +240,8 @@ namespace TriangleNet.Tools
         /// </remarks>
         public int Bandwidth()
         {
+            if (mesh == null) return 0;
+
             // Lower and upper bandwidth of the matrix
             int ml = 0, mu = 0;
 

Triangle.NET/Triangle/Triangle.csproj

@@ -82,6 +82,7 @@
     <Compile Include="Sampler.cs" />
     <Compile Include="Smoothing\ISmoother.cs" />
     <Compile Include="Tools\BoundedVoronoi.cs" />
+    <Compile Include="Tools\CuthillMcKee.cs" />
     <Compile Include="Tools\QualityMeasure.cs" />
     <Compile Include="Tools\Statistic.cs" />
     <Compile Include="Algorithm\SweepLine.cs" />