vsa-lib.hpp

/*
  This file is part of LOOS.
 
  LOOS (Lightweight Object-Oriented Structure library)
  Copyright (c) 2010 Tod D. Romo
  Department of Biochemistry and Biophysics
  School of Medicine & Dentistry, University of Rochester
 
  This package (LOOS) is free software: you can redistribute it and/or modify
  it under the terms of the GNU General Public License as published by
  the Free Software Foundation under version 3 of the License.
 
  This package 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 General Public License for more details.
 
  You should have received a copy of the GNU General Public License
  along with this program.  If not, see <http://www.gnu.org/licenses/>.
 
*/
 
 
 
#if !defined(LOOS_VSA_LIB_HPP)
#define LOOS_VSA_LIB_HPP
 
#include <loos.hpp>
#include "enm-lib.hpp"
 
#if defined(__linux__) || defined(__CYGWIN__) || defined(__FreeBSD__)
extern "C" {
  void dsygvx_(int*, char*, char*, char*, int*, double*, int*, double*, int*, double*, double*, int*, int*, double*, int*, double*, double*, int*, double*, int*, int*, int*, int*);
  void dpotrf_(char*, int*, double*, int*, int*);
  void dtrmm_(char*, char*, char*, char*, int*, int*, double*, double*, int*, double*, int*);
}
#endif
 
 
namespace ENM {
 
 
 
  class VSA : public ElasticNetworkModel {
  public:
 
 
    VSA(SuperBlock* blocker, const uint subn) : 
      ElasticNetworkModel(blocker),
      subset_size_(subn)
    { prefix_ = "vsa"; }
 
 
    VSA(SuperBlock* blocker, const uint subn, const loos::DoubleMatrix& M) :
      ElasticNetworkModel(blocker),
      subset_size_(subn),
      masses_(M)
    { prefix_ = "vsa"; }
 
 
 
    void solve();
  
 
    void setMasses(const loos::DoubleMatrix& M) {
      masses_ = M;
    };
 
 
    void free() {
      masses_.reset();
      Msp_.reset();
      Hssp_.reset();
    }
 
 
  private:
    boost::tuple<loos::DoubleMatrix, loos::DoubleMatrix> eigenDecomp(loos::DoubleMatrix& A, loos::DoubleMatrix& B);
    loos::DoubleMatrix massWeight(loos::DoubleMatrix& U, loos::DoubleMatrix& M);
 
  
  private:
    uint subset_size_;
    loos::DoubleMatrix masses_;
 
    loos::DoubleMatrix Msp_;
    loos::DoubleMatrix Hssp_;
  };
 
 
};
 
 
#endif