1 December 2009

© Copyright 1973, 2009 Larry Andrews. All rights reserved

## LGPL NOTICESThis library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation; either version 2.1 of the License, or (at your option) any later version. This library is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU* Lesser General Public License for more details. You should have received a copy of the GNU Lesser General Public License along with this library; if not, write to the Free Software Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA |

Originally conceived in 1973 to support computer graphics, Vector is a modest linear algebra package specifically written for 3-dimensions. Although most of the algorithms are not restrictive, no higher space is supported. The original version was Fortran, here included with little modification. The C++ version is largely a translation of the Fortran version. 3-D vector algebra is well supported, and Vector has been widely used in the crystallographic community for programs to support graphics, analysis, and structure construction.

This release provides a vector_3d.h and vector_3d.cpp, and a FORTRAN library, and Vector_3_Test.cpp, with example/test programs

The Vector package is available at www.sourceforge.net/projects/vector. A source zip is available at downloads.sourceforge.net/vector/vector-0.1.3.zip. Later tarballs may be available.

The C++ code includes 4 classes. Vector_3 and Matrix_3x3 comprise most of the code. Line_3 and Plane_3 are small, rather incomplete classes for use with the other two.

Included also is a group of random numbers functions derived from code originally written by Rob Harrison. The function RandomOrientation, which returns a random rotation matrix, is the likely main one for use with Vector.

- Installation
- vector_3d.h, The declarations for Vector
- vector_3d.cpp, The definitions for Vector
- TestVector_3d.cpp, Test code for Vector

The header vector_3d.h and definitions source code vector_3d.cpp should be in the directory of the program to be compiled.

The provided interface is:

