Changeset df793a for src


Ignore:
Timestamp:
Oct 9, 2008, 6:27:56 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:
ad85c2
Parents:
ba4432
Message:

BUGFIXES: CyclicStructureAnalysis() now compatible to disconnected subgraphs, AssignKeySetsToFragment() and FillBondStructureFromReference() memory cleanup corrected

+ molecule::DepthFirstSearchAnalysis() now just returns BackEdgeStack, not MinimumRingSize. CyclicStructureAnalysis() is called during FragmentMolecule(), after subgraphs bonds list have been filled by FillBondStructureFromReference().
+ new function molecule::PickLocalBackEdges(), as the BackEdgeStack returned by DepthFirstSearchAnalysis() co
ntains only global bonds, not the local ones for the subgraph, we have to step through it and pick the right
ones out.
+ molecule::FragmentMolecule() now calls molecule::CyclicStructureAnalysis() separately for each subgraph, along with a BackEdgeStack filled by PickLocalBackEdges(), and allocates&initialises MinimumRingSize array. Als
o AssignKeySetsToFragment() frees the LocalListOfAtoms now (FreeList=true), now longer after the following wh
ile
+ molecule::CyclicStructureAnalysis() takes a local BackEdgeStack and analysis the subgraphs cycles, returnin
g minimum ring size
+ MoleculeLeafClass::AssignKeySetsToFragment() now frees memory for ListOfLocalAtoms when FreeList is set. BUGFIX: test of first key was testing against ..->nr != -1. However, LocalListOfAtoms was not even initialised correctly to NULL, hence ...->nr pointed in some cases to nowhere. Now it test atom* against NULL.
+ MoleculeLeafClass::FillBondStructureFromReference() now frees memory for ListOfLocalAtoms when FreeList is set correctly (only free initial pointer when FragmentCounter == 0, also it was decreased not before but after freeing, hence we free'd the wrong list). Also, father replaced by GetTrueFather() (makes the function moregenerally useable, was not a bug).
+ ParseCommandLineOptions() option 'D': adapted to changes in DepthFirstSearchAnalysis() in a similar manner
to FragmentMolecule()
+ molecule::CountCyclicBonds() adapted but does not perform CyclicStructureAnalysis()
+ molecule::CreateAdjacencyList() counts the bonds that could not be brought to covalently corrected degree (i.e. the remaining ionic atoms)
+ molecule::CreateListOfBondsPerAtom() prints atom names and number, which is helpful as name contains global

and number contains local number (helped in finding above bugs)

+ CreateFatherLookupTable(): BUGFIX: LookupTable was not initialised to NULL (see above)

Location:
src
Files:
5 edited

Legend:

Unmodified
Added
Removed

  • src/builder.cpp

    rba4432 rdf793a  
    953953                cout << Verbose(1) << "Depth-First-Search Analysis." << endl;
    954954                MoleculeLeafClass *Subgraphs = NULL;      // list of subgraphs from DFS analysis
    955                 int *MinimumRingSize = NULL;
     955                int *MinimumRingSize = new int[mol->AtomCount];
     956                atom ***ListOfLocalAtoms = NULL;
     957                int FragmentCounter = 0;
     958                class StackClass<bond *> *BackEdgeStack = NULL;
     959                class StackClass<bond *> *LocalBackEdgeStack = NULL;
    956960                mol->CreateAdjacencyList((ofstream *)&cout, atof(argv[argptr]), configuration.GetIsAngstroem());
    957961                mol->CreateListOfBondsPerAtom((ofstream *)&cout);
    958                 Subgraphs = mol->DepthFirstSearchAnalysis((ofstream *)&cout, MinimumRingSize);
    959                 delete[](MinimumRingSize);
     962                Subgraphs = mol->DepthFirstSearchAnalysis((ofstream *)&cout, BackEdgeStack);
    960963                if (Subgraphs != NULL) {
     964                  Subgraphs->next->FillBondStructureFromReference((ofstream *)&cout, mol, (FragmentCounter = 0), ListOfLocalAtoms, false);  // we want to keep the created ListOfLocalAtoms
    961965                  while (Subgraphs->next != NULL) {
    962966                    Subgraphs = Subgraphs->next;
     967                    LocalBackEdgeStack = new StackClass<bond *> (Subgraphs->Leaf->BondCount);
     968                    Subgraphs->Leaf->PickLocalBackEdges((ofstream *)&cout, ListOfLocalAtoms[FragmentCounter++], BackEdgeStack, LocalBackEdgeStack);
     969                    Subgraphs->Leaf->CyclicStructureAnalysis((ofstream *)&cout, BackEdgeStack, MinimumRingSize);
     970                    delete(LocalBackEdgeStack);
    963971                    delete(Subgraphs->previous);
    964972                  }
    965973                  delete(Subgraphs);
     974                  for (int i=0;i<FragmentCounter;i++)
     975                    Free((void **)&ListOfLocalAtoms[FragmentCounter], "ParseCommandLineOptions: **ListOfLocalAtoms[]");
     976                  Free((void **)&ListOfLocalAtoms, "ParseCommandLineOptions: ***ListOfLocalAtoms");
    966977                }
     978                delete(BackEdgeStack);
     979                delete[](MinimumRingSize);
    967980              }
    968981              argptr+=1;
    … …  
    12551268  int Z;
    12561269  int j, axis, count, faktor;
    1257   int *MinimumRingSize = NULL;
     1270  clock_t start,end;
     1271//  int *MinimumRingSize = NULL;
    12581272  MoleculeLeafClass *Subgraphs = NULL;
    1259   clock_t start,end;
     1273//  class StackClass<bond *> *BackEdgeStack = NULL;
    12601274  element **Elements;
    12611275  Vector **vectors;
    … …  
    14231437       
    14241438      case 'o': // create the connection matrix
    1425         cout << Verbose(0) << "What's the maximum bond distance: ";
    1426         cin >> tmp1;
    1427         start = clock();
    1428         mol->CreateAdjacencyList((ofstream *)&cout, tmp1, configuration.GetIsAngstroem());
    1429         //mol->CreateListOfBondsPerAtom((ofstream *)&cout);
    1430         Subgraphs = mol->DepthFirstSearchAnalysis((ofstream *)&cout, MinimumRingSize);
    1431         while (Subgraphs->next != NULL) {
    1432           Subgraphs = Subgraphs->next;
    1433           delete(Subgraphs->previous);
     1439        {
     1440          cout << Verbose(0) << "What's the maximum bond distance: ";
     1441          cin >> tmp1;
     1442          start = clock();
     1443          mol->CreateAdjacencyList((ofstream *)&cout, tmp1, configuration.GetIsAngstroem());
     1444          mol->CreateListOfBondsPerAtom((ofstream *)&cout);
     1445//          Subgraphs = mol->DepthFirstSearchAnalysis((ofstream *)&cout, BackEdgeStack);
     1446//          while (Subgraphs->next != NULL) {
     1447//            Subgraphs = Subgraphs->next;
     1448//            Subgraphs->Leaf->CyclicStructureAnalysis((ofstream *)&cout, BackEdgeStack, MinimumRingSize);
     1449//            delete(Subgraphs->previous);
     1450//          }
     1451//          delete(Subgraphs);    // we don't need the list here, so free everything
     1452//          delete[](MinimumRingSize);
     1453//          Subgraphs = NULL;
     1454          end = clock();
     1455          cout << Verbose(0) << "Clocks for this operation: " << (end-start) << ", time: " << ((double)(end-start)/CLOCKS_PER_SEC) << "s." << endl;
    14341456        }
    1435         delete(Subgraphs);    // we don't need the list here, so free everything
    1436         delete[](MinimumRingSize);
    1437         Subgraphs = NULL;
    1438         end = clock();
    1439         cout << Verbose(0) << "Clocks for this operation: " << (end-start) << ", time: " << ((double)(end-start)/CLOCKS_PER_SEC) << "s." << endl;
    14401457        break;
    14411458       

  • src/helpers.hpp

    rba4432 rdf793a  
    100100    status = false;
    101101  } else {
     102    for (int i=0;i<count;i++)
     103      LookupTable[i] = NULL;
    102104    Walker = start;
    103105    while (Walker->next != end) { // create a lookup table (Atom::nr -> atom) used as a marker table lateron

  • src/moleculelist.cpp

    rba4432 rdf793a  
    589589    while (Walker->next != Leaf->end) {
    590590      Walker = Walker->next;
    591       AtomNo = Walker->father->nr;  // global id of the current walker
     591      AtomNo = Walker->GetTrueFather()->nr;  // global id of the current walker
    592592      for(int i=0;i<reference->NumberOfBondsPerAtom[AtomNo];i++) { // go through father's bonds and copy them all
    593593        Binder = reference->ListOfBondsPerAtom[AtomNo][i];
    594         OtherWalker = ListOfLocalAtoms[FragmentCounter][Binder->GetOtherAtom(Walker->father)->nr];    // local copy of current bond partner of walker
     594        OtherWalker = ListOfLocalAtoms[FragmentCounter][Binder->GetOtherAtom(Walker->GetTrueFather())->nr];    // local copy of current bond partner of walker
    595595        if (OtherWalker != NULL) {
    596596          if (OtherWalker->nr > Walker->nr)
    597597          Leaf->AddBond(Walker, OtherWalker, Binder->BondDegree);
    598598        } else {
    599           *out << Verbose(1) << "OtherWalker = ListOfLocalAtoms[" << FragmentCounter << "][" << Binder->GetOtherAtom(Walker->father)->nr << "] is NULL!" << endl;
     599          *out << Verbose(1) << "OtherWalker = ListOfLocalAtoms[" << FragmentCounter << "][" << Binder->GetOtherAtom(Walker->GetTrueFather())->nr << "] is NULL!" << endl;
    600600          status = false;
    601601        }
    … …  
    608608  }
    609609 
     610  FragmentCounter--;
    610611  if (FreeList) {
    611612    // free the index lookup list
    612613    Free((void **)&ListOfLocalAtoms[FragmentCounter], "MoleculeLeafClass::FillBondStructureFromReference - **ListOfLocalAtoms[]");
    613     if (ListOfLocalAtoms[FragmentCounter] == NULL)
     614    if (FragmentCounter == 0) // first fragments frees the initial pointer to list
    614615      Free((void **)&ListOfLocalAtoms, "MoleculeLeafClass::FillBondStructureFromReference - ***ListOfLocalAtoms");
    615616  }
    … …  
    661662 * \param *out output stream fro debugging
    662663 * \param ***ListOfLocalAtoms Lookup table for each subgraph and index of each atom in global molecule, may be NULL on start, then it is filled
    663  * \param &FragmentCounter counts the fragments as we move along the list
     664 * \param FragmentCounter counts the fragments as we move along the list
    664665 * \param GlobalAtomCount number of atoms in the complete molecule
    665666 * \param &FreeList true - ***ListOfLocalAtoms is free'd before return, false - it is not
    666667 * \return true - succes, false - failure
    667668 */
    668 bool MoleculeLeafClass::FillListOfLocalAtoms(ofstream *out, atom ***&ListOfLocalAtoms, int &FragmentCounter, int GlobalAtomCount, bool &FreeList)
     669bool MoleculeLeafClass::FillListOfLocalAtoms(ofstream *out, atom ***&ListOfLocalAtoms, const int FragmentCounter, const int GlobalAtomCount, bool &FreeList)
    669670{
    670671  bool status = true;
    … …  
    726727    KeySet *TempSet = new KeySet;
    727728    for(Graph::iterator runner = KeySetList->begin();runner != KeySetList->end(); runner++) { // key sets contain global numbers!
    728       if ( ListOfLocalAtoms[FragmentCounter][reference->FindAtom(*((*runner).first.begin()))->nr]->nr != -1) {// as we may assume that that bond structure is unchanged, we only test the first key in each set
     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
    729730        // translate keyset to local numbers
    730731        for(KeySet::iterator sprinter = (*runner).first.begin(); sprinter != (*runner).first.end(); sprinter++)
    … …  
    744745    if (next != NULL)
    745746      next->AssignKeySetsToFragment(out, reference, KeySetList, ListOfLocalAtoms, FragmentList, FragmentCounter, FreeList);
    746     FragmentCounter--;
    747747  } else
    748748    *out << Verbose(1) << "KeySetList is NULL or empty." << endl;
    749749 
     750  FragmentCounter--;
     751  if (FreeList) {
     752    // free the index lookup list
     753    Free((void **)&ListOfLocalAtoms[FragmentCounter], "MoleculeLeafClass::AssignKeySetsToFragment - **ListOfLocalAtoms[]");
     754    if (FragmentCounter == 0) // first fragments frees the initial pointer to list
     755      Free((void **)&ListOfLocalAtoms, "MoleculeLeafClass::AssignKeySetsToFragment - ***ListOfLocalAtoms");
     756  }
    750757  return status;
    751758};

  • src/molecules.cpp

    rba4432 rdf793a  
    22442244  int *MinimumRingSize = NULL;
    22452245  MoleculeLeafClass *Subgraphs = NULL;
     2246  class StackClass<bond *> *BackEdgeStack = NULL;
    22462247  bond *Binder = first;
    22472248  if ((Binder->next != last) && (Binder->next->Type == Undetermined)) {
    22482249    *out << Verbose(0) << "No Depth-First-Search analysis performed so far, calling ..." << endl;
    2249     Subgraphs = DepthFirstSearchAnalysis(out, MinimumRingSize);
     2250    Subgraphs = DepthFirstSearchAnalysis(out, BackEdgeStack);
    22502251    while (Subgraphs->next != NULL) {
    22512252      Subgraphs = Subgraphs->next;
    … …  
    22602261      No++;   
    22612262  }
     2263  delete(BackEdgeStack);
    22622264  return No;
    22632265};
    … …  
    23412343  double *matrix = ReturnFullMatrixforSymmetric(cell_size);
    23422344  Vector x;
     2345  int FalseBondDegree = 0;
    23432346 
    23442347  BondDistance = bonddistance; // * ((IsAngstroem) ? 1. : 1./AtomicLengthToAngstroem);
    … …  
    23592362      j += i+1;
    23602363      divisor[i] = (int)floor(cell_size[j]/bonddistance); // take smaller value such that size of linked cell is at least bonddistance
    2361       *out << Verbose(1) << "divisor[" << i << "]  = " << divisor[i] << "." << endl;
     2364      //*out << Verbose(1) << "divisor[" << i << "]  = " << divisor[i] << "." << endl;
    23622365    }
    23632366    // 2a. allocate memory for the cell list
    … …  
    23842387      }
    23852388      index = n[2] + (n[1] + n[0] * divisor[1]) * divisor[2];
    2386       *out << Verbose(1) << "Atom " << *Walker << " goes into cell number [" << n[0] << "," << n[1] << "," << n[2] << "] = " << index << "." << endl;
     2389      //*out << Verbose(1) << "Atom " << *Walker << " goes into cell number [" << n[0] << "," << n[1] << "," << n[2] << "] = " << index << "." << endl;
    23872390      // add copy atom to this cell
    23882391      if (CellList[index] == NULL)  // allocate molecule if not done
    … …  
    24262429                        distance =   OtherWalker->x.PeriodicDistance(&(Walker->x), cell_size);
    24272430                        if ((OtherWalker->father->nr > Walker->father->nr) && (distance <= MaxDistance*MaxDistance) && (distance >= MinDistance*MinDistance)) { // create bond if distance is smaller
    2428                           *out << Verbose(0) << "Adding Bond between " << *Walker << " and " << *OtherWalker << "." << endl;
     2431                          //*out << Verbose(0) << "Adding Bond between " << *Walker << " and " << *OtherWalker << "." << endl;
    24292432                          AddBond(Walker->father, OtherWalker->father, 1);  // also increases molecule::BondCount
    24302433                          BondCount++;
    … …  
    24852488              ListOfBondsPerAtom[Walker->nr][CandidateBondNo]->BondDegree++;
    24862489              *out << Verbose(2) << "Increased bond degree for bond " << *ListOfBondsPerAtom[Walker->nr][CandidateBondNo] << "." << endl;
    2487             }
     2490            } else
     2491              *out << Verbose(2) << "Could not find correct degree for bond " << *ListOfBondsPerAtom[Walker->nr][CandidateBondNo] << "." << endl;
     2492              FalseBondDegree++;
    24882493          }
    24892494                    }
    … …  
    24922497                } else
    24932498                        *out << Verbose(1) << "BondCount is " << BondCount << ", no bonds between any of the " << AtomCount << " atoms." << endl;
    2494           *out << Verbose(1) << "I detected " << BondCount << " bonds in the molecule with distance " << bonddistance << "." << endl;
     2499          *out << Verbose(1) << "I detected " << BondCount << " bonds in the molecule with distance " << bonddistance << ", " << FalseBondDegree << " bonds could not be corrected." << endl;
    24952500               
    24962501          // output bonds for debugging (if bond chain list was correctly installed)
    … …  
    25132518 * We use the algorithm from [Even, Graph Algorithms, p.62].
    25142519 * \param *out output stream for debugging
    2515  * \param *&MinimumRingSize contains smallest ring size in molecular structure on return or -1 if no rings were found
     2520 * \param *&BackEdgeStack NULL pointer to StackClass with all the found back edges, allocated and filled on return
    25162521 * \return list of each disconnected subgraph as an individual molecule class structure
    25172522 */
    2518 MoleculeLeafClass * molecule::DepthFirstSearchAnalysis(ofstream *out, int *&MinimumRingSize)
    2519 {
    2520   class StackClass<atom *> *AtomStack;
    2521   AtomStack = new StackClass<atom *>(AtomCount);
    2522   class StackClass<bond *> *BackEdgeStack = new StackClass<bond *> (BondCount);
     2523MoleculeLeafClass * molecule::DepthFirstSearchAnalysis(ofstream *out, class StackClass<bond *> *&BackEdgeStack)
     2524{
     2525  class StackClass<atom *> *AtomStack = new StackClass<atom *>(AtomCount);
     2526  BackEdgeStack = new StackClass<bond *> (BondCount);
    25232527  MoleculeLeafClass *SubGraphs = new MoleculeLeafClass(NULL);
    25242528  MoleculeLeafClass *LeafWalker = SubGraphs;
    … …  
    26522656    LeafWalker->Leaf->Output(out);
    26532657    *out << endl;
    2654    
     2658
    26552659    // step on to next root
    26562660    while ((Root != end) && (Root->GraphNr != -1)) {
    … …  
    26702674  }
    26712675
    2672   // analysis of the cycles (print rings, get minimum cycle length)
    2673   CyclicStructureAnalysis(out, BackEdgeStack, MinimumRingSize);
    26742676
    26752677  *out << Verbose(1) << "Final graph info for each atom is:" << endl;
    … …  
    27092711 * as cyclic and print out the cycles.
    27102712 * \param *out output stream for debugging
    2711  * \param *BackEdgeStack stack with all back edges found during DFS scan
     2713 * \param *BackEdgeStack stack with all back edges found during DFS scan. Beware: This stack contains the bonds from the total molecule, not from the subgraph!
    27122714 * \param *&MinimumRingSize contains smallest ring size in molecular structure on return or -1 if no rings were found, if set is maximum search distance
    27132715 * \todo BFS from the not-same-LP to find back to starting point of tributary cycle over more than one bond
    … …  
    27302732    ColorList[i] = white;
    27312733  }
    2732   MinimumRingSize = new int[AtomCount];
    2733   for(int i=AtomCount;i--;)
    2734     MinimumRingSize[i] = AtomCount;
    2735 
    27362734 
    27372735  *out << Verbose(1) << "Back edge list - ";
    … …  
    35353533        MoleculeListClass *BondFragments = NULL;
    35363534  int *SortIndex = NULL;
    3537   int *MinimumRingSize = NULL;
     3535  int *MinimumRingSize = new int[AtomCount];
    35383536  int FragmentCounter;
    35393537  MoleculeLeafClass *MolecularWalker = NULL;
    … …  
    35413539  fstream File;
    35423540  bool FragmentationToDo = true;
     3541  class StackClass<bond *> *BackEdgeStack = NULL, *LocalBackEdgeStack = NULL;
    35433542  bool CheckOrder = false;
    35443543  Graph **FragmentList = NULL;
    … …  
    35723571
    35733572  // ===== 2. perform a DFS analysis to gather info on cyclic structure and a list of disconnected subgraphs =====
    3574   Subgraphs = DepthFirstSearchAnalysis(out, MinimumRingSize);
     3573  Subgraphs = DepthFirstSearchAnalysis(out, BackEdgeStack);
    35753574  // fill the bond structure of the individually stored subgraphs
    35763575  Subgraphs->next->FillBondStructureFromReference(out, this, (FragmentCounter = 0), ListOfLocalAtoms, false);  // we want to keep the created ListOfLocalAtoms
    3577 
     3576  // analysis of the cycles (print rings, get minimum cycle length) for each subgraph
     3577  for(int i=AtomCount;i--;)
     3578    MinimumRingSize[i] = AtomCount;
     3579  MolecularWalker = Subgraphs;
     3580  FragmentCounter = 0;
     3581  while (MolecularWalker->next != NULL) {
     3582    MolecularWalker = MolecularWalker->next;
     3583    *out << Verbose(0) << "Analysing the cycles of subgraph " << MolecularWalker->Leaf << " with nr. " << FragmentCounter << "." << endl;
     3584    LocalBackEdgeStack = new StackClass<bond *> (MolecularWalker->Leaf->BondCount);
     3585//    // check the list of local atoms for debugging
     3586//    *out << Verbose(0) << "ListOfLocalAtoms for this subgraph is:" << endl;
     3587//    for (int i=0;i<AtomCount;i++)
     3588//      if (ListOfLocalAtoms[FragmentCounter][i] == NULL)
     3589//        *out << "\tNULL";
     3590//      else
     3591//        *out << "\t" << ListOfLocalAtoms[FragmentCounter][i]->Name;
     3592    MolecularWalker->Leaf->PickLocalBackEdges(out, ListOfLocalAtoms[FragmentCounter++], BackEdgeStack, LocalBackEdgeStack);
     3593    MolecularWalker->Leaf->CyclicStructureAnalysis(out, LocalBackEdgeStack, MinimumRingSize);
     3594    delete(LocalBackEdgeStack);
     3595  }
     3596 
    35783597  // ===== 3. if structure still valid, parse key set file and others =====
    35793598  FragmentationToDo = FragmentationToDo && ParseKeySetFile(out, configuration->configpath, ParsedFragmentList);
    … …  
    35843603  // =================================== Begin of FRAGMENTATION ===============================
    35853604  // ===== 6a. assign each keyset to its respective subgraph =====
    3586   Subgraphs->next->AssignKeySetsToFragment(out, this, ParsedFragmentList, ListOfLocalAtoms, FragmentList, (FragmentCounter = 0), false);
     3605  Subgraphs->next->AssignKeySetsToFragment(out, this, ParsedFragmentList, ListOfLocalAtoms, FragmentList, (FragmentCounter = 0), true);
    35873606
    35883607  // ===== 6b. prepare and go into the adaptive (Order<0), single-step (Order==0) or incremental (Order>0) cycle
    … …  
    36213640  delete[](MinimumRingSize);
    36223641 
    3623   // free the index lookup list
    3624   for (int i=FragmentCounter;i--;)
    3625     Free((void **)&ListOfLocalAtoms[i], "molecule::FragmentMolecule - *ListOfLocalAtoms[]");
    3626   Free((void **)&ListOfLocalAtoms, "molecule::FragmentMolecule - **ListOfLocalAtoms");
    36273642
    36283643  // ==================================== End of FRAGMENTATION ============================================
    … …  
    36983713};
    36993714
     3715
     3716/** Picks from a global stack with all back edges the ones in the fragment.
     3717 * \param *out output stream for debugging
     3718 * \param **ListOfLocalAtoms array of father atom::nr to local atom::nr (reverse of atom::father)
     3719 * \param *ReferenceStack stack with all the back egdes
     3720 * \param *LocalStack stack to be filled
     3721 * \return true - everything ok, false - ReferenceStack was empty
     3722 */
     3723bool molecule::PickLocalBackEdges(ofstream *out, atom **ListOfLocalAtoms, class StackClass<bond *> *&ReferenceStack, class StackClass<bond *> *&LocalStack)
     3724{
     3725  bool status = true;
     3726  if (ReferenceStack->IsEmpty()) {
     3727    cerr << "ReferenceStack is empty!" << endl;
     3728    return false;
     3729  }
     3730  bond *Binder = ReferenceStack->PopFirst();
     3731  bond *FirstBond = Binder;   // mark the first bond, so that we don't loop through the stack indefinitely
     3732  atom *Walker = NULL, *OtherAtom = NULL;
     3733  ReferenceStack->Push(Binder);
     3734 
     3735  do {  // go through all bonds and push local ones
     3736    Walker = ListOfLocalAtoms[Binder->leftatom->nr];  // get one atom in the reference molecule
     3737    if (Walker == NULL) // if this Walker exists in the subgraph ...
     3738      continue;
     3739    for(int i=0;i<NumberOfBondsPerAtom[Walker->nr];i++) {    // go through the local list of bonds
     3740      OtherAtom = ListOfBondsPerAtom[Walker->nr][i]->GetOtherAtom(Walker);
     3741      if (OtherAtom == ListOfLocalAtoms[Binder->rightatom->nr]) { // found the bond
     3742        LocalStack->Push(ListOfBondsPerAtom[Walker->nr][i]);
     3743        break;
     3744      }
     3745    }
     3746    Binder = ReferenceStack->PopFirst();  // loop the stack for next item
     3747    ReferenceStack->Push(Binder);
     3748  } while (FirstBond != Binder);
     3749 
     3750  return status;
     3751};
     3752
    37003753/** Stores pairs (Atom::nr, Atom::AdaptiveOrder) into file.
    37013754 * Atoms not present in the file get "-1".
    … …  
    38403893  while (Walker->next != end) {
    38413894    Walker = Walker->next;
    3842     *out << Verbose(4) << "Atom " << Walker->Name << " with " << NumberOfBondsPerAtom[Walker->nr] << " bonds: ";
     3895    *out << Verbose(4) << "Atom " << Walker->Name << "/" << Walker->nr << " with " << NumberOfBondsPerAtom[Walker->nr] << " bonds: ";
    38433896    TotalDegree = 0;
    38443897    for (int j=0;j<NumberOfBondsPerAtom[Walker->nr];j++) {

  • src/molecules.hpp

    rba4432 rdf793a  
    284284 
    285285  // Graph analysis
    286   MoleculeLeafClass * DepthFirstSearchAnalysis(ofstream *out, int *&MinimumRingSize);
     286  MoleculeLeafClass * DepthFirstSearchAnalysis(ofstream *out, class StackClass<bond *> *&BackEdgeStack);
    287287  void CyclicStructureAnalysis(ofstream *out, class StackClass<bond *> *BackEdgeStack, int *&MinimumRingSize);
     288  bool PickLocalBackEdges(ofstream *out, atom **ListOfLocalAtoms, class StackClass<bond *> *&ReferenceStack, class StackClass<bond *> *&LocalStack);
    288289  bond * FindNextUnused(atom *vertex);
    289290  void SetNextComponentNumber(atom *vertex, int nr);
    … …  
    378379  bool FillRootStackForSubgraphs(ofstream *out, KeyStack *&RootStack, bool *AtomMask, int &FragmentCounter);
    379380  bool AssignKeySetsToFragment(ofstream *out, molecule *reference, Graph *KeySetList, atom ***&ListOfLocalAtoms, Graph **&FragmentList, int &FragmentCounter, bool FreeList = false);
    380   bool FillListOfLocalAtoms(ofstream *out, atom ***&ListOfLocalAtoms, int &FragmentCounter, int GlobalAtomCount, bool &FreeList);
     381  bool FillListOfLocalAtoms(ofstream *out, atom ***&ListOfLocalAtoms, const int FragmentCounter, const int GlobalAtomCount, bool &FreeList);
    381382  void TranslateIndicesToGlobalIDs(ofstream *out, Graph **FragmentList, int &FragmentCounter, int &TotalNumberOfKeySets, Graph &TotalGraph);
    382383  int Count() const;
Note: See TracChangeset for help on using the changeset viewer.