Index: src/boundary.cpp
===================================================================
--- src/boundary.cpp	(revision 5c7bf82e2f560799e5659f7449a8e27e32e57f91)
+++ src/boundary.cpp	(revision 08ef350a55813b0a4db675d13b4659553f889216)
@@ -476,5 +476,6 @@
   for (int axis = 0; axis < NDIM; axis++)
     for (Boundaries::iterator runner = BoundaryPoints[axis].begin(); runner != BoundaryPoints[axis].end(); runner++)
-        mol->TesselStruct->AddPoint(runner->second.second);
+        if (!mol->TesselStruct->AddPoint(runner->second.second, 0))
+          *out << Verbose(3) << "WARNING: Point " << *(runner->second.second) << " is already present!" << endl;
 
   *out << Verbose(2) << "I found " << mol->TesselStruct->PointsOnBoundaryCount << " points on the convex boundary." << endl;
@@ -560,4 +561,123 @@
 };
 
+/** Creates a convex envelope from a given non-convex one.
+ * -# We go through all PointsOnBoundary.
+ * -# We look at each of its lines and CheckConvexityCriterion().
+ * -# If any returns concave, this Point cannot be on the final convex envelope.
+ * -# Hence, we remove it and re-create all its triangles from its getCircleOfConnectedPoints()
+ * -# We calculate the additional volume
+ * -# We go over the points until none yields a concave spot anymore.
+ * 
+ * Note: This routine - for free - calculates the difference in volume between convex and
+ * non-convex envelope, as the former is easy to calculate - VolumeOfConvexEnvelope() - it
+ * can be used to compute volumes of arbitrary shapes.
+ * \param *out output stream for debugging
+ * \param *TesselStruct non-convex envelope, is changed in return!
+ * \param *configuration contains IsAngstroem()
+ * \return volume difference between the non- and the created convex envelope
+ */
+double ConvexizeNonconvexEnvelope(ofstream *out, class Tesselation *TesselStruct)
+{
+  double volume = 0;
+  class BoundaryPointSet *point = NULL;
+  class BoundaryLineSet *line = NULL;
+  class BoundaryTriangleSet *triangle = NULL;
+  Vector OldPoint, TetraederVector[3];
+  *out << Verbose(0) << "Begin of ConvexizeNonconvexEnvelope" << endl;
+
+  if (TesselStruct == NULL) {
+    *out << Verbose(1) << "ERROR: TesselStruct is empty!" << endl;
+    return volume;
+  }
+
+  bool AllConvex = true;
+  bool Convexity = false;
+  do {
+    AllConvex = true;
+    for (PointMap::iterator PointRunner = TesselStruct->PointsOnBoundary.begin(); PointRunner != TesselStruct->PointsOnBoundary.end(); PointRunner++) {
+      point = PointRunner->second;
+      *out << Verbose(1) << "INFO: Current point is " << *point << "." << endl;
+      Convexity = true;
+      for (LineMap::iterator LineRunner = point->lines.begin(); LineRunner != point->lines.end(); LineRunner++) {
+        line = LineRunner->second;
+        *out << Verbose(2) << "INFO: Current line is " << *line << "." << endl;
+        if (line->CheckConvexityCriterion(out)) { 
+          // convex
+        } else { // concave
+          Convexity = false; // we have to go again through all ...
+          break;  // no need to check the others too
+        }
+      }
+      AllConvex = AllConvex && Convexity;
+      if (!Convexity) {
+        *out << Verbose(1) << "... point cannot be on convex envelope." << endl;
+        OldPoint.CopyVector(point->node->node);
+        // get list of connected points
+        BoundaryPointSet *Otherpoint = line->GetOtherEndpoint(point);
+        list<TesselPoint*> *CircleofPoints = TesselStruct->getCircleOfConnectedPoints(out, point->node, Otherpoint->node->node);
+        // remove all triangles
+        int *numbers = NULL;
+        int count = 0;
+        int i=0;
+        for (LineMap::iterator LineRunner = point->lines.begin(); LineRunner != point->lines.end(); LineRunner++)
+          count+=LineRunner->second->triangles.size();
+        numbers = new int[count];
+        for (LineMap::iterator LineRunner = point->lines.begin(); LineRunner != point->lines.end(); LineRunner++) {
+          line = LineRunner->second;
+          for (TriangleMap::iterator TriangleRunner = line->triangles.begin(); TriangleRunner != line->triangles.end(); TriangleRunner++) {
+            triangle = TriangleRunner->second;
+            *out << Verbose(2) << "Erasing triangle " << *triangle << "." << endl;
+            numbers[i] = triangle->Nr;
+            TesselStruct->TrianglesOnBoundary.erase(triangle->Nr);
+            delete(triangle);
+            i++;
+          }
+        }
+        *out << Verbose(1) << i << " triangles were removed." << endl;
+
+        // re-create all triangles by going through connected points list
+        list<TesselPoint*>::iterator OtherCircleRunner = CircleofPoints->end();
+        OtherCircleRunner--;
+        i=0;
+        for (list<TesselPoint*>::iterator CircleRunner = CircleofPoints->begin(); CircleRunner != CircleofPoints->end(); CircleRunner++) {
+          *out << Verbose(2) << "Adding new triangle points."<< endl;
+          TesselStruct->AddPoint(point->node, 0);
+          TesselStruct->AddPoint(*CircleRunner, 1);
+          TesselStruct->AddPoint(*OtherCircleRunner, 2);
+          *out << Verbose(2) << "Adding new triangle lines."<< endl;
+          TesselStruct->AddTriangleLine(TesselStruct->BPS[0], TesselStruct->BPS[1], 0); 
+          TesselStruct->AddTriangleLine(TesselStruct->BPS[0], TesselStruct->BPS[2], 1); 
+          TesselStruct->AddTriangleLine(TesselStruct->BPS[1], TesselStruct->BPS[2], 2); 
+          TesselStruct->BTS = new class BoundaryTriangleSet(TesselStruct->BLS, numbers[i]);
+          TesselStruct->TrianglesOnBoundary.insert(TrianglePair(numbers[i], triangle));
+          *out << Verbose(2) << "Created triangle " << *TesselStruct->BTS << "." << endl;
+          // calculate volume summand as a general tetraeder
+          for (int j=0;j<3;j++) {
+            TetraederVector[j].CopyVector(TesselStruct->BPS[j]->node->node);
+            TetraederVector[j].SubtractVector(&OldPoint);
+          }
+          OldPoint.CopyVector(&TetraederVector[0]);
+          OldPoint.VectorProduct(&TetraederVector[1]);
+          volume += 1./6. * fabs(OldPoint.ScalarProduct(&TetraederVector[2]));
+          // advance other shank
+          OtherCircleRunner++;
+          if (OtherCircleRunner == CircleofPoints->end())
+            OtherCircleRunner = CircleofPoints->begin();
+          // advance number
+          i++;
+          if (i > count)
+            *out << Verbose(2) << "WARNING: Maximum of numbers reached!" << endl;
+        }
+        *out << Verbose(1) << i << " triangles were created." << endl;
+        
+        delete[](numbers);
+      }
+    }
+  } while (!AllConvex);
+  *out << Verbose(0) << "End of ConvexizeNonconvexEnvelope" << endl;
+
+  // end
+  return volume;
+};
 
 /** Determines the volume of a cluster.
Index: src/boundary.hpp
===================================================================
--- src/boundary.hpp	(revision 5c7bf82e2f560799e5659f7449a8e27e32e57f91)
+++ src/boundary.hpp	(revision 08ef350a55813b0a4db675d13b4659553f889216)
@@ -39,4 +39,5 @@
 void Find_convex_border(ofstream *out, molecule* mol, class LinkedCell *LCList, const char *filename);
 void Find_non_convex_border(ofstream *out, molecule* mol, class LinkedCell *LC, const char *tempbasename, const double RADIUS);
+double ConvexizeNonconvexEnvelope(ofstream *out, class Tesselation *TesselStruct);
 void Find_next_suitable_point(class BoundaryTriangleSet *BaseTriangle, class BoundaryLineSet *BaseLine, atom*& OptCandidate, Vector *OptCandidateCenter, double *ShortestAngle, const double RADIUS, LinkedCell *LC);
 Boundaries *GetBoundaryPoints(ofstream *out, molecule *mol);
Index: src/builder.cpp
===================================================================
--- src/builder.cpp	(revision 5c7bf82e2f560799e5659f7449a8e27e32e57f91)
+++ src/builder.cpp	(revision 08ef350a55813b0a4db675d13b4659553f889216)
@@ -1871,9 +1871,11 @@
                 cout << Verbose(1) << "Storing tecplot data in " << argv[argptr] << "." << endl;
                 LinkedCell LCList(mol, 10.);
-                class Tesselation *TesselStruct = NULL;
-                Find_convex_border((ofstream *)&cout, mol, &LCList, argv[argptr]);
-                double clustervolume = VolumeOfConvexEnvelope((ofstream *)&cout, TesselStruct, &configuration);
-                cout << Verbose(0) << "The tesselated surface area is " << clustervolume << "." << endl;
-                delete(TesselStruct);
+                //Find_convex_border((ofstream *)&cout, mol, &LCList, argv[argptr]);
+                Find_non_convex_border((ofstream *)&cout, mol, &LCList, argv[argptr], 10.);
+
+		double volumedifference = ConvexizeNonconvexEnvelope((ofstream *)&cout, mol->TesselStruct);
+                double clustervolume = VolumeOfConvexEnvelope((ofstream *)&cout, mol->TesselStruct, &configuration);
+                cout << Verbose(0) << "The tesselated surface area is " << clustervolume << " " << (configuration.GetIsAngstroem() ? "angstrom" : "atomiclength") << "^3." << endl;
+                cout << Verbose(0) << "The non-convex tesselated surface area is " << clustervolume-volumedifference << " " << (configuration.GetIsAngstroem() ? "angstrom" : "atomiclength") << "^3." << endl;
                 argptr+=1;
               }
Index: src/tesselation.cpp
===================================================================
--- src/tesselation.cpp	(revision 5c7bf82e2f560799e5659f7449a8e27e32e57f91)
+++ src/tesselation.cpp	(revision 08ef350a55813b0a4db675d13b4659553f889216)
@@ -208,5 +208,5 @@
  * \return NULL - if endpoint not contained in BoundaryLineSet, or pointer to BoundaryPointSet otherwise
  */
