Changeset fee69b

Timestamp:

Jan 11, 2010, 9:16:36 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:

241485

Parents:

fcad4b

git-author:

Frederik Heber <heber@…> (01/11/10 09:13:49)

git-committer:

Frederik Heber <heber@…> (01/11/10 09:16:36)

Message:

Fixed Tesselation::IsInnerPoint()

• Both unit tests TesselationInsideOutsideTest and VectorUnitTest are working (and also all others btw);
• FIX: BoundaryTriangleSet::GetIntersectionInsideTriangle() - checks whether intersection is one endpoint and performs the GetIntersectionOfTwoLinesOnPlane() right away to check whether we are inside of triangle
• FIX: Vector::GetIntersectionOfTwoLinesOnPlane() - parallelity is now checked correctly and in accordance with interface description fixed in VectorUnitTest.

Signed-off-by: Frederik Heber <heber@tabletINS.(none)>

Location:

src

Files:

• tesselation.cpp (modified) (2 diffs)
• vector.cpp (modified) (5 diffs)

src/tesselation.cpp

@@ -402 +402 @@
 bool BoundaryTriangleSet::GetIntersectionInsideTriangle(Vector *MolCenter, Vector *x, Vector *Intersection)
 {
-        Info FunctionInfo(__func__);
+  Info FunctionInfo(__func__);
   Vector CrossPoint;
   Vector helper;
@@ -415 +415 @@
   Log() << Verbose(1) << "INFO: Intersection is " << *Intersection << "." << endl;
 
+  if (Intersection->DistanceSquared(endpoints[0]->node->node) < MYEPSILON) {
+    Log() << Verbose(1) << "Intersection coindices with first endpoint." << endl;
+    return true;
+  }   else if (Intersection->DistanceSquared(endpoints[1]->node->node) < MYEPSILON) {
+    Log() << Verbose(1) << "Intersection coindices with second endpoint." << endl;
+    return true;
+  }   else if (Intersection->DistanceSquared(endpoints[2]->node->node) < MYEPSILON) {
+    Log() << Verbose(1) << "Intersection coindices with third endpoint." << endl;
+    return true;
+  }
   // Calculate cross point between one baseline and the line from the third endpoint to intersection
   int i=0;
   do {
     if (CrossPoint.GetIntersectionOfTwoLinesOnPlane(endpoints[i%3]->node->node, endpoints[(i+1)%3]->node->node, endpoints[(i+2)%3]->node->node, Intersection, &NormalVector)) {
+      helper.CopyVector(endpoints[(i+1)%3]->node->node);
+      helper.SubtractVector(endpoints[i%3]->node->node);
       CrossPoint.SubtractVector(endpoints[i%3]->node->node);  // cross point was returned as absolute vector
+      const double s = CrossPoint.ScalarProduct(&helper)/helper.NormSquared();
+      Log() << Verbose(1) << "INFO: Factor s is " << s << "." << endl;
+      if ((s < -MYEPSILON) || ((s-1.) > MYEPSILON)) {
+        Log() << Verbose(1) << "INFO: Crosspoint " << CrossPoint << "outside of triangle." << endl;
+        i=4;
+        break;
+      }
       i++;
     } else 

src/vector.cpp

@@ -306 +306 @@
 {
   Info FunctionInfo(__func__);
-  Vector a;
-  Vector b;
-  Vector c;
 
   GSLMatrix *M = new GSLMatrix(4,4);
@@ -319 +316 @@
     M->Set(3, i, Line2b->x[i]);
   }
-  Log() << Verbose(1) << "Coefficent matrix is:" << endl;
-  for (int i=0;i<4;i++) {
-    for (int j=0;j<4;j++)
-      cout << "\t" << M->Get(i,j);
-    cout << endl;
-  }
+  
+  //Log() << Verbose(1) << "Coefficent matrix is:" << endl;
+  //for (int i=0;i<4;i++) {
+  //  for (int j=0;j<4;j++)
+  //    cout << "\t" << M->Get(i,j);
+  //  cout << endl;
+  //}
   if (fabs(M->Determinant()) > MYEPSILON) {
     Log() << Verbose(1) << "Determinant of coefficient matrix is NOT zero." << endl;
@@ -333 +331 @@
 
   // constuct a,b,c
+  Vector a;
+  Vector b;
+  Vector c;
+  Vector d;
   a.CopyVector(Line1b);
   a.SubtractVector(Line1a);
@@ -339 +341 @@
   c.CopyVector(Line2a);
   c.SubtractVector(Line1a);
+  d.CopyVector(Line2b);
+  d.SubtractVector(Line1b);
   Log() << Verbose(1) << "INFO: a = " << a << ", b = " << b << ", c = " << c << "." << endl;
+  if ((a.NormSquared() < MYEPSILON) || (b.NormSquared() < MYEPSILON)) {
+   Zero();
+   Log() << Verbose(1) << "At least one of the lines is ill-defined, i.e. offset equals second vector." << endl;
+   return false;
+  }
 
   // check for parallelity
   Vector parallel;
-  parallel.CopyVector(&a);
-  parallel.SubtractVector(&b);
-  if (parallel.NormSquared() < MYEPSILON) {
+  double factor = 0.;
+  double pfactor = 0.;
+  if (fabs(a.ScalarProduct(&b)*a.ScalarProduct(&b)/a.NormSquared()/b.NormSquared() - 1.) < MYEPSILON) {
+    parallel.CopyVector(Line1a);
+    parallel.SubtractVector(Line2a);
+    factor = parallel.ScalarProduct(&a)/a.Norm();
+    if ((factor >= -MYEPSILON) && (factor - 1. < MYEPSILON)) {
+      CopyVector(Line2a);
+      Log() << Verbose(1) << "Lines conincide." << endl;
+      return true;
+    } else {
+      parallel.CopyVector(Line1a);
+      parallel.SubtractVector(Line2b);
+      factor = parallel.ScalarProduct(&a)/a.Norm();
+      if ((factor >= -MYEPSILON) && (factor - 1. < MYEPSILON)) {
+        CopyVector(Line2b);
+        Log() << Verbose(1) << "Lines conincide." << endl;
+        return true;
+      }
+    }
     Log() << Verbose(1) << "Lines are parallel." << endl;
+    Zero();
     return false;
   }
@@ -370 +397 @@
   Log() << Verbose(1) << "Intersection is at " << *this << "." << endl;
 
-  if ((s >=0 ) && (s<=1))
-    return true;
-  else
-    return false;
+  return true;
 };