/* * Project: MoleCuilder * Description: creates and alters molecular systems * Copyright (C) 2012 University of Bonn. All rights reserved. * Copyright (C) 2013 Frederik Heber. All rights reserved. * Please see the COPYING file or "Copyright notice" in builder.cpp for details. * * * This file is part of MoleCuilder. * * MoleCuilder 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, either version 2 of the License, or * (at your option) any later version. * * MoleCuilder 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 MoleCuilder. If not, see . */ /* * Extractors.cpp * * Created on: 15.10.2012 * Author: heber */ // include config.h #ifdef HAVE_CONFIG_H #include #endif #include "CodePatterns/MemDebug.hpp" #include #include #include #include #include #include #include "CodePatterns/Assert.hpp" #include "CodePatterns/IteratorAdaptors.hpp" #include "CodePatterns/Log.hpp" #include "CodePatterns/toString.hpp" #include "LinearAlgebra/Vector.hpp" #include "FunctionApproximation/Extractors.hpp" #include "FunctionApproximation/FunctionArgument.hpp" using namespace boost::assign; FunctionModel::arguments_t Extractors::gatherAllDistanceArguments( const Fragment::positions_t& positions, const Fragment::charges_t& charges, const size_t globalid) { FunctionModel::arguments_t result; // go through current configuration and gather all other distances Fragment::positions_t::const_iterator firstpositer = positions.begin(); for (;firstpositer != positions.end(); ++firstpositer) { Fragment::positions_t::const_iterator secondpositer = positions.begin();//firstpositer; for (; secondpositer != positions.end(); ++secondpositer) { if (firstpositer == secondpositer) continue; argument_t arg; const Vector firsttemp((*firstpositer)[0],(*firstpositer)[1],(*firstpositer)[2]); const Vector secondtemp((*secondpositer)[0],(*secondpositer)[1],(*secondpositer)[2]); arg.distance = firsttemp.distance(secondtemp); arg.types = std::make_pair( charges[ std::distance(positions.begin(), firstpositer) ], charges[ std::distance(positions.begin(), secondpositer) ] ); arg.indices = std::make_pair( std::distance( positions.begin(), firstpositer), std::distance( positions.begin(), secondpositer) ); arg.globalid = globalid; result.push_back(arg); } } return result; } FunctionModel::arguments_t Extractors::gatherAllSymmetricDistanceArguments( const Fragment::positions_t& positions, const Fragment::charges_t& charges, const size_t globalid) { FunctionModel::arguments_t result; // go through current configuration and gather all other distances Fragment::positions_t::const_iterator firstpositer = positions.begin(); for (;firstpositer != positions.end(); ++firstpositer) { Fragment::positions_t::const_iterator secondpositer = firstpositer; for (; secondpositer != positions.end(); ++secondpositer) { if (firstpositer == secondpositer) continue; argument_t arg; const Vector firsttemp((*firstpositer)[0],(*firstpositer)[1],(*firstpositer)[2]); const Vector secondtemp((*secondpositer)[0],(*secondpositer)[1],(*secondpositer)[2]); arg.distance = firsttemp.distance(secondtemp); arg.types = std::make_pair( charges[ std::distance(positions.begin(), firstpositer) ], charges[ std::distance(positions.begin(), secondpositer) ] ); arg.indices = std::make_pair( std::distance( positions.begin(), firstpositer), std::distance( positions.begin(), secondpositer) ); arg.globalid = globalid; LOG(3, "DEBUG: Created argument " << arg << "."); result.push_back(arg); } } return result; } Fragment::positions_t Extractors::_detail::gatherPositionsFromTargets( const Fragment::positions_t& positions, const Fragment::charges_t& charges, const chargeiters_t &targets ) { Fragment::positions_t filtered_positions; for (chargeiters_t::const_iterator firstpairiter = targets.begin(); firstpairiter != targets.end(); ++firstpairiter) { Fragment::positions_t::const_iterator positer = positions.begin(); const size_t steps = std::distance(charges.begin(), *firstpairiter); std::advance(positer, steps); filtered_positions.push_back(*positer); } return filtered_positions; } FunctionModel::arguments_t Extractors::_detail::gatherDistancesFromTargets( const Fragment::positions_t& positions, const Fragment::charges_t& charges, const chargeiters_t &targets, const size_t globalid ) { Fragment::positions_t filtered_positions; Fragment::charges_t filtered_charges; for (chargeiters_t::const_iterator firstpairiter = targets.begin(); firstpairiter != targets.end(); ++firstpairiter) { Fragment::positions_t::const_iterator positer = positions.begin(); const size_t steps = std::distance(charges.begin(), *firstpairiter); std::advance(positer, steps); filtered_positions.push_back(*positer); filtered_charges.push_back(**firstpairiter); } return Extractors::gatherAllSymmetricDistanceArguments( filtered_positions, filtered_charges, globalid); } Extractors::elementcounts_t Extractors::_detail::getElementCounts( const Fragment::charges_t elements ) { elementcounts_t elementcounts; for (Fragment::charges_t::const_iterator elementiter = elements.begin(); elementiter != elements.end(); ++elementiter) { // insert new element std::pair< elementcounts_t::iterator, bool> inserter = elementcounts.insert( std::make_pair( *elementiter, 1) ); // if already present, just increase its count if (!inserter.second) ++(inserter.first->second); } return elementcounts; } Extractors::elementtargets_t Extractors::_detail::convertElementcountsToTargets( const Fragment::charges_t &charges, const elementcounts_t &elementcounts ) { elementtargets_t elementtargets; for (elementcounts_t::const_iterator countiter = elementcounts.begin(); countiter != elementcounts.end(); ++countiter) { chargeiter_t chargeiter = charges.begin(); const element_t &element = countiter->first; const count_t &count = countiter->second; for (count_t i = 0; i < count; ++i) { chargeiter_t tempiter = std::find(chargeiter, charges.end(), element); if (tempiter != charges.end()) { // try to insert new list std::pair< elementtargets_t::iterator, bool> inserter = elementtargets.insert( std::make_pair( countiter->first, chargeiters_t(1, tempiter)) ); // if already present, append to it if (!inserter.second) { inserter.first->second.push_back(tempiter); } else { // if created, increase vector's reserve to known size inserter.first->second.reserve(countiter->second); } // search from this element onwards then chargeiter = ++tempiter; } else { ELOG(1, "Could not find desired number " << count << " of element " << element << " in fragment with " << charges << "."); return Extractors::elementtargets_t(); } } } return elementtargets; } Extractors::elementtargets_t Extractors::_detail::convertChargesToTargetMap( const Fragment::charges_t& charges, Fragment::charges_t elements ) { // place each charge into a map elementtargets_t completeelementtargets; for (chargeiter_t chargeiter = charges.begin(); chargeiter != charges.end(); ++chargeiter) { std::pair< elementtargets_t::iterator, bool> inserter = completeelementtargets.insert( std::make_pair( *chargeiter, chargeiters_t(1, chargeiter)) ); // if already present, append to it if (!inserter.second) { inserter.first->second.push_back(chargeiter); } } // pick out desired charges only std::sort(elements.begin(), elements.end()); Fragment::charges_t::iterator eraseiter = std::unique(elements.begin(), elements.end()); elements.erase(eraseiter, elements.end()); elementtargets_t elementtargets; for(Fragment::charges_t::const_iterator iter = elements.begin(); iter != elements.end(); ++iter) { elementtargets_t::const_iterator finditer = completeelementtargets.find(*iter); ASSERT( finditer != completeelementtargets.end(), "Extractors::_detail::convertChargesToTargetMap() - no element "+toString(*iter)+" present?"); std::pair< elementtargets_t::iterator, bool> inserter = elementtargets.insert( std::make_pair( finditer->first, finditer->second) ); ASSERT( inserter.second, "Extractors::_detail::convertChargesToTargetMap() - key twice?"); } return elementtargets; } Extractors::chargeiters_t Extractors::_detail::realignElementtargets( const elementtargets_t &elementtargets, const Fragment::charges_t elements, const elementcounts_t &elementcounts ) { chargeiters_t targets; elementcounts_t counts; // how many chargeiters of this element have been used if (!elements.empty()) { // skip if no elements given targets.reserve(elements.size()); for (Fragment::charges_t::const_iterator elementiter = elements.begin(); elementiter != elements.end(); ++elementiter) { const element_t &element = *elementiter; count_t &count = counts[element]; // if not present, std::map creates instances with default of 0 #ifndef NDEBUG { elementcounts_t::const_iterator testiter = elementcounts.find(element); ASSERT( (testiter != elementcounts.end()) && (count < testiter->second), "Extractors::_detail::realignElementTargets() - we want to use more chargeiters for element " +toString(element)+" than we counted initially."); } #endif elementtargets_t::const_iterator targetiter = elementtargets.find(element); if (targetiter != elementtargets.end()) { const chargeiters_t &chargeiters = targetiter->second; const chargeiter_t &chargeiter = chargeiters[count++]; targets.push_back(chargeiter); } } } return targets; } FunctionModel::arguments_t Extractors::gatherAllDistancesFromFragment( const Fragment::positions_t& positions, const Fragment::charges_t& charges, const Fragment::charges_t elements, const size_t globalid ) { /// The main problem here is that we have to know how many same /// elements (but different atoms!) we are required to find. Hence, /// we first have to count same elements, then get different targets /// for each and then associated them in correct order back again. // 0. if no elements given, we return empty arguments if (elements.empty()) return FunctionModel::arguments_t(); // 1. we have to place each charge into a map as unique chargeiter, i.e. map elementtargets_t elementtargets = Extractors::_detail::convertChargesToTargetMap( charges, elements); // 2. now we have to combine elementcounts out of elementtargets per desired element // in a combinatorial fashion targets_per_combination_t combinations = Extractors::_detail::CombineChargesAndTargets( elements, elementtargets); // 3. finally, convert chargeiters into argument list FunctionModel::arguments_t args = Extractors::_detail::convertTargetsToArguments( positions, charges, combinations, globalid); return args; } Extractors::targets_per_combination_t Extractors::_detail::CombineChargesAndTargets( const Fragment::charges_t& elements, const elementtargets_t& elementtargets ) { // recursively create all correct combinations of targets targets_per_combination_t combinations; chargeiters_t currenttargets; boost::function addFunction = boost::bind(&targets_per_combination_t::push_back, boost::ref(combinations), _1); pickLastElementAsTarget(elements, elementtargets, currenttargets, addFunction); return combinations; } const Fragment::position_t& getPositionToChargeIter( const Fragment::positions_t& positions, const Fragment::charges_t& charges, const Extractors::chargeiter_t &iter ) { Fragment::positions_t::const_iterator positer = positions.begin(); std::advance(positer, std::distance(charges.begin(), iter)); const Fragment::position_t &position = *positer; return position; } FunctionModel::arguments_t Extractors::_detail::convertTargetsToArguments( const Fragment::positions_t& positions, const Fragment::charges_t& charges, const targets_per_combination_t combinations, const size_t globalid ) { FunctionModel::arguments_t args; // create arguments from each combination. We cannot use // gatherallSymmetricDistanceArguments() because it would not create the // correct indices. for (targets_per_combination_t::const_iterator iter = combinations.begin(); iter != combinations.end(); ++iter) { for(chargeiters_t::const_iterator firstiter = iter->begin(); firstiter != iter->end(); ++firstiter) { const Fragment::position_t &firstpos = getPositionToChargeIter(positions, charges, *firstiter); const Vector firsttemp(firstpos[0],firstpos[1],firstpos[2]); for(chargeiters_t::const_iterator seconditer = firstiter; seconditer != iter->end(); ++seconditer) { if (firstiter == seconditer) continue; const Fragment::position_t &secondpos = getPositionToChargeIter(positions, charges, *seconditer); const Vector secondtemp(secondpos[0],secondpos[1],secondpos[2]); argument_t arg; arg.distance = firsttemp.distance(secondtemp); arg.indices.first = std::distance(charges.begin(), *firstiter); arg.indices.second = std::distance(charges.begin(), *seconditer); arg.types.first = **firstiter; arg.types.second = **seconditer; args.push_back( arg ); } } } return args; } void Extractors::_detail::pickLastElementAsTarget( Fragment::charges_t elements, elementtargets_t elementtargets, chargeiters_t ¤ttargets, boost::function &addFunction ) { // get last element from charges ASSERT( !elements.empty(), "Extractors::_detail::pickLastElementAsTarget() - no elements given to pick targets for."); const Fragment::charge_t charge = elements.back(); elements.pop_back(); elementtargets_t::iterator iter = elementtargets.find(charge); if (iter == elementtargets.end()) return; bool NotEmpty = !iter->second.empty(); while (NotEmpty) { // get last target from the vector of chargeiters chargeiter_t target = iter->second.back(); iter->second.pop_back(); // remove this key if empty if (iter->second.empty()) { elementtargets.erase(iter); NotEmpty = false; } currenttargets.push_back(target); if (elements.empty()) { // call add function { std::stringstream targetstream; BOOST_FOREACH( chargeiter_t target, currenttargets ) { targetstream << " " << *target; } LOG(3, "DEBUG: Adding set" << targetstream.str() << "."); } addFunction(currenttargets); } else { // if not, call us recursively pickLastElementAsTarget(elements, elementtargets, currenttargets, addFunction); } // pop last in currenset again currenttargets.pop_back(); } } Extractors::chargeiters_t Extractors::_detail::gatherTargetsFromFragment( const Fragment::charges_t& charges, const Fragment::charges_t elements ) { /// The main problem here is that we have to know how many same /// elements (but different atoms!) we are required to find. Hence, /// we first have to count same elements, then get different targets /// for each and then associated them in correct order back again. // 1. we have to make elements unique with counts, hence convert to map elementcounts_t elementcounts = Extractors::_detail::getElementCounts(elements); // 2. then for each element we need as many targets (chargeiters) as counts elementtargets_t elementtargets = Extractors::_detail::convertElementcountsToTargets(charges, elementcounts); // 3. we go again through elements and use one found target for each count // in that order chargeiters_t targets = Extractors::_detail::realignElementtargets(elementtargets, elements, elementcounts); #ifndef NDEBUG // check all for debugging for (chargeiters_t::const_iterator chargeiter = targets.begin(); chargeiter != targets.end(); ++chargeiter) ASSERT( *chargeiter != charges.end(), "Extractors::gatherTargetsFromFragment() - we have not found enough targets?!"); #endif return targets; } Fragment::positions_t Extractors::gatherPositionsFromFragment( const Fragment::positions_t positions, const Fragment::charges_t charges, const Fragment::charges_t& elements ) { // 1.-3. gather correct charge positions chargeiters_t targets = Extractors::_detail::gatherTargetsFromFragment(charges, elements); // 4. convert position_t to Vector return Extractors::_detail::gatherPositionsFromTargets( positions, charges, targets); } FunctionModel::arguments_t Extractors::gatherDistancesFromFragment( const Fragment::positions_t positions, const Fragment::charges_t charges, const Fragment::charges_t& elements, const size_t globalid ) { // 1.-3. gather correct charge positions chargeiters_t targets = Extractors::_detail::gatherTargetsFromFragment(charges, elements); // 4. convert position_t to Vector return Extractors::_detail::gatherDistancesFromTargets( positions, charges, targets, globalid); } FunctionModel::list_of_arguments_t Extractors::reorderArgumentsByIncreasingDistance( const FunctionModel::list_of_arguments_t &listargs ) { FunctionModel::list_of_arguments_t returnargs; for (FunctionModel::list_of_arguments_t::const_iterator iter = listargs.begin(); iter != listargs.end(); ++iter) { const FunctionModel::arguments_t &args = *iter; FunctionModel::arguments_t sortedargs(args); std::sort(sortedargs.begin(), sortedargs.end(), argument_t::DistanceComparator); returnargs.push_back(sortedargs); } return returnargs; } struct ParticleTypesComparator { bool operator()(const argument_t::types_t &a, const argument_t::types_t &b) { if (a.first < a.second) { if (b.first < b.second) { if (a.first < b.first) return true; else if (a.first > b.first) return false; else return (a.second < b.second); } else { if (a.first < b.second) return true; else if (a.first > b.second) return false; else return (a.second < b.first); } } else { if (b.first < b.second) { if (a.second < b.first) return true; else if (a.second > b.first) return false; else return (a.first < b.second); } else { if (a.second < b.second) return true; else if (a.second > b.second) return false; else return (a.first < b.first); } } } }; std::ostream& operator<<(std::ostream &out, const argument_t::types_t &a) { out << "[" << a.first << "," << a.second << "]"; return out; } FunctionModel::list_of_arguments_t Extractors::reorderArgumentsByParticleTypes( const FunctionModel::list_of_arguments_t &listargs, const ParticleTypes_t &_types ) { FunctionModel::list_of_arguments_t returnargs; for (FunctionModel::list_of_arguments_t::const_iterator iter = listargs.begin(); iter != listargs.end(); ++iter) { const FunctionModel::arguments_t &args = *iter; /// We place all arguments into multimap according to particle type pair. // here, we need a special comparator such that types in key pair are always // properly ordered. typedef std::multimap< argument_t::types_t, argument_t, ParticleTypesComparator> TypePair_Argument_Map_t; TypePair_Argument_Map_t argument_map; for(FunctionModel::arguments_t::const_iterator iter = args.begin(); iter != args.end(); ++iter) { argument_map.insert( std::make_pair(iter->types, *iter) ); } LOG(4, "DEBUG: particle_type map is " << argument_map << "."); /// Then, we create the desired unique keys typedef std::vector UniqueTypes_t; UniqueTypes_t UniqueTypes; for (ParticleTypes_t::const_iterator firstiter = _types.begin(); firstiter != _types.end(); ++firstiter) { for (ParticleTypes_t::const_iterator seconditer = firstiter; seconditer != _types.end(); ++seconditer) { if (seconditer == firstiter) continue; UniqueTypes.push_back( std::make_pair(*firstiter, *seconditer) ); } } LOG(4, "DEBUG: Created unique types as keys " << UniqueTypes << "."); /// Finally, we use the unique key list to pick corresponding arguments from the map FunctionModel::arguments_t sortedargs; sortedargs.reserve(args.size()); while (!argument_map.empty()) { // note that particle_types_t may be flipped, i.e. 1,8 is equal to 8,1, but we // must maintain the correct order in indices in accordance with the order // in _types, i.e. 1,8,1 must match with e.g. ids 1,0,2 where 1 has type 1, // 0 has type 8, and 2 has type 2. // In other words: We do not want to flip/modify arguments such that they match // with the specific type pair we seek but then this comes at the price that we // have flip indices when the types in a pair are flipped. typedef std::vector indices_t; //!> here, we gather the indices as we discover them indices_t indices; indices.resize(_types.size(), (size_t)-1); // these are two iterators that create index pairs in the same way as we have // created type pairs. If a -1 is still present in indices, then the index is // still arbitrary but is then set by the next found index indices_t::iterator firstindex = indices.begin(); indices_t::iterator secondindex = firstindex+1; //!> here, we gather the current bunch of arguments as we find them FunctionModel::arguments_t argumentbunch; argumentbunch.reserve(UniqueTypes.size()); for (UniqueTypes_t::const_iterator typeiter = UniqueTypes.begin(); typeiter != UniqueTypes.end(); ++typeiter) { // have all arguments to same type pair as list within the found range std::pair< TypePair_Argument_Map_t::iterator, TypePair_Argument_Map_t::iterator> range_t = argument_map.equal_range(*typeiter); LOG(4, "DEBUG: Set of arguments to current key [" << typeiter->first << "," << typeiter->second << "] is " << std::list( MapValueIterator(range_t.first), MapValueIterator(range_t.second) ) << "."); // the first key is always easy and is pivot which the rest has to be associated to if (typeiter == UniqueTypes.begin()) { const argument_t & arg = range_t.first->second; if ((typeiter->first == arg.types.first) && (typeiter->second == arg.types.second)) { // store in correct order *firstindex = arg.indices.first; *secondindex = arg.indices.second; } else { // store in flipped order *firstindex = arg.indices.second; *secondindex = arg.indices.first; } argumentbunch.push_back(arg); argument_map.erase(range_t.first); LOG(4, "DEBUG: Gathered first argument " << arg << "."); } else { // go through the range and pick the first argument matching the index constraints for (TypePair_Argument_Map_t::iterator argiter = range_t.first; argiter != range_t.second; ++argiter) { // seconditer may be -1 still const argument_t &arg = argiter->second; if (arg.indices.first == *firstindex) { if ((arg.indices.second == *secondindex) || (*secondindex == (size_t)-1)) { if (*secondindex == (size_t)-1) *secondindex = arg.indices.second; argumentbunch.push_back(arg); argument_map.erase(argiter); LOG(4, "DEBUG: Gathered another argument " << arg << "."); break; } } else if ((arg.indices.first == *secondindex) || (*secondindex == (size_t)-1)) { if (arg.indices.second == *firstindex) { if (*secondindex == (size_t)-1) *secondindex = arg.indices.first; argumentbunch.push_back(arg); argument_map.erase(argiter); LOG(4, "DEBUG: Gathered another (flipped) argument " << arg << "."); break; } } } } // move along in indices and check bounds ++secondindex; if (secondindex == indices.end()) { ++firstindex; if (firstindex != indices.end()-1) secondindex = firstindex+1; } } ASSERT( (firstindex == indices.end()-1) && (secondindex == indices.end()), "Extractors::reorderArgumentsByParticleTypes() - we have not gathered enough arguments."); ASSERT( argumentbunch.size() == UniqueTypes.size(), "Extractors::reorderArgumentsByParticleTypes() - we have not gathered enough arguments."); // place bunch of arguments in return args LOG(3, "DEBUG: Given types " << _types << " and found indices " << indices << "."); LOG(3, "DEBUG: Final bunch of arguments is " << argumentbunch << "."); sortedargs.insert(sortedargs.end(), argumentbunch.begin(), argumentbunch.end()); } returnargs.push_back(sortedargs); } return returnargs; } FunctionModel::list_of_arguments_t Extractors::filterArgumentsByParticleTypes( const FunctionModel::arguments_t &args, const ParticleTypes_t &_types ) { typedef std::list< argument_t > ListArguments_t; ListArguments_t availableList(args.begin(), args.end()); LOG(2, "DEBUG: Initial list of args is " << args << "."); // TODO: fill a lookup map such that we don't have O(M^3) scaling, if M is number // of types (and we always must have M(M-1)/2 args) but O(M^2 log(M)). However, as // M is very small (<=3), this is not necessary fruitful now. // typedef ParticleTypes_t firsttype; // typedef ParticleTypes_t secondtype; // typedef std::map< firsttype, std::map< secondtype, boost::ref(args) > > ArgsLookup_t; // ArgsLookup_t ArgsLookup; // basically, we have two choose any two pairs out of types but only those // where the first is less than the latter. Hence, we start the second // iterator at the current position of the first one and skip the equal case. FunctionModel::arguments_t allargs; allargs.reserve(args.size()); for (ParticleTypes_t::const_iterator firstiter = _types.begin(); firstiter != _types.end(); ++firstiter) { for (ParticleTypes_t::const_iterator seconditer = firstiter; seconditer != _types.end(); ++seconditer) { if (seconditer == firstiter) continue; LOG(3, "DEBUG: Looking for (" << *firstiter << "," << *seconditer << ") in all args."); // search the right one in _args (we might allow switching places of // firstiter and seconditer, as distance is symmetric). ListArguments_t::iterator iter = availableList.begin(); while (iter != availableList.end()) { LOG(4, "DEBUG: Current args is " << *iter << "."); if ((iter->types.first == *firstiter) && (iter->types.second == *seconditer)) { allargs.push_back( *iter ); iter = availableList.erase(iter); LOG(4, "DEBUG: Accepted argument."); } else if ((iter->types.first == *seconditer) && (iter->types.second == *firstiter)) { allargs.push_back( *iter ); iter = availableList.erase(iter); LOG(4, "DEBUG: Accepted (flipped) argument."); } else { ++iter; LOG(4, "DEBUG: Rejected argument."); } } } } LOG(2, "DEBUG: Final list of args is " << allargs << "."); // first, we bring together tuples of distances that belong together FunctionModel::list_of_arguments_t singlelist_allargs; singlelist_allargs.push_back(allargs); FunctionModel::list_of_arguments_t sortedargs = reorderArgumentsByParticleTypes(singlelist_allargs, _types); ASSERT( sortedargs.size() == (size_t)1, "Extractors::filterArgumentsByParticleTypes() - reordering did not generate a single list."); // then we split up the tuples of arguments and place each into single list FunctionModel::list_of_arguments_t returnargs; FunctionModel::arguments_t::const_iterator argiter = sortedargs.begin()->begin(); const size_t num_types = _types.size(); const size_t args_per_tuple = num_types * (num_types-1) / 2; while (argiter != sortedargs.begin()->end()) { FunctionModel::arguments_t currenttuple(args_per_tuple); const FunctionModel::arguments_t::const_iterator startiter = argiter; std::advance(argiter, args_per_tuple); #ifndef NDEBUG FunctionModel::arguments_t::const_iterator endoutiter = #endif std::copy(startiter, argiter, currenttuple.begin()); ASSERT( endoutiter == currenttuple.end(), "Extractors::filterArgumentsByParticleTypes() - currenttuple not initialized to right size."); returnargs.push_back(currenttuple); } LOG(2, "DEBUG: We have generated " << returnargs.size() << " tuples of distances."); return returnargs; } FunctionModel::arguments_t Extractors::combineArguments( const FunctionModel::arguments_t &firstargs, const FunctionModel::arguments_t &secondargs) { FunctionModel::arguments_t args = concatenateArguments(firstargs, secondargs); std::sort(args.begin(), args.end(), boost::bind(&argument_t::operator<, _1, _2)); FunctionModel::arguments_t::iterator iter = std::unique(args.begin(), args.end(), boost::bind(&argument_t::operator==, _1, _2)); args.erase(iter, args.end()); return args; } FunctionModel::arguments_t Extractors::concatenateArguments( const FunctionModel::arguments_t &firstargs, const FunctionModel::arguments_t &secondargs) { FunctionModel::arguments_t args(firstargs); args.insert(args.end(), secondargs.begin(), secondargs.end()); return args; } FunctionModel::list_of_arguments_t Extractors::concatenateListOfArguments( const FunctionModel::list_of_arguments_t &firstlistargs, const FunctionModel::list_of_arguments_t &secondlistargs) { FunctionModel::list_of_arguments_t listargs(firstlistargs); listargs.insert(listargs.end(), secondlistargs.begin(), secondlistargs.end()); return listargs; }