Changeset ad37ab

Timestamp:

Oct 27, 2009, 11:51:52 AM (15 years ago)

Author:

Frederik Heber <heber@…>

Branches:

Action_Thermostats, Add_AtomRandomPerturbation, Add_FitFragmentPartialChargesAction, Add_RotateAroundBondAction, Add_SelectAtomByNameAction, Added_ParseSaveFragmentResults, AddingActions_SaveParseParticleParameters, Adding_Graph_to_ChangeBondActions, Adding_MD_integration_tests, Adding_ParticleName_to_Atom, Adding_StructOpt_integration_tests, AtomFragments, Automaking_mpqc_open, AutomationFragmentation_failures, Candidate_v1.5.4, Candidate_v1.6.0, Candidate_v1.6.1, ChangeBugEmailaddress, ChangingTestPorts, ChemicalSpaceEvaluator, CombiningParticlePotentialParsing, Combining_Subpackages, Debian_Package_split, Debian_package_split_molecuildergui_only, Disabling_MemDebug, Docu_Python_wait, EmpiricalPotential_contain_HomologyGraph, EmpiricalPotential_contain_HomologyGraph_documentation, Enable_parallel_make_install, Enhance_userguide, Enhanced_StructuralOptimization, Enhanced_StructuralOptimization_continued, Example_ManyWaysToTranslateAtom, Exclude_Hydrogens_annealWithBondGraph, FitPartialCharges_GlobalError, Fix_BoundInBox_CenterInBox_MoleculeActions, Fix_ChargeSampling_PBC, Fix_ChronosMutex, Fix_FitPartialCharges, Fix_FitPotential_needs_atomicnumbers, Fix_ForceAnnealing, Fix_IndependentFragmentGrids, Fix_ParseParticles, Fix_ParseParticles_split_forward_backward_Actions, Fix_PopActions, Fix_QtFragmentList_sorted_selection, Fix_Restrictedkeyset_FragmentMolecule, Fix_StatusMsg, Fix_StepWorldTime_single_argument, Fix_Verbose_Codepatterns, Fix_fitting_potentials, Fixes, ForceAnnealing_goodresults, ForceAnnealing_oldresults, ForceAnnealing_tocheck, ForceAnnealing_with_BondGraph, ForceAnnealing_with_BondGraph_continued, ForceAnnealing_with_BondGraph_continued_betteresults, ForceAnnealing_with_BondGraph_contraction-expansion, FragmentAction_writes_AtomFragments, FragmentMolecule_checks_bonddegrees, GeometryObjects, Gui_Fixes, Gui_displays_atomic_force_velocity, ImplicitCharges, IndependentFragmentGrids, IndependentFragmentGrids_IndividualZeroInstances, IndependentFragmentGrids_IntegrationTest, IndependentFragmentGrids_Sole_NN_Calculation, JobMarket_RobustOnKillsSegFaults, JobMarket_StableWorkerPool, JobMarket_unresolvable_hostname_fix, MoreRobust_FragmentAutomation, ODR_violation_mpqc_open, PartialCharges_OrthogonalSummation, PdbParser_setsAtomName, PythonUI_with_named_parameters, QtGui_reactivate_TimeChanged_changes, Recreated_GuiChecks, Rewrite_FitPartialCharges, RotateToPrincipalAxisSystem_UndoRedo, SaturateAtoms_findBestMatching, SaturateAtoms_singleDegree, StoppableMakroAction, Subpackage_CodePatterns, Subpackage_JobMarket, Subpackage_LinearAlgebra, Subpackage_levmar, Subpackage_mpqc_open, Subpackage_vmg, Switchable_LogView, ThirdParty_MPQC_rebuilt_buildsystem, TrajectoryDependenant_MaxOrder, TremoloParser_IncreasedPrecision, TremoloParser_MultipleTimesteps, TremoloParser_setsAtomName, Ubuntu_1604_changes, stable

Children:

776b64

Parents:

70ff32

Message:

analyzer.cpp and boundary.cpp refactored.

• ticket #39 is done with respect to these files
• both do not contain member functions (anymore), hence ticket #4 and ticket #38 are both fulfilled.

Location:

src

Files:

• analyzer.cpp (modified) (1 diff)
• boundary.cpp (modified) (19 diffs)
• boundary.hpp (modified) (1 diff)
• molecule.hpp (modified) (1 diff)

src/analyzer.cpp

@@ -59 +59 @@
   bool NoTime = false;
   bool NoHCorrection = true;
