Examples/MeshGenToolbox/mgt/Mesh.h

+/*****************************************************************************/
+/*  PMGT                                                                     */
+/*                                                                           */
+/*  January, 2007                                                            */
+/*                                                                           */
+/*  Copyright 2007                                                           */
+/*  Wen Yu                                                                   */
+/*  Department of Computing and Software                                     */
+/*  McMaster University                                                      */
+/*  yuw4@mcmaster.ca                                                         */
+/*                                                                           */
+/*  This program may be freely redistributed under the condition that the    */
+/*    copyright notices are not removed. You may distribute modified versions*/
+/*    of this code UNDER THE CONDITION THAT THIS CODE AND ANY MODIFICATIONS  */
+/*    MADE TO IT IN THE SAME FILE REMAIN UNDER COPYRIGHT OF THE ORIGINAL     */
+/*    AUTHOR, BOTH SOURCE AND OBJECT CODE ARE MADE FREELY AVAILABLE WITHOUT  */
+/*    CHARGE, AND CLEAR NOTICE IS GIVEN OF THE MODIFICATIONS.                */
+/*                                                                           */
+/*  This is the serial version of PMGT                                       */
+/*  This file includes classes for the mesh module.                          */
+/*****************************************************************************/
+
+#ifndef MESH
+#define MESH
+#endif
+
+#ifndef VERTEX
+#include "Vertex.h"
+#endif
+
+#ifndef EDGE
+#include "Edge.h"
+#endif
+
+#ifndef CELL
+#include "Cell.h"
+#endif
+
+namespace PMGT{
+//namespace Mesh{
+
+//using namespace Vertex;
+//using namespace Edge;
+//using namespace Cell;
+
+
+typedef vector<int> Ids;
+class Mesh{
+public:
+	//constructor/destructor
+	Mesh(){}
+
+	// ***** size *****
+	inline int n_of_v(){ return vertices_.size(); }
+	inline int n_of_e(){ return edges_.size();}
+	inline int n_of_c(){ return cells_.size(); }
+
+	// ***** Id validation check *****
+	inline bool is_valid_vId(int vId){
+		return (vId < vertices_.size()) && (vId >= 0);
+	}
+
+	inline bool is_valid_hId(int hId){
+		return (hId < (edges_.size()<<1)) && (hId >= 0);
+	}
+
+	inline bool is_valid_eId(int eId){
+		return eId < edges_.size() && eId >= 0;
+	}
+
+	inline bool is_valid_cId(int cId){
+		return cId < cells_.size() && cId >= 0;
+	}
+
+
+	// ***** get entity from id *****
+	inline Vertex& vertex(int vId){
+		assert(is_valid_vId(vId));
+		return vertices_[vId];
+	}
+
+	inline Halfedge& halfedge(int hId){
+		assert(is_valid_hId(hId));
+		return edges_[hId>>1].pair()[hId & 1];
+	}
+
+	inline Edge& edge(int eId){
+		assert(is_valid_eId(eId));
+		return edges_[eId];
+	}
+
+	inline Cell& cell(int cId){
+		assert(is_valid_cId(cId));
+		return cells_[cId];
+	}
+
+	// ***** boundary check *****
+	inline bool is_b_h(int hId){
+		assert(is_valid_hId(hId));
+		return halfedge(hId).IcId() == -1;
+	}
+
+	inline bool is_b_v(int vId){
+		assert(is_valid_vId(vId));
+		if(vertex(vId).OhId() == -1) { return true;}
+		return is_b_h(ops_hId(vertex(vId).OhId()));
+	}
+
+	inline bool is_b_e(int eId){
+		assert(is_valid_eId(eId));
+		return is_b_h(eId << 1) || is_b_h((eId << 1) + 1);
+	}
+	inline bool is_b_c(int cId){
+		assert(is_valid_cId(cId));
+		int hId = cell(cId).LhId();
+		return is_b_h(ops_hId(hId)) ||
+			is_b_h(ops_hId(next_hId(hId))) ||
+			is_b_h(ops_hId(prev_hId(hId)));
+	}
+
+	// ***** connectivity *****
+
+	// ----- eId <==> hId -----
+
+	inline int eId_from_hId(int hId){
+		assert(is_valid_hId(hId));
+		return (hId >> 1);
+	}
+	inline int hId_from_eId(int eId, int index){
+		assert(is_valid_eId(eId));
+		assert( index < 2 && index >= 0 );
+		return (eId << 1) + index;
+	}
+
+	// ----- halfedge-halfedge -----
+
+	inline int ops_hId(int hId){
+		assert(is_valid_hId(hId));
+		 return ((hId & 1) ? hId - 1: hId + 1);
+	}
+	inline int next_hId(int hId){
+		assert(is_valid_hId(hId));
+		return halfedge(hId).NhId();
+	}
+	inline int prev_hId(int hId){
+		assert(is_valid_hId(hId));
+		return halfedge(hId).PhId();
+	}
+	// ccw rotated halfedge
+	inline int ccw_hId(int hId){
+		assert(is_valid_hId(hId));
+		if (is_b_h(hId)) return -1;
+		return ops_hId(prev_hId(hId));
+	}
+	// cw rotated halfedge
+	inline int cw_hId(int hId){
+		assert(is_valid_hId(hId));
+		if (is_b_h(ops_hId(hId))) return -1;
+		return next_hId(ops_hId(hId));
+	}
+
+	// ----- vertex -----
+
+	// first outgoing halfedge
+	inline int OhId_from_vId(int vId){
+		assert(is_valid_vId(vId));
+		return vertex(vId).OhId();
+	}
+
+	// ----- halfedge -----
+
+	// targate vertex
+	inline int tvId_from_hId(int hId){
+		assert(is_valid_hId(hId));
+		return halfedge(hId).TvId();
+	}
+
+	// start vertex
+	inline int svId_from_hId(int hId){
+		assert(is_valid_hId(hId));
+		return tvId_from_hId(ops_hId(hId));
+	}
+
+	// its cell
+	inline int IcId_from_hId(int hId){
+		assert(is_valid_hId(hId));
+		return halfedge(hId).IcId();
+	}
+
+	// ----- edge -----
+	// vertices
+	inline vector<int> vIds_from_eId(int eId){
+		assert(is_valid_eId(eId));
+		vector<int> vIds;
+		vIds.push_back(svId_from_hId(hId_from_eId(eId, 0)));
+		vIds.push_back(tvId_from_hId(hId_from_eId(eId, 0)));
+		return vIds;
+	}
+
+	// cells
+	inline vector<int> cIds_from_eId(int eId){
+		assert(is_valid_eId(eId));
+		vector<int> cIds;
+		cIds.push_back(IcId_from_hId(hId_from_eId(eId, 0)));
+		cIds.push_back(IcId_from_hId(hId_from_eId(eId, 1)));
+		return cIds;
+	}
+
+	// ----- cell -----
+
+	// its longest halfedge
+	inline int LhId_from_cId(int cId){
+		assert(is_valid_cId(cId));
+		return cell(cId).LhId();
+	}
+
+	// three vertices
+	inline vector<int> vIds_from_cId(int cId){
+		assert(is_valid_cId(cId));
+		vector<int> vIds;
+		vIds.push_back(tvId_from_hId(cell(cId).LhId()));
+		vIds.push_back(tvId_from_hId(next_hId(cell(cId).LhId())));
+		vIds.push_back(svId_from_hId(cell(cId).LhId()));
+		return vIds;
+	}
+
+	// ***** isolation test *****
+
+	// is the vertex is isolated
+	inline bool is_isolated_v(int vId){
+		assert(is_valid_vId(vId));
+		return !is_valid_hId(OhId_from_vId(vId));
+	}
+
+	// ***** circulators *****
+
+	// ----- vertex - vertex (ccw) one ring circulator -----
+
+	inline int vv_first(int vId){
+		assert(is_valid_vId(vId));
+		return tvId_from_hId(OhId_from_vId(vId));
+	}
+
+	// Returns -1 if there is no next
+	inline int vv_next(int svId, int tvId){
+		assert(is_valid_vId(svId));
+		assert(is_valid_vId(tvId));
+		int hId = voh_first(svId);
+		while (tvId_from_hId(hId) != tvId && hId != -1){
+			hId = voh_next(svId, hId);
+		}
+		if (hId != -1){
+			if (is_b_h(hId)) return -1;
+			return tvId_from_hId(next_hId(hId));
+		}
+		assert(false);
+	}
+
+	// ----- vertex - outgoing-halfedge (ccw) one ring circulator -----
+	inline int voh_first(int vId){
+		assert(is_valid_vId(vId));
+		return OhId_from_vId(vId);
+	}
+
+	// Returns -1 if there is no next
+	inline int voh_next(int vId, int hId){
+		assert(is_valid_vId(vId));
+		assert(is_valid_hId(hId));
+		assert(vId == svId_from_hId(hId));
+		return ccw_hId(hId);
+	}
+
+	// ----- vertex - incoming-halfedge (ccw) one ring circulator -----
+	inline int vih_first(int vId){
+		assert(is_valid_vId(vId));
+		if ( is_isolated_v(vId)) return -1;
+		return ops_hId(voh_first(vId));
+	}
+
+	// Returns -1 if there is no next
+	inline int vih_next(int vId, int hId){
+		assert(is_valid_vId(vId));
+		assert(is_valid_hId(hId));
+		assert(vId == tvId_from_hId(hId));
+		int ops = ops_hId(hId);
+		if (is_b_h(ops)) return -1;
+		//return prev_hId(ops);
+		return  ops_hId(voh_next(vId, ops));
+	}
+
+	// ----- vertex - edge (ccw) one ring circulator -----
+	inline int ve_first(int vId){
+		assert(is_valid_vId(vId));
+		if (is_isolated_v(vId)) return -1;
+		return eId_from_hId(voh_first(vId));
+	}
+
+	// Return -1 if there is no next
+	inline int ve_next(int vId, int eId){
+		assert(is_valid_vId(vId));
+		assert(is_valid_eId(eId));
+		int hId = hId_from_eId(eId,0);
+		if (vId == svId_from_hId(hId)){ // hId is an outgoing halfedge
+			int voNext = voh_next(vId, hId);
+			if (voNext == -1) return -1;
+			return eId_from_hId(voNext);
+		}
+		//hId is an incoming halfedge
+		assert(vId == tvId_from_hId(hId));
+		int viNext = vih_next(vId, hId);
+		if (viNext == -1) return -1;
+		return eId_from_hId(viNext);
+	}
+
+	// ----- vertex - cell (cww) one ring circulator -----
+
+	inline int vc_first(int vId){
+		assert(is_valid_vId(vId));
+		if (is_isolated_v(vId)) return -1;
+		return IcId_from_hId(voh_first(vId));
+	}
+
+	// Returns -1 if there is no next
+	inline int vc_next(int vId, int cId){
+		assert(is_valid_vId(vId));
+		assert(is_valid_cId(cId));
+		// find the halfedge whose targate vertex is vId
+		int hId = LhId_from_cId(cId);
+		if (tvId_from_hId(hId) != vId){
+			if (tvId_from_hId(prev_hId(hId)) == vId)
+				hId = prev_hId(hId);
+			else if (tvId_from_hId(next_hId(hId)) == vId)
+				hId = next_hId(hId);
+			else assert(false);
+		}
+		if (is_b_h(ops_hId(hId))) return -1;
+		return IcId_from_hId(ops_hId(hId));
+	}
+
+	// ----- cell - vertex (cww) one ring circulator -----
+	// The first vertex is the tartage vertex of the longest halfedge
+	inline int cv_first(int cId){
+		assert(is_valid_cId(cId));
+		return tvId_from_hId(LhId_from_cId(cId));
+	}
+
+	inline int cv_next(int cId, int vId){
+		assert(is_valid_cId(cId));
+		assert(is_valid_vId(vId));
+		// find the halfedge whose start vertex is vId
+		int hId = LhId_from_cId(cId);
+		if (svId_from_hId(hId) != vId){
+			if (svId_from_hId(prev_hId(hId)) == vId)
+				hId = prev_hId(hId);
+			else if (svId_from_hId(next_hId(hId)) == vId)
+				hId = next_hId(hId);
+			else assert(false);
+		}
+		return tvId_from_hId(hId);
+	}
+
+	// ----- cell - halfedge (cww) one ring circulator -----
+	// The first halfedge is the longest halfedge
+	inline int ch_first(int cId){
+		assert(is_valid_cId(cId));
+		return LhId_from_cId(cId);
+	}
+	inline int ch_next(int cId, int hId){
+		assert(is_valid_cId(cId));
+		assert(is_valid_hId(hId));
+		assert(IcId_from_hId(hId) == cId);
+		return next_hId(hId);
+	}
+
+	// ----- cell - edge (cww) one ring circulator -----
+	// The first edge is the longest edge
+	inline int ce_first(int cId){
+		assert(is_valid_cId(cId));
+		return eId_from_hId(ch_first(cId));
+	}
+	// next should change
+	inline int ce_next(int cId, int eId){
+		assert(is_valid_cId(cId));
+		assert(is_valid_eId(eId));
+		int hId = hId_from_eId(eId, 0);
+		if (IcId_from_hId(hId) == cId)
+			return eId_from_hId(ch_next(cId,hId));
+		else if (IcId_from_hId(ops_hId(hId)) == cId)
+			return eId_from_hId(ch_next(cId, ops_hId(hId)));
+		else assert(false);
+	}
+
+	// ----- cell - cell (cww) one ring circulator -----
+
+	// The first cell is the cell adjacent the longest edge
+	// if the given cell is not a boundary cell.
+	// If one edge is boundary edge, returns the cell adjacent to the next edge (ccw).
+	// If two edges are boundary edges, refturns the only adjacent cell.
+	// Returns -1 if there is no adjacent cell.
+	inline int cc_first(int cId){
+		assert(is_valid_cId(cId));
+		int hId = LhId_from_cId(cId);
+		if (!is_b_c(cId)) return IcId_from_hId(ops_hId(hId));
+
+		if (!is_b_h(ops_hId(hId))){
+			if ( is_b_h(ops_hId(next_hId(hId))) && !is_b_h(ops_hId(prev_hId(hId))) )
+				return IcId_from_hId(ops_hId(prev_hId(hId)));
+			else return IcId_from_hId(ops_hId(hId));
+		}
+		else {
+			if ( !is_b_h(ops_hId(next_hId(hId))) )
+				return IcId_from_hId(ops_hId(next_hId(hId)));
+			else return IcId_from_hId(ops_hId(prev_hId(hId)));
+		}
+	}
+
+	// Returns -1 if there is no next
+	inline int cc_next(int scId, int tcId){
+		assert(is_valid_cId(tcId));
+		assert(is_valid_cId(scId));
+		int hId = LhId_from_cId(scId);
+
+		if (IcId_from_hId(ops_hId(hId)) != tcId){
+			if (IcId_from_hId(ops_hId(next_hId(hId))) == tcId)
+				hId = next_hId(hId);
+			else if (IcId_from_hId(ops_hId(prev_hId(hId))) == tcId)
+				hId = prev_hId(hId);
+			else assert(false);
+		}
+		return IcId_from_hId(ops_hId(next_hId(hId)));
+	}
+
+
+
+	// ***** find an entity *****
+
+	// Returns the id of the vertex that has given coordinate.
+	// Returns -1 if the vertex does not exist.
+	int find_vertex(Point point);
+
+	// Returns the id of the edge that has given points.
+	// Returns -1 if the edge does not exist.
+	int find_edge(int svId, int tvId);
+
+	// Returns the id of the halfedge that has given start and targate points.
+	// Returns -1 if the halfedge does not exist.
+	int find_halfedge(int svId, int tvId);
+
+	// Returns the id of the cell that has given vertices.
+	// Returns -1 if the cell does not exist.
+	int find_cell(vector<int> vIds);
+
+
+	// ***** Geometric properties *****
+
+	// Returns the distance of two vertex
+	inline double dist(int svId, int tvId){
+		assert(is_valid_vId(svId) && is_valid_vId(tvId));
+		Point 	ps = vertex(svId).point(),
+			pt = vertex(tvId).point();
+		return sqrt((ps[0]-pt[0])*(ps[0]-pt[0]) +
+				(ps[1]-pt[1])*(ps[1]-pt[1]));
+	}
+
+	inline double dist(int hId){
+		assert(is_valid_hId(hId));
+		return dist(svId_from_hId(hId), tvId_from_hId(hId));
+	}
+
+	// ----- mid point of a halfedge
+	inline Point mid_point(int hId){
+		assert(is_valid_hId(hId));
+		double sx = vertex(svId_from_hId(hId)).x();
+		double sy = vertex(svId_from_hId(hId)).y();
+		double tx = vertex(tvId_from_hId(hId)).x();
+		double ty = vertex(tvId_from_hId(hId)).y();
+		Point mid;
+		mid.push_back((sx+tx)/2);
+		mid.push_back((sy+ty)/2);
+		return mid;
+	}
+
+	// ----- mid point of two points
+	inline Point mid_point(int svId, int tvId){
+		double sx = vertex(svId).x();
+		double sy = vertex(svId).y();
+		double tx = vertex(tvId).x();
+		double ty = vertex(tvId).y();
+		Point mid;
+		mid.push_back((sx+tx)/2);
+		mid.push_back((sy+ty)/2);
+		return mid;
+	}
+
+	// ----- mid point of the longest halfedge of a cell
+	inline Point longest_mid(int cId){
+		assert(is_valid_cId(cId));
+		return mid_point(LhId_from_cId(cId));
+	}
+
+	// ----- centroid point of a cell -----
+	inline Point centroid(int cId){
+		assert(is_valid_cId(cId));
+		Point point = longest_mid(cId);
+		// find the opposite point to the longest halfedge
+		Vertex& v = vertex(tvId_from_hId(next_hId(LhId_from_cId(cId))));
+		point[0] = point[0]*2/3 + v.x()/3;
+		point[1] = point[1]*2/3 + v.y()/3;
+
+		return point;
+	}
+
+
+	// ----- size of a cell ----- //
+	// The size of a cell if the length of longest halfedge of the cell
+	inline double size(int cId){
+		assert(is_valid_cId(cId));
+		return dist(LhId_from_cId(cId));
+	}
+
+	// ----- Test for counter clockwise ----- //
+	inline bool is_ccw(Point pa, Point pb, Point pc){
+		double detleft = (pa[0] - pc[0]) * (pb[1] - pc[1]);
+		double detright = (pa[1] - pc[1]) * (pb[0] - pc[0]);
+  		return detleft - detright > DBL_EPSILON*100;
+	}
+
+	inline bool is_ccw(vector<int> vIds){
+		assert(vIds.size() == 3);
+		assert(is_valid_vId(vIds[0]));
+		assert(is_valid_vId(vIds[1]));
+		assert(is_valid_vId(vIds[2]));
+		return is_ccw(vertex(vIds[0]).point(),
+						vertex(vIds[1]).point(),
+						vertex(vIds[2]).point());
+	}
+
+/*	// ----- Test for collinear ------ //
+	inline bool is_collinear(Point pa, Point pb, Point pc){
+		double detleft = (pa[0] - pc[0]) * (pb[1] - pc[1]);
+		double detright = (pa[1] - pc[1]) * (pb[0] - pc[0]);
+  		return detleft - detright == 0;
+	}
+*/
+	// ----- Test for the same side ----- //
+
+	inline bool is_same_side(int hId, Point pa, Point pb){
+		Point ps = vertex(svId_from_hId(hId)).point();
+		Point pt = vertex(tvId_from_hId(hId)).point();
+		double deta = (pa[1]-ps[1])*(pt[0]-ps[0]) -
+					(pa[0]-ps[0])*(pt[1]-ps[1]);
+		double detb = (pb[1]-ps[1])*(pt[0]-ps[0]) -
+					(pb[0]-ps[0])*(pt[1]-ps[1]);
+		return deta*detb > 0;
+	}
+
+/*	inline bool is_in_halfedge(int hId, Point p){
+		Point ps = vertex(svId_from_hId(hId)).point();
+		Point pt = vertex(tvId_from_hId(hId)).point();
+		return is_collinear(ps, pt, p);
+	}
+*/
+	/** \ Add an entity
+	*/
+
+	// Create a vertex and add it to the end of the vertex list.
+	// Returns false if the vertex already exists.
+	bool add_vertex(Point point);
+
+
+	// Create an edge and add it to the end of the edge list.
+	// Returns false if the edge already exists.
+	bool add_edge(int svId, int tvId);
+
+
+	// Create a cell and add it to the end of the cell list.
+	// Returns false if the cell already exists.
+	bool add_cell(vector<int> vIds);
+
+	/** \ Update an entity
+	*/
+
+	// update a cell using three vertices.
+	// Returns false if the new cell already exists.
+	bool update_cell(vector<int> vIds, int cId);
+
+
+	// split 2 cells by split their common longest edge
+
+	bool longest_ref2(int hId, int o_hId, int next_hId);
+
+	// ****** delete an entity *****
+
+	// Delete vertices that marked to be deleted from the list.
+	bool delete_vertices();
+
+	// Delete edges that marked to be deleted from the list.
+	bool delete_edges();
+
+	// Delete cells that marked to be deleted from the list.
+	bool delete_cells();
+
+	// ***** create a mesh ******
+	// create a mesh using two input variables.
+	// One is the coordinates of points
+	// The other is the vertex list of cells
+	bool create(vector<Point>, vector<Ids>);
+
+	// ***** write a mesh *****
+	// write the mesh to two variables
+	// One is the coordinates of points
+	// The other is the vertex list of cells
+	void write_var(vector<Point>&, vector<Ids>&);
+
+
+
+	// display the internal information of a mesh on the screen
+	void display_info();
+
+
+private:
+	vector<Vertex> 		vertices_;
+	vector<Edge>		edges_;
+	vector<Cell>		cells_;
+
+
+};
+
+
+//}
+
+}

