source: src/Fragmentation/KeySetsContainer.cpp@ 336da8

Action_Thermostats Add_AtomRandomPerturbation Add_FitFragmentPartialChargesAction Add_RotateAroundBondAction Add_SelectAtomByNameAction Adding_Graph_to_ChangeBondActions Adding_MD_integration_tests Adding_StructOpt_integration_tests Automaking_mpqc_open AutomationFragmentation_failures Candidate_v1.5.4 Candidate_v1.6.0 Candidate_v1.6.1 ChangeBugEmailaddress ChangingTestPorts ChemicalSpaceEvaluator Combining_Subpackages Debian_Package_split Debian_package_split_molecuildergui_only Disabling_MemDebug Docu_Python_wait EmpiricalPotential_contain_HomologyGraph EmpiricalPotential_contain_HomologyGraph_documentation Enable_parallel_make_install Enhance_userguide Enhanced_StructuralOptimization Enhanced_StructuralOptimization_continued Example_ManyWaysToTranslateAtom Exclude_Hydrogens_annealWithBondGraph FitPartialCharges_GlobalError Fix_ChargeSampling_PBC Fix_ChronosMutex Fix_FitPartialCharges Fix_FitPotential_needs_atomicnumbers Fix_ForceAnnealing Fix_IndependentFragmentGrids Fix_ParseParticles Fix_ParseParticles_split_forward_backward_Actions Fix_StatusMsg Fix_StepWorldTime_single_argument Fix_Verbose_Codepatterns ForceAnnealing_goodresults ForceAnnealing_oldresults ForceAnnealing_tocheck ForceAnnealing_with_BondGraph ForceAnnealing_with_BondGraph_continued ForceAnnealing_with_BondGraph_continued_betteresults ForceAnnealing_with_BondGraph_contraction-expansion GeometryObjects Gui_displays_atomic_force_velocity IndependentFragmentGrids IndependentFragmentGrids_IndividualZeroInstances IndependentFragmentGrids_IntegrationTest IndependentFragmentGrids_Sole_NN_Calculation JobMarket_RobustOnKillsSegFaults JobMarket_StableWorkerPool JobMarket_unresolvable_hostname_fix ODR_violation_mpqc_open PartialCharges_OrthogonalSummation PythonUI_with_named_parameters QtGui_reactivate_TimeChanged_changes Recreated_GuiChecks RotateToPrincipalAxisSystem_UndoRedo StoppableMakroAction Subpackage_CodePatterns Subpackage_JobMarket Subpackage_LinearAlgebra Subpackage_levmar Subpackage_mpqc_open Subpackage_vmg ThirdParty_MPQC_rebuilt_buildsystem TrajectoryDependenant_MaxOrder TremoloParser_IncreasedPrecision TremoloParser_MultipleTimesteps Ubuntu_1604_changes stable
Last change on this file since 336da8 was 34af97, checked in by Frederik Heber <heber@…>, 8 years ago

FIX: KeySetsContainer::insert(KeySetsContainer&) was not fully working.

  • sets were simply appended and according to their index sets.
  • this works for normal index sets but not for forceindexsets, which have the additional (excluded) hydrogens that take not part when counting the bond order.
  • As FragmentationResultsContainer relied on index and forceindex to have the same fragment order, this was no longer true.
  • Added extensive unit test on KeySetsContainer::insert().
  • Property mode set to 100644
