Ignore:
Timestamp:
Mar 2, 2011, 9:53:08 PM (14 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:
1a4d4fe
Parents:
2e4105
git-author:
Frederik Heber <heber@…> (02/25/11 20:52:50)
git-committer:
Frederik Heber <heber@…> (03/02/11 21:53:08)
Message:

Moved CyclicStructureAnalysis from molecule_graph.cpp into own functor in Graph/.

File:
1 edited

Legend:

Unmodified
Added
Removed

  • src/molecule_graph.cpp

    r2e4105 re73ad9a  
    4444#define MAXBONDS 8
    4545
    46 struct BFSAccounting
    47 {
    48   atom **PredecessorList;
    49   int *ShortestPathList;
    50   enum GraphEdge::Shading *ColorList;
    51   std::deque<atom *> *BFSStack;
    52   std::deque<atom *> *TouchedStack;
    53   int AtomCount;
    54   int BondOrder;
    55   atom *Root;
    56   bool BackStepping;
    57   int CurrentGraphNr;
    58   int ComponentNr;
    59 };
    6046
    6147/** Accounting data for Depth First Search.
    … …  
    200186
    201187
    202 /** Sets atom::GraphNr and atom::LowpointNr to BFSAccounting::CurrentGraphNr.
     188/** Sets atom::GraphNr and atom::LowpointNr to DFSAccounting::CurrentGraphNr.
    203189 * \param *Walker current node
    204190 * \param &BFS structure with accounting data for BFS
    … …  
    518504}
    519505;
    520 
    521 /** Initialise each vertex as white with no predecessor, empty queue, color Root lightgray.
    522  * \param &BFS accounting structure
    523  * \param AtomCount number of entries in the array to allocate
    524  */
    525 void InitializeBFSAccounting(struct BFSAccounting &BFS, int AtomCount)
    526 {
    527   BFS.AtomCount = AtomCount;
    528   BFS.PredecessorList = new atom*[AtomCount];
    529   BFS.ShortestPathList = new int[AtomCount];
    530   BFS.ColorList = new enum GraphEdge::Shading[AtomCount];
    531   BFS.BFSStack = new std::deque<atom *> (AtomCount);
    532   BFS.TouchedStack = new std::deque<atom *> (AtomCount);
    533 
    534   for (int i = AtomCount; i--;) {
    535     BFS.ShortestPathList[i] = -1;
    536     BFS.PredecessorList[i] = 0;
    537     BFS.ColorList[i] = GraphEdge::white;
    538   }
    539 };
    540 
    541 /** Free's accounting structure.
    542  * \param &BFS accounting structure
    543  */
    544 void FinalizeBFSAccounting(struct BFSAccounting &BFS)
    545 {
    546   delete[](BFS.PredecessorList);
    547   delete[](BFS.ShortestPathList);
    548   delete[](BFS.ColorList);
    549   delete (BFS.BFSStack);
    550   delete (BFS.TouchedStack);
    551   BFS.AtomCount = 0;
    552 };
    553 
    554 /** Clean the accounting structure.
    555  * \param &BFS accounting structure
    556  */
    557 void CleanBFSAccounting(struct BFSAccounting &BFS)
    558 {
    559   atom *Walker = NULL;
    560   while (!BFS.TouchedStack->empty()) {
    561     Walker = BFS.TouchedStack->front();
    562     BFS.TouchedStack->pop_front();
    563     BFS.PredecessorList[Walker->getNr()] = NULL;
    564     BFS.ShortestPathList[Walker->getNr()] = -1;
    565     BFS.ColorList[Walker->getNr()] = GraphEdge::white;
    566   }
    567 };
    568 
    569 /** Resets shortest path list and BFSStack.
    570  * \param *&Walker current node, pushed onto BFSAccounting::BFSStack and BFSAccounting::TouchedStack
    571  * \param &BFS accounting structure
    572  */
    573 void ResetBFSAccounting(atom *&Walker, struct BFSAccounting &BFS)
    574 {
    575   BFS.ShortestPathList[Walker->getNr()] = 0;
    576   BFS.BFSStack->clear(); // start with empty BFS stack
    577   BFS.BFSStack->push_front(Walker);
    578   BFS.TouchedStack->push_front(Walker);
    579 };
    580 
    581 /** Performs a BFS from \a *Root, trying to find the same node and hence a cycle.
    582  * \param *&BackEdge the edge from root that we don't want to move along
    583  * \param &BFS accounting structure
    584  */
    585 void CyclicStructureAnalysis_CyclicBFSFromRootToRoot(bond *&BackEdge, struct BFSAccounting &BFS)
    586 {
    587   atom *Walker = NULL;
    588   atom *OtherAtom = NULL;
    589   do { // look for Root
    590     ASSERT(!BFS.BFSStack->empty(), "CyclicStructureAnalysis_CyclicBFSFromRootToRoot() - BFS.BFSStack is empty!");
    591     Walker = BFS.BFSStack->front();
    592     BFS.BFSStack->pop_front();
    593     DoLog(2) && (Log() << Verbose(2) << "Current Walker is " << *Walker << ", we look for SP to Root " << *BFS.Root << "." << endl);
    594     const BondList& ListOfBonds = Walker->getListOfBonds();
    595     for (BondList::const_iterator Runner = ListOfBonds.begin();
    596         Runner != ListOfBonds.end();
    597         ++Runner) {
    598       if ((*Runner) != BackEdge) { // only walk along DFS spanning tree (otherwise we always find SP of one being backedge Binder)
    599         OtherAtom = (*Runner)->GetOtherAtom(Walker);
    600 #ifdef ADDHYDROGEN
    601         if (OtherAtom->getType()->getAtomicNumber() != 1) {
    602 #endif
    603         DoLog(2) && (Log() << Verbose(2) << "Current OtherAtom is: " << OtherAtom->getName() << " for bond " << *(*Runner) << "." << endl);
    604         if (BFS.ColorList[OtherAtom->getNr()] == GraphEdge::white) {
    605           BFS.TouchedStack->push_front(OtherAtom);
    606           BFS.ColorList[OtherAtom->getNr()] = GraphEdge::lightgray;
    607           BFS.PredecessorList[OtherAtom->getNr()] = Walker; // Walker is the predecessor
    608           BFS.ShortestPathList[OtherAtom->getNr()] = BFS.ShortestPathList[Walker->getNr()] + 1;
    609           DoLog(2) && (Log() << Verbose(2) << "Coloring OtherAtom " << OtherAtom->getName() << " lightgray, its predecessor is " << Walker->getName() << " and its Shortest Path is " << BFS.ShortestPathList[OtherAtom->getNr()] << " egde(s) long." << endl);
    610           //if (BFS.ShortestPathList[OtherAtom->getNr()] < MinimumRingSize[Walker->GetTrueFather()->nr]) { // Check for maximum distance
    611           DoLog(3) && (Log() << Verbose(3) << "Putting OtherAtom into queue." << endl);
    612           BFS.BFSStack->push_front(OtherAtom);
    613           //}
    614         } else {
    615           DoLog(3) && (Log() << Verbose(3) << "Not Adding, has already been visited." << endl);
    616         }
    617         if (OtherAtom == BFS.Root)
    618           break;
    619 #ifdef ADDHYDROGEN
    620       } else {
    621         DoLog(2) && (Log() << Verbose(2) << "Skipping hydrogen atom " << *OtherAtom << "." << endl);
    622         BFS.ColorList[OtherAtom->getNr()] = GraphEdge::black;
    623       }
    624 #endif
    625       } else {
    626         DoLog(2) && (Log() << Verbose(2) << "Bond " << *(*Runner) << " not Visiting, is the back edge." << endl);
    627       }
    628     }
    629     BFS.ColorList[Walker->getNr()] = GraphEdge::black;
    630     DoLog(1) && (Log() << Verbose(1) << "Coloring Walker " << Walker->getName() << " black." << endl);
    631     if (OtherAtom == BFS.Root) { // if we have found the root, check whether this cycle wasn't already found beforehand
    632       // step through predecessor list
    633       while (OtherAtom != BackEdge->rightatom) {
    634         if (!OtherAtom->GetTrueFather()->IsCyclic) // if one bond in the loop is not marked as cyclic, we haven't found this cycle yet
    635           break;
    636         else
    637           OtherAtom = BFS.PredecessorList[OtherAtom->getNr()];
    638       }
    639       if (OtherAtom == BackEdge->rightatom) { // if each atom in found cycle is cyclic, loop's been found before already
    640         DoLog(3) && (Log() << Verbose(3) << "This cycle was already found before, skipping and removing seeker from search." << endl);
    641         do {
    642           ASSERT(!BFS.TouchedStack->empty(), "CyclicStructureAnalysis_CyclicBFSFromRootToRoot() - BFS.TouchedStack is empty!");
    643           OtherAtom = BFS.TouchedStack->front();
    644           BFS.TouchedStack->pop_front();
    645           if (BFS.PredecessorList[OtherAtom->getNr()] == Walker) {
    646             DoLog(4) && (Log() << Verbose(4) << "Removing " << *OtherAtom << " from lists and stacks." << endl);
    647             BFS.PredecessorList[OtherAtom->getNr()] = NULL;
    648             BFS.ShortestPathList[OtherAtom->getNr()] = -1;
    649             BFS.ColorList[OtherAtom->getNr()] = GraphEdge::white;
    650             // rats ... deque has no find()
    651             std::deque<atom *>::iterator iter = find(
    652                 BFS.BFSStack->begin(),
    653                 BFS.BFSStack->end(),
    654                 OtherAtom);
    655             ASSERT(iter != BFS.BFSStack->end(),
    656                 "CyclicStructureAnalysis_CyclicBFSFromRootToRoot() - can't find "+toString(*OtherAtom)+" on stack!");
    657             BFS.BFSStack->erase(iter);
    658           }
    659         } while ((!BFS.TouchedStack->empty()) && (BFS.PredecessorList[OtherAtom->getNr()] == NULL));
    660         BFS.TouchedStack->push_front(OtherAtom); // last was wrongly popped
    661         OtherAtom = BackEdge->rightatom; // set to not Root
    662       } else
    663         OtherAtom = BFS.Root;
    664     }
    665   } while ((!BFS.BFSStack->empty()) && (OtherAtom != BFS.Root) && (OtherAtom != NULL)); // || (ShortestPathList[OtherAtom->getNr()] < MinimumRingSize[Walker->GetTrueFather()->getNr()])));
    666 };
    667 
    668 /** Climb back the BFSAccounting::PredecessorList and find cycle members.
    669  * \param *&OtherAtom
    670  * \param *&BackEdge denotes the edge we did not want to travel along when doing CyclicBFSFromRootToRoot()
    671  * \param &BFS accounting structure
    672  * \param *&MinimumRingSize minimum distance from this node possible without encountering oneself, set on return for each atom
    673  * \param &MinRingSize global minimum distance from one node without encountering oneself, set on return
    674  */
    675 void CyclicStructureAnalysis_RetrieveCycleMembers(atom *&OtherAtom, bond *&BackEdge, struct BFSAccounting &BFS, int *&MinimumRingSize, int &MinRingSize)
    676 {
    677   atom *Walker = NULL;
    678   int NumCycles = 0;
    679   int RingSize = -1;
    680 
    681   if (OtherAtom == BFS.Root) {
    682     // now climb back the predecessor list and thus find the cycle members
    683     NumCycles++;
    684     RingSize = 1;
    685     BFS.Root->GetTrueFather()->IsCyclic = true;
    686     DoLog(1) && (Log() << Verbose(1) << "Found ring contains: ");
    687     Walker = BFS.Root;
    688     while (Walker != BackEdge->rightatom) {
    689       DoLog(0) && (Log() << Verbose(0) << Walker->getName() << " <-> ");
    690       Walker = BFS.PredecessorList[Walker->getNr()];
    691       Walker->GetTrueFather()->IsCyclic = true;
    692       RingSize++;
    693     }
    694     DoLog(0) && (Log() << Verbose(0) << Walker->getName() << "  with a length of " << RingSize << "." << endl << endl);
    695     // walk through all and set MinimumRingSize
    696     Walker = BFS.Root;
    697     MinimumRingSize[Walker->GetTrueFather()->getNr()] = RingSize;
    698     while (Walker != BackEdge->rightatom) {
    699       Walker = BFS.PredecessorList[Walker->getNr()];
    700       if (RingSize < MinimumRingSize[Walker->GetTrueFather()->getNr()])
    701         MinimumRingSize[Walker->GetTrueFather()->getNr()] = RingSize;
    702     }
    703     if ((RingSize < MinRingSize) || (MinRingSize == -1))
    704       MinRingSize = RingSize;
    705   } else {
    706     DoLog(1) && (Log() << Verbose(1) << "No ring containing " << *BFS.Root << " with length equal to or smaller than " << MinimumRingSize[BFS.Root->GetTrueFather()->getNr()] << " found." << endl);
    707   }
    708 };
    709 
    710 /** From a given node performs a BFS to touch the next cycle, for whose nodes \a *&MinimumRingSize is set and set it accordingly.
    711  * \param *&Root node to look for closest cycle from, i.e. \a *&MinimumRingSize is set for this node
    712  * \param *&MinimumRingSize minimum distance from this node possible without encountering oneself, set on return for each atom
    713  * \param AtomCount number of nodes in graph
    714  */
    715 void CyclicStructureAnalysis_BFSToNextCycle(atom *&Root, atom *&Walker, int *&MinimumRingSize, int AtomCount)
    716 {
    717   struct BFSAccounting BFS;
    718   atom *OtherAtom = Walker;
    719 
    720   InitializeBFSAccounting(BFS, AtomCount);
    721 
    722   ResetBFSAccounting(Walker, BFS);
    723   while (OtherAtom != NULL) { // look for Root
    724     ASSERT(!BFS.BFSStack->empty(), "CyclicStructureAnalysis_BFSToNextCycle() - BFS.BFSStack is empty!");
    725     Walker = BFS.BFSStack->front();
    726     BFS.BFSStack->pop_front();
    727     //Log() << Verbose(2) << "Current Walker is " << *Walker << ", we look for SP to Root " << *Root << "." << endl;
    728     const BondList& ListOfBonds = Walker->getListOfBonds();
    729     for (BondList::const_iterator Runner = ListOfBonds.begin();
    730         Runner != ListOfBonds.end();
    731         ++Runner) {
    732       // "removed (*Runner) != BackEdge) || " from next if, is u
    733       if ((ListOfBonds.size() == 1)) { // only walk along DFS spanning tree (otherwise we always find SP of 1 being backedge Binder), but terminal hydrogens may be connected via backedge, hence extra check
    734         OtherAtom = (*Runner)->GetOtherAtom(Walker);
    735         //Log() << Verbose(2) << "Current OtherAtom is: " << OtherAtom->Name << " for bond " << *Binder << "." << endl;
    736         if (BFS.ColorList[OtherAtom->getNr()] == GraphEdge::white) {
    737           BFS.TouchedStack->push_front(OtherAtom);
    738           BFS.ColorList[OtherAtom->getNr()] = GraphEdge::lightgray;
    739           BFS.PredecessorList[OtherAtom->getNr()] = Walker; // Walker is the predecessor
    740           BFS.ShortestPathList[OtherAtom->getNr()] = BFS.ShortestPathList[Walker->getNr()] + 1;
    741           //Log() << Verbose(2) << "Coloring OtherAtom " << OtherAtom->Name << " lightgray, its predecessor is " << Walker->Name << " and its Shortest Path is " << ShortestPathList[OtherAtom->getNr()] << " egde(s) long." << endl;
    742           if (OtherAtom->GetTrueFather()->IsCyclic) { // if the other atom is connected to a ring
    743             MinimumRingSize[Root->GetTrueFather()->getNr()] = BFS.ShortestPathList[OtherAtom->getNr()] + MinimumRingSize[OtherAtom->GetTrueFather()->getNr()];
    744             OtherAtom = NULL; //break;
    745             break;
    746           } else
    747             BFS.BFSStack->push_front(OtherAtom);
    748         } else {
    749           //Log() << Verbose(3) << "Not Adding, has already been visited." << endl;
    750         }
    751       } else {
    752         //Log() << Verbose(3) << "Not Visiting, is a back edge." << endl;
    753       }
    754     }
    755     BFS.ColorList[Walker->getNr()] = GraphEdge::black;
    756     //Log() << Verbose(1) << "Coloring Walker " << Walker->Name << " black." << endl;
    757   }
    758   //CleanAccountingLists(TouchedStack, PredecessorList, ShortestPathList, ColorList);
    759 
    760   FinalizeBFSAccounting(BFS);
    761 }
    762 ;
    763 
    764 /** All nodes that are not in cycles get assigned a \a *&MinimumRingSizeby BFS to next cycle.
    765  * \param *&MinimumRingSize array with minimum distance without encountering onself for each atom
    766  * \param &MinRingSize global minium distance
    767  * \param &NumCyles number of cycles in graph
    768  * \param *mol molecule with atoms
    769  */
    770 void CyclicStructureAnalysis_AssignRingSizetoNonCycleMembers(int *&MinimumRingSize, int &MinRingSize, int &NumCycles, const molecule * const mol)
    771 {
    772   atom *Root = NULL;
    773   atom *Walker = NULL;
    774   if (MinRingSize != -1) { // if rings are present
    775     // go over all atoms
    776     for (molecule::const_iterator iter = mol->begin(); iter != mol->end(); ++iter) {
    777       Root = *iter;
    778 
    779       if (MinimumRingSize[Root->GetTrueFather()->getNr()] == mol->getAtomCount()) { // check whether MinimumRingSize is set, if not BFS to next where it is
    780         Walker = Root;
    781 
    782         //Log() << Verbose(1) << "---------------------------------------------------------------------------------------------------------" << endl;
    783         CyclicStructureAnalysis_BFSToNextCycle(Root, Walker, MinimumRingSize, mol->getAtomCount());
    784 
    785       }
    786       DoLog(1) && (Log() << Verbose(1) << "Minimum ring size of " << *Root << " is " << MinimumRingSize[Root->GetTrueFather()->getNr()] << "." << endl);
    787     }
    788     DoLog(1) && (Log() << Verbose(1) << "Minimum ring size is " << MinRingSize << ", over " << NumCycles << " cycles total." << endl);
    789   } else
    790     DoLog(1) && (Log() << Verbose(1) << "No rings were detected in the molecular structure." << endl);
    791 }
    792 ;
    793 
    794 /** Analyses the cycles found and returns minimum of all cycle lengths.
    795  * We begin with a list of Back edges found during DepthFirstSearchAnalysis(). We go through this list - one end is the Root,
    796  * the other our initial Walker - and do a Breadth First Search for the Root. We mark down each Predecessor and as soon as
    797  * we have found the Root via BFS, we may climb back the closed cycle via the Predecessors. Thereby we mark atoms and bonds
    798  * as cyclic and print out the cycles.
    799  * \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!
    800  * \param *&MinimumRingSize contains smallest ring size in molecular structure on return or -1 if no rings were found, if set is maximum search distance
    801  * \todo BFS from the not-same-LP to find back to starting point of tributary cycle over more than one bond
    802  */
    803 void molecule::CyclicStructureAnalysis(
    804     std::deque<bond *> * BackEdgeStack,
    805     int *&MinimumRingSize
    806     ) const
    807 {
    808   struct BFSAccounting BFS;
    809   atom *Walker = NULL;
    810   atom *OtherAtom = NULL;
    811   bond *BackEdge = NULL;
    812   int NumCycles = 0;
    813   int MinRingSize = -1;
    814 
    815   InitializeBFSAccounting(BFS, getAtomCount());
    816 
    817   //Log() << Verbose(1) << "Back edge list - ";
    818   //BackEdgeStack->Output(out);
    819 
    820   DoLog(1) && (Log() << Verbose(1) << "Analysing cycles ... " << endl);
    821   NumCycles = 0;
    822   while (!BackEdgeStack->empty()) {
    823     BackEdge = BackEdgeStack->front();
    824     BackEdgeStack->pop_front();
    825     // this is the target
    826     BFS.Root = BackEdge->leftatom;
    827     // this is the source point
    828     Walker = BackEdge->rightatom;
    829 
    830     ResetBFSAccounting(Walker, BFS);
    831 
    832     DoLog(1) && (Log() << Verbose(1) << "---------------------------------------------------------------------------------------------------------" << endl);
    833     OtherAtom = NULL;
    834     CyclicStructureAnalysis_CyclicBFSFromRootToRoot(BackEdge, BFS);
    835 
    836     CyclicStructureAnalysis_RetrieveCycleMembers(OtherAtom, BackEdge, BFS, MinimumRingSize, MinRingSize);
    837 
    838     CleanBFSAccounting(BFS);
    839   }
    840   FinalizeBFSAccounting(BFS);
    841 
    842   CyclicStructureAnalysis_AssignRingSizetoNonCycleMembers(MinimumRingSize, MinRingSize, NumCycles, this);
    843 };
    844506
    845507/** Sets the next component number.
Note: See TracChangeset for help on using the changeset viewer.