-  int counter;
+  int counter = 0;
 
   cout << "ANOVA Analyzer" << endl;

src/boundary.cpp

@@ -34 +34 @@
   bool BoundaryFreeFlag = false;
   Boundaries *BoundaryPoints = BoundaryPtr;
+  double OldComponent = 0.;
+  double tmp = 0.;
+  double w1 = 0.;
+  double w2 = 0.;
+  Vector DistanceVector;
+  Vector OtherVector;
+  int component = 0;
+  int Othercomponent = 0;
+  Boundaries::iterator Neighbour;
+  Boundaries::iterator OtherNeighbour;
+  double *GreatestDiameter = new double[NDIM];
+
   if (BoundaryPoints == NULL) {
     BoundaryFreeFlag = true;
@@ -41 +53 @@
   }
   // determine biggest "diameter" of cluster for each axis
-  Boundaries::iterator Neighbour, OtherNeighbour;
-  double *GreatestDiameter = new double[NDIM];
   for (int i = 0; i < NDIM; i++)
     GreatestDiameter[i] = 0.;
-  double OldComponent, tmp, w1, w2;
-  Vector DistanceVector, OtherVector;
-  int component, Othercomponent;
   for (int axis = 0; axis < NDIM; axis++)
     { // regard each projected plane
@@ -135 +142 @@
   Vector *MolCenter = mol->DetermineCenterOfAll(out);
   Vector helper;
+  BoundariesTestPair BoundaryTestPair;
+  Vector AxisVector;
+  Vector AngleReferenceVector;
+  Vector AngleReferenceNormalVector;
+  Vector ProjectedVector;
+  Boundaries *BoundaryPoints = new Boundaries[NDIM]; // first is alpha, second is (r, nr)
+  double radius = 0.;
+  double angle = 0.;
 
   *out << Verbose(1) << "Finding all boundary points." << endl;
-  Boundaries *BoundaryPoints = new Boundaries[NDIM]; // first is alpha, second is (r, nr)
-  BoundariesTestPair BoundaryTestPair;
-  Vector AxisVector, AngleReferenceVector, AngleReferenceNormalVector;
-  double radius, angle;
   // 3a. Go through every axis
   for (int axis = 0; axis < NDIM; axis++) {
@@ -156 +167 @@
     while (Walker->next != mol->end) {
       Walker = Walker->next;
-      Vector ProjectedVector;
       ProjectedVector.CopyVector(&Walker->x);
       ProjectedVector.SubtractVector(MolCenter);
@@ -267 +277 @@
 
           // calculate each length
-          double a = SideA.Norm();
-          //double b = SideB.Norm();
-          //double c = SideC.Norm();
-          double h = SideH.Norm();
+          const double a = SideA.Norm();
+          //const double b = SideB.Norm();
+          //const double c = SideC.Norm();
+          const double h = SideH.Norm();
           // calculate the angles
-          double alpha = SideA.Angle(&SideH);
-          double beta = SideA.Angle(&SideC);
-          double gamma = SideB.Angle(&SideH);
-          double delta = SideC.Angle(&SideH);
-          double MinDistance = a * sin(beta) / (sin(delta)) * (((alpha < M_PI / 2.) || (gamma < M_PI / 2.)) ? 1. : -1.);
+          const double alpha = SideA.Angle(&SideH);
+          const double beta = SideA.Angle(&SideC);
+          const double gamma = SideB.Angle(&SideH);
+          const double delta = SideC.Angle(&SideH);
+          const double MinDistance = a * sin(beta) / (sin(delta)) * (((alpha < M_PI / 2.) || (gamma < M_PI / 2.)) ? 1. : -1.);
           //*out << Verbose(2) << " I calculated: a = " << a << ", h = " << h << ", beta(" << left->second.second->Name << "," << left->second.second->Name << "-" << right->second.second->Name << ") = " << beta << ", delta(" << left->second.second->Name << "," << runner->second.second->Name << ") = " << delta << ", Min = " << MinDistance << "." << endl;
           *out << Verbose(1) << "Checking CoG distance of runner " << *runner->second.second << " " << h << " against triangle's side length spanned by (" << *left->second.second << "," << *right->second.second << ") of " << MinDistance << "." << endl;
@@ -386 +396 @@
 
   // 3d. check all baselines whether the peaks of the two adjacent triangles with respect to center of baseline are convex, if not, make the baseline between the two peaks and baseline endpoints become the new peaks
-  bool AllConvex;
+  bool AllConvex = true;
   class BoundaryLineSet *line = NULL;
   do {
@@ -506 +516 @@
   class BoundaryPointSet *point = NULL;
   class BoundaryLineSet *line = NULL;
-  bool Concavity;
+  bool Concavity = false;
   char dummy[MAXSTRINGSIZE];
-  PointMap::iterator PointRunner, PointAdvance;
-  LineMap::iterator LineRunner, LineAdvance;
-  TriangleMap::iterator TriangleRunner, TriangleAdvance;
+  PointMap::iterator PointRunner;
+  PointMap::iterator PointAdvance;
+  LineMap::iterator LineRunner;
+  LineMap::iterator LineAdvance;
+  TriangleMap::iterator TriangleRunner;
+  TriangleMap::iterator TriangleAdvance;
+  int run = 0;
 
   *out << Verbose(0) << "Begin of ConvexizeNonconvexEnvelope" << endl;
@@ -521 +535 @@
 
   // First step: RemovePointFromTesselatedSurface
-  int run = 0;
-  double tmp;
   do {
     Concavity = false;
@@ -564 +576 @@
       // take highest of both lines
       if (TesselStruct->IsConvexRectangle(out, line) == NULL) {
-        tmp = TesselStruct->PickFarthestofTwoBaselines(out, line);
+        const double tmp = TesselStruct->PickFarthestofTwoBaselines(out, line);
         volume += tmp;
-        if (tmp != 0) {
+        if (tmp != 0.) {
           mol->TesselStruct->FlipBaseline(out, line);
           Concavity = true;
@@ -619 +631 @@
   bool IsAngstroem = configuration->GetIsAngstroem();
   double volume = 0.;
-  double PyramidVolume = 0.;
-  double G, h;
-  Vector x, y;
-  double a, b, c;
+  Vector x;
+  Vector y;
 
   // 6a. Every triangle forms a pyramid with the center of gravity as its peak, sum up the volumes
@@ -634 +644 @@
       y.CopyVector(runner->second->endpoints[0]->node->node);
       y.SubtractVector(runner->second->endpoints[2]->node->node);
-      a = sqrt(runner->second->endpoints[0]->node->node->DistanceSquared(runner->second->endpoints[1]->node->node));
-      b = sqrt(runner->second->endpoints[0]->node->node->DistanceSquared(runner->second->endpoints[2]->node->node));
-      c = sqrt(runner->second->endpoints[2]->node->node->DistanceSquared(runner->second->endpoints[1]->node->node));
-      G = sqrt(((a + b + c) * (a + b + c) - 2 * (a * a + b * b + c * c)) / 16.); // area of tesselated triangle
+      const double a = sqrt(runner->second->endpoints[0]->node->node->DistanceSquared(runner->second->endpoints[1]->node->node));
+      const double b = sqrt(runner->second->endpoints[0]->node->node->DistanceSquared(runner->second->endpoints[2]->node->node));
+      const double c = sqrt(runner->second->endpoints[2]->node->node->DistanceSquared(runner->second->endpoints[1]->node->node));
+      const double G = sqrt(((a + b + c) * (a + b + c) - 2 * (a * a + b * b + c * c)) / 16.); // area of tesselated triangle
       x.MakeNormalVector(runner->second->endpoints[0]->node->node, runner->second->endpoints[1]->node->node, runner->second->endpoints[2]->node->node);
       x.Scale(runner->second->endpoints[1]->node->node->ScalarProduct(&x));
-      h = x.Norm(); // distance of CoG to triangle
-      PyramidVolume = (1. / 3.) * G * h; // this formula holds for _all_ pyramids (independent of n-edge base or (not) centered peak)
+      const double h = x.Norm(); // distance of CoG to triangle
+      const double PyramidVolume = (1. / 3.) * G * h; // this formula holds for _all_ pyramids (independent of n-edge base or (not) centered peak)
       *out << Verbose(2) << "Area of triangle is " << G << " "
           << (IsAngstroem ? "angstrom" : "atomiclength") << "^2, height is "
@@ -694 +704 @@
  * \param celldensity desired average density in final cell
  */
-void
-PrepareClustersinWater(ofstream *out, config *configuration, molecule *mol,
-    double ClusterVolume, double celldensity)
+void PrepareClustersinWater(ofstream *out, config *configuration, molecule *mol, double ClusterVolume, double celldensity)
 {
+  bool IsAngstroem = NULL;
+  double *GreatestDiameter = NULL;
+  Boundaries *BoundaryPoints = NULL;
+  class Tesselation *TesselStruct = NULL;
+  Vector BoxLengths;
+  int repetition[NDIM] = { 1, 1, 1 };
+  int TotalNoClusters = 1;
+  atom *Walker = NULL;
+  double totalmass = 0.;
+  double clustervolume = 0.;
+  double cellvolume = 0.;
+
   // transform to PAS
   mol->PrincipalAxisSystem(out, true);
 
-  // some preparations beforehand
-  bool IsAngstroem = configuration->GetIsAngstroem();
-  Boundaries *BoundaryPoints = GetBoundaryPoints(out, mol);
-  class Tesselation *TesselStruct = NULL;
+  IsAngstroem = configuration->GetIsAngstroem();
+  GreatestDiameter = GetDiametersOfCluster(out, BoundaryPoints, mol, IsAngstroem);
+  BoundaryPoints = GetBoundaryPoints(out, mol);
   LinkedCell LCList(mol, 10.);
   FindConvexBorder(out, mol, &LCList, NULL);
-  double clustervolume;
+
+  // some preparations beforehand
   if (ClusterVolume == 0)
     clustervolume = VolumeOfConvexEnvelope(out, TesselStruct, configuration);
   else
     clustervolume = ClusterVolume;
-  double *GreatestDiameter = GetDiametersOfCluster(out, BoundaryPoints, mol, IsAngstroem);
-  Vector BoxLengths;
-  int repetition[NDIM] =
-    { 1, 1, 1 };
-  int TotalNoClusters = 1;
+
   for (int i = 0; i < NDIM; i++)
     TotalNoClusters *= repetition[i];
 
   // sum up the atomic masses
-  double totalmass = 0.;
-  atom *Walker = mol->start;
-  while (Walker->next != mol->end)
-    {
+  Walker = mol->start;
+  while (Walker->next != mol->end) {
       Walker = Walker->next;
       totalmass += Walker->type->mass;
-    }
-  *out << Verbose(0) << "RESULT: The summed mass is " << setprecision(10)
-      << totalmass << " atomicmassunit." << endl;
-
-  *out << Verbose(0) << "RESULT: The average density is " << setprecision(10)
-      << totalmass / clustervolume << " atomicmassunit/"
-      << (IsAngstroem ? "angstrom" : "atomiclength") << "^3." << endl;
+  }
+  *out << Verbose(0) << "RESULT: The summed mass is " << setprecision(10) << totalmass << " atomicmassunit." << endl;
+  *out << Verbose(0) << "RESULT: The average density is " << setprecision(10) << totalmass / clustervolume << " atomicmassunit/" << (IsAngstroem ? "angstrom" : "atomiclength") << "^3." << endl;
 
   // solve cubic polynomial
-  *out << Verbose(1) << "Solving equidistant suspension in water problem ..."
-      << endl;
-  double cellvolume;
+  *out << Verbose(1) << "Solving equidistant suspension in water problem ..." << endl;
   if (IsAngstroem)
-    cellvolume = (TotalNoClusters * totalmass / SOLVENTDENSITY_A - (totalmass
-        / clustervolume)) / (celldensity - 1);
+    cellvolume = (TotalNoClusters * totalmass / SOLVENTDENSITY_A - (totalmass / clustervolume)) / (celldensity - 1);
   else
-    cellvolume = (TotalNoClusters * totalmass / SOLVENTDENSITY_a0 - (totalmass
-        / clustervolume)) / (celldensity - 1);
-  *out << Verbose(1) << "Cellvolume needed for a density of " << celldensity
-      << " g/cm^3 is " << cellvolume << " " << (IsAngstroem ? "angstrom"
-      : "atomiclength") << "^3." << endl;
-
-  double minimumvolume = TotalNoClusters * (GreatestDiameter[0]
-      * GreatestDiameter[1] * GreatestDiameter[2]);
-  *out << Verbose(1)
-      << "Minimum volume of the convex envelope contained in a rectangular box is "
-      << minimumvolume << " atomicmassunit/" << (IsAngstroem ? "angstrom"
-      : "atomiclength") << "^3." << endl;
-  if (minimumvolume > cellvolume)
-    {
-      cerr << Verbose(0)
-          << "ERROR: the containing box already has a greater volume than the envisaged cell volume!"
-          << endl;
-      cout << Verbose(0)
-          << "Setting Box dimensions to minimum possible, the greatest diameters."
-          << endl;
-      for (int i = 0; i < NDIM; i++)
-        BoxLengths.x[i] = GreatestDiameter[i];
-      mol->CenterEdge(out, &BoxLengths);
-    }
-  else
-    {
-      BoxLengths.x[0] = (repetition[0] * GreatestDiameter[0] + repetition[1]
-          * GreatestDiameter[1] + repetition[2] * GreatestDiameter[2]);
-      BoxLengths.x[1] = (repetition[0] * repetition[1] * GreatestDiameter[0]
-          * GreatestDiameter[1] + repetition[0] * repetition[2]
-          * GreatestDiameter[0] * GreatestDiameter[2] + repetition[1]
-          * repetition[2] * GreatestDiameter[1] * GreatestDiameter[2]);
-      BoxLengths.x[2] = minimumvolume - cellvolume;
-      double x0 = 0., x1 = 0., x2 = 0.;
-      if (gsl_poly_solve_cubic(BoxLengths.x[0], BoxLengths.x[1],
-          BoxLengths.x[2], &x0, &x1, &x2) == 1) // either 1 or 3 on return
-        *out << Verbose(0) << "RESULT: The resulting spacing is: " << x0
-            << " ." << endl;
-      else
-        {
-          *out << Verbose(0) << "RESULT: The resulting spacings are: " << x0
-              << " and " << x1 << " and " << x2 << " ." << endl;
-          x0 = x2; // sorted in ascending order
-        }
-
-      cellvolume = 1;
-      for (int i = 0; i < NDIM; i++)
-        {
-          BoxLengths.x[i] = repetition[i] * (x0 + GreatestDiameter[i]);
-          cellvolume *= BoxLengths.x[i];
-        }
-
-      // set new box dimensions
-      *out << Verbose(0) << "Translating to box with these boundaries." << endl;
-      mol->SetBoxDimension(&BoxLengths);
-      mol->CenterInBox((ofstream *) &cout);
-    }
+    cellvolume = (TotalNoClusters * totalmass / SOLVENTDENSITY_a0 - (totalmass / clustervolume)) / (celldensity - 1);
+  *out << Verbose(1) << "Cellvolume needed for a density of " << celldensity << " g/cm^3 is " << cellvolume << " " << (IsAngstroem ? "angstrom" : "atomiclength") << "^3." << endl;
+
+  double minimumvolume = TotalNoClusters * (GreatestDiameter[0] * GreatestDiameter[1] * GreatestDiameter[2]);
+  *out << Verbose(1) << "Minimum volume of the convex envelope contained in a rectangular box is " << minimumvolume << " atomicmassunit/" << (IsAngstroem ? "angstrom" : "atomiclength") << "^3." << endl;
+  if (minimumvolume > cellvolume) {
+    cerr << Verbose(0) << "ERROR: the containing box already has a greater volume than the envisaged cell volume!" << endl;
+    cout << Verbose(0) << "Setting Box dimensions to minimum possible, the greatest diameters." << endl;
+    for (int i = 0; i < NDIM; i++)
+      BoxLengths.x[i] = GreatestDiameter[i];
+    mol->CenterEdge(out, &BoxLengths);
+  } else {
+    BoxLengths.x[0] = (repetition[0] * GreatestDiameter[0] + repetition[1] * GreatestDiameter[1] + repetition[2] * GreatestDiameter[2]);
+    BoxLengths.x[1] = (repetition[0] * repetition[1] * GreatestDiameter[0] * GreatestDiameter[1] + repetition[0] * repetition[2] * GreatestDiameter[0] * GreatestDiameter[2] + repetition[1] * repetition[2] * GreatestDiameter[1] * GreatestDiameter[2]);
+    BoxLengths.x[2] = minimumvolume - cellvolume;
+    double x0 = 0.;
+    double x1 = 0.;
+    double x2 = 0.;
+    if (gsl_poly_solve_cubic(BoxLengths.x[0], BoxLengths.x[1], BoxLengths.x[2], &x0, &x1, &x2) == 1) // either 1 or 3 on return
+      *out << Verbose(0) << "RESULT: The resulting spacing is: " << x0 << " ." << endl;
+    else {
+      *out << Verbose(0) << "RESULT: The resulting spacings are: " << x0 << " and " << x1 << " and " << x2 << " ." << endl;
+      x0 = x2; // sorted in ascending order
+    }
+
+    cellvolume = 1.;
+    for (int i = 0; i < NDIM; i++) {
+      BoxLengths.x[i] = repetition[i] * (x0 + GreatestDiameter[i]);
+      cellvolume *= BoxLengths.x[i];
+    }
+
+    // set new box dimensions
+    *out << Verbose(0) << "Translating to box with these boundaries." << endl;
+    mol->SetBoxDimension(&BoxLengths);
+    mol->CenterInBox((ofstream *) &cout);
+  }
   // update Box of atoms by boundary
   mol->SetBoxDimension(&BoxLengths);
-  *out << Verbose(0) << "RESULT: The resulting cell dimensions are: "
-      << BoxLengths.x[0] << " and " << BoxLengths.x[1] << " and "
-      << BoxLengths.x[2] << " with total volume of " << cellvolume << " "
-      << (IsAngstroem ? "angstrom" : "atomiclength") << "^3." << endl;
-}
-;
+  *out << Verbose(0) << "RESULT: The resulting cell dimensions are: " << BoxLengths.x[0] << " and " << BoxLengths.x[1] << " and " << BoxLengths.x[2] << " with total volume of " << cellvolume << " " << (IsAngstroem ? "angstrom" : "atomiclength") << "^3." << endl;
+};
 
 
@@ -836 +818 @@
   atom *Walker = NULL;
   bond *Binder = NULL;
-  int i;
+  int i = 0;
   LinkedCell *LCList[List->ListOfMolecules.size()];
+  double phi[NDIM];
 
   *out << Verbose(0) << "Begin of FillBoxWithMolecule" << endl;
@@ -858 +841 @@
   filler->CountAtoms(out);
   atom * CopyAtoms[filler->AtomCount];
-  int nr = 0;
 
   // calculate filler grid in [0,1]^3
@@ -903 +885 @@
 
           // go through all atoms
-          nr=0;
           Walker = filler->start;
           while (Walker->next != filler->end) {
@@ -916 +897 @@
             // ... and rotation matrix
             if (DoRandomRotation) {
-              double phi[NDIM];
               for (int i=0;i<NDIM;i++) {
                 phi[i] = rand()/(RAND_MAX/(2.*M_PI));
@@ -975 +955 @@
 {
   bool freeLC = false;
+  LineMap::iterator baseline;
+  LineMap::iterator testline;
+  bool OneLoopWithoutSuccessFlag = false;  // marks whether we went once through all baselines without finding any without two triangles
+  bool TesselationFailFlag = false;
 
   *out << Verbose(1) << "Entering search for non convex hull. " << endl;
@@ -985 +969 @@
     mol->TesselStruct = new Tesselation;
   }
-  LineMap::iterator baseline;
-  LineMap::iterator testline;
+
   *out << Verbose(0) << "Begin of FindNonConvexBorder\n";
-  bool OneLoopWithoutSuccessFlag = false;  // marks whether we went once through all baselines without finding any without two triangles
-  bool TesselationFailFlag = false;
 
   // initialise Linked Cell
@@ -1066 +1047 @@
 
 
-/** Finds a hole of sufficient size in \a this molecule to embed \a *srcmol into it.
+/** Finds a hole of sufficient size in \a *mols to embed \a *srcmol into it.
  * \param *out output stream for debugging
- * \param *srcmol molecule to embed into
+ * \param *mols molecules in the domain to embed in between
+ * \param *srcmol embedding molecule
  * \return *Vector new center of \a *srcmol for embedding relative to \a this
  */
-Vector* molecule::FindEmbeddingHole(ofstream *out, molecule *srcmol)
+Vector* FindEmbeddingHole(ofstream *out, MoleculeListClass *mols, molecule *srcmol)
 {
   Vector *Center = new Vector;

src/boundary.hpp

@@ -59 +59 @@
 void StoreTrianglesinFile(ofstream *out, molecule *mol, const char *filename, const char *extraSuffix);
 bool RemoveAllBoundaryPoints(ofstream *out, class Tesselation *TesselStruct, molecule *mol, char *filename);
+Vector* FindEmbeddingHole(ofstream *out, MoleculeListClass *mols, molecule *srcmol);
 
 

src/molecule.hpp

@@ -228 +228 @@
   void PrincipalAxisSystem(ofstream *out, bool DoRotate);
   double VolumeOfConvexEnvelope(ofstream *out, bool IsAngstroem);
-  Vector* FindEmbeddingHole(ofstream *out, molecule *srcmol);
-
 
   double ConstrainedPotential(ofstream *out, struct EvaluatePotential &Params);