Changeset cb80d4 for src/Dynamics/ForceAnnealing.hpp
- Timestamp:
- Jul 20, 2017, 9:38:38 AM (8 years ago)
- Branches:
- ForceAnnealing_with_BondGraph_continued
- Children:
- d3d964
- Parents:
- 461a7f
- git-author:
- Frederik Heber <frederik.heber@…> (06/29/17 14:53:38)
- git-committer:
- Frederik Heber <frederik.heber@…> (07/20/17 09:38:38)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
src/Dynamics/ForceAnnealing.hpp
r461a7f rcb80d4 278 278 279 279 // for each atom we need to gather weights and then project the gradient 280 typedef std::deque<double> weights_t; 281 typedef std::map<atomId_t, weights_t > weights_per_atom_t; 280 typedef std::map<atomId_t, BondVectors::weights_t > weights_per_atom_t; 282 281 std::vector<weights_per_atom_t> weights_per_atom(2); 283 282 for (size_t timestep = 0; timestep <= 1; ++timestep) { … … 297 296 298 297 // gather subset of BondVectors for the current atom 299 std::vector<Vector> BondVectors = bv.getAtomsBondVectorsAtStep(walker, CurrentStep); 298 const std::vector<Vector> BondVectors = 299 bv.getAtomsBondVectorsAtStep(walker, CurrentStep); 300 300 301 301 // go through all its bonds and calculate what magnitude is represented 302 302 // by the others i.e. sum of scalar products against other bonds 303 std::pair<weights_per_atom_t::iterator, bool> inserter =303 const std::pair<weights_per_atom_t::iterator, bool> inserter = 304 304 weights_per_atom[timestep].insert( 305 std::make_pair(walker.getId(), weights_t()) ); 305 std::make_pair(walker.getId(), 306 bv.getWeightsForAtomAtStep(walker, CurrentStep)) ); 306 307 ASSERT( inserter.second, 307 308 "ForceAnnealing::operator() - weight map for atom "+toString(walker) 308 309 +" and time step "+toString(timestep)+" already filled?"); 309 weights_t &weights = inserter.first->second; 310 for (std::vector<Vector>::const_iterator iter = BondVectors.begin(); 311 iter != BondVectors.end(); ++iter) { 312 std::vector<double> scps; 313 scps.reserve(BondVectors.size()); 314 std::transform( 315 BondVectors.begin(), BondVectors.end(), 316 std::back_inserter(scps), 317 boost::bind(static_cast< double (*)(double) >(&fabs), 318 boost::bind(&Vector::ScalarProduct, boost::cref(*iter), _1)) 319 ); 320 const double scp_sum = std::accumulate(scps.begin(), scps.end(), 0.); 321 ASSERT( (scp_sum-1.) > -MYEPSILON, 322 "ForceAnnealing() - sum of weights must be equal or larger one but is " 323 +toString(scp_sum)); 324 weights.push_back( 1./scp_sum ); 325 } 326 LOG(4, "DEBUG: Weights for atom #" << walker.getId() << ": " << weights); 327 328 // for testing we check whether all weighted scalar products now come out as 1. 329 #ifndef NDEBUG 330 for (std::vector<Vector>::const_iterator iter = BondVectors.begin(); 331 iter != BondVectors.end(); ++iter) { 332 std::vector<double> scps; 333 scps.reserve(BondVectors.size()); 334 std::transform( 335 BondVectors.begin(), BondVectors.end(), 336 weights.begin(), 337 std::back_inserter(scps), 338 boost::bind(static_cast< double (*)(double) >(&fabs), 339 boost::bind(std::multiplies<double>(), 340 boost::bind(&Vector::ScalarProduct, boost::cref(*iter), _1), 341 _2)) 342 ); 343 const double scp_sum = std::accumulate(scps.begin(), scps.end(), 0.); 344 ASSERT( fabs(scp_sum - 1.) < MYEPSILON, 345 "ForceAnnealing::operator() - for BondVector "+toString(*iter) 346 +" we have weighted scalar product of "+toString(scp_sum)+" != 1."); 347 } 348 #endif 310 BondVectors::weights_t &weights = inserter.first->second; 311 ASSERT( weights.size() == ListOfBonds.size(), 312 "ForceAnnealing::operator() - number of weights " 313 +toString(weights.size())+" does not match number of bonds " 314 +toString(ListOfBonds.size())+", error in calculation?"); 349 315 350 316 // projected gradient over all bonds and place in one of projected_forces 351 317 // using the obtained weights 352 318 { 353 weights_t::const_iterator weightiter = weights.begin();319 BondVectors::weights_t::const_iterator weightiter = weights.begin(); 354 320 std::vector<Vector>::const_iterator vectoriter = BondVectors.begin(); 355 321 Vector forcesum_check;
Note:
See TracChangeset
for help on using the changeset viewer.