Examples/MeshGenToolbox/mgt/Output.cpp

+/*****************************************************************************/
+/*  PMGT                                                                     */
+/*                                                                           */
+/*  January, 2007                                                            */
+/*                                                                           */
+/*  Copyright 2007                                                           */
+/*  Wen Yu                                                                   */
+/*  Department of Computing and Software                                     */
+/*  McMaster University                                                      */
+/*  yuw4@mcmaster.ca                                                         */
+/*                                                                           */
+/*  This program may be freely redistributed under the condition that the    */
+/*    copyright notices are not removed. You may distribute modified versions*/
+/*    of this code UNDER THE CONDITION THAT THIS CODE AND ANY MODIFICATIONS  */
+/*    MADE TO IT IN THE SAME FILE REMAIN UNDER COPYRIGHT OF THE ORIGINAL     */
+/*    AUTHOR, BOTH SOURCE AND OBJECT CODE ARE MADE FREELY AVAILABLE WITHOUT  */
+/*    CHARGE, AND CLEAR NOTICE IS GIVEN OF THE MODIFICATIONS.                */
+/*                                                                           */
+/*  This is the serial version of PMGT                                       */
+/*****************************************************************************/
+
+#include "Output.h"
+
+namespace PMGT{
+	// write the mesh to two files
+	// One is the coordinates of points
+	// The other is the vertex list of cells
+	bool Output::write_file(Mesh& mesh, char* vertexFile, char* cellFile){
+		FILE *v;
+		FILE *c;
+		int i;
+		vector<int> vIds(3);
+
+  		/* Open files.  If NULL is returned there was an error */
+  		if((v = fopen(vertexFile, "w")) == NULL) {
+    			printf("Error Opening Vertex File.\n");
+    			exit(1);
+  		}
+ 		if((c = fopen(cellFile, "w")) == NULL) {
+    			printf("Error Opening Cell File.\n");
+    			exit(1);
+  		}
+
+		// write vertex file
+		for (i=0; i<mesh.n_of_v(); i++){
+			fprintf(v, "%.6g %.6g\n",
+				mesh.vertex(i).x(), mesh.vertex(i).y());
+		}
+		// write cell file
+		for (i=0; i<mesh.n_of_c(); i++){
+			vIds[0] = mesh.tvId_from_hId(mesh.LhId_from_cId(i));
+			vIds[1] = mesh.tvId_from_hId(mesh.next_hId(mesh.LhId_from_cId(i)));
+			vIds[2] = mesh.tvId_from_hId(mesh.prev_hId(mesh.LhId_from_cId(i)));
+			fprintf(c, "%d %d %d\n", vIds[0], vIds[1], vIds[2]);
+		}
+		fclose(v);
+		fclose(c);
+	}
+
+}

