Index: /anuga_core/source/anuga/alpha_shape/alpha_shape.py
===================================================================
--- /anuga_core/source/anuga/alpha_shape/alpha_shape.py	(revision 4138)
+++ /anuga_core/source/anuga/alpha_shape/alpha_shape.py	(revision 4139)
@@ -179,5 +179,5 @@
 
         ## Build Delaunay triangulation
-        import anuga.mesh_engine.triang as triang
+        from anuga.mesh_engine.mesh_engine import generate_mesh
         points = []
         seglist = []
@@ -186,5 +186,4 @@
         pointattlist = []
         segattlist = []
-        trilist = []
 
         points = [(pt[0], pt[1]) for pt in self.points]
@@ -192,6 +191,6 @@
         mode = "Qzcn"
         #print "computing delaunay triangulation ... \n" 
-        tridata = triang.genMesh(points,seglist,holelist,regionlist,
-                                 pointattlist,segattlist,trilist,mode)
+        tridata = generate_mesh(points,seglist,holelist,regionlist,
+                                 pointattlist,segattlist,mode)
         #print tridata
         #print "point attlist: ", tridata['generatedpointattributelist'],"\n"
Index: /anuga_core/source/anuga/mesh_engine/compile.py
===================================================================
--- /anuga_core/source/anuga/mesh_engine/compile.py	(revision 4138)
+++ /anuga_core/source/anuga/mesh_engine/compile.py	(revision 4139)
@@ -319,8 +319,7 @@
               print "********** Manually doing dependencies **************"
               compile(['triang.c','triangle.c'])
-            elif  filename == 'mesh_engine.c': # not part of ANUGA
-              continue
+            elif  filename == 'mesh_engine_c_layer.c': 
               #print "********** Manually doing dependencies **************"
-              #compile(['mesh_engine.c','triangle.c'])
+              compile(['mesh_engine_c_layer.c','triangle.c'])
               
             else:
Index: /anuga_core/source/anuga/mesh_engine/me_compile.py
===================================================================
--- /anuga_core/source/anuga/mesh_engine/me_compile.py	(revision 4138)
+++ /anuga_core/source/anuga/mesh_engine/me_compile.py	(revision 4139)
@@ -315,7 +315,7 @@
           print '--------------- Trying to compile c-extension %s' %filename
           try:
-            if filename == 'mesh_engine.c':
+            if filename == 'mesh_engine_c_layer.c':
               print "********** Manually doing dependencies **************"
-              compile(['mesh_engine.c','triangle.c'])
+              compile(['mesh_engine_c_layer.c','triangle.c'])
             else:
               compile(filename)
