Changeset d30402 for src/boundary.cpp

Timestamp:

Jul 27, 2009, 8:16:03 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:

a20e6a

Parents:

36ec71

Message:

Split VolumeOfConvexEnvelope() into find_convex_border() and remainder.

• Convex envelope is not working anymore, in the state of fixing it and trying to refactor code a bit.
• planned to have VolumeOfConvexEnvelope() to be replaced by VolumeOfEnvelope(), which can also do Non-Convex-Envelopes

File:

• src/boundary.cpp (modified) (8 diffs)

src/boundary.cpp

@@ -760 +760 @@
 }
 
-/** Determines the volume of a cluster.
- * Determines first the convex envelope, then tesselates it and calculates its volume.
+/** Tesselates the convex boundary by finding all boundary points.
  * \param *out output stream for debugging
+ * \param *mol molecule structure with Atom's and Bond's
+ * \param *TesselStruct Tesselation filled with points, lines and triangles on boundary on return
+ * \param *LCList atoms in LinkedCell list
  * \param *filename filename prefix for output of vertex data
- * \param *configuration needed for path to store convex envelope file
- * \param *BoundaryPoints NDIM set of boundary points on the projected plane per axis, on return if desired
- * \param *mol molecule structure representing the cluster
- * \return determined volume of the cluster in cubed config:GetIsAngstroem()
- */
-double
-VolumeOfConvexEnvelope(ofstream *out, const char *filename, config *configuration,
-    Boundaries *BoundaryPtr, molecule *mol)
-{
-  bool IsAngstroem = configuration->GetIsAngstroem();
+ * \return *TesselStruct is filled with convex boundary and tesselation is stored under \a *filename.
+ */
+void Find_convex_border(ofstream *out, molecule* mol, class Tesselation *&TesselStruct, class LinkedCell *LCList, const char *filename)
+{
   atom *Walker = NULL;
-  struct Tesselation *TesselStruct = new Tesselation;
   bool BoundaryFreeFlag = false;
-  Boundaries *BoundaryPoints = BoundaryPtr;
-  double volume = 0.;
-  double PyramidVolume = 0.;
-  double G, h;
-  Vector x, y;
-  double a, b, c;
-
-  //Find_non_convex_border(out, tecplot, *TesselStruct, mol); // Is now called from command line.
+  Boundaries *BoundaryPoints = NULL;
+
+  if (TesselStruct != NULL) // free if allocated
+    delete(TesselStruct);
+  TesselStruct = new class Tesselation;
 
   // 1. calculate center of gravity
@@ -793 +785 @@
   *out << Verbose(1) << "Translating system to Center of Gravity." << endl;
   Walker = mol->start;
-  while (Walker->next != mol->end)
-    {
-      Walker = Walker->next;
-      Walker->x.Translate(CenterOfGravity);
-    }
+  while (Walker->next != mol->end) {
+    Walker = Walker->next;
+    Walker->x.Translate(CenterOfGravity);
+  }
 
   // 3. Find all points on the boundary
-  if (BoundaryPoints == NULL)
-    {
+  if (BoundaryPoints == NULL) {
       BoundaryFreeFlag = true;
       BoundaryPoints = GetBoundaryPoints(out, mol);
-    }
-  else
-    {
+  } else {
       *out << Verbose(1) << "Using given boundary points set." << endl;
-    }
+  }
 
   // 4. fill the boundary point list
   for (int axis = 0; axis < NDIM; axis++)
-    for (Boundaries::iterator runner = BoundaryPoints[axis].begin(); runner
-        != BoundaryPoints[axis].end(); runner++)
-      {
+    for (Boundaries::iterator runner = BoundaryPoints[axis].begin(); runner != BoundaryPoints[axis].end(); runner++)
         TesselStruct->AddPoint(runner->second.second);
-      }
 
   *out << Verbose(2) << "I found " << TesselStruct->PointsOnBoundaryCount
@@ -833 +818 @@
 
   // 5b. go through all lines, that are not yet part of two triangles (only of one so far)
-  TesselStruct->TesselateOnBoundary(out, configuration, mol);
+  TesselStruct->TesselateOnBoundary(out, mol);
 
   *out << Verbose(2) << "I created " << TesselStruct->TrianglesOnBoundaryCount
@@ -839 +824 @@
       << " lines and " << TesselStruct->PointsOnBoundaryCount << " points."
       << endl;
+
+  // 6. translate all points back from CoG
+  *out << Verbose(1) << "Translating system back from Center of Gravity."
+      << endl;
+  CenterOfGravity->Scale(-1);
+  Walker = mol->start;
+  while (Walker->next != mol->end)
+    {
+      Walker = Walker->next;
+      Walker->x.Translate(CenterOfGravity);
+    }
+
+  // 7. Store triangles in tecplot file
+  if (filename != NULL) {
+    string OutputName(filename);
+    OutputName.append(TecplotSuffix);
+    ofstream *tecplot = new ofstream(OutputName.c_str());
+    write_tecplot_file(out, tecplot, TesselStruct, mol, 0);
+    tecplot->close();
+    delete(tecplot);
+    ofstream *rasterplot = new ofstream(OutputName.c_str());
+    write_raster3d_file(out, rasterplot, TesselStruct, mol);
+    rasterplot->close();
+    delete(rasterplot);
+  }
+
+  // free reference lists
+  if (BoundaryFreeFlag)
+    delete[] (BoundaryPoints);
+
+};
+
+
+/** Determines the volume of a cluster.
+ * Determines first the convex envelope, then tesselates it and calculates its volume.
+ * \param *out output stream for debugging
+ * \param *TesselStruct Tesselation filled with points, lines and triangles on boundary on return
+ * \param *configuration needed for path to store convex envelope file
+ * \return determined volume of the cluster in cubed config:GetIsAngstroem()
+ */
+double
+VolumeOfConvexEnvelope(ofstream *out, class Tesselation *TesselStruct, class config *configuration)
+{
+  bool IsAngstroem = configuration->GetIsAngstroem();
+  double volume = 0.;
+  double PyramidVolume = 0.;
+  double G, h;
+  Vector x, y;
+  double a, b, c;
 
   // 6a. Every triangle forms a pyramid with the center of gravity as its peak, sum up the volumes
@@ -874 +908 @@
       << endl;
 
-  // 7. translate all points back from CoG
-  *out << Verbose(1) << "Translating system back from Center of Gravity."
-      << endl;
-  CenterOfGravity->Scale(-1);
-  Walker = mol->start;
-  while (Walker->next != mol->end)
-    {
-      Walker = Walker->next;
-      Walker->x.Translate(CenterOfGravity);
-    }
-
-  // 8. Store triangles in tecplot file
-  string OutputName(filename);
-  OutputName.append(TecplotSuffix);
-  ofstream *tecplot = new ofstream(OutputName.c_str());
-  write_tecplot_file(out, tecplot, TesselStruct, mol, 0);
-  tecplot->close();
-  delete(tecplot);
-
-  // free reference lists
-  if (BoundaryFreeFlag)
-    delete[] (BoundaryPoints);
-
   return volume;
 }
@@ -919 +930 @@
   bool IsAngstroem = configuration->GetIsAngstroem();
   Boundaries *BoundaryPoints = GetBoundaryPoints(out, mol);
+  class Tesselation *TesselStruct = NULL;
+  LinkedCell LCList(mol, 10.);
+  Find_convex_border(out, mol, TesselStruct, &LCList, NULL);
   double clustervolume;
   if (ClusterVolume == 0)
-    clustervolume = VolumeOfConvexEnvelope(out, NULL, configuration,
-        BoundaryPoints, mol);
+    clustervolume = VolumeOfConvexEnvelope(out, TesselStruct, configuration);
   else
     clustervolume = ClusterVolume;
-  double *GreatestDiameter = GetDiametersOfCluster(out, BoundaryPoints, mol,
-      IsAngstroem);
+  double *GreatestDiameter = GetDiametersOfCluster(out, BoundaryPoints, mol, IsAngstroem);
   Vector BoxLengths;
   int repetition[NDIM] =
@@ -1229 +1241 @@
  */
 void
-Tesselation::TesselateOnBoundary(ofstream *out, config *configuration,
-    molecule *mol)
+Tesselation::TesselateOnBoundary(ofstream *out, molecule *mol)
 {
   bool flag;
@@ -2792 +2803 @@
  * \param *mol molecule structure with Atom's and Bond's
  * \param *Tess Tesselation filled with points, lines and triangles on boundary on return
+ * \param *LCList atoms in LinkedCell list
  * \param *filename filename prefix for output of vertex data
  * \para RADIUS radius of the virtual sphere