Changeset 5de3c9 for src/molecules.cpp

Timestamp:

May 8, 2008, 5:31:04 PM (17 years ago)

Author:

Frederik Heber <heber@…>

Branches:

Action_Thermostats, Add_AtomRandomPerturbation, Add_FitFragmentPartialChargesAction, Add_RotateAroundBondAction, Add_SelectAtomByNameAction, Added_ParseSaveFragmentResults, AddingActions_SaveParseParticleParameters, Adding_Graph_to_ChangeBondActions, Adding_MD_integration_tests, Adding_ParticleName_to_Atom, Adding_StructOpt_integration_tests, AtomFragments, Automaking_mpqc_open, AutomationFragmentation_failures, Candidate_v1.5.4, Candidate_v1.6.0, Candidate_v1.6.1, ChangeBugEmailaddress, ChangingTestPorts, ChemicalSpaceEvaluator, CombiningParticlePotentialParsing, 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_BoundInBox_CenterInBox_MoleculeActions, Fix_ChargeSampling_PBC, Fix_ChronosMutex, Fix_FitPartialCharges, Fix_FitPotential_needs_atomicnumbers, Fix_ForceAnnealing, Fix_IndependentFragmentGrids, Fix_ParseParticles, Fix_ParseParticles_split_forward_backward_Actions, Fix_PopActions, Fix_QtFragmentList_sorted_selection, Fix_Restrictedkeyset_FragmentMolecule, Fix_StatusMsg, Fix_StepWorldTime_single_argument, Fix_Verbose_Codepatterns, Fix_fitting_potentials, Fixes, ForceAnnealing_goodresults, ForceAnnealing_oldresults, ForceAnnealing_tocheck, ForceAnnealing_with_BondGraph, ForceAnnealing_with_BondGraph_continued, ForceAnnealing_with_BondGraph_continued_betteresults, ForceAnnealing_with_BondGraph_contraction-expansion, FragmentAction_writes_AtomFragments, FragmentMolecule_checks_bonddegrees, GeometryObjects, Gui_Fixes, Gui_displays_atomic_force_velocity, ImplicitCharges, IndependentFragmentGrids, IndependentFragmentGrids_IndividualZeroInstances, IndependentFragmentGrids_IntegrationTest, IndependentFragmentGrids_Sole_NN_Calculation, JobMarket_RobustOnKillsSegFaults, JobMarket_StableWorkerPool, JobMarket_unresolvable_hostname_fix, MoreRobust_FragmentAutomation, ODR_violation_mpqc_open, PartialCharges_OrthogonalSummation, PdbParser_setsAtomName, PythonUI_with_named_parameters, QtGui_reactivate_TimeChanged_changes, Recreated_GuiChecks, Rewrite_FitPartialCharges, RotateToPrincipalAxisSystem_UndoRedo, SaturateAtoms_findBestMatching, SaturateAtoms_singleDegree, StoppableMakroAction, Subpackage_CodePatterns, Subpackage_JobMarket, Subpackage_LinearAlgebra, Subpackage_levmar, Subpackage_mpqc_open, Subpackage_vmg, Switchable_LogView, ThirdParty_MPQC_rebuilt_buildsystem, TrajectoryDependenant_MaxOrder, TremoloParser_IncreasedPrecision, TremoloParser_MultipleTimesteps, TremoloParser_setsAtomName, Ubuntu_1604_changes, stable

Children:

6f6a8e

Parents:

a529de

Message:

molecule::CheckOrderAtSite() now interprets negative Orders as adaptive increase by parsing EnergyPerFragment.da
t

positive Order is global increase till all sites have at least this order, Order 0 means single global increase step, negative Order is the exponent in 10^{{-Order} to give the threshold value: ENERGYPERFRAGMENT is scanned in
to a map (sorted by values) and all fragments still above the threshold are taken into AtomMask for increase. Jo
iner has to write this file with (Root id of Fragment, energy contribution) pairs.}

File:

• src/molecules.cpp (modified) (5 diffs)

src/molecules.cpp

@@ -2017 +2017 @@
     *out << Verbose(1) << "done." << endl;
   } else {
-    *out << Verbose(1) << "failed to open file" << line << "." << endl;
+    *out << Verbose(1) << "failed to open " << line << "." << endl;
     status = false;
   }
@@ -2137 +2137 @@
 /** Checks whether the OrderAtSite is still bewloe \a Order at some site.
  * \param *out output stream for debugging
- * \param &FragmentationToDo return boolean
- * \param Order desired Order
+ * \param *AtomMask defines true/false per global Atom::nr to mask in/out each nuclear site, used to activate given number of site to increment order adaptively
+ * \param *GlobalKeySetList list of keysets with global ids (valid in "this" molecule) needed for adaptive increase
+ * \param Order desired Order if positive, desired exponent in threshold criteria if negative (0 is single-step)
+ * \param *path path to ENERGYPERFRAGMENT file (may be NULL if Order is non-negative)
  * \return true - needs further fragmentation, false - does not need fragmentation
  */
