gentoo-overlay/sci-chemistry/ball/files/ball-1.4.1-Added-MAX_PENALTY-option-to-bond-order-assignment.patch

975 lines
39 KiB
Diff

From: Anna Dehof <anne@bioinf.uni-sb.de>
Date: Mon, 26 Mar 2012 19:11:39 +0200
Subject: Added MAX_PENALTY option to bond order assignment
---
include/BALL/DATATYPE/GRAPH/treeWidth.h | 28 +--
.../STRUCTURE/BONDORDERS/FPTBondOrderStrategy.h | 16 +-
include/BALL/STRUCTURE/assignBondOrderProcessor.h | 30 ++-
source/APPLICATIONS/UTILITIES/assign_bond_orders.C | 34 +++-
source/DATATYPE/GRAPH/treeWidth.C | 2 +-
source/STRUCTURE/BONDORDERS/FPTBondOrderStrategy.C | 199 +++++++++-----------
source/STRUCTURE/assignBondOrderProcessor.C | 44 +++--
7 files changed, 186 insertions(+), 167 deletions(-)
diff --git a/include/BALL/DATATYPE/GRAPH/treeWidth.h b/include/BALL/DATATYPE/GRAPH/treeWidth.h
index f7687ef..53ecffc 100644
--- a/include/BALL/DATATYPE/GRAPH/treeWidth.h
+++ b/include/BALL/DATATYPE/GRAPH/treeWidth.h
@@ -100,14 +100,14 @@ namespace BALL
typedef boost::adjacency_list<boost::vecS, boost::vecS, boost::directedS,
boost::property<boost::vertex_bag_content_t, std::set<OriginalVertexType>,
- boost::property<boost::vertex_bag_special_t, OriginalVertexType,
- boost::property<boost::vertex_bag_type_t, int> > >,
+ boost::property<boost::vertex_bag_special_t, OriginalVertexType,
+ boost::property<boost::vertex_bag_type_t, int> > >,
boost::no_property> TreeDecompositionGraph;
typedef typename boost::graph_traits<TreeDecompositionGraph>::vertex_descriptor TreeDecompositionBag;
typedef boost::iterator_property_map<typename std::vector<TreeDecompositionBag>::iterator,
- typename boost::property_map<TreeDecompositionGraph, boost::vertex_index_t>::type>
+ typename boost::property_map<TreeDecompositionGraph, boost::vertex_index_t>::type>
TreeDecompositionParentMap;
typedef boost::graph_as_tree<TreeDecompositionGraph, TreeDecompositionParentMap> TreeDecomposition;
@@ -131,7 +131,7 @@ namespace BALL
class ComponentFilter_
{
public:
- ComponentFilter_(ComponentMap cm, Position i)
+ ComponentFilter_(ComponentMap cm, Position i)
: cm_(cm),
component_(i)
{ }
@@ -144,7 +144,7 @@ namespace BALL
protected:
ComponentMap cm_;
- Position component_;
+ Position component_;
};
/** PropertyWriter for graphviz output.
@@ -167,7 +167,7 @@ namespace BALL
// TODO: would UndirectedGraph suffice here?
MolecularGraph const* input_;
- std::vector<boost::shared_ptr<EditableGraph> > components_;
+ std::vector<boost::shared_ptr<EditableGraph> > components_;
std::vector<boost::shared_ptr<TreeDecomposition> > nice_tree_decompositions_;
std::vector<boost::shared_ptr<TreeDecompositionGraph> > nice_tree_decomposition_graphs_;
@@ -201,7 +201,7 @@ namespace BALL
* @tparam Reducer the reducer which removes a vertex from the graph to reduce it's size
*/
template<class Criterion, class Reducer>
- class GeneralLowerBoundAlgorithm
+ class GeneralLowerBoundAlgorithm
: public UnaryFunctor<UndirectedGraph, Size>
{
public:
@@ -284,7 +284,7 @@ namespace BALL
* @throw BALL::GRAPH::UnconnectedGraphException if called on unconnected graphs
*/
template<class Criterion>
- class GreedyX
+ class GreedyX
: public UnaryFunctor<UndirectedGraph, typename std::pair<
std::vector<boost::graph_traits<typename UndirectedGraph::vertex_descriptor> >, Size> >
{
@@ -296,7 +296,7 @@ namespace BALL
* A criterium for GreedyFillIn which search for a vertex with
* minimum number of additional edges after eliminating
*/
- struct FillInHeuristic
+ struct FillInHeuristic
{
VertexType& operator() (UndirectedGraph& graph);
@@ -323,8 +323,8 @@ namespace BALL
*/
enum SIMPLICIAL_TYPE
{
- NOT_SIMPLICIAL,
- ALMOST_SIMPLICIAL,
+ NOT_SIMPLICIAL,
+ ALMOST_SIMPLICIAL,
IS_SIMPLICIAL
};
@@ -484,7 +484,7 @@ namespace BALL
TreeDecompositionBag buildJoin_(TreeDecompositionBag node, TreeDecompositionBag left,
TreeDecompositionBag right, bool do_forget);
- TreeDecompositionBag buildSingle_(TreeDecompositionBag node, int node_type,
+ TreeDecompositionBag buildSingle_(TreeDecompositionBag node, int node_type,
TreeDecompositionBag child);
TreeDecompositionBag buildLinkage_(TreeDecompositionBag node, TreeDecompositionBag child);
@@ -492,8 +492,8 @@ namespace BALL
TreeDecompositionBag linkWithIntroduceNodes_(TreeDecompositionContent parent_set, TreeDecompositionBag child);
TreeDecompositionBag linkWithForgetNodes_ (TreeDecompositionContent parent_set, TreeDecompositionBag child);
- TreeDecompositionBag branch_(TreeDecompositionBag node, int node_type,
- typename std::vector<TreeDecompositionBag>::iterator begin,
+ TreeDecompositionBag branch_(TreeDecompositionBag node, int node_type,
+ typename std::vector<TreeDecompositionBag>::iterator begin,
typename std::vector<TreeDecompositionBag>::iterator end);
boost::shared_ptr<TreeDecomposition> tree_;
diff --git a/include/BALL/STRUCTURE/BONDORDERS/FPTBondOrderStrategy.h b/include/BALL/STRUCTURE/BONDORDERS/FPTBondOrderStrategy.h
index 55cef31..ba129b0 100644
--- a/include/BALL/STRUCTURE/BONDORDERS/FPTBondOrderStrategy.h
+++ b/include/BALL/STRUCTURE/BONDORDERS/FPTBondOrderStrategy.h
@@ -1117,6 +1117,11 @@ namespace BALL
Size max_heap_size_;
/**
+ * The number of solutions produced so far.
+ */
+ Size num_computed_solutions_;
+
+ /**
* current upperbound. This algorithm will just iterate solutions which are better than this upperbound;
*/
Penalty upper_bound_;
@@ -1248,14 +1253,6 @@ namespace BALL
Size solution_number, Penalty upper_bound = infinite_penalty);
/**
- * Construct a DPBackTrackingCombiner with the given FPTBondOrder and the number of solutions
- * @param bondAssignments vector with the bond assignments. Call #compute before constructing
- * @param solutionNumber the maximum number of solutions you want to backtrack
- */
- DPBackTrackingCombiner_(std::vector<FPTBondOrderAssignment_>& bond_assignments,
- Size solution_number, Penalty upper_bound = infinite_penalty);
-
- /**
* Copy constructor
*/
DPBackTrackingCombiner_(DPBackTrackingCombiner_ const& copy);
@@ -1407,9 +1404,6 @@ namespace BALL
*/
std::vector<std::vector<int> > const* atom_to_block_;
- /// upper bound on the penalty values for which we compute solutions
- Penalty upper_bound_;
-
/**
* A shared pointer to the computing data, so that you can copy this instance without copying
* all the computing data
diff --git a/include/BALL/STRUCTURE/assignBondOrderProcessor.h b/include/BALL/STRUCTURE/assignBondOrderProcessor.h
index f8065d6..0d14af8 100644
--- a/include/BALL/STRUCTURE/assignBondOrderProcessor.h
+++ b/include/BALL/STRUCTURE/assignBondOrderProcessor.h
@@ -183,17 +183,35 @@ namespace BALL
/** the maximal number of solutions to compute
*
- * If set to zero all optimal solutions will be computed.
+ * If set to zero all optimal or all up to \link MAX_PENALTY MAX_PENALTY \endlink
+ * solutions will be computed.
*
- * @see Option::COMPUTE_ALSO_NON_OPTIMAL_SOLUTIONS
+ * @see Option::COMPUTE_ALSO_NON_OPTIMAL_SOLUTIONS
+ * @see Option::MAX_PENALTY
*/
static const char* MAX_NUMBER_OF_SOLUTIONS;
- /** compute also non-optimal solutions but not more than
- * \link MAX_NUMBER_OF_SOLUTIONS MAX_NUMBER_OF_SOLUTIONS \endlink solutions.
+ /** the maximal penalty score allowed
+ *
+ * This option respects option \link MAX_NUMBER_OF_SOLUTIONS MAX_NUMBER_OF_SOLUTIONS \endlink
+ * if specified.
*
+ * If set to -1 this option will be ignored.
+ *
+ * @see Option::COMPUTE_ALSO_NON_OPTIMAL_SOLUTIONS
* @see Option::MAX_NUMBER_OF_SOLUTIONS
*/
+ static const char* MAX_PENALTY;
+
+ /** compute also sub-optimal solutions but not more than
+ * \link MAX_NUMBER_OF_SOLUTIONS MAX_NUMBER_OF_SOLUTIONS \endlink solutions.
+ *
+ * Alternatively \link Option::MAX_PENALTY Option::MAX_PENALTY \endlink
+ * allows to specify a maximal penalty.
+ *
+ * @see Option::MAX_NUMBER_OF_SOLUTIONS
+ * @see Option::MAX_PENALTY
+ */
static const char* COMPUTE_ALSO_NON_OPTIMAL_SOLUTIONS;
/** weighting of bond length penalties wrt valence penalties.
@@ -227,6 +245,7 @@ namespace BALL
static const String INIFile;
static const int MAX_BOND_ORDER;
static const int MAX_NUMBER_OF_SOLUTIONS;
+ static const int MAX_PENALTY;
static const bool COMPUTE_ALSO_NON_OPTIMAL_SOLUTIONS;
static const float BOND_LENGTH_WEIGHTING;
static const bool APPLY_FIRST_SOLUTION;
@@ -716,6 +735,9 @@ namespace BALL
// the max number of solutions to compute
int max_number_of_solutions_;
+ // the max penalty score
+ int max_penalty_;
+
// flag to indicate, whether also non-optimal solutions should be computed
bool compute_also_non_optimal_solutions_;
diff --git a/source/APPLICATIONS/UTILITIES/assign_bond_orders.C b/source/APPLICATIONS/UTILITIES/assign_bond_orders.C
index cefef9b..7053baa 100644
--- a/source/APPLICATIONS/UTILITIES/assign_bond_orders.C
+++ b/source/APPLICATIONS/UTILITIES/assign_bond_orders.C
@@ -16,10 +16,10 @@ using namespace BALL;
int main(int argc, char** argv)
{
- if ((argc < 3) || (argc > 4))
+ if (argc < 3)
// if (argc != 3)
{
- Log << "Usage:" << argv[0] << " <Molecule infile> <Molecule outfile> [num_sol]" << endl;
+ Log << "Usage:" << argv[0] << " <Molecule infile> <Molecule outfile> [--num-sol=n] [--strategy=(fpt|a_star|ilp)] [--max-penalty=m]" << endl;
return 1;
}
@@ -49,24 +49,40 @@ int main(int argc, char** argv)
// set the options:
//
// the solution strategy (A*, ILP or FPT)
- abop.options.set(AssignBondOrderProcessor::Option::ALGORITHM, AssignBondOrderProcessor::Algorithm::A_STAR);
+ abop.options.set(AssignBondOrderProcessor::Option::ALGORITHM, AssignBondOrderProcessor::Algorithm::FPT);
// specify the inifile with the atomic valence penalties
abop.options.set(AssignBondOrderProcessor::Option::INIFile, AssignBondOrderProcessor::Default::INIFile);
// options for considering bond length as well
- abop.options.setReal(AssignBondOrderProcessor::Option::BOND_LENGTH_WEIGHTING, 0);
- abop.options.setReal(AssignBondOrderProcessor::Option::USE_FINE_PENALTY, true);
+// abop.options.setReal(AssignBondOrderProcessor::Option::BOND_LENGTH_WEIGHTING, 0);
+// abop.options.setReal(AssignBondOrderProcessor::Option::USE_FINE_PENALTY, true);
// the combination of the following two options causes the computation of all optimal solutions
abop.options.setInteger(AssignBondOrderProcessor::Option::MAX_NUMBER_OF_SOLUTIONS, 0);
abop.options.setBool(AssignBondOrderProcessor::Option::COMPUTE_ALSO_NON_OPTIMAL_SOLUTIONS, false);
- if (argc == 4)
+ if (argc > 3)
{
- Log << " with compute fixed number (" << argv[3] << ") solutions" << endl;
- abop.options.setInteger(AssignBondOrderProcessor::Option::MAX_NUMBER_OF_SOLUTIONS, String(argv[3]).toInt());
- abop.options.setBool(AssignBondOrderProcessor::Option::COMPUTE_ALSO_NON_OPTIMAL_SOLUTIONS, true);
+ for (Index i=3; i<argc; ++i)
+ {
+ if (String(argv[i]).hasPrefix("--strategy"))
+ {
+ abop.options.set(AssignBondOrderProcessor::Option::ALGORITHM, String(argv[i]).getField(1, "="));
+ }
+ else if (String(argv[i]).hasPrefix("--num-sol"))
+ {
+ Log << " with compute fixed number (" << argv[i] << ") solutions" << endl;
+ abop.options.setInteger(AssignBondOrderProcessor::Option::MAX_NUMBER_OF_SOLUTIONS, String(argv[i]).getField(1, "=").toInt());
+ abop.options.setBool(AssignBondOrderProcessor::Option::COMPUTE_ALSO_NON_OPTIMAL_SOLUTIONS, true);
+ }
+ else if (String(argv[i]).hasPrefix("--max-penalty"))
+ {
+ abop.options.setInteger(AssignBondOrderProcessor::Option::MAX_PENALTY, String(argv[i]).getField(1, "=").toInt());
+ abop.options.setBool(AssignBondOrderProcessor::Option::COMPUTE_ALSO_NON_OPTIMAL_SOLUTIONS, true);
+ Log << " with maximum penalty (" << abop.options.getInteger(AssignBondOrderProcessor::Option::MAX_PENALTY) << ")" << endl;
+ }
+ }
}
// define input and output properties
diff --git a/source/DATATYPE/GRAPH/treeWidth.C b/source/DATATYPE/GRAPH/treeWidth.C
index 7485305..cd0a879 100644
--- a/source/DATATYPE/GRAPH/treeWidth.C
+++ b/source/DATATYPE/GRAPH/treeWidth.C
@@ -12,7 +12,7 @@ namespace BALL
// find all vertices in the current bag
TreeDecompositionContent content = boost::get(boost::vertex_bag_content, *td_, v);
- for (typename TreeDecompositionContent::const_iterator tdc_it = content.begin(); tdc_it != content.end(); ++tdc_it)
+ for (TreeDecompositionContent::const_iterator tdc_it = content.begin(); tdc_it != content.end(); ++tdc_it)
{
TreeWidth<MolecularGraph>::OriginalVertexType ov = *tdc_it;
Atom const* atom = boost::get(boost::vertex_atom_ptr, *original_graph_, ov);
diff --git a/source/STRUCTURE/BONDORDERS/FPTBondOrderStrategy.C b/source/STRUCTURE/BONDORDERS/FPTBondOrderStrategy.C
index 42040e5..107823e 100644
--- a/source/STRUCTURE/BONDORDERS/FPTBondOrderStrategy.C
+++ b/source/STRUCTURE/BONDORDERS/FPTBondOrderStrategy.C
@@ -12,10 +12,6 @@ namespace BALL
const FPTBondOrderStrategy::Penalty FPTBondOrderStrategy::infinite_penalty = 1e5;
const FPTBondOrderStrategy::Valence FPTBondOrderStrategy::max_valence = 8;
- String FPTBondOrderStrategy::Option::UPPER_PENALTY_BOUND = "upper_penalty_bound";
-
- FPTBondOrderStrategy::Penalty FPTBondOrderStrategy::Default::UPPER_PENALTY_BOUND = infinite_penalty;
-
FPTBondOrderStrategy::FPTBondOrderStrategy(AssignBondOrderProcessor* parent)
: BondOrderAssignmentStrategy(parent)
{
@@ -33,15 +29,11 @@ namespace BALL
bool FPTBondOrderStrategy::readOptions(const Options& options)
{
- upper_bound_ = options.get(Option::UPPER_PENALTY_BOUND).toFloat();
-
return true;
}
void FPTBondOrderStrategy::setDefaultOptions()
{
- abop->options.setDefault(FPTBondOrderStrategy::Option::UPPER_PENALTY_BOUND,
- FPTBondOrderStrategy::Default::UPPER_PENALTY_BOUND);
}
void FPTBondOrderStrategy::init()
@@ -64,16 +56,22 @@ namespace BALL
std::vector<boost::shared_ptr<TreeDecomposition> > & ntds = computing_data_->tw->getNiceTreeDecompositions();
std::vector<FPTBondOrderAssignment_*>& bond_assignments = computing_data_->bond_assignments;
+ int max_penalty = abop->max_penalty_;
+ if (max_penalty == -1)
+ max_penalty = infinite_penalty;
+ else
+ max_penalty += 1;
+
bond_assignments.reserve(ntds.size());
for (Position i = 0; i < ntds.size(); ++i)
{
- bond_assignments.push_back(new FPTBondOrderAssignment_(*this, ntds[i], upper_bound_));
- bond_assignments[i]->compute();
+ bond_assignments.push_back(new FPTBondOrderAssignment_(*this, ntds[i], max_penalty));
+ Penalty result = bond_assignments[i]->compute();
}
// initialize backtracking
combiner_ = boost::shared_ptr<DPBackTrackingCombiner_>(
- new DPBackTrackingCombiner_(bond_assignments, abop->max_number_of_solutions_, upper_bound_));
+ new DPBackTrackingCombiner_(bond_assignments, abop->max_number_of_solutions_, max_penalty));
}
boost::shared_ptr<BondOrderAssignment> FPTBondOrderStrategy::computeNextSolution()
@@ -144,7 +142,7 @@ namespace BALL
int penalty = (*penalties_)[start_index + position];
return (penalty < 0) ? infinite_penalty : static_cast<float>(penalty);
- }
+ }
else
{
return infinite_penalty;
@@ -155,20 +153,20 @@ namespace BALL
//* DPConfig_ *
//*****************************************************************************************************
- FPTBondOrderStrategy::DPConfig_::DPConfig_()
- : consumed_valences(0),
+ FPTBondOrderStrategy::DPConfig_::DPConfig_()
+ : consumed_valences(0),
bond_assignments(0)
{
}
- FPTBondOrderStrategy::DPConfig_::DPConfig_(Size atoms, Size bonds)
- : consumed_valences(atoms, 0),
+ FPTBondOrderStrategy::DPConfig_::DPConfig_(Size atoms, Size bonds)
+ : consumed_valences(atoms, 0),
bond_assignments(bonds, 0)
{
}
- FPTBondOrderStrategy::DPConfig_::DPConfig_(std::vector<Valence> const& v, std::vector<BondOrder> const& bo)
- : consumed_valences(v),
+ FPTBondOrderStrategy::DPConfig_::DPConfig_(std::vector<Valence> const& v, std::vector<BondOrder> const& bo)
+ : consumed_valences(v),
bond_assignments(bo)
{
}
@@ -258,12 +256,12 @@ namespace BALL
//* DPTable_ *
//*****************************************************************************************************
- FPTBondOrderStrategy::DPTable_::DPTable_()
+ FPTBondOrderStrategy::DPTable_::DPTable_()
: table()
{
}
- FPTBondOrderStrategy::DPTable_::DPTable_(DPTable_ const& dptab)
+ FPTBondOrderStrategy::DPTable_::DPTable_(DPTable_ const& dptab)
: table(dptab.table)
{
}
@@ -345,7 +343,7 @@ namespace BALL
if (insertion.second)
{
return true;
- }
+ }
else
{
if (insertion.first->second > penalty)
@@ -353,7 +351,7 @@ namespace BALL
insertion.first->second = penalty;
return true;
- }
+ }
else
{
return false;
@@ -365,14 +363,14 @@ namespace BALL
//* AdditionalBagProperties_ *
//*****************************************************************************************************
- FPTBondOrderStrategy::AdditionalBagProperties_::AdditionalBagProperties_()
- : bonds(),
- table(new DPTable_())
+ FPTBondOrderStrategy::AdditionalBagProperties_::AdditionalBagProperties_()
+ : bonds(),
+ table(new DPTable_())
{
}
- FPTBondOrderStrategy::AdditionalBagProperties_::AdditionalBagProperties_(AdditionalBagProperties_ const& copy)
- : bonds(copy.bonds),
+ FPTBondOrderStrategy::AdditionalBagProperties_::AdditionalBagProperties_(AdditionalBagProperties_ const& copy)
+ : bonds(copy.bonds),
table(new DPTable_(*copy.table))
{
}
@@ -384,7 +382,7 @@ namespace BALL
FPTBondOrderStrategy::AdditionalBagProperties_& FPTBondOrderStrategy::AdditionalBagProperties_::operator=(AdditionalBagProperties_ const& copy)
{
- if (&copy != this)
+ if (&copy != this)
{
bonds = copy.bonds;
table = new DPTable_(*copy.table);
@@ -447,7 +445,7 @@ namespace BALL
{
throw Exception::NullPointer(__FILE__, __LINE__);
}
- case TreeWidth<MolecularGraph>::LEAF_BAG:
+ case TreeWidth<MolecularGraph>::LEAF_BAG:
break;
case TreeWidth<MolecularGraph>::END_BAG:
if (boost::num_vertices(*molecule_) == 0 && (bag == root(*ntd_)))
@@ -455,7 +453,7 @@ namespace BALL
// empty molecule -> empty table
bag_properties.table->insert(DPConfig_(), 0);
return bag_properties.table;
- }
+ }
else
{
// else nice tree decomposition is damaged
@@ -603,7 +601,7 @@ namespace BALL
return bonds;
}
- void FPTBondOrderStrategy::FPTBondOrderAssignment_::computeIntroduceBag(TreeDecompositionBag& bag, DPTable_& child,
+ void FPTBondOrderStrategy::FPTBondOrderAssignment_::computeIntroduceBag(TreeDecompositionBag& bag, DPTable_& child,
AdditionalBagProperties_& property)
{
typedef TreeWidth<MolecularGraph>::TreeDecompositionContent TreeDecompositionContent;
@@ -614,7 +612,7 @@ namespace BALL
TreeDecompositionContent vertices = boost::get(boost::vertex_bag_content, *ntd_, bag);
Size num_valences = vertices.size();
-
+
std::vector<MolecularGraphTraits::EdgeType>& bonds = property.bonds;
Size num_bonds = bonds.size();
@@ -654,7 +652,7 @@ namespace BALL
{
// remember the indizes of the new introduced bonds
indices.push_back(vindex);
- }
+ }
else
{
conf.bond_assignments[vindex] = child_entry.first.get().bond_assignments[cindex++];
@@ -662,7 +660,7 @@ namespace BALL
}
// if there are any introduced bonds we have to fill them with each possible value
- if (indices.size() > 0)
+ if (indices.size() > 0)
{
std::vector<int> bond_values(indices.size());
@@ -685,7 +683,7 @@ namespace BALL
++current_index;
}
}
- }
+ }
else
{
table.insert(conf, child_entry.second);
@@ -768,7 +766,7 @@ namespace BALL
return parent_->getPenaltyFor_(forgotten_vertex, forgotten_valence) + child_row.second;
}
- void FPTBondOrderStrategy::FPTBondOrderAssignment_::computeForgetBag(TreeDecompositionBag& bag, DPTable_& child,
+ void FPTBondOrderStrategy::FPTBondOrderAssignment_::computeForgetBag(TreeDecompositionBag& bag, DPTable_& child,
AdditionalBagProperties_& property)
{
typedef TreeWidth<MolecularGraph>::TreeDecompositionContent TreeDecompositionContent;
@@ -795,7 +793,7 @@ namespace BALL
}
}
- void FPTBondOrderStrategy::FPTBondOrderAssignment_::computeRootBag(TreeDecompositionBag& bag,
+ void FPTBondOrderStrategy::FPTBondOrderAssignment_::computeRootBag(TreeDecompositionBag& bag,
DPTable_& child, AdditionalBagProperties_& bag_properties)
{
DPConfig_ empty(0, 0);
@@ -816,7 +814,8 @@ namespace BALL
}
void FPTBondOrderStrategy::FPTBondOrderAssignment_::computeJoinBag(TreeDecompositionBag& bag,
- DPTable_& left_child, DPTable_& right_child, AdditionalBagProperties_& property)
+ DPTable_& left_child, DPTable_& right_child,
+ AdditionalBagProperties_& property)
{
typedef TreeWidth<MolecularGraph>::TreeDecompositionContent TreeDecompositionContent;
@@ -1268,22 +1267,23 @@ namespace BALL
//* DPBackTracking_ *
//*****************************************************************************************************
- FPTBondOrderStrategy::DPBackTracking_::DPBackTracking_(FPTBondOrderAssignment_& bond_assignment, Size max_num_solutions,
- std::vector<MolecularGraphTraits::EdgeType> const& bonds, Penalty upper_bound)
- : bond_assignment_(&bond_assignment),
- current_state_(NULL),
- queue_(),
- max_num_solutions_(max_num_solutions),
+ FPTBondOrderStrategy::DPBackTracking_::DPBackTracking_(FPTBondOrderAssignment_& bond_assignment, Size max_num_solutions,
+ std::vector<MolecularGraphTraits::EdgeType> const& bonds, Penalty upper_bound)
+ : bond_assignment_(&bond_assignment),
+ current_state_(NULL),
+ queue_(),
+ max_num_solutions_(max_num_solutions),
bonds_(&bonds),
- bags_(boost::shared_ptr<std::vector<TreeDecompositionBag> >(new std::vector<TreeDecompositionBag>)),
- max_heap_size_(max_num_solutions),
+ bags_(boost::shared_ptr<std::vector<TreeDecompositionBag> >(new std::vector<TreeDecompositionBag>)),
+ max_heap_size_(1e6),
+ num_computed_solutions_(0),
upper_bound_(upper_bound)
{
// order bags in preorder
if (bond_assignment.ntd_ == NULL)
{
throw Exception::NullPointer(__FILE__, __LINE__);
- }
+ }
else
{
bags_->reserve(boost::num_vertices(bond_assignment.ntd_->_g));
@@ -1303,14 +1303,15 @@ namespace BALL
}
}
- FPTBondOrderStrategy::DPBackTracking_::DPBackTracking_(DPBackTracking_ const& copy)
- : bond_assignment_(copy.bond_assignment_),
- current_state_(NULL),
+ FPTBondOrderStrategy::DPBackTracking_::DPBackTracking_(DPBackTracking_ const& copy)
+ : bond_assignment_(copy.bond_assignment_),
+ current_state_(NULL),
queue_(),
- max_num_solutions_(copy.max_num_solutions_),
- bonds_(copy.bonds_),
- bags_(copy.bags_),
+ max_num_solutions_(copy.max_num_solutions_),
+ bonds_(copy.bonds_),
+ bags_(copy.bags_),
max_heap_size_(copy.max_heap_size_),
+ num_computed_solutions_(copy.num_computed_solutions_),
upper_bound_(copy.upper_bound_)
{
if (copy.current_state_ != NULL)
@@ -1345,6 +1346,7 @@ namespace BALL
bonds_ = copy.bonds_;
bags_ = copy.bags_;
max_heap_size_ = copy.max_heap_size_;
+ num_computed_solutions_ = copy.num_computed_solutions_;
upper_bound_ = copy.upper_bound_;
}
@@ -1365,6 +1367,7 @@ namespace BALL
std::vector<BackTrackingState_*> copy(queue_.begin(), queue_.end());
queue_.clear();
+ num_computed_solutions_ = 0;
for (std::vector<BackTrackingState_*>::iterator iter = copy.begin(); iter != copy.end(); ++iter)
{
@@ -1372,7 +1375,7 @@ namespace BALL
}
}
- bool FPTBondOrderStrategy::DPBackTracking_::StateComparator_::operator() (BackTrackingState_ const * left,
+ bool FPTBondOrderStrategy::DPBackTracking_::StateComparator_::operator() (BackTrackingState_ const * left,
BackTrackingState_ const * right) const
{
return *left < *right;
@@ -1397,7 +1400,7 @@ namespace BALL
if (order < properties.size())
{
return properties[order];
- }
+ }
else
{
throw Exception::IndexOverflow(__FILE__, __LINE__, static_cast<Index>(order), properties.size());
@@ -1426,12 +1429,12 @@ namespace BALL
bool FPTBondOrderStrategy::DPBackTracking_::hasMoreSolutions() const
{
- return (!queue_.empty() && (max_heap_size_ > 0));
+ return (!queue_.empty() && (!max_num_solutions_ || (num_computed_solutions_ <= max_num_solutions_)));
}
void FPTBondOrderStrategy::DPBackTracking_::nextSolution()
{
- if (queue_.empty() || max_heap_size_ == 0)
+ if (queue_.empty() || max_heap_size_ == 0 || ((max_num_solutions_ > 0) && (num_computed_solutions_ > max_num_solutions_)))
{
throw Exception::OutOfRange(__FILE__, __LINE__);
}
@@ -1444,7 +1447,7 @@ namespace BALL
std::multiset<BackTrackingState_*, StateComparator_>::iterator first = queue_.begin();
current_state_ = *first;
queue_.erase(first);
- --max_heap_size_;
+ ++num_computed_solutions_;
TreeDecomposition& tree = *bond_assignment_->ntd_;
@@ -1469,7 +1472,7 @@ namespace BALL
{
TreeDecomposition::children_iterator c_2 = c_i;
++c_2;
- visitJoin(*current_state_, bag, getTable(boost::get(boost::vertex_index, tree, *c_i)),
+ visitJoin(*current_state_, bag, getTable(boost::get(boost::vertex_index, tree, *c_i)),
getTable(boost::get(boost::vertex_index, tree, *c_2)));
break;
}
@@ -1769,7 +1772,7 @@ namespace BALL
if (d < bonds_->size())
{
return d;
- }
+ }
else
{
throw Exception::IndexOverflow(__FILE__, __LINE__, static_cast<Index>(d), bonds_->size());
@@ -1795,7 +1798,7 @@ namespace BALL
bool FPTBondOrderStrategy::DPBackTracking_::isSolutionNeeded(Penalty penalty)
{
- if (max_heap_size_ == 0) {return false;}
+ if (max_heap_size_ == 0 || ((max_num_solutions_ > 0) && (num_computed_solutions_ > max_num_solutions_))) {return false;}
if (queue_.size() >= max_heap_size_) {return Maths::isLess(penalty, upper_bound_);}
return Maths::isLessOrEqual(penalty, upper_bound_);
}
@@ -1846,7 +1849,7 @@ namespace BALL
if (queue_.empty() || max_heap_size_ == 0)
{
return FPTBondOrderStrategy::infinite_penalty;
- }
+ }
else
{
return (*queue_.begin())->assignment.penalty;
@@ -1857,25 +1860,25 @@ namespace BALL
//* DPBackTrackingCombiner_ *
//*****************************************************************************************************
- FPTBondOrderStrategy::DPBackTrackingCombiner_::DPBackTrackingCombiner_(DPBackTrackingCombiner_ const& copy)
- : backtrackers_(copy.deepCopyOfBacktrackers_()),
+ FPTBondOrderStrategy::DPBackTrackingCombiner_::DPBackTrackingCombiner_(DPBackTrackingCombiner_ const& copy)
+ : backtrackers_(copy.deepCopyOfBacktrackers_()),
priority_queue_(copy.priority_queue_),
- component_solutions_(copy.component_solutions_),
+ component_solutions_(copy.component_solutions_),
assignment_(copy.assignment_),
- solution_number_(copy.solution_number_),
- optimum_(copy.optimum_),
+ solution_number_(copy.solution_number_),
+ optimum_(copy.optimum_),
upper_bound_(copy.upper_bound_)
{
}
- FPTBondOrderStrategy::DPBackTrackingCombiner_::DPBackTrackingCombiner_(std::vector<FPTBondOrderAssignment_*>& bond_assignments,
- Size solution_number, Penalty upper_bound)
- : backtrackers_(),
- priority_queue_(),
+ FPTBondOrderStrategy::DPBackTrackingCombiner_::DPBackTrackingCombiner_(std::vector<FPTBondOrderAssignment_*>& bond_assignments,
+ Size solution_number, Penalty upper_bound)
+ : backtrackers_(),
+ priority_queue_(),
component_solutions_(bond_assignments.size()),
- assignment_(),
- solution_number_(solution_number),
- optimum_(FPTBondOrderStrategy::infinite_penalty),
+ assignment_(),
+ solution_number_(solution_number),
+ optimum_(FPTBondOrderStrategy::infinite_penalty),
upper_bound_(upper_bound)
{
if (bond_assignments.empty())
@@ -1897,7 +1900,7 @@ namespace BALL
FPTBondOrderStrategy::EdgeComparator_ ec(&graph);
std::sort(sorted_edges.begin(), sorted_edges.end(), ec);
- for (std::vector<FPTBondOrderAssignment_*>::const_iterator iter = bond_assignments.begin();
+ for (std::vector<FPTBondOrderAssignment_*>::const_iterator iter = bond_assignments.begin();
iter != bond_assignments.end(); ++iter)
{
backtrackers_.push_back(new DPBackTracking_(**iter, solution_number_, sorted_edges, upper_bound_));
@@ -1906,39 +1909,6 @@ namespace BALL
initialize_();
}
- FPTBondOrderStrategy::DPBackTrackingCombiner_::DPBackTrackingCombiner_(std::vector<FPTBondOrderAssignment_>& bond_assignments,
- Size solution_number, Penalty upper_bound)
- : backtrackers_(),
- priority_queue_(),
- component_solutions_(bond_assignments.size()),
- assignment_(),
- solution_number_(solution_number),
- optimum_(FPTBondOrderStrategy::infinite_penalty),
- upper_bound_(upper_bound)
- {
- backtrackers_.reserve(bond_assignments.size());
-
- MolecularGraph& graph = *bond_assignments[0].molecule_;
- sorted_edges.reserve(boost::num_edges(graph));
-
- BGL_FORALL_EDGES(edge_it, graph, MolecularGraph)
- {
- sorted_edges.push_back(edge_it);
- }
-
- // sort bonds - the second vertex could be in false order
- FPTBondOrderStrategy::EdgeComparator_ ec(&graph);
- std::sort(sorted_edges.begin(), sorted_edges.end(), ec);
-
- for (std::vector<FPTBondOrderAssignment_>::iterator iter = bond_assignments.begin();
- iter != bond_assignments.end(); ++iter)
- {
- backtrackers_.push_back(new DPBackTracking_(*iter, solution_number_, sorted_edges, upper_bound_));
- }
-
- initialize_();
- }
-
FPTBondOrderStrategy::DPBackTrackingCombiner_::~DPBackTrackingCombiner_()
{
clear();
@@ -1946,14 +1916,14 @@ namespace BALL
void FPTBondOrderStrategy::DPBackTrackingCombiner_::clear()
{
- for (std::vector<DPBackTracking_*>::iterator iter = backtrackers_.begin();
+ for (std::vector<DPBackTracking_*>::iterator iter = backtrackers_.begin();
iter != backtrackers_.end(); ++iter)
{
delete *iter;
}
}
- FPTBondOrderStrategy::DPBackTrackingCombiner_&
+ FPTBondOrderStrategy::DPBackTrackingCombiner_&
FPTBondOrderStrategy::DPBackTrackingCombiner_::operator = (DPBackTrackingCombiner_ const& copy)
{
if (this != &copy)
@@ -1972,7 +1942,7 @@ namespace BALL
return *this;
}
- std::pair<Size, FPTBondOrderStrategy::Penalty>
+ std::pair<Size, FPTBondOrderStrategy::Penalty>
FPTBondOrderStrategy::DPBackTrackingCombiner_::getNextMinimumBackTracker_() const
{
if (backtrackers_.size() == 1)
@@ -2034,7 +2004,7 @@ namespace BALL
assignment_ = ass;
return;
- }
+ }
else
{
throw Exception::OutOfRange(__FILE__, __LINE__);
@@ -2061,7 +2031,7 @@ namespace BALL
assignment_ = priority_queue_.top();
priority_queue_.pop();
- }
+ }
else
{
throw Exception::OutOfRange(__FILE__, __LINE__);
@@ -2151,7 +2121,7 @@ namespace BALL
return assignment_;
}
- std::vector<FPTBondOrderStrategy::DPBackTracking_*>
+ std::vector<FPTBondOrderStrategy::DPBackTracking_*>
FPTBondOrderStrategy::DPBackTrackingCombiner_::deepCopyOfBacktrackers_() const
{
std::vector<DPBackTracking_*> ts;
@@ -2166,5 +2136,4 @@ namespace BALL
return ts;
}
-
}
diff --git a/source/STRUCTURE/assignBondOrderProcessor.C b/source/STRUCTURE/assignBondOrderProcessor.C
index f8acc4f..9863a54 100644
--- a/source/STRUCTURE/assignBondOrderProcessor.C
+++ b/source/STRUCTURE/assignBondOrderProcessor.C
@@ -108,6 +108,9 @@ namespace BALL
const char* AssignBondOrderProcessor::Option::MAX_NUMBER_OF_SOLUTIONS = "max_number_of_solutions";
const int AssignBondOrderProcessor::Default::MAX_NUMBER_OF_SOLUTIONS = 10;
+ const char* AssignBondOrderProcessor::Option::MAX_PENALTY = "max_penalty_score";
+ const int AssignBondOrderProcessor::Default::MAX_PENALTY = -1;
+
const char* AssignBondOrderProcessor::Option::COMPUTE_ALSO_NON_OPTIMAL_SOLUTIONS = "compute_also_non_optimal_solutions";
const bool AssignBondOrderProcessor::Default::COMPUTE_ALSO_NON_OPTIMAL_SOLUTIONS = false;
@@ -148,6 +151,7 @@ namespace BALL
max_bond_order_(),
alpha_(),
max_number_of_solutions_(),
+ max_penalty_(-1),
compute_also_non_optimal_solutions_(),
add_missing_hydrogens_(),
compute_also_connectivity_(),
@@ -237,13 +241,14 @@ namespace BALL
ret &= s_it->second->readOptions(options);
}
- max_bond_order_ = options.getInteger(Option::MAX_BOND_ORDER);
- alpha_ = options.getReal(Option::BOND_LENGTH_WEIGHTING);
- max_number_of_solutions_ = options.getInteger(Option::MAX_NUMBER_OF_SOLUTIONS);
+ max_bond_order_ = options.getInteger(Option::MAX_BOND_ORDER);
+ alpha_ = options.getReal(Option::BOND_LENGTH_WEIGHTING);
+ max_number_of_solutions_ = options.getInteger(Option::MAX_NUMBER_OF_SOLUTIONS);
+ max_penalty_ = options.getInteger(Option::MAX_PENALTY);
compute_also_non_optimal_solutions_ = options.getBool(Option::COMPUTE_ALSO_NON_OPTIMAL_SOLUTIONS);
- add_missing_hydrogens_ = options.getBool(Option::ADD_HYDROGENS);
+ add_missing_hydrogens_ = options.getBool(Option::ADD_HYDROGENS);
compute_also_connectivity_ = options.getBool(Option::COMPUTE_ALSO_CONNECTIVITY);
- use_fine_penalty_ = options.getBool(Option::USE_FINE_PENALTY);
+ use_fine_penalty_ = options.getBool(Option::USE_FINE_PENALTY);
if (max_bond_order_ <= 0)
{
@@ -257,6 +262,12 @@ namespace BALL
<< " : Error in options! Please check the option Option::MAX_NUMBER_OF_SOLUTIONS." << endl;
ret = false;
}
+ if (max_penalty_ < -1)
+ {
+ Log.error() << __FILE__ << " " << __LINE__
+ << " : Error in options! Please check the option Option::MAX_PENALTY." << endl;
+ ret = false;
+ }
if ((alpha_ < 0) || ((alpha_ > 1)))
{
Log.error() << __FILE__ << " " << __LINE__
@@ -280,14 +291,13 @@ namespace BALL
|| (options.getBool(Option::OVERWRITE_SELECTED_BONDS) == true)
|| (options.getBool(Option::OVERWRITE_SINGLE_BOND_ORDERS) == false)
|| (options.getBool(Option::OVERWRITE_DOUBLE_BOND_ORDERS) == false)
- || (options.getBool(Option::OVERWRITE_TRIPLE_BOND_ORDERS) == false)
- || (options.getInteger(Option::MAX_NUMBER_OF_SOLUTIONS) == 0))
+ || (options.getBool(Option::OVERWRITE_TRIPLE_BOND_ORDERS) == false) )
)
{
Log.error() << __FILE__ << " " << __LINE__
<< " : Error in options! FPT cannot be used with these options." << endl
- << " Switch to solution strategy ASTAR by setting Option::ALGORITHM to Algorithm::ASTAR." << endl;
-
+ << " Consider switch to solution strategy ASTAR by setting Option::ALGORITHM to Algorithm::ASTAR." << endl
+ << " Abort." << endl;
ret = false;
}
@@ -328,6 +338,7 @@ cout << " \t Penalty file " << options[Option::Option::INIFile] << endl;
cout << " \t alpha: " << options[Option::BOND_LENGTH_WEIGHTING] << endl;
cout << " \t max bond order: " << options[Option::MAX_BOND_ORDER] << endl;
cout << " \t max number of solutions " << options[Option::MAX_NUMBER_OF_SOLUTIONS] << endl;
+cout << " \t max penalty " << options[Option::MAX_PENALTY] << endl;
cout << " \t compute also non-optimal solutions: " << options.getBool(Option::COMPUTE_ALSO_NON_OPTIMAL_SOLUTIONS) << endl;
cout << " \t compute also connectivity: " << options.getBool(Option::COMPUTE_ALSO_CONNECTIVITY) << endl;
cout << " \t connectivity cutoff: " << options[Option::CONNECTIVITY_CUTOFF] << endl;
@@ -354,7 +365,7 @@ cout << endl;
// check the options
// What kind of composite do we have?
- // Do we have a molecule? (Nothingelse is allowed)
+ // Do we have a molecule? (Nothing else is allowed)
if (RTTI::isKindOf<Molecule>(ac))
{
// Store the AtomContainer
@@ -517,7 +528,7 @@ cout << "preassignPenaltyClasses_:" << preassignPenaltyClasses_() << " precomput
#ifdef BALL_HAS_LPSOLVE
strategy = strategies_["ILP"].get();
#else
- Log.error() << "Error: BALL was configured without lpsolve support! Try A_STAR instead!" <<
+ Log.error() << "Error: BALL was configured without lpsolve support! Try A_STAR or FPT instead!" <<
__FILE__ << " " << __LINE__<< std::endl;
return Processor::ABORT;
@@ -552,16 +563,20 @@ cout << "preassignPenaltyClasses_:" << preassignPenaltyClasses_() << " precomput
solutions_.push_back(*solution);
// Do we have to find more solutions?
- bool found_another = true;
- bool last_sol_is_optimal = true;
+ bool found_another = true;
+ bool last_sol_is_optimal = true;
+ bool consider_max_penalty = (max_penalty_ > -1);
+ double last_penalty = getTotalPenalty(0);
while ( found_another
&& ((getNumberOfComputedSolutions() < max_number_of_solutions_) || (!max_number_of_solutions_))
&& (last_sol_is_optimal || (compute_also_non_optimal_solutions_))
+ && ((last_penalty <= max_penalty_) || (!consider_max_penalty))
)
{
found_another = computeNextSolution(options.getBool(Option::APPLY_FIRST_SOLUTION));
last_sol_is_optimal &= (fabs(getTotalPenalty(0) - getTotalPenalty(solutions_.size()-1)) < 1.e-4);
+ last_penalty = getTotalPenalty(solutions_.size()-1);
}
#if defined DEBUG_TIMER
timer_.stop();
@@ -1206,6 +1221,9 @@ cout << " ~~~~~~~~ added hydrogen dump ~~~~~~~~~~~~~~~~" << endl;
options.setDefaultInteger(AssignBondOrderProcessor::Option::MAX_NUMBER_OF_SOLUTIONS,
AssignBondOrderProcessor::Default::MAX_NUMBER_OF_SOLUTIONS);
+ options.setDefaultInteger(AssignBondOrderProcessor::Option::MAX_PENALTY,
+ AssignBondOrderProcessor::Default::MAX_PENALTY);
+
options.setDefaultBool(AssignBondOrderProcessor::Option::COMPUTE_ALSO_NON_OPTIMAL_SOLUTIONS,
AssignBondOrderProcessor::Default::COMPUTE_ALSO_NON_OPTIMAL_SOLUTIONS);