Index: uga_core/source/anuga/mesh_engine/me_test.py
===================================================================
--- /anuga_core/source/anuga/mesh_engine/me_test.py	(revision 4138)
+++ 	(revision )
@@ -1,378 +1,0 @@
-#!/usr/bin/env python
-#
-
-import sys
-
-import unittest
-import anuga.mesh_engine.mesh_engine as triang
-
-from Numeric import array
-
-class triangTestCase(unittest.TestCase):
-    def setUp(self):
-        pass
-
-    def tearDown(self):
-        pass
-
-    def metestrectangle(self):
-
-        points = []
-        seglist = []
-        segattlist = []
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
-	# This will cause a seg fault
-        #points = [[0.0,0.0],[0.0,10.0],[3.0,0.0],[3.0,10.0]]
-
-        # automatically generate pointattlist
-	pointattlist = []
-        for point in points:
-            pointattlist.append([])
-        #print "pointattlist",pointattlist
-	#pointattlist = [[],[],[],[]]
-        trilist = []
-        mode = "Qzcn"
-
-        
-        test = array([1.0, 2.0, 3.0])
-        print "test", test
-        
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist,trilist, mode,test)
-        self.failUnless(data['generatedtrianglelist'] ==[(1, 0, 2), (2, 3, 1)],
-                        'generatedtrianglelist is wrong!')
-        #print "data['generatedsegmentlist']",data['generatedsegmentlist']
-        self.failUnless(data['generatedsegmentlist'] ==[(2, 0), (3, 2), (1, 3),
-                                                        (0, 1)]
-              ,          'generatedsegmentlist is wrong!')
-        self.failUnless(data['generatedpointlist'] ==[(0.0, 0.0), (0.0, 10.0),
-                                                      (3.0, 0.0), (3.0, 10.0)],
-                        ' is wrong!')
-
-        #print "data['generatedtriangleneighborlist']",
-        #data['generatedtriangleneighborlist']
-        self.failUnless(data['generatedtriangleneighborlist'] ==[(-1, 1, -1),
-                                                                 (-1, 0, -1)],
-                        'generatedtriangleneighborlist is wrong!')
-
-
-    def testrectangleII(self):
-
-        points = []
-        seglist = []
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
-	pointattlist = [[],[],[],[]]
-	regionlist.append( (1.2,1.2,5.0) )
-	seglist = [(0,1),(1,3),(3,2),(2,0)]
-        segattlist = [0,0,0,0]
-        trilist = []
-        mode = "Qzp"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist,trilist, mode)
-
-        self.failUnless(data['generatedtrianglelist'] ==[(1, 0, 2), (2, 3, 1)],
-                        'trianglelist is wrong!')
-        self.failUnless(data['generatedsegmentlist'] ==[(0, 1), (1, 3),
-                                                        (3, 2), (2, 0)],
-                        'segmentlist is wrong!')
-        self.failUnless(data['generatedpointlist'] ==[(0.0, 0.0), (0.0, 10.0),
-                                                      (3.0, 0.0), (3.0, 10.0)],
-                        ' is wrong!')
-
-    def testsegmarker(self):
-
-        points = []
-        seglist = []
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
-	pointattlist = [[],[],[],[]]
-	regionlist.append( (1.2,1.2,5.0) )
-	seglist = [(0,1),(1,3),(3,2),(2,0)]
-        segattlist = [1.0,2.0,3.0,4.0]
-        trilist = []
-        mode = "Qzp"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist,trilist, mode)
-
-        self.failUnless(data['generatedtrianglelist'] ==[(1, 0, 2), (2, 3, 1)],
-                        'trianglelist is wrong!')
-        self.failUnless(data['generatedsegmentlist'] ==[(0, 1), (1, 3),
-                                                        (3, 2), (2, 0)],
-                        'segmentlist is wrong!')
-        self.failUnless(data['generatedpointlist'] ==[(0.0, 0.0), (0.0, 10.0),
-                                                      (3.0, 0.0), (3.0, 10.0)],
-                        ' is wrong!')
-        self.failUnless(data['generatedsegmentmarkerlist'] ==[1,2,3,4],
-                        ' is wrong!')
-
-    def testbad_region(self):
-
-        points = []
-        seglist = []
-        holelist = []
-        regionlist = [(1.2,1.2)]
-
-	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
-	pointattlist = [[],[],[],[]]
-	seglist = [(0,1),(1,3),(3,2),(2,0)]
-        segattlist = [0,0,0,0]
-        trilist = []
-        mode = "Qzpn"
-        try:
-            data = triang.genMesh(points,seglist,holelist,
-                                  regionlist,pointattlist,segattlist,
-                                  trilist, mode)
-
-        except TypeError:
-            pass
-        else:
-            self.failUnless(0 ==1,
-                        'bad region list did not raise error!')
-
-    def testregion_with_maxarea(self):
-
-        points = []
-        seglist = []
-        holelist = []
-        regionlist = [(3,1,1.0)]
-
-	points = [(0.0,0.0),(6.0,0.0),(6.0,6.0),(0.0,6.0)]
-	pointattlist = [[],[],[],[]]
-	seglist = [(0,1),(1,2),(3,2),(3,0),(0,2)]
-        segattlist = [0,0,0,0,0]
-        trilist = []
-        mode = "Qzpna36a"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist, trilist, mode)
-
-        self.failUnless(len(data['generatedtrianglelist']) == 2,
-                        'testregion_with_maxarea 1: # of tris is wrong!')
-        ## Another test case
-        regionlist = [(3,1,1.0),(1,3,1.0,8.0)]
-        mode = "Qzp21na36a"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist, trilist, mode)
-        #print "len(data['generatedtrianglelist']",len(data['generatedtrianglelist'])
-        # for Duncan On unix this returns a 7 triangle result.
-        # for Duncan on Windows returns a 6 triangle result.
-        # for Ole on nautilus this returns 6
-        # for Duncan on nautilus this returns 7
-        # ??, it seems to be the results from triangle that is
-        # causing the different results, and we are treating
-        # triangle as a back box.
-
-        self.failUnless(len(data['generatedtrianglelist']) >= 6,
-                        'testregion_with_maxarea 2: # of tris is wrong!')
-        ## Another test case
-        regionlist = [(3,1,1.0,8.0),(1,3,1.0,8.0)]
-        mode = "Qzpna36a"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist, trilist, mode)
-        #print "len(data['generatedtrianglelist']",len(data['generatedtrianglelist'])
-        # On unix this returns a 10 triangle result.
-        # Windows returns a 8 triangle result.
-        self.failUnless(len(data['generatedtrianglelist']) >= 8,
-                        'testregion_with_maxarea 3: # of tris is wrong!')
-
-        ## Another test case
-        regionlist = [(3,1,1.0),(1,3,1.0,8.0)]
-        mode = "Qzpna8a"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist, trilist, mode)
-        #print "len(data['generatedtrianglelist']",len(data['generatedtrianglelist'])
-        # On unix this returns a 10 triangle result.
-        # Windows returns a 8 triangle result.
-
-        self.failUnless(len(data['generatedtrianglelist']) >= 8,
-                        'testregion_with_maxarea 4: # of tris is wrong!')
-
-    def testbad_point(self):
-
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0)]
-	pointattlist = [[],[],[],[]]
-	seglist = [(0,1),(1,3),(3,2),(2,0)]
-        segattlist = [0,0,0,0]
-        trilist = []
-        mode = "Qzpn"
-        try:
-            data = triang.genMesh(points,seglist,holelist,regionlist,
-                                  pointattlist,segattlist, trilist, mode)
-
-        except TypeError:
-            pass
-        else:
-            self.failUnless(0 ==1,
-                        'bad point list did not raise error!')
-
-    def testbad_hole(self):
-
-        holelist = [(9.0)]
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
-	pointattlist = [[],[],[],[]]
-	seglist = [(0,1),(1,3),(3,2),(2,0)]
-        segattlist = [0,0,0,0]
-        trilist = []
-        mode = "Qzpn"
-        try:
-            data = triang.genMesh(points,seglist,holelist,regionlist,
-                                  pointattlist,segattlist, trilist, mode)
-
-        except TypeError:
-            pass
-        else:
-            self.failUnless(0 ==1,
-                        'bad hole list did not raise error!')
-
-    def testbad_segment(self):
-
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
-	pointattlist = [[],[],[],[]]
-	seglist = [(0,1),(1,3),(3,2),(2)]
-        segattlist = [0,0,0,0]
-        trilist = []
-        mode = "Qzpn"
-        try:
-            data = triang.genMesh(points,seglist,holelist,regionlist,
-                                  pointattlist,segattlist, trilist, mode)
-
-        except TypeError:
-            pass
-        else:
-            self.failUnless(0 ==1,
-                        'bad segment list did not raise error!')
-
-    def testbad_segattlist(self):
-
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
-	pointattlist = [[],[],[],[]]
-	seglist = [(0,1),(1,3),(3,2),(2,0)]
-        segattlist = [0,0]
-        trilist = []
-        mode = "Qzpn"
-        try:
-            data = triang.genMesh(points,seglist,holelist,regionlist,
-                                  pointattlist,segattlist, trilist, mode)
-
-            self.failUnless(data['trianglelist'] ==[(1, 0, 2), (2, 3, 1)],
-                        'trianglelist is wrong!')
-        except TypeError:
-            pass
-        else:
-            self.failUnless(0 ==1,
-                        'bad segment attribute list did not raise error!')
-
-    def testrectangle_regions(self):
-
-        points = []
-        seglist = []
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
-	pointattlist = [[],[],[],[]]
-        #it seems that
-        #triangle only associates one region with a triangle
-	regionlist.append( (1.3,1.3,88.33) )
-	regionlist.append( (1.2,1.2,77,55) )
-	seglist = [(0,1),(1,3),(3,2),(2,0)]
-        segattlist = [0,0,0,0]
-        trilist = []
-        mode = "QAzpq"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist, trilist, mode)
-
-        self.failUnless(data['generatedtriangleattributelist'] ==[[77.0], [77.0], [77.0], [77.0]],
-                        'triangleattributelist is wrong!')
-
-
-    def BAD_SEG_FAULT_test_lone_verts(self):
-
-        points = []
-        seglist = []
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0),(0.0,10.0)]
-	pointattlist = [[],[],[],[],[]]
-	regionlist.append( (1.2,1.2,5.0) )
-	seglist = [(0,4),(4,1),(1,3),(3,2),(2,0)]
-        segattlist = [0,0,0,0,0]
-        trilist = []
-        mode = "QpznAa2000.1a"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist,trilist, mode)
-        #print "data['generatedtrianglelist']", data['generatedtrianglelist']
-        self.failUnless(data['generatedtrianglelist'] ==[(4, 0, 2), (2, 3, 4)],
-                        'trianglelist is wrong!')
-        self.failUnless(data['generatedsegmentlist'] ==[(0, 4), (4, 3),
-                                                        (3, 2), (2, 0)],
-                        'segmentlist is wrong!')
-        self.failUnless(data['generatedpointlist'] ==[(0.0, 0.0), (0.0, 10.0),
-                                                      (3.0, 0.0), (3.0, 10.0),
-                                                      (0.0,10.0)],
-                        ' is wrong!')
-        self.failUnless(data['lonepointlist'] ==[1],
-                        'lonepointlist is wrong!')
-        #print "", data['lonepointlist']
-
-    def BAD_SEG_FAULT_test_lone_vertsII(self):
-
-        points = []
-        seglist = []
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,0.0),(0.0,10.0),(0.0,10.0),(10.0,10.0),
-                  (10.0,10.0),(0.0,10.0),(10.0,0.0)]
-
-	pointattlist = []
-        for point in points:
-            pointattlist.append([])
-	seglist = [(0,1),(1,2),(2,3),(3,4),(4,5),(5,7),(7,0)]
-	segattlist = []
-        for seg in seglist:
-            segattlist.append(0)
-        trilist = []
-        mode = "pznAa2000.1a"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist,trilist, mode)
-        #print "data['generatedtrianglelist']", data['generatedtrianglelist']
-        self.failUnless(data['generatedtrianglelist'] ==[(6, 1, 7), (7, 5, 6)],
-                        'trianglelist is wrong!')
-        self.failUnless(data['generatedsegmentlist'] ==[(1, 6), (6, 5),
-                                                        (5, 7), (7, 1)],
-                        'segmentlist is wrong!')
-        self.failUnless(data['generatedpointlist'] ==[(0.0,0.0),(0.0,0.0),
-                                                      (0.0,10.0),(0.0,10.0),
-                                                      (10.0,10.0),
-                                                      (10.0,10.0),(0.0,10.0),
-                                                      (10.0,0.0)],
-                        ' is wrong!')
-        self.failUnless(data['lonepointlist'] ==[0,2,3,4],
-                        'lonepointlist is wrong!')
-
-if __name__ == "__main__":
-
-    #suite = unittest.makeSuite(triangTestCase,'test')
-    suite = unittest.makeSuite(triangTestCase,'metest')
-    runner = unittest.TextTestRunner()
-    runner.run(suite)
Index: uga_core/source/anuga/mesh_engine/mesh_engine.c
===================================================================
--- /anuga_core/source/anuga/mesh_engine/mesh_engine.c	(revision 4138)
+++ 	(revision )
@@ -1,631 +1,0 @@
-/* this appears to be necessary to run on Linux */
-
-#undef _STDLIB_H_
-#define _STDLIB_H_
-
-
-/*
-
-    This code is testing passing numeric structures into the
-    triangle interface, rather than python lists.
-    It is not operational code.
-    
-    This code is based on code by A. Pletzer
-
-    This code interfaces pmesh directly with "triangle", a general        
-    purpose triangulation code. In doing so, Python data structures          
-    are passed to C for  use by "triangle" and the results are then          
-    converted back. This was accomplished using the Python/C API.            
-                                                                           
-    I. Lists of points,edges,holes, and regions must normally be 
-     created for proper triangulation. However, if there are no holes,
-     the  hole list need not be specified (however, a 
-     measure of control over the shape of the convex hull is lost; 
-     triangle decides). 
-     
-     points list: [(x1,y1),(x2,y2),...] (Tuples of doubles)                  
-     edge list: [(Ea,Eb),(Ef,Eg),...] (Tuples of integers)                   
-     hole list: [(x1,y1),...](Tuples of doubles, one inside each hole region)
-     regionlist: [ (x1,y1,index,area),...] (Tuple of doubles)                
-     pointattributelist:[ (a1,a2..),...] (Tuple of lists) 
-     mode: String - the commands sent to triangle
-          
-     This is all that is needed to accomplish an intial triangulation.       
-     A dictionary of the Triangle output is returned.                                   
-     Precondition
-     End list in the pointattributelist has to have the same length
-     
-    Compilation can be done as follows                                       
-                    
-    c++ -I/usr/local/include/python1.5 -shared -o triangmodule.so            
-    triangmodule.c triangle.o                                                */
-
-
-#ifdef SINGLE  
-#define REAL float 
-#else  
-#define REAL double 
-#endif 
-
-#include "triangle.h" 
-
-#ifndef _STDLIB_H_ 
-#ifndef _WIN32
-extern "C" void *malloc(); 
-extern "C" void free(); 
-#endif
-#endif  
-
-#include "Python.h" 
-
-#include "Numeric/arrayobject.h"
-
-
- static PyObject *triang_genMesh(PyObject *self, PyObject *args){
-     
-  /* Mesh generation routine
-     pre-conditions:
-     mode must have a 'z' switch, to avoid segmentation faults
-     
-     shallow_water_ext.c gravity function
-  */
-    
-  struct triangulateio in,out;
-
-  PyObject *hull;
-  PyObject *seglist;
-  PyObject *holelist;
-  PyObject *regionlist;
-  PyObject *trianglelist;
-  PyObject *mode;
-  PyObject *firstPointAttr;
-  PyObject *pointattributelist;
-  PyObject *segmarkerlist;
-  PyObject *Attrlist;
-  PyArrayObject *test;
-  REAL Attr;
-  int i, j, iatt, n, write_here;
-  int a,b,c;
-  int marker;
-  int tuplesize;
-    
-  int index = 0;
-  double x,y;
-  PyObject *elems;
-  PyObject *pair;
-     
-  char *mod;
-  PyObject *holder, *holderlist, *attributelist;
-  int listsize,attsize ;
-  PyObject  *ii;
-  
-  int *points_connected;
-  
-  /* used for testing numeric arrays*/
-  int n_test;     
-  if(!PyArg_ParseTuple(args,(char *)"OOOOOOOOO",&hull,&seglist,&holelist,&regionlist,&pointattributelist,&segmarkerlist,&trianglelist,&mode,&test)){
-    return NULL;
-  }
-  if(!PyList_Check(hull)){
-    PyErr_SetString(PyExc_TypeError, 
-		    "incorrect first argument for 'genMesh': list required.");
-    return NULL;
-  }
-  if(!PyList_Check(hull)){
-    PyErr_SetString(PyExc_TypeError,
-		    "incorrect second argument for 'genMesh': list required.");
-    return NULL;
-  }
-  if(!PyList_Check(hull)){
-    PyErr_SetString(PyExc_TypeError,
-		    "incorrect third argument for 'genMesh': list required.");
-    return NULL;
-  }
-  if(!PyList_Check(hull)){
-    PyErr_SetString(PyExc_TypeError,
-		    "incorrect fourth argument for 'genMesh': list required.");
-    return NULL;
-  }
-  if(!PyList_Check(hull)){
-    PyErr_SetString(PyExc_TypeError,
-		    "incorrect fifth argument for 'genMesh': list required.");
-    return NULL;
-  }
-  if(!PyList_Check(hull)){
-    PyErr_SetString(PyExc_TypeError,
-		    "incorrect sixth argument for 'genMesh': list required.");
-    return NULL;
-  }
-  if(!PyList_Check(hull)){
-    PyErr_SetString(PyExc_TypeError,
-		    "incorrect seventh argument for 'genMesh': list required.");
-    return NULL;
-  }
-  if(!PyList_Check(hull)){
-    PyErr_SetString(PyExc_TypeError,
-		    "incorrect eighth argument for 'genMesh': list required.");
-    return NULL;
-  }
-     
-  /* printf("starting\n"); */    
-  n_test=test -> dimensions[0];
-  printf("test -> dimensions[0] %i\n", n); 
-     
-  /* Initialize  points */
-  in.numberofpoints =  PyList_Size(hull); 
-  in.pointlist = (REAL *)malloc(in.numberofpoints*2*sizeof(REAL));
-  in.pointmarkerlist = (int *)NULL; /* malloc(in.numberofpoints*sizeof(int)); */
-     
-  /* Initialize and fill up region list */
-  in.numberofregions = PyList_Size(regionlist);
-  
-  /*printf(" numberofregions %i\n",in.numberofregions );*/
-  in.regionlist = (REAL *)malloc(in.numberofregions*4*sizeof(REAL));
-  for(i=0; i<in.numberofregions;i++){
-    elems = PyList_GetItem(regionlist,i);
-    tuplesize = PySequence_Length(elems);
-    /*printf("tuplesize %i\n",tuplesize);*/
-    if (tuplesize>0){
-      in.regionlist[4*i] = PyFloat_AsDouble(PyTuple_GetItem(elems,0));
-    } else {
-      PyErr_SetString(PyExc_TypeError,
-		      "data sent to trianglemodule error: A region has no x value.");
-      return NULL;
-    }
-    if (tuplesize>1){
-      in.regionlist[4*i+1] = PyFloat_AsDouble(PyTuple_GetItem(elems,1));
-    } else {
-      PyErr_SetString(PyExc_TypeError,
-		    "data sent to trianglemodule error: A region has no y value.");
-      return NULL;
-    }
-    if (tuplesize>2){
-      /* 2 is attribute*/
-      in.regionlist[4*i+2] = PyFloat_AsDouble(PyTuple_GetItem(elems,2));
-    } else {
-      PyErr_SetString(PyExc_TypeError,
-		    "data sent to trianglemodule error: A region has no attribute value.");
-      return NULL;
-    }
-    if (tuplesize>3){
-      in.regionlist[4*i+3] = PyFloat_AsDouble(PyTuple_GetItem(elems,3));
-    } else {
-      /* 0.0 is interpreted as no local max area by triangle  */
-      in.regionlist[4*i+3] = (REAL)0.0;
-    } 
-  }
-  
-  /*Initialize segments */
-  in.numberofsegments = PyList_Size(seglist);
-  in.segmentlist = (int *)malloc(in.numberofsegments*2*sizeof(int));
-  in.segmentmarkerlist = (int *)malloc(in.numberofsegments*sizeof(int));
-
-    
-  /*Initialize triangles */
-  in.numberofcorners = 3;
-  in.numberoftriangles = PyList_Size(trianglelist);
-  if (in.numberoftriangles  != 0) {
-    in.trianglelist = (int *)malloc(in.numberofsegments*3*sizeof(int));
-  } else {
-    in.trianglelist = (int *)NULL;
-  }
-  in.numberoftriangleattributes = 0;  
-  in.triangleattributelist  = (REAL *)NULL; 
-     
-     
-     
-  /*Initialize holes */
-  in.numberofholes = PyList_Size(holelist);
-  if (in.numberofholes != 0) {
-    in.holelist = (REAL *)malloc(in.numberofholes*2*sizeof(REAL));
-  } else {
-    in.holelist = (REAL *)NULL;
-  }     
-     
-  /* Fill up triangle list */
-  index = 0;
-  for(i = 0; i < in.numberoftriangles; i++) {
-    pair =  PyList_GetItem(trianglelist,i);
-    if (3 != PySequence_Length(pair)) {
-      PyErr_SetString(PyExc_TypeError,
-		      "data sent to trianglemodule error:Found a triangle without 3 vertices.");
-      return NULL;
-    }
-    a = (int)PyInt_AsLong(PyTuple_GetItem(pair,0));
-    in.trianglelist[index] = a;
-    index++;
-    b = (int)PyInt_AsLong(PyTuple_GetItem(pair,1));
-    in.trianglelist[index] = b;
-    index++;
-    c = (int)PyInt_AsLong(PyTuple_GetItem(pair,2));
-    in.trianglelist[index] = c;
-    index++;
-  }
-     
-     
-  /*  print the triangle list DEBUGGING
-      printf ("numberoftriangles %i\n", in.numberoftriangles);
-      for(i = 0; i < in.numberoftriangles; i++) {
-      printf("(%i,%i,%i)",in.trianglelist[i* 3+ 0],in.trianglelist[i* 3+ 1], in.trianglelist[i* 3+ 2]);
-      }*/
-     
-  /* Fill up segment list */
-  index = 0;
-  if (in.numberofsegments != PyList_Size(segmarkerlist)) {
-    PyErr_SetString(PyExc_TypeError,
-		    "data sent to trianglemodule error:number of segments doesn't = number of segment attributes.");
-    return NULL;
-  }
-    
-  for(i = 0; i < in.numberofsegments; i++) {
-    pair =  PyList_GetItem(seglist,i);
-    if (2 != PySequence_Length(pair)) {
-      PyErr_SetString(PyExc_TypeError,
-		      "data sent to trianglemodule error:Found a segment without 2 vertices.");
-      return NULL;
-    }
-    a = (int)PyInt_AsLong(PyTuple_GetItem(pair,0));
-    in.segmentlist[index] = a;
-    index++;
-    b = (int)PyInt_AsLong(PyTuple_GetItem(pair,1));
-    in.segmentlist[index] = b;
-    index++;
-    /* Fill up segment marker list */
-    marker =  (int)PyInt_AsLong(PyList_GetItem(segmarkerlist,i));
-    in.segmentmarkerlist[i] = marker;
-  }
-  /* Fill up hole list */
-  index = 0;
-  for(i = 0; i < in.numberofholes; i++) {
-    pair =  PyList_GetItem(holelist,i);
-    if (2 != PySequence_Length(pair)) {
-      PyErr_SetString(PyExc_TypeError,
-		      "data sent to trianglemodule error:Found a pair without 2 vertices.");
-      return NULL;
-    }
-    x = PyFloat_AsDouble(PyTuple_GetItem(pair,0));
-    in.holelist[index] = x;
-    index++;
-    y = PyFloat_AsDouble(PyTuple_GetItem(pair,1));
-    in.holelist[index] = y;
-    index++;
-  }
-     
-  /* Fill up point list */
-  index = 0;
-  for(i = 0; i < in.numberofpoints; i++) {
-    pair =  PyList_GetItem(hull,i);
-    if (2 != PySequence_Length(pair)) {
-      PyErr_SetString(PyExc_TypeError,
-		      "data sent to trianglemodule error:Found a vertex without an x,y values.");
-      return NULL;
-    }
-    x = PyFloat_AsDouble(PyTuple_GetItem(pair,0));
-    in.pointlist[index] = x;
-    index++;
-    y = PyFloat_AsDouble(PyTuple_GetItem(pair,1));
-    in.pointlist[index] = y;
-    index++;
-  }
-
-  /* Initialize  point attributes */
-/*      firstPointAttr = PyList_GetItem(pointattributelist,0); */
-/*      if (firstPointAttr != NULL) { */
-/*        //in.numberofpointattributes = PyList_Size(firstPointAttr); */
-/*        printf("in.numberofpointattributes =%i",PyList_Size(firstPointAttr)); */
-/*        //printf("in.numberofpointattributes =%i",in.numberofpointattributes); */
-/*      } else { */
-/*        printf("(firstPointAttr is NULL\n"); */
-/*        in.pointattributelist = NULL; */
-/*        in.numberofpointattributes = 0; */
-/*      } */
-  n = PyList_Size(pointattributelist);
-  if (n < 0) {
-    in.numberofpointattributes = 0; 
-    in.pointattributelist = (REAL *)NULL;
-  } else {
-    /* Get the # of attributes based on the size of the first point*/
-    firstPointAttr = PyList_GetItem(pointattributelist,0);
-    if (firstPointAttr == NULL) {
-      in.numberofpointattributes = 0; 
-      in.pointattributelist = (REAL *)NULL;
-    } else {
-      in.numberofpointattributes = PyList_Size(firstPointAttr);
-      in.pointattributelist = (REAL *)malloc(in.numberofpoints
-					     *in.numberofpointattributes
-					     *sizeof(REAL)); 
-     
-     
-      /* fill the point attribute list */
-      if (in.numberofpointattributes != 0) {
-	for(j = 0; j < in.numberofpoints; j++){
-	  Attrlist = PyList_GetItem(pointattributelist,j);
-	  if (in.numberofpointattributes != PySequence_Length(Attrlist)) {
-	    PyErr_SetString(PyExc_TypeError,
-			    "data sent to trianglemodule error:List size of attributes is different from vertex to vertex.");
-	    return NULL;
-	  }
-	  for(i = 0; i < in.numberofpointattributes; i++){
-	    Attr = PyFloat_AsDouble(PyList_GetItem(Attrlist,i));
-	    in.pointattributelist[in.numberofpointattributes*j+i]= (REAL)Attr;
-	    /* 	 printf("i =%i\n",i);  */
-	    /* 	 printf("j =%i\n",j);  */
-	    /* 	 printf("in.pointattributelist[4*j+i]
-		 =%f\n",in.pointattributelist[4*j+i]); */
-	  }
-	}
-      } 
-    }
-  }
-    
-     
-  /* this half from triangulate*/
-  /* set up the switch arguments to the triangulation routine */
-
-  mod = PyString_AS_STRING(mode);
-     
-         
-  out.pointlist = (REAL *)NULL;
-  out.pointmarkerlist = (int *)NULL;
-  out.numberofpointattributes = in.numberofpointattributes;
-  out.pointattributelist = (REAL *)NULL;
-  
-  out.trianglelist = (int *)NULL;
-  out.triangleattributelist = (REAL *)NULL;                   
-  out.trianglearealist = (REAL *)NULL;                        
-  out.neighborlist = (int *)NULL;
-     
-  out.segmentlist = (int *)NULL;
-  out.segmentmarkerlist = (int *)NULL;
-     
-  out.edgelist = (int *)NULL;
-  out.edgemarkerlist = (int *)NULL;
-     
-  out.holelist = (REAL *)NULL;
-  out.regionlist = (REAL *)NULL;
-    
-  
-  /*printf("\n\nTriangulate input args: %s \n\n", mod); */
-  triangulate(mod, &in, &out, (struct triangulateio *)NULL );
-
-  /*
-    ------- Pass point numbers,coordinates and neighbors back to Python ------
-    we send back a dictionary:                                               
-    { index : [ coordinates, [connections], Attribute ] } 
-  */
-  holder = PyDict_New();
-     
-  /* list of int's, used to keep track of which verts are connected to
-     triangles. */
-  points_connected = (int *)malloc(out.numberofpoints*sizeof(int));
-  
-  /* Add pointlist */
-  listsize = out.numberofpoints;
-  holderlist = PyList_New(listsize);
-     
-  for(i=0; i<listsize;i++){
-    PyObject *mlist = Py_BuildValue((char *)"(d,d)", 
-				    out.pointlist[i*2  ], out.pointlist[i*2+1]);  
-    PyList_SetItem(holderlist,i, mlist);
-    points_connected[i] = 0;
-  }  
-  ii=PyString_FromString("generatedpointlist");
-  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);
-    
-  /* Add point marker list */
-  listsize = out.numberofpoints;
-  holderlist = PyList_New(listsize);
-  
-  for(i=0; i<listsize;i++){
-    PyObject *mlist = Py_BuildValue((char *)"d", 
-				    out.pointmarkerlist[i]);     
-    PyList_SetItem(holderlist,i, mlist); 
-  }  
-  ii=PyString_FromString("generatedpointmarkerlist");
-  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);
-    
-  /* Add point attribute list */
-  listsize = out.numberofpoints;
-  holderlist = PyList_New(listsize);
-  
-  for(i=0; i<listsize;i++){
-    attsize = out.numberofpointattributes;
-    attributelist = PyList_New(attsize);
-    for(iatt=0; iatt<attsize;iatt++){
-      PyObject *mlist = Py_BuildValue((char *)"d", 
-				      out.pointattributelist[i*attsize + iatt]);     
-      PyList_SetItem(attributelist,iatt, mlist);
-    }      
-    PyList_SetItem(holderlist,i, attributelist);
-  }  
-  ii=PyString_FromString("generatedpointattributelist");
-  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);  
- 
-  
-  /* Add triangle list */
-  listsize = out.numberoftriangles;
-  holderlist = PyList_New(listsize);
-  for(i=0; i<listsize;i++){
-    PyObject *mlist = Py_BuildValue((char *)"(i,i,i)", 
-				    out.trianglelist[i*3  ], out.trianglelist [i*3+1], out.trianglelist [i*3+2]);    
-    PyList_SetItem(holderlist,i, mlist);
-    /* printf(" A vert index %i\n",out.trianglelist[i*3] );
-    printf(" A vert index %i\n",out.trianglelist[i*3+1] );
-    printf(" A vert index %i\n",out.trianglelist[i*3+2] ); */
-    points_connected[out.trianglelist[i*3]] = 1;
-    points_connected[out.trianglelist[i*3+1]] = 1;
-    points_connected[out.trianglelist[i*3+2]] = 1;
-  }    
-  ii=PyString_FromString("generatedtrianglelist");
-  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);
-  
-  /* convert the points_connected vector from a true(1) false(0) vector, where
-     index is the vert, to a vector of the lone verts, at the beggining
-     of the vector. */
-  write_here = 0;   
-  for(i=0; i<out.numberofpoints;i++){
-    if (0 == points_connected[i]) {
-      points_connected[write_here] = i;
-      write_here ++;
-    }
-  }   
-  /* printf(" ******************** \n" );
-  for(i=0; i<write_here;i++){
-    printf(" A vert index %i\n",points_connected[i] );
-    } */
-  
-  listsize = write_here;
-  holderlist = PyList_New(listsize);
-  for(i=0; i<listsize;i++){
-   PyObject *mlist = Py_BuildValue((char *)"i", points_connected[i]);    
-    PyList_SetItem(holderlist,i, mlist);
-  }
-  ii=PyString_FromString("lonepointlist");
-  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);  
-     
-  /* Add triangle attribute list */
-  listsize = out.numberoftriangles;
-  holderlist = PyList_New(listsize);
-     
-  for(i=0; i<listsize; i++){
-    attsize = out.numberoftriangleattributes;
-    attributelist = PyList_New(attsize);       
-    for(iatt=0; iatt<attsize;iatt++){
-      PyObject *mlist = Py_BuildValue((char *)"d",out.triangleattributelist[i*attsize + iatt]);     
-      PyList_SetItem(attributelist,iatt, mlist);
-    }      
-    PyList_SetItem(holderlist,i, attributelist);
-  }  
-  ii=PyString_FromString("generatedtriangleattributelist");
-  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);
-  
-  /* Add segment list */
-  listsize = out.numberofsegments;
-  holderlist = PyList_New(listsize);
-  for(i=0; i<listsize;i++){
-    PyObject *mlist = Py_BuildValue((char *)"(i,i)", 
-				    out.segmentlist[i*2  ],
-				    out.segmentlist [i*2+1]);
-    PyList_SetItem(holderlist,i, mlist);
-  }    
-  ii=PyString_FromString("generatedsegmentlist");
-  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);  
-  
-  /* Add segment marker list */
-  listsize = out.numberofsegments;
-  holderlist = PyList_New(listsize);
-  for(i=0; i<listsize;i++){
-    PyObject *mlist = Py_BuildValue((char *)"i", 
-				    out.segmentmarkerlist[i]);
-    PyList_SetItem(holderlist,i, mlist);
-  }    
-  ii=PyString_FromString("generatedsegmentmarkerlist");
-  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);  
-  /* Add triangle neighbor list */
-  if (out.neighborlist != NULL) {
-    listsize = out.numberoftriangles;
-    holderlist = PyList_New(listsize);
-    for(i=0; i<listsize;i++){
-      PyObject *mlist = Py_BuildValue((char *)"(i,i,i)", 
-				      out.neighborlist[i*3  ],
-				      out.neighborlist [i*3+1],
-				      out.neighborlist [i*3+2]);    
-      PyList_SetItem(holderlist,i, mlist);
-    }    
-    ii=PyString_FromString("generatedtriangleneighborlist");
-    PyDict_SetItem(holder, ii, holderlist);
-    Py_DECREF(ii); Py_DECREF(holderlist);
-  }   
-     
-  /* Free in/out structure memory */
-  
-  /* OUT */
-
-  if(!out.pointlist){
-    free(out.pointlist);  out.pointlist=NULL;
-  }
-  if(!out.pointmarkerlist){    
-    free(out.pointmarkerlist); out.pointmarkerlist=NULL;
-  }
-  if(!out.pointattributelist){    
-    free(out.pointattributelist); out.pointattributelist=NULL;
-  }
-    
- 
-  if(!out.trianglelist){    
-    free(out.trianglelist); out.trianglelist=NULL;
-  }
-  if(!out.triangleattributelist){    
-    free(out.triangleattributelist); out.triangleattributelist=NULL;
-  }
-  if(!out.trianglearealist){    
-    free(out.trianglearealist); out.trianglearealist=NULL;
-  }
-  if(!out.neighborlist){    
-    free(out.neighborlist); out.neighborlist=NULL;
-  }
-  if(!out.segmentlist){
-    free(out.segmentlist); out.segmentlist =NULL;
-  }
-  if(!out.segmentmarkerlist){
-    free(out.segmentmarkerlist); out.segmentmarkerlist  =NULL;
-  }
-  if(!out.edgelist){
-    free(out.edgelist);  out.edgelist=NULL;
-  }
-  if(!out.edgemarkerlist){
-    free(out.edgemarkerlist);  out.edgemarkerlist=NULL;
-  }
-  if(!out.holelist){
-    free(out.holelist); out.holelist=NULL;
-  }
-  if(!out.regionlist){
-    free(out.regionlist); out.regionlist=NULL;
-  }
-  if(!points_connected ){
-    free(points_connected ); points_connected =NULL;
-  }
-    
-  
-  
-  /*  IN */
-  
-  if(!in.pointlist){
-    free(in.pointlist);  in.pointlist   =NULL;             
-  }
-  if(!in.pointattributelist){
-    free(in.pointattributelist); in.pointattributelist =NULL;       
-  }
-  if(!in.pointmarkerlist){
-    free(in.pointmarkerlist); in.pointmarkerlist    =NULL;        
-  }
-  if(!in.segmentlist){
-    free(in.segmentlist);  in.segmentlist    =NULL;         
-  }
-  if(!in.segmentmarkerlist){
-    free(in.segmentmarkerlist); in.segmentmarkerlist    =NULL;   
-  }
-  if(!in.regionlist){
-    free(in.regionlist); in.regionlist  =NULL;            
-  }
-  if(!in.holelist){
-    free(in.holelist); in.holelist=NULL;
-  }
-
-  //return Py_BuildValue((char *)"O", holder);
-  return Py_BuildValue((char *)"O", holder);
-}
-
-/* end of my function*/
-
-
-   
-    
-static PyMethodDef triang_methods[] = {
-  {(char *)"genMesh", triang_genMesh, 1, (char *)"Passes hull points to triangle.All arguments are lists.(<hull>,<seglist>,<holelist>,<regionlist>)-><hull>"}, 
-  {NULL,NULL}
-};
-
-void initmesh_engine(){
-  Py_InitModule((char *)"mesh_engine",triang_methods);
-}    
Index: /anuga_core/source/anuga/mesh_engine/mesh_engine.py
===================================================================
--- /anuga_core/source/anuga/mesh_engine/mesh_engine.py	(revision 4139)
+++ /anuga_core/source/anuga/mesh_engine/mesh_engine.py	(revision 4139)
@@ -0,0 +1,133 @@
+#!/usr/bin/env python
+
+import sys
+
+from types import ListType, TupleType
+
+import anuga.mesh_engine.mesh_engine_c_layer as triang
+
+from Numeric import array, Float, Int
+
+from anuga.utilities.numerical_tools import ensure_numeric
+from anuga.utilities.anuga_exceptions import ANUGAError
+    
+def generate_mesh(points=None,
+                  segments=None,holes=None,regions=None,
+                  pointatts=None,segatts=None,
+                  mode=None, dummy_test=None):
+    """
+    
+    """
+    #FIXME (DSG-DSG): Catch parameters that are lists,
+    #instead of lists of lists
+    # check shape[1] is 2 etc
+
+    if points is None:
+        points = []
+
+    if segments is None:
+        segments = []
+
+    if holes is None:
+        holes = []
+        
+    if regions is None:
+        regions = []
+
+    if dummy_test is None:
+        dummy_test  = []
+        
+    try:
+        points =  ensure_numeric(points, Float)
+    except ValueError:
+        msg = 'ERROR: Inconsistent points array.'
+        raise ANUGAError, msg
+    if points.shape[1] <>2:
+        msg = 'ERROR: Bad shape points array.'
+        raise ANUGAError, msg
+
+    # This is after points is numeric
+    if pointatts is None or pointatts == []:
+        pointatts = [[] for x in range(points.shape[0])]
+        
+    try:
+        segments = ensure_numeric(segments, Int)
+    except ValueError:
+        msg = 'ERROR: Inconsistent segments array.'
+        raise ANUGAError, msg
+    
+    # This is after segments is numeric
+    if segatts is None or segatts == []:
+        segatts = [0 for x in range(segments.shape[0])]
+        
+    try:
+        holes = ensure_numeric(holes, Float)
+    except ValueError:
+        msg = 'ERROR: Inconsistent holess array.'
+        raise ANUGAError, msg
+
+   
+    regions = add_area_tag(regions)
+    try:
+        regions = ensure_numeric(regions, Float)
+    except ValueError:
+        msg = 'ERROR: Inconsistent regions array.'
+        raise ANUGAError, msg
+        
+    if not regions.shape[0] == 0 and regions.shape[1] <= 2:
+        msg = 'ERROR: Bad shape points array.'
+        raise ANUGAError, msg
+    
+    try:
+        pointatts = ensure_numeric(pointatts, Float)
+    except ValueError:
+        msg = 'ERROR: Inconsistent point attributes array.'
+        raise ANUGAError, msg
+
+    if pointatts.shape[0] <> points.shape[0]:
+        msg = """ERROR: Point attributes array not the same shape as
+        point array."""
+        raise ANUGAError, msg
+    
+    try:
+        segatts = ensure_numeric(segatts, Int)
+    except ValueError:
+        msg = 'ERROR: Inconsistent point attributes array.'
+        raise ANUGAError, msg
+    if segatts.shape[0] <> segments.shape[0]:
+        msg = """ERROR: Segment attributes array not the same shape as
+        segment array."""
+        raise ANUGAError, msg
+    return triang.genMesh(points,segments,holes,regions,
+                          pointatts,segatts, mode, points)
+
+def add_area_tag(regions):
+    """
+    So, what is the format?
+    A list with
+    [x,y,region_tag,area] OR [x,y,region_tag]
+    if it's [x,y,region_tag], add a 4th element, value of 0.0.
+    """
+    if isinstance(regions, ListType):
+        for i, region in enumerate(regions):
+            if len(region) == 3:
+                if isinstance(region, TupleType):
+                    #FIXME: How do you convert a tuple to a list?
+                    # I can do it a stupid way..
+                    tuple = region[:]
+                    regions[i] = []
+                    for j in tuple:
+                        regions[i].append(j)
+                    regions[i].append(0.0)
+                else:
+                    regions[i].append(0.0)
+                    
+                # let ensure numeric catch this..    
+                #len(region) <= 2:
+                #msg = 'ERROR: Inconsistent regions array.'
+                #raise msg
+            #elif
+    return regions
+
+if __name__ == "__main__":
+    pass 
Index: /anuga_core/source/anuga/mesh_engine/mesh_engine_c_layer.c
===================================================================
--- /anuga_core/source/anuga/mesh_engine/mesh_engine_c_layer.c	(revision 4139)
+++ /anuga_core/source/anuga/mesh_engine/mesh_engine_c_layer.c	(revision 4139)
@@ -0,0 +1,643 @@
+/* this appears to be necessary to run on Linux */
+
+#undef _STDLIB_H_
+#define _STDLIB_H_
+
+
+/*
+
+    This code interfaces pmesh directly with "triangle", a general        
+    purpose triangulation code. In doing so, Python Numeric data structures          
+    are passed to C for  use by "triangle" and the results are then          
+    converted back. This was accomplished using the Python/C API.            
+                                                                           
+    I. Arrays of points,edges,holes, and regions must normally be 
+     created for proper triangulation. However, if there are no holes,
+     the  hole list need not be specified (however, a 
+     measure of control over the shape of the convex hull is lost; 
+     triangle decides). 
+     
+     points array: [(x1,y1),(x2,y2),...] ( doubles)                  
+     edge array: [(Ea,Eb),(Ef,Eg),...] ( integers)                   
+     hole array: [(x1,y1),...]( doubles, one inside each hole region)
+     region array: [ (x1,y1,index,area),...] ( doubles)                
+     pointattribute array:[ (a1,a2..),...] (doubles) 
+     mode: String - the commands sent to triangle
+          
+     This is all that is needed to accomplish an intial triangulation.       
+     A dictionary of the Triangle output is returned.                                   
+     Precondition
+     End list in the pointattributelist has to have the same length
+     
+                                                 */
+
+
+#ifdef SINGLE  
+#define REAL float 
+#else  
+#define REAL double 
+#endif 
+
+#include "triangle.h" 
+
+#ifndef _STDLIB_H_ 
+#ifndef _WIN32
+extern "C" void *malloc(); 
+extern "C" void free(); 
+#endif
+#endif  
+
+#include "Python.h" 
+
+#include "Numeric/arrayobject.h"
+
+
+ static PyObject *triang_genMesh(PyObject *self, PyObject *args){
+     
+  /* Mesh generation routine
+     pre-conditions:
+     mode must have a 'z' switch, to avoid segmentation faults
+     
+     shallow_water_ext.c gravity function
+  */
+    
+  struct triangulateio in,out;
+
+  struct triangulateio in_test;
+  
+  PyArrayObject *pointlist;
+  PyArrayObject *seglist;
+  PyArrayObject *holelist;
+  PyArrayObject *regionlist;
+  PyObject *mode;
+  PyArrayObject *pointattributelist;
+  PyArrayObject *segmarkerlist;
+  PyArrayObject *test;
+  REAL Attr;
+  int i, j, iatt, n, write_here;
+  int a,b,c;
+  int marker;
+  int tuplesize;
+    
+  int index = 0;
+  double x,y;
+  PyObject *elems;
+  PyObject *pair;
+     
+  char *mod;
+  PyObject *holder, *holderlist, *attributelist;
+  int listsize,attsize ;
+  PyObject  *ii;
+  
+  int *points_connected;
+  int *lone_verts;
+  
+  /* used for testing numeric arrays*/
+  int n_test;     
+  if(!PyArg_ParseTuple(args,(char *)"OOOOOOOO",&pointlist,&seglist,&holelist,&regionlist,&pointattributelist,&segmarkerlist,&mode,&test)){
+    return NULL;
+  }
+  /*if(!PyList_Check(hull)){
+    PyErr_SetString(PyExc_TypeError, 
+		    "incorrect first argument for 'genMesh': list required.");
+    return NULL;
+    } */
+     
+  /* printf("starting\n");   
+  n_test=test -> dimensions[0];
+  printf("test -> dimensions[0] %i\n", n_test); 
+  printf("test -> data %f\n", ((double *) test -> data)[0]); /*returns 1.0*/
+     
+  /* Initialize  points 
+  in.numberofpoints =  PyList_Size(hull); 
+  in.pointlist = (REAL *)malloc(in.numberofpoints*2*sizeof(REAL));
+  in.pointmarkerlist = (int *)NULL; /* malloc(in.numberofpoints*sizeof(int)); */
+  
+  /* Initialize  points */
+  in.numberofpoints =  pointlist-> dimensions[0]; 
+  in.pointlist = (double *) pointlist -> data;
+  in.pointmarkerlist = (int *)NULL; 
+  
+  /* Initialize and fill up region list */
+  in.numberofregions = regionlist-> dimensions[0]; 
+  in.regionlist = (double *) regionlist -> data;
+  
+  /* Initialize and fill up region list 
+  in.numberofregions = PyList_Size(regionlist);
+  
+  /*printf(" numberofregions %i\n",in.numberofregions );
+  in.regionlist = (REAL *)malloc(in.numberofregions*4*sizeof(REAL));
+  for(i=0; i<in.numberofregions;i++){
+    elems = PyList_GetItem(regionlist,i);
+    tuplesize = PySequence_Length(elems);
+    /*printf("tuplesize %i\n",tuplesize);
+    if (tuplesize>0){
+      in.regionlist[4*i] = PyFloat_AsDouble(PyTuple_GetItem(elems,0));
+    } else {
+      PyErr_SetString(PyExc_TypeError,
+		      "data sent to trianglemodule error: A region has no x value.");
+      return NULL;
+    }
+    if (tuplesize>1){
+      in.regionlist[4*i+1] = PyFloat_AsDouble(PyTuple_GetItem(elems,1));
+    } else {
+      PyErr_SetString(PyExc_TypeError,
+		    "data sent to trianglemodule error: A region has no y value.");
+      return NULL;
+    }
+    if (tuplesize>2){
+      /* 2 is attribute
+      in.regionlist[4*i+2] = PyFloat_AsDouble(PyTuple_GetItem(elems,2));
+    } else {
+      PyErr_SetString(PyExc_TypeError,
+		    "data sent to trianglemodule error: A region has no attribute value.");
+      return NULL;
+    }
+    if (tuplesize>3){
+      in.regionlist[4*i+3] = PyFloat_AsDouble(PyTuple_GetItem(elems,3));
+    } else {
+      /* 0.0 is interpreted as no local max area by triangle  
+      in.regionlist[4*i+3] = (REAL)0.0;
+    } 
+  }*/
+  
+  /*Initialize segments and segment markers */
+  in.numberofsegments = seglist -> dimensions[0]; 
+  in.segmentlist = (int *) seglist -> data;
+  in.segmentmarkerlist = (int *) segmarkerlist -> data;
+  
+  /*Initialize segments 
+  in.numberofsegments = PyList_Size(seglist);
+  in.segmentlist = (int *)malloc(in.numberofsegments*2*sizeof(int));
+  in.segmentmarkerlist = (int *)malloc(in.numberofsegments*sizeof(int)); */
+
+    
+  /*Initialize triangles */
+  in.numberoftriangles = 0;
+  in.trianglelist = (int *)NULL;
+  in.numberoftriangleattributes = 0;  
+  in.triangleattributelist  = (REAL *)NULL; 
+     
+  /*Initialize holes
+  in.numberofholes = PyList_Size(holelist);
+  if (in.numberofholes != 0) {
+    in.holelist = (REAL *)malloc(in.numberofholes*2*sizeof(REAL));
+  } else {
+    in.holelist = (REAL *)NULL;
+  }      */
+     
+  /*Initialize holes */
+  in.numberofholes = holelist  -> dimensions[0];
+  if (in.numberofholes != 0) {
+    in.holelist =  (double *) holelist -> data; 
+  } else {
+    in.holelist = (REAL *)NULL;
+  }     
+  
+      
+  /* Fill up segment list 
+  index = 0;
+  if (in.numberofsegments != PyList_Size(segmarkerlist)) {
+    PyErr_SetString(PyExc_TypeError,
+		    "data sent to trianglemodule error:number of segments doesn't = number of segment attributes.");
+    return NULL;
+  }
+    
+  for(i = 0; i < in.numberofsegments; i++) { */
+    /*
+    pair =  PyList_GetItem(seglist,i);
+    if (2 != PySequence_Length(pair)) {
+      PyErr_SetString(PyExc_TypeError,
+		      "data sent to trianglemodule error:Found a segment without 2 vertices.");
+      return NULL;
+    }
+    a = (int)PyInt_AsLong(PyTuple_GetItem(pair,0));
+    in.segmentlist[index] = a;
+    index++;
+    b = (int)PyInt_AsLong(PyTuple_GetItem(pair,1));
+    in.segmentlist[index] = b;
+    index++;
+    */
+    /* Fill up segment marker list 
+    marker =  (int)PyInt_AsLong(PyList_GetItem(segmarkerlist,i));
+    in.segmentmarkerlist[i] = marker; 
+  } */
+  /* Fill up hole list 
+  index = 0;
+  for(i = 0; i < in.numberofholes; i++) {
+    pair =  PyList_GetItem(holelist,i);
+    if (2 != PySequence_Length(pair)) {
+      PyErr_SetString(PyExc_TypeError,
+		      "data sent to trianglemodule error:Found a pair without 2 vertices.");
+      return NULL;
+    }
+    x = PyFloat_AsDouble(PyTuple_GetItem(pair,0));
+    in.holelist[index] = x;
+    index++;
+    y = PyFloat_AsDouble(PyTuple_GetItem(pair,1));
+    in.holelist[index] = y;
+    index++;
+  } */
+     
+  /* Fill up point list
+  index = 0;
+  for(i = 0; i < in.numberofpoints; i++) {
+    pair =  PyList_GetItem(hull,i);
+    if (2 != PySequence_Length(pair)) {
+      PyErr_SetString(PyExc_TypeError,
+		      "data sent to trianglemodule error:Found a vertex without an x,y values.");
+      return NULL;
+    }
+    x = PyFloat_AsDouble(PyTuple_GetItem(pair,0));
+    in.pointlist[index] = x;
+    index++;
+    y = PyFloat_AsDouble(PyTuple_GetItem(pair,1));
+    in.pointlist[index] = y;
+    index++;
+  }  */
+  
+  /* for(i = 0; i < in.numberofpoints * 2; i++) {
+    printf("in.pointlist[index] %f\n", in.pointlist[i]);  
+    } */
+  /*
+  for(i = 0; i < in.numberofpoints * 2; i++) {
+    printf("in_test.pointlist[index] %f\n", in_test.pointlist[i]);  
+  }
+  
+  printf("in_test.numberofpoints %i\n", in_test.numberofpoints); 
+  printf("test -> data %f\n", ((double *) test -> data)[1,1]); /*returns 1.0*/
+
+  
+  /* Initialize  points
+  in.numberofpoints = in_test.numberofpoints;
+  in.pointlist = in_test.pointlist;
+  in.pointmarkerlist =  in_test.pointmarkerlist;  */
+  
+  /* Initialize  point attributes */
+/*      firstPointAttr = PyList_GetItem(pointattributelist,0); */
+/*      if (firstPointAttr != NULL) { */
+/*        //in.numberofpointattributes = PyList_Size(firstPointAttr); */
+/*        printf("in.numberofpointattributes =%i",PyList_Size(firstPointAttr)); */
+/*        //printf("in.numberofpointattributes =%i",in.numberofpointattributes); */
+/*      } else { */
+/*        printf("(firstPointAttr is NULL\n"); */
+/*        in.pointattributelist = NULL; */
+/*        in.numberofpointattributes = 0; */
+/*      } */
+
+ /* Point attribute list */
+  if (0 == pointattributelist -> dimensions[0]) {
+    in.numberofpointattributes = 0;
+    in.pointattributelist =  NULL;
+  } else {
+    if (0 == pointattributelist -> dimensions[1]) {
+    in.numberofpointattributes = 0;
+    in.pointattributelist =  NULL;
+    } else {
+      in.numberofpointattributes = pointattributelist -> dimensions[1];
+      in.pointattributelist =  (double *) pointattributelist -> data;
+    }
+  }
+   
+  /* 
+  n = PyList_Size(pointattributelist);
+  if (n < 0) {
+    in.numberofpointattributes = 0; 
+    in.pointattributelist = (REAL *)NULL;
+  } else {
+    /* Get the # of attributes based on the size of the first point
+    firstPointAttr = PyList_GetItem(pointattributelist,0);
+    if (firstPointAttr == NULL) {
+      in.numberofpointattributes = 0; 
+      in.pointattributelist = (REAL *)NULL;
+    } else {
+      in.numberofpointattributes = PyList_Size(firstPointAttr);
+      in.pointattributelist = (REAL *)malloc(in.numberofpoints
+					     *in.numberofpointattributes
+					     *sizeof(REAL)); 
+     
+     
+      /* fill the point attribute list 
+      if (in.numberofpointattributes != 0) {
+	for(j = 0; j < in.numberofpoints; j++){
+	  Attrlist = PyList_GetItem(pointattributelist,j);
+	  if (in.numberofpointattributes != PySequence_Length(Attrlist)) {
+	    PyErr_SetString(PyExc_TypeError,
+			    "data sent to trianglemodule error:List size of attributes is different from vertex to vertex.");
+	    return NULL;
+	  }
+	  for(i = 0; i < in.numberofpointattributes; i++){
+	    Attr = PyFloat_AsDouble(PyList_GetItem(Attrlist,i));
+	    in.pointattributelist[in.numberofpointattributes*j+i]= (REAL)Attr;
+	    /* 	 printf("i =%i\n",i);  
+	    /* 	 printf("j =%i\n",j);  
+	    /* 	 printf("in.pointattributelist[4*j+i]
+		 =%f\n",in.pointattributelist[4*j+i]); 
+	  }
+	}
+      } 
+    }
+    } */
+    
+     
+  /* this half from triangulate*/
+  /* set up the switch arguments to the triangulation routine */
+
+  /* DEBUGGING
+  printf(" ***  hello world\n" );  
+  
+  printf ("numberofholess %i\n", in.numberofholes);
+  for(i = 0; i < in.numberofholes; i++) {
+    printf("(%f,%f)\n",in.holelist[i* 2+ 0],in.holelist[i* 2+ 1]);
+  }
+     
+      printf ("numberoftriangles %i\n", in.numberoftriangles);
+      for(i = 0; i < in.numberoftriangles; i++) {
+      printf("(%i,%i,%i)",in.trianglelist[i* 3+ 0],in.trianglelist[i* 3+ 1], in.trianglelist[i* 3+ 2]);
+      }
+ 
+  printf(" ***  see ya world\n" );       */
+ 
+      
+  mod = PyString_AS_STRING(mode);
+     
+         
+  out.pointlist = (REAL *)NULL;
+  out.pointmarkerlist = (int *)NULL;
+  out.numberofpointattributes = in.numberofpointattributes;
+  out.pointattributelist = (REAL *)NULL;
+  
+  out.trianglelist = (int *)NULL;
+  out.triangleattributelist = (REAL *)NULL;                   
+  out.trianglearealist = (REAL *)NULL;                        
+  out.neighborlist = (int *)NULL;
+     
+  out.segmentlist = (int *)NULL;
+  out.segmentmarkerlist = (int *)NULL;
+     
+  out.edgelist = (int *)NULL;
+  out.edgemarkerlist = (int *)NULL;
+     
+  out.holelist = (REAL *)NULL;
+  out.regionlist = (REAL *)NULL;
+    
+  
+  /*printf("\n\nTriangulate input args: %s \n\n", mod); */
+  triangulate(mod, &in, &out, (struct triangulateio *)NULL );
+
+  /*
+    ------- Pass point numbers,coordinates and neighbors back to Python ------
+    we send back a dictionary:                                               
+    { index : [ coordinates, [connections], Attribute ] } 
+  */
+  holder = PyDict_New();
+     
+  /* list of int's, used to keep track of which verts are connected to
+     triangles. */
+  points_connected = (int *)malloc(out.numberofpoints*sizeof(int));
+  /* lone_verts = (int *)malloc(out.numberofpoints*sizeof(int)); */ 
+  
+  /* Initialise lone vert list */
+  for(i=0; i<out.numberofpoints;i++){
+    points_connected[i] = 0;
+    /* lone_verts[i] = 0; */
+  }  
+  
+  /* Add triangle list */
+  listsize = out.numberoftriangles;
+  holderlist = PyList_New(listsize);
+  for(i=0; i<listsize;i++){
+    PyObject *mlist = Py_BuildValue((char *)"(i,i,i)", 
+				    out.trianglelist[i*3  ], out.trianglelist [i*3+1], out.trianglelist [i*3+2]);    
+    PyList_SetItem(holderlist,i, mlist);
+    /* printf(" A vert index %i\n",out.trianglelist[i*3] );
+    printf(" A vert index %i\n",out.trianglelist[i*3+1] );
+    printf(" A vert index %i\n",out.trianglelist[i*3+2] ); */
+    points_connected[out.trianglelist[i*3]] = 1;
+    points_connected[out.trianglelist[i*3+1]] = 1;
+    points_connected[out.trianglelist[i*3+2]] = 1;
+    /* lone_verts[out.trianglelist[i*3]] = 1;
+    lone_verts[out.trianglelist[i*3+1]] = 1;
+    lone_verts[out.trianglelist[i*3+2]] = 1; */
+    /*  printf(" Add triangle list \n");*/
+  }    
+  ii=PyString_FromString("generatedtrianglelist");
+  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);
+  
+  /* convert the points_connected vector from a true(1) false(0) vector, where
+     index is the vert, to a vector of the lone verts, at the beggining
+     of the vector. */
+  write_here = 0;   
+  for(i=0; i<out.numberofpoints;i++){
+    /* lone_verts[i-write_here] = lone_verts[i]; */
+    if (0 == points_connected[i]) {
+      points_connected[write_here] = i;
+      write_here ++;
+    }
+  }   
+  /* printf(" ******************** \n" );
+  for(i=0; i<write_here;i++){
+    printf(" A vert index %i\n",points_connected[i] );
+    } */
+  
+  listsize = write_here;
+  holderlist = PyList_New(listsize);
+  for(i=0; i<listsize;i++){
+   PyObject *mlist = Py_BuildValue((char *)"i", points_connected[i]);    
+    PyList_SetItem(holderlist,i, mlist);
+  }
+  ii=PyString_FromString("lonepointlist");
+  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist); 
+     
+  
+  /* Add pointlist */
+  listsize = out.numberofpoints;
+  holderlist = PyList_New(listsize);
+     
+  for(i=0; i<out.numberofpoints;i++){
+    PyObject *mlist = Py_BuildValue((char *)"(d,d)", 
+				    out.pointlist[i*2  ], out.pointlist[i*2+1]);  
+    PyList_SetItem(holderlist,i, mlist);
+  }  
+  ii=PyString_FromString("generatedpointlist");
+  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);
+  
+  /* Add point marker list */
+  listsize = out.numberofpoints;
+  holderlist = PyList_New(listsize);
+  
+  for(i=0; i<listsize;i++){
+    PyObject *mlist = Py_BuildValue((char *)"d", 
+				    out.pointmarkerlist[i]);     
+    PyList_SetItem(holderlist,i, mlist); 
+  }  
+  ii=PyString_FromString("generatedpointmarkerlist");
+  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);
+    
+  /* Add point attribute list */
+  listsize = out.numberofpoints;
+  holderlist = PyList_New(listsize);
+  
+  for(i=0; i<listsize;i++){
+    attsize = out.numberofpointattributes;
+    attributelist = PyList_New(attsize);
+    for(iatt=0; iatt<attsize;iatt++){
+      PyObject *mlist = Py_BuildValue((char *)"d", 
+				      out.pointattributelist[i*attsize + iatt]);     
+      PyList_SetItem(attributelist,iatt, mlist);
+    }      
+    PyList_SetItem(holderlist,i, attributelist);
+  }  
+  ii=PyString_FromString("generatedpointattributelist");
+  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);  
+ 
+  
+  /* Add triangle attribute list */
+  listsize = out.numberoftriangles;
+  holderlist = PyList_New(listsize);
+     
+  for(i=0; i<listsize; i++){
+    attsize = out.numberoftriangleattributes;
+    attributelist = PyList_New(attsize);       
+    for(iatt=0; iatt<attsize;iatt++){
+      PyObject *mlist = Py_BuildValue((char *)"d",out.triangleattributelist[i*attsize + iatt]);  
+      PyList_SetItem(attributelist,iatt, mlist);
+    }      
+    PyList_SetItem(holderlist,i, attributelist);
+  }  
+  ii=PyString_FromString("generatedtriangleattributelist");
+  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);
+  
+  /* Add segment list */
+  listsize = out.numberofsegments;
+  holderlist = PyList_New(listsize);
+  for(i=0; i<listsize;i++){
+    PyObject *mlist = Py_BuildValue((char *)"(i,i)", 
+				    out.segmentlist[i*2  ],
+				    out.segmentlist [i*2+1]);
+    PyList_SetItem(holderlist,i, mlist);
+  }    
+  ii=PyString_FromString("generatedsegmentlist");
+  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);  
+  
+  /* Add segment marker list */
+  listsize = out.numberofsegments;
+  holderlist = PyList_New(listsize);
+  for(i=0; i<listsize;i++){
+    PyObject *mlist = Py_BuildValue((char *)"i", 
+				    out.segmentmarkerlist[i]);
+    PyList_SetItem(holderlist,i, mlist);
+  }    
+  ii=PyString_FromString("generatedsegmentmarkerlist");
+  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);  
+  /* Add triangle neighbor list */
+  if (out.neighborlist != NULL) {
+    listsize = out.numberoftriangles;
+    holderlist = PyList_New(listsize);
+    for(i=0; i<listsize;i++){
+      PyObject *mlist = Py_BuildValue((char *)"(i,i,i)", 
+				      out.neighborlist[i*3  ],
+				      out.neighborlist [i*3+1],
+				      out.neighborlist [i*3+2]);    
+      PyList_SetItem(holderlist,i, mlist);
+    }    
+    ii=PyString_FromString("generatedtriangleneighborlist");
+    PyDict_SetItem(holder, ii, holderlist);
+    Py_DECREF(ii); Py_DECREF(holderlist);
+  }   
+     
+  /* Free in/out structure memory */
+  
+  /* OUT */
+
+  if(!out.pointlist){
+    free(out.pointlist);  out.pointlist=NULL;
+  }
+  if(!out.pointmarkerlist){    
+    free(out.pointmarkerlist); out.pointmarkerlist=NULL;
+  }
+  if(!out.pointattributelist){    
+    free(out.pointattributelist); out.pointattributelist=NULL;
+  }
+    
+ 
+  if(!out.trianglelist){    
+    free(out.trianglelist); out.trianglelist=NULL;
+  }
+  if(!out.triangleattributelist){    
+    free(out.triangleattributelist); out.triangleattributelist=NULL;
+  }
+  if(!out.trianglearealist){    
+    free(out.trianglearealist); out.trianglearealist=NULL;
+  }
+  if(!out.neighborlist){    
+    free(out.neighborlist); out.neighborlist=NULL;
+  }
+  if(!out.segmentlist){
+    free(out.segmentlist); out.segmentlist =NULL;
+  }
+  if(!out.segmentmarkerlist){
+    free(out.segmentmarkerlist); out.segmentmarkerlist  =NULL;
+  }
+  if(!out.edgelist){
+    free(out.edgelist);  out.edgelist=NULL;
+  }
+  if(!out.edgemarkerlist){
+    free(out.edgemarkerlist);  out.edgemarkerlist=NULL;
+  }
+  if(!out.holelist){
+    free(out.holelist); out.holelist=NULL;
+  }
+  if(!out.regionlist){
+    free(out.regionlist); out.regionlist=NULL;
+  }
+  if(!points_connected ){
+    free(points_connected ); points_connected =NULL;
+  }
+    
+  
+  
+  /*  IN */
+  /* These are all numeric structures, so I shouldn't free them,
+  but should I decrease the ref?
+  Based on examples of this stuff I don't have to decrease
+  the ref
+    
+  if(!in.pointlist){
+    free(in.pointlist);  in.pointlist   =NULL;             
+  }
+  if(!in.pointattributelist){
+    free(in.pointattributelist); in.pointattributelist =NULL;       
+  }
+  if(!in.pointmarkerlist){
+    free(in.pointmarkerlist); in.pointmarkerlist    =NULL;        
+  }
+  if(!in.segmentlist){
+    free(in.segmentlist);  in.segmentlist    =NULL;         
+  }
+  if(!in.segmentmarkerlist){
+    free(in.segmentmarkerlist); in.segmentmarkerlist    =NULL;   
+  }
+  if(!in.regionlist){
+    free(in.regionlist); in.regionlist  =NULL;            
+  }
+  if(!in.holelist){
+    free(in.holelist); in.holelist=NULL;
+    } */
+
+  return Py_BuildValue((char *)"O", holder);
+}
+
+/* end of my function*/
+
+
+   
+    
+static PyMethodDef triang_methods[] = {
+  {(char *)"genMesh", triang_genMesh, 1, (char *)"Passes hull points to triangle.All arguments are lists.(<hull>,<seglist>,<holelist>,<regionlist>)-><hull>"}, 
+  {NULL,NULL}
+};
+
+void initmesh_engine_c_layer(){
+  Py_InitModule((char *)"mesh_engine_c_layer",triang_methods);
+}    
Index: /anuga_core/source/anuga/mesh_engine/test_generate_mesh.py
===================================================================
--- /anuga_core/source/anuga/mesh_engine/test_generate_mesh.py	(revision 4139)
+++ /anuga_core/source/anuga/mesh_engine/test_generate_mesh.py	(revision 4139)
@@ -0,0 +1,526 @@
+#!/usr/bin/env python
+#
+
+import sys
+
+import unittest
+from anuga.mesh_engine.mesh_engine import generate_mesh
+
+from Numeric import array, Float, Int
+
+from anuga.utilities.numerical_tools import ensure_numeric
+
+from anuga.utilities.anuga_exceptions import ANUGAError
+
+class triangTestCase(unittest.TestCase):
+    def setUp(self):
+        pass
+
+    def tearDown(self):
+        pass
+
+    def metestrectangleII(self):
+
+        points = []
+        seglist = []
+        holelist = []
+        regionlist = []
+
+	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
+	pointattlist = [[],[],[],[]]
+	regionlist.append( (1.2,1.2,5.0) )
+	seglist = [(0,1),(1,3),(3,2),(2,0)]
+        segattlist = [0,0,0,0]
+        
+        mode = "Qzp"
+        
+        points =  ensure_numeric(points, Float)
+        seglist = test = ensure_numeric(seglist, Int)
+        data = triang.genMesh(points,seglist,holelist,regionlist,
+                              pointattlist,segattlist, mode, test)
+
+        self.failUnless(data['generatedtrianglelist'] ==[(1, 0, 2), (2, 3, 1)],
+                        'trianglelist is wrong!')
+        self.failUnless(data['generatedsegmentlist'] ==[(0, 1), (1, 3),
+                                                        (3, 2), (2, 0)],
+                        'segmentlist is wrong!')
+        self.failUnless(data['generatedpointlist'] ==[(0.0, 0.0), (0.0, 10.0),
+                                                      (3.0, 0.0), (3.0, 10.0)],
+                        ' is wrong!')
+
+
+    def testrectangleIIb(self):
+        # segattlist = []
+        points = []
+        seglist = []
+        holelist = []
+        regionlist = []
+
+	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
+	pointattlist = []
+	regionlist.append( (1.2,1.2,5.0) )
+	seglist = [(0,1),(1,3),(3,2),(2,0)]
+        segattlist = []
+        
+        mode = "Qzp"
+        
+        data = generate_mesh(points,seglist,holelist,regionlist,
+                              pointattlist,segattlist, mode, points)
+
+        self.failUnless(data['generatedtrianglelist'] ==[(1, 0, 2), (2, 3, 1)],
+                        'trianglelist is wrong!')
+        self.failUnless(data['generatedsegmentlist'] ==[(0, 1), (1, 3),
+                                                        (3, 2), (2, 0)],
+                        'segmentlist is wrong!')
+        self.failUnless(data['generatedpointlist'] ==[(0.0, 0.0), (0.0, 10.0),
+                                                      (3.0, 0.0), (3.0, 10.0)],
+                        ' is wrong!')
+     
+    def testrectangleIIc(self):
+        # segattlist = None
+        points = []
+        seglist = []
+        holelist = []
+        regionlist = []
+
+	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
+	pointattlist = None
+	regionlist.append( (1.2,1.2,5.0) )
+	seglist = [(0,1),(1,3),(3,2),(2,0)]
+        segattlist = None
+        
+        mode = "Qzp"
+        
+        data = generate_mesh(points,seglist,holelist,regionlist,
+                              pointattlist,segattlist, mode, points)
+
+        self.failUnless(data['generatedtrianglelist'] ==[(1, 0, 2), (2, 3, 1)],
+                        'trianglelist is wrong!')
+        self.failUnless(data['generatedsegmentlist'] ==[(0, 1), (1, 3),
+                                                        (3, 2), (2, 0)],
+                        'segmentlist is wrong!')
+        self.failUnless(data['generatedpointlist'] ==[(0.0, 0.0), (0.0, 10.0),
+                                                      (3.0, 0.0), (3.0, 10.0)],
+                        ' is wrong!')
+        
+    def test_bad_point_atts(self):
+
+        points = []
+        seglist = []
+        holelist = []
+        regionlist = []
+
+	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
+	pointattlist = [[],[],[]]
+	regionlist.append( (1.2,1.2,5.0) )
+	seglist = [(0,1),(1,3),(3,2),(2,0)]
+        segattlist = [0,0,0,0]
+        
+        mode = "Qzp"
+        
+        try:
+            data = generate_mesh(points,seglist,holelist,
+                                  regionlist,pointattlist,segattlist,
+                                   mode, points)
+        except ANUGAError:
+            pass
+        else:
+            self.failUnless(0 ==1,
+                        'bad list did not raise error!')
+            
+        
+    def test_bad_point_attsII(self):
+
+        points = []
+        seglist = []
+        holelist = []
+        regionlist = []
+
+	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
+	pointattlist = [[1],[2],[3],[4,8]]
+	regionlist.append( (1.2,1.2,5.0) )
+	seglist = [(0,1),(1,3),(3,2),(2,0)]
+        segattlist = [0,0,0,0]
+        
+        mode = "Qzp"
+        
+        try:
+            data = generate_mesh(points,seglist,holelist,
+                                  regionlist,pointattlist,segattlist,
+                                   mode, points)
+        except ANUGAError:
+            pass
+        else:
+            self.failUnless(0 ==1,
+                        'bad list did not raise error!')
+            
+    def testsegmarker(self):
+
+        points = []
+        seglist = []
+        holelist = []
+        regionlist = []
+
+	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
+	pointattlist = [[],[],[],[]]
+	regionlist.append( (1.2,1.2,5.0) )
+	seglist = [(0,1),(1,3),(3,2),(2,0)]
+        segattlist = [1.0,2.0,3.0,4.0]
+        
+        mode = "Qzp"
+        data = generate_mesh(points,seglist,holelist,regionlist,
+                              pointattlist,segattlist, mode, points)
+
+        self.failUnless(data['generatedtrianglelist'] ==[(1, 0, 2), (2, 3, 1)],
+                        'trianglelist is wrong!')
+        self.failUnless(data['generatedsegmentlist'] ==[(0, 1), (1, 3),
+                                                        (3, 2), (2, 0)],
+                        'segmentlist is wrong!')
+        self.failUnless(data['generatedpointlist'] ==[(0.0, 0.0), (0.0, 10.0),
+                                                      (3.0, 0.0), (3.0, 10.0)],
+                        'generatedpointlist is wrong!')
+        self.failUnless(data['generatedsegmentmarkerlist'] ==[1,2,3,4],
+                        'generatedsegmentmarkerlist is wrong!')
+
+    def testbad_region(self):
+
+        points = []
+        seglist = []
+        holelist = []
+        regionlist = [(1.2,1.2)]
+
+	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
+	pointattlist = [[],[],[],[]]
+	seglist = [(0,1),(1,3),(3,2),(2,0)]
+        segattlist = [0,0,0,0]
+        
+        mode = "Qzpn"
+        try:
+            data = generate_mesh(points,seglist,holelist,
+                                  regionlist,pointattlist,segattlist,
+                                   mode, points)
+        except ANUGAError:
+            pass
+        else:
+            self.failUnless(0 ==1,
+                        'bad region list did not raise error!')
+
+    def testbad_regionII(self):
+
+        points = []
+        seglist = []
+        holelist = []
+        regionlist = [(1.2,1.2), (1.2,1.25,1.0)]
+
+	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
+	pointattlist = [[],[],[],[]]
+	seglist = [(0,1),(1,3),(3,2),(2,0)]
+        segattlist = [0,0,0,0]
+        
+        mode = "Qzpn"
+        try:
+            data = generate_mesh(points,seglist,holelist,
+                                  regionlist,pointattlist,segattlist,
+                                   mode, points)
+        except ANUGAError:
+            pass
+        else:
+            self.failUnless(0==1,
+                        'bad region list did not raise error!')
+
+    def testregion_with_maxarea(self):
+
+        points = []
+        seglist = []
+        holelist = []
+        regionlist = [(3,1,1.0)]
+
+	points = [(0.0,0.0),(6.0,0.0),(6.0,6.0),(0.0,6.0)]
+	pointattlist = [[],[],[],[]]
+	seglist = [(0,1),(1,2),(3,2),(3,0),(0,2)]
+        segattlist = [0,0,0,0,0]
+        
+        mode = "Qzpna36a"
+        data = generate_mesh(points,seglist,holelist,regionlist,
+                              pointattlist,segattlist,  mode, points)
+
+        self.failUnless(len(data['generatedtrianglelist']) == 2,
+                        'testregion_with_maxarea 1: # of tris is wrong!')
+        ## Another test case
+        regionlist = [(3,1,1.0),(1,3,1.0,8.0)]
+        mode = "Qzp21na36a"
+        data = generate_mesh(points,seglist,holelist,regionlist,
+                              pointattlist,segattlist,  mode, points)
+        #print "len(data['generatedtrianglelist']",len(data['generatedtrianglelist'])
+        # for Duncan On unix this returns a 7 triangle result.
+        # for Duncan on Windows returns a 6 triangle result.
+        # for Ole on nautilus this returns 6
+        # for Duncan on nautilus this returns 7
+        # ??, it seems to be the results from triangle that is
+        # causing the different results, and we are treating
+        # triangle as a back box.
+
+        self.failUnless(len(data['generatedtrianglelist']) >= 6,
+                        'testregion_with_maxarea 2: # of tris is wrong!')
+        ## Another test case
+        regionlist = [(3,1,1.0,8.0),(1,3,1.0,8.0)]
+        mode = "Qzpna36a"
+        data = generate_mesh(points,seglist,holelist,regionlist,
+                              pointattlist,segattlist,  mode, points)
+        #print "len(data['generatedtrianglelist']",len(data['generatedtrianglelist'])
+        # On unix this returns a 10 triangle result.
+        # Windows returns a 8 triangle result.
+        self.failUnless(len(data['generatedtrianglelist']) >= 8,
+                        'testregion_with_maxarea 3: # of tris is wrong!')
+
+        ## Another test case
+        regionlist = [(3,1,1.0),(1,3,1.0,8.0)]
+        mode = "Qzpna8a"
+        data = generate_mesh(points,seglist,holelist,regionlist,
+                              pointattlist,segattlist,  mode, points)
+        #print "len(data['generatedtrianglelist']",len(data['generatedtrianglelist'])
+        # On unix this returns a 10 triangle result.
+        # Windows returns a 8 triangle result.
+
+        self.failUnless(len(data['generatedtrianglelist']) >= 8,
+                        'testregion_with_maxarea 4: # of tris is wrong!')
+
+
+    def FIXME_testbad_hole(self):
+
+        holelist = [(9.0)]
+        regionlist = []
+
+	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
+	pointattlist = [[],[],[],[]]
+	seglist = [(0,1),(1,3),(3,2),(2,0)]
+        segattlist = [0,0,0,0]
+        
+        mode = "Qzpn"
+        try:
+            data = generate_mesh(points,seglist,holelist,regionlist,
+                                  pointattlist,segattlist,  mode, points)
+
+        except TypeError:
+            pass
+        else:
+            self.failUnless(0 ==1,
+                        'bad hole list did not raise error!')
+
+
+    def testbad_segattlist(self):
+
+        holelist = []
+        regionlist = []
+
+	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
+	pointattlist = [[],[],[],[]]
+	seglist = [(0,1),(1,3),(3,2),(2,0)]
+        segattlist = [0,0]
+        
+        mode = "Qzpn"
+        try:
+            data = generate_mesh(points,seglist,holelist,regionlist,
+                                  pointattlist,segattlist, mode, points)
+            print "data",data
+            self.failUnless(data['generatedtrianglelist'] ==[(1, 0, 2), (2, 3, 1)],
+                        'trianglelist is wrong!')
+        except ANUGAError:
+            pass
+        else:
+            self.failUnless(0 ==1,
+                        'bad segment attribute list did not raise error!')
+
+    def testrectangle_regions(self):
+
+        points = []
+        seglist = []
+        holelist = []
+        regionlist = []
+
+	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
+	pointattlist = [[],[],[],[]]
+        #it seems that
+        #triangle only associates one region with a triangle
+	regionlist.append( (1.3,1.3,88.33) )
+	regionlist.append( (1.2,1.2,77,55) )
+	seglist = [(0,1),(1,3),(3,2),(2,0)]
+        segattlist = [0,0,0,0]
+        mode = "QAzpq"
+        data = generate_mesh(points,seglist,holelist,regionlist,
+                              pointattlist,segattlist, mode, points)
+
+        self.failUnless(data['generatedtriangleattributelist'] ==[[77.0], [77.0], [77.0], [77.0]],
+                        'triangleattributelist is wrong!')
+
+    def testrectangle_regionsII(self):
+
+        points = []
+        seglist = []
+        holelist = []
+        regionlist = []
+
+	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
+	pointattlist = [[],[],[],[]]
+        #it seems that
+        #triangle only associates one region with a triangle
+	regionlist.append( [1.3,1.3,88.33] )
+	regionlist.append( [1.2,1.2,77,55] )
+	seglist = [(0,1),(1,3),(3,2),(2,0)]
+        segattlist = [0,0,0,0]
+        
+        mode = "QAzpq"
+        data = generate_mesh(points,seglist,holelist,regionlist,
+                              pointattlist,segattlist, mode, points)
+
+        self.failUnless(data['generatedtriangleattributelist'] ==[[77.0], [77.0], [77.0], [77.0]],
+                        'triangleattributelist is wrong!')
+
+
+    def test_lone_verts(self):
+        points = []
+        seglist = []
+        holelist = []
+        regionlist = []
+
+	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0),(0.0,10.0)]
+	pointattlist = [[],[],[],[],[]]
+	regionlist.append( (1.2,1.2,5.0) )
+	seglist = [(0,4),(4,1),(1,3),(3,2),(2,0)]
+        segattlist = [0,0,0,0,0]
+        
+        mode = "QpznAa2000.1a"
+        data = generate_mesh(points,seglist,holelist,regionlist,
+                              pointattlist,segattlist, mode, points)
+        #print "data['generatedtrianglelist']", data['generatedtrianglelist']
+        self.failUnless(data['generatedtrianglelist'] ==[(4, 0, 2), (2, 3, 4)],
+                        'trianglelist is wrong!')
+        self.failUnless(data['generatedsegmentlist'] ==[(0, 4), (4, 3),
+                                                        (3, 2), (2, 0)],
+                        'segmentlist is wrong!')
+        self.failUnless(data['generatedpointlist'] ==[(0.0, 0.0), (0.0, 10.0),
+                                                      (3.0, 0.0), (3.0, 10.0),
+                                                      (0.0,10.0)],
+                        ' is wrong!')
+        #print "data['lonepointlist']", data['lonepointlist']
+        self.failUnless(data['lonepointlist'] ==[1],
+                        'lonepointlist is wrong!')
+
+    def test_lone_vertsII(self):
+
+        points = []
+        seglist = []
+        holelist = []
+        regionlist = []
+
+	points = [(0.0,0.0),(0.0,0.0),(0.0,10.0),(0.0,10.0),(10.0,10.0),
+                  (10.0,10.0),(0.0,10.0),(10.0,0.0)]
+
+	pointattlist = []
+        for point in points:
+            pointattlist.append([])
+	seglist = [(0,1),(1,2),(2,3),(3,4),(4,5),(5,7),(7,0)]
+	segattlist = []
+        for seg in seglist:
+            segattlist.append(0)
+        
+        mode = "QpznAa2000.1a"
+        data = generate_mesh(points,seglist,holelist,regionlist,
+                              pointattlist,segattlist, mode, points)
+        #print "data['generatedtrianglelist']", data['generatedtrianglelist']
+        self.failUnless(data['generatedtrianglelist'] ==[(6, 1, 7), (7, 5, 6)],
+                        'trianglelist is wrong!')
+        self.failUnless(data['generatedsegmentlist'] ==[(1, 6), (6, 5),
+                                                        (5, 7), (7, 1)],
+                        'segmentlist is wrong!')
+        self.failUnless(data['generatedpointlist'] ==[(0.0,0.0),(0.0,0.0),
+                                                      (0.0,10.0),(0.0,10.0),
+                                                      (10.0,10.0),
+                                                      (10.0,10.0),(0.0,10.0),
+                                                      (10.0,0.0)],
+                        ' is wrong!')
+        self.failUnless(data['lonepointlist'] ==[0,2,3,4],
+                        'lonepointlist is wrong!')
+
+    def test_lone_vertsIII(self):
+
+        points = []
+        seglist = []
+        holelist = []
+        regionlist = []
+
+	points = [(-3.,-5.),(0.0,0.0),(0.0,0.0),
+                  (-5.,-8.),(0.0,10.0),(0.0,10.0),
+                  (-3.,-5.),(10.0,10.0), (10.0,10.0),
+                  (0.0,10.0),
+                  (10.0,0.0),
+                  (-12.,45.)]
+
+	pointattlist = []
+        for point in points:
+            pointattlist.append([])
+	seglist = [(1,2),(2,4),(4,5),(5,7),(7,8),(8,10),(10,1)]
+	segattlist = []
+        for seg in seglist:
+            segattlist.append(0)
+            """
+            0 1
+            1 2
+            2 4
+            3 5
+            4 7
+            5 8
+            6 9
+            7 10
+            """
+        
+        mode = "QpznAa2000.1a"
+        data = generate_mesh(points,seglist,holelist,regionlist,
+                              pointattlist,segattlist, mode, points)
+        #print "data['generatedtrianglelist']", data['generatedtrianglelist']
+        #self.failUnless(data['generatedtrianglelist'] ==[(9, 2, 10),
+         #                                                (10, 8, 9)],
+          #              'trianglelist is wrong!')
+        self.failUnless(data['generatedtrianglelist'] ==[(2, 10, 8),
+                                                         (8, 9, 2)],
+                        'trianglelist is wrong!')
+        #print "data['generatedsegmentlist']",data['generatedsegmentlist']
+        self.failUnless(data['generatedsegmentlist'] ==[(9, 2), (9, 8),
+                                                        (8, 10), (2, 10)],
+                        'segmentlist is wrong!')
+        self.failUnless(data['generatedpointlist'] ==points,
+                        ' is wrong!')
+        self.failUnless(data['lonepointlist'] ==[0,1,3,4,5,6,7,11],
+                        'lonepointlist is wrong!')
+        
+    def test_lone_verts4(self):
+        points = []
+        seglist = []
+        holelist = []
+        regionlist = []
+
+	points = [(-56.,-78),(0.0,0.0),(0.0,10.0),(10.,0.)]
+	pointattlist = []
+	seglist = [(1,2),(2,3),(3,1)]
+        segattlist = []
+        
+        mode = "QpznAa2000.1a"
+        data = generate_mesh(points,seglist,holelist,regionlist,
+                              pointattlist,segattlist, mode, points)
+        #print "data['generatedtrianglelist']", data['generatedtrianglelist']
+        #self.failUnless(data['generatedtrianglelist'] ==[(4, 0, 2),(2, 3, 4)],
+         #               'trianglelist is wrong!')
+        self.failUnless(data['generatedsegmentlist'] ==seglist,
+                        'segmentlist is wrong!')
+        self.failUnless(data['generatedpointlist'] ==points,
+                        ' is wrong!')
+        #print "data['lonepointlist']", data['lonepointlist']
+        self.failUnless(data['lonepointlist'] ==[0],
+                        'lonepointlist is wrong!')
+if __name__ == "__main__":
+
+    suite = unittest.makeSuite(triangTestCase,'test')
+    #suite = unittest.makeSuite(triangTestCase,'test_lone_verts4')
+    ###suite = unittest.makeSuite(triangTestCase,'testrectangleII')
+    runner = unittest.TextTestRunner() #verbosity=2)
+    runner.run(suite)
Index: uga_core/source/anuga/mesh_engine/test_triangmodule.py
===================================================================
--- /anuga_core/source/anuga/mesh_engine/test_triangmodule.py	(revision 4138)
+++ 	(revision )
@@ -1,371 +1,0 @@
-#!/usr/bin/env python
-#
-
-import sys
-
-import unittest
-import anuga.mesh_engine.triang as triang
-
-
-class triangTestCase(unittest.TestCase):
-    def setUp(self):
-        pass
-
-    def tearDown(self):
-        pass
-
-    def testrectangle(self):
-
-        points = []
-        seglist = []
-        segattlist = []
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
-	# This will cause a seg fault
-        #points = [[0.0,0.0],[0.0,10.0],[3.0,0.0],[3.0,10.0]]
-
-        # automatically generate pointattlist
-	pointattlist = []
-        for point in points:
-            pointattlist.append([])
-        #print "pointattlist",pointattlist
-	#pointattlist = [[],[],[],[]]
-        trilist = []
-        mode = "Qzcn"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist,trilist, mode)
-        self.failUnless(data['generatedtrianglelist'] ==[(1, 0, 2), (2, 3, 1)],
-                        'generatedtrianglelist is wrong!')
-        #print "data['generatedsegmentlist']",data['generatedsegmentlist']
-        self.failUnless(data['generatedsegmentlist'] ==[(2, 0), (3, 2), (1, 3),
-                                                        (0, 1)]
-              ,          'generatedsegmentlist is wrong!')
-        self.failUnless(data['generatedpointlist'] ==[(0.0, 0.0), (0.0, 10.0),
-                                                      (3.0, 0.0), (3.0, 10.0)],
-                        ' is wrong!')
-
-        #print "data['generatedtriangleneighborlist']",
-        #data['generatedtriangleneighborlist']
-        self.failUnless(data['generatedtriangleneighborlist'] ==[(-1, 1, -1),
-                                                                 (-1, 0, -1)],
-                        'generatedtriangleneighborlist is wrong!')
-
-
-    def testrectangleII(self):
-
-        points = []
-        seglist = []
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
-	pointattlist = [[],[],[],[]]
-	regionlist.append( (1.2,1.2,5.0) )
-	seglist = [(0,1),(1,3),(3,2),(2,0)]
-        segattlist = [0,0,0,0]
-        trilist = []
-        mode = "Qzp"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist,trilist, mode)
-
-        self.failUnless(data['generatedtrianglelist'] ==[(1, 0, 2), (2, 3, 1)],
-                        'trianglelist is wrong!')
-        self.failUnless(data['generatedsegmentlist'] ==[(0, 1), (1, 3),
-                                                        (3, 2), (2, 0)],
-                        'segmentlist is wrong!')
-        self.failUnless(data['generatedpointlist'] ==[(0.0, 0.0), (0.0, 10.0),
-                                                      (3.0, 0.0), (3.0, 10.0)],
-                        ' is wrong!')
-
-    def testsegmarker(self):
-
-        points = []
-        seglist = []
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
-	pointattlist = [[],[],[],[]]
-	regionlist.append( (1.2,1.2,5.0) )
-	seglist = [(0,1),(1,3),(3,2),(2,0)]
-        segattlist = [1.0,2.0,3.0,4.0]
-        trilist = []
-        mode = "Qzp"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist,trilist, mode)
-
-        self.failUnless(data['generatedtrianglelist'] ==[(1, 0, 2), (2, 3, 1)],
-                        'trianglelist is wrong!')
-        self.failUnless(data['generatedsegmentlist'] ==[(0, 1), (1, 3),
-                                                        (3, 2), (2, 0)],
-                        'segmentlist is wrong!')
-        self.failUnless(data['generatedpointlist'] ==[(0.0, 0.0), (0.0, 10.0),
-                                                      (3.0, 0.0), (3.0, 10.0)],
-                        ' is wrong!')
-        self.failUnless(data['generatedsegmentmarkerlist'] ==[1,2,3,4],
-                        ' is wrong!')
-
-    def testbad_region(self):
-
-        points = []
-        seglist = []
-        holelist = []
-        regionlist = [(1.2,1.2)]
-
-	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
-	pointattlist = [[],[],[],[]]
-	seglist = [(0,1),(1,3),(3,2),(2,0)]
-        segattlist = [0,0,0,0]
-        trilist = []
-        mode = "Qzpn"
-        try:
-            data = triang.genMesh(points,seglist,holelist,
-                                  regionlist,pointattlist,segattlist,
-                                  trilist, mode)
-
-        except TypeError:
-            pass
-        else:
-            self.failUnless(0 ==1,
-                        'bad region list did not raise error!')
-
-    def testregion_with_maxarea(self):
-
-        points = []
-        seglist = []
-        holelist = []
-        regionlist = [(3,1,1.0)]
-
-	points = [(0.0,0.0),(6.0,0.0),(6.0,6.0),(0.0,6.0)]
-	pointattlist = [[],[],[],[]]
-	seglist = [(0,1),(1,2),(3,2),(3,0),(0,2)]
-        segattlist = [0,0,0,0,0]
-        trilist = []
-        mode = "Qzpna36a"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist, trilist, mode)
-
-        self.failUnless(len(data['generatedtrianglelist']) == 2,
-                        'testregion_with_maxarea 1: # of tris is wrong!')
-        ## Another test case
-        regionlist = [(3,1,1.0),(1,3,1.0,8.0)]
-        mode = "Qzp21na36a"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist, trilist, mode)
-        #print "len(data['generatedtrianglelist']",len(data['generatedtrianglelist'])
-        # for Duncan On unix this returns a 7 triangle result.
-        # for Duncan on Windows returns a 6 triangle result.
-        # for Ole on nautilus this returns 6
-        # for Duncan on nautilus this returns 7
-        # ??, it seems to be the results from triangle that is
-        # causing the different results, and we are treating
-        # triangle as a back box.
-
-        self.failUnless(len(data['generatedtrianglelist']) >= 6,
-                        'testregion_with_maxarea 2: # of tris is wrong!')
-        ## Another test case
-        regionlist = [(3,1,1.0,8.0),(1,3,1.0,8.0)]
-        mode = "Qzpna36a"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist, trilist, mode)
-        #print "len(data['generatedtrianglelist']",len(data['generatedtrianglelist'])
-        # On unix this returns a 10 triangle result.
-        # Windows returns a 8 triangle result.
-        self.failUnless(len(data['generatedtrianglelist']) >= 8,
-                        'testregion_with_maxarea 3: # of tris is wrong!')
-
-        ## Another test case
-        regionlist = [(3,1,1.0),(1,3,1.0,8.0)]
-        mode = "Qzpna8a"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist, trilist, mode)
-        #print "len(data['generatedtrianglelist']",len(data['generatedtrianglelist'])
-        # On unix this returns a 10 triangle result.
-        # Windows returns a 8 triangle result.
-
-        self.failUnless(len(data['generatedtrianglelist']) >= 8,
-                        'testregion_with_maxarea 4: # of tris is wrong!')
-
-    def testbad_point(self):
-
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0)]
-	pointattlist = [[],[],[],[]]
-	seglist = [(0,1),(1,3),(3,2),(2,0)]
-        segattlist = [0,0,0,0]
-        trilist = []
-        mode = "Qzpn"
-        try:
-            data = triang.genMesh(points,seglist,holelist,regionlist,
-                                  pointattlist,segattlist, trilist, mode)
-
-        except TypeError:
-            pass
-        else:
-            self.failUnless(0 ==1,
-                        'bad point list did not raise error!')
-
-    def testbad_hole(self):
-        
-        holelist = [(9.0)]
-        regionlist = []
-        
-        points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
-        pointattlist = [[],[],[],[]]
-        seglist = [(0,1),(1,3),(3,2),(2,0)]
-        segattlist = [0,0,0,0]
-        trilist = []
-        mode = "Qzpn"
-        try:
-            data = triang.genMesh(points,seglist,holelist,regionlist,
-                                  pointattlist,segattlist, trilist, mode)
-
-        except TypeError:
-            pass
-        else:
-            self.failUnless(0 ==1,
-                        'bad hole list did not raise error!')
-
-    def testbad_segment(self):
-
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
-	pointattlist = [[],[],[],[]]
-	seglist = [(0,1),(1,3),(3,2),(2)]
-        segattlist = [0,0,0,0]
-        trilist = []
-        mode = "Qzpn"
-        try:
-            data = triang.genMesh(points,seglist,holelist,regionlist,
-                                  pointattlist,segattlist, trilist, mode)
-
-        except TypeError:
-            pass
-        else:
-            self.failUnless(0 ==1,
-                        'bad segment list did not raise error!')
-
-    def testbad_segattlist(self):
-
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
-	pointattlist = [[],[],[],[]]
-	seglist = [(0,1),(1,3),(3,2),(2,0)]
-        segattlist = [0,0]
-        trilist = []
-        mode = "Qzpn"
-        try:
-            data = triang.genMesh(points,seglist,holelist,regionlist,
-                                  pointattlist,segattlist, trilist, mode)
-
-            self.failUnless(data['trianglelist'] ==[(1, 0, 2), (2, 3, 1)],
-                        'trianglelist is wrong!')
-        except TypeError:
-            pass
-        else:
-            self.failUnless(0 ==1,
-                        'bad segment attribute list did not raise error!')
-
-    def testrectangle_regions(self):
-
-        points = []
-        seglist = []
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
-	pointattlist = [[],[],[],[]]
-        #it seems that
-        #triangle only associates one region with a triangle
-	regionlist.append( (1.3,1.3,88.33) )
-	regionlist.append( (1.2,1.2,77,55) )
-	seglist = [(0,1),(1,3),(3,2),(2,0)]
-        segattlist = [0,0,0,0]
-        trilist = []
-        mode = "QAzpq"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist, trilist, mode)
-
-        self.failUnless(data['generatedtriangleattributelist'] ==[[77.0], [77.0], [77.0], [77.0]],
-                        'triangleattributelist is wrong!')
-
-
-    def BAD_SEG_FAULT_test_lone_verts(self):
-
-        points = []
-        seglist = []
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0),(0.0,10.0)]
-	pointattlist = [[],[],[],[],[]]
-	regionlist.append( (1.2,1.2,5.0) )
-	seglist = [(0,4),(4,1),(1,3),(3,2),(2,0)]
-        segattlist = [0,0,0,0,0]
-        trilist = []
-        mode = "QpznAa2000.1a"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist,trilist, mode)
-        #print "data['generatedtrianglelist']", data['generatedtrianglelist']
-        self.failUnless(data['generatedtrianglelist'] ==[(4, 0, 2), (2, 3, 4)],
-                        'trianglelist is wrong!')
-        self.failUnless(data['generatedsegmentlist'] ==[(0, 4), (4, 3),
-                                                        (3, 2), (2, 0)],
-                        'segmentlist is wrong!')
-        self.failUnless(data['generatedpointlist'] ==[(0.0, 0.0), (0.0, 10.0),
-                                                      (3.0, 0.0), (3.0, 10.0),
-                                                      (0.0,10.0)],
-                        ' is wrong!')
-        self.failUnless(data['lonepointlist'] ==[1],
-                        'lonepointlist is wrong!')
-        #print "", data['lonepointlist']
-
-    def BAD_SEG_FAULT_test_lone_vertsII(self):
-
-        points = []
-        seglist = []
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,0.0),(0.0,10.0),(0.0,10.0),(10.0,10.0),
-                  (10.0,10.0),(0.0,10.0),(10.0,0.0)]
-
-	pointattlist = []
-        for point in points:
-            pointattlist.append([])
-	seglist = [(0,1),(1,2),(2,3),(3,4),(4,5),(5,7),(7,0)]
-	segattlist = []
-        for seg in seglist:
-            segattlist.append(0)
-        trilist = []
-        mode = "pznAa2000.1a"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist,trilist, mode)
-        #print "data['generatedtrianglelist']", data['generatedtrianglelist']
-        self.failUnless(data['generatedtrianglelist'] ==[(6, 1, 7), (7, 5, 6)],
-                        'trianglelist is wrong!')
-        self.failUnless(data['generatedsegmentlist'] ==[(1, 6), (6, 5),
-                                                        (5, 7), (7, 1)],
-                        'segmentlist is wrong!')
-        self.failUnless(data['generatedpointlist'] ==[(0.0,0.0),(0.0,0.0),
-                                                      (0.0,10.0),(0.0,10.0),
-                                                      (10.0,10.0),
-                                                      (10.0,10.0),(0.0,10.0),
-                                                      (10.0,0.0)],
-                        ' is wrong!')
-        self.failUnless(data['lonepointlist'] ==[0,2,3,4],
-                        'lonepointlist is wrong!')
-
-if __name__ == "__main__":
-
-    suite = unittest.makeSuite(triangTestCase,'test')
-    runner = unittest.TextTestRunner()
-    runner.run(suite)
Index: uga_core/source/anuga/mesh_engine/test_triangmoduleII.py
===================================================================
--- /anuga_core/source/anuga/mesh_engine/test_triangmoduleII.py	(revision 4138)
+++ 	(revision )
@@ -1,85 +1,0 @@
-#!/usr/bin/env python
-#
-
-import sys
-
-import unittest
-import anuga.mesh_engine.triang as triang
-
-"""
-Note, this test is seperated from the other tests because if it is added
-the triangle code seg faults.  Don't know why
-"""
-class triangTestCase(unittest.TestCase):
-
-    def setUp(self):
-        pass
-
-    def tearDown(self):
-        pass
-
-    def BAD_test_lone_vertsII(self):
-
-        points = []
-        seglist = []
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,0.0),(0.0,10.0),(0.0,10.0),(10.0,10.0),
-                  (10.0,10.0),(0.0,10.0),(10.0,0.0)]
-
-	pointattlist = []
-        for point in points:
-            pointattlist.append([])
-	seglist = [(0,1),(1,2),(2,3),(3,4),(4,5),(5,7),(7,0)]
-	segattlist = []
-        for seg in seglist:
-            segattlist.append(0)
-        trilist = []
-        mode = "QpznAa2000.1a"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist,trilist, mode)
-        #print "data['generatedtrianglelist']", data['generatedtrianglelist']
-        self.failUnless(data['generatedtrianglelist'] ==[(6, 1, 7), (7, 5, 6)],
-                        'trianglelist is wrong!')
-        self.failUnless(data['generatedsegmentlist'] ==[(1, 6), (6, 5),
-                                                        (5, 7), (7, 1)],
-                        'segmentlist is wrong!')
-        self.failUnless(data['generatedpointlist'] ==[(0.0,0.0),(0.0,0.0),
-                                                      (0.0,10.0),(0.0,10.0),
-                                                      (10.0,10.0),
-                                                      (10.0,10.0),(0.0,10.0),
-                                                      (10.0,0.0)],
-                        ' is wrong!')
-        self.failUnless(data['lonepointlist'] ==[0,2,3,4],
-                        'lonepointlist is wrong!')
-
-    def test_lone_vertsIII(self):
-
-        points = []
-        seglist = []
-        holelist = []
-        regionlist = []
-
-	points = [(0.0,0.0),(0.0,0.0),(0.0,10.0),(0.0,10.0),(10.0,10.0),
-                  (10.0,10.0),(0.0,10.0),(10.0,0.0)]
-
-	pointattlist = []
-        for point in points:
-            pointattlist.append([])
-	seglist = []
-	segattlist = []
-        for seg in seglist:
-            segattlist.append(0)
-        trilist = []
-        mode = "Qzp"
-        data = triang.genMesh(points,seglist,holelist,regionlist,
-                              pointattlist,segattlist,trilist, mode)
-
-        self.failUnless(data['lonepointlist'] ==[0,1,2,3,4,5,6,7],
-                        'lonepointlist is wrong!')
-if __name__ == "__main__":
-
-    suite = unittest.makeSuite(triangTestCase,'test')
-    runner = unittest.TextTestRunner()
-    runner.run(suite)
Index: uga_core/source/anuga/mesh_engine/triang.c
===================================================================
--- /anuga_core/source/anuga/mesh_engine/triang.c	(revision 4138)
+++ 	(revision )
@@ -1,618 +1,0 @@
-/* this appears to be necessary to run on Linux */
-
-#undef _STDLIB_H_
-#define _STDLIB_H_
-
-
-/*
-    This code is based on code by A. Pletzer
-
-    This code interfaces pmesh directly with "triangle", a general        
-    purpose triangulation code. In doing so, Python data structures          
-    are passed to C for  use by "triangle" and the results are then          
-    converted back. This was accomplished using the Python/C API.            
-                                                                           
-    I. Lists of points,edges,holes, and regions must normally be 
-     created for proper triangulation. However, if there are no holes,
-     the  hole list need not be specified (however, a 
-     measure of control over the shape of the convex hull is lost; 
-     triangle decides). 
-     
-     points list: [(x1,y1),(x2,y2),...] (Tuples of doubles)                  
-     edge list: [(Ea,Eb),(Ef,Eg),...] (Tuples of integers)                   
-     hole list: [(x1,y1),...](Tuples of doubles, one inside each hole region)
-     regionlist: [ (x1,y1,index,area),...] (Tuple of doubles)                
-     pointattributelist:[ (a1,a2..),...] (Tuple of lists) 
-     mode: String - the commands sent to triangle
-          
-     This is all that is needed to accomplish an intial triangulation.       
-     A dictionary of the Triangle output is returned.                                   
-     Precondition
-     End list in the pointattributelist has to have the same length
-     
-    Compilation can be done as follows                                       
-                    
-    c++ -I/usr/local/include/python1.5 -shared -o triangmodule.so            
-    triangmodule.c triangle.o                                                */
-
-
-#ifdef SINGLE  
-#define REAL float 
-#else  
-#define REAL double 
-#endif 
-
-#include "triangle.h" 
-
-#ifndef _STDLIB_H_ 
-#ifndef _WIN32
-extern "C" void *malloc(); 
-extern "C" void free(); 
-#endif
-#endif  
-
-#include "Python.h" 
-
-
-
- static PyObject *triang_genMesh(PyObject *self, PyObject *args){
-     
-  /* Mesh generation routine
-     pre-conditions:
-     mode must have a 'z' switch, to avoid segmentation faults
-  */
-    
-  struct triangulateio in,out;
-
-  PyObject *hull;
-  PyObject *seglist;
-  PyObject *holelist;
-  PyObject *regionlist;
-  PyObject *trianglelist;
-  PyObject *mode;
-  PyObject *firstPointAttr;
-  PyObject *pointattributelist;
-  PyObject *segmarkerlist;
-  PyObject *Attrlist;
-  REAL Attr;
-  int i, j, iatt, n, write_here;
-  int a,b,c;
-  int marker;
-  int tuplesize;
-    
-  int index = 0;
-  double x,y;
-  PyObject *elems;
-  PyObject *pair;
-     
-  char *mod;
-  PyObject *holder, *holderlist, *attributelist;
-  int listsize,attsize ;
-  PyObject  *ii;
-  
-  int *points_connected;
-       
-  if(!PyArg_ParseTuple(args,(char *)"OOOOOOOO",&hull,&seglist,&holelist,&regionlist,&pointattributelist,&segmarkerlist,&trianglelist,&mode)){
-    return NULL;
-  }
-  if(!PyList_Check(hull)){
-    PyErr_SetString(PyExc_TypeError, 
-		    "incorrect first argument for 'genMesh': list required.");
-    return NULL;
-  }
-  if(!PyList_Check(hull)){
-    PyErr_SetString(PyExc_TypeError,
-		    "incorrect second argument for 'genMesh': list required.");
-    return NULL;
-  }
-  if(!PyList_Check(hull)){
-    PyErr_SetString(PyExc_TypeError,
-		    "incorrect third argument for 'genMesh': list required.");
-    return NULL;
-  }
-  if(!PyList_Check(hull)){
-    PyErr_SetString(PyExc_TypeError,
-		    "incorrect fourth argument for 'genMesh': list required.");
-    return NULL;
-  }
-  if(!PyList_Check(hull)){
-    PyErr_SetString(PyExc_TypeError,
-		    "incorrect fifth argument for 'genMesh': list required.");
-    return NULL;
-  }
-  if(!PyList_Check(hull)){
-    PyErr_SetString(PyExc_TypeError,
-		    "incorrect sixth argument for 'genMesh': list required.");
-    return NULL;
-  }
-  if(!PyList_Check(hull)){
-    PyErr_SetString(PyExc_TypeError,
-		    "incorrect seventh argument for 'genMesh': list required.");
-    return NULL;
-  }
-  if(!PyList_Check(hull)){
-    PyErr_SetString(PyExc_TypeError,
-		    "incorrect eighth argument for 'genMesh': list required.");
-    return NULL;
-  }
-     
-  /* printf("starting\n"); */    
-     
-  /* Initialize  points */
-  in.numberofpoints =  PyList_Size(hull); 
-  in.pointlist = (REAL *)malloc(in.numberofpoints*2*sizeof(REAL));
-  in.pointmarkerlist = (int *)NULL; /* malloc(in.numberofpoints*sizeof(int)); */
-     
-  /* Initialize and fill up region list */
-  in.numberofregions = PyList_Size(regionlist);
-  
-  /*printf(" numberofregions %i\n",in.numberofregions );*/
-  in.regionlist = (REAL *)malloc(in.numberofregions*4*sizeof(REAL));
-  for(i=0; i<in.numberofregions;i++){
-    elems = PyList_GetItem(regionlist,i);
-    tuplesize = PySequence_Length(elems);
-    /*printf("tuplesize %i\n",tuplesize);*/
-    if (tuplesize>0){
-      in.regionlist[4*i] = PyFloat_AsDouble(PyTuple_GetItem(elems,0));
-    } else {
-      PyErr_SetString(PyExc_TypeError,
-		      "data sent to trianglemodule error: A region has no x value.");
-      return NULL;
-    }
-    if (tuplesize>1){
-      in.regionlist[4*i+1] = PyFloat_AsDouble(PyTuple_GetItem(elems,1));
-    } else {
-      PyErr_SetString(PyExc_TypeError,
-		    "data sent to trianglemodule error: A region has no y value.");
-      return NULL;
-    }
-    if (tuplesize>2){
-      /* 2 is attribute*/
-      in.regionlist[4*i+2] = PyFloat_AsDouble(PyTuple_GetItem(elems,2));
-    } else {
-      PyErr_SetString(PyExc_TypeError,
-		    "data sent to trianglemodule error: A region has no attribute value.");
-      return NULL;
-    }
-    if (tuplesize>3){
-      in.regionlist[4*i+3] = PyFloat_AsDouble(PyTuple_GetItem(elems,3));
-    } else {
-      /* 0.0 is interpreted as no local max area by triangle  */
-      in.regionlist[4*i+3] = (REAL)0.0;
-    } 
-  }
-  
-  /*Initialize segments */
-  in.numberofsegments = PyList_Size(seglist);
-  in.segmentlist = (int *)malloc(in.numberofsegments*2*sizeof(int));
-  in.segmentmarkerlist = (int *)malloc(in.numberofsegments*sizeof(int));
-
-    
-  /*Initialize triangles */
-  in.numberofcorners = 3;
-  in.numberoftriangles = PyList_Size(trianglelist);
-  if (in.numberoftriangles  != 0) {
-    in.trianglelist = (int *)malloc(in.numberofsegments*3*sizeof(int));
-  } else {
-    in.trianglelist = (int *)NULL;
-  }
-  in.numberoftriangleattributes = 0;  
-  in.triangleattributelist  = (REAL *)NULL; 
-     
-     
-     
-  /*Initialize holes */
-  in.numberofholes = PyList_Size(holelist);
-  if (in.numberofholes != 0) {
-    in.holelist = (REAL *)malloc(in.numberofholes*2*sizeof(REAL));
-  } else {
-    in.holelist = (REAL *)NULL;
-  }     
-     
-  /* Fill up triangle list */
-  index = 0;
-  for(i = 0; i < in.numberoftriangles; i++) {
-    pair =  PyList_GetItem(trianglelist,i);
-    if (3 != PySequence_Length(pair)) {
-      PyErr_SetString(PyExc_TypeError,
-		      "data sent to trianglemodule error:Found a triangle without 3 vertices.");
-      return NULL;
-    }
-    a = (int)PyInt_AsLong(PyTuple_GetItem(pair,0));
-    in.trianglelist[index] = a;
-    index++;
-    b = (int)PyInt_AsLong(PyTuple_GetItem(pair,1));
-    in.trianglelist[index] = b;
-    index++;
-    c = (int)PyInt_AsLong(PyTuple_GetItem(pair,2));
-    in.trianglelist[index] = c;
-    index++;
-  }
-     
-     
-  /*  print the triangle list DEBUGGING
-      printf ("numberoftriangles %i\n", in.numberoftriangles);
-      for(i = 0; i < in.numberoftriangles; i++) {
-      printf("(%i,%i,%i)",in.trianglelist[i* 3+ 0],in.trianglelist[i* 3+ 1], in.trianglelist[i* 3+ 2]);
-      }*/
-     
-  /* Fill up segment list */
-  index = 0;
-  if (in.numberofsegments != PyList_Size(segmarkerlist)) {
-    PyErr_SetString(PyExc_TypeError,
-		    "data sent to trianglemodule error:number of segments doesn't = number of segment attributes.");
-    return NULL;
-  }
-    
-  for(i = 0; i < in.numberofsegments; i++) {
-    pair =  PyList_GetItem(seglist,i);
-    if (2 != PySequence_Length(pair)) {
-      PyErr_SetString(PyExc_TypeError,
-		      "data sent to trianglemodule error:Found a segment without 2 vertices.");
-      return NULL;
-    }
-    a = (int)PyInt_AsLong(PyTuple_GetItem(pair,0));
-    in.segmentlist[index] = a;
-    index++;
-    b = (int)PyInt_AsLong(PyTuple_GetItem(pair,1));
-    in.segmentlist[index] = b;
-    index++;
-    /* Fill up segment marker list */
-    marker =  (int)PyInt_AsLong(PyList_GetItem(segmarkerlist,i));
-    in.segmentmarkerlist[i] = marker;
-  }
-  /* Fill up hole list */
-  index = 0;
-  for(i = 0; i < in.numberofholes; i++) {
-    pair =  PyList_GetItem(holelist,i);
-    if (2 != PySequence_Length(pair)) {
-      PyErr_SetString(PyExc_TypeError,
-		      "data sent to trianglemodule error:Found a pair without 2 vertices.");
-      return NULL;
-    }
-    x = PyFloat_AsDouble(PyTuple_GetItem(pair,0));
-    in.holelist[index] = x;
-    index++;
-    y = PyFloat_AsDouble(PyTuple_GetItem(pair,1));
-    in.holelist[index] = y;
-    index++;
-  }
-     
-  /* Fill up point list */
-  index = 0;
-  for(i = 0; i < in.numberofpoints; i++) {
-    pair =  PyList_GetItem(hull,i);
-    if (2 != PySequence_Length(pair)) {
-      PyErr_SetString(PyExc_TypeError,
-		      "data sent to trianglemodule error:Found a vertex without an x,y values.");
-      return NULL;
-    }
-    x = PyFloat_AsDouble(PyTuple_GetItem(pair,0));
-    in.pointlist[index] = x;
-    index++;
-    y = PyFloat_AsDouble(PyTuple_GetItem(pair,1));
-    in.pointlist[index] = y;
-    index++;
-  }
-
-  /* Initialize  point attributes */
-/*      firstPointAttr = PyList_GetItem(pointattributelist,0); */
-/*      if (firstPointAttr != NULL) { */
-/*        //in.numberofpointattributes = PyList_Size(firstPointAttr); */
-/*        printf("in.numberofpointattributes =%i",PyList_Size(firstPointAttr)); */
-/*        //printf("in.numberofpointattributes =%i",in.numberofpointattributes); */
-/*      } else { */
-/*        printf("(firstPointAttr is NULL\n"); */
-/*        in.pointattributelist = NULL; */
-/*        in.numberofpointattributes = 0; */
-/*      } */
-  n = PyList_Size(pointattributelist);
-  if (n < 0) {
-    in.numberofpointattributes = 0; 
-    in.pointattributelist = (REAL *)NULL;
-  } else {
-    /* Get the # of attributes based on the size of the first point*/
-    firstPointAttr = PyList_GetItem(pointattributelist,0);
-    if (firstPointAttr == NULL) {
-      in.numberofpointattributes = 0; 
-      in.pointattributelist = (REAL *)NULL;
-    } else {
-      in.numberofpointattributes = PyList_Size(firstPointAttr);
-      in.pointattributelist = (REAL *)malloc(in.numberofpoints
-					     *in.numberofpointattributes
-					     *sizeof(REAL)); 
-     
-     
-      /* fill the point attribute list */
-      if (in.numberofpointattributes != 0) {
-	for(j = 0; j < in.numberofpoints; j++){
-	  Attrlist = PyList_GetItem(pointattributelist,j);
-	  if (in.numberofpointattributes != PySequence_Length(Attrlist)) {
-	    PyErr_SetString(PyExc_TypeError,
-			    "data sent to trianglemodule error:List size of attributes is different from vertex to vertex.");
-	    return NULL;
-	  }
-	  for(i = 0; i < in.numberofpointattributes; i++){
-	    Attr = PyFloat_AsDouble(PyList_GetItem(Attrlist,i));
-	    in.pointattributelist[in.numberofpointattributes*j+i]= (REAL)Attr;
-	    /* 	 printf("i =%i\n",i);  */
-	    /* 	 printf("j =%i\n",j);  */
-	    /* 	 printf("in.pointattributelist[4*j+i]
-		 =%f\n",in.pointattributelist[4*j+i]); */
-	  }
-	}
-      } 
-    }
-  }
-    
-     
-  /* this half from triangulate*/
-  /* set up the switch arguments to the triangulation routine */
-
-  mod = PyString_AS_STRING(mode);
-     
-         
-  out.pointlist = (REAL *)NULL;
-  out.pointmarkerlist = (int *)NULL;
-  out.numberofpointattributes = in.numberofpointattributes;
-  out.pointattributelist = (REAL *)NULL;
-  
-  out.trianglelist = (int *)NULL;
-  out.triangleattributelist = (REAL *)NULL;                   
-  out.trianglearealist = (REAL *)NULL;                        
-  out.neighborlist = (int *)NULL;
-     
-  out.segmentlist = (int *)NULL;
-  out.segmentmarkerlist = (int *)NULL;
-     
-  out.edgelist = (int *)NULL;
-  out.edgemarkerlist = (int *)NULL;
-     
-  out.holelist = (REAL *)NULL;
-  out.regionlist = (REAL *)NULL;
-    
-  
-  /*printf("\n\nTriangulate input args: %s \n\n", mod); */
-  triangulate(mod, &in, &out, (struct triangulateio *)NULL );
-
-  /*
-    ------- Pass point numbers,coordinates and neighbors back to Python ------
-    we send back a dictionary:                                               
-    { index : [ coordinates, [connections], Attribute ] } 
-  */
-  holder = PyDict_New();
-     
-  /* list of int's, used to keep track of which verts are connected to
-     triangles. */
-  points_connected = (int *)malloc(out.numberofpoints*sizeof(int));
-  
-  /* Add pointlist */
-  listsize = out.numberofpoints;
-  holderlist = PyList_New(listsize);
-     
-  for(i=0; i<listsize;i++){
-    PyObject *mlist = Py_BuildValue((char *)"(d,d)", 
-				    out.pointlist[i*2  ], out.pointlist[i*2+1]);  
-    PyList_SetItem(holderlist,i, mlist);
-    points_connected[i] = 0;
-  }  
-  ii=PyString_FromString("generatedpointlist");
-  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);
-    
-  /* Add point marker list */
-  listsize = out.numberofpoints;
-  holderlist = PyList_New(listsize);
-  
-  for(i=0; i<listsize;i++){
-    PyObject *mlist = Py_BuildValue((char *)"d", 
-				    out.pointmarkerlist[i]);     
-    PyList_SetItem(holderlist,i, mlist); 
-  }  
-  ii=PyString_FromString("generatedpointmarkerlist");
-  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);
-    
-  /* Add point attribute list */
-  listsize = out.numberofpoints;
-  holderlist = PyList_New(listsize);
-  
-  for(i=0; i<listsize;i++){
-    attsize = out.numberofpointattributes;
-    attributelist = PyList_New(attsize);
-    for(iatt=0; iatt<attsize;iatt++){
-      PyObject *mlist = Py_BuildValue((char *)"d", 
-				      out.pointattributelist[i*attsize + iatt]);     
-      PyList_SetItem(attributelist,iatt, mlist);
-    }      
-    PyList_SetItem(holderlist,i, attributelist);
-  }  
-  ii=PyString_FromString("generatedpointattributelist");
-  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);  
- 
-  
-  /* Add triangle list */
-  listsize = out.numberoftriangles;
-  holderlist = PyList_New(listsize);
-  for(i=0; i<listsize;i++){
-    PyObject *mlist = Py_BuildValue((char *)"(i,i,i)", 
-				    out.trianglelist[i*3  ], out.trianglelist [i*3+1], out.trianglelist [i*3+2]);    
-    PyList_SetItem(holderlist,i, mlist);
-    /* printf(" A vert index %i\n",out.trianglelist[i*3] );
-    printf(" A vert index %i\n",out.trianglelist[i*3+1] );
-    printf(" A vert index %i\n",out.trianglelist[i*3+2] ); */
-    points_connected[out.trianglelist[i*3]] = 1;
-    points_connected[out.trianglelist[i*3+1]] = 1;
-    points_connected[out.trianglelist[i*3+2]] = 1;
-  }    
-  ii=PyString_FromString("generatedtrianglelist");
-  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);
-  
-  /* convert the points_connected vector from a true(1) false(0) vector, where
-     index is the vert, to a vector of the lone verts, at the beggining
-     of the vector. */
-  write_here = 0;   
-  for(i=0; i<out.numberofpoints;i++){
-    if (0 == points_connected[i]) {
-      points_connected[write_here] = i;
-      write_here ++;
-    }
-  }   
-  /* printf(" ******************** \n" );
-  for(i=0; i<write_here;i++){
-    printf(" A vert index %i\n",points_connected[i] );
-    } */
-  
-  listsize = write_here;
-  holderlist = PyList_New(listsize);
-  for(i=0; i<listsize;i++){
-   PyObject *mlist = Py_BuildValue((char *)"i", points_connected[i]);    
-    PyList_SetItem(holderlist,i, mlist);
-  }
-  ii=PyString_FromString("lonepointlist");
-  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);  
-     
-  /* Add triangle attribute list */
-  listsize = out.numberoftriangles;
-  holderlist = PyList_New(listsize);
-     
-  for(i=0; i<listsize; i++){
-    attsize = out.numberoftriangleattributes;
-    attributelist = PyList_New(attsize);       
-    for(iatt=0; iatt<attsize;iatt++){
-      PyObject *mlist = Py_BuildValue((char *)"d",out.triangleattributelist[i*attsize + iatt]);     
-      PyList_SetItem(attributelist,iatt, mlist);
-    }      
-    PyList_SetItem(holderlist,i, attributelist);
-  }  
-  ii=PyString_FromString("generatedtriangleattributelist");
-  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);
-  
-  /* Add segment list */
-  listsize = out.numberofsegments;
-  holderlist = PyList_New(listsize);
-  for(i=0; i<listsize;i++){
-    PyObject *mlist = Py_BuildValue((char *)"(i,i)", 
-				    out.segmentlist[i*2  ],
-				    out.segmentlist [i*2+1]);
-    PyList_SetItem(holderlist,i, mlist);
-  }    
-  ii=PyString_FromString("generatedsegmentlist");
-  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);  
-  
-  /* Add segment marker list */
-  listsize = out.numberofsegments;
-  holderlist = PyList_New(listsize);
-  for(i=0; i<listsize;i++){
-    PyObject *mlist = Py_BuildValue((char *)"i", 
-				    out.segmentmarkerlist[i]);
-    PyList_SetItem(holderlist,i, mlist);
-  }    
-  ii=PyString_FromString("generatedsegmentmarkerlist");
-  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);  
-  /* Add triangle neighbor list */
-  if (out.neighborlist != NULL) {
-    listsize = out.numberoftriangles;
-    holderlist = PyList_New(listsize);
-    for(i=0; i<listsize;i++){
-      PyObject *mlist = Py_BuildValue((char *)"(i,i,i)", 
-				      out.neighborlist[i*3  ],
-				      out.neighborlist [i*3+1],
-				      out.neighborlist [i*3+2]);    
-      PyList_SetItem(holderlist,i, mlist);
-    }    
-    ii=PyString_FromString("generatedtriangleneighborlist");
-    PyDict_SetItem(holder, ii, holderlist);
-    Py_DECREF(ii); Py_DECREF(holderlist);
-  }   
-     
-  /* Free in/out structure memory */
-  
-  /* OUT */
-
-  if(!out.pointlist){
-    free(out.pointlist);  out.pointlist=NULL;
-  }
-  if(!out.pointmarkerlist){    
-    free(out.pointmarkerlist); out.pointmarkerlist=NULL;
-  }
-  if(!out.pointattributelist){    
-    free(out.pointattributelist); out.pointattributelist=NULL;
-  }
-    
- 
-  if(!out.trianglelist){    
-    free(out.trianglelist); out.trianglelist=NULL;
-  }
-  if(!out.triangleattributelist){    
-    free(out.triangleattributelist); out.triangleattributelist=NULL;
-  }
-  if(!out.trianglearealist){    
-    free(out.trianglearealist); out.trianglearealist=NULL;
-  }
-  if(!out.neighborlist){    
-    free(out.neighborlist); out.neighborlist=NULL;
-  }
-  if(!out.segmentlist){
-    free(out.segmentlist); out.segmentlist =NULL;
-  }
-  if(!out.segmentmarkerlist){
-    free(out.segmentmarkerlist); out.segmentmarkerlist  =NULL;
-  }
-  if(!out.edgelist){
-    free(out.edgelist);  out.edgelist=NULL;
-  }
-  if(!out.edgemarkerlist){
-    free(out.edgemarkerlist);  out.edgemarkerlist=NULL;
-  }
-  if(!out.holelist){
-    free(out.holelist); out.holelist=NULL;
-  }
-  if(!out.regionlist){
-    free(out.regionlist); out.regionlist=NULL;
-  }
-  if(!points_connected ){
-    free(points_connected ); points_connected =NULL;
-  }
-    
-  
-  
-  /*  IN */
-  
-  if(!in.pointlist){
-    free(in.pointlist);  in.pointlist   =NULL;             
-  }
-  if(!in.pointattributelist){
-    free(in.pointattributelist); in.pointattributelist =NULL;       
-  }
-  if(!in.pointmarkerlist){
-    free(in.pointmarkerlist); in.pointmarkerlist    =NULL;        
-  }
-  if(!in.segmentlist){
-    free(in.segmentlist);  in.segmentlist    =NULL;         
-  }
-  if(!in.segmentmarkerlist){
-    free(in.segmentmarkerlist); in.segmentmarkerlist    =NULL;   
-  }
-  if(!in.regionlist){
-    free(in.regionlist); in.regionlist  =NULL;            
-  }
-  if(!in.holelist){
-    free(in.holelist); in.holelist=NULL;
-  }
-
-  //return Py_BuildValue((char *)"O", holder);
-  return Py_BuildValue((char *)"O", holder);
-}
-
-/* end of my function*/
-
-
-   
-    
-static PyMethodDef triang_methods[] = {
-  {(char *)"genMesh", triang_genMesh, 1, (char *)"Passes hull points to triangle.All arguments are lists.(<hull>,<seglist>,<holelist>,<regionlist>)-><hull>"}, 
-  {NULL,NULL}
-};
-
-void inittriang(){
-  Py_InitModule((char *)"triang",triang_methods);
-}    
Index: /anuga_core/source/anuga/pmesh/mesh.py
===================================================================
--- /anuga_core/source/anuga/pmesh/mesh.py	(revision 4138)
+++ /anuga_core/source/anuga/pmesh/mesh.py	(revision 4139)
@@ -5,5 +5,5 @@
    Note: A .index attribute is added to objects such as vertices and
    segments, often when reading and writing to files.  This information
