--- trunk/src/applications/hydrodynamics/RoughShell.cpp	2006/03/15 15:51:44	898
+++ trunk/src/applications/hydrodynamics/RoughShell.cpp	2006/05/24 18:31:12	972
@@ -40,81 +40,93 @@
  */
 
 #include "applications/hydrodynamics/RoughShell.hpp"
-
+#include "applications/hydrodynamics/ShapeBuilder.hpp"
+#include "brains/SimInfo.hpp"
 namespace oopse {
-
-RoughShell::RoughShell(StuntDouble* sd, const DynamicProperty& extraParams) : HydrodynamicsModel(sd, extraParams), sdShape_(sd){
-    DynamicProperty::const_iterator iter = extraParams.find("Sigma");
-    if (iter != extraParams.end()) {
-        boost::any param = iter->second;
-        sigma_ = boost::any_cast<double>(param);
+  
+  RoughShell::RoughShell(StuntDouble* sd, SimInfo* info) : ApproximationModel(sd, info){
+    shape_=ShapeBuilder::createShape(sd);
+    Globals* simParams = info->getSimParams();
+    if (simParams->haveBeadSize()) {
+      sigma_ = simParams->getBeadSize();
     }else {
-        std::cout << "RoughShell Model Error\n" ;
+      
     }
-}
-
-struct BeadLattice {
+  }
+  
+  struct BeadLattice {
     Vector3d origin;
-    double radius;
+    RealType radius;
     bool interior;
-};
-
-struct ExteriorFunctor : public std::unary_function<BeadLattice, bool>{
-
+  };
+  
+  struct ExteriorFunctor : public std::unary_function<BeadLattice, bool>{
+    
     bool operator() (const BeadLattice& bead) {
-        return !bead.interior;
+      return !bead.interior;
     }
-
-};
-
-struct InteriorFunctor  : public std::unary_function<BeadLattice, bool>{
-
+    
+  };
+  
+  struct InteriorFunctor  : public std::unary_function<BeadLattice, bool>{
+    
     bool operator() (const BeadLattice& bead) {
-        return bead.interior;
+      return bead.interior;
     }
-
-};
-bool RoughShell::createBeads(std::vector<BeadParam>& beads) {
-    std::pair<Vector3d, Vector3d> boxBoundary = sdShape_.getBox();
-    double len = boxBoundary.second[0] - boxBoundary.first[0];
-    int numLattices = static_cast<int>(len/sigma_) + 1;
+    
+  };
+  bool RoughShell::createBeads(std::vector<BeadParam>& beads) {
+    std::pair<Vector3d, Vector3d> boxBoundary = shape_->getBoundingBox();
+    RealType firstMin = std::min(std::min(boxBoundary.first[0], boxBoundary.first[1]), boxBoundary.first[2]);
+    RealType secondMax = std::max(std::max(boxBoundary.second[0], boxBoundary.second[1]), boxBoundary.second[2]);
+    RealType len = secondMax - firstMin;
+    int numLattices = static_cast<int>(len/sigma_) + 2;
     Grid3D<BeadLattice>  grid(numLattices, numLattices, numLattices);
-
+    
     //fill beads
     for (int i = 0; i < numLattices; ++i) {
-        for (int j = 0; j < numLattices; ++j) {
-            for (int k = 0; k < numLattices; ++k) {
-                BeadLattice& currentBead = grid(i, j, k);
-                currentBead.origin = Vector3d(i*sigma_ + boxBoundary.first[0], j *sigma_ + boxBoundary.first[1], k*sigma_+ boxBoundary.first[2]);
-                currentBead.radius = sigma_;
-                currentBead.interior = sdShape_.isInterior(grid(i, j, k).origin);                
-            }
+      for (int j = 0; j < numLattices; ++j) {
+        for (int k = 0; k < numLattices; ++k) {
+          BeadLattice& currentBead = grid(i, j, k);
+          currentBead.origin = Vector3d((i-1)*sigma_ + boxBoundary.first[0], (j-1) *sigma_ + boxBoundary.first[1], (k-1)*sigma_+ boxBoundary.first[2]);
+          currentBead.radius = sigma_;
+          currentBead.interior = shape_->isInterior(grid(i, j, k).origin);                
         }
+      }
     }
-
+    
     //remove embedded beads 
     for (int i = 0; i < numLattices; ++i) {
-        for (int j = 0; j < numLattices; ++j) {
-            for (int k = 0; k < numLattices; ++k) {
-                 std::vector<BeadLattice> neighborCells = grid.getAllNeighbors(i, j, k);
-                 //if one of its neighbor cells is exterior, current cell is on the surface
-
-                 std::vector<BeadLattice>::iterator ei = std::find_if(neighborCells.begin(), neighborCells.end(), ExteriorFunctor());                
-                 std::vector<BeadLattice>::iterator ii = std::find_if(neighborCells.begin(), neighborCells.end(), InteriorFunctor());                
-                 
-                  if (ei != neighborCells.end() && ii != neighborCells.end()) {
-                      BeadParam surfaceBead;
-                      surfaceBead.atomName = "Bead";
-                      surfaceBead.pos = grid(i, j, k).origin;
-                      surfaceBead.radius = grid(i, j, k).radius;
-                      beads.push_back(surfaceBead);
-                  }
-
+      for (int j = 0; j < numLattices; ++j) {
+        for (int k = 0; k < numLattices; ++k) {
+          std::vector<BeadLattice> neighborCells = grid.getAllNeighbors(i, j, k);
+          //if one of its neighbor cells is exterior, current cell is on the surface
+          
+          if (grid(i, j, k).interior){
+            
+            bool allNeighBorIsInterior = true;
+            for (std::vector<BeadLattice>::iterator l = neighborCells.begin(); l != neighborCells.end(); ++l) {
+              if (!l->interior) {
+                allNeighBorIsInterior = false;
+                break;
+              }
             }
+            
+            if (allNeighBorIsInterior)
+              continue;
+            
+            BeadParam surfaceBead;
+            surfaceBead.atomName = "H";
+            surfaceBead.pos = grid(i, j, k).origin;
+            surfaceBead.radius = grid(i, j, k).radius;
+            beads.push_back(surfaceBead);                    
+            
+          }
         }
+      }
     }
-     
+    
     return true;
+  }
+  
 }
-
-}