This version works for the explicit FEA code provided by Professor Paolino.
authorIsaac Dooley <idooley2@illinois.edu>
Sat, 3 Mar 2007 03:21:57 +0000 (03:21 +0000)
committerIsaac Dooley <idooley2@illinois.edu>
Sat, 3 Mar 2007 03:21:57 +0000 (03:21 +0000)
src/libs/ck-libs/ParFUM-Tops/Makefile
src/libs/ck-libs/ParFUM-Tops/ParFUM_TOPS.C
src/libs/ck-libs/ParFUM-Tops/ParFUM_TOPS.h
src/libs/ck-libs/ParFUM-Tops/ParFUM_TOPS_Iterators.C

index 515c54dfe666a3d9116aa8a0609a73579f5bbc61..90f51666eca9f1f441521ad27069f6eb4e660585 100644 (file)
@@ -1,5 +1,5 @@
 CDIR=../../../..
-CHARMC=$(CDIR)/bin/charmc $(OPTS) -g
+CHARMC=$(CDIR)/bin/charmc $(OPTS)
 
 #Headers to be copied to include directory so application users can see them
 HEADERS= ParFUM_TOPS.h  $(INTERNALHEADERS)
index 9600f14ae7207d377122b10c163cc7d2a234903a..938f577408cdcd439762181c1f7ec7ff5fe8a581 100644 (file)
@@ -8,15 +8,22 @@
    @todo add code to generate ghost layers
    @todo Support multiple models
 
+   @note FEM_DATA+0 holds the elemAttr or nodeAtt data
+   @note FEM_DATA+1 holds the id 
+   @note FEM_CONN holds the element connectivity
+   @note FEM_COORD holds the nodal coordinates
+
 */
 
 #include "ParFUM_TOPS.h"
 #include "ParFUM.decl.h"