File size: 8.8 KB
Line 
1/*
2 * Project: MoleCuilder
3 * Description: creates and alters molecular systems
4 * Copyright (C) 2010-2012 University of Bonn. All rights reserved.
5 * Copyright (C) 2013 Frederik Heber. All rights reserved.
6 *
7 *
8 * This file is part of MoleCuilder.
9 *
10 * MoleCuilder is free software: you can redistribute it and/or modify
11 * it under the terms of the GNU General Public License as published by
12 * the Free Software Foundation, either version 2 of the License, or
13 * (at your option) any later version.
14 *
15 * MoleCuilder is distributed in the hope that it will be useful,
16 * but WITHOUT ANY WARRANTY; without even the implied warranty of
17 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18 * GNU General Public License for more details.
19 *
20 * You should have received a copy of the GNU General Public License
21 * along with MoleCuilder. If not, see <http://www.gnu.org/licenses/>.
22 */
23
24/*
25 * KeySetsContainer.cpp
26 *
27 * Created on: Sep 15, 2011
28 * Author: heber
29 */
30
31// include config.h
32#ifdef HAVE_CONFIG_H
33#include <config.h>
34#endif
35
36#include "CodePatterns/MemDebug.hpp"
37
38#include <fstream>
39#include <sstream>
40
41#include <boost/lexical_cast.hpp>
42#include <boost/tokenizer.hpp>
43
44#include "CodePatterns/Assert.hpp"
45#include "CodePatterns/Log.hpp"
46
47#include "Fragmentation/defs.hpp"
48#include "Fragmentation/helpers.hpp"
49#include "Helpers/defs.hpp"
50#include "KeySetsContainer.hpp"
51
52/** Constructor of KeySetsContainer class.
53 */
54KeySetsContainer::KeySetsContainer() :
55 FragmentCounter(0),
56 Order(0),
57 FragmentsPerOrder(0)
58{};
59
60/** Destructor of KeySetsContainer class.
61 */
62KeySetsContainer::~KeySetsContainer() {
63// for(int i=FragmentCounter;i--;)
64// delete[](KeySets[i]);
65// for(int i=Order;i--;)
66// delete[](OrderSet[i]);
67// delete[](KeySets);
68// delete[](OrderSet);
69// delete[](AtomCounter);
70// delete[](FragmentsPerOrder);
71};
72
73/** Parsing KeySets into array.
74 * \param *name directory with keyset file
75 * \param *ACounter number of atoms per fragment
76 * \param FCounter number of fragments
77 * \return parsing succesful
78 */
79bool KeySetsContainer::ParseKeySets(const std::string path, const std::string name, const int FCounter) {
80 ifstream input;
81 //char *FragmentNumber = NULL;
82 stringstream file;
83 char filename[1023];
84
85 FragmentCounter = FCounter;
86 KeySets.resize(FragmentCounter);
87 AtomCounter.resize(FragmentCounter);
88
89 // open file
90 file << path << "/" << name;
91 input.open(file.str().c_str(), ios::in);
92 if (input.fail()) {
93 LOG(0, endl << "KeySetsContainer::ParseKeySets: Unable to open " << file.str() << ", is the directory correct?");
94 return false;
95 }
96
97 // parse each line, i.e. get the index set per fragment
98 LOG(0, "Parsing key sets.");
99 for(int i=0;(i<FragmentCounter) && (!input.eof());i++) {
100 input.getline(filename, 1023);
101 string line(filename);
102// LOG(2, "DEBUG: line is " << line << ".");
103 std::stringstream set_output;
104 typedef boost::tokenizer<boost::char_separator<char> > tokenizer;
105 boost::char_separator<char> sep(" \t");
106 tokenizer tok(line, sep);
107 for(tokenizer::iterator beg=tok.begin();
108 beg != tok.end();++beg){
109 const int tempvalue = boost::lexical_cast<int>(*beg);
110 KeySets[i].push_back(tempvalue);
111// set_output << " " << KeySets[i].back();
112 }
113// LOG(2, "DEBUG: Scanned keys are '" << set_output.str() << "'.");
114 AtomCounter[i] = KeySets[i].size();
115// {
116// std::stringstream output;
117// FragmentNumber = FixedDigitNumber(FragmentCounter, i);
118// output << FRAGMENTPREFIX << FragmentNumber << "[" << AtomCounter[i] << "]:" << set_output.str();
119// delete[](FragmentNumber);
120// LOG(0, output.str());
121// }
122 }
123 input.close();
124 return true;
125};
126
127/** Parse many body terms, associating each fragment to a certain bond order.
128 * \return parsing succesful
129 */
130bool KeySetsContainer::ParseManyBodyTerms()
131{
132 int Counter;
133
134 LOG(0, "Creating Fragment terms.");
135 // scan through all to determine maximum order
136 Order=0;
137 for(int i=FragmentCounter;i--;) {
138 Counter=0;
139 for(int j=AtomCounter[i];j--;)
140 if (KeySets[i][j] != -1)
141 Counter++;
142 if (Counter > Order)
143 Order = Counter;
144 }
145 LOG(0, "Found Order is " << Order << ".");
146
147 // scan through all to determine fragments per order
148 FragmentsPerOrder.resize(Order);
149// for(int i=Order;i--;)
150// FragmentsPerOrder[i] = 0;
151 for(int i=FragmentCounter;i--;) {
152 Counter=0;
153 for(int j=AtomCounter[i];j--;)
154 if (KeySets[i][j] != -1)
155 Counter++;
156 FragmentsPerOrder[Counter-1]++;
157 }
158 for(int i=0;i<Order;i++)
159 LOG(0, "Found No. of Fragments of Order " << i+1 << " is " << FragmentsPerOrder[i] << ".");
160
161 // scan through all to gather indices to each order set
162 OrderSet.resize(Order);
163 for(int i=Order;i--;)
164 OrderSet[i].resize(FragmentsPerOrder[i]);
165 for(int i=Order;i--;)
166 FragmentsPerOrder[i] = 0;
167 for(int i=FragmentCounter;i--;) {
168 Counter=0;
169 for(int j=AtomCounter[i];j--;)
170 if (KeySets[i][j] != -1)
171 Counter++;
172 OrderSet[Counter-1][FragmentsPerOrder[Counter-1]] = i;
173 FragmentsPerOrder[Counter-1]++;
174 }
175 std::stringstream output;
176 output << "Printing OrderSet: " << std::endl;
177 for(int i=0;i<Order;i++) {
178 for (int j=0;j<FragmentsPerOrder[i];j++) {
179 output << " " << OrderSet[i][j];
180 }
181 output << std::endl;
182 }
183 LOG(0, output.str());
184
185 return true;
186};
187
188/** Compares each entry in \a *SmallerSet if it is containted in \a *GreaterSet.
189 * \param GreaterSet index to greater set
190 * \param SmallerSet index to smaller set
191 * \return true if all keys of SmallerSet contained in GreaterSet
192 */
193bool KeySetsContainer::Contains(const int GreaterSet, const int SmallerSet)
194{
195 bool result = true;
196 bool intermediate;
197 if ((GreaterSet < 0) || (SmallerSet < 0) || (GreaterSet > FragmentCounter) || (SmallerSet > FragmentCounter)) // index out of bounds
198 return false;
199 for(int i=AtomCounter[SmallerSet];i--;) {
200 intermediate = false;
201 for (int j=AtomCounter[GreaterSet];j--;)
202 intermediate = (intermediate || ((KeySets[SmallerSet][i] == KeySets[GreaterSet][j]) || (KeySets[SmallerSet][i] == -1)));
203 result = result && intermediate;
204 }
205
206 return result;
207};
208
209/** Comparison operator for class KeySetsContainer.
210 *
211 * @param other instance to compare to
212 * @return true - both instances are the same in each member variable.
213 */
214bool KeySetsContainer::operator==(const KeySetsContainer &other) const
215{
216 // compare KeySets
217 if (KeySets != other.KeySets)
218 return false;
219 // compare OrderSet
220 if (OrderSet != other.OrderSet)
221 return false;
222 // compare AtomCounter
223 if (AtomCounter != other.AtomCounter)
224 return false;
225 // compare FragmentCounter
226 if (FragmentCounter != other.FragmentCounter)
227 return false;
228 // compare Order
229 if (Order != other.Order)
230 return false;
231 // compare FragmentsPerOrder
232 if (FragmentsPerOrder != other.FragmentsPerOrder)
233 return false;
234
235 return true;
236}
237
238void KeySetsContainer::clear()
239{
240 KeySets.clear();
241 OrderSet.clear();
242 AtomCounter.clear();
243 FragmentCounter = 0;
244 Order = 0;
245 FragmentsPerOrder.clear();
246}
247
248void KeySetsContainer::insert(const KeySetsContainer &other)
249{
250 // append sets and their atom count
251 KeySets.insert(KeySets.end(), other.KeySets.begin(), other.KeySets.end());
252 AtomCounter.insert(AtomCounter.end(), other.AtomCounter.begin(), other.AtomCounter.end());
253
254 // ASSUME: that the fragments inside other are number consecutively
255
256 // find the smallest number
257 int smallest_index = std::numeric_limits<int>::max();
258 for(size_t i=0;i<other.OrderSet.size();++i)
259 for(IntVector::const_iterator iter = other.OrderSet[i].begin();
260 iter != other.OrderSet[i].end(); ++iter)
261 smallest_index = std::min(smallest_index, (*iter));
262
263 // renumber incoming fragments and insert into OrderSets
264 int nr_fragments = 0;
265 for(size_t i=0;i<other.OrderSet.size();++i) {
266 for(IntVector::const_iterator iter = other.OrderSet[i].begin();
267 iter != other.OrderSet[i].end(); ++iter) {
268 if (i >= OrderSet.size() )
269 OrderSet.resize(i+1);
270 OrderSet[i].push_back((*iter)+FragmentCounter-smallest_index);
271 if (i >= FragmentsPerOrder.size())
272 FragmentsPerOrder.resize(i+1);
273 ++(FragmentsPerOrder[i]);
274 }
275 nr_fragments += other.OrderSet[i].size();
276 }
277 Order = std::max(Order, other.Order);
278 // the following assert must not necessarily be true because of cycles:
279 // Order 4 with added cycles is still order 4.
280 //ASSERT( Order == (int)OrderSet.size(),
281 // "KeySetsContainer::insert() - order not max of either container?");
282 FragmentCounter += nr_fragments;
283}
284
285void KeySetsContainer::insert(const IntVector &indices, const size_t order)
286{
287 KeySets.push_back(indices);
288 AtomCounter.push_back(indices.size());
289 if (order > OrderSet.size() )
290 OrderSet.resize(order);
291 OrderSet[order-1].push_back(FragmentCounter);
292 if (order > FragmentsPerOrder.size())
293 FragmentsPerOrder.resize(order);
294 ++(FragmentsPerOrder[order-1]);
295 ++FragmentCounter;
296}
Note: See TracBrowser for help on using the repository browser.