Changes in src/vector.cpp [29812d:d09ff7]

File:

• src/vector.cpp (modified) (4 diffs)

src/vector.cpp

@@ -225 +225 @@
   Direction.CopyVector(LineVector);
   Direction.SubtractVector(Origin);
+  Direction.Normalize();
   //*out << Verbose(4) << "INFO: Direction is " << Direction << "." << endl;
   factor = Direction.ScalarProduct(PlaneNormal);
@@ -234 +235 @@
   helper.SubtractVector(Origin);
   factor = helper.ScalarProduct(PlaneNormal)/factor;
+  if (factor < MYEPSILON) { // Origin is in-plane
+    //*out << Verbose(2) << "Origin of line is in-plane, simple." << endl;
+    CopyVector(Origin);
+    return true;
+  }
   //factor = Origin->ScalarProduct(PlaneNormal)*(-PlaneOffset->ScalarProduct(PlaneNormal))/(Direction.ScalarProduct(PlaneNormal));
   Direction.Scale(factor);
@@ -655 +661 @@
 };
 
+/** Given a box by its matrix \a *M and its inverse *Minv the vector is made to point within that box.
+ * \param *M matrix of box
+ * \param *Minv inverse matrix
+ */
+void Vector::WrapPeriodically(const double *M, const double *Minv)
+{
+  MatrixMultiplication(Minv);
+  // truncate to [0,1] for each axis
+  for (int i=0;i<NDIM;i++) {
+    x[i] += 0.5;  // set to center of box
+    while (x[i] >= 1.)
+      x[i] -= 1.;
+    while (x[i] < 0.)
+      x[i] += 1.;
+  }
+  MatrixMultiplication(M);
+};
+
 /** Do a matrix multiplication.
  * \param *matrix NDIM_NDIM array
  */
-void Vector::MatrixMultiplication(double *M)
+void Vector::MatrixMultiplication(const double *M)
 {
   Vector C;
@@ -700 +724 @@
  * \param *matrix NDIM_NDIM array
  */
-void Vector::InverseMatrixMultiplication(double *A)
+void Vector::InverseMatrixMultiplication(const double *A)
 {
   Vector C;