#include "vector_3d.h"FUNCTIONS IN Vector_3Vector_3(void); // constructor Default constructor - x,y,z are all initialized to DBL_MAXexplicit Vector_3(const doublearg); // constructorargis used to initialize each of x,y,zVector_3(const doubleax, const doubleay, const doubleaz); // constructor x,y,z are initialized toax,ay,az, respectivelyVector_3(const Vector_3&input_vector); // copy constructor the copy constructor returns a copy ofinput_vectorvirtual ~Vector_3(void); // destructorSCALAR OPERATIONS ON VECTORSVector_3 operator*(const doublearg)const; returns a vector that is the input multiplied by the floating point scaling factor,argVector_3 operator/(const doublearg)const; returns a vector that is the input divided by the floating point scaling factor,argOPERATIONS ONLY ON VECTORSdouble Dot(const Vector_3&)const; returns the dot product of two vectorsVector_3 DyadProduct(const Vector_3&)const; returns a vector x(result)=x1*x2, y(result)=y1*y2, z(result)=s1*z2Vector_3 Cross(const Vector_3&)const; typical cross product of two vectors (also known as "vector product")Matrix_3x3 TensorProduct(const Vector_3&)const; where the dot product isv1timesv2-transpose, this returns v1-transpose timesv2, populating a matrix with all pairs of x,y,z(1) times x,y,z(2) (also known as "outer product")Vector_3 operator+(const Vector_3&)const; returns the vector sum of two vectorsVector_3 operator-(const Vector_3&)const; returns the vector difference of two vectorsVector_3 operator-( void )const; unary minus operator for a vectorVector_3& operator+=(const Vector_3&); in place vector sum of two vectorsVector_3& operator-=(const Vector_3&); in place subtraction of two vectorsdouble Norm( void )const; returns the length of a vectordouble SquaredLength( void )const; returns the squared length of a vector (avoids the square root computation required by Norm(vector))Vector_3 UnitV( void )const; returns a unit vector with the same direction (therefore same direction cosines) as the input vectorVector_3& operator=(const Vector_3&); assignment operator of vectorsstatic Vector_3GetCenterOfMass(const Container&); returns the vector to the center of mass of in input group of vectorsstatic Matrix_3x3Pair(const Vector_3&v1,const Vector_3&v2,const Vector_3&v3,const Vector_3&v4); returns a matrix that will rotatev1ontov3and placev2as close as possible tov4while maintaining the alignment ofv1,v3static Vector_3GetXAxis( void ); returns the vector (1,0,0)static Vector_3GetYAxis( void ); returns the vector (0,1,0)static Vector_3GetZAxis( void ); returns the vector (0,0,1)static Vector_3GetZeroVector(void ); returns the vector (0,0,0)bool operator==(const Vector_3&)const;bool operator!=(const Vector_3&)const;static double Angle(const Vector_3&v1;,const Vector_3&v2,const Vector_3&v3); returns the angle between the vector fromv2tov1and the vector fromv2tov3. Torsion and Angle differ from most of the functions in Vector_3 in that they use the ends of vectorsstatic double Torsion(const Vector_3&v1,const Vector_3&v2,const Vector_3&v3,const Vector_3&v4); returns the torsion angle from the sequence to the ends ofv1,v2,v3,v4, respectively. Torsion and Angle differ from most of the functions in Vector_3 in that they use the ends of vectorsdouble DistanceFromPlane(const VPlane&)const; returns the distance of a point from a 3-space planedouble DistanceFromLine(const VLine&)const; returns the distance of a point from a 3-space lineOPERATIONS USING MATRICESVector_3 MV(const Matrix_3x3&)const; premultiply a vector by a matrixCOMPLEX OPERATIONSMatrix_3x3 Rotmat(const doubleangle)const; returns the matrix that will rotate points about the input vector by the specifiedangleMatrix_3x3 InertiaTensorForPoint(const doubleweight)const; returns the inertia tensor for a point (used to compute best line and plane)weightis used to weight the importance of the input pointsACCESS FUNCTIONSdouble operator[](const int&n )const; for input values of 0,1,2,operator[]returns x,y,z, respectivelydouble at(const int&)const; for input values of 0,1,2,atreturns x,y,z, respectivelyfriend functionsVector_3 operator*(const double&value, const Vector_3&v); returns a vector that isvaluetimes the input vectorVPlane_3 BestPlane(const Container&); returns the plane that best approximates the list of points in the input container (which can be any container class, such as std::vector, that has proper iterators)VLine_3 BestLine(const Container&); returns the line that best approximates the list of points in the input container (which can be any container class, such as std::vector, that has proper iterators)FUNCTIONS IN Matrix_3x3Matrix_3x3( void ); // constructor constructs a Matrix_3x3 with all entries initialized to DBL_MAXMatrix_3x3(const double,const double,const double,const double,const double,const double,const double,const double,const double); // constructor returns a Matrix_3x3 with the elements equal to the 9 input values, in row orderMatrix_3x3(const Matrix_3x3&); // copy constructorvirtual ~Matrix_3x3( void ); // destructor MATRIX/VECTOR OPERATIONSVector_3 MV(const Vector_3&) const; returns a vector that is the input vector premultiplied by the input matrix OPERATIONS ON MATRICES ONLYMatrix_3x3 operator*(const Matrix_3x3&)const; returns the product of the two input matricesMatrix_3x3 operator+(const Matrix_3x3&)const; returns the sum of the two input matricesMatrix_3x3 operator-(const Matrix_3x3&)const; returns the difference of the two input matricesMatrix_3x3 operator*(const doubled)const; returns a matrix that is multiplied bydMatrix_3x3 operator/(const doubled)const; returns a matrix that is divided bydMatrix_3x3& operator+=(const Matrix_3x3&); in place sum of two matricesMatrix_3x3& operator-=(const Matrix_3x3&); in place difference of two matricesMatrix_3x3& operator*=(const Matrix_3x3&); in place product of two matricesMatrix_3x3& operator=(const Matrix_3x3&); // copy constructorMatrix_3x3 Transpose( void )const; returns the transpose of the input matrixMatrix_3x3 Inver( void )const; returns the inverse of the input matrixvoid UnitMatrix( void ); in place convesion of the input matrix to a unit matrixdouble Det( void )const; returns the determinant of the input matrixvoid Zero( void ); converts the input matrix to a zero matrixbool Eigen1(double&eigenvalue,Vector_3&eigenvector)const; returns the largest eigenvalue of the input matrix and its corresponding eigenvectorbool Eigen3(std::vector& eigenvalues,std::vector& eigenvectors)const; returns all 3eigenvaluesof the input matrix and their correspondingeigenvectorsdouble operator[](const int&)const; returns then-th element of the input matrix (in row order)Matrix_3x3 Cart(const doublea,const doubleb,const doublec,const doublealpha,const doublebeta,const doublegamma)Cartis used to generate a matrix that will convert from coordinates in a non-orthogonal coordinate system to an orthogonal basis. This is the solution for the common problem in crystallography, where coordinates are often reported in coordinates that are measured using the unit cell dimensions.Cartwill return a matrix that when multiplied by a vector expressed in fractional coordinates will generate the corresponding position in an orthonormal system. NOTE: VERY IMPORTANT: The coordinates created by that operation may NOT be comparable to those generated by some other conversion software. Each system assumes a particular relationship between the two system, and they are ALL correct. Basically, the conversion matrix can be multiplied by ANY rotation matrix at all, and the result will still generate correct orthogonal coordinates, just different ones. The inverse of the matrix thatCartgenerates will convert from orthogonal coordinates back to fractional ONLY IF the orginal transformation was made with the same convention thatCartuses. The convention used inCartis that the a-crystallographic axis will be aligned parallel to x, the b-axis as close as possible to y, and the c-axis by construction of a right handed-coordinate system.OTHER FUNCTIONSstatic std::vectorAB_BestRotation const ContainerType&v1,const ContainerType&v2,Vector_3¢erOfMass1,Vector_3¢erOfMass2,Matrix_3x3&rotationMatrixFrom2To1,double&rmsd,const intpairCount=0 ) A common problem in crystallography is superimposing two similar structures. AB_BestRotation is a new (2009) algorithm for solving this problem. Currently, this is a somewhat preliminary version of the method. The transformations that it produces are nearly identical to the standard Kabsch algorithm's results. Interest in it centers on its abilities to use some interesting weighting schemes and possibly substructure determination. This early version does not yet implement these. Its behavior for small numbers of input points has not been fully investigated, but it seems quite stable for the number of atoms found in a typical protein structure.RETURN: the returned vector is a quaternion of four elements: cos(theta/2), sin(theta/2)*x, sin(theta/2)*y, sin(theta/2)*z), where x,y,z are the direction cosines of the axis of rotation.const ContainerType&v1v1is the list of input points for the first (reference) position.ContainerTypecan be any STL container or other containers (e.g. NearTree) that have a resonable iterator and a size() function.const ContainerType&v2Likev1, but this is the structure to be rotated ontov1.v2andv1must have exactly the same number of points in corresponding positions.Vector_3¢erOfMass1The returned center of mass of the points inv1.Vector_3¢erOfMass2The returned center of mass of the points inv2.Matrix_3x3&rotationMatrixFrom2To1The returned rotation matrix that will rotatev2ontov1(as well as possible).double&rmsdThe returned root mean square distance between the corresponding points ofv1andv2.const intpairCountpairCountcontrols the nature of the calculation. IfpairCountis zero (the default value), then every pair of atoms (that is in its original and rotated positions) is used with every other pair in the computation (an O(n^2) process). IfpairCountis greater than zero, then each is compared withpairCountother pairs. IfpairCountis one, then the process is exactly O(n). Larger values have the potential to give more representative output, but preliminary results seem to indicate that for a protein, this would not make much difference. For small molecules, the exhaustive search is probably safer, and is quite fast. For values ofpairCountlarger than 0, the other pairs are chosen by random selection.FUNCTIONS IN Line_3Line_3( void ); // constructorLine_3(const Vector_3&; direction,const Vector_3&CenterOfMass ); constructs a Line_3 given a direction of the line and a point on the lineVector_3 LineAxis( void )const; returns the direction vector of a Line_3Vector_3 BasePoint( void )const; returns the point on a Line_3 that is closest to the coordinate system originvirtual ~Line_3( void ); // destructordouble DistanceFromLine(const Vector_3&v)const; for a given point, returns the distance of the point from the input lineFUNCTIONS IN Plane_3Plane_3( void ); // constructorPlane_3::Plane_3(const Vector_3&direction,const Vector_3&CenterOfMass); constructs a plane given a vector perpendicular to the plane and a point in the planeVector_3 NormalVector( void )const; returns the normal to the planeVector_3 BasePoint( void )const; returns the point in the input plane that is closest to the originvirtual ~Plane_3( void ); // destructordouble DistanceFromPlane(const Vector_3&v)const; returns the distance of the input pointvfrom the input planeOTHER FUNCTIONSMatrix_3x3 RandomOrientation( void ) returns a rotation matrix from a well-sampled distributionVector_3 RandomVector( void ) returns a vector that points in a random direction

Updated 1 December 2009