+#include "ParFUM_internals.h"
 
 int elem_attr_size, node_attr_size;
 
 TopModel* topModel_Create_Init(int elem_attr_sz, int node_attr_sz){
-
+  CkAssert(elem_attr_sz > 0);
+  CkAssert(node_attr_sz > 0);
   elem_attr_size = elem_attr_sz;
   node_attr_size = node_attr_sz;
 
@@ -27,27 +34,40 @@ TopModel* topModel_Create_Init(int elem_attr_sz, int node_attr_sz){
   char* temp_array = new char[16]; // just some junk array to use below
   
   // Allocate element connectivity
-  FEM_Mesh_data(which_mesh,FEM_ELEM+0,FEM_CONN,temp_array, 0, 0, FEM_INDEX_0, 3);
-  // Allocate node coords
-  FEM_Mesh_data(which_mesh,FEM_NODE,FEM_COORD,temp_array, 0, 0, FEM_DOUBLE, 3);
+  FEM_Mesh_data(which_mesh,FEM_ELEM+0,FEM_CONN,temp_array, 0, 0, FEM_INDEX_0, 4);
   // Allocate element attributes
   FEM_Mesh_data(which_mesh,FEM_ELEM+0,FEM_DATA+0,temp_array, 0, 0, FEM_BYTE, elem_attr_size);
+  // Allocate element Id array
+  FEM_Mesh_data(which_mesh,FEM_ELEM+0,FEM_DATA+1,temp_array, 0, 0, FEM_INT, 1);
+
+  // Allocate node coords
+  FEM_Mesh_data(which_mesh,FEM_NODE,FEM_COORD,temp_array, 0, 0, FEM_DOUBLE, 3);
   // Allocate node attributes
   FEM_Mesh_data(which_mesh,FEM_NODE+0,FEM_DATA+0,temp_array, 0, 0, FEM_BYTE, node_attr_size);
+  // Allocate node Id array
+  FEM_Mesh_data(which_mesh,FEM_NODE+0,FEM_DATA+1,temp_array, 0, 0, FEM_INT, 1);
 
   delete[] temp_array;
 
   // Don't Allocate the Global Number attribute for the elements and nodes  
   // It will be automatically created upon calls to void FEM_Entity::setGlobalno(int r,int g) {
-  
+
   FEM_Mesh_allocate_valid_attr(which_mesh, FEM_NODE);
   FEM_Mesh_allocate_valid_attr(which_mesh, FEM_ELEM+0);
   
+
+  // Setup ghost layers, assume Tet4
+  // const int tet4vertices[4] = {0,1,2,3};
+  //  FEM_Add_ghost_layer(1,1);
+  // FEM_Add_ghost_elem(0,4,tet4vertices);
+
   return mesh;
 }
 
 TopModel* topModel_Create_Driver(int elem_attr_sz, int node_attr_sz){
   // This only uses a single mesh, so don't create multiple TopModels of these
+  CkAssert(elem_attr_sz > 0);
+  CkAssert(node_attr_sz > 0);
   elem_attr_size = elem_attr_sz;
   node_attr_size = node_attr_sz;
   int which_mesh=FEM_Mesh_default_read();
@@ -66,15 +86,18 @@ TopNode topModel_InsertNode(TopModel* m, double x, double y, double z){
 }
 
 
-/** Set id of a node */
+/** Set id of a node 
+@todo Make this work with ghosts
+*/
 void topNode_SetId(TopModel* m, TopNode n, TopID id){
-       // just set the nodal id attribute to id
-  m->node.setGlobalno(n,id);
+  CkAssert(n>=0);
+  FEM_DataAttribute * at = (FEM_DataAttribute*) m->node.lookup(FEM_DATA+1,"topNode_SetId");
+  at->getInt()(n,0)=id;
 }
 
 /** Insert an element */
 TopElement topModel_InsertElem(TopModel*m, TopElemType type, TopNode* nodes){
-  assert(type ==  TOP_ELEMENT_TET4);
+  CkAssert(type ==  TOP_ELEMENT_TET4);
   int conn[4];
   conn[0] = nodes[0];
   conn[1] = nodes[1];
@@ -84,9 +107,13 @@ TopElement topModel_InsertElem(TopModel*m, TopElemType type, TopNode* nodes){
   return newEl;
 }
 
-/** Set id of an element */
+/** Set id of an element 
+@todo Make this work with ghosts
+*/
 void topElement_SetId(TopModel* m, TopElement e, TopID id){
-  m->elem[0].setGlobalno(e,id);
+  CkAssert(e>=0);
+  FEM_DataAttribute * at = (FEM_DataAttribute*) m->elem[0].lookup(FEM_DATA+1,"topElement_SetId");
+  at->getInt()(e,0)=id;
 }
 
 
@@ -102,6 +129,7 @@ void topElement_SetId(TopModel* m, TopElement e, TopID id){
 
 */
 void topNode_SetAttrib(TopModel* m, TopNode n, void* d){
+  CkAssert(n>=0);
   FEM_DataAttribute * at = (FEM_DataAttribute*) m->node.lookup(FEM_DATA+0,"topNode_SetAttrib");
   AllocTable2d<unsigned char> &dataTable  = at->getChar();
   unsigned char *data = dataTable.getData();
@@ -112,6 +140,7 @@ void topNode_SetAttrib(TopModel* m, TopNode n, void* d){
 See topNode_SetAttrib() for description
 */
 void topElement_SetAttrib(TopModel* m, TopElement e, void* d){
+  CkAssert(e>=0);
   FEM_DataAttribute * at = (FEM_DataAttribute*) m->elem[0].lookup(FEM_DATA+0,"topElem_SetAttrib");
   AllocTable2d<unsigned char> &dataTable  = at->getChar();
   unsigned char *data = dataTable.getData();
@@ -123,6 +152,9 @@ void topElement_SetAttrib(TopModel* m, TopElement e, void* d){
 See topNode_SetAttrib() for description
 */
 void* topElement_GetAttrib(TopModel* m, TopElement e){
+  if(! m->elem[0].is_valid_any_idx(e))
+       return NULL;
+
   FEM_DataAttribute * at = (FEM_DataAttribute*) m->elem[0].lookup(FEM_DATA+0,"topElem_GetAttrib");
   AllocTable2d<unsigned char> &dataTable  = at->getChar();
   unsigned char *data = dataTable.getData();
@@ -133,6 +165,9 @@ void* topElement_GetAttrib(TopModel* m, TopElement e){
 See topNode_SetAttrib() for description
 */
 void* topNode_GetAttrib(TopModel* m, TopNode n){
+  if(! m->node.is_valid_any_idx(n))
+       return NULL;
+
   FEM_DataAttribute * at = (FEM_DataAttribute*) m->node.lookup(FEM_DATA+0,"topNode_GetAttrib");
   AllocTable2d<unsigned char> &dataTable  = at->getChar();
   unsigned char *data = dataTable.getData();
@@ -143,42 +178,117 @@ void* topNode_GetAttrib(TopModel* m, TopNode n){
 
 /** 
        Get node via id 
-       @todo Implement this function
+       @todo Re-implement this function with some hash to make it fast. 
+       @note In the explicit FEA example, this is just used during initialization, so speed is not too important.
+       @todo Does not work with ghosts yet.
 */
-TopNode topModel_GetNodeAtId(TopModel*,TopID){
+TopNode topModel_GetNodeAtId(TopModel* m, TopID id){
   // lookup node via global ID
-  assert(0);
+  FEM_DataAttribute * at = (FEM_DataAttribute*) m->node.lookup(FEM_DATA+1,"topModel_GetNodeAtId");
+  for(int i=0;i<at->getInt().size();++i){
+       if(at->getInt()(i,0)==id){
+         return i;
+       }
+  }
+  return -1;
 }
 
 /** 
        Get elem via id
        @todo Implement this function
+       @note Is this function even supposed to exist?
  */
-TopElement topModel_GetElemAtId(TopModel*,TopID){
-  assert(0);
+TopElement topModel_GetElemAtId(TopModel*m,TopID id){
+  CkPrintf("Error: topModel_GetElemAtId() not yet implemented");
+  CkExit();
 }
 
 
 
 TopNode topElement_GetNode(TopModel* m,TopElement e,int idx){
-  assert(0);
-}
+  CkAssert(e>=0);
+  const AllocTable2d<int> &conn = m->elem[0].getConn();
+  CkAssert(idx>=0 && idx<conn.width());
+  CkAssert(idx<conn.size());
 
+  int node = conn(e,idx);
+
+  return conn(e,idx);
+}
 
 int topNode_GetId(TopModel* m, TopNode n){
-  assert(0);
+  CkAssert(n>=0);
+  FEM_DataAttribute * at = (FEM_DataAttribute*) m->node.lookup(FEM_DATA+1,"topNode_SetId");
+  return at->getInt()(n,0);
 }
 
+
+/** @todo handle ghost nodes as appropriate */
 int topModel_GetNNodes(TopModel *model){
-  assert(0);
+  return model->node.count_valid();
 }
 
+/** @todo Fix to return the width of the conn array */
 int topElement_GetNNodes(TopModel* model, TopElement elem){
-  assert(0);
+  return 4;
 }
 
+/** @todo make sure we are in a getting mesh */
 void topNode_GetPosition(TopModel*model, TopNode node,double*x,double*y,double*z){
-  assert(0);
+  CkAssert(node>=0);
+  double coord[3]={7.01,7.01,7.01};
+
+  // lookup node via global ID
+  FEM_DataAttribute * at = (FEM_DataAttribute*) model->node.lookup(FEM_COORD,"topModel_GetNodeAtId");
+  CkAssert(at->getDouble().width()==3);
+  //  CkPrintf("node=%d, size=%d\n",node,at->getDouble().size() );
+  CkAssert(node<at->getDouble().size());
+  *x = at->getDouble()(node,0);
+  *y = at->getDouble()(node,1);
+  *z = at->getDouble()(node,2);
+}
+
+void topModel_Sync(TopModel*m){
+  MPI_Barrier(MPI_COMM_WORLD);
+
+  //  CkPrintf("%d: %d local, %d ghost elements\n", FEM_My_partition(), m->elem[0].size(),m->elem[0].ghost->size() );
+  //  CkPrintf("%d: %d local, %d ghost valid elements\n", FEM_My_partition(), m->elem[0].count_valid(),m->elem[0].ghost->count_valid() );
+
+}
+
+void topModel_TestIterators(TopModel*m){
+  CkAssert(m->elem[0].ghost!=NULL);
+  CkAssert(m->node.ghost!=NULL);
+
+  int expected_elem_count = m->elem[0].count_valid() + m->elem[0].ghost->count_valid(); 
+  int iterated_elem_count = 0;
+
+  int expected_node_count = m->node.count_valid() + m->node.ghost->count_valid(); 
+  int iterated_node_count = 0;
+
+  int myId = FEM_My_partition();
+
+
+  TopNodeItr* itr = topModel_CreateNodeItr(m);
+  for(topNodeItr_Begin(itr);topNodeItr_IsValid(itr);topNodeItr_Next(itr)){
+       iterated_node_count++;
+       TopNode node = topNodeItr_GetCurr(itr);
+       void* na = topNode_GetAttrib(m,node);
+       CkAssert(na != NULL);
+  }
+    
+  TopElemItr* e_itr = topModel_CreateElemItr(m);
+  for(topElemItr_Begin(e_itr);topElemItr_IsValid(e_itr);topElemItr_Next(e_itr)){
+       iterated_elem_count++;
+       TopElement elem = topElemItr_GetCurr(e_itr);
+       void* ea = topElement_GetAttrib(m,elem);
+       CkAssert(ea != NULL);
+  }
+  
+  CkAssert(iterated_node_count == expected_node_count);
+  CkAssert(iterated_elem_count==expected_elem_count);
+  
 }
 
 
index 9d68fd137b4238cb2b9b878652e8671bee721ffa..d2cec40105217d6ad504543f18ce2dc5e1422074 100644 (file)
@@ -61,9 +61,9 @@ Sample usage:
 typedef FEM_Mesh TopModel;
 
 /** Tops uses some bit patterns for these, but we just use TopNode as a signed value to represent the corresponding ParFUM node. A non-negative value is a local node, while a negative value is a ghost. */
-typedef unsigned TopNode;
+typedef long TopNode;
 /** See notes for ::TopNode */
-typedef unsigned TopElement;
+typedef long TopElement;
 
 
 enum {
@@ -123,6 +123,9 @@ TopModel* topModel_Create_Driver(int elem_attr_sz, int node_attr_sz);
 /** Cleanup a model. Currently does nothing */
 void topModel_Destroy(TopModel* m);
 
+/** Synchronize element and node data for ghost layers. Should be called every timestep */
+void topModel_Sync(TopModel*m);
+
 /** Insert a node */
 TopNode topModel_InsertNode(TopModel*, double x, double y, double z);
 
@@ -144,9 +147,6 @@ void topElement_SetAttrib(TopModel*, TopElement, void*);
 /** Get node via id */
 TopNode topModel_GetNodeAtId(TopModel*,TopID);
 
-/** Get elem via id */
-TopElement topModel_GetElemAtId(TopModel*,TopID);
-
 /** Get nodal attribute */
 void* topNode_GetAttrib(TopModel*, TopNode);
 
@@ -158,7 +158,6 @@ TopNode topElement_GetNode(TopModel*,TopElement,int idx);
 
 int topNode_GetId(TopModel* m, TopNode n);
 
-
 int topModel_GetNNodes(TopModel *model);
 
 int topElement_GetNNodes(TopModel* model, TopElement elem);
@@ -202,5 +201,8 @@ void topElemItr_Next(TopElemItr*);
 /** Get TopElement associated with the iterator */
 TopElement topElemItr_GetCurr(TopElemItr*);
 
+/** Perform sanity check on iterators. This checks to make sure that the count of the itereated elements and nodes matches that returned by ParFUM's countValid() */
+void topModel_TestIterators(TopModel*m);
+
 
 #endif
index b13f5a2624796f09b87df0a9617ceb3d835e9a50..05d2628059e52799a007df803a8b76337bc88390 100644 (file)
@@ -1,5 +1,8 @@
 #include "ParFUM_TOPS.h"
 
+
+
+
 /**************************************************************************
  *     Iterator for nodes 
  */
@@ -15,54 +18,51 @@ void topNodeItr_Destroy(TopNodeItr* itr){
 }
 
 void topNodeItr_Begin(TopNodeItr* itr){
-    itr->parfum_index = 0;
+  if(itr->model->node.ghost != NULL){
+       itr->parfum_index =  FEM_To_ghost_index(itr->model->node.ghost->size());
+  }
+  else{
+       itr->parfum_index =  0;
+  }  
+
+  // Make sure we start with a valid one:
+  while((!itr->model->node.is_valid_any_idx(itr->parfum_index)) &&
+               (itr->parfum_index < itr->model->node.size()))
+       itr->parfum_index++;
+  
+  if(itr->parfum_index==itr->model->node.size()){
+       itr->parfum_index = itr->model->node.size()+1000; // way past the end
+  }  
+
+#ifdef PTOPS_ITERATOR_PRINT
+  CkPrintf("Initializing Node Iterator to %d\n", itr->parfum_index);
+#endif
 }
 
 bool topNodeItr_IsValid(TopNodeItr*itr){
-     return itr->model->node.is_valid_any_idx(itr->parfum_index);
+       return itr->model->node.is_valid_any_idx(itr->parfum_index);
 }
 
 void topNodeItr_Next(TopNodeItr* itr){
-
-    if(!topNodeItr_IsValid(itr))
-        return;
-
-    // advance index until we hit a valid index
-    itr->parfum_index++;
-
-    if(itr->parfum_index > 0) {// local nodes
-
-        while ((! itr->model->node.is_valid_any_idx(itr->parfum_index)) &&
-                  (itr->parfum_index<itr->model->node.size()))
-        {
-            itr->parfum_index++;
-        }
-
-        if(itr->model->node.is_valid_any_idx(itr->parfum_index)) {
-            return;
-        } else {
-            // cycle to most negative index possible for ghosts
-                 printf("Node iterator switched to ghosts\n");
-            itr->parfum_index = FEM_To_ghost_index(itr->model->node.ghost->size());
-        }
-    }
-
-    // just go through ghost nodes
-    
-    while ( (! itr->model->node.ghost->
-                is_valid_any_idx(FEM_To_ghost_index(itr->parfum_index)))
-                &&
-                itr->parfum_index<0)
-    {
-        itr->parfum_index++;
-    }
-
-    if(itr->parfum_index==0){
-        itr->parfum_index = itr->model->node.size()+1000; // way past the end
-    }
-
+  CkAssert(topNodeItr_IsValid(itr));
+  
+  // advance index until we hit a valid index
+  itr->parfum_index++;
+
+  while((!itr->model->node.is_valid_any_idx(itr->parfum_index)) &&
+               (itr->parfum_index < itr->model->node.size()))
+       itr->parfum_index++;
+  
+  if(itr->parfum_index==itr->model->node.size()){
+       itr->parfum_index = itr->model->node.size()+1000; // way past the end
+  }  
+
+#ifdef PTOPS_ITERATOR_PRINT
+  CkPrintf("Advancing Node Iterator to %d\n", itr->parfum_index);
+#endif
 }
 
+
 TopNode topNodeItr_GetCurr(TopNodeItr*itr){
        return itr->parfum_index;
 }
@@ -83,53 +83,53 @@ void topElemItr_Destroy(TopElemItr* itr){
 }
 
 void topElemItr_Begin(TopElemItr* itr){
-    itr->parfum_index = 0;
+  if(itr->model->elem[0].ghost != NULL){
+       itr->parfum_index =  FEM_To_ghost_index(itr->model->elem[0].ghost->size());
+  }
+  else{
+       itr->parfum_index =  0;
+  }  
+
+  // Make sure we start with a valid one:
+  while((!itr->model->elem[0].is_valid_any_idx(itr->parfum_index)) &&
+               (itr->parfum_index < itr->model->elem[0].size()))
+       itr->parfum_index++;
+  
+  if(itr->parfum_index==itr->model->elem[0].size()){
+       itr->parfum_index = itr->model->elem[0].size()+1000; // way past the end
+  }  
+
+#ifdef PTOPS_ITERATOR_PRINT
+  CkPrintf("Initializing Elem[0] Iterator to %d\n", itr->parfum_index);
+#endif
+
 }
 
 bool topElemItr_IsValid(TopElemItr*itr){
-     return itr->model->elem[0].is_valid_any_idx(itr->parfum_index);
+  return itr->model->elem[0].is_valid_any_idx(itr->parfum_index);
 }
 
 void topElemItr_Next(TopElemItr* itr){
-
-    if(!topElemItr_IsValid(itr))
-        return;
-
-    // advance index until we hit a valid index
-    itr->parfum_index++;
-
-    if(itr->parfum_index > 0) {// non-ghosts
-
-        while ((! itr->model->elem[0].is_valid_any_idx(itr->parfum_index)) &&
-                  (itr->parfum_index<itr->model->elem[0].size()))
-        {
-            itr->parfum_index++;
-        }
-
-        if(itr->model->elem[0].is_valid_any_idx(itr->parfum_index)) {
-            return;
-        } else {
-            // cycle to most negative index possible for ghosts
-                 printf("Elem iterator switched to ghosts\n");
-            itr->parfum_index = FEM_To_ghost_index(itr->model->elem[0].ghost->size());
-        }
-    }
-
-    // just go through ghosts    
-    while ( (! itr->model->elem[0].ghost->
-                is_valid_any_idx(FEM_To_ghost_index(itr->parfum_index)))
-                &&
-                itr->parfum_index<0)
-    {
-        itr->parfum_index++;
-    }
-
-    if(itr->parfum_index==0){
-        itr->parfum_index = itr->model->elem[0].size()+1000; // way past the end
-    }
-
+  CkAssert(topElemItr_IsValid(itr));
+  
+  // advance index until we hit a valid index
+  itr->parfum_index++;
+
+  while((!itr->model->elem[0].is_valid_any_idx(itr->parfum_index)) &&
+               (itr->parfum_index < itr->model->elem[0].size()))
+       itr->parfum_index++;
+
+  
+  if(itr->parfum_index==itr->model->elem[0].size()){
+       itr->parfum_index = itr->model->elem[0].size()+1000; // way past the end
+  }  
+
+#ifdef PTOPS_ITERATOR_PRINT
+  CkPrintf("Advancing Elem Iterator to %d\n", itr->parfum_index);
+#endif
 }
 
+
 TopNode topElemItr_GetCurr(TopElemItr*itr){
        return itr->parfum_index;
 }