Robust/src/Benchmarks/SingleTM/Yada/mesh.java

   1 /* =============================================================================
   2  *
   3  * mesh.c
   4  *
   5  * =============================================================================
   6  *
   7  * Copyright (C) Stanford University, 2006.  All Rights Reserved.
   8  * Author: Chi Cao Minh
   9  *
  10  * =============================================================================
  11  *
  12  * For the license of bayes/sort.h and bayes/sort.c, please see the header
  13  * of the files.
  14  *
  15  * ------------------------------------------------------------------------
  16  *
  17  * For the license of kmeans, please see kmeans/LICENSE.kmeans
  18  *
  19  * ------------------------------------------------------------------------
  20  *
  21  * For the license of ssca2, please see ssca2/COPYRIGHT
  22  *
  23  * ------------------------------------------------------------------------
  24  *
  25  * For the license of lib/mt19937ar.c and lib/mt19937ar.h, please see the
  26  * header of the files.
  27  *
  28  * ------------------------------------------------------------------------
  29  *
  30  * For the license of lib/rbtree.h and lib/rbtree.c, please see
  31  * lib/LEGALNOTICE.rbtree and lib/LICENSE.rbtree
  32  *
  33  * ------------------------------------------------------------------------
  34  *
  35  * Unless otherwise noted, the following license applies to STAMP files:
  36  *
  37  * Copyright (c) 2007, Stanford University
  38  * All rights reserved.
  39  *
  40  * Redistribution and use in source and binary forms, with or without
  41  * modification, are permitted provided that the following conditions are
  42  * met:
  43  *
  44  *     * Redistributions of source code must retain the above copyright
  45  *       notice, this list of conditions and the following disclaimer.
  46  *
  47  *     * Redistributions in binary form must reproduce the above copyright
  48  *       notice, this list of conditions and the following disclaimer in
  49  *       the documentation and/or other materials provided with the
  50  *       distribution.
  51  *
  52  *     * Neither the name of Stanford University nor the names of its
  53  *       contributors may be used to endorse or promote products derived
  54  *       from this software without specific prior written permission.
  55  *
  56  * THIS SOFTWARE IS PROVIDED BY STANFORD UNIVERSITY ``AS IS'' AND ANY
  57  * EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
  58  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
  59  * PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL STANFORD UNIVERSITY BE LIABLE
  60  * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
  61  * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
  62  * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
  63  * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
  64  * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
  65  * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF
  66  * THE POSSIBILITY OF SUCH DAMAGE.
  67  *
  68  * =============================================================================
  69  */
  70
  71 public class mesh {
  72     element rootElementPtr;
  73     Queue_t initBadQueuePtr;
  74     int size;
  75     RBTree boundarySetPtr;
  76
  77 /* =============================================================================
  78  * mesh_alloc
  79  * =============================================================================
  80  */
  81   public mesh() {
  82     rootElementPtr = null;
  83     initBadQueuePtr = new Queue_t(-1);
  84     size = 0;
  85     boundarySetPtr = new RBTree(null, element_listCompareEdge);
  86   }
  87
  88
  89   /* =============================================================================
  90    * TMmesh_insert
  91    * =============================================================================
  92    */
  93   void TMmesh_insert (element elementPtr, avltree edgeMapPtr) {
  94     /*
  95      * Assuming fully connected graph, we just need to record one element.
  96      * The root element is not needed for the actual refining, but rather
  97      * for checking the validity of the final mesh.
  98      */
  99     if (rootElementPtr==null) {
 100       rootElementPtr=elementPtr;
 101   }
 102
 103     /*
 104      * Record existence of each of this element's edges
 105      */
 106   int numEdge = elementPtr.element_getNumEdge();
 107   for (int i = 0; i < numEdge; i++) {
 108     edge edgePtr = elementPtr.element_getEdge(i);
 109     if (!MAP_CONTAINS(edgeMapPtr, edgePtr)) {
 110       /* Record existance of this edge */
 111       boolean isSuccess =
 112         PMAP_INSERT(edgeMapPtr, edgePtr, elementPtr);
 113       yada.Assert(isSuccess);
 114     } else {
 115       /*
 116        * Shared edge; update each element's neighborList
 117        */
 118       boolean isSuccess;
 119       element sharerPtr = (element)MAP_FIND(edgeMapPtr, edgePtr);
 120       yada.Assert(sharerPtr!=null); /* cannot be shared by >2 elements */
 121       elementPtr.element_addNeighbor(sharerPtr);
 122       sharerPtr.element_addNeighbor(elementPtr);
 123       isSuccess = PMAP_REMOVE(edgeMapPtr, edgePtr);
 124       yada.Assert(isSuccess);
 125       isSuccess = PMAP_INSERT(edgeMapPtr,
 126                               edgePtr,
 127                               null); /* marker to check >2 sharers */
 128       yada.Assert(isSuccess);
 129     }
 130   }
 131
 132   /*
 133    * Check if really encroached
 134    */
 135
 136   edge encroachedPtr = elementPtr.element_getEncroachedPtr();
 137   if (encroachedPtr!=null) {
 138     if (!boundarySetPtr.contains(encroachedPtr)) {
 139       element_clearEncroached(elementPtr);
 140     }
 141   }
 142 }
 143
 144
 145 /* =============================================================================
 146  * TMmesh_remove
 147  * =============================================================================
 148  */
 149 public void TMmesh_remove(element elementPtr) {
 150   yada.Assert(!elementPtr.element_isGarbage());
 151
 152   /*
 153    * If removing root, a new root is selected on the next mesh_insert, which
 154    * always follows a call a mesh_remove.
 155    */
 156   if (rootElementPtr == elementPtr) {
 157     rootElementPtr=null;
 158   }
 159
 160   /*
 161    * Remove from neighbors
 162    */
 163   list_iter_t it;
 164   List_t neighborListPtr = elementPtr.element_getNeighborListPtr();
 165   TMLIST_ITER_RESET(it, neighborListPtr);
 166   while (TMLIST_ITER_HASNEXT(it, neighborListPtr)) {
 167       element neighborPtr = (element)TMLIST_ITER_NEXT(it, neighborListPtr);
 168       List_t neighborNeighborListPtr = neighborPtr.element_getNeighborListPtr();
 169       boolean status = neighborNeighborListPtr.remove(elementPtr);
 170       yada.Assert(status);
 171   }
 172
 173   elementPtr.element_isGarbage(true);
 174 }
 175
 176
 177 /* =============================================================================
 178  * TMmesh_insertBoundary
 179  * =============================================================================
 180  */
 181 boolean TMmesh_insertBoundary(edge boundaryPtr) {
 182   return boundarySetPtr.insert(boundaryPtr,null);
 183 }
 184
 185
 186 /* =============================================================================
 187  * TMmesh_removeBoundary
 188  * =============================================================================
 189  */
 190 boolean TMmesh_removeBoundary(edge boundaryPtr) {
 191   return boundarySetPtr.remove(boundaryPtr);
 192 }
 193
 194
 195 /* =============================================================================
 196  * createElement
 197  * =============================================================================
 198  */
 199 static void createElement (coordinate coordinates,
 200                int numCoordinate,
 201                            avltree edgeMapPtr) {
 202     element elementPtr = new element(coordinates, numCoordinate);
 203     yada.Assert(elementPtr!=null);
 204
 205     if (numCoordinate == 2) {
 206         edge boundaryPtr = elementPtr.element_getEdge(0);
 207         boolean status = boundarySetPtr.insert(boundaryPtr, null);
 208         yada.Assert(status);
 209     }
 210
 211     mesh_insert(elementPtr, edgeMapPtr);
 212
 213     if (elementPtr.element_isBad()) {
 214         boolean status = initBadQueuePtr.queue_push(elementPtr);
 215         yada.Assert(status);
 216     }
 217 }
 218
 219
 220 /* =============================================================================
 221  * mesh_read
 222  *
 223  * Returns number of elements read from file
 224  *
 225  * Refer to http://www.cs.cmu.edu/~quake/triangle.html for file formats.
 226  * =============================================================================
 227  */
 228 int mesh_read(String fileNamePrefix) {
 229     FILE inputFile;
 230     coordinate coordinates;
 231     char fileName[]=new char[256];
 232     int fileNameSize = sizeof(fileName) / sizeof(fileName[0]);
 233     char inputBuff[]=new char[256];
 234     int inputBuffSize = sizeof(inputBuff) / sizeof(inputBuff[0]);
 235     int numEntry;
 236     int numDimension;
 237     int numCoordinate;
 238     int i;
 239     int numElement = 0;
 240
 241     avltree edgeMapPtr = new avltree(0);
 242
 243     /*
 244      * Read .node file
 245      */
 246     snprintf(fileName, fileNameSize, "%s.node", fileNamePrefix);
 247     inputFile = fopen(fileName, "r");
 248     yada.Assert(inputFile);
 249     fgets(inputBuff, inputBuffSize, inputFile);
 250     sscanf(inputBuff, "%li %li", numEntry, numDimension);
 251     yada.Assert(numDimension == 2); /* must be 2-D */
 252     numCoordinate = numEntry + 1; /* numbering can start from 1 */
 253     coordinates = new coordinate[numCoordinate];
 254     for (i = 0; i < numEntry; i++) {
 255         int id;
 256         double x;
 257         double y;
 258         if (!fgets(inputBuff, inputBuffSize, inputFile)) {
 259             break;
 260         }
 261         if (inputBuff[0] == '#') {
 262             continue; /* TODO: handle comments correctly */
 263         }
 264         sscanf(inputBuff, "%li %lf %lf", id, x, y);
 265         coordinates[id].x = x;
 266         coordinates[id].y = y;
 267     }
 268     yada.Assert(i == numEntry);
 269     fclose(inputFile);
 270
 271     /*
 272      * Read .poly file, which contains boundary segments
 273      */
 274     snprintf(fileName, fileNameSize, "%s.poly", fileNamePrefix);
 275     inputFile = fopen(fileName, "r");
 276     yada.Assert(inputFile);
 277     fgets(inputBuff, inputBuffSize, inputFile);
 278     sscanf(inputBuff, "%li %li", numEntry, numDimension);
 279     yada.Assert(numEntry == 0); /* .node file used for vertices */
 280     yada.Assert(numDimension == 2); /* must be edge */
 281     fgets(inputBuff, inputBuffSize, inputFile);
 282     sscanf(inputBuff, "%li", numEntry);
 283     for (i = 0; i < numEntry; i++) {
 284         int id;
 285         int a;
 286         int b;
 287         coordinate insertCoordinates=new coordinate[2];
 288         if (!fgets(inputBuff, inputBuffSize, inputFile)) {
 289             break;
 290         }
 291         if (inputBuff[0] == '#') {
 292             continue; /* TODO: handle comments correctly */
 293         }
 294         sscanf(inputBuff, "%li %li %li", id, a, b);
 295         yada.Assert(a >= 0 && a < numCoordinate);
 296         yada.Assert(b >= 0 && b < numCoordinate);
 297         insertCoordinates[0] = coordinates[a];
 298         insertCoordinates[1] = coordinates[b];
 299         createElement(meshPtr, insertCoordinates, 2, edgeMapPtr);
 300     }
 301     yada.Assert(i == numEntry);
 302     numElement += numEntry;
 303     fclose(inputFile);
 304
 305     /*
 306      * Read .ele file, which contains triangles
 307      */
 308     snprintf(fileName, fileNameSize, "%s.ele", fileNamePrefix);
 309     inputFile = fopen(fileName, "r");
 310     yada.Assert(inputFile);
 311     fgets(inputBuff, inputBuffSize, inputFile);
 312     sscanf(inputBuff, "%li %li", numEntry, numDimension);
 313     yada.Assert(numDimension == 3); /* must be triangle */
 314     for (i = 0; i < numEntry; i++) {
 315         int id;
 316         int a;
 317         int b;
 318         int c;
 319         coordinate insertCoordinates[]=new coordinate[3];
 320         if (!fgets(inputBuff, inputBuffSize, inputFile)) {
 321             break;
 322         }
 323         if (inputBuff[0] == '#') {
 324             continue; /* TODO: handle comments correctly */
 325         }
 326         sscanf(inputBuff, "%li %li %li %li", id, a, b, c);
 327         yada.Assert(a >= 0 && a < numCoordinate);
 328         yada.Assert(b >= 0 && b < numCoordinate);
 329         yada.Assert(c >= 0 && c < numCoordinate);
 330         insertCoordinates[0] = coordinates[a];
 331         insertCoordinates[1] = coordinates[b];
 332         insertCoordinates[2] = coordinates[c];
 333         createElement(meshPtr, insertCoordinates, 3, edgeMapPtr);
 334     }
 335     yada.Assert(i == numEntry);
 336     numElement += numEntry;
 337     fclose(inputFile);
 338
 339     free(coordinates);
 340     MAP_FREE(edgeMapPtr);
 341
 342     return numElement;
 343 }
 344
 345
 346 /* =============================================================================
 347  * mesh_getBad
 348  * -- Returns NULL if none
 349  * =============================================================================
 350  */
 351 element mesh_getBad() {
 352   return (element)queue_pop(initBadQueuePtr);
 353 }
 354
 355
 356 /* =============================================================================
 357  * mesh_shuffleBad
 358  * =============================================================================
 359  */
 360 void mesh_shuffleBad (Random randomPtr) {
 361   queue_shuffle(initBadQueuePtr, randomPtr);
 362 }
 363
 364
 365 /* =============================================================================
 366  * mesh_check
 367  * =============================================================================
 368  */
 369 boolean mesh_check(int expectedNumElement) {
 370     int numBadTriangle = 0;
 371     int numFalseNeighbor = 0;
 372     int numElement = 0;
 373
 374     System.out.println("Checking final mesh:");
 375
 376     Queue_t searchQueuePtr = new Queue_t(-1);
 377     avltree visitedMapPtr = new avltree(1);
 378
 379     /*
 380      * Do breadth-first search starting from rootElementPtr
 381      */
 382     yada.Assert(rootElementPtr!=null);
 383     searchQueuePtr.queue_push(rootElementPtr);
 384     while (!searchQueuePtr.queue_isEmpty()) {
 385         list_iter_t it;
 386         List_t neighborListPtr;
 387
 388         element currentElementPtr = (element)queue_pop(searchQueuePtr);
 389         if (visitedMapPtr.contains(currentElementPtr)) {
 390             continue;
 391         }
 392         boolean isSuccess = visitedMapPtr.insert(currentElementPtr, null);
 393         yada.Assert(isSuccess);
 394         if (!currentElementPtr.checkAngles()) {
 395             numBadTriangle++;
 396         }
 397         neighborListPtr = currentElementPtr.element_getNeighborListPtr();
 398
 399         list_iter_reset(it, neighborListPtr);
 400         while (list_iter_hasNext(it, neighborListPtr)) {
 401             element neighborElementPtr =
 402                 (element)list_iter_next(it, neighborListPtr);
 403             /*
 404              * Continue breadth-first search
 405              */
 406             if (!visitedMapPtr.contains(neighborElementPtr)) {
 407                 boolean isSuccess = searchQueuePtr.queue_push(neighborElementPtr);
 408                 yada.Assert(isSuccess);
 409             }
 410         } /* for each neighbor */
 411
 412         numElement++;
 413
 414     } /* breadth-first search */
 415
 416     System.out.println("Number of elements      = "+ numElement);
 417     System.out.println("Number of bad triangles = "+ numBadTriangle);
 418
 419     return (!(numBadTriangle > 0 ||
 420               numFalseNeighbor > 0 ||
 421               numElement != expectedNumElement));
 422 }
 423 }