diff --git a/doc/tests/scripts/sidechain_steps.py b/doc/tests/scripts/sidechain_steps.py index 34f0bef26684f53b20106ace122045674e14132b..67d890d6c348c2e25daf714686e5072f4e94fba7 100644 --- a/doc/tests/scripts/sidechain_steps.py +++ b/doc/tests/scripts/sidechain_steps.py @@ -66,7 +66,7 @@ for i,r in enumerate(prot.residues): # buildup a graph given the rotamer groups -graph = sidechain.Graph.CreateFromFRMList(rotamer_groups) +graph = sidechain.RotamerGraph.CreateFromFRMList(rotamer_groups) # and get a solution out of it solution = graph.Solve() diff --git a/sidechain/doc/graph.rst b/sidechain/doc/graph.rst index f7d4802311d89a3786fefc360c2134db35fd1e11..e4228af58af089758566e665f2766234e4a1264d 100644 --- a/sidechain/doc/graph.rst +++ b/sidechain/doc/graph.rst @@ -12,10 +12,10 @@ decomposition and dead end elimination, decomposes the graph in a tree form and finally comes back with a solution. -The Interaction Graph +The Rotamer Graph -------------------------------------------------------------------------------- -.. class:: Graph +.. class:: RotamerGraph The Graph object has no constructor exported to python. It is meant to be constructed using the static create functions. @@ -40,9 +40,9 @@ The Interaction Graph .. method:: Solve([max_complexity=inf, initial_epsilon=0.02]) - The method solves a graph by iteratively performing dead end elimination - and edge decomposition until the number of possible permutations, that - have to be enumerated in the resulting tree, fall below + The method solves a rotamer graph by iteratively performing dead end + elimination and edge decomposition until the number of possible + permutations, that have to be enumerated in the resulting tree, fall below **max_complecity**. In every iteration the epsilon value gets doubled to enforce more edges to be approximated by edge decomposition. This function assumes all self energies of the rotamers to be calculated and properly set! diff --git a/sidechain/doc/index.rst b/sidechain/doc/index.rst index 98c570f303b9a79bdea34dd6a866864acb93ce05..e8528c16e10237efb6e3efc84d2f5247a9f9c875 100644 --- a/sidechain/doc/index.rst +++ b/sidechain/doc/index.rst @@ -44,7 +44,7 @@ Contents: Rotamers <rotamer> Frame <frame> Rotamer Library <rotamer_lib> - Graph <graph> + RotamerGraph <graph> The Settings - Control Things... <sidechain_settings> Disulfid Bond Evaluation <disulfid> Loading Libraries <loading> diff --git a/sidechain/pymod/_reconstruct_sidechains.py b/sidechain/pymod/_reconstruct_sidechains.py index 589b00e758e165c7828c2f62e81c656404f942eb..145c60673774781a17cfc834afd4826373f505cf 100644 --- a/sidechain/pymod/_reconstruct_sidechains.py +++ b/sidechain/pymod/_reconstruct_sidechains.py @@ -496,9 +496,9 @@ def Reconstruct(ent, keep_sidechains=False, build_disulfids=True, # set up graph and solve to get best rotamers if use_frm: - graph = sidechain.Graph.CreateFromFRMList(rotamer_groups) + graph = sidechain.RotamerGraph.CreateFromFRMList(rotamer_groups) else: - graph = sidechain.Graph.CreateFromRRMList(rotamer_groups) + graph = sidechain.RotamerGraph.CreateFromRRMList(rotamer_groups) solution = graph.Solve(100000000,0.02) diff --git a/sidechain/pymod/export_graph.cc b/sidechain/pymod/export_graph.cc index 37c0201b5807190494f6cebd2d667d9a7edd08e6..87f647daa1a20ada25889b68cf8de204c3d42f8e 100644 --- a/sidechain/pymod/export_graph.cc +++ b/sidechain/pymod/export_graph.cc @@ -9,25 +9,25 @@ using namespace promod3; namespace { -GraphPtr WrapRRMList(boost::python::list& rotamer_groups, - boost::python::list& rt_operators) { +RotamerGraphPtr WrapRRMList(boost::python::list& rotamer_groups, + boost::python::list& rt_operators) { std::vector<RRMRotamerGroupPtr> v_rotamer_groups; std::vector<geom::Transform> v_rt_operators; core::ConvertListToVector(rotamer_groups, v_rotamer_groups); core::ConvertListToVector(rt_operators, v_rt_operators); - return Graph::CreateFromRRMList(v_rotamer_groups,v_rt_operators); + return RotamerGraph::CreateFromRRMList(v_rotamer_groups,v_rt_operators); } -GraphPtr WrapFRMList(boost::python::list& rotamer_groups, - boost::python::list& rt_operators) { +RotamerGraphPtr WrapFRMList(boost::python::list& rotamer_groups, + boost::python::list& rt_operators) { std::vector<FRMRotamerGroupPtr> v_rotamer_groups; std::vector<geom::Transform> v_rt_operators; core::ConvertListToVector(rotamer_groups, v_rotamer_groups); core::ConvertListToVector(rt_operators, v_rt_operators); - return Graph::CreateFromFRMList(v_rotamer_groups,v_rt_operators); + return RotamerGraph::CreateFromFRMList(v_rotamer_groups,v_rt_operators); } -boost::python::list WrapSolve(GraphPtr graph, uint64_t max_complexity, +boost::python::list WrapSolve(RotamerGraphPtr graph, uint64_t max_complexity, Real initial_epsilon) { std::vector<int> solution = graph->Solve(max_complexity,initial_epsilon); boost::python::list return_list; @@ -40,19 +40,19 @@ boost::python::list WrapSolve(GraphPtr graph, uint64_t max_complexity, void export_Graph() { - class_<Graph,boost::noncopyable>("Graph", no_init) + class_<RotamerGraph,boost::noncopyable>("RotamerGraph", no_init) .def("CreateFromRRMList",&WrapRRMList,(arg("rotamer_groups"), arg("rt_operators")=boost::python::list())).staticmethod("CreateFromRRMList") .def("CreateFromFRMList",&WrapFRMList,(arg("rotamer_groups"), arg("rt_operators")=boost::python::list())).staticmethod("CreateFromFRMList") - .def("Prune",&Graph::Prune,(arg("epsilon"),arg("consider_all_nodes")=false)) - .def("Reset", &Graph::Reset) - .def("GetNumNodes", &Graph::GetNumNodes) - .def("GetNumEdges", &Graph::GetNumEdges) - .def("GetNumActiveNodes", &Graph::GetNumActiveNodes) - .def("GetNumActiveEdges", &Graph::GetNumActiveEdges) + .def("Prune",&RotamerGraph::Prune,(arg("epsilon"),arg("consider_all_nodes")=false)) + .def("Reset", &RotamerGraph::Reset) + .def("GetNumNodes", &RotamerGraph::GetNumNodes) + .def("GetNumEdges", &RotamerGraph::GetNumEdges) + .def("GetNumActiveNodes", &RotamerGraph::GetNumActiveNodes) + .def("GetNumActiveEdges", &RotamerGraph::GetNumActiveEdges) .def("Solve",&WrapSolve,(arg("max_complexity")=std::numeric_limits<uint64_t>::max(),arg("initial_epsilon")=0.02)) ; - register_ptr_to_python<GraphPtr>(); + register_ptr_to_python<RotamerGraphPtr>(); } diff --git a/sidechain/src/graph.cc b/sidechain/src/graph.cc index ddc2ab4d99a5624d64855e9eed65c98b5ac87cff..d70e75613a98d671562c3758ecfa55bc82df91ee 100644 --- a/sidechain/src/graph.cc +++ b/sidechain/src/graph.cc @@ -12,31 +12,36 @@ namespace{ namespace promod3{ namespace sidechain{ -Edge::Edge(uint index, Node* node1, Node* node2, uint num_i, uint num_j, - Real* pairwise_energies): index_(index), - active_(true), - node1_(node1), node2_(node2), - maxi_(num_i), maxj_(num_j), - pairwise_energies_(pairwise_energies){ +RotamerGraphEdge::RotamerGraphEdge(uint index, + RotamerGraphNode* node1, + RotamerGraphNode* node2, + uint num_i, uint num_j, + Real* pairwise_energies): index_(index), + active_(true), + node1_(node1), + node2_(node2), + maxi_(num_i), + maxj_(num_j), + pairwise_energies_(pairwise_energies){ } -Edge::~Edge(){ +RotamerGraphEdge::~RotamerGraphEdge(){ delete [] pairwise_energies_; //Nodes don't get deleted here! //It is assumed, that this happens in the graph object } -void Edge::Reset(){ +void RotamerGraphEdge::Reset(){ active_ = true; } -void Edge::Deactivate(){ +void RotamerGraphEdge::Deactivate(){ if(!active_) return; active_ = false; } -bool Edge::Decompose(Real thresh){ +bool RotamerGraphEdge::Decompose(Real thresh){ if(!active_) return false; @@ -46,10 +51,10 @@ bool Edge::Decompose(Real thresh){ uint n_total = node2_->GetNumRotamers(); int temp; - for(Node::iterator i = node1_->active_rotamers_begin(); + for(RotamerGraphNode::iterator i = node1_->active_rotamers_begin(); i != node1_->active_rotamers_end(); ++i){ temp = (*i) * n_total; - for(Node::iterator j = node2_->active_rotamers_begin(); + for(RotamerGraphNode::iterator j = node2_->active_rotamers_begin(); j != node2_->active_rotamers_end(); ++j){ summed_e += pairwise_energies_[temp + (*j)]; } @@ -61,11 +66,11 @@ bool Edge::Decompose(Real thresh){ Real bl[n]; uint m_counter = 0; - for(Node::iterator i = node1_->active_rotamers_begin(); + for(RotamerGraphNode::iterator i = node1_->active_rotamers_begin(); i != node1_->active_rotamers_end(); ++i){ summed_e = 0.0; int temp = (*i) * n_total; - for(Node::iterator j = node2_->active_rotamers_begin(); + for(RotamerGraphNode::iterator j = node2_->active_rotamers_begin(); j != node2_->active_rotamers_end(); ++j){ summed_e += pairwise_energies_[temp + (*j)]; } @@ -75,19 +80,20 @@ bool Edge::Decompose(Real thresh){ Real e; uint n_counter = 0; - for(Node::iterator i = node2_->active_rotamers_begin(); + for(RotamerGraphNode::iterator i = node2_->active_rotamers_begin(); i != node2_->active_rotamers_end(); ++i){ summed_e = 0.0; - for(Node::iterator j = node1_->active_rotamers_begin(); + for(RotamerGraphNode::iterator j = node1_->active_rotamers_begin(); j != node1_->active_rotamers_end(); ++j){ summed_e += pairwise_energies_[(*j) * node2_->GetNumRotamers() + (*i)]; } bl[n_counter] = summed_e/m - avg_value; m_counter = 0; - for(Node::iterator j = node1_->active_rotamers_begin(); + for(RotamerGraphNode::iterator j = node1_->active_rotamers_begin(); j != node1_->active_rotamers_end(); ++j){ - e = std::abs(pairwise_energies_[(*j)*n_total+(*i)]-ak[m_counter]-bl[n_counter]); + e = std::abs(pairwise_energies_[(*j)*n_total+(*i)]- + ak[m_counter]-bl[n_counter]); if(e > thresh) return false; ++m_counter; } @@ -95,14 +101,14 @@ bool Edge::Decompose(Real thresh){ } for(uint i = 0; i < m; ++i){ - for(Node::iterator j = node1_->active_rotamers_begin(); + for(RotamerGraphNode::iterator j = node1_->active_rotamers_begin(); j != node1_->active_rotamers_end(); ++j){ node1_->AddValue(*j,ak[i]); } } for(uint i = 0; i < n; ++i){ - for(Node::iterator j = node2_->active_rotamers_begin(); + for(RotamerGraphNode::iterator j = node2_->active_rotamers_begin(); j != node2_->active_rotamers_end(); ++j){ node2_->AddValue(*j,bl[i]); } @@ -113,33 +119,38 @@ bool Edge::Decompose(Real thresh){ return true; } -Real Edge::EMinDiff(const Node* node_ptr, uint idx1, uint idx2) const{ +Real RotamerGraphEdge::EMinDiff(const RotamerGraphNode* node_ptr, + uint idx1, uint idx2) const{ Real min = std::numeric_limits<Real>::max(); if(node_ptr == node1_){ int temp1 = idx1 * node2_->GetNumRotamers(); int temp2 = idx2 * node2_->GetNumRotamers(); - for(Node::const_iterator i = node2_->active_rotamers_begin(); + for(RotamerGraphNode::const_iterator i = node2_->active_rotamers_begin(); i != node2_->active_rotamers_end(); ++i){ - min = std::min(min,pairwise_energies_[temp1 + (*i)] - pairwise_energies_[temp2 + (*i)]); + min = std::min(min,pairwise_energies_[temp1 + (*i)] - + pairwise_energies_[temp2 + (*i)]); } return min; } else if(node_ptr == node2_){ int num = node2_->GetNumRotamers(); - for(Node::const_iterator i = node1_->active_rotamers_begin(); + for(RotamerGraphNode::const_iterator i = node1_->active_rotamers_begin(); i != node1_->active_rotamers_end(); ++i){ - min = std::min(min,pairwise_energies_[(*i)*num + idx1] - pairwise_energies_[(*i)*num + idx2]); + min = std::min(min,pairwise_energies_[(*i)*num + idx1] - + pairwise_energies_[(*i)*num + idx2]); } return min; } - throw promod3::Error("Cannot evaluate Node which is unconnected to the edge of interest!"); + String err = "Cannot evaluate Node which is unconnected to the edge of "; + err += "interest!"; + throw promod3::Error(err); } -void Edge::FillActiveEnergies(Real* pairwise_energies) const{ - for(Node::const_iterator i = node1_->active_rotamers_begin(); +void RotamerGraphEdge::FillActiveEnergies(Real* pairwise_energies) const{ + for(RotamerGraphNode::const_iterator i = node1_->active_rotamers_begin(); i != node1_->active_rotamers_end(); ++i){ - for(Node::const_iterator j = node2_->active_rotamers_begin(); + for(RotamerGraphNode::const_iterator j = node2_->active_rotamers_begin(); j != node2_->active_rotamers_end(); ++j){ *pairwise_energies = pairwise_energies_[*i*node2_->GetNumRotamers() + *j]; ++pairwise_energies; @@ -147,13 +158,13 @@ void Edge::FillActiveEnergies(Real* pairwise_energies) const{ } } -Node* Edge::GetOtherNode(Node* node){ +RotamerGraphNode* RotamerGraphEdge::GetOtherNode(RotamerGraphNode* node){ if(node == node1_) return node2_; if(node == node2_) return node1_; throw promod3::Error("Node does not belong to edge!"); } -Node::Node(uint index, +RotamerGraphNode::RotamerGraphNode(uint index, uint num_rotamers, Real* internal_energies, Real* frame_energies): active_(true),index_(index), @@ -169,16 +180,17 @@ Node::Node(uint index, sorted_active_rotamers_.push_back(std::make_pair(i,self_energies_[i])); } - std::sort(sorted_active_rotamers_.begin(), sorted_active_rotamers_.end(),active_rotamer_sorter); + std::sort(sorted_active_rotamers_.begin(), sorted_active_rotamers_.end(), + active_rotamer_sorter); } -Node::~Node(){ +RotamerGraphNode::~RotamerGraphNode(){ delete [] internal_energies_; delete [] frame_energies_; delete [] self_energies_; } -void Node::Reset(){ +void RotamerGraphNode::Reset(){ //clear everything active_rotamers_.clear(); @@ -192,16 +204,20 @@ void Node::Reset(){ self_energies_[i] = internal_energies_[i] + frame_energies_[i]; sorted_active_rotamers_.push_back(std::make_pair(i,self_energies_[i])); } - std::sort(sorted_active_rotamers_.begin(), sorted_active_rotamers_.end(),active_rotamer_sorter); + std::sort(sorted_active_rotamers_.begin(), sorted_active_rotamers_.end(), + active_rotamer_sorter); for(uint i = 0; i < edges_.size(); ++i){ active_edges_.push_back(i); } - for(std::vector<Edge*>::iterator i = edges_.begin(); i != edges_.end(); ++i) (*i)->Reset(); + for(std::vector<RotamerGraphEdge*>::iterator i = edges_.begin(); + i != edges_.end(); ++i){ + (*i)->Reset(); + } } -void Node::Deactivate(){ +void RotamerGraphNode::Deactivate(){ if(!active_) return; active_ = false; //lets first figure out which is the lowest energy @@ -218,35 +234,37 @@ void Node::Deactivate(){ active_rotamers_[0] = index; //deactivate all edges connected to this node - for(iterator i = this->active_edges_begin(); i != this->active_edges_end(); ++i){ - Node* other = edges_[*i]->GetOtherNode(this); + for(iterator i = this->active_edges_begin(); + i != this->active_edges_end(); ++i){ + RotamerGraphNode* other = edges_[*i]->GetOtherNode(this); other->RemoveFromActiveEdges(edges_[*i]); edges_[*i]->Deactivate(); } active_edges_.clear(); } -void Node::AddEdge(Edge* edge){ +void RotamerGraphNode::AddEdge(RotamerGraphEdge* edge){ edges_.push_back(edge); active_edges_.push_back(edges_.size()-1); } -void Node::AddValue(Real value){ +void RotamerGraphNode::AddValue(Real value){ for(std::vector<uint>::iterator i = active_rotamers_.begin(); i != active_rotamers_.end(); ++i){ self_energies_[*i] += value; } } -void Node::AddValue(uint node_idx, Real value){ +void RotamerGraphNode::AddValue(uint node_idx, Real value){ self_energies_[node_idx] += value; } -void Node::RemoveFromActiveEdges(Edge* edge){ +void RotamerGraphNode::RemoveFromActiveEdges(RotamerGraphEdge* edge){ uint counter = 0; for(; counter < edges_.size(); ++counter){ if(edge == edges_[counter]){ - iterator i = std::find(active_edges_.begin(), active_edges_.end(), counter); + iterator i = std::find(active_edges_.begin(), active_edges_.end(), + counter); if(i != active_edges_.end()){ active_edges_.erase(i); } @@ -256,7 +274,7 @@ void Node::RemoveFromActiveEdges(Edge* edge){ } } -bool Node::DEE(){ +bool RotamerGraphNode::DEE(){ if(active_rotamers_.size() <= 1) return false; @@ -268,8 +286,10 @@ bool Node::DEE(){ if(i == rotamer_to_eliminate) continue; sum = self_energies_[sorted_active_rotamers_[rotamer_to_eliminate].first] - self_energies_[sorted_active_rotamers_[i].first]; - for(iterator j = this->active_edges_.begin(); j != this->active_edges_.end(); ++j){ - sum += edges_[*j]->EMinDiff(this,sorted_active_rotamers_[rotamer_to_eliminate].first, + for(iterator j = this->active_edges_.begin(); + j != this->active_edges_.end(); ++j){ + sum += edges_[*j]->EMinDiff(this, + sorted_active_rotamers_[rotamer_to_eliminate].first, sorted_active_rotamers_[i].first); } if(sum >= 0){ @@ -278,7 +298,8 @@ bool Node::DEE(){ active_rotamers_.end(), sorted_active_rotamers_[rotamer_to_eliminate].first); active_rotamers_.erase(it); - sorted_active_rotamers_.erase(sorted_active_rotamers_.begin() + rotamer_to_eliminate); + sorted_active_rotamers_.erase(sorted_active_rotamers_.begin() + + rotamer_to_eliminate); break; } } @@ -288,7 +309,7 @@ bool Node::DEE(){ return something_happened; } -void Node::FillActiveSelfEnergies(Real* self_energies) const{ +void RotamerGraphNode::FillActiveSelfEnergies(Real* self_energies) const{ for(const_iterator i = this->active_rotamers_begin(); i != this->active_rotamers_end(); ++i){ *self_energies = self_energies_[*i]; @@ -296,26 +317,26 @@ void Node::FillActiveSelfEnergies(Real* self_energies) const{ } } -Graph::~Graph(){ - for(std::vector<Node*>::iterator i = nodes_.begin(); +RotamerGraph::~RotamerGraph(){ + for(std::vector<RotamerGraphNode*>::iterator i = nodes_.begin(); i != nodes_.end(); ++i){ delete *i; } - for(std::vector<Edge*>::iterator i = edges_.begin(); + for(std::vector<RotamerGraphEdge*>::iterator i = edges_.begin(); i != edges_.end(); ++i){ delete *i; } } -GraphPtr Graph::CreateFromRRMList( +RotamerGraphPtr RotamerGraph::CreateFromRRMList( const std::vector<RRMRotamerGroupPtr>& rotamer_groups, const std::vector<geom::Transform>& rt_operators) { core::ScopedTimerPtr prof = core::StaticRuntimeProfiler::StartScoped( "Graph::CreateFromRRMList", 2); - GraphPtr graph(new Graph); + RotamerGraphPtr graph(new RotamerGraph); //let's create the nodes graph->nodes_.resize(rotamer_groups.size()); @@ -330,16 +351,20 @@ GraphPtr Graph::CreateFromRRMList( internal_energies[j] = (**i)[j]->GetInternalEnergy(); frame_energies[j] = (**i)[j]->GetFrameEnergy(); } - graph->nodes_[counter] = new Node(counter,(*i)->size(),internal_energies,frame_energies); + graph->nodes_[counter] = new RotamerGraphNode(counter, + (*i)->size(), + internal_energies, + frame_energies); } //let's create the edges and directly add them to the corresponding nodes for(uint i = 0; i < rotamer_groups.size(); ++i){ for(uint j = i+1; j < rotamer_groups.size(); ++j){ - Real* pairwise_energies = rotamer_groups[i]->CalculatePairwiseEnergies(rotamer_groups[j],0.01); + Real* pairwise_energies = + rotamer_groups[i]->CalculatePairwiseEnergies(rotamer_groups[j],0.01); if(pairwise_energies!=NULL){ - Edge* edge = new Edge(graph->edges_.size(), + RotamerGraphEdge* edge = new RotamerGraphEdge(graph->edges_.size(), graph->nodes_[i],graph->nodes_[j], rotamer_groups[i]->size(), rotamer_groups[j]->size(), @@ -363,7 +388,8 @@ GraphPtr Graph::CreateFromRRMList( for(uint i = 0; i < rotamer_groups.size(); ++i){ for(uint j = 0; j < transformed_groups.size(); ++j){ - Real* pairwise_energies = rotamer_groups[i]->CalculatePairwiseEnergies(transformed_groups[j],0.01); + Real* pairwise_energies = + rotamer_groups[i]->CalculatePairwiseEnergies(transformed_groups[j],0.01); if(pairwise_energies!=NULL){ if(i > j){ int i_size = rotamer_groups[i]->size(); @@ -371,7 +397,7 @@ GraphPtr Graph::CreateFromRRMList( Real* inverted_energies = new Real[i_size * j_size]; Invert(pairwise_energies,inverted_energies,i_size,j_size); delete [] pairwise_energies; - Edge* edge = new Edge(graph->edges_.size(), + RotamerGraphEdge* edge = new RotamerGraphEdge(graph->edges_.size(), graph->nodes_[j], graph->nodes_[i], rotamer_groups[j]->size(), @@ -383,7 +409,7 @@ GraphPtr Graph::CreateFromRRMList( continue; } - Edge* edge = new Edge(graph->edges_.size(), + RotamerGraphEdge* edge = new RotamerGraphEdge(graph->edges_.size(), graph->nodes_[i], graph->nodes_[j], rotamer_groups[i]->size(), @@ -401,14 +427,14 @@ GraphPtr Graph::CreateFromRRMList( } -GraphPtr Graph::CreateFromFRMList( +RotamerGraphPtr RotamerGraph::CreateFromFRMList( const std::vector<FRMRotamerGroupPtr>& rotamer_groups, const std::vector<geom::Transform>& rt_operators) { core::ScopedTimerPtr prof = core::StaticRuntimeProfiler::StartScoped( "Graph::CreateFromFRMList", 2); - GraphPtr graph(new Graph); + RotamerGraphPtr graph(new RotamerGraph); //let's create the nodes graph->nodes_.resize(rotamer_groups.size()); @@ -423,15 +449,19 @@ GraphPtr Graph::CreateFromFRMList( internal_energies[j] = (**i)[j]->GetInternalEnergy(); frame_energies[j] = (**i)[j]->GetFrameEnergy(); } - graph->nodes_[counter] = new Node(counter,(*i)->size(),internal_energies,frame_energies); + graph->nodes_[counter] = new RotamerGraphNode(counter, + (*i)->size(), + internal_energies, + frame_energies); } //let's create the edges and directly add them to the corresponding nodes for(uint i = 0; i < rotamer_groups.size(); ++i){ for(uint j = i+1; j < rotamer_groups.size(); ++j){ - Real* pairwise_energies = rotamer_groups[i]->CalculatePairwiseEnergies(rotamer_groups[j],0.01); + Real* pairwise_energies = + rotamer_groups[i]->CalculatePairwiseEnergies(rotamer_groups[j],0.01); if(pairwise_energies != NULL){ - Edge* edge = new Edge(graph->edges_.size(), + RotamerGraphEdge* edge = new RotamerGraphEdge(graph->edges_.size(), graph->nodes_[i], graph->nodes_[j], rotamer_groups[i]->size(), @@ -456,7 +486,8 @@ GraphPtr Graph::CreateFromFRMList( for(uint i = 0; i < rotamer_groups.size(); ++i){ for(uint j = 0; j < transformed_groups.size(); ++j){ - Real* pairwise_energies = rotamer_groups[i]->CalculatePairwiseEnergies(transformed_groups[j],0.01); + Real* pairwise_energies = + rotamer_groups[i]->CalculatePairwiseEnergies(transformed_groups[j],0.01); if(pairwise_energies!=NULL){ if(i > j){ int i_size = rotamer_groups[i]->size(); @@ -464,7 +495,7 @@ GraphPtr Graph::CreateFromFRMList( Real* inverted_energies = new Real[i_size * j_size]; Invert(pairwise_energies,inverted_energies,i_size,j_size); delete [] pairwise_energies; - Edge* edge = new Edge(graph->edges_.size(), + RotamerGraphEdge* edge = new RotamerGraphEdge(graph->edges_.size(), graph->nodes_[j], graph->nodes_[i], rotamer_groups[j]->size(), @@ -476,7 +507,7 @@ GraphPtr Graph::CreateFromFRMList( continue; } - Edge* edge = new Edge(graph->edges_.size(), + RotamerGraphEdge* edge = new RotamerGraphEdge(graph->edges_.size(), graph->nodes_[i], graph->nodes_[j], rotamer_groups[i]->size(), @@ -493,54 +524,58 @@ GraphPtr Graph::CreateFromFRMList( return graph; } -size_t Graph::GetNumActiveRotamers() const{ +size_t RotamerGraph::GetNumActiveRotamers() const{ size_t num_rotamers = 0; - for(std::vector<Node*>::const_iterator i = nodes_.begin(); i != nodes_.end(); ++i){ + for(std::vector<RotamerGraphNode*>::const_iterator i = nodes_.begin(); + i != nodes_.end(); ++i){ if(!(*i)->IsActive()) continue; num_rotamers += (*i)->GetNumActiveRotamers(); } return num_rotamers; } -size_t Graph::GetNumActiveEdges() const{ +size_t RotamerGraph::GetNumActiveEdges() const{ size_t num_edges = 0; - for(std::vector<Edge*>::const_iterator i = edges_.begin(); i != edges_.end(); ++i){ + for(std::vector<RotamerGraphEdge*>::const_iterator i = edges_.begin(); + i != edges_.end(); ++i){ if((*i)->IsActive()) ++num_edges; } return num_edges; } -size_t Graph::GetNumActiveNodes() const{ +size_t RotamerGraph::GetNumActiveNodes() const{ size_t num_nodes = 0; - for(std::vector<Node*>::const_iterator i = nodes_.begin(); i != nodes_.end(); ++i){ + for(std::vector<RotamerGraphNode*>::const_iterator i = nodes_.begin(); + i != nodes_.end(); ++i){ if((*i)->IsActive()) ++num_nodes; } return num_nodes; } -size_t Graph::GetNumRotamers() const{ +size_t RotamerGraph::GetNumRotamers() const{ size_t num_rotamers = 0; - for(std::vector<Node*>::const_iterator i = nodes_.begin(); i != nodes_.end(); ++i){ + for(std::vector<RotamerGraphNode*>::const_iterator i = nodes_.begin(); + i != nodes_.end(); ++i){ num_rotamers += (*i)->GetNumRotamers(); } return num_rotamers; } -void Graph::GetNumRotamers(std::vector<int>& num_rotamers) const{ +void RotamerGraph::GetNumRotamers(std::vector<int>& num_rotamers) const{ num_rotamers.resize(nodes_.size()); for(uint i = 0; i < nodes_.size(); ++i){ num_rotamers[i] = nodes_[i]->GetNumRotamers(); } } -void Graph::GetNumActiveRotamers(std::vector<int>& num_rotamers) const{ +void RotamerGraph::GetNumActiveRotamers(std::vector<int>& num_rotamers) const{ num_rotamers.resize(nodes_.size()); for(uint i = 0; i < nodes_.size(); ++i){ num_rotamers[i] = nodes_[i]->GetNumActiveRotamers(); } } -void Graph::Prune(Real epsilon, bool consider_all_nodes){ +void RotamerGraph::Prune(Real epsilon, bool consider_all_nodes){ //nodes/edges, where something happened std::vector<bool> hot_edges(edges_.size(),true); @@ -571,9 +606,9 @@ void Graph::Prune(Real epsilon, bool consider_all_nodes){ if(nodes_[i]->DEE()){ //let's set the status hot to the edges, as well as other //nodes connected to that node - for(Node::iterator j = nodes_[i]->active_edges_begin(); + for(RotamerGraphNode::iterator j = nodes_[i]->active_edges_begin(); j != nodes_[i]->active_edges_end(); ++j){ - Edge* current_edge = nodes_[i]->GetEdge(*j); + RotamerGraphEdge* current_edge = nodes_[i]->GetEdge(*j); hot_edges[current_edge->GetIndex()] = true; //one of the following connected node is the checked node //itself... will be set to false again... @@ -589,7 +624,7 @@ void Graph::Prune(Real epsilon, bool consider_all_nodes){ } } -void Graph::Reset(){ +void RotamerGraph::Reset(){ for(node_iterator i = nodes_begin(); i != nodes_end(); ++i){ (*i)->Reset(); } @@ -598,7 +633,8 @@ void Graph::Reset(){ } } -std::vector<int> Graph::Solve(uint64_t max_complexity, Real initial_epsilon) { +std::vector<int> RotamerGraph::Solve(uint64_t max_complexity, + Real initial_epsilon) { core::ScopedTimerPtr prof = core::StaticRuntimeProfiler::StartScoped( "Graph::Solve", 2); @@ -621,13 +657,15 @@ std::vector<int> Graph::Solve(uint64_t max_complexity, Real initial_epsilon) { //the indices are based on the active rotamers in the rotamer groups //we have to map that int actual_pos = 0;; - for(node_iterator i = this->nodes_begin(); i != this->nodes_end(); ++i, ++actual_pos){ + for(node_iterator i = this->nodes_begin(); + i != this->nodes_end(); ++i, ++actual_pos){ return_vec[actual_pos] = (*i)->GetOverallIndex(return_vec[actual_pos]); } return return_vec; } -void Graph::Invert(Real* energies, Real* inverted_energies, int i_size, int j_size){ +void RotamerGraph::Invert(Real* energies, Real* inverted_energies, + int i_size, int j_size){ Real* e_ptr = energies; for(int i = 0; i < i_size; ++i){ for(int j = 0; j < j_size; ++j){ diff --git a/sidechain/src/graph.hh b/sidechain/src/graph.hh index 8eb80fae7e4b4fe5b76ba2276c006f57bff0b7bb..31d83a4e7d6d6d5d9cb1071618f8a5cc0f0ac63e 100644 --- a/sidechain/src/graph.hh +++ b/sidechain/src/graph.hh @@ -1,5 +1,5 @@ -#ifndef PROMOD3_GRAPH_COMPONENTS_HH -#define PROMOD3_GRAPH_COMPONENTS_HH +#ifndef PROMOD3_ROTAMER_GRAPH_COMPONENTS_HH +#define PROMOD3_ROTAMER_GRAPH_COMPONENTS_HH #include <math.h> #include <vector> @@ -19,18 +19,19 @@ namespace promod3{ namespace sidechain{ -class Edge; -class Node; -class Graph; +class RotamerGraphEdge; +class RotamerGraphNode; +class RotamerGraph; -typedef boost::shared_ptr<Graph> GraphPtr; +typedef boost::shared_ptr<RotamerGraph> RotamerGraphPtr; -class Edge{ +class RotamerGraphEdge{ public: - Edge(uint index, Node* node1, Node* node2, uint num_i, uint num_j, Real* pairwise_energies_); + RotamerGraphEdge(uint index, RotamerGraphNode* node1, RotamerGraphNode* node2, + uint num_i, uint num_j, Real* pairwise_energies_); - ~Edge(); + ~RotamerGraphEdge(); void Reset(); @@ -44,38 +45,38 @@ public: Real GetPairwiseEnergy(uint i, uint j) const { return pairwise_energies_[i * maxj_ + j]; } - Node* GetNode1() const { return node1_; } + RotamerGraphNode* GetNode1() const { return node1_; } - Node* GetNode2() const { return node2_; } + RotamerGraphNode* GetNode2() const { return node2_; } - Node* GetOtherNode(Node* node); + RotamerGraphNode* GetOtherNode(RotamerGraphNode* node); //resolves the expression min_x(EPair(idx1,x)-EPair(idx2,x)) //in the DEE formalism, where x only considers the active //rotamers from the OTHER node than node_ptr - Real EMinDiff(const Node* node_ptr, uint idx1, uint idx2) const; + Real EMinDiff(const RotamerGraphNode* node_ptr, uint idx1, uint idx2) const; void FillActiveEnergies(Real* pairwise_energies) const; private: uint index_; bool active_; - Node* node1_; - Node* node2_; + RotamerGraphNode* node1_; + RotamerGraphNode* node2_; uint maxi_; uint maxj_; Real* pairwise_energies_; }; -class Node{ +class RotamerGraphNode{ public: - Node(uint index, + RotamerGraphNode(uint index, uint num_rotamers, Real* internal_energies, Real* frame_energies); - ~Node(); + ~RotamerGraphNode(); void Reset(); @@ -85,11 +86,11 @@ public: uint GetIndex() { return index_; } - void AddEdge(Edge* edge); + void AddEdge(RotamerGraphEdge* edge); - Edge* GetEdge(uint index) const { return edges_.at(index); } + RotamerGraphEdge* GetEdge(uint index) const { return edges_.at(index); } - void RemoveFromActiveEdges(Edge* edge); + void RemoveFromActiveEdges(RotamerGraphEdge* edge); void AddValue(Real val); @@ -135,20 +136,20 @@ private: std::vector<uint> active_rotamers_; std::vector<uint> active_edges_; std::vector<std::pair<uint,Real> > sorted_active_rotamers_; - std::vector<Edge*> edges_; + std::vector<RotamerGraphEdge*> edges_; }; -class Graph{ +class RotamerGraph{ public: - ~Graph(); + ~RotamerGraph(); - static GraphPtr CreateFromRRMList(const std::vector<RRMRotamerGroupPtr>& rotamer_groups, - const std::vector<geom::Transform>& rt_operators = std::vector<geom::Transform>()); + static RotamerGraphPtr CreateFromRRMList(const std::vector<RRMRotamerGroupPtr>& rotamer_groups, + const std::vector<geom::Transform>& rt_operators = std::vector<geom::Transform>()); - static GraphPtr CreateFromFRMList(const std::vector<FRMRotamerGroupPtr>& rotamer_groups, - const std::vector<geom::Transform>& rt_operators = std::vector<geom::Transform>()); + static RotamerGraphPtr CreateFromFRMList(const std::vector<FRMRotamerGroupPtr>& rotamer_groups, + const std::vector<geom::Transform>& rt_operators = std::vector<geom::Transform>()); size_t GetNumActiveRotamers() const; @@ -173,15 +174,15 @@ public: std::vector<int> Solve(uint64_t max_complexity = std::numeric_limits<uint64_t>::max(), Real intial_epsilon = 0.02); - Node* GetNode(uint index) { return nodes_.at(index); } + RotamerGraphNode* GetNode(uint index) { return nodes_.at(index); } - Edge* GetEdge(uint index) { return edges_.at(index); } + RotamerGraphEdge* GetEdge(uint index) { return edges_.at(index); } - typedef std::vector<Node*>::const_iterator const_node_iterator; - typedef std::vector<Node*>::iterator node_iterator; + typedef std::vector<RotamerGraphNode*>::const_iterator const_node_iterator; + typedef std::vector<RotamerGraphNode*>::iterator node_iterator; - typedef std::vector<Edge*>::const_iterator const_edge_iterator; - typedef std::vector<Edge*>::iterator edge_iterator; + typedef std::vector<RotamerGraphEdge*>::const_iterator const_edge_iterator; + typedef std::vector<RotamerGraphEdge*>::iterator edge_iterator; node_iterator nodes_begin() { return nodes_.begin(); } node_iterator nodes_end() { return nodes_.end(); } @@ -194,15 +195,15 @@ public: const_edge_iterator edges_end() const{ return edges_.end(); } private: - Graph() { } + RotamerGraph() { } - Graph(const Graph& other) { } + RotamerGraph(const RotamerGraph& other) { } static void Invert(Real* energies, Real* inverted_energies, int i_size , int j_size); - std::vector<Node*> nodes_; - std::vector<Edge*> edges_; + std::vector<RotamerGraphNode*> nodes_; + std::vector<RotamerGraphEdge*> edges_; }; }}//ns diff --git a/sidechain/src/tree.cc b/sidechain/src/tree.cc index 466efbf0752818f19f8c0652da07ae600f94f46f..e348b64c46465d1464309479ff6d05766c43bfa8 100644 --- a/sidechain/src/tree.cc +++ b/sidechain/src/tree.cc @@ -312,7 +312,7 @@ Tree::~Tree(){ //they only contain pointers to the memory adresses just deleted } -uint64_t Tree::GenerateTree(Graph* graph){ +uint64_t Tree::GenerateTree(RotamerGraph* graph){ if(graph->GetNumNodes() != num_nodes_){ std::stringstream ss; ss << "Cannot do the tree decomposition, given graph object differs in number of nodes "; @@ -336,16 +336,16 @@ uint64_t Tree::GenerateTree(Graph* graph){ //to generate the matrix we need a ptr/index mapper... //could be implemented a bit more efficient - std::map<Node*,int> ptr_index_map; + std::map<RotamerGraphNode*,int> ptr_index_map; int counter = 0; - for(Graph::node_iterator i = graph->nodes_begin(); + for(RotamerGraph::node_iterator i = graph->nodes_begin(); i != graph->nodes_end(); ++i){ ptr_index_map[*i] = counter; ++counter; } int a,b; - for(Graph::edge_iterator i = graph->edges_begin(); + for(RotamerGraph::edge_iterator i = graph->edges_begin(); i != graph->edges_end(); ++i){ if(!(*i)->IsActive()) continue; a = ptr_index_map[(*i)->GetNode1()]; @@ -415,7 +415,7 @@ uint64_t Tree::GenerateTree(Graph* graph){ return complexity; } -void Tree::Solve(Graph* graph, std::vector<int>& overall_solution){ +void Tree::Solve(RotamerGraph* graph, std::vector<int>& overall_solution){ if(graph->GetNumNodes() != num_nodes_){ std::stringstream ss; ss << "Cannot get complexity, given graph object differs in number of nodes "; @@ -431,7 +431,7 @@ void Tree::Solve(Graph* graph, std::vector<int>& overall_solution){ //fill up self energies uint pos = 0; - for(Graph::const_node_iterator i = graph->nodes_begin(); + for(RotamerGraph::const_node_iterator i = graph->nodes_begin(); i != graph->nodes_end(); ++i, ++pos){ Real* self_energies = new Real[num_rotamers[pos]]; (*i)->FillActiveSelfEnergies(self_energies); @@ -442,15 +442,15 @@ void Tree::Solve(Graph* graph, std::vector<int>& overall_solution){ for(uint i = 0; i < num_nodes_; ++i){ pairwise_energies_[i].clear(); } - std::map<Node*,int> ptr_index_map; + std::map<RotamerGraphNode*,int> ptr_index_map; int counter = 0; - for(Graph::node_iterator i = graph->nodes_begin(); + for(RotamerGraph::node_iterator i = graph->nodes_begin(); i != graph->nodes_end(); ++i){ ptr_index_map[*i] = counter; ++counter; } int a,b; - for(Graph::edge_iterator i = graph->edges_begin(); + for(RotamerGraph::edge_iterator i = graph->edges_begin(); i != graph->edges_end(); ++i){ if(!(*i)->IsActive()) continue; a = ptr_index_map[(*i)->GetNode1()]; diff --git a/sidechain/src/tree.hh b/sidechain/src/tree.hh index 187d060565e17d8a426d86c23cb7309de894995c..5cd4405e0dcc94ac85db41694bd0eebf037407c6 100644 --- a/sidechain/src/tree.hh +++ b/sidechain/src/tree.hh @@ -24,15 +24,15 @@ typedef boost::shared_ptr<Bag> BagPtr; class Tree{ public: - Tree(Graph* graph): num_nodes_(graph->GetNumNodes()) { } + Tree(RotamerGraph* graph): num_nodes_(graph->GetNumNodes()) { } ~Tree(); size_t GetNumNodes() const { return num_nodes_; } - uint64_t GenerateTree(Graph* graph); + uint64_t GenerateTree(RotamerGraph* graph); - void Solve(Graph* graph, std::vector<int>& overall_solution); + void Solve(RotamerGraph* graph, std::vector<int>& overall_solution); private: uint num_nodes_;