linear_algebra.h

Go to the documentation of this file.

// -*- mode: c++; -*-
//
//The Biomolecule Toolkit (BTK) is a C++ library for use in the
//modeling, analysis, and design of biological macromolecules.
//Copyright (C) 2001, Tim Robertson <kid50@users.sourceforge.net>
//
//This program 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 program 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 program; if not, write to the Free Software Foundation, 
//Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.


#ifndef BTK_LINEAR_ALGEBRA_H
#define BTK_LINEAR_ALGEBRA_H


#include <cassert>
#include <cmath>

#include <iosfwd>
#include <vector>

#include <boost/numeric/ublas/vector.hpp>
#include <boost/numeric/ublas/matrix.hpp>
#include <boost/numeric/ublas/symmetric.hpp>
#include <boost/numeric/ublas/io.hpp>

namespace uBLAS = boost::numeric::ublas;

namespace BTK {


const double PI = 3.1415926535897323;
const double DEG2RAD = PI/180.0;
const double RAD2DEG = 180.0/PI;

typedef uBLAS::vector<double> BTKVector;
typedef uBLAS::matrix<double> BTKMatrix;
typedef uBLAS::symmetric_matrix<double> BTKSymmetricMatrix;

// BTKVector convenience functions
//
inline 
BTKVector 
new_vector(double x=0.,
           double y=0.,
           double z=0.)
{
  BTKVector temp(3);
  temp[0]=x; temp[1]=y; temp[2]=z;
  return temp;
}

inline double length(const BTKVector& v) { return norm_2(v); }

inline 
BTKVector 
cross(const BTKVector& lhs, 
      const BTKVector& rhs) 
{
  BTKVector temp(3);
  
  temp[0] = lhs[1] * rhs[2] - lhs[2] * rhs[1];
  temp[1] = lhs[2] * rhs[0] - lhs[0] * rhs[2];
  temp[2] = lhs[0] * rhs[1] - lhs[1] * rhs[0];
  
  return temp;
}

BTKVector 
project(const BTKVector& a,
        const BTKVector& b);

BTKVector 
project_normal(const BTKVector& a,
               const BTKVector& b);

bool 
within_sqr_dist(const BTKVector& a,
                const BTKVector& b,
                double r_squared);

inline 
BTKVector& 
normalize(BTKVector & v) 
{
  v /= norm_2(v);
  return v;
}

inline 
BTKVector 
normalize(const BTKVector& v) 
{
  BTKVector tmp = v / norm_2(v);
  return tmp;
}

inline 
double 
cosine_vector_angle(const BTKVector& v1,
                    const BTKVector& v2) 
{
  return (prec_inner_prod(v1,v2) / (norm_2(v1) * norm_2(v2)));
}

inline 
double 
vector_angle(const BTKVector& v1,
             const BTKVector& v2) 
{
  return acos(cosine_vector_angle(v1,v2));
}

inline
double
point_angle(const BTKVector& v1,
            const BTKVector& v2,
            const BTKVector& v3){
  return vector_angle(v1-v2,v3-v2);
}


double 
vector_dihedral(const BTKVector& v1,
                const BTKVector& v2,
                const BTKVector& v3);


double 
point_dihedral(const BTKVector& a,
               const BTKVector& b,
               const BTKVector& c,
               const BTKVector& d);

BTKVector
set_vector_from_dihedral(const BTKVector& v3,
                         const BTKVector& v2,
                         const BTKVector& v1,
                         double len34,
                         double ang234,
                         double dih1234);

BTKVector
set_vector_from_two_angles(const BTKVector& v3,
                           const BTKVector& v2,
                           const BTKVector& v1,
                           double len34,
                           double ang234,
                           double ang134);


struct EigenState {
  EigenState() {}
  EigenState(double e,BTKVector& v) :
    eigenvalue(e), eigenvector(v) {}

  bool 
  operator<(const EigenState& right) const 
  {
      return eigenvalue < right.eigenvalue;
  }

  double eigenvalue;
  BTKVector eigenvector;
};


std::ostream& 
operator<<(std::ostream& os, 
           const std::vector<EigenState>& es);

std::vector<EigenState>
diagonalize_symmetric(const BTKMatrix& m);





BTKMatrix 
create_rotation_matrix(const BTKVector& axis, 
                       double theta);
void 
initialize_rotation_matrix(const BTKVector& axis, 
                           double theta,
                           BTKMatrix& rm);


BTKMatrix 
create_rotation_matrix(double phi, 
                       double theta, 
                       double psi);

void 
initialize_rotation_matrix(double phi, 
                           double theta, 
                           double psi,
                           BTKMatrix& rm);

}

#endif // LINEAR_ALGEBRA_H

---