-bool molecule::CheckOrderAtSite(ofstream *out, int Order)
-{
-  atom *Walker = NULL;
+bool molecule::CheckOrderAtSite(ofstream *out, bool *AtomMask, Graph *GlobalKeySetList, int Order, char *path)
+{
+  atom *Walker = start;
   bool status = false;
-  
-  Walker = start;
-  while (Walker->next != end) {
-    Walker = Walker->next;
-#ifdef ADDHYDROGEN
-    if (Walker->type->Z != 1) // skip hydrogen
-#endif
-      if (Walker->AdaptiveOrder < Order)
-        status = true;
-  }
-  if (!status)
-    *out << Verbose(1) << "Order at every site is already equal or above desired order " << Order << "." << endl;
+  ifstream InputFile;
+  stringstream line;
+  
+  if (Order < 0) { // adaptive increase of BondOrder per site
+    for(int i=0;i<AtomCount;i++)
+      AtomMask[i] = false;
+    // parse the EnergyPerFragment file
+    line.str(path);
+    line << "/" << FRAGMENTPREFIX << ENERGYPERFRAGMENT;
+    InputFile.open(line.str().c_str(), ios::in);
+    if (InputFile != NULL) {
+      char *buffer = (char *) Malloc(sizeof(char)*MAXSTRINGSIZE, "molecule::CheckOrderAtSite: *buffer");
+      int lines = 0;
+      // count the number of lines, i.e. the number of fragments
+      while(!InputFile.eof()) {
+        InputFile.getline(buffer, MAXSTRINGSIZE);
+        lines++;
+      }
+      InputFile.clear();
+      InputFile.seekg(ios::beg);
+      map<double, int> FragmentEnergies;
+      int No;
+      int Value;
+      // each line represents a fragment root (Atom::nr) id and its energy contribution
+      while(!InputFile.eof()) {
+        InputFile.getline(buffer, MAXSTRINGSIZE);
+        InputFile >> No;
+        InputFile >> Value;
+        FragmentEnergies.insert( make_pair(Value, No) );
+      }
+      for(map<double, int>::iterator runner = FragmentEnergies.lower_bound(pow(1.,Order)); runner != FragmentEnergies.begin(); runner++) {
+        // as the smallest number in each set has always been the root (we use global id to keep the doubles away), seek smallest and insert into AtomMask
+        AtomMask[(*runner).second] = true;
+      }
+      // close and done
+      InputFile.close();
+      InputFile.clear();
+      Free((void **)&buffer, "molecule::CheckOrderAtSite: *buffer");
+    } else {
+      cerr << "Unable to parse " << FRAGMENTPREFIX << ENERGYPERFRAGMENT << " file." << endl;
+      return false;
+    }
+    // pick a given number of highest values and set AtomMask
+  } else if (Order > 0) { // global increase of Bond Order
+    while (Walker->next != end) {
+      Walker = Walker->next;
+  #ifdef ADDHYDROGEN
+      if (Walker->type->Z != 1) // skip hydrogen
+  #endif
+      {
+        AtomMask[Walker->nr] = true;  // include all (non-hydrogen) atoms
+        if (Walker->AdaptiveOrder < Order)
+          status = true;
+      }
+    }
+    if (!status)
+      *out << Verbose(1) << "Order at every site is already equal or above desired order " << Order << "." << endl;
+  } else {  // single stepping, just check
+    if (AtomMask[AtomCount] != true)  // if true, we would've done a step already
+      for(int i=0;i<=AtomCount;i++)
+        AtomMask[i] = true;
+    else
+      status = false;
+  }
  
   return status;
@@ -2226 +2280 @@
   atom **ListOfAtoms = NULL;
   atom ***ListOfLocalAtoms = NULL;
+  bool *AtomMask = NULL; 
   
   *out << endl;
@@ -2266 +2321 @@
     // ===== 6a. assign each keyset to its respective subgraph ===== 
     Subgraphs->next->AssignKeySetsToFragment(out, this, ParsedFragmentList, ListOfLocalAtoms, FragmentList, (FragmentCounter = 0), false);
-    delete(ParsedFragmentList);
 
     KeyStack *RootStack = new KeyStack[Subgraphs->next->Count()];
-    while (CheckOrderAtSite(out, Order)) {
+    AtomMask = new bool[AtomCount+1];
+    for (int i=0;i<AtomCount;i++) // need to set to false to recognize in single-stepping (one global step set them all to true)
+      AtomMask[i] = false;
+    while (CheckOrderAtSite(out, AtomMask, ParsedFragmentList, Order)) {
+      AtomMask[AtomCount] = true;   // last plus one entry is used as marker that we have been through this loop once already in CheckOrderAtSite()
       // ===== 6b. fill RootStack for each subgraph (second adaptivity check) ===== 
-      Subgraphs->next->FillRootStackForSubgraphs(out, RootStack, Order, (FragmentCounter = 0));
+      Subgraphs->next->FillRootStackForSubgraphs(out, RootStack, AtomMask, Order, (FragmentCounter = 0));
 
       // ===== 7. fill the bond fragment list =====
@@ -2293 +2351 @@
     }
     delete[](RootStack);
-  }
+    delete[](AtomMask);
+  }
+  delete(ParsedFragmentList);
   
   // free the index lookup list