Changeset 1907a7 for src/molecules.cpp

Timestamp:

Apr 2, 2009, 4:12:54 PM (16 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:

ca3ccc

Parents:

d8b94a

Message:

Basic implementation of Multiple molecules.

builder.cpp:

• MoleculeListClass *molecules instead of molecule *mol
• in menu replaced list of elements and atoms by list of molecules
• new option: SetActiveMolecule
• new submenus/functions that all work on MoleculeListClass:
  • EditMolecule()
  • ManipulateMolecules() (global operations on all atoms in a molecule)
  • MergeMolecules()
  • ManipulateAtoms() (local operations on single atoms)
• SaveConfig(): SimpleMultiAdd() of all molecules into a single new one that is stored with the present configuration.

molecules.hpp:

• new STL list MoleculeList
• new variables in MoleculeListClass: list of molecules as STL list
• new variables in molecule: ActiveFlag, name of molecule, Center of molecule (used on merging)

moleculelist.cpp:

• replaced listofmolecules array by STL list everywhere (only smaller changes necessary)
• new merging function: SimpleMerge, SimpleAdd, SimpleMultiMerge, SimpleMultiAdd, (EmbedMerge, ScatterMerge ... both not finished). Add does not while merge does delete the src molecules.
• new function: Enumerate(). Output of all molecules with number of atoms and chemical formula
• new function: NumberOfActiveMolecules(). Counts the number of molecules in the list with ActiveFlag set.
• new function: insert(). Inserts a molecule into the list with a unique index

molecules.cpp:

• new function: SetNameFromFilename. Takes basename of a filename and sets name accordingly.
• new function: UnlinkAtom. Only removes atom from list, does not delete it from memory.

atom.cpp:

• Output() also accepts specific comment instead of "# molecule nr ..."

File:

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

src/molecules.cpp

@@ -62 +62 @@
         cell_size[0] = cell_size[2] = cell_size[5]= 20.;
         cell_size[1] = cell_size[3] = cell_size[4]= 0.;
+        strcpy(name,"none");
 };
 
@@ -596 +597 @@
         cerr << Verbose(1) << "molecule::RemoveBond: Function not implemented yet." << endl;
         return false;
+};
+
+/** Set molecule::name from the basename without suffix in the given \a *filename.
+ * \param *filename filename
+ */
+void molecule::SetNameFromFilename(char *filename)
+{
+  int length = 0;
+  char *molname = strrchr(filename, '/')+sizeof(char);  // search for filename without dirs
+  char *endname = strrchr(filename, '.');
+  if ((endname == NULL) || (endname < molname))
+    length = strlen(molname);
+  else
+    length = strlen(molname) - strlen(endname);
+  strncpy(name, molname, length);
 };
 
@@ -1187 +1203 @@
 };
 
-/** Removes atom from molecule list.
+/** Removes atom from molecule list and deletes it.
  * \param *pointer atom to be removed
  * \return true - succeeded, false - atom not found in list
@@ -1201 +1217 @@
         Trajectories.erase(pointer);
         return remove(pointer, start, end);
+};
+
+/** Removes atom from molecule list, but does not delete it.
+ * \param *pointer atom to be removed
+ * \return true - succeeded, false - atom not found in list
+ */
+bool molecule::UnlinkAtom(atom *pointer)
+{
+  if (pointer == NULL)
+    return false;
+  if (ElementsInMolecule[pointer->type->Z] != 0)  // this would indicate an error
+    ElementsInMolecule[pointer->type->Z]--; // decrease number of atom of this element
+  else
+    cerr << "ERROR: Atom " << pointer->Name << " is of element " << pointer->type->Z << " but the entry in the table of the molecule is 0!" << endl;
+  if (ElementsInMolecule[pointer->type->Z] == 0)  // was last atom of this element?
+    ElementCount--;
+  Trajectories.erase(pointer);
+  unlink(pointer);
+  return true;
 };
 
@@ -1727 +1762 @@
                                                 };
                                                 AddBond(Walker, OtherWalker); //Add the bond between the two atoms with respective indices.
