Ignore:
Timestamp:
Jul 23, 2009, 1:45:24 PM (15 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:
51c910
Parents:
ce5ac3 (diff), 437922 (diff)
Note: this is a merge changeset, the changes displayed below correspond to the merge itself.
Use the (diff) links above to see all the changes relative to each parent.
git-author:
Frederik Heber <heber@…> (07/23/09 12:34:47)
git-committer:
Frederik Heber <heber@…> (07/23/09 13:45:24)
Message:

Merge branch 'MultipleMolecules'

Conflicts:

molecuilder/src/analyzer.cpp
molecuilder/src/atom.cpp
molecuilder/src/boundary.cpp
molecuilder/src/boundary.hpp
molecuilder/src/builder.cpp
molecuilder/src/config.cpp
molecuilder/src/datacreator.hpp
molecuilder/src/helpers.cpp
molecuilder/src/joiner.cpp
molecuilder/src/moleculelist.cpp
molecuilder/src/molecules.cpp
molecuilder/src/molecules.hpp
molecuilder/src/parser.cpp
molecuilder/src/parser.hpp
molecuilder/src/vector.cpp
molecuilder/src/verbose.cpp

merges:

compilation fixes:

File:
1 edited

Legend:

Unmodified
Added
Removed
  • src/moleculelist.cpp

    rce5ac3 rd067d45  
    11/** \file MoleculeListClass.cpp
    2  * 
     2 *
    33 * Function implementations for the class MoleculeListClass.
    4  * 
     4 *
    55 */
    66
     
    1313MoleculeListClass::MoleculeListClass()
    1414{
    15 };
    16 
    17 /** constructor for MoleculeListClass.
    18  * \param NumMolecules number of molecules to allocate for
    19  * \param NumAtoms number of atoms to allocate for
    20  */
    21 MoleculeListClass::MoleculeListClass(int NumMolecules, int NumAtoms)
    22 {
    23   ListOfMolecules = (molecule **) Malloc(sizeof(molecule *)*NumMolecules, "MoleculeListClass:MoleculeListClass: **ListOfMolecules");
    24   for (int i=NumMolecules;i--;)
    25     ListOfMolecules[i] = NULL;
    26   NumberOfMolecules = NumMolecules;
    27   NumberOfTopAtoms = NumAtoms;
    28 };
    29 
     15  // empty lists
     16  ListOfMolecules.clear();
     17  MaxIndex = 1;
     18};
    3019
    3120/** Destructor for MoleculeListClass.
     
    3423{
    3524  cout << Verbose(3) << this << ": Freeing ListOfMolcules." << endl;
    36   for (int i=NumberOfMolecules;i--;) {
    37     if (ListOfMolecules[i] != NULL) { // if NULL don't free
    38       cout << Verbose(4) << "ListOfMolecules: Freeing " << ListOfMolecules[i] << "." << endl;
    39       delete(ListOfMolecules[i]);
    40       ListOfMolecules[i] = NULL;
    41     }
     25  for (MoleculeList::iterator ListRunner = ListOfMolecules.begin(); ListRunner != ListOfMolecules.end(); ListRunner++) {
     26    cout << Verbose(4) << "ListOfMolecules: Freeing " << *ListRunner << "." << endl;
     27    delete (*ListRunner);
    4228  }
    4329  cout << Verbose(4) << "Freeing ListOfMolecules." << endl;
    44   Free((void **)&ListOfMolecules, "MoleculeListClass:MoleculeListClass: **ListOfMolecules");
     30  ListOfMolecules.clear(); // empty list
     31};
     32
     33/** Insert a new molecule into the list and set its number.
     34 * \param *mol molecule to add to list.
     35 * \return true - add successful
     36 */
     37void MoleculeListClass::insert(molecule *mol)
     38{
     39  mol->IndexNr = MaxIndex++;
     40  ListOfMolecules.push_back(mol);
    4541};
    4642
     
    5753  atom *aWalker = NULL;
    5854  atom *bWalker = NULL;
    59  
     55
    6056  // sort each atom list and put the numbers into a list, then go through
    6157  //cout << "Comparing fragment no. " << *(molecule **)a << " to " << *(molecule **)b << "." << endl;
    62   if ( (**(molecule **)a).AtomCount < (**(molecule **)b).AtomCount ) {
     58  if ((**(molecule **) a).AtomCount < (**(molecule **) b).AtomCount) {
    6359    return -1;
    64   } else { if ((**(molecule **)a).AtomCount > (**(molecule **)b).AtomCount)
    65     return +1;
     60  } else {
     61    if ((**(molecule **) a).AtomCount > (**(molecule **) b).AtomCount)
     62      return +1;
    6663    else {
    67       Count = (**(molecule **)a).AtomCount;
     64      Count = (**(molecule **) a).AtomCount;
    6865      aList = new int[Count];
    6966      bList = new int[Count];
    70  
     67
    7168      // fill the lists
    72       aWalker = (**(molecule **)a).start;
    73       bWalker = (**(molecule **)b).start;
     69      aWalker = (**(molecule **) a).start;
     70      bWalker = (**(molecule **) b).start;
    7471      Counter = 0;
    7572      aCounter = 0;
    7673      bCounter = 0;
    77       while ((aWalker->next != (**(molecule **)a).end) && (bWalker->next != (**(molecule **)b).end)) {
     74      while ((aWalker->next != (**(molecule **) a).end) && (bWalker->next != (**(molecule **) b).end)) {
    7875        aWalker = aWalker->next;
    7976        bWalker = bWalker->next;
     
    9087      // check if AtomCount was for real
    9188      flag = 0;
    92       if ((aWalker->next == (**(molecule **)a).end) && (bWalker->next != (**(molecule **)b).end)) {
     89      if ((aWalker->next == (**(molecule **) a).end) && (bWalker->next != (**(molecule **) b).end)) {
    9390        flag = -1;
    9491      } else {
    95         if ((aWalker->next != (**(molecule **)a).end) && (bWalker->next == (**(molecule **)b).end))
     92        if ((aWalker->next != (**(molecule **) a).end) && (bWalker->next == (**(molecule **) b).end))
    9693          flag = 1;
    9794      }
    9895      if (flag == 0) {
    9996        // sort the lists
    100         gsl_heapsort(aList,Count, sizeof(int), CompareDoubles);
    101         gsl_heapsort(bList,Count, sizeof(int), CompareDoubles);
    102         // compare the lists 
    103        
     97        gsl_heapsort(aList, Count, sizeof(int), CompareDoubles);
     98        gsl_heapsort(bList, Count, sizeof(int), CompareDoubles);
     99        // compare the lists
     100
    104101        flag = 0;
    105         for(int i=0;i<Count;i++) {
     102        for (int i = 0; i < Count; i++) {
    106103          if (aList[i] < bList[i]) {
    107104            flag = -1;
     
    114111        }
    115112      }
    116       delete[](aList);
    117       delete[](bList);
     113      delete[] (aList);
     114      delete[] (bList);
    118115      return flag;
    119116    }
    120117  }
    121   return  -1;
     118  return -1;
     119};
     120
     121/** Output of a list of all molecules.
     122 * \param *out output stream
     123 */
     124void MoleculeListClass::Enumerate(ofstream *out)
     125{
     126  element* Elemental = NULL;
     127  atom *Walker = NULL;
     128  int Counts[MAX_ELEMENTS];
     129  double size=0;
     130  Vector Origin;
     131
     132  // header
     133  *out << "Index\tName\t\tAtoms\tFormula\tCenter\tSize" << endl;
     134  cout << Verbose(0) << "-----------------------------------------------" << endl;
     135  if (ListOfMolecules.size() == 0)
     136    *out << "\tNone" << endl;
     137  else {
     138    Origin.Zero();
     139    for (MoleculeList::iterator ListRunner = ListOfMolecules.begin(); ListRunner != ListOfMolecules.end(); ListRunner++) {
     140      // reset element counts
     141      for (int j = 0; j<MAX_ELEMENTS;j++)
     142        Counts[j] = 0;
     143      // count atoms per element and determine size of bounding sphere
     144      size=0.;
     145      Walker = (*ListRunner)->start;
     146      while (Walker->next != (*ListRunner)->end) {
     147        Walker = Walker->next;
     148        Counts[Walker->type->Z]++;
     149        if (Walker->x.DistanceSquared(&Origin) > size)
     150          size = Walker->x.DistanceSquared(&Origin);
     151      }
     152      // output Index, Name, number of atoms, chemical formula
     153      *out << ((*ListRunner)->ActiveFlag ? "*" : " ") << (*ListRunner)->IndexNr << "\t" << (*ListRunner)->name << "\t\t" << (*ListRunner)->AtomCount << "\t";
     154      Elemental = (*ListRunner)->elemente->end;
     155      while(Elemental->previous != (*ListRunner)->elemente->start) {
     156        Elemental = Elemental->previous;
     157        if (Counts[Elemental->Z] != 0)
     158          *out << Elemental->symbol << Counts[Elemental->Z];
     159      }
     160      // Center and size
     161      *out << "\t" << (*ListRunner)->Center << "\t" << sqrt(size) << endl;
     162    }
     163  }
     164};
     165
     166/** Returns the molecule with the given index \a index.
     167 * \param index index of the desired molecule
     168 * \return pointer to molecule structure, NULL if not found
     169 */
     170molecule * MoleculeListClass::ReturnIndex(int index)
     171{
     172  for(MoleculeList::iterator ListRunner = ListOfMolecules.begin(); ListRunner != ListOfMolecules.end(); ListRunner++)
     173    if ((*ListRunner)->IndexNr == index)
     174      return (*ListRunner);
     175  return NULL;
     176};
     177
     178/** Simple merge of two molecules into one.
     179 * \param *mol destination molecule
     180 * \param *srcmol source molecule
     181 * \return true - merge successful, false - merge failed (probably due to non-existant indices
     182 */
     183bool MoleculeListClass::SimpleMerge(molecule *mol, molecule *srcmol)
     184{
     185  if (srcmol == NULL)
     186    return false;
     187
     188  // put all molecules of src into mol
     189  atom *Walker = srcmol->start;
     190  atom *NextAtom = Walker->next;
     191  while (NextAtom != srcmol->end) {
     192    Walker = NextAtom;
     193    NextAtom = Walker->next;
     194    srcmol->UnlinkAtom(Walker);
     195    mol->AddAtom(Walker);
     196  }
     197
     198  // remove src
     199  ListOfMolecules.remove(srcmol);
     200  delete(srcmol);
     201  return true;
     202};
     203
     204/** Simple add of one molecules into another.
     205 * \param *mol destination molecule
     206 * \param *srcmol source molecule
     207 * \return true - merge successful, false - merge failed (probably due to non-existant indices
     208 */
     209bool MoleculeListClass::SimpleAdd(molecule *mol, molecule *srcmol)
     210{
     211  if (srcmol == NULL)
     212    return false;
     213
     214  // put all molecules of src into mol
     215  atom *Walker = srcmol->start;
     216  atom *NextAtom = Walker->next;
     217  while (NextAtom != srcmol->end) {
     218    Walker = NextAtom;
     219    NextAtom = Walker->next;
     220    Walker = mol->AddCopyAtom(Walker);
     221    Walker->father = Walker;
     222  }
     223
     224  return true;
     225};
     226
     227/** Simple merge of a given set of molecules into one.
     228 * \param *mol destination molecule
     229 * \param *src index of set of source molecule
     230 * \param N number of source molecules
     231 * \return true - merge successful, false - some merges failed (probably due to non-existant indices)
     232 */
     233bool MoleculeListClass::SimpleMultiMerge(molecule *mol, int *src, int N)
     234{
     235  bool status = true;
     236  // check presence of all source molecules
     237  for (int i=0;i<N;i++) {
     238    molecule *srcmol = ReturnIndex(src[i]);
     239    status = status && SimpleMerge(mol, srcmol);
     240  }
     241  return status;
     242};
     243
     244/** Simple add of a given set of molecules into one.
     245 * \param *mol destination molecule
     246 * \param *src index of set of source molecule
     247 * \param N number of source molecules
     248 * \return true - merge successful, false - some merges failed (probably due to non-existant indices)
     249 */
     250bool MoleculeListClass::SimpleMultiAdd(molecule *mol, int *src, int N)
     251{
     252  bool status = true;
     253  // check presence of all source molecules
     254  for (int i=0;i<N;i++) {
     255    molecule *srcmol = ReturnIndex(src[i]);
     256    status = status && SimpleAdd(mol, srcmol);
     257  }
     258  return status;
     259};
     260
     261/** Scatter merge of a given set of molecules into one.
     262 * Scatter merge distributes the molecules in such a manner that they don't overlap.
     263 * \param *mol destination molecule
     264 * \param *src index of set of source molecule
     265 * \param N number of source molecules
     266 * \return true - merge successful, false - merge failed (probably due to non-existant indices
     267 * \TODO find scatter center for each src molecule
     268 */
     269bool MoleculeListClass::ScatterMerge(molecule *mol, int *src, int N)
     270{
     271  // check presence of all source molecules
     272  for (int i=0;i<N;i++) {
     273    // get pointer to src molecule
     274    molecule *srcmol = ReturnIndex(src[i]);
     275    if (srcmol == NULL)
     276      return false;
     277  }
     278  // adapt each Center
     279  for (int i=0;i<N;i++) {
     280    // get pointer to src molecule
     281    molecule *srcmol = ReturnIndex(src[i]);
     282    //srcmol->Center.Zero();
     283    srcmol->Translate(&srcmol->Center);
     284  }
     285  // perform a simple multi merge
     286  SimpleMultiMerge(mol, src, N);
     287  return true;
     288};
     289
     290/** Embedding merge of a given set of molecules into one.
     291 * Embedding merge inserts one molecule into the other.
     292 * \param *mol destination molecule
     293 * \param *srcmol source molecule
     294 * \return true - merge successful, false - merge failed (probably due to non-existant indices
     295 * \TODO find embedding center
     296 */
     297bool MoleculeListClass::EmbedMerge(molecule *mol, molecule *srcmol)
     298{
     299  if (srcmol == NULL)
     300    return false;
     301
     302  // calculate center for merge
     303  srcmol->Center.CopyVector(mol->FindEmbeddingHole((ofstream *)&cout, srcmol));
     304  srcmol->Center.Zero();
     305
     306  // perform simple merge
     307  SimpleMerge(mol, srcmol);
     308  return true;
    122309};
    123310
     
    127314void MoleculeListClass::Output(ofstream *out)
    128315{
    129   *out<< Verbose(1) << "MoleculeList: ";
    130   for (int i=0;i<NumberOfMolecules;i++)
    131     *out << ListOfMolecules[i] << "\t";
     316  *out << Verbose(1) << "MoleculeList: ";
     317  for (MoleculeList::iterator ListRunner = ListOfMolecules.begin(); ListRunner != ListOfMolecules.end(); ListRunner++)
     318    *out << *ListRunner << "\t";
    132319  *out << endl;
    133320};
     
    145332  atom *Runner = NULL;
    146333  double ***FitConstant = NULL, **correction = NULL;
    147   int a,b;
     334  int a, b;
    148335  ofstream output;
    149336  ifstream input;
     
    165352  input.open(line.c_str());
    166353  if (input == NULL) {
    167     cerr << endl << "Unable to open " << line << ", is the directory correct?" << endl;
     354    cerr << endl << "Unable to open " << line << ", is the directory correct?"
     355        << endl;
    168356    return false;
    169357  }
    170   a=0;
    171   b=-1; // we overcount by one
     358  a = 0;
     359  b = -1; // we overcount by one
    172360  while (!input.eof()) {
    173361    input.getline(ParsedLine, 1023);
    174362    zeile.str(ParsedLine);
    175     int i=0;
     363    int i = 0;
    176364    while (!zeile.eof()) {
    177365      zeile >> distance;
    178366      i++;
    179367    }
    180     if (i > a) 
    181       a = i; 
     368    if (i > a)
     369      a = i;
    182370    b++;
    183371  }
    184372  cout << "I recognized " << a << " columns and " << b << " rows, ";
    185373  input.close();
    186  
     374
    187375  // 0b. allocate memory for constants
    188   FitConstant = (double ***) Malloc(sizeof(double **)*3, "MoleculeListClass::AddHydrogenCorrection: ***FitConstant");
    189   for (int k=0;k<3;k++) {
    190     FitConstant[k] = (double **) Malloc(sizeof(double *)*a, "MoleculeListClass::AddHydrogenCorrection: **FitConstant[]");
    191     for (int i=a;i--;) {
    192       FitConstant[k][i] = (double *) Malloc(sizeof(double)*b, "MoleculeListClass::AddHydrogenCorrection: *FitConstant[][]");
     376  FitConstant = (double ***) Malloc(sizeof(double **) * 3, "MoleculeListClass::AddHydrogenCorrection: ***FitConstant");
     377  for (int k = 0; k < 3; k++) {
     378    FitConstant[k] = (double **) Malloc(sizeof(double *) * a, "MoleculeListClass::AddHydrogenCorrection: **FitConstant[]");
     379    for (int i = a; i--;) {
     380      FitConstant[k][i] = (double *) Malloc(sizeof(double) * b, "MoleculeListClass::AddHydrogenCorrection: *FitConstant[][]");
    193381    }
    194382  }
    195383  // 0c. parse in constants
    196   for (int i=0;i<3;i++) {
     384  for (int i = 0; i < 3; i++) {
    197385    line = path;
    198386    line.append("/");
    199387    line += FRAGMENTPREFIX;
    200     sprintf(ParsedLine, "%d", i+1);
     388    sprintf(ParsedLine, "%d", i + 1);
    201389    line += ParsedLine;
    202390    line += FITCONSTANTSUFFIX;
     
    206394      return false;
    207395    }
    208     int k = 0,l;
     396    int k = 0, l;
    209397    while ((!input.eof()) && (k < b)) {
    210398      input.getline(ParsedLine, 1023);
     
    223411    input.close();
    224412  }
    225   for(int k=0;k<3;k++) {
     413  for (int k = 0; k < 3; k++) {
    226414    cout << "Constants " << k << ":" << endl;
    227     for (int j=0;j<b;j++) {
    228       for (int i=0;i<a;i++) {
     415    for (int j = 0; j < b; j++) {
     416      for (int i = 0; i < a; i++) {
    229417        cout << FitConstant[k][i][j] << "\t";
    230418      }
     
    233421    cout << endl;
    234422  }
    235  
     423
    236424  // 0d. allocate final correction matrix
    237   correction = (double **) Malloc(sizeof(double *)*a, "MoleculeListClass::AddHydrogenCorrection: **correction");
    238   for (int i=a;i--;)
    239     correction[i] = (double *) Malloc(sizeof(double)*b, "MoleculeListClass::AddHydrogenCorrection: *correction[]");
    240        
     425  correction = (double **) Malloc(sizeof(double *) * a, "MoleculeListClass::AddHydrogenCorrection: **correction");
     426  for (int i = a; i--;)
     427    correction[i] = (double *) Malloc(sizeof(double) * b, "MoleculeListClass::AddHydrogenCorrection: *correction[]");
     428
    241429  // 1a. go through every molecule in the list
    242   for(int i=NumberOfMolecules;i--;) {
     430  for (MoleculeList::iterator ListRunner = ListOfMolecules.begin(); ListRunner != ListOfMolecules.end(); ListRunner++) {
    243431    // 1b. zero final correction matrix
    244     for (int k=a;k--;)
    245       for (int j=b;j--;)
     432    for (int k = a; k--;)
     433      for (int j = b; j--;)
    246434        correction[k][j] = 0.;
    247435    // 2. take every hydrogen that is a saturated one
    248     Walker = ListOfMolecules[i]->start;
    249     while (Walker->next != ListOfMolecules[i]->end) {
     436    Walker = (*ListRunner)->start;
     437    while (Walker->next != (*ListRunner)->end) {
    250438      Walker = Walker->next;
    251       //cout << Verbose(1) << "Walker: " << *Walker << " with first bond " << *ListOfMolecules[i]->ListOfBondsPerAtom[Walker->nr][0] << "." << endl;
    252       if ((Walker->type->Z == 1) && ((Walker->father == NULL) || (Walker->father->type->Z != 1))) { // if it's a hydrogen
    253         Runner = ListOfMolecules[i]->start;
    254         while (Runner->next != ListOfMolecules[i]->end) {
     439      //cout << Verbose(1) << "Walker: " << *Walker << " with first bond " << *(*Runner)->ListOfBondsPerAtom[Walker->nr][0] << "." << endl;
     440      if ((Walker->type->Z == 1) && ((Walker->father == NULL)
     441          || (Walker->father->type->Z != 1))) { // if it's a hydrogen
     442        Runner = (*ListRunner)->start;
     443        while (Runner->next != (*ListRunner)->end) {
    255444          Runner = Runner->next;
    256           //cout << Verbose(2) << "Runner: " << *Runner << " with first bond " << *ListOfMolecules[i]->ListOfBondsPerAtom[Runner->nr][0] << "." << endl;
     445          //cout << Verbose(2) << "Runner: " << *Runner << " with first bond " << *(*Runner)->ListOfBondsPerAtom[Runner->nr][0] << "." << endl;
    257446          // 3. take every other hydrogen that is the not the first and not bound to same bonding partner
    258           if ((Runner->type->Z == 1) && (Runner->nr > Walker->nr) && (ListOfMolecules[i]->ListOfBondsPerAtom[Runner->nr][0]->GetOtherAtom(Runner) != ListOfMolecules[i]->ListOfBondsPerAtom[Walker->nr][0]->GetOtherAtom(Walker))) { // (hydrogens have only one bonding partner!)
     447          if ((Runner->type->Z == 1) && (Runner->nr > Walker->nr) && ((*ListRunner)->ListOfBondsPerAtom[Runner->nr][0]->GetOtherAtom(Runner) != (*ListRunner)->ListOfBondsPerAtom[Walker->nr][0]->GetOtherAtom(Walker))) { // (hydrogens have only one bonding partner!)
    259448            // 4. evaluate the morse potential for each matrix component and add up
    260             distance = sqrt(Runner->x.Distance(&Walker->x));
    261             //cout << "Fragment " << i << ": " << *Runner << "<= " << distance << "=>" << *Walker << ":" << endl;
    262             for(int k=0;k<a;k++) {
    263               for (int j=0;j<b;j++) {
    264                 switch(k) {
    265                 case 1:
    266                 case 7:
    267                 case 11:
    268                   tmp = pow(FitConstant[0][k][j] * ( 1. - exp(-FitConstant[1][k][j] * (distance - FitConstant[2][k][j]) ) ),2);
    269                   break;
    270                 default:
    271                   tmp = FitConstant[0][k][j] * pow( distance, FitConstant[1][k][j]) + FitConstant[2][k][j];
     449            distance = Runner->x.Distance(&Walker->x);
     450            //cout << "Fragment " << (*ListRunner)->name << ": " << *Runner << "<= " << distance << "=>" << *Walker << ":" << endl;
     451            for (int k = 0; k < a; k++) {
     452              for (int j = 0; j < b; j++) {
     453                switch (k) {
     454                  case 1:
     455                  case 7:
     456                  case 11:
     457                    tmp = pow(FitConstant[0][k][j] * (1. - exp(-FitConstant[1][k][j] * (distance - FitConstant[2][k][j]))), 2);
     458                    break;
     459                  default:
     460                    tmp = FitConstant[0][k][j] * pow(distance, FitConstant[1][k][j]) + FitConstant[2][k][j];
    272461                };
    273                 correction[k][j] -= tmp;    // ground state is actually lower (disturbed by additional interaction)
     462                correction[k][j] -= tmp; // ground state is actually lower (disturbed by additional interaction)
    274463                //cout << tmp << "\t";
    275464              }
     
    285474    line.append("/");
    286475    line += FRAGMENTPREFIX;
    287     FragmentNumber = FixedDigitNumber(NumberOfMolecules, i);
     476    FragmentNumber = FixedDigitNumber(ListOfMolecules.size(), (*ListRunner)->IndexNr);
    288477    line += FragmentNumber;
    289     delete(FragmentNumber);
     478    delete (FragmentNumber);
    290479    line += HCORRECTIONSUFFIX;
    291480    output.open(line.c_str());
    292481    output << "Time\t\tTotal\t\tKinetic\t\tNonLocal\tCorrelation\tExchange\tPseudo\t\tHartree\t\t-Gauss\t\tEwald\t\tIonKin\t\tETotal" << endl;
    293     for (int j=0;j<b;j++) {
    294       for(int i=0;i<a;i++)
     482    for (int j = 0; j < b; j++) {
     483      for (int i = 0; i < a; i++)
    295484        output << correction[i][j] << "\t";
    296485      output << endl;
     
    303492  output.open(line.c_str());
    304493  output << "Time\t\tTotal\t\tKinetic\t\tNonLocal\tCorrelation\tExchange\tPseudo\t\tHartree\t\t-Gauss\t\tEwald\t\tIonKin\t\tETotal" << endl;
    305   for (int j=0;j<b;j++) {
    306     for(int i=0;i<a;i++)
     494  for (int j = 0; j < b; j++) {
     495    for (int i = 0; i < a; i++)
    307496      output << 0 << "\t";
    308497    output << endl;
     
    310499  output.close();
    311500  // 6. free memory of parsed matrices
    312   FitConstant = (double ***) Malloc(sizeof(double **)*a, "MoleculeListClass::AddHydrogenCorrection: ***FitConstant");
    313   for (int k=0;k<3;k++) {
    314     FitConstant[k] = (double **) Malloc(sizeof(double *)*a, "MoleculeListClass::AddHydrogenCorrection: **FitConstant[]");
    315     for (int i=a;i--;) {
    316       FitConstant[k][i] = (double *) Malloc(sizeof(double)*b, "MoleculeListClass::AddHydrogenCorrection: *FitConstant[][]");
     501  FitConstant = (double ***) Malloc(sizeof(double **) * a, "MoleculeListClass::AddHydrogenCorrection: ***FitConstant");
     502  for (int k = 0; k < 3; k++) {
     503    FitConstant[k] = (double **) Malloc(sizeof(double *) * a, "MoleculeListClass::AddHydrogenCorrection: **FitConstant[]");
     504    for (int i = a; i--;) {
     505      FitConstant[k][i] = (double *) Malloc(sizeof(double) * b, "MoleculeListClass::AddHydrogenCorrection: *FitConstant[][]");
    317506    }
    318507  }
     
    327516 * \return true - file written successfully, false - writing failed
    328517 */
    329 bool MoleculeListClass::StoreForcesFile(ofstream *out, char *path, int *SortIndex)
     518bool MoleculeListClass::StoreForcesFile(ofstream *out, char *path,
     519    int *SortIndex)
    330520{
    331521  bool status = true;
     
    342532    //cout << Verbose(1) << "Final AtomicForcesList: ";
    343533    //output << prefix << "Forces" << endl;
    344     for(int j=0;j<NumberOfMolecules;j++) {
    345       //if (TEList[j] != 0) {
    346       runner = ListOfMolecules[j]->elemente->start;
    347       while (runner->next != ListOfMolecules[j]->elemente->end) { // go through every element
     534    for (MoleculeList::iterator ListRunner = ListOfMolecules.begin(); ListRunner != ListOfMolecules.end(); ListRunner++) {
     535      runner = (*ListRunner)->elemente->start;
     536      while (runner->next != (*ListRunner)->elemente->end) { // go through every element
    348537        runner = runner->next;
    349         if (ListOfMolecules[j]->ElementsInMolecule[runner->Z]) { // if this element got atoms
    350           Walker = ListOfMolecules[j]->start;
    351           while (Walker->next != ListOfMolecules[j]->end) { // go through every atom of this element
     538        if ((*ListRunner)->ElementsInMolecule[runner->Z]) { // if this element got atoms
     539          Walker = (*ListRunner)->start;
     540          while (Walker->next != (*ListRunner)->end) { // go through every atom of this element
    352541            Walker = Walker->next;
    353542            if (Walker->type->Z == runner->Z) {
    354543              if ((Walker->GetTrueFather() != NULL) && (Walker->GetTrueFather() != Walker)) {// if there is a rea
    355544                //cout << "Walker is " << *Walker << " with true father " << *( Walker->GetTrueFather()) << ", it
    356                 ForcesFile <<  SortIndex[Walker->GetTrueFather()->nr] << "\t";
    357                 } else  // otherwise a -1 to indicate an added saturation hydrogen
    358                   ForcesFile << "-1\t";
    359               }
     545                ForcesFile << SortIndex[Walker->GetTrueFather()->nr] << "\t";
     546              } else
     547                // otherwise a -1 to indicate an added saturation hydrogen
     548                ForcesFile << "-1\t";
    360549            }
    361550          }
     551        }
    362552      }
    363553      ForcesFile << endl;
     
    370560  }
    371561  ForcesFile.close();
    372  
     562
    373563  return status;
    374564};
     
    382572 * \return true - success (each file was written), false - something went wrong.
    383573 */
    384 bool MoleculeListClass::OutputConfigForListOfFragments(ofstream *out, const char *fragmentprefix, config *configuration, int *SortIndex, bool DoPeriodic, bool DoCentering)
     574bool MoleculeListClass::OutputConfigForListOfFragments(ofstream *out, config *configuration, int *SortIndex)
    385575{
    386576  ofstream outputFragment;
     
    395585  int FragmentCounter = 0;
    396586  ofstream output;
    397   string basis("3-21G");
    398  
     587
    399588  // store the fragments as config and as xyz
    400   for(int i=0;i<NumberOfMolecules;i++) {
     589  for (MoleculeList::iterator ListRunner = ListOfMolecules.begin(); ListRunner != ListOfMolecules.end(); ListRunner++) {
    401590    // save default path as it is changed for each fragment
    402591    path = configuration->GetDefaultPath();
     
    407596
    408597    // correct periodic
    409     if (DoPeriodic)
    410       ListOfMolecules[i]->ScanForPeriodicCorrection(out);
     598    (*ListRunner)->ScanForPeriodicCorrection(out);
    411599
    412600    // output xyz file
    413     FragmentNumber = FixedDigitNumber(NumberOfMolecules, FragmentCounter++);
    414     sprintf(FragmentName, "%s/%s%s.conf.xyz", configuration->configpath, fragmentprefix, FragmentNumber);
     601    FragmentNumber = FixedDigitNumber(ListOfMolecules.size(), FragmentCounter++);
     602    sprintf(FragmentName, "%s/%s%s.conf.xyz", configuration->configpath, FRAGMENTPREFIX, FragmentNumber);
    415603    outputFragment.open(FragmentName, ios::out);
    416     *out << Verbose(2) << "Saving " << fragmentprefix << " No. " << FragmentNumber << "/" << FragmentCounter-1 << " as XYZ ...";
    417     if ((intermediateResult = ListOfMolecules[i]->OutputXYZ(&outputFragment)))
     604    *out << Verbose(2) << "Saving bond fragment No. " << FragmentNumber << "/" << FragmentCounter - 1 << " as XYZ ...";
     605    if ((intermediateResult = (*ListRunner)->OutputXYZ(&outputFragment)))
    418606      *out << " done." << endl;
    419607    else
     
    425613    // list atoms in fragment for debugging
    426614    *out << Verbose(2) << "Contained atoms: ";
    427     Walker = ListOfMolecules[i]->start;
    428     while (Walker->next != ListOfMolecules[i]->end) {
     615    Walker = (*ListRunner)->start;
     616    while (Walker->next != (*ListRunner)->end) {
    429617      Walker = Walker->next;
    430618      *out << Walker->Name << " ";
    431619    }
    432620    *out << endl;
    433    
     621
    434622    // center on edge
    435     if (DoCentering) {
    436       ListOfMolecules[i]->CenterEdge(out, &BoxDimension);
    437       ListOfMolecules[i]->SetBoxDimension(&BoxDimension);  // update Box of atoms by boundary
    438       int j = -1;
    439       for (int k=0;k<NDIM;k++) {
    440         j += k+1;
    441         BoxDimension.x[k] = 2.5* (configuration->GetIsAngstroem() ? 1. : 1./AtomicLengthToAngstroem);
    442         ListOfMolecules[i]->cell_size[j] += BoxDimension.x[k]*2.;
    443       }
    444       ListOfMolecules[i]->Translate(&BoxDimension);
    445     }
     623    (*ListRunner)->CenterEdge(out, &BoxDimension);
     624    (*ListRunner)->SetBoxDimension(&BoxDimension); // update Box of atoms by boundary
     625    int j = -1;
     626    for (int k = 0; k < NDIM; k++) {
     627      j += k + 1;
     628      BoxDimension.x[k] = 2.5 * (configuration->GetIsAngstroem() ? 1. : 1. / AtomicLengthToAngstroem);
     629      (*ListRunner)->cell_size[j] += BoxDimension.x[k] * 2.;
     630    }
     631    (*ListRunner)->Translate(&BoxDimension);
    446632
    447633    // also calculate necessary orbitals
    448     ListOfMolecules[i]->CountElements();  // this is a bugfix, atoms should should actually be added correctly to this fragment
    449     ListOfMolecules[i]->CalculateOrbitals(*configuration);
    450    
     634    (*ListRunner)->CountElements(); // this is a bugfix, atoms should shoulds actually be added correctly to this fragment
     635    (*ListRunner)->CalculateOrbitals(*configuration);
     636
    451637    // change path in config
    452638    //strcpy(PathBackup, configuration->configpath);
    453     sprintf(FragmentName, "%s/%s%s/", PathBackup, fragmentprefix, FragmentNumber);
     639    sprintf(FragmentName, "%s/%s%s/", PathBackup, FRAGMENTPREFIX, FragmentNumber);
    454640    configuration->SetDefaultPath(FragmentName);
    455    
     641
    456642    // and save as config
    457     sprintf(FragmentName, "%s/%s%s.conf", configuration->configpath, fragmentprefix, FragmentNumber);
    458     *out << Verbose(2) << "Saving " << fragmentprefix << " No. " << FragmentNumber << "/" << FragmentCounter-1 << " as config ...";
    459     if ((intermediateResult = configuration->Save(FragmentName, ListOfMolecules[i]->elemente, ListOfMolecules[i])))
     643    sprintf(FragmentName, "%s/%s%s.conf", configuration->configpath, FRAGMENTPREFIX, FragmentNumber);
     644    *out << Verbose(2) << "Saving bond fragment No. " << FragmentNumber << "/" << FragmentCounter - 1 << " as config ...";
     645    if ((intermediateResult = configuration->Save(FragmentName, (*ListRunner)->elemente, (*ListRunner))))
    460646      *out << " done." << endl;
    461647    else
     
    466652    configuration->SetDefaultPath(PathBackup);
    467653
    468 
    469654    // and save as mpqc input file
    470     sprintf(FragmentName, "%s/%s%s.in", configuration->configpath, fragmentprefix, FragmentNumber);
    471     *out << Verbose(2) << "Saving " << fragmentprefix << " No. " << FragmentNumber << "/" << FragmentCounter-1 << " as mpqc input ...";
    472     if ((intermediateResult = configuration->SaveMPQC(FragmentName, ListOfMolecules[i])))
     655    sprintf(FragmentName, "%s/%s%s.conf", configuration->configpath, FRAGMENTPREFIX, FragmentNumber);
     656    *out << Verbose(2) << "Saving bond fragment No. " << FragmentNumber << "/" << FragmentCounter - 1 << " as mpqc input ...";
     657    if ((intermediateResult = configuration->SaveMPQC(FragmentName, (*ListRunner))))
    473658      *out << " done." << endl;
    474659    else
    475660      *out << " failed." << endl;
    476        
     661
    477662    result = result && intermediateResult;
    478663    //outputFragment.close();
    479664    //outputFragment.clear();
    480     delete(FragmentNumber);
     665    delete (FragmentNumber);
    481666    //Free((void **)&FragmentNumber, "MoleculeListClass::OutputConfigForListOfFragments: *FragmentNumber");
    482667  }
    483668  cout << " done." << endl;
    484  
     669
    485670  // printing final number
    486671  *out << "Final number of fragments: " << FragmentCounter << "." << endl;
    487      
     672
    488673  return result;
    489674};
     675
     676/** Counts the number of molecules with the molecule::ActiveFlag set.
     677 * \return number of molecules with ActiveFlag set to true.
     678 */
     679int MoleculeListClass::NumberOfActiveMolecules()
     680{
     681  int count = 0;
     682  for (MoleculeList::iterator ListRunner = ListOfMolecules.begin(); ListRunner != ListOfMolecules.end(); ListRunner++)
     683    count += ((*ListRunner)->ActiveFlag ? 1 : 0);
     684  return count;
     685};
     686
    490687
    491688/******************************************* Class MoleculeLeafClass ************************************************/
     
    498695MoleculeLeafClass::MoleculeLeafClass(MoleculeLeafClass *PreviousLeaf = NULL)
    499696{
    500 //  if (Up != NULL)
    501 //    if (Up->DownLeaf == NULL) // are we the first down leaf for the upper leaf?
    502 //      Up->DownLeaf = this;
    503 //  UpLeaf = Up;
    504 //  DownLeaf = NULL;
     697  //  if (Up != NULL)
     698  //    if (Up->DownLeaf == NULL) // are we the first down leaf for the upper leaf?
     699  //      Up->DownLeaf = this;
     700  //  UpLeaf = Up;
     701  //  DownLeaf = NULL;
    505702  Leaf = NULL;
    506703  previous = PreviousLeaf;
     
    518715MoleculeLeafClass::~MoleculeLeafClass()
    519716{
    520 //  if (DownLeaf != NULL) {// drop leaves further down
    521 //    MoleculeLeafClass *Walker = DownLeaf;
    522 //    MoleculeLeafClass *Next;
    523 //    do {
    524 //      Next = Walker->NextLeaf;
    525 //      delete(Walker);
    526 //      Walker = Next;
    527 //    } while (Walker != NULL);
    528 //    // Last Walker sets DownLeaf automatically to NULL
    529 //  }
     717  //  if (DownLeaf != NULL) {// drop leaves further down
     718  //    MoleculeLeafClass *Walker = DownLeaf;
     719  //    MoleculeLeafClass *Next;
     720  //    do {
     721  //      Next = Walker->NextLeaf;
     722  //      delete(Walker);
     723  //      Walker = Next;
     724  //    } while (Walker != NULL);
     725  //    // Last Walker sets DownLeaf automatically to NULL
     726  //  }
    530727  // remove the leaf itself
    531728  if (Leaf != NULL) {
    532     delete(Leaf);
     729    delete (Leaf);
    533730    Leaf = NULL;
    534731  }
    535732  // remove this Leaf from level list
    536   if (previous != NULL)   
     733  if (previous != NULL)
    537734    previous->next = next;
    538 //  } else { // we are first in list (connects to UpLeaf->DownLeaf)
    539 //    if ((NextLeaf != NULL) && (NextLeaf->UpLeaf == NULL))
    540 //      NextLeaf->UpLeaf = UpLeaf;  // either null as we are top level or the upleaf of the first node
    541 //    if (UpLeaf != NULL)
    542 //      UpLeaf->DownLeaf = NextLeaf;  // either null as we are only leaf or NextLeaf if we are just the first
    543 //  }
    544 //  UpLeaf = NULL;
     735  //  } else { // we are first in list (connects to UpLeaf->DownLeaf)
     736  //    if ((NextLeaf != NULL) && (NextLeaf->UpLeaf == NULL))
     737  //      NextLeaf->UpLeaf = UpLeaf;  // either null as we are top level or the upleaf of the first node
     738  //    if (UpLeaf != NULL)
     739  //      UpLeaf->DownLeaf = NextLeaf;  // either null as we are only leaf or NextLeaf if we are just the first
     740  //  }
     741  //  UpLeaf = NULL;
    545742  if (next != NULL) // are we last in list
    546743    next->previous = previous;
     
    581778    return false;
    582779  }
    583  
     780
    584781  if (status) {
    585     *out << Verbose(1) << "Creating adjacency list for subgraph " << this << "." << endl;
     782    *out << Verbose(1) << "Creating adjacency list for subgraph " << this
     783        << "." << endl;
    586784    Walker = Leaf->start;
    587785    while (Walker->next != Leaf->end) {
    588786      Walker = Walker->next;
    589       AtomNo = Walker->GetTrueFather()->nr;  // global id of the current walker
    590       for(int i=0;i<reference->NumberOfBondsPerAtom[AtomNo];i++) { // go through father's bonds and copy them all
     787      AtomNo = Walker->GetTrueFather()->nr; // global id of the current walker
     788      for (int i = 0; i < reference->NumberOfBondsPerAtom[AtomNo]; i++) { // go through father's bonds and copy them all
    591789        Binder = reference->ListOfBondsPerAtom[AtomNo][i];
    592         OtherWalker = ListOfLocalAtoms[FragmentCounter][Binder->GetOtherAtom(Walker->GetTrueFather())->nr];    // local copy of current bond partner of walker
     790        OtherWalker = ListOfLocalAtoms[FragmentCounter][Binder->GetOtherAtom(Walker->GetTrueFather())->nr]; // local copy of current bond partner of walker
    593791        if (OtherWalker != NULL) {
    594792          if (OtherWalker->nr > Walker->nr)
    595           Leaf->AddBond(Walker, OtherWalker, Binder->BondDegree);
     793            Leaf->AddBond(Walker, OtherWalker, Binder->BondDegree);
    596794        } else {
    597795          *out << Verbose(1) << "OtherWalker = ListOfLocalAtoms[" << FragmentCounter << "][" << Binder->GetOtherAtom(Walker->GetTrueFather())->nr << "] is NULL!" << endl;
     
    606804    FragmentCounter--;
    607805  }
    608  
     806
    609807  if ((FreeList) && (ListOfLocalAtoms != NULL)) {
    610808    // free the index lookup list
    611     Free((void **)&ListOfLocalAtoms[FragmentCounter], "MoleculeLeafClass::FillBondStructureFromReference - **ListOfLocalAtoms[]");
     809    Free((void **) &ListOfLocalAtoms[FragmentCounter], "MoleculeLeafClass::FillBondStructureFromReference - **ListOfLocalAtoms[]");
    612810    if (FragmentCounter == 0) // first fragments frees the initial pointer to list
    613       Free((void **)&ListOfLocalAtoms, "MoleculeLeafClass::FillBondStructureFromReference - ***ListOfLocalAtoms");
     811      Free((void **) &ListOfLocalAtoms, "MoleculeLeafClass::FillBondStructureFromReference - ***ListOfLocalAtoms");
    614812  }
    615813  FragmentCounter--;
     
    626824 * \return true - stack is non-empty, fragmentation necessary, false - stack is empty, no more sites to update
    627825 */
    628 bool MoleculeLeafClass::FillRootStackForSubgraphs(ofstream *out, KeyStack *&RootStack, bool *AtomMask, int &FragmentCounter)
     826bool MoleculeLeafClass::FillRootStackForSubgraphs(ofstream *out,
     827    KeyStack *&RootStack, bool *AtomMask, int &FragmentCounter)
    629828{
    630829  atom *Walker = NULL, *Father = NULL;
    631830
    632831  if (RootStack != NULL) {
    633     // find first root candidates 
     832    // find first root candidates
    634833    if (&(RootStack[FragmentCounter]) != NULL) {
    635       RootStack[FragmentCounter].clear(); 
     834      RootStack[FragmentCounter].clear();
    636835      Walker = Leaf->start;
    637836      while (Walker->next != Leaf->end) { // go through all (non-hydrogen) atoms
     
    639838        Father = Walker->GetTrueFather();
    640839        if (AtomMask[Father->nr]) // apply mask
    641     #ifdef ADDHYDROGEN
     840#ifdef ADDHYDROGEN
    642841          if (Walker->type->Z != 1) // skip hydrogen
    643     #endif
    644             RootStack[FragmentCounter].push_front(Walker->nr);
     842#endif
     843          RootStack[FragmentCounter].push_front(Walker->nr);
    645844      }
    646845      if (next != NULL)
    647846        next->FillRootStackForSubgraphs(out, RootStack, AtomMask, ++FragmentCounter);
    648     }  else {
    649       *out << Verbose(1) << "Rootstack[" << FragmentCounter  << "] is NULL." << endl;
     847    } else {
     848      *out << Verbose(1) << "Rootstack[" << FragmentCounter << "] is NULL." << endl;
    650849      return false;
    651850    }
     
    669868{
    670869  bool status = true;
    671  
     870
    672871  int Counter = Count();
    673872  if (ListOfLocalAtoms == NULL) { // allocated initial pointer
    674873    // allocate and set each field to NULL
    675     ListOfLocalAtoms = (atom ***) Malloc(sizeof(atom **)*Counter, "MoleculeLeafClass::FillBondStructureFromReference - ***ListOfLocalAtoms");
     874    ListOfLocalAtoms = (atom ***) Malloc(sizeof(atom **) * Counter, "MoleculeLeafClass::FillBondStructureFromReference - ***ListOfLocalAtoms");
    676875    if (ListOfLocalAtoms != NULL) {
    677       for (int i=Counter;i--;)
     876      for (int i = Counter; i--;)
    678877        ListOfLocalAtoms[i] = NULL;
    679878      FreeList = FreeList && true;
     
    681880      status = false;
    682881  }
    683  
     882
    684883  if ((ListOfLocalAtoms != NULL) && (ListOfLocalAtoms[FragmentCounter] == NULL)) { // allocate and fill list of this fragment/subgraph
    685884    status = status && CreateFatherLookupTable(out, Leaf->start, Leaf->end, ListOfLocalAtoms[FragmentCounter], GlobalAtomCount);
    686885    FreeList = FreeList && true;
    687886  }
    688  
     887
    689888  return status;
    690889};
     
    700899 * \retuen true - success, false - failure
    701900 */
    702 bool MoleculeLeafClass::AssignKeySetsToFragment(ofstream *out, molecule *reference, Graph *KeySetList, atom ***&ListOfLocalAtoms, Graph **&FragmentList, int &FragmentCounter, bool FreeList)
     901bool MoleculeLeafClass::AssignKeySetsToFragment(ofstream *out,
     902    molecule *reference, Graph *KeySetList, atom ***&ListOfLocalAtoms,
     903    Graph **&FragmentList, int &FragmentCounter, bool FreeList)
    703904{
    704905  bool status = true;
    705906  int KeySetCounter = 0;
    706  
     907
    707908  *out << Verbose(1) << "Begin of AssignKeySetsToFragment." << endl;
    708909  // fill ListOfLocalAtoms if NULL was given
     
    715916  if (FragmentList == NULL) {
    716917    KeySetCounter = Count();
    717     FragmentList = (Graph **) Malloc(sizeof(Graph *)*KeySetCounter, "MoleculeLeafClass::AssignKeySetsToFragment - **FragmentList");
    718     for(int i=KeySetCounter;i--;)
     918    FragmentList = (Graph **) Malloc(sizeof(Graph *) * KeySetCounter, "MoleculeLeafClass::AssignKeySetsToFragment - **FragmentList");
     919    for (int i = KeySetCounter; i--;)
    719920      FragmentList[i] = NULL;
    720921    KeySetCounter = 0;
    721922  }
    722  
     923
    723924  if ((KeySetList != NULL) && (KeySetList->size() != 0)) { // if there are some scanned keysets at all
    724925    // assign scanned keysets
     
    726927      FragmentList[FragmentCounter] = new Graph;
    727928    KeySet *TempSet = new KeySet;
    728     for(Graph::iterator runner = KeySetList->begin();runner != KeySetList->end(); runner++) { // key sets contain global numbers!
    729       if ( ListOfLocalAtoms[FragmentCounter][reference->FindAtom(*((*runner).first.begin()))->nr] != NULL) {// as we may assume that that bond structure is unchanged, we only test the first key in each set
     929    for (Graph::iterator runner = KeySetList->begin(); runner != KeySetList->end(); runner++) { // key sets contain global numbers!
     930      if (ListOfLocalAtoms[FragmentCounter][reference->FindAtom(*((*runner).first.begin()))->nr] != NULL) {// as we may assume that that bond structure is unchanged, we only test the first key in each set
    730931        // translate keyset to local numbers
    731         for(KeySet::iterator sprinter = (*runner).first.begin(); sprinter != (*runner).first.end(); sprinter++)
     932        for (KeySet::iterator sprinter = (*runner).first.begin(); sprinter != (*runner).first.end(); sprinter++)
    732933          TempSet->insert(ListOfLocalAtoms[FragmentCounter][reference->FindAtom(*sprinter)->nr]->nr);
    733934        // insert into FragmentList
    734         FragmentList[FragmentCounter]->insert(GraphPair (*TempSet, pair<int,double>(KeySetCounter++, (*runner).second.second)));
     935        FragmentList[FragmentCounter]->insert(GraphPair(*TempSet, pair<int, double> (KeySetCounter++, (*runner).second.second)));
    735936      }
    736937      TempSet->clear();
    737938    }
    738     delete(TempSet);
     939    delete (TempSet);
    739940    if (KeySetCounter == 0) {// if there are no keysets, delete the list
    740941      *out << Verbose(1) << "KeySetCounter is zero, deleting FragmentList." << endl;
    741       delete(FragmentList[FragmentCounter]);
     942      delete (FragmentList[FragmentCounter]);
    742943    } else
    743944      *out << Verbose(1) << KeySetCounter << " keysets were assigned to subgraph " << FragmentCounter << "." << endl;
     
    748949  } else
    749950    *out << Verbose(1) << "KeySetList is NULL or empty." << endl;
    750  
     951
    751952  if ((FreeList) && (ListOfLocalAtoms != NULL)) {
    752953    // free the index lookup list
    753     Free((void **)&ListOfLocalAtoms[FragmentCounter], "MoleculeLeafClass::AssignKeySetsToFragment - **ListOfLocalAtoms[]");
     954    Free((void **) &ListOfLocalAtoms[FragmentCounter], "MoleculeLeafClass::AssignKeySetsToFragment - **ListOfLocalAtoms[]");
    754955    if (FragmentCounter == 0) // first fragments frees the initial pointer to list
    755       Free((void **)&ListOfLocalAtoms, "MoleculeLeafClass::AssignKeySetsToFragment - ***ListOfLocalAtoms");
     956      Free((void **) &ListOfLocalAtoms, "MoleculeLeafClass::AssignKeySetsToFragment - ***ListOfLocalAtoms");
    756957  }
    757958  *out << Verbose(1) << "End of AssignKeySetsToFragment." << endl;
     
    765966 * \param &TotalNumberOfKeySets global key set counter
    766967 * \param &TotalGraph Graph to be filled with global numbers
    767  */
    768 void MoleculeLeafClass::TranslateIndicesToGlobalIDs(ofstream *out, Graph **FragmentList, int &FragmentCounter, int &TotalNumberOfKeySets, Graph &TotalGraph)
     968 */
     969void MoleculeLeafClass::TranslateIndicesToGlobalIDs(ofstream *out,
     970    Graph **FragmentList, int &FragmentCounter, int &TotalNumberOfKeySets,
     971    Graph &TotalGraph)
    769972{
    770973  *out << Verbose(1) << "Begin of TranslateIndicesToGlobalIDs." << endl;
    771974  KeySet *TempSet = new KeySet;
    772975  if (FragmentList[FragmentCounter] != NULL) {
    773     for(Graph::iterator runner = FragmentList[FragmentCounter]->begin(); runner != FragmentList[FragmentCounter]->end(); runner++) {
    774       for(KeySet::iterator sprinter = (*runner).first.begin(); sprinter != (*runner).first.end(); sprinter++)
     976    for (Graph::iterator runner = FragmentList[FragmentCounter]->begin(); runner != FragmentList[FragmentCounter]->end(); runner++) {
     977      for (KeySet::iterator sprinter = (*runner).first.begin(); sprinter != (*runner).first.end(); sprinter++)
    775978        TempSet->insert((Leaf->FindAtom(*sprinter))->GetTrueFather()->nr);
    776       TotalGraph.insert(GraphPair(*TempSet, pair<int,double>(TotalNumberOfKeySets++, (*runner).second.second)));
     979      TotalGraph.insert(GraphPair(*TempSet, pair<int, double> (TotalNumberOfKeySets++, (*runner).second.second)));
    777980      TempSet->clear();
    778981    }
    779     delete(TempSet);
     982    delete (TempSet);
    780983  } else {
    781984    *out << Verbose(1) << "FragmentList is NULL." << endl;
     
    793996{
    794997  if (next != NULL)
    795     return next->Count()+1;
     998    return next->Count() + 1;
    796999  else
    797     return 1;
    798 };
     1000    return 1;
     1001};
     1002
Note: See TracChangeset for help on using the changeset viewer.