Examples/MeshGenToolbox/mgt/Output.h

+/*****************************************************************************/
+/*  PMGT                                                                     */
+/*                                                                           */
+/*  January, 2007                                                            */
+/*                                                                           */
+/*  Copyright 2007                                                           */
+/*  Wen Yu                                                                   */
+/*  Department of Computing and Software                                     */
+/*  McMaster University                                                      */
+/*  yuw4@mcmaster.ca                                                         */
+/*                                                                           */
+/*  This program may be freely redistributed under the condition that the    */
+/*    copyright notices are not removed. You may distribute modified versions*/
+/*    of this code UNDER THE CONDITION THAT THIS CODE AND ANY MODIFICATIONS  */
+/*    MADE TO IT IN THE SAME FILE REMAIN UNDER COPYRIGHT OF THE ORIGINAL     */
+/*    AUTHOR, BOTH SOURCE AND OBJECT CODE ARE MADE FREELY AVAILABLE WITHOUT  */
+/*    CHARGE, AND CLEAR NOTICE IS GIVEN OF THE MODIFICATIONS.                */
+/*                                                                           */
+/*  This is the serial version of PMGT                                       */
+/*  This file includes classes for the output module.                        */
+/*****************************************************************************/
+
+#ifndef OUTPUT
+#define OUTPUT
+#endif
+
+#ifndef MESH
+#include "Mesh.h"
+#endif
+
+namespace PMGT{
+class Output{
+public:
+	// write the mesh to two files
+	// One is the coordinates of points
+	// The other is the vertex list of cells
+	bool write_file(Mesh& mesh, char* vertexFile, char* cellFile);
+};
+}
+
+