-   should not be used as percistant information.  It is not the 'index' of
+   should not be used as persistant information.  It is not the 'index' of
    an element in a list.
 
@@ -16,5 +16,4 @@
 import sys
 import math
-import anuga.mesh_engine.triang as triang
 import re
 import os
@@ -23,11 +22,16 @@
 import types
 import exceptions
-
-import load_mesh
+from Numeric import array, Float, Int
+
+
+#import load_mesh
 from anuga.coordinate_transforms.geo_reference import Geo_reference,DEFAULT_ZONE
 from anuga.utilities.polygon import point_in_polygon 
 import load_mesh.loadASCII
 import anuga.alpha_shape.alpha_shape
-from anuga.geospatial_data.geospatial_data import Geospatial_data, ensure_geospatial, ensure_absolute
+from anuga.geospatial_data.geospatial_data import Geospatial_data, ensure_geospatial, ensure_absolute, ensure_numeric
+from anuga.mesh_engine.mesh_engine import generate_mesh
+
+#import anuga.mesh_engine_b.mesh_engine as triang
 
 try:  
@@ -467,4 +471,55 @@
 Segment.set_default_tag("")       
 
+class Rigid_triangulation:
+    """ 
+    This is a triangulation that can't have triangles added or taken away.
+
+    It just represents the triangulation, not the mesh outline needed to
+    build the triangulation.
+
+      This is the guts of the data structure;
+        generated vertex list: [(x1,y1),(x2,y2),...] (Tuples of doubles)
+        generated segment list: [(point1,point2),(p3,p4),...]
+            (Tuples of integers) 
+        generated segment tag list: [tag,tag,...] list of strings
+
+        generated triangle list: [(p1,p2,p3), (p4,p5,p6),....] tuple of points
+
+        generated triangle attribute list: [s1,s2,...] list of strings
+
+        generated triangle neighbor list: [(t1,t2,t3), (t4,t5,t6),....]
+            tuple of triangles
+
+            Should I do the basic structure like triangle, general
+            mesh or something else? How about like triangle, since
+            that's where the info is from, and it'll fit easier into
+            this file..
+
+             Removing the loners is difficult, since all the vert's
+             after it must be removed.
+
+             This happens in set_triangulation.
+            
+    """
+
+    def __init__(self,
+                 triangles,
+                 segments,
+                 vertices,
+                 triangle_tags,
+                 triangle_neighbor,
+                 segment_tags,
+                 ):
+       
+        self.triangles=ensure_numeric(triangles) 
+        self.triangle_neighbor=ensure_numeric(triangle_neighbor)
+        self.triangle_tags=triangle_tags
+        self.segments=ensure_numeric(segments) 
+        
+        self.segment_tags=segment_tags # string
+        self.vertices=ensure_numeric(vertices)
+
+        
+        
 class Mesh:
     """
@@ -1054,4 +1109,7 @@
         #print "mesh#generateMesh# self.mode",self.mode  
         meshDict = self.Mesh2triangList()
+
+        #FIXME (DSG-DSG)  move below section into generate_mesh.py
+        #                  & 4 functions eg segment_strings2ints
         #print "*************************!@!@ This is going to triangle   !@!@"
         #print meshDict
@@ -1067,5 +1125,5 @@
         #print "%%%%%%%%%%%%%%%%%%%%%%%%%%%regionlist",meshDict['regionlist']
         #print "meshDict['segmenttaglist']", meshDict['segmenttaglist'
-        generatedMesh = triang.genMesh(
+        generatedMesh = generate_mesh(
                               meshDict['pointlist'],
                               meshDict['segmentlist'],
@@ -1074,6 +1132,6 @@
                               meshDict['pointattributelist'],
                               meshDict['segmenttaglist'],
-                              [],  # since the trianglelist isn't used
-                              self.mode)
+                              self.mode,
+                              meshDict['pointlist'])
         #print "%%%%%%%%%%%%%%%%%%%%%%%%%%%generated",generatedMesh
         generatedMesh['generatedsegmentmarkerlist'] = \
@@ -1081,8 +1139,10 @@
                                   segconverter)
         #print "processed gen",generatedMesh['generatedsegmentmarkerlist']
+        #print "generatedtriangleattributelist",generatedMesh['generatedtriangleattributelist']
         generatedMesh['generatedtriangleattributelist'] = \
          region_ints2strings(generatedMesh['generatedtriangleattributelist'],
                                   regionconverter)
 
+        #FIXME (DSG-DSG)  move above section into generate_mesh.py
 
         if len(generatedMesh['generatedpointattributelist'][0])==0:
@@ -1400,4 +1460,6 @@
             -1 means there's no triangle neighbor
         triangle attribute list: [(T1att), (T2att), ...]
+
+
             (list of a list of strings)
         lone point list:[point1, ...] (list of integers)
Index: /anuga_core/source/anuga/pmesh/test_mesh.py
===================================================================
--- /anuga_core/source/anuga/pmesh/test_mesh.py	(revision 4138)
+++ /anuga_core/source/anuga/pmesh/test_mesh.py	(revision 4139)
@@ -1068,5 +1068,66 @@
                         ,
                         'exported Ascii xya file is wrong')
-           
+     
+    def test_lone_vert_in_mesh_gen_c_layer(self):
+        # currently just a copy of the above test
+        a = Vertex (0,0)
+        b = Vertex (0,3)
+        c = Vertex (3,3)
+        d = Vertex (1,2)
+        e = Vertex (3,1)
+        #f = Vertex (3,1)
+      
+        s1 = Segment(a,b, tag = 50)
+        s2 = Segment(b,c, tag = 40)
+        s3 = Segment(c,a, tag = 30)
+     
+        r1 = Region(2, 1,tag = 1.3)
+        h1 = Hole(1,4)
+        # Warning mesh can't produce this type of data structure its self
+        m = Mesh(userVertices=[a,b,c,d,e],
+                 userSegments=[s1,s2,s3],
+                 regions=[r1],
+                 holes = [h1])
+        
+        fileName = tempfile.mktemp(".xya")
+        #fileName = 't.xya'
+        #os.remove(fileName)
+        m.exportPointsFile(fileName)
+        file = open(fileName)
+        lFile = file.read().split('\n')
+        file.close()
+
+        os.remove(fileName)
+        self.failUnless(lFile[0] == "" and
+                        lFile[1] == "0,0" and
+                        lFile[2] == "0,3" and
+                        lFile[3] == "3,3" 
+                        ,
+                        'exported Ascii xya file is wrong')
+        self.failUnless(lFile[4] == "1,2" and
+                        lFile[5] == "3,1" 
+                        ,
+                        'exported Ascii xya file is wrong')
+        
+        # vertex e is outside of the outline, so
+        # it is a loner and it is removed.
+        m.generateMesh("Q", maxArea = 2.1)
+        fileName = tempfile.mktemp(".xya")
+        #fileName = 't.xya'
+        #m.export_mesh_file('m.tsh')
+        m.exportPointsFile(fileName)
+        file = open(fileName)
+        lFile = file.read().split('\n')
+        file.close()
+        os.remove(fileName)
+        
+        self.failUnless(lFile[0] == "" and
+                        lFile[1] == "0.0,0.0" and
+                        lFile[2] == "0.0,3.0" and
+                        lFile[3] == "3.0,3.0" and
+                        lFile[4] == "1.0,2.0"
+                        ,
+                        'exported Ascii xya file is wrong')
+        
     def test_exportPointsFilefile2(self):
         m = Mesh()
@@ -2276,5 +2337,6 @@
 if __name__ == "__main__":
     suite = unittest.makeSuite(meshTestCase,'test')
-    #suite = unittest.makeSuite(meshTestCase,'test_add_points_and_segments')
+    #suite = unittest.makeSuite(meshTestCase,'test_exportASCIIsegmentoutlinefile')
+    #suite = unittest.makeSuite(meshTestCase,'test_e')
     runner = unittest.TextTestRunner() #verbosity=2)
     runner.run(suite)
