Changeset 4898


Ignore:
Timestamp:
Dec 31, 2007, 7:37:59 AM (17 years ago)
Author:
duncan
Message:

Outputting numeric arrays from mesh_engine. they aren't used yet though.

Location:
anuga_core/source/anuga
Files:
1 added
4 edited

Legend:

Unmodified
Added
Removed

  • anuga_core/source/anuga/mesh_engine/mesh_engine.py

    r4894 r4898  
    88#import anuga.mesh_engine.list_dic as triang
    99
    10 from Numeric import array, Float, Int32
     10from Numeric import array, Float, Int32, reshape
    1111
    1212from anuga.utilities.numerical_tools import ensure_numeric
    … …  
    1818                  mode=None, dummy_test=None):
    1919    """
    20    
     20    pointatts can be a list of lists.
     21
     22    generatedtriangleattributelist is used to represent tagged regions.
     23    #FIXME (DSG-DSG): add comments
    2124    """
    2225    #FIXME (DSG-DSG): Catch parameters that are lists,
    … …  
    5255    if pointatts is None or pointatts == []:
    5356        pointatts = [[] for x in range(points.shape[0])]
    54        
    5557    try:
    5658        # If Int is used, instead of Int32, it fails in Linux
    5759        segments = ensure_numeric(segments, Int32)
     60       
    5861    except ValueError:
    5962        msg = 'ERROR: Inconsistent segments array.'
    … …  
    6366    if segatts is None or segatts == []:
    6467        segatts = [0 for x in range(segments.shape[0])]
    65        
    6668    try:
    6769        holes = ensure_numeric(holes, Float)
    … …  
    8385   
    8486    try:
    85         #print "pointatts",pointatts
    8687        pointatts = ensure_numeric(pointatts, Float)
    8788    except (ValueError, TypeError):
    … …  
    9394        point array."""
    9495        raise ANUGAError, msg
     96    if len(pointatts.shape) == 1:
     97        pointatts = reshape(pointatts,(pointatts.shape[0],1))
    9598   
    9699    try:
    … …  
    113116        #    input vertices and vertices ``eaten'' by holes).  - output a
    114117        #    list of neighboring triangles
     118        # EG handles lone verts!
    115119           
    116120    #print "points",points
    … …  
    121125    #print "pointatts", pointatts
    122126    #print "segatts", segatts
    123     #print "mode", mode
     127    #XSprint "mode", mode
    124128    #print "yeah"
    125     mesh_dict, r_test = triang.genMesh(points,segments,holes,regions,
     129    mesh_dict, trianglelist, pointlist, pointmarkerlist, pointattributelist, triangleattributelist, segmentlist, segmentmarkerlist, neighborlist = triang.genMesh(points,segments,holes,regions,
    126130                          pointatts,segatts, mode, segments.flat)
     131    # the values as arrays
     132    mesh_dict['trianglelist'] = trianglelist
     133    mesh_dict['pointlist'] = pointlist
     134
     135    # WARNING - pointmarkerlist IS UNTESTED
     136    mesh_dict['pointmarkerlist'] = pointmarkerlist
     137    mesh_dict['pointattributelist'] = pointattributelist
     138    mesh_dict['triangleattributelist'] = triangleattributelist
     139    mesh_dict['segmentlist'] = segmentlist
     140    mesh_dict['segmentmarkerlist'] =  segmentmarkerlist
     141    mesh_dict['triangleneighborlist'] = neighborlist
    127142    mesh_dict['qaz'] = 1 #debugging
    128143    ##print "r_test", r_test
     144   
    129145    return mesh_dict
    130146

  • anuga_core/source/anuga/mesh_engine/mesh_engine_c_layer.c

    r4894 r4898  
    3636     shape and mesh_engine.
    3737
    38      to return numeric arrays, check how it is done in
    39      quantity_ext.c compute_gradients
    4038         
    4139     Precondition
    … …  
    9290 
    9391 
    94   PyArrayObject *r_test;
     92  PyArrayObject *gentrianglelist;
     93  PyArrayObject *genpointlist;
     94  PyArrayObject *genseglist;
     95  PyArrayObject *genpointmarkerlist;
     96  PyArrayObject *genpointattributelist;
     97  PyArrayObject *gentriangleattributelist;
     98  PyArrayObject *gensegmentlist;
     99  PyArrayObject *gensegmentmarkerlist;
     100  PyArrayObject *genneighborlist;
    95101  PyObject *R;
    96102
    97   int dimensions[1];
     103  int dimensions[2];
    98104   
    99105  REAL Attr;
    … …  
    149155 
    150156  /* Point attribute list */
     157   /*printf ("in.pointattributelist -> dimensions[0] %i\n", pointattributelist -> dimensions[0]);
     158  printf ("in.pointattributelist -> dimensions[1] %i\n", pointattributelist -> dimensions[1]); */
    151159  if (0 == pointattributelist -> dimensions[0]) {
    152160    in.numberofpointattributes = 0;
    … …  
    220228 
    221229 
     230 
     231  /* printf(" ***  back from triangulate\n" );    */
    222232  /*
     233    ------- Pass point numbers,coordinates and neighbors back to Python ------
     234    we send back a dictionary:                                               
     235    { index : [ coordinates, [connections], Attribute ] }
     236  */
     237  holder = PyDict_New();   
     238 
     239  /*
     240 
     241     to return numeric arrays, check how it is done in
     242     abs quantity_ext.c compute_gradients
     243     
    223244  PyArray_FromDims allolws you to create a Numeric array with unitialized data.
    224245   The first argument is the size of the second argument (
    … …  
    230251  */
    231252 
    232   //Py_Initialize();
    233   // Testing passing a numeric array out
    234   dimensions[0] = 4;
    235   // Allocate space for return vectors a and b (don't DECREF)
    236   r_test = (PyArrayObject *) PyArray_FromDims(1, dimensions, PyArray_DOUBLE);
    237  
    238   /* printf(" ***  back from triangulate\n" );    */
    239   /*
    240     ------- Pass point numbers,coordinates and neighbors back to Python ------
    241     we send back a dictionary:                                               
    242     { index : [ coordinates, [connections], Attribute ] }
    243   */
    244   holder = PyDict_New();   
    245  
     253 
     254  /* Add triangle list */
     255  dimensions[0] = out.numberoftriangles;
     256  dimensions[1] = 3;   
     257  gentrianglelist = (PyArrayObject *) PyArray_FromDims(2,
     258                                                    dimensions,
     259                                                    PyArray_INT);
     260  gentrianglelist -> data = out.trianglelist;
     261   
    246262  /* Add triangle list */
    247263  listsize = out.numberoftriangles;
    … …  
    255271  ii=PyString_FromString("generatedtrianglelist");
    256272  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);
     273 
     274 
     275  /* Add pointlist */
     276  dimensions[0] = out.numberofpoints;
     277  dimensions[1] = 2;   
     278  genpointlist = (PyArrayObject *) PyArray_FromDims(2,
     279                                                 dimensions,
     280                                                 PyArray_DOUBLE);
     281  /*                                             
     282  (double*) genpointlist -> data = out.pointlist;               
     283  ((double*) genpointlist -> data) = out.pointlist;
     284  ( genpointlist -> data) = (double*) out.pointlist;
     285 
     286  */
     287  genpointlist -> data = out.pointlist;
    257288     
    258289  /* Add pointlist */
    … …  
    269300 
    270301  /* Add point marker list */
     302  dimensions[0] = out.numberofpoints;
     303  genpointmarkerlist = (PyArrayObject *) PyArray_FromDims(1,
     304                                                    dimensions,
     305                                                    PyArray_INT);
     306  genpointmarkerlist -> data = out.pointmarkerlist;
     307 
     308  /* Add point marker list */
    271309  listsize = out.numberofpoints;
    272310  holderlist = PyList_New(listsize);
    273311 
    274312  for(i=0; i<listsize;i++){
    275     PyObject *mlist = Py_BuildValue((char *)"d",
     313    PyObject *mlist = Py_BuildValue((char *)"i",
    276314                                    out.pointmarkerlist[i]);     
    277315    PyList_SetItem(holderlist,i, mlist);
    … …  
    280318  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist);
    281319   
     320  /* Add point attribute list */
     321  dimensions[0] = out.numberofpoints;
     322  dimensions[1] = out.numberofpointattributes;   
     323  genpointattributelist = (PyArrayObject *) PyArray_FromDims(2,
     324                                                 dimensions,
     325                                                 PyArray_DOUBLE);
     326  genpointattributelist -> data = out.pointattributelist;       
     327 
     328 
    282329  /* Add point attribute list */
    283330  listsize = out.numberofpoints;
    … …  
    297344  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist); 
    298345 
     346  /* Add triangle attribute list */
     347  dimensions[0] = out.numberoftriangles;
     348  dimensions[1] = out.numberoftriangleattributes;   
     349  gentriangleattributelist = (PyArrayObject *) PyArray_FromDims(2,
     350                                                 dimensions,
     351                                                 PyArray_DOUBLE);
     352  gentriangleattributelist -> data = out.triangleattributelist;
    299353 
    300354  /* Add triangle attribute list */
    … …  
    315369 
    316370  /* Add segment list */
     371  dimensions[0] = out.numberofsegments;
     372  dimensions[1] = 2;   
     373  gensegmentlist = (PyArrayObject *) PyArray_FromDims(2,
     374                                                    dimensions,
     375                                                    PyArray_INT);
     376  gensegmentlist -> data = out.segmentlist;
     377 
     378 
     379  /* Add segment list */
    317380  listsize = out.numberofsegments;
    318381  holderlist = PyList_New(listsize);
    … …  
    327390 
    328391  /* Add segment marker list */
     392  dimensions[0] = out.numberofsegments;
     393  gensegmentmarkerlist = (PyArrayObject *) PyArray_FromDims(1,
     394                                                    dimensions,
     395                                                    PyArray_INT);
     396  gensegmentmarkerlist -> data = out.segmentmarkerlist;
     397 
     398 
     399  /* Add segment marker list */
    329400  listsize = out.numberofsegments;
    330401  holderlist = PyList_New(listsize);
    … …  
    335406  }   
    336407  ii=PyString_FromString("generatedsegmentmarkerlist");
    337   PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii); Py_DECREF(holderlist); 
     408  PyDict_SetItem(holder, ii, holderlist); Py_DECREF(ii);
     409  Py_DECREF(holderlist); 
     410 
     411  /* Add triangle neighbor list */
     412  dimensions[0] = out.numberoftriangles;
     413  dimensions[1] = 3;   
     414  genneighborlist = (PyArrayObject *) PyArray_FromDims(2,
     415                                                    dimensions,
     416                                                    PyArray_INT);
     417  genneighborlist -> data = out.neighborlist;
     418 
     419 
    338420  /* Add triangle neighbor list */
    339421  if (out.neighborlist != NULL) {
    … …  
    399481 
    400482  /* R = Py_BuildValue((char *)"O", holder); */
    401   R = Py_BuildValue((char *)"OO", holder, PyArray_Return(r_test));
     483  R = Py_BuildValue((char *)"OOOOOOOOO", holder
     484                    ,PyArray_Return(gentrianglelist)
     485                    ,PyArray_Return(genpointlist)
     486                    ,PyArray_Return(genpointmarkerlist)
     487                    ,PyArray_Return(genpointattributelist)
     488                    ,PyArray_Return(gentriangleattributelist)
     489                    ,PyArray_Return(gensegmentlist)
     490                    ,PyArray_Return(gensegmentmarkerlist)
     491                    ,PyArray_Return(genneighborlist)
     492                    );
    402493  Py_DECREF(holder); /** This fixed a  memory problem ticket#189 */
    403   Py_DECREF(r_test);
     494  Py_DECREF(gentrianglelist);
     495  Py_DECREF(genpointlist);
     496  Py_DECREF(genpointmarkerlist);
     497  Py_DECREF(genpointattributelist);
     498  Py_DECREF(gentriangleattributelist);
     499  Py_DECREF(gensegmentlist);
     500  Py_DECREF(gensegmentmarkerlist);
     501  Py_DECREF(genneighborlist);
    404502  return R;
    405503}

  • anuga_core/source/anuga/mesh_engine/test_generate_mesh.py

    r4893 r4898  
    159159    def testsegmarker(self):
    160160
    161         points = []
    162         seglist = []
    163161        holelist = []
    164162        regionlist = []
    … …  
    259257        # for Ole on nautilus this returns 6
    260258        # for Duncan on nautilus this returns 7
    261         # ??, it seems to be the results from triangle that is
     259        # It seems to be the results from triangle that is
    262260        # causing the different results, and we are treating
    263261        # triangle as a back box.
    … …  
    526524        self.failUnless(data['lonepointlist'] ==[0],
    527525                        'lonepointlist is wrong!')
     526
     527
     528
     529    def test_transition_to_arrays(self):
     530        # segattlist = []
     531        points = []
     532        seglist = []
     533        holelist = []
     534        regionlist = []
     535
     536        points = [(0.0,0.0),(0.0,10.0),(3.0,0.0),(3.0,10.0)]
     537        pointattlist = []
     538        # 5.0 is the region tag, 99.0 is the max area
     539        # REgion tag isn't working though
     540        regionlist.append( [0.2,0.2,5.0,99.0] )
     541        seglist = [(0,1),(1,3),(3,2),(2,0)]
     542        segattlist = [21,22,23,24]
     543         #The 'A' has to be there to get the region marker stuff working
     544        mode = "QzpnA"
     545        #mode = "jQpznAa2000.1a"
     546        data = generate_mesh(points,seglist,holelist,regionlist,
     547                              pointattlist,segattlist, mode, points)
     548        #print "data", data
     549
     550        self.failUnless(data['generatedtrianglelist'] ==[(1, 0, 2), (2, 3, 1)],
     551                        'trianglelist is wrong!')
     552        #print "data['generatedtrianglelist']",data['generatedtrianglelist']
     553        #print "data['trianglelist']", data['trianglelist']
     554        self.failUnless(data['generatedtrianglelist'] ==data['trianglelist'],
     555                        'trianglelist is wrong!')
     556        self.failUnless(data['generatedsegmentlist'] ==[(0, 1), (1, 3),
     557                                                        (3, 2), (2, 0)],
     558                        'segmentlist is wrong!')
     559        self.failUnless(data['generatedpointlist'] ==[(0.0, 0.0), (0.0, 10.0),
     560                                                      (3.0, 0.0), (3.0, 10.0)],
     561                        ' is wrong!')
     562        self.failUnless(data['generatedpointlist'] ==data['pointlist'],
     563                        ' is wrong!')
     564        self.failUnless(data['generatedtriangleattributelist'] == [[5.0],
     565                                                                   [5.0]],
     566                        ' is wrong!')
     567       
     568        self.failUnless(data['generatedtriangleattributelist'] == \
     569                        data['triangleattributelist'],
     570                        ' is wrong!')
     571        self.failUnless(data['generatedsegmentlist'] == seglist,
     572                        ' is wrong!')
     573        self.failUnless(data['generatedsegmentlist'] == data['segmentlist'],
     574                        ' is wrong!')
     575        self.failUnless(data['generatedsegmentmarkerlist'] == segattlist,
     576                        ' is wrong!')
     577        self.failUnless(data['generatedsegmentmarkerlist'] == \
     578                        data['segmentmarkerlist'],
     579                        ' is wrong!')
     580        # I copied these answers from the output, so bad test..
     581        self.failUnless(data['generatedtriangleneighborlist'] == \
     582                        [(-1, 1, -1), (-1, 0, -1)],
     583                        ' is wrong!')
     584        self.failUnless(data['generatedtriangleneighborlist'] == \
     585                        data['triangleneighborlist'],
     586                        ' is wrong!')
     587                       
     588    def test_pointattlist(self):
     589        # segattlist = []
     590        points = []
     591        seglist = []
     592        holelist = []
     593        regionlist = []
     594
     595        points = [(0.0,0.0),(0.0,4.0),(4.0,2.0),(2.0,0.0)]
     596        pointattlist = [0.,0.,10.,10.]
     597        regionlist.append( [0.2,0.2,2.1, 99.] )
     598        seglist = [(0,1),(1,2),(2,3),(3,0)]
     599        segattlist = [11,12,13,14]
     600        mode = "Qzp"
     601        data = generate_mesh(points,seglist,holelist,regionlist,
     602                              pointattlist,segattlist, mode, points)
     603        self.failUnless(data['generatedpointattributelist'] == [[0.0],[0.0],
     604                                                                [10],[10]],
     605                        ' is wrong!')
     606        self.failUnless(data['generatedpointattributelist'] == \
     607                        data['pointattributelist'],
     608                        ' is wrong!')
     609       
     610       
     611        pointattlist = [[0.],[0.],[10.],[10.]]
     612        mode = "Qzp"       
     613        data = generate_mesh(points,seglist,holelist,regionlist,
     614                              pointattlist,segattlist, mode, points)
     615        self.failUnless(data['generatedpointattributelist'] == [[0.0],[0.0],
     616                                                                [10],[10]],
     617                        ' is wrong!')
     618        self.failUnless(data['generatedpointattributelist'] == \
     619                        data['pointattributelist'],
     620                        ' is wrong!')
     621        pointattlist = [[0.,1],[0.,1],[10.,20],[10.,20]]
     622        mode = "Qzp"       
     623        data = generate_mesh(points,seglist,holelist,regionlist,
     624                              pointattlist,segattlist, mode, points)
     625        #print "data", data
     626        self.failUnless(data['generatedpointattributelist'] == pointattlist,
     627                        ' is wrong!')
     628        self.failUnless(data['generatedpointattributelist'] == \
     629                        data['pointattributelist'],
     630                        ' is wrong!')
     631     
    528632if __name__ == "__main__":
    529633
    530634    suite = unittest.makeSuite(triangTestCase,'test')
    531635    #suite = unittest.makeSuite(triangTestCase,'test_lone_verts4')
    532     #suite = unittest.makeSuite(triangTestCase,'testrectangleIIb')
     636    #suite = unittest.makeSuite(triangTestCase,'test_transition_to_arrays')
    533637    runner = unittest.TextTestRunner() #verbosity=2)
    534638    runner.run(suite)

  • anuga_core/source/anuga/pmesh/mesh.py

    r4894 r4898  
    11631163        #FIXME (DSG-DSG)  move below section into generate_mesh.py
    11641164        #                  & 4 functions eg segment_strings2ints
     1165        # Actually, because of region_list.append((1.0,2.0,""))
     1166        # don't move it, without careful thought
    11651167        #print "*************************!@!@ This is going to triangle   !@!@"
    11661168        #print meshDict
Note: See TracChangeset for help on using the changeset viewer.