Examples/MeshGenToolbox/mgt/Refining.cpp

+/*****************************************************************************/
+/*  PMGT                                                                     */
+/*                                                                           */
+/*  January, 2007                                                            */
+/*                                                                           */
+/*  Copyright 2007                                                           */
+/*  Wen Yu                                                                   */
+/*  Department of Computing and Software                                     */
+/*  McMaster University                                                      */
+/*  yuw4@mcmaster.ca                                                         */
+/*                                                                           */
+/*  This program may be freely redistributed under the condition that the    */
+/*    copyright notices are not removed. You may distribute modified versions*/
+/*    of this code UNDER THE CONDITION THAT THIS CODE AND ANY MODIFICATIONS  */
+/*    MADE TO IT IN THE SAME FILE REMAIN UNDER COPYRIGHT OF THE ORIGINAL     */
+/*    AUTHOR, BOTH SOURCE AND OBJECT CODE ARE MADE FREELY AVAILABLE WITHOUT  */
+/*    CHARGE, AND CLEAR NOTICE IS GIVEN OF THE MODIFICATIONS.                */
+/*                                                                           */
+/*  This is the serial version of PMGT                                       */
+/*****************************************************************************/
+#ifndef REFINING
+#include "Refining.h"
+#endif
+
+namespace PMGT{
+	// split cells at their centroid points
+	void Refining::split(Mesh& mesh){
+		vector<int> vIds;
+		vector<int> newIds;
+		int i,j, newvId;
+		Point p;
+
+		for (i=0; i<mesh.n_of_c(); i++){
+			if (mesh.cell(i).mFlag() == REFINE){
+				// get three vertices
+				vIds.clear();
+				vIds.push_back(mesh.cv_first(i));
+				vIds.push_back(mesh.cv_next(i, vIds[0]));
+				vIds.push_back(mesh.cv_next(i, vIds[1]));
+				// add the centroid point
+				Point np = mesh.centroid(i);
+				mesh.add_vertex(np);
+				newvId = mesh.find_vertex(np);
+				assert(newvId !=-1);
+				// add three new edges
+				for (j=0; j<3; j++){
+					mesh.add_edge(newvId, vIds[j]);
+				}
+				// add/update three cells
+				for (j=0; j<3; j++){
+					newIds.clear();
+					newIds.push_back(newvId);
+					newIds.push_back(vIds[j]);
+					newIds.push_back(vIds[(j+1)%3]);
+					if (j==0){
+						mesh.update_cell(newIds, i);
+					}
+					else
+						mesh.add_cell(newIds);
+				}
+			}// if
+		}
+
+	}
+
+
+
+	void Refining::longest_edge(Mesh& mesh){
+
+		for(int i=0; i<mesh.n_of_c(); i++){
+			if (mesh.cell(i).mFlag() == REFINE){
+				int LhId = mesh.cell(i).LhId();
+				int ops_LhId = mesh.ops_hId(LhId);
+				int next_cell = mesh.IcId_from_hId(ops_LhId);
+
+				vector<int> lepp;
+				lepp.clear();
+
+				while(1){
+					lepp.push_back(LhId);
+					lepp.push_back(ops_LhId);
+					if (next_cell == -1){
+						break;
+					}
+					if (mesh.LhId_from_cId(next_cell) == ops_LhId){
+						break;
+					}
+
+					LhId = mesh.cell(next_cell).LhId();;
+					ops_LhId = mesh.ops_hId(LhId);
+					next_cell = mesh.IcId_from_hId(ops_LhId);
+				} // while(1)
+
+				assert(lepp.size()%2 == 0);
+
+				// refine LEPP
+				while (lepp.size() > 0){
+					int hId = lepp[lepp.size()-1];
+					lepp.pop_back();
+					assert(lepp.size() > 0);
+					int ops_hId = lepp[lepp.size()-1];
+					lepp.pop_back();
+
+					assert(ops_hId == mesh.ops_hId(hId));
+					int next_hId = -1;
+					if (lepp.size() > 0){
+						next_hId = lepp[lepp.size()-1];
+					}
+
+					mesh.longest_ref2(hId, ops_hId, next_hId);
+				}// while
+
+			} // if REFINE
+		}// for i
+
+	}
+
+
+
+}

Examples/MeshGenToolbox/mgt/Refining.h