-inline class BoundaryPointSet *BoundaryLineSet::GetOtherEndpoint(class BoundaryPointSet *point)
+class BoundaryPointSet *BoundaryLineSet::GetOtherEndpoint(class BoundaryPointSet *point)
 {
   if (endpoints[0] == point)
@@ -1044,6 +1044,5 @@
         // add Walker to boundary points
         *out << Verbose(2) << "Adding " << *Walker << " to BoundaryPoints." << endl;
-        AddPoint(Walker);
-        if (BPS[0] != NULL)
+        if (AddPoint(Walker,0))
           NewPoint = BPS[0];
         else
@@ -1100,16 +1099,20 @@
 /** Adds an point to the tesselation::PointsOnBoundary list.
  * \param *Walker point to add
- */
-void
-Tesselation::AddPoint(TesselPoint *Walker)
+ * \param n TesselStruct::BPS index to put pointer into
+ * \return true - new point was added, false - point already present
+ */
+bool
+Tesselation::AddPoint(TesselPoint *Walker, int n)
 {
   PointTestPair InsertUnique;
-  BPS[0] = new class BoundaryPointSet(Walker);
-  InsertUnique = PointsOnBoundary.insert(PointPair(Walker->nr, BPS[0]));
-  if (InsertUnique.second) // if new point was not present before, increase counter
+  BPS[n] = new class BoundaryPointSet(Walker);
+  InsertUnique = PointsOnBoundary.insert(PointPair(Walker->nr, BPS[n]));
+  if (InsertUnique.second) { // if new point was not present before, increase counter
     PointsOnBoundaryCount++;
-  else {
-    delete(BPS[0]);
-    BPS[0] = NULL;
+    return true;
+  } else {
+    delete(BPS[n]);
+    BPS[n] = InsertUnique.first->second;
+    return false;
   }
 }
Index: src/tesselation.hpp
===================================================================
--- src/tesselation.hpp	(revision 5c7bf82e2f560799e5659f7449a8e27e32e57f91)
+++ src/tesselation.hpp	(revision 08ef350a55813b0a4db675d13b4659553f889216)
@@ -190,5 +190,5 @@
     virtual ~Tesselation();
 
-    void AddPoint(TesselPoint *Walker);
+    bool AddPoint(TesselPoint *Walker, int n);
     void AddTrianglePoint(TesselPoint* Candidate, int n);
     void AddTriangleLine(class BoundaryPointSet *a, class BoundaryPointSet *b, int n);
@@ -236,9 +236,10 @@
     virtual bool IsEnd();
 
-  private:
-    class BoundaryPointSet *TPS[3]; //this is a Storage for pointers to triangle points, this and BPS[2] needed due to AddLine restrictions
     class BoundaryPointSet *BPS[2];
     class BoundaryLineSet *BLS[3];
     class BoundaryTriangleSet *BTS;
+
+  private:
+    class BoundaryPointSet *TPS[3]; //this is a Storage for pointers to triangle points, this and BPS[2] needed due to AddLine restrictions
 
     PointMap::iterator InternalPointer;
