Changeset 0a4f7f for src/boundary.cpp

Timestamp:

Apr 7, 2010, 3:45:38 PM (15 years ago)

Author:

Tillmann Crueger <crueger@…>

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:

72e7fa

Parents:

f9352d

Message:

Made data internal data-structure of vector class private

• Replaced occurences of access to internals with operator
• moved Vector-class into LinAlg-Module
• Reworked Vector to allow clean modularization
• Added Plane class to describe arbitrary planes in 3d space

File:

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

src/boundary.cpp

@@ -18 +18 @@
 #include "tesselation.hpp"
 #include "tesselationhelpers.hpp"
+#include "Plane.hpp"
 
 #include<gsl/gsl_poly.h>
@@ -80 +81 @@
               DistanceVector.SubtractVector(&Neighbour->second.second->x);
               do { // seek for neighbour pair where it flips
-                  OldComponent = DistanceVector.x[Othercomponent];
+                  OldComponent = DistanceVector[Othercomponent];
                   Neighbour++;
                   if (Neighbour == BoundaryPoints[axis].end()) // make it wrap around
@@ -88 +89 @@
                   //Log() << Verbose(2) << "OldComponent is " << OldComponent << ", new one is " << DistanceVector.x[Othercomponent] << "." << endl;
                 } while ((runner != Neighbour) && (fabs(OldComponent / fabs(
-                  OldComponent) - DistanceVector.x[Othercomponent] / fabs(
-                  DistanceVector.x[Othercomponent])) < MYEPSILON)); // as long as sign does not flip
+                  OldComponent) - DistanceVector[Othercomponent] / fabs(
+                  DistanceVector[Othercomponent])) < MYEPSILON)); // as long as sign does not flip
               if (runner != Neighbour) {
                   OtherNeighbour = Neighbour;
@@ -102 +103 @@
                   //Log() << Verbose(1) << "OtherComponents to Neighbour and OtherNeighbour are " << DistanceVector.x[Othercomponent] << " and " << OtherVector.x[Othercomponent] << "." << endl;
                   // do linear interpolation between points (is exact) to extract exact intersection between Neighbour and OtherNeighbour
-                  w1 = fabs(OtherVector.x[Othercomponent]);
-                  w2 = fabs(DistanceVector.x[Othercomponent]);
-                  tmp = fabs((w1 * DistanceVector.x[component] + w2
-                      * OtherVector.x[component]) / (w1 + w2));
+                  w1 = fabs(OtherVector[Othercomponent]);
+                  w2 = fabs(DistanceVector[Othercomponent]);
+                  tmp = fabs((w1 * DistanceVector[component] + w2
+                      * OtherVector[component]) / (w1 + w2));
                   // mark if it has greater diameter
                   //Log() << Verbose(1) << "Comparing current greatest " << GreatestDiameter[component] << " to new " << tmp << "." << endl;
@@ -159 +160 @@
     AngleReferenceVector.Zero();
     AngleReferenceNormalVector.Zero();
-    AxisVector.x[axis] = 1.;
-    AngleReferenceVector.x[(axis + 1) % NDIM] = 1.;
-    AngleReferenceNormalVector.x[(axis + 2) % NDIM] = 1.;
+    AxisVector[axis] = 1.;
+    AngleReferenceVector[(axis + 1) % NDIM] = 1.;
+    AngleReferenceNormalVector[(axis + 2) % NDIM] = 1.;
 
     Log() << Verbose(1) << "Axisvector is " << AxisVector << " and AngleReferenceVector is " << AngleReferenceVector << ", and AngleReferenceNormalVector is " << AngleReferenceNormalVector << "." << endl;
@@ -636 +637 @@
       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 = Plane(*(runner->second->endpoints[0]->node->node),
+                *(runner->second->endpoints[1]->node->node),
+                *(runner->second->endpoints[2]->node->node)).getNormal();
       x.Scale(runner->second->endpoints[1]->node->node->ScalarProduct(&x));
       const double h = x.Norm(); // distance of CoG to triangle
@@ -751 +754 @@
     Log() << Verbose(0) << "Setting Box dimensions to minimum possible, the greatest diameters." << endl;
     for (int i = 0; i < NDIM; i++)
-      BoxLengths.x[i] = GreatestDiameter[i];
+      BoxLengths[i] = GreatestDiameter[i];
     mol->CenterEdge(&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;
+    BoxLengths[0] = (repetition[0] * GreatestDiameter[0] + repetition[1] * GreatestDiameter[1] + repetition[2] * GreatestDiameter[2]);
+    BoxLengths[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[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
+    if (gsl_poly_solve_cubic(BoxLengths[0], BoxLengths[1], BoxLengths[2], &x0, &x1, &x2) == 1) // either 1 or 3 on return
       Log() << Verbose(0) << "RESULT: The resulting spacing is: " << x0 << " ." << endl;
     else {
@@ -769 +772 @@
     cellvolume = 1.;
     for (int i = 0; i < NDIM; i++) {
-      BoxLengths.x[i] = repetition[i] * (x0 + GreatestDiameter[i]);
-      cellvolume *= BoxLengths.x[i];
+      BoxLengths[i] = repetition[i] * (x0 + GreatestDiameter[i]);
+      cellvolume *= BoxLengths[i];
     }
 
@@ -780 +783 @@
   // update Box of atoms by boundary
   mol->SetBoxDimension(&BoxLengths);
-  Log() << 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;
+  Log() << Verbose(0) << "RESULT: The resulting cell dimensions are: " << BoxLengths[0] << " and " << BoxLengths[1] << " and " << BoxLengths[2] << " with total volume of " << cellvolume << " " << (IsAngstroem ? "angstrom" : "atomiclength") << "^3." << endl;
 };
 
@@ -839 +842 @@
   FillerDistance.InverseMatrixMultiplication(M);
   for(int i=0;i<NDIM;i++)
-    N[i] = (int) ceil(1./FillerDistance.x[i]);
+    N[i] = static_cast<int>(ceil(1./FillerDistance[i]));
   Log() << Verbose(1) << "INFO: Grid steps are " << N[0] << ", " << N[1] << ", " << N[2] << "." << endl;
 
@@ -880 +883 @@
           // create molecule random translation vector ...
           for (int i=0;i<NDIM;i++)
-            FillerTranslations.x[i] = RandomMolDisplacement*(rand()/(RAND_MAX/2.) - 1.);
+            FillerTranslations[i] = RandomMolDisplacement*(rand()/(RAND_MAX/2.) - 1.);
           Log() << Verbose(2) << "INFO: Translating this filler by " << FillerTranslations << "." << endl;
 
@@ -892 +895 @@
             // create atomic random translation vector ...
             for (int i=0;i<NDIM;i++)
-              AtomTranslations.x[i] = RandomAtomDisplacement*(rand()/(RAND_MAX/2.) - 1.);
+              AtomTranslations[i] = RandomAtomDisplacement*(rand()/(RAND_MAX/2.) - 1.);
 
             // ... and rotation matrix