Changeset a439e5

Timestamp:

Aug 6, 2010, 12:40:05 PM (14 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:

6e5084

Parents:

311d688

git-author:

Frederik Heber <heber@…> (08/05/10 14:27:34)

git-committer:

Frederik Heber <heber@…> (08/06/10 12:40:05)

Message:

New functions vor class Vector and Matrix.

Matrix:

• transpose() (returning copy and on itself)
• zero() (initialize to zero)
• transformToEigenbasis() (columns are eigenvectors, return value has eigenvalues)

Vector:

• GreatestComponent() (yields index of greatest component)
• SmallestComponent() (same for smallest)

Location:

src/LinearAlgebra

Files:

• Line.cpp (modified) (1 diff)
• Matrix.cpp (modified) (4 diffs)
• Matrix.hpp (modified) (2 diffs)
• Vector.cpp (modified) (2 diffs)
• Vector.hpp (modified) (1 diff)

src/LinearAlgebra/Line.cpp

@@ -11 +11 @@
 
 #include <cmath>
+#include <iostream>
 
 #include "LinearAlgebra/Vector.hpp"

src/LinearAlgebra/Matrix.cpp

@@ -16 +16 @@
 
 #include <gsl/gsl_blas.h>
+#include <gsl/gsl_eigen.h>
+#include <gsl/gsl_matrix.h>
+#include <gsl/gsl_multimin.h>
+#include <gsl/gsl_vector.h>
 #include <cmath>
 #include <iostream>
@@ -92 +96 @@
     for(int j=NDIM;j--;){
       set(i,j,i==j);
+    }
+  }
+}
+
+void Matrix::zero(){
+  for(int i=NDIM;i--;){
+    for(int j=NDIM;j--;){
+      set(i,j,0.);
     }
   }
@@ -192 +204 @@
        - at(2,2)*at(1,0)*at(0,1);
 }
+
+Matrix Matrix::transpose() const
+{
+  Matrix copy(*this);
+  copy.transpose();
+  return copy;
+}
+
+
+void Matrix::transpose()
+{
+  double tmp;
+  for (int i=0;i<NDIM;i++)
+    for (int j=i+1;j<NDIM;j++) {
+      tmp = at(j,i);
+      at(i,j) = tmp;
+      at(j,i) = tmp;
+    }
+}
+
 
 Matrix Matrix::invert() const{
@@ -213 +245 @@
 }
 
+Vector Matrix::transformToEigenbasis()
+{
+  gsl_eigen_symmv_workspace *T = gsl_eigen_symmv_alloc(NDIM);
+  gsl_vector *eval = gsl_vector_alloc(NDIM);
+  gsl_matrix *evec = gsl_matrix_alloc(NDIM, NDIM);
+  gsl_eigen_symmv(content->content, eval, evec, T);
+  gsl_eigen_symmv_free(T);
+  gsl_matrix_memcpy(content->content, evec);
+  Vector evalues(gsl_vector_get(eval,0), gsl_vector_get(eval,1), gsl_vector_get(eval,2));
+  return evalues;
+}
+
 Matrix &Matrix::operator*=(const double factor){
   gsl_matrix_scale(content->content, factor);

src/LinearAlgebra/Matrix.hpp

@@ -51 +51 @@
 
   /**
+   * Set all matrix entries to zero.
+   */
+  void zero();
+
+  /**
    * Access the matrix at index (i,j)
    */
@@ -95 +100 @@
    */
   Matrix invert() const;
+
+  /**
+   * Diagonalizes a matrix and sets its rows to the resulting eigenvalues.
+   * The eigenvalues are returned as a vector.
+   *
+   * Rather costly, so use precomputation as often as possible.
+   */
+  Vector transformToEigenbasis();
+
+  /**
+   * Calculate the transpose of the matrix.
+   */
+  Matrix transpose() const;
+  void transpose();
 
   // operators

src/LinearAlgebra/Vector.cpp

@@ -17 +17 @@
 #include <iostream>
 #include <gsl/gsl_blas.h>
+#include <gsl/gsl_vector.h>
 
 
@@ -92 +93 @@
   return (sqrt(DistanceSquared(y)));
 };
+
+size_t Vector::GreatestComponent() const
+{
+  int greatest = 0;
+  for (int i=1;i<NDIM;i++) {
+    if (at(i) > at(greatest))
+      greatest = i;
+  }
+  return greatest;
+}
+
+size_t Vector::SmallestComponent() const
+{
+  int smallest = 0;
+  for (int i=1;i<NDIM;i++) {
+    if (at(i) < at(smallest))
+      smallest = i;
+  }
+  return smallest;
+}
+
 
 Vector Vector::getClosestPoint(const Vector &point) const{

src/LinearAlgebra/Vector.hpp

@@ -63 +63 @@
   std::pair<Vector,Vector> partition(const Vector&) const;
   std::pair<pointset,Vector> partition(const pointset&) const;
+  size_t GreatestComponent() const;
+  size_t SmallestComponent() const;
 
   // Accessors ussually come in pairs... and sometimes even more than that