-                
+
                                         }
 
@@ -2529 +2564 @@
                                         cerr << "KeySet file must be corrupt as there are two equal key sets therein!" << endl;
                                 }
-                                //FragmentList->ListOfMolecules[NumberOfFragments++] = StoreFragmentFromKeySet(out, CurrentSet, IsAngstroem);
                         }
                 }
@@ -3088 +3122 @@
         //if (FragmentationToDo) {              // we should always store the fragments again as coordination might have changed slightly without changing bond structure
                 // allocate memory for the pointer array and transmorph graphs into full molecular fragments
-                BondFragments = new MoleculeListClass(TotalGraph.size(), AtomCount);
+                BondFragments = new MoleculeListClass();
                 int k=0;
                 for(Graph::iterator runner = TotalGraph.begin(); runner != TotalGraph.end(); runner++) {
                         KeySet test = (*runner).first;
                         *out << "Fragment No." << (*runner).second.first << " with TEFactor " << (*runner).second.second << "." << endl;
-                        BondFragments->ListOfMolecules[k] = StoreFragmentFromKeySet(out, test, configuration);
+                        BondFragments->insert(StoreFragmentFromKeySet(out, test, configuration));
                         k++;
                 }
-                *out << k << "/" << BondFragments->NumberOfMolecules << " fragments generated from the keysets." << endl;
+                *out << k << "/" << BondFragments->ListOfMolecules.size() << " fragments generated from the keysets." << endl;
 
                 // ===== 9. Save fragments' configuration and keyset files et al to disk ===
-                if (BondFragments->NumberOfMolecules != 0) {
+                if (BondFragments->ListOfMolecules.size() != 0) {
                         // create the SortIndex from BFS labels to order in the config file
                         CreateMappingLabelsToConfigSequence(out, SortIndex);
 
-                        *out << Verbose(1) << "Writing " << BondFragments->NumberOfMolecules << " possible bond fragmentation configs" << endl;
+                        *out << Verbose(1) << "Writing " << BondFragments->ListOfMolecules.size() << " possible bond fragmentation configs" << endl;
                         if (BondFragments->OutputConfigForListOfFragments(out, configuration, SortIndex))
                                 *out << Verbose(1) << "All configs written." << endl;
@@ -3629 +3663 @@
 };
 
-/** Creates \a MoleculeListClass of all unique fragments of the \a molecule containing \a Order atoms or vertices.
- * The picture to have in mind is that of a DFS "snake" of a certain length \a Order, i.e. as in the infamous
- * computer game, that winds through the connected graph representing the molecule. Color (white,
- * lightgray, darkgray, black) indicates whether a vertex has been discovered so far or not. Labels will help in
- * creating only unique fragments and not additional ones with vertices simply in different sequence.
- * The Predecessor is always the one that came before in discovering, needed on backstepping. And
- * finally, the ShortestPath is needed for removing vertices from the snake stack during the back-
- * stepping.
- * \param *out output stream for debugging
- * \param Order number of atoms in each fragment
- * \param *configuration configuration for writing config files for each fragment
- * \return List of all unique fragments with \a Order atoms
- */
-/*
-MoleculeListClass * molecule::CreateListOfUniqueFragmentsOfOrder(ofstream *out, int Order, config *configuration)
-{
-        atom **PredecessorList = (atom **) Malloc(sizeof(atom *)*AtomCount, "molecule::CreateListOfUniqueFragmentsOfOrder: **PredecessorList");
-        int *ShortestPathList = (int *) Malloc(sizeof(int)*AtomCount, "molecule::CreateListOfUniqueFragmentsOfOrder: *ShortestPathList");
-        int *Labels = (int *) Malloc(sizeof(int)*AtomCount, "molecule::CreateListOfUniqueFragmentsOfOrder: *Labels");
-        enum Shading *ColorVertexList = (enum Shading *) Malloc(sizeof(enum Shading)*AtomCount, "molecule::CreateListOfUniqueFragmentsOfOrder: *ColorList");
-        enum Shading *ColorEdgeList = (enum Shading *) Malloc(sizeof(enum Shading)*BondCount, "molecule::CreateListOfUniqueFragmentsOfOrder: *ColorBondList");
-        StackClass<atom *> *RootStack = new StackClass<atom *>(AtomCount);
-        StackClass<atom *> *TouchedStack = new StackClass<atom *>((int)pow(4,Order)+2); // number of atoms reached from one with maximal 4 bonds plus Root itself
-        StackClass<atom *> *SnakeStack = new StackClass<atom *>(Order+1); // equal to Order is not possible, as then the StackClass<atom *> cannot discern between full and empty stack!
-        MoleculeLeafClass *Leaflet = NULL, *TempLeaf = NULL;
-        MoleculeListClass *FragmentList = NULL;
-        atom *Walker = NULL, *OtherAtom = NULL, *Root = NULL, *Removal = NULL;
-        bond *Binder = NULL;
-        int RunningIndex = 0, FragmentCounter = 0;
-
-        *out << Verbose(1) << "Begin of CreateListOfUniqueFragmentsOfOrder." << endl;
-
-        // reset parent list
-        *out << Verbose(3) << "Resetting labels, parent, predecessor, color and shortest path lists." << endl;
-        for (int i=0;i<AtomCount;i++) { // reset all atom labels
-                // initialise each vertex as white with no predecessor, empty queue, color lightgray, not labelled, no sons
-                Labels[i] = -1;
-                SonList[i] = NULL;
-                PredecessorList[i] = NULL;
-                ColorVertexList[i] = white;
-                ShortestPathList[i] = -1;
-        }
-        for (int i=0;i<BondCount;i++)
-                ColorEdgeList[i] = white;
-        RootStack->ClearStack();        // clearstack and push first atom if exists
-        TouchedStack->ClearStack();
-        Walker = start->next;
-        while ((Walker != end)
-#ifdef ADDHYDROGEN
-         && (Walker->type->Z == 1)
-#endif
-                                                                                                ) { // search for first non-hydrogen atom
-                *out << Verbose(4) << "Current Root candidate is " << Walker->Name << "." << endl;
-                Walker = Walker->next;
-        }
-        if (Walker != end)
-                RootStack->Push(Walker);
-        else
-                *out << Verbose(0) << "ERROR: Could not find an appropriate Root atom!" << endl;
-        *out << Verbose(3) << "Root " << Walker->Name << " is on AtomStack, beginning loop through all vertices ..." << endl;
-
-        ///// OUTER LOOP ////////////
-        while (!RootStack->IsEmpty()) {
-                // get new root vertex from atom stack
-                Root = RootStack->PopFirst();
-                ShortestPathList[Root->nr] = 0;
-                if (Labels[Root->nr] == -1)
-                        Labels[Root->nr] = RunningIndex++; // prevent it from getting again on AtomStack
-                PredecessorList[Root->nr] = Root;
-                TouchedStack->Push(Root);
-                *out << Verbose(0) << "Root for this loop is: " << Root->Name << ".\n";
-
-                // clear snake stack
-                SnakeStack->ClearStack();
-                //SnakeStack->TestImplementation(out, start->next);
-
-                ///// INNER LOOP ////////////
-                // Problems:
-                // - what about cyclic bonds?
-                Walker = Root;
-                do {
-                        *out << Verbose(1) << "Current Walker is: " << Walker->Name;
-                        // initial setting of the new Walker: label, color, shortest path and put on stacks
-                        if (Labels[Walker->nr] == -1)   {       // give atom a unique, monotonely increasing number
-                                Labels[Walker->nr] = RunningIndex++;
-                                RootStack->Push(Walker);
-                        }
-                        *out << ", has label " << Labels[Walker->nr];
-                        if ((ColorVertexList[Walker->nr] == white) || ((Binder != NULL) && (ColorEdgeList[Binder->nr] == white))) {     // color it if newly discovered and push on stacks (and if within reach!)
-                                if ((Binder != NULL) && (ColorEdgeList[Binder->nr] == white)) {
-                                        // Binder ought to be set still from last neighbour search
-                                        *out << ", coloring bond " << *Binder << " black";
-                                        ColorEdgeList[Binder->nr] = black; // mark this bond as used
-                                }
-                                if (ShortestPathList[Walker->nr] == -1) {
-                                        ShortestPathList[Walker->nr] = ShortestPathList[PredecessorList[Walker->nr]->nr]+1;
-                                        TouchedStack->Push(Walker); // mark every atom for lists cleanup later, whose shortest path has been changed
-                                }
-                                if ((ShortestPathList[Walker->nr] < Order) && (ColorVertexList[Walker->nr] != darkgray)) {      // if not already on snake stack
-                                        SnakeStack->Push(Walker);
-                                        ColorVertexList[Walker->nr] = darkgray; // mark as dark gray of on snake stack
-                                }
-                        }
-                        *out << ", SP of " << ShortestPathList[Walker->nr]      << " and its color is " << GetColor(ColorVertexList[Walker->nr]) << "." << endl;
-
-                        // then check the stack for a newly stumbled upon fragment
-                        if (SnakeStack->ItemCount() == Order) { // is stack full?
-                                // store the fragment if it is one and get a removal candidate
-                                Removal = StoreFragmentFromStack(out, Root, Walker, Leaflet, SnakeStack, ShortestPathList, SonList, Labels, &FragmentCounter, configuration);
-                                // remove the candidate if one was found
-                                if (Removal != NULL) {
-                                        *out << Verbose(2) << "Removing item " << Removal->Name << " with SP of " << ShortestPathList[Removal->nr] << " from snake stack." << endl;
-                                        SnakeStack->RemoveItem(Removal);
-                                        ColorVertexList[Removal->nr] = lightgray; // return back to not on snake stack but explored marking
-                                        if (Walker == Removal) { // if the current atom is to be removed, we also have to take a step back
-                                                Walker = PredecessorList[Removal->nr];
-                                                *out << Verbose(2) << "Stepping back to " << Walker->Name << "." << endl;
-                                        }
-                                }
-                        } else
-                                Removal = NULL;
-
-                        // finally, look for a white neighbour as the next Walker
-                        Binder = NULL;
-                        if ((Removal == NULL) || (Walker != PredecessorList[Removal->nr])) {    // don't look, if a new walker has been set above
-                                *out << Verbose(2) << "Snake has currently " << SnakeStack->ItemCount() << " item(s)." << endl;
-                                OtherAtom = NULL; // this is actually not needed, every atom has at least one neighbour
-                                if (ShortestPathList[Walker->nr] < Order) {
-                                        for(int i=0;i<NumberOfBondsPerAtom[Walker->nr];i++) {
-                                                Binder = ListOfBondsPerAtom[Walker->nr][i];
-                                                *out << Verbose(2) << "Current bond is " << *Binder << ": ";
-                                                OtherAtom = Binder->GetOtherAtom(Walker);
-                                                if ((Labels[OtherAtom->nr] != -1) && (Labels[OtherAtom->nr] < Labels[Root->nr])) { // we don't step up to labels bigger than us
-                                                        *out << "Label " << Labels[OtherAtom->nr] << " is smaller than Root's " << Labels[Root->nr] << "." << endl;
-                                                        //ColorVertexList[OtherAtom->nr] = lightgray;           // mark as explored
-                                                } else { // otherwise check its colour and element
-                                                        if (
-#ifdef ADDHYDROGEN
-                                                        (OtherAtom->type->Z != 1) &&
-#endif
-                                                                                (ColorEdgeList[Binder->nr] == white)) { // skip hydrogen, look for unexplored vertices
-                                                                *out << "Moving along " << GetColor(ColorEdgeList[Binder->nr]) << " bond " << Binder << " to " << ((ColorVertexList[OtherAtom->nr] == white) ? "unexplored" : "explored") << " item: " << OtherAtom->Name << "." << endl;
-                                                                // i find it currently rather sensible to always set the predecessor in order to find one's way back
-                                                                //if (PredecessorList[OtherAtom->nr] == NULL) {
-                                                                PredecessorList[OtherAtom->nr] = Walker;
-                                                                *out << Verbose(3) << "Setting Predecessor of " << OtherAtom->Name << " to " << PredecessorList[OtherAtom->nr]->Name << "." << endl;
-                                                                //} else {
-                                                                //      *out << Verbose(3) << "Predecessor of " << OtherAtom->Name << " is " << PredecessorList[OtherAtom->nr]->Name << "." << endl;
-                                                                //}
-                                                                Walker = OtherAtom;
-                                                                break;
-                                                        } else {
-                                                                if (OtherAtom->type->Z == 1)
-                                                                        *out << "Links to a hydrogen atom." << endl;
-                                                                else
-                                                                        *out << "Bond has not white but " << GetColor(ColorEdgeList[Binder->nr]) << " color." << endl;
-                                                        }
-                                                }
-                                        }
-                                } else {        // means we have stepped beyond the horizon: Return!
-                                        Walker = PredecessorList[Walker->nr];
-                                        OtherAtom = Walker;
-                                        *out << Verbose(3) << "We have gone too far, stepping back to " << Walker->Name << "." << endl;
-                                }
-                                if (Walker != OtherAtom) {      // if no white neighbours anymore, color it black
-                                        *out << Verbose(2) << "Coloring " << Walker->Name << " black." << endl;
-                                        ColorVertexList[Walker->nr] = black;
-                                        Walker = PredecessorList[Walker->nr];
-                                }
-                        }
-                } while ((Walker != Root) || (ColorVertexList[Root->nr] != black));
-                *out << Verbose(2) << "Inner Looping is finished." << endl;
-
-                // if we reset all AtomCount atoms, we have again technically O(N^2) ...
-                *out << Verbose(2) << "Resetting lists." << endl;
-                Walker = NULL;
-                Binder = NULL;
-                while (!TouchedStack->IsEmpty()) {
-                        Walker = TouchedStack->PopLast();
-                        *out << Verbose(3) << "Re-initialising entries of " << *Walker << "." << endl;
-                        for(int i=0;i<NumberOfBondsPerAtom[Walker->nr];i++)
-                                ColorEdgeList[ListOfBondsPerAtom[Walker->nr][i]->nr] = white;
-                        PredecessorList[Walker->nr] = NULL;
-                        ColorVertexList[Walker->nr] = white;
-                        ShortestPathList[Walker->nr] = -1;
-                }
-        }
-        *out << Verbose(1) << "Outer Looping over all vertices is done." << endl;
-
-        // copy together
-        *out << Verbose(1) << "Copying all fragments into MoleculeList structure." << endl;
-        FragmentList = new MoleculeListClass(FragmentCounter, AtomCount);
-        RunningIndex = 0;
-        while ((Leaflet != NULL) && (RunningIndex < FragmentCounter))   {
-                FragmentList->ListOfMolecules[RunningIndex++] = Leaflet->Leaf;
-                Leaflet->Leaf = NULL; // prevent molecule from being removed
-                TempLeaf = Leaflet;
-                Leaflet = Leaflet->previous;
-                delete(TempLeaf);
-        };
-
-        // free memory and exit
-        Free((void **)&PredecessorList, "molecule::CreateListOfUniqueFragmentsOfOrder: **PredecessorList");
-        Free((void **)&ShortestPathList, "molecule::CreateListOfUniqueFragmentsOfOrder: *ShortestPathList");
-        Free((void **)&Labels, "molecule::CreateListOfUniqueFragmentsOfOrder: *Labels");
-        Free((void **)&ColorVertexList, "molecule::CreateListOfUniqueFragmentsOfOrder: *ColorList");
-        delete(RootStack);
-        delete(TouchedStack);
-        delete(SnakeStack);
-
-        *out << Verbose(1) << "End of CreateListOfUniqueFragmentsOfOrder." << endl;
-        return FragmentList;
-};
-*/
-
 /** Structure containing all values in power set combination generation.
  */