+/*****************************************************************************/
+/*  PMGT                                                                     */
+/*                                                                           */
+/*  January, 2007                                                            */
+/*                                                                           */
+/*  Copyright 2007                                                           */
+/*  Wen Yu                                                                   */
+/*  Department of Computing and Software                                     */
+/*  McMaster University                                                      */
+/*  yuw4@mcmaster.ca                                                         */
+/*                                                                           */
+/*  This program may be freely redistributed under the condition that the    */
+/*    copyright notices are not removed. You may distribute modified versions*/
+/*    of this code UNDER THE CONDITION THAT THIS CODE AND ANY MODIFICATIONS  */
+/*    MADE TO IT IN THE SAME FILE REMAIN UNDER COPYRIGHT OF THE ORIGINAL     */
+/*    AUTHOR, BOTH SOURCE AND OBJECT CODE ARE MADE FREELY AVAILABLE WITHOUT  */
+/*    CHARGE, AND CLEAR NOTICE IS GIVEN OF THE MODIFICATIONS.                */
+/*                                                                           */
+/*  This is the serial version of PMGT                                       */
+/*  This file contains the class for refine module.                          */
+/*  There are two algorithms:                                                */
+/*    split -- the point insertion algorithm for refining                    */
+/*    longest_edge -- the longest side bisection algorithm for refining      */
+/*****************************************************************************/
+
+#ifndef REFINING
+#define REFINING
+#endif
+
+#ifndef MESH
+#include "Mesh.h"
+#endif
+
+namespace PMGT{
+class Refining{
+public:
+	// split cells at their centroid points
+	void split(Mesh& mesh);
+
+	// refine cells at their longest edges
+	void longest_edge(Mesh& mesh);
+
+};
+}

Examples/MeshGenToolbox/mgt/Service.cpp

+/*****************************************************************************/
+/*  PMGT                                                                     */
+/*                                                                           */
+/*  January, 2007                                                            */
+/*                                                                           */
+/*  Copyright 2007                                                           */
+/*  Wen Yu                                                                   */
+/*  Department of Computing and Software                                     */
+/*  McMaster University                                                      */
+/*  yuw4@mcmaster.ca                                                         */
+/*                                                                           */
+/*  This program may be freely redistributed under the condition that the    */
+/*    copyright notices are not removed. You may distribute modified versions*/
+/*    of this code UNDER THE CONDITION THAT THIS CODE AND ANY MODIFICATIONS  */
+/*    MADE TO IT IN THE SAME FILE REMAIN UNDER COPYRIGHT OF THE ORIGINAL     */
+/*    AUTHOR, BOTH SOURCE AND OBJECT CODE ARE MADE FREELY AVAILABLE WITHOUT  */
+/*    CHARGE, AND CLEAR NOTICE IS GIVEN OF THE MODIFICATIONS.                */
+/*                                                                           */
+/*  This program used for correctness test of PMGT.                          */
+/*****************************************************************************/
+
+#include "Service.h"
+
+namespace PMGT{
+
+// auxiliary functions
+vector<int> Service::boundary_list(Mesh& mesh){
+	vector<int> list;
+	vector<int> b_edges;
+	vector<int>::iterator tempIterator;
+
+	for(int i=0; i<mesh.n_of_e(); i++){
+		if (mesh.is_b_e(i)){
+			b_edges.push_back(i);
+		}
+	}
+	vector<int> vIds;
+	vIds = mesh.vIds_from_eId(b_edges[0]);
+	list.push_back(vIds[0]);
+	//list.push_back(vIds[1]);
+	int first_vId = vIds[0];
+	int current = vIds[1];
+	b_edges[0] = b_edges[b_edges.size()-1];
+	b_edges.pop_back();
+	bool error = false, found;
+
+	int next_vertex, i;
+	while(b_edges.size()!=0){
+		found = false;
+		for( i=0; i<b_edges.size(); i++){
+			vIds = mesh.vIds_from_eId(b_edges[i]);
+			for(int j=0; j<2; j++){
+				if(vIds[j] == current){ // find next boundary edge
+					list.push_back(vIds[j]);
+					b_edges[i] = b_edges[b_edges.size()-1];
+					b_edges.pop_back();
+					current = vIds[(j+1)%2];
+					found = true;
+					break;
+				}
+			}
+			if(found) break;
+		}
+		if (!found){ // the boundary is not close
+			error = true;
+			break;
+		}
+		else if(first_vId==current){// find a close boundary
+			if(b_edges.size()!=0){ // there are boundary edges left
+				error = true;
+			}
+			break; //find an unique close boundary
+		}
+	}
+
+	if (error){
+		list.clear();
+	}
+	return list;
+}
+
+
+bool Service::in_boundary(Mesh& mesh1, int vId1, Mesh& mesh2, vector<int> b_list2){
+	double x = mesh1.vertex(vId1).x();
+	double y = mesh1.vertex(vId1).x();
+	for (int i=0; i<b_list2.size()-2; i++){
+		double xa = mesh2.vertex(b_list2[i]).x();
+		double ya = mesh2.vertex(b_list2[i]).y();
+		double xb = mesh2.vertex(b_list2[i+1]).x();
+		double yb = mesh2.vertex(b_list2[i+1]).y();
+		double xc = mesh2.vertex(b_list2[i+2]).x();
+		double yc = mesh2.vertex(b_list2[i+2]).y();
+
+		double detAB = (y-ya)*(xb-xa) - (x-xa)*(yb-ya);
+		double detBC = (y-yb)*(xc-xb) - (x-xb)*(yc-yb);
+
+		if (detAB*detBC < 0) return false;
+	}
+	return true;
+}
+
+
+double Service::area(Mesh& mesh, int vId1, int vId2, int vId3){
+	double x1 = mesh.vertex(vId1).x();
+	double y1 = mesh.vertex(vId1).y();
+	double x2 = mesh.vertex(vId2).x();
+	double y2 = mesh.vertex(vId2).y();
+	double x3 = mesh.vertex(vId3).x();
+	double y3 = mesh.vertex(vId3).y();
+
+	return abs((x1*y2-x2*y1+x1*y3-x3*y1+x2*y3-x3*y2)/2.0);
+}
+
+
+bool Service::inside(Mesh& mesh, int vId, int cId){
+	double x = mesh.vertex(vId).x();
+	double y = mesh.vertex(vId).y();
+
+	vector<int> vIds = mesh.vIds_from_cId(cId);
+	double x1 = mesh.vertex(vIds[0]).x();
+	double y1 = mesh.vertex(vIds[0]).y();
+	double x2 = mesh.vertex(vIds[1]).x();
+	double y2 = mesh.vertex(vIds[1]).y();
+	double x3 = mesh.vertex(vIds[2]).x();
+	double y3 = mesh.vertex(vIds[2]).y();
+
+	double det12 = (y-y1)*(x2-x1) - (x-x1)*(y2-y1);
+	double det23 = (y-y2)*(x3-x2) - (x-x2)*(y3-y2);
+	double det31 = (y-y3)*(x1-x3) - (x-x3)*(y1-y3);
+	//(x,y) in the same side of 12, 23, 31
+	return (det12*det23>0)&&(det23*det31>0);
+
+}
+
+// http://local.wasp.uwa.edu.au/~pbourke/geometry/lineline2d
+bool Service::across(Mesh& mesh, int eId1, int eId2){
+
+	vector<int> vIds1 = mesh.vIds_from_eId(eId1);
+	double x1,y1,x2,y2,x3,y3,x4,y4;
+	x1 = mesh.vertex(vIds1[0]).x();
+	y1 = mesh.vertex(vIds1[0]).y();
+	x2 = mesh.vertex(vIds1[1]).x();
+	y2 = mesh.vertex(vIds1[1]).y();
+
+	vector<int> vIds2 = mesh.vIds_from_eId(eId2);
+	x3 = mesh.vertex(vIds2[0]).x();
+	y3 = mesh.vertex(vIds2[0]).y();
+	x4 = mesh.vertex(vIds2[1]).x();
+	y4 = mesh.vertex(vIds2[1]).y();
+
+	double det312 = (y3-y1)*(x2-x1) - (x3-x1)*(y2-y1);
+	double det412 = (y4-y1)*(x2-x1) - (x4-x1)*(y2-y1);
+
+	double det134 = (y1-y3)*(x4-x3) - (x1-x3)*(y4-y3);
+	double det234 = (y2-y3)*(x4-x3) - (x2-x3)*(y4-y3);
+
+
+	// 3 and 4 in the different side of 12 and
+	// 1 and 2 in the different side of 34
+	return (det312*det412<0 && det134*det234<0);
+}
+
+
+// functions test correctness
+bool Service::length_gt0(Mesh& mesh){
+	for (int i=0; i<mesh.n_of_e(); i++){
+		if (mesh.dist(mesh.hId_from_eId(i,0)) < 0){
+			cout << "length_gt0 is false." << endl;
+			return false;
+		}
+	}
+	return true;
+}
+
+bool Service::areas_gt0(Mesh& mesh){
+	for (int i=0; i<mesh.n_of_c(); i++){
+		vector<int> vIds = mesh.vIds_from_cId(i);
+		if (area(mesh, vIds[0], vIds[1], vIds[2]) <= 0){
+			cout << "areas_gt0 is false." << endl;
+			return false;
+		}
+	}
+	return true;
+}
+
+bool Service::bounded(Mesh& mesh){
+	return (boundary_list(mesh).size()!=0);
+}
+
+
+bool Service::no_interior_intersection(Mesh& mesh){
+	for (int i=0; i<mesh.n_of_v(); i++){
+		for (int j=0; j<mesh.n_of_c(); j++){
+			if (inside(mesh, i, j)){
+				cout << "no_interior_intersection is false." << endl;
+				return false;
+			}
+		}
+	}
+
+	for (int i=0; i<mesh.n_of_e(); i++){
+		for (int j=0; j<mesh.n_of_e(); j++){
+			if (across(mesh, i, j)){
+				cout << "no_interior_intersection is false." << endl;
+				return false;
+			}
+		}
+	}
+	return true;
+}
+
+
+bool Service::conformal(Mesh& mesh){
+	for (int i=0; i<mesh.n_of_v(); i++){
+		for (int j=0; j<mesh.n_of_e(); j++){
+			vector<int> vIds = mesh.vIds_from_eId(j);
+			double x = mesh.vertex(i).x();
+			double y = mesh.vertex(i).y();
+			double x1 = mesh.vertex(vIds[0]).x();
+			double y1 = mesh.vertex(vIds[0]).y();
+			double x2 = mesh.vertex(vIds[1]).x();
+			double y2 = mesh.vertex(vIds[1]).y();
+			if((x<x1&&y<y1 && x>x2&&y>y2) || (x>x1&&y>y1 && x<x2&&y<y2)){
+				if (area(mesh,i,vIds[0],vIds[1])<10E-10){
+					cout << "conformal is false." << endl;
+					return false;
+				}
+			}
+		}
+	}
+	return true;
+
+}
+
+
+bool Service::counterclockwise(Mesh& mesh){
+	for (int i=0; i<mesh.n_of_c(); i++){
+		vector<int> vIds = mesh.vIds_from_cId(i);
+		Vertex v1 = mesh.vertex(vIds[0]);
+		Vertex v2 = mesh.vertex(vIds[1]);
+		Vertex v3 = mesh.vertex(vIds[2]);
+		double detleft = (v1.x()-v3.x())*(v2.y()-v3.y());
+		double detright = (v1.y()-v3.y())*(v2.x()-v3.x());
+		//return detleft - detright > DBL_EPSILON*100;
+		if ( detleft < detright + DBL_EPSILON*100){
+			cout << "counterclockwise is false."<< endl;
+			return false;
+		}
+	}
+	return true;
+}
+
+bool Service::covering_up(Mesh& mesh1, Mesh& mesh2){
+	vector<int> b_list = boundary_list(mesh1);
+	for (int i=0; i<mesh2.n_of_v(); i++){
+		if (!in_boundary(mesh2, i, mesh1, b_list)){
+			return false;
+		}
+	}
+	b_list = boundary_list(mesh2);
+	for (int i=0; i<mesh1.n_of_v(); i++){
+		if (!in_boundary(mesh1, i, mesh2, b_list)){
+			return false;
+		}
+	}
+
+	return true;
+}
+
+
+}
+

Examples/MeshGenToolbox/mgt/Service.h

+/*****************************************************************************/
+/*  PMGT                                                                     */
+/*                                                                           */
+/*  January, 2007                                                            */
+/*                                                                           */
+/*  Copyright 2007                                                           */
+/*  Wen Yu                                                                   */
+/*  Department of Computing and Software                                     */
+/*  McMaster University                                                      */
+/*  yuw4@mcmaster.ca                                                         */
+/*                                                                           */
+/*  This program may be freely redistributed under the condition that the    */
+/*    copyright notices are not removed. You may distribute modified versions*/
+/*    of this code UNDER THE CONDITION THAT THIS CODE AND ANY MODIFICATIONS  */
+/*    MADE TO IT IN THE SAME FILE REMAIN UNDER COPYRIGHT OF THE ORIGINAL     */
+/*    AUTHOR, BOTH SOURCE AND OBJECT CODE ARE MADE FREELY AVAILABLE WITHOUT  */
+/*    CHARGE, AND CLEAR NOTICE IS GIVEN OF THE MODIFICATIONS.                */
+/*                                                                           */
+/*  This program is used for correctness test of PMGT.                       */
+/*  This file includes classes for the service module.                       */
+/*****************************************************************************/
+
+#ifndef SERVICE
+#define SERVICE
+#endif
+
+#ifndef VERTEX
+#include "Vertex.h"
+#endif
+
+#ifndef EDGE
+#include "Edge.h"
+#endif
+
+#ifndef CELL
+#include "Cell.h"
+#endif
+
+#ifndef MESH
+#include "Mesh.h"
+#endif
+
+namespace PMGT{
+class Service{
+public:
+
+	// auxiliary functions
+	vector<int> boundary_list(Mesh& mesh);
+	bool in_boundary(Mesh& mesh1, int vId1, Mesh& mesh2, vector<int> b_list2);
+	double area(Mesh& mesh, int vId1, int vId2, int vId3);
+	bool inside(Mesh& mesh, int vId, int cId);
+	bool across(Mesh& mesh, int eId1, int eId2);
+
+
+	// functions test correctness
+	inline bool euler(Mesh& mesh){
+		return (mesh.n_of_v() +
+				mesh.n_of_c() -
+				mesh.n_of_e() == 1);
+	}
+	bool length_gt0(Mesh& mesh);
+	bool areas_gt0(Mesh& mesh);
+	bool no_interior_intersection(Mesh& mesh);
+	bool conformal(Mesh& mesh);
+	bool bounded(Mesh& mesh);
+	bool counterclockwise(Mesh& mesh);
+
+	// functions in the service module
+	inline bool is_valid_mesh(Mesh& mesh){
+		return areas_gt0(mesh) && bounded(mesh) &&
+		conformal(mesh) && no_interior_intersection(mesh) &&
+		length_gt0(mesh) && counterclockwise(mesh) &&
+		euler(mesh);
+
+	}
+	bool covering_up(Mesh& mesh1, Mesh& mesh2);
+
+
+};
+}
+

Examples/MeshGenToolbox/mgt/Vertex.h

+/*****************************************************************************/
+/*  PMGT                                                                     */
+/*                                                                           */
+/*  January, 2007                                                            */
+/*                                                                           */
+/*  Copyright 2007                                                           */
+/*  Wen Yu                                                                   */
+/*  Department of Computing and Software                                     */
+/*  McMaster University                                                      */
+/*  yuw4@mcmaster.ca                                                         */
+/*                                                                           */
+/*  This program may be freely redistributed under the condition that the    */
+/*    copyright notices are not removed. You may distribute modified versions*/
+/*    of this code UNDER THE CONDITION THAT THIS CODE AND ANY MODIFICATIONS  */
+/*    MADE TO IT IN THE SAME FILE REMAIN UNDER COPYRIGHT OF THE ORIGINAL     */
+/*    AUTHOR, BOTH SOURCE AND OBJECT CODE ARE MADE FREELY AVAILABLE WITHOUT  */
+/*    CHARGE, AND CLEAR NOTICE IS GIVEN OF THE MODIFICATIONS.                */
+/*                                                                           */
+/*  This is the serial version of PMGT                                       */
+/*  This file includes classes for the vertex module.                        */
+/*****************************************************************************/
+
+#ifndef VERTEX
+#define VERTEX
+#endif
+
+#include <vector>
+#include <cfloat>
+#include <cmath>
+#include <iostream>
+#include <cstdlib>
+using namespace std;
+
+namespace PMGT{
+//namespace Vertex{
+// type definition
+
+// Point defines a coordinate of a point.
+// The first is x coordinate,
+// and the second is y coordinate.
+typedef vector<double>	Point;
+
+// MFLAG is a flag that indicates what should be done to a cell.
+enum MFLAG {COARSEN = -1, UNCHANGE = 0, REFINE = 1};
+
+class Vertex{
+public:
+
+	//constractor
+	Vertex(Point point, int vId=-1, int OhId=-1, bool del=false){
+	point_ = point;
+	vId_= vId;
+	OhId_= OhId;
+	del_= del;
+	}
+
+	Vertex(double x=0, double y=0, int vId=-1, int OhId=-1, int del=false){
+		point_.push_back(x);
+		point_.push_back(y);
+		vId_ = vId;
+		OhId_ = OhId;
+		del_ = del;
+	}
+
+	inline double x(){ return point_[0]; }
+	inline double y(){ return point_[1]; }
+	inline Point& point(){ return point_;}
+	inline int vId(){ return vId_;}
+	inline int OhId(){ return OhId_; }
+	inline bool del(){ return del_; }
+
+
+	inline void setPoint(double x, double y){
+		point_[0] = x;
+		point_[1] = y;
+	}
+	inline void setPoint(Point point){
+		point_ = point;
+	}
+	inline void setVId(int vId){ vId_ = vId; }
+	inline void setOhId(int OhId){ OhId_ = OhId; }
+	inline void setDel(bool del) { del_ = del; }
+
+	inline bool is_isolated_v(){ return OhId_ == -1; }
+	inline void set_isolated_v(){ OhId_ = -1; }
+
+	inline bool same(double x, double y){
+		return (x==point_[0] && y==point_[1]);
+	}
+	inline bool same(Point point){
+		return point_ == point;
+	}
+	inline void vcopy(Vertex& other){
+		point_[0] = other.point_[0];
+		point_[1] = other.point_[1];
+		vId_ = other.vId_;
+		OhId_ = other.OhId_;
+		del_ = other.del_;
+	}
+	
+	inline double distance(Point& p){
+		return sqrt((p[0]-point_[0])*(p[0]-point_[0])+(p[1]-point_[1])*(p[1]-point_[1]));
+	}
+
+private:
+	Point	point_;	// coordinate
+	int		vId_;	// id
+	int		OhId_;	// the first (ccw) outgoing halfedge.
+	bool 	del_;	// a flag to indicate if the vertex
+					// is marked to be deleted
+}; // end of class Vertex
+
+//} // end of namespace Vertex
+} // end of namespace PMGT

Examples/MeshGenToolbox/mgt/makefile

+OBJECTS = Service.o Coarsening.o Refining.o Output.o Input.o Mesh.o
+	
+SOURCES = Service.cpp Coarsening.cpp Refining.cpp Output.cpp Input.cpp Mesh.cpp
+
+HEADERS = Service.h Coarsening.h Refining.h Input.h Input.h Mesh.h Cell.h Edge.h Vertex.h
+
+# all testing excutable
+all: $(OBJECTS)
+	
+# all object files of PMGT
+$(OBJECTS): $(SOURCES) $(HEADERS)
+	g++ -c $(SOURCES)
+	
+# remove generated files
+clean:
+	rm -f *.o
\ No newline at end of file

Examples/MeshGenToolbox/testdriver/5_test_sr_arc.cpp

+/*****************************************************************************/
+/*  PMGT                                                                     */
+/*                                                                           */
+/*  January, 2007                                                            */
+/*                                                                           */
+/*  Copyright 2007                                                           */
+/*  Wen Yu                                                                   */
+/*  Department of Computing and Software                                     */
+/*  McMaster University                                                      */
+/*  yuw4@mcmaster.ca                                                         */
+/*                                                                           */
+/*  This program may be freely redistributed under the condition that the    */
+/*    copyright notices are not removed. You may distribute modified versions*/
+/*    of this code UNDER THE CONDITION THAT THIS CODE AND ANY MODIFICATIONS  */
+/*    MADE TO IT IN THE SAME FILE REMAIN UNDER COPYRIGHT OF THE ORIGINAL     */
+/*    AUTHOR, BOTH SOURCE AND OBJECT CODE ARE MADE FREELY AVAILABLE WITHOUT  */
+/*    CHARGE, AND CLEAR NOTICE IS GIVEN OF THE MODIFICATIONS.                */
+/*                                                                           */
+/*  This is the test case TC5.                                               */
+/*  This test case is design for test correctness of PMGT.                   */
+/*****************************************************************************/
+#include <iostream>
+
+#ifndef SERVICE
+#include "Service.h"
+#endif
+
+#ifndef COARSENING
+#include "Coarsening.h"
+#endif
+
+#ifndef REFINING
+#include "Refining.h"
+#endif
+
+#ifndef OUTPUT
+#include "Output.h"
+#endif
+
+#ifndef INPUT
+#include "Input.h"
+#endif
+
+#ifndef MESH
+#include "Mesh.h"
+#endif
+
+using namespace PMGT;
+int main(){
+
+	Mesh mesh;
+	Input input;
+	Refining ref;
+	Output output;
+	
+	input.read_file(mesh, "oldVertices.dat", "oldCells.dat");
+
+
+	double x = 0.4;
+	int count = 0;
+	char vfile[32];
+	char* vp = vfile;
+	char cfile[32];
+	char* cp = cfile;
+	char number[8];
+	while(1){
+
+		// refine
+		count++;
+		sprintf(number, "%i", count);
+
+		vp = strcpy(vp, "newVertices");
+		vp = strcat(vp, number);
+		vp = strcat(vp,".dat");
+
+		cp = strcpy(cp, "newCells");
+		cp = strcat(cp, number);
+		cp = strcat(cp,".dat");
+
+		int n = mesh.n_of_c();
+		vector<int> vIds;
+		Point zero; 
+		Vertex vtx[3];
+		
+		zero.clear();
+		zero.push_back(0);
+		zero.push_back(0);
+
+		for (int i=0; i<n; i++){
+			vIds.clear();
+			vIds = mesh.vIds_from_cId(i);
+			for (int j=0; j<3; j++){
+				vtx[j] = mesh.vertex(vIds[j]);
+			}
+
+			double first = vtx[0].distance(zero)-0.7;
+			double second = vtx[1].distance(zero)-0.7;
+			if (first*second<DBL_EPSILON){
+				mesh.cell(i).setMflag(REFINE);
+			}
+			else{
+				double third = vtx[2].distance(zero)-0.7;
+				if (first*third < DBL_EPSILON){
+					mesh.cell(i).setMflag(REFINE);
+				}
+			}
+		}
+
+		ref.longest_edge(mesh);
+		output.write_file(mesh,vp, cp);
+
+		if (count==15){
+			break;
+		}
+	}
+
+	cout << "There are " << count << " meshes generated." << endl;
+	
+	return 0;
+}

Examples/MeshGenToolbox/testdriver/ReadMe.txt

+This is the serial version of PMGT.
+
+There is only one test case (TC5).
+
+To compile this test file, just type "make"
+
+To remove generated object files and the executable file, type "make cleanAll"
+
+To remove generated output files, type "make cleanOutput"
\ No newline at end of file

Examples/MeshGenToolbox/testdriver/makefile

+OBJECTS = ../mgt/Service.o ../mgt/Coarsening.o ../mgt/Refining.o ../mgt/Output.o ../mgt/Input.o ../mgt/Mesh.o
+TEST_OBJ = 5_test_sr_arc.o
+	
+TEST_SRC = 5_test_sr_arc.cpp
+
+ALL_TEST = 5_test_sr_arc	
+
+# all testing executable
+all: $(TEST_OBJ)
+	@cd ../mgt; $(MAKE) $@
+	g++ -o 5_test_sr_arc $(OBJECTS) 5_test_sr_arc.o
+
+$(TEST_OBJ): $(TEST_SRC)
+	g++ -c $(TEST_SRC) -I../mgt
+
+# remove generated testing and PMGT files
+clean:
+	@cd ../mgt; $(MAKE) $@
+	rm -f *.o $(ALL_TEST)
+
+# remove generated testing files.
+# keep the PMGT files
+cleanTest:
+	rm -f $(TEST_OBJ) $(ALL_TEST)
+
+# remove generated output files
+cleanOutput:
+	rm newVertices* newCells*

Examples/MeshGenToolbox/testdriver/makefile2.txt

+OBJECTS = Service.o Coarsening.o Refining.o Output.o Input.o Mesh.o
+TEST_OBJ = 5_test_sr_arc.o
+	
+SOURCES = Service.cpp Coarsening.cpp Refining.cpp Output.cpp Input.cpp Mesh.cpp
+TEST_SRC = 5_test_sr_arc.cpp
+
+HEADERS = Service.h Coarsening.h Refining.h Input.h Input.h Mesh.h Cell.h Edge.h Vertex.h
+
+ALL_TEST = 5_test_sr_arc	
+
+# all testing executable
+all: $(TEST_OBJ)
+	@cd ../mgt; $(MAKE) $@
+	g++ -o 5_test_sr_arc $(OBJECTS) 5_test_sr_arc.o -L../mgt
+	
+5_test_sr_arc:5_test_sr_arc.o
+	g++ -o 5_test_sr_arc $(OBJECTS) 5_test_sr_arc.o
+	
+$(TEST_OBJ): $(TEST_SRC)
+	g++ -c $(TEST_SRC)
+
+# remove generated testing and PMGT files
+clean:
+	@cd src; $(MAKE) $@
+	rm *.o $(ALL_TEST)
+	
+# remove generated testing files.
+# keep the PMGT files
+cleanTest:
+	rm $(TEST_OBJ) $(ALL_TEST)
+	
+# remove generated output files
+cleanOutput:
+	rm newVertices* newCells*
\ No newline at end of file

Examples/MeshGenToolbox/testdriver/oldCells.dat

+0 1 3
+0 3 2

Examples/MeshGenToolbox/testdriver/oldVertices.dat

+0.0 0.0
+1.0 0.0
+0.0 1.0
+1.0 1.0

Examples/MeshGenToolbox/testdriver/show_serial.m

+function show_serial(n)
+
+% plot original mesh
+load oldVertices.dat; load oldCells.dat;
+oldCells(:,:) = oldCells(:,:)+1;
+triplot(oldCells, oldVertices(:,1), oldVertices(:,2))
+title('Original Mesh')
+axis equal
+
+% plot new mesh
+for (i=1:n)
+    vfile = strcat('newVertices',num2str(i),'.dat');
+    vfid = fopen(vfile, 'r');
+    v = fscanf(vfid,'%f %f');
+    fclose(vfid);
+
+    [nv,nnv] = size(v);
+    vv = zeros(nv/2, 2);
+    [mv,mmv] = size(vv);
+    for(j=1:mv)
+        vv(j,1) = v((j-1)*2+1);
+        vv(j,2) = v((j-1)*2+2);
+    end
+    
+    
+    cfile = strcat('newCells',num2str(i),'.dat');
+    cfid = fopen(cfile, 'r');
+    c = fscanf(cfid,'%g %g %g');
+    fclose(cfid);
+
+    [nc,nnc] = size(c);
+    cc = zeros(nc/3, 3);
+    [mc,mmc] = size(cc);
+    for(j=1:mc)
+        cc(j,1) = c((j-1)*3+1);
+        cc(j,2) = c((j-1)*3+2);
+        cc(j,3) = c((j-1)*3+3);
+    end
+    cc(:,:) = cc(:,:)+1;
+    % plot
+    figure
+    triplot(cc, vv(:,1), vv(:,2))
+    t = strcat('Mesh for "', vfile,'" and "', cfile,'"');
+    title(t)
+    axis equal
+end
+
+
+ 
+% cells2(:,:) = cells2(:,:)+1;
+% cells3(:,:) = cells3(:,:)+1;
+%     
+% % plot
+% triplot(cells1, vertices1(:,1), vertices1(:,2))
+% title('Test Creation')
+% axis equal
+% 
+% figure
+% triplot(cells2, vertices2(:,1), vertices2(:,2))
+% title('Test Refinement')
+% axis equal
+% 
+%  figure
+%  triplot(cells3, vertices3(:,1), vertices3(:,2))
+%  title('Test Coarsening')
+%  axis equal
\ No newline at end of file

Examples/MeshGenToolbox/testdriver/testingCode.cpp

+// testing code
+
+#include <iostream>
+
+
+#ifndef MESH
+#include "Mesh.h"
+#endif
+
+#ifndef INPUT
+#include "Input.h"
+#endif
+
+#ifndef OUTPUT
+#include "Output.h"
+#endif
+
+
+#ifndef REFINING
+#include "Refining.h"
+#endif
+
+#ifndef COARSENING
+#include "Coarsening.h"
+#endif
+
+#ifndef SERVICE
+#include "Service.h"
+#endif
+
+
+//using namespace PMGT::Vertex;
+//using namespace PMGT::Edge;
+//using namespace PMGT::Cell;
+//using namespace PMGT::Mesh;
+
+using namespace PMGT;
+
+int main(){
+
+	char* vertexFile = "oldVertices.dat";
+	char* cellFile = "oldCells.dat";
+
+	int i,nc;
+
+	Mesh mesh;
+	Input input;
+	Output output;
+	Refining ref;
+	Coarsening cos;
+	Service sev;
+
+	input.read_file(mesh, vertexFile, cellFile);
+	if (sev.is_valid_mesh(mesh)){
+		cout << "Mesh creatation is ok. The # of cells is " << mesh.n_of_c() << endl;
+	}
+	else{
+		cout << "Mesh creatation is wrong." << endl;
+	}
+
+	nc=mesh.n_of_c();
+	for (i=0; i<nc; i++){
+		mesh.cell(i).setMflag(REFINE);
+	}
+	//ref.split(mesh);
+	ref.longest_edge(mesh);
+	if (sev.is_valid_mesh(mesh)){
+		cout << "Mesh refinement is ok. The # of cells is " << mesh.n_of_c()  << endl;
+	}
+	else{
+		cout << "Mesh refinement is wrong." << endl;
+	}
+
+
+	nc=mesh.n_of_c();
+	for (i=0; i<nc; i++){
+			mesh.cell(i).setMflag(COARSEN);
+	}
+	cos.edge_collapse(mesh);
+	if (sev.is_valid_mesh(mesh)){
+		cout << "Mesh coarsening is ok. The # of cells is " << mesh.n_of_c()  << endl;
+	}
+	else{
+		cout << "Mesh coarsening is wrong." << endl;
+	}
+
+	output.write_file(mesh, "newVertices.dat", "newCells.dat");
+	return 0;
+}