Changes in / [ce4487:04488a]

Files:

• src/Actions/ActionRegistry.cpp (modified) (6 diffs)
• src/Actions/ActionRegistry.hpp (modified) (3 diffs)
• src/Actions/CmdAction/BondLengthTableAction.cpp (modified) (1 diff)
• src/Actions/FragmentationAction/DepthFirstSearchAction.cpp (modified) (1 diff)
• src/Actions/FragmentationAction/FragmentationAction.cpp (modified) (3 diffs)
• src/Actions/Makefile.am (modified) (2 diffs)
• src/Actions/MapOfActions.cpp (modified) (3 diffs)
• src/Actions/MoleculeAction/FillWithMoleculeAction.cpp (modified) (1 diff)
• src/Actions/MoleculeAction/LinearInterpolationofTrajectoriesAction.cpp (modified) (1 diff)
• src/Actions/MoleculeAction/SaveAdjacencyAction.cpp (modified) (1 diff)
• src/Actions/MoleculeAction/SaveBondsAction.cpp (modified) (1 diff)
• src/Actions/ParserAction/LoadXyzAction.cpp (modified) (4 diffs)
• src/Actions/WorldAction/InputAction.cpp (added)
• src/Actions/WorldAction/InputAction.hpp (added)
• src/Actions/WorldAction/OutputAction.cpp (added)
• src/Actions/WorldAction/OutputAction.hpp (added)
• src/CommandLineParser.cpp (modified) (1 diff)
• src/ConfigFileBuffer.cpp (added)
• src/ConfigFileBuffer.hpp (added)
• src/Descriptors/AtomDescriptor.cpp (modified) (1 diff)
• src/Descriptors/MoleculeDescriptor.cpp (modified) (1 diff)
• src/Helpers/MemDebug.cpp (modified) (1 diff)
• src/Legacy/oldmenu.cpp (modified) (2 diffs)
• src/Makefile.am (modified) (10 diffs)
• src/Parser/ChangeTracker.cpp (modified) (3 diffs)
• src/Parser/ChangeTracker.hpp (modified) (2 diffs)
• src/Parser/FormatParser.cpp (modified) (2 diffs)
• src/Parser/FormatParserStorage.cpp (added)
• src/Parser/FormatParserStorage.hpp (added)
• src/Parser/MpqcParser.cpp (added)
• src/Parser/MpqcParser.hpp (added)
• src/Parser/PcpParser.cpp (added)
• src/Parser/PcpParser.hpp (added)
• src/Parser/TremoloParser.cpp (modified) (2 diffs)
• src/Parser/TremoloParser.hpp (modified) (1 diff)
• src/Parser/XyzParser.cpp (modified) (2 diffs)
• src/Parser/XyzParser.hpp (modified) (1 diff)
• src/Patterns/Singleton.hpp (modified) (3 diffs)
• src/Patterns/Singleton_impl.hpp (modified) (2 diffs)
• src/ThermoStatContainer.cpp (added)
• src/ThermoStatContainer.hpp (added)
• src/UIElements/CommandLineUI/CommandLineWindow.cpp (modified) (3 diffs)
• src/UIElements/Dialog.cpp (modified) (1 diff)
• src/World.cpp (modified) (9 diffs)
• src/World.hpp (modified) (3 diffs)
• src/atom.cpp (modified) (2 diffs)
• src/atom.hpp (modified) (1 diff)
• src/atom_trajectoryparticle.cpp (modified) (4 diffs)
• src/builder.cpp (modified) (7 diffs)
• src/config.cpp (modified) (20 diffs)
• src/config.hpp (modified) (6 diffs)
• src/molecule.cpp (modified) (5 diffs)
• src/molecule.hpp (modified) (5 diffs)
• src/molecule_dynamics.cpp (modified) (6 diffs)
• src/molecule_fragmentation.cpp (modified) (20 diffs)
• src/molecule_graph.cpp (modified) (7 diffs)
• src/moleculelist.cpp (modified) (12 diffs)
• src/periodentafel.cpp (modified) (2 diffs)
• src/tesselation.cpp (modified) (3 diffs)
• src/unittests/LinkedCellUnitTest.cpp (modified) (2 diffs)
• src/unittests/Makefile.am (modified) (3 diffs)
• src/unittests/ObserverTest.cpp (modified) (2 diffs)
• src/unittests/ParserUnitTest.cpp (modified) (11 diffs)
• src/unittests/ParserUnitTest.hpp (modified) (2 diffs)
• src/unittests/manipulateAtomsTest.cpp (modified) (3 diffs)
• tests/regression/Domain/6/post/test-x.conf.xyz (modified) (1 diff)
• tests/regression/Domain/6/post/test-y.conf.xyz (modified) (1 diff)
• tests/regression/Domain/6/post/test-z.conf.xyz (modified) (1 diff)
• tests/regression/Domain/6/post/test.conf.xyz (modified) (1 diff)
• tests/regression/Filling/1/post/test.conf (modified) (3 diffs)
• tests/regression/Simple_configuration/2/post/test.conf.xyz (modified) (1 diff)
• tests/regression/Simple_configuration/3/post/test.conf.xyz (modified) (1 diff)
• tests/regression/Simple_configuration/4/post/test.conf (modified) (3 diffs)
• tests/regression/Simple_configuration/4/post/test.conf.xyz (modified) (1 diff)
• tests/regression/Simple_configuration/5/post/test.conf (modified) (3 diffs)
• tests/regression/Simple_configuration/8/post/test.conf.xyz (modified) (1 diff)
• tests/regression/testsuite-fragmentation.at (modified) (2 diffs)
• tests/regression/testsuite-simple_configuration.at (modified) (2 diffs)
• tests/regression/testsuite-standard_options.at (modified) (1 diff)

src/Actions/ActionRegistry.cpp

@@ -19 +19 @@
 using namespace std;
 
+/** Constructor for class ActionRegistry.
+ */
 ActionRegistry::ActionRegistry()
 {
 }
 
+/** Destructor for class ActionRegistry.
+ */
 ActionRegistry::~ActionRegistry()
 {
@@ -32 +36 @@
 }
 
+/** Returns pointer to an action named by \a name.
+ * \param name name of action
+ * \return pointer to Action
+ */
 Action* ActionRegistry::getActionByName(const std::string name){
   map<const string,Action*>::iterator iter;
@@ -39 +47 @@
 }
 
+/** States whether action is present or not.
+ * \note This iss needed as ActionRegistry::getActionByName() ASSERT()s that action is in map.
+ * \param name name of action
+ * \return true - Action present, false - Action absent
+ */
 bool ActionRegistry::isActionByNamePresent(const std::string name){
   map<const string,Action*>::iterator iter;
@@ -45 +58 @@
 }
 
+/** Registers an Action with the ActionRegistry.
+ * \param *action pointer to Action.
+ */
 void ActionRegistry::registerAction(Action* action){
   pair<map<const string,Action*>::iterator,bool> ret;
+  //cout << "Trying to register action with name " << action->getName() << "." << endl;
   ret = actionMap.insert(pair<const string,Action*>(action->getName(),action));
   ASSERT(ret.second,"Two actions with the same name added to registry");
 }
 
+/** Unregisters an Action.
+ * \param *action pointer to Action.
+ */
 void ActionRegistry::unregisterAction(Action* action){
+  //cout << "Unregistering action with name " << action->getName() << "." << endl;
   actionMap.erase(action->getName());
 }
 
+/** Returns an iterator pointing to the start of the map of Action's.
+ * \return begin iterator
+ */
 std::map<const std::string,Action*>::iterator ActionRegistry::getBeginIter()
 {
@@ -60 +84 @@
 }
 
+/** Returns an iterator pointing to the end of the map of Action's.
+ * \return end iterator
+ */
 std::map<const std::string,Action*>::iterator ActionRegistry::getEndIter()
 {
@@ -65 +92 @@
 }
 
+/** Returns a const iterator pointing to the start of the map of Action's.
+ * \return constant begin iterator
+ */
+std::map<const std::string,Action*>::const_iterator ActionRegistry::getBeginIter() const
+{
+  return actionMap.begin();
+}
+
+/** Returns a const iterator pointing to the end of the map of Action's.
+ * \return constant end iterator
+ */
+std::map<const std::string,Action*>::const_iterator ActionRegistry::getEndIter() const
+{
+  return actionMap.end();
+}
+
+/** Prints the contents of the ActionRegistry \a &m to \a &ost.
+ * \param &ost output stream
+ * \param &m reference to ActionRegistry
+ * \return reference to the above out stream for concatenation
+ */
+ostream& operator<<(ostream& ost, const ActionRegistry& m)
+{
+  ost << "ActionRegistry contains:" << endl;
+  for (std::map<const std::string,Action*>::const_iterator iter = m.getBeginIter(); iter != m.getEndIter(); ++iter) {
+    ost << "\t" << iter->first << " with pointer " << iter->second << endl;
+  }
+  return ost;
+};
+
+
+
 CONSTRUCT_SINGLETON(ActionRegistry)

src/Actions/ActionRegistry.hpp

@@ -9 +9 @@
 #define ACTIONREGISTRY_HPP_
 
+#include <iostream>
 #include <string>
 #include <map>
@@ -26 +27 @@
 
   std::map<const std::string,Action*>::iterator getBeginIter();
+  std::map<const std::string,Action*>::const_iterator getBeginIter() const;
   std::map<const std::string,Action*>::iterator getEndIter();
+  std::map<const std::string,Action*>::const_iterator getEndIter() const;
 
 private:
@@ -36 +39 @@
 };
 
+std::ostream& operator<<(std::ostream& ost, const ActionRegistry& m);
+
 #endif /* ACTIONREGISTRY_HPP_ */

src/Actions/CmdAction/BondLengthTableAction.cpp

@@ -9 +9 @@
 
 #include "Actions/CmdAction/BondLengthTableAction.hpp"
+#include "bondgraph.hpp"
 #include "config.hpp"
 #include "log.hpp"

src/Actions/FragmentationAction/DepthFirstSearchAction.cpp

@@ -10 +10 @@
 #include "Actions/FragmentationAction/DepthFirstSearchAction.hpp"
 #include "atom.hpp"
+#include "bondgraph.hpp"
 #include "config.hpp"
 #include "log.hpp"

src/Actions/FragmentationAction/FragmentationAction.cpp

@@ -10 +10 @@
 #include "Actions/FragmentationAction/FragmentationAction.hpp"
 #include "atom.hpp"
+#include "bondgraph.hpp"
 #include "config.hpp"
 #include "log.hpp"
@@ -41 +42 @@
   double distance = -1.;
   int order = 0;
+  std::string path;
   config *configuration = World::getInstance().getConfig();
   int ExitFlag = 0;
 
   cout << "pre-dialog"<< endl;
-  dialog->queryMolecule(NAME, &mol, MapOfActions::getInstance().getDescription(NAME));
+  dialog->queryString(NAME, &path, MapOfActions::getInstance().getDescription(NAME));
+  dialog->queryMolecule("molecule-by-id", &mol, MapOfActions::getInstance().getDescription("molecule-by-id"));
   dialog->queryDouble("distance", &distance, MapOfActions::getInstance().getDescription("distance"));
   dialog->queryInt("order", &order, MapOfActions::getInstance().getDescription("order"));
@@ -58 +61 @@
     DoLog(0) && (Log() << Verbose(0) << "Fragmenting molecule with current connection matrix ..." << endl);
     if (mol->hasBondStructure()) {
-      ExitFlag = mol->FragmentMolecule(order, configuration);
+      ExitFlag = mol->FragmentMolecule(order, path);
     }
     World::getInstance().setExitFlag(ExitFlag);

src/Actions/Makefile.am

@@ -120 +120 @@
   WorldAction/CenterOnEdgeAction.cpp \
   WorldAction/ChangeBoxAction.cpp \
+  WorldAction/InputAction.cpp \
+  WorldAction/OutputAction.cpp \
   WorldAction/RemoveSphereOfAtomsAction.cpp \
   WorldAction/RepeatBoxAction.cpp \
@@ -131 +133 @@
   WorldAction/CenterOnEdgeAction.hpp \
   WorldAction/ChangeBoxAction.hpp \
+  WorldAction/InputAction.hpp \
+  WorldAction/OutputAction.hpp \
   WorldAction/RemoveSphereOfAtomsAction.hpp \
   WorldAction/RepeatBoxAction.hpp \

src/Actions/MapOfActions.cpp

@@ -83 +83 @@
   DescriptionMap["fragment-mol"] = "create for a given molecule into fragments up to given order";
   DescriptionMap["help"] = "Give this help screen";
+  DescriptionMap["input"] = "specify input files";
   DescriptionMap["linear-interpolate"] = "linear interpolation in discrete steps between start and end position of a molecule";
   DescriptionMap["nonconvex-envelope"] = "create the non-convex envelope for a molecule";
   DescriptionMap["molecular-volume"] = "calculate the volume of a given molecule";
+  DescriptionMap["output"] = "specify output formats";
   DescriptionMap["pair-correlation"] = "pair correlation analysis between two elements, element and point or element and surface";
   DescriptionMap["parse-xyz"] = "parse xyz file into World";
@@ -185 +187 @@
   TypeMap["fastparsing"] = Boolean;
   TypeMap["fill-molecule"] = String;
-  TypeMap["fragment-mol"] = Molecule;
+  TypeMap["fragment-mol"] = String;
   TypeMap["input"] = String;
   TypeMap["linear-interpolate"] = String;
   TypeMap["molecular-volume"] = Molecule;
   TypeMap["nonconvex-envelope"] = Molecule;
+  TypeMap["output"] = String;
   TypeMap["parse-xyz"] = String;
   TypeMap["pair-correlation"] = String;
@@ -262 +265 @@
   generic.insert("fragment-mol");
   generic.insert("help");
-        generic.insert("linear-interpolate");
+  generic.insert("input");
+  generic.insert("linear-interpolate");
 //  generic.insert("molecular-volume");
   generic.insert("nonconvex-envelope");
+  generic.insert("output");
         generic.insert("pair-correlation");
-//      generic.insert("parse-xyz");
+        generic.insert("parse-xyz");
 //  generic.insert("principal-axis-system");
   generic.insert("remove-atom");

src/Actions/MoleculeAction/FillWithMoleculeAction.cpp

@@ -102 +102 @@
       World::getInstance().getMolecules()->insert(Filling);
     }
+    for (molecule::iterator iter = filler->begin(); !filler->empty(); iter = filler->begin()) {
+      atom *Walker = *iter;
+      filler->erase(iter);
+      World::getInstance().destroyAtom(Walker);
+    }
     World::getInstance().destroyMolecule(filler);
 

src/Actions/MoleculeAction/LinearInterpolationofTrajectoriesAction.cpp

@@ -69 +69 @@
     if (IdMapping)
       DoLog(1) && (Log() << Verbose(1) << "Using Identity for the permutation map." << endl);
-    char outputname[MAXSTRINGSIZE];
-    strcpy(outputname, filename.c_str());
-    // TODO: LinearInterpolationBetweenConfiguration should use stream, not the filename directly! (better for unit test)
-    if (!mol->LinearInterpolationBetweenConfiguration(start, end, outputname, *(World::getInstance().getConfig()), IdMapping))
-      DoLog(2) && (Log() << Verbose(2) << "Could not store " << outputname << " files." << endl);
+    if (!mol->LinearInterpolationBetweenConfiguration(start, end, filename, *(World::getInstance().getConfig()), IdMapping))
+      DoLog(2) && (Log() << Verbose(2) << "Could not store " << filename << " files." << endl);
     else
       DoLog(2) && (Log() << Verbose(2) << "Steps created and " << filename << " files stored." << endl);

src/Actions/MoleculeAction/SaveAdjacencyAction.cpp

@@ -65 +65 @@
     World::getInstance().getConfig()->BG->ConstructBondGraph(mol);
     // TODO: sollte stream nicht filename benutzen, besser fuer unit test
-    char outputname[MAXSTRINGSIZE];
-    strcpy(outputname, filename.c_str());
-    mol->StoreAdjacencyToFile(NULL, outputname);
+    mol->StoreAdjacencyToFile(filename);
     delete dialog;
     return Action::success;

src/Actions/MoleculeAction/SaveBondsAction.cpp

@@ -64 +64 @@
     DoLog(0) && (Log() << Verbose(0) << "Storing bonds to path " << filename << "." << endl);
     World::getInstance().getConfig()->BG->ConstructBondGraph(mol);
-    // TODO: sollte stream, nicht filenamen direkt nutzen, beser fuer unit tests
-    char outputname[MAXSTRINGSIZE];
-    strcpy(outputname, filename.c_str());
-    mol->StoreBondsToFile(NULL, outputname);
+    // TODO: sollte stream, nicht filenamen direkt nutzen, besser fuer unit tests
+    mol->StoreBondsToFile(filename);
     delete dialog;
     return Action::success;

src/Actions/ParserAction/LoadXyzAction.cpp

@@ -8 +8 @@
 #include "Helpers/MemDebug.hpp"
 
+using namespace std;
+
 #include "Actions/ParserAction/LoadXyzAction.hpp"
 #include "Parser/XyzParser.hpp"
 
 #include <iostream>
+#include <set>
 #include <string>
+#include <vector>
 
-using namespace std;
 
 #include "UIElements/UIFactory.hpp"
@@ -21 +24 @@
 
 #include "atom.hpp"
+#include "log.hpp"
 #include "molecule.hpp"
+#include "verbose.hpp"
+#include "World.hpp"
 
 /****** ParserLoadXyzAction *****/
@@ -37 +43 @@
 //};
 
-const char ParserLoadXyzAction::NAME[] = "LoadXyz";
+const char ParserLoadXyzAction::NAME[] = "parse-xyz";
 
 ParserLoadXyzAction::ParserLoadXyzAction() :
@@ -48 +54 @@
 Action::state_ptr ParserLoadXyzAction::performCall() {
   string filename;
-  XyzParser parser;
   Dialog *dialog = UIFactory::getInstance().makeDialog();
 
-  dialog->queryString("filename",&filename, "Filename of the xyz file");
+  dialog->queryString(NAME,&filename, NAME);
 
   if(dialog->display()) {
+    DoLog(1) && (Log() << Verbose(1) << "Parsing xyz file for new atoms." << endl);
     // parse xyz file
     ifstream input;
     input.open(filename.c_str());
-    if (!input.fail())
+    if (!input.fail()) {
+      // TODO: Remove the insertion into molecule when saving does not depend on them anymore. Also, remove molecule.hpp include
+      set <atom*> UniqueList;
+      {
+        vector<atom *> ListBefore = World::getInstance().getAllAtoms();
+        for (vector<atom *>::iterator runner = ListBefore.begin();runner != ListBefore.end(); ++runner)
+          UniqueList.insert(*runner);
+      }
+      XyzParser parser; // briefly instantiate a parser which is removed at end of focus
       parser.load(&input);
+      {
+        vector<atom *> ListAfter = World::getInstance().getAllAtoms();
+        pair< set<atom *>::iterator, bool > Inserter;
+        if (UniqueList.size() != ListAfter.size()) { // only create if new atoms have been parsed
+          MoleculeListClass *molecules = World::getInstance().getMolecules();
+          molecule *mol= NULL;
+          if (molecules->ListOfMolecules.empty()) {
+            mol = World::getInstance().createMolecule();
+            molecules->insert(mol);
+          } else {
+            mol = *(molecules->ListOfMolecules.begin());
+          }
+          for (vector<atom *>::iterator runner = ListAfter.begin(); runner != ListAfter.end(); ++runner) {
+            Inserter = UniqueList.insert(*runner);
+            if (Inserter.second) { // if not present, then new (just parsed) atom, add ...
+              cout << "Adding new atom " << **runner << " to new mol." << endl;
+              mol->AddAtom(*runner);
+            }
+          }
+          mol->doCountAtoms();
+        } else {
+          cout << "No atoms parsed?" << endl;
+        }
+      }
+    } else {
+      DoeLog(1) && (eLog() << Verbose(1) << "Could not open file " << filename << "." << endl);
+    }
     input.close();
   }

src/CommandLineParser.cpp

@@ -56 +56 @@
 {
   // go through all arguments
+  cout << Verbose(1) << "By default putting input into the sequence." << endl;
+  // TODO: This may be bad, because const string "input" is destroyed at end of function
+  SequenceOfActions.push_back("input");
   for (int i=1;i<argc;i++) {
     (cout << Verbose(1) << "Checking on " << argv[i] << endl);

src/Descriptors/AtomDescriptor.cpp

@@ -12 +12 @@
 
 #include "World.hpp"
-
 #include "atom.hpp"
 
 #include <boost/bind.hpp>
-#include <cassert>
+
 #include <iostream>
 

src/Descriptors/MoleculeDescriptor.cpp

@@ -16 +16 @@
 
 #include <boost/bind.hpp>
-#include <cassert>
 #include <iostream>
 

src/Helpers/MemDebug.cpp

@@ -6 +6 @@
  */
 
-#ifndef NDBEGUG
+#ifndef NDEBUG
 #ifndef NO_MEMDEBUG
 

src/Legacy/oldmenu.cpp

@@ -609 +609 @@
 {
   int Order1;
+  std::string path;
   clock_t start, end;
 
@@ -614 +615 @@
   Log() << Verbose(0) << "What's the desired bond order: ";
   cin >> Order1;
+  DoLog(0) && (Log() << Verbose(0) << "What's the output path and prefix [e.g. /home/foo/BondFragment]: ");
+  cin >> path;
   if (mol->hasBondStructure()) {
     start = clock();
-    mol->FragmentMolecule(Order1, configuration);
+    mol->FragmentMolecule(Order1, path);
     end = clock();
     Log() << Verbose(0) << "Clocks for this operation: " << (end-start) << ", time: " << ((double)(end-start)/CLOCKS_PER_SEC) << "s." << endl;

src/Makefile.am

@@ -63 +63 @@
 
 ACTIONSHEADER = \
-  ${ANALYSISACTIONHEADER} \
-  ${ATOMACTIONHEADER} \
-  ${CMDACTIONHEADER} \
-  ${FRAGMENTATIONACTIONHEADER} \
-  ${MOLECULEACTIONHEADER} \
-  ${PARSERACTIONHEADER} \
-  ${TESSELATIONACTIONHEADER} \
-  ${WORLDACTIONHEADER} \
   Actions/Action.hpp \
   Actions/ActionHistory.hpp \
@@ -88 +80 @@
   Parser/ChangeTracker.cpp \
   Parser/FormatParser.cpp \
+  Parser/FormatParserStorage.cpp \
+  Parser/MpqcParser.cpp \
+  Parser/PcpParser.cpp \
   Parser/TremoloParser.cpp \
   Parser/XyzParser.cpp
+
 PARSERHEADER = \
   Parser/ChangeTracker.hpp \
   Parser/FormatParser.hpp \
+  Parser/FormatParserStorage.hpp \
+  Parser/MpqcParser.hpp \
+  Parser/PcpParser.hpp \
   Parser/TremoloParser.hpp \
   Parser/XyzParser.hpp
@@ -183 +182 @@
   CommandLineParser.cpp \
   config.cpp \
+  ConfigFileBuffer.cpp \
   element.cpp \
   elements_db.cpp \
@@ -208 +208 @@
   tesselation.cpp \
   tesselationhelpers.cpp \
+  ThermoStatContainer.cpp \
   triangleintersectionlist.cpp \
   vector.cpp \
@@ -228 +229 @@
   CommandLineParser.hpp \
   config.hpp \
+  ConfigFileBuffer.hpp \
   defs.hpp \
   element.hpp \
@@ -250 +252 @@
   tesselation.hpp \
   tesselationhelpers.hpp \
+  ThermoStatContainer.hpp \
   triangleintersectionlist.hpp \
   verbose.hpp \
@@ -265 +268 @@
 INCLUDES = -I$(top_srcdir)/src/unittests -I$(top_srcdir)/src/Actions -I$(top_srcdir)/src/UIElements
 
-noinst_LIBRARIES = libmolecuilder.a libgslwrapper.a libmenu.a
-bin_PROGRAMS = molecuilder molecuildergui joiner analyzer 
+noinst_LIBRARIES = libmolecuilder.a libgslwrapper.a libmenu.a libparser.a
+bin_PROGRAMS = molecuilder molecuildergui joiner analyzer
+
 molecuilderdir = ${bindir}
 
@@ -272 +276 @@
 
 libmenu_a_SOURCES = ${UISOURCE} ${UIHEADER} 
-
+libparser_a_SOURCES = ${PARSERSOURCE} ${PARSERHEADER}
 libgslwrapper_a_SOURCES = ${LINALGSOURCE} ${LINALGHEADER}
 
@@ -282 +286 @@
 molecuilder_SOURCES = ${LEGACYSOURCE} builder.cpp
 molecuilder_SOURCES += $(srcdir)/version.c
-molecuilder_LDADD = UIElements/libMolecuilderUI.a Actions/libMolecuilderActions.a libmolecuilder.a libgslwrapper.a $(BOOST_LIB) ${BOOST_THREAD_LIB} ${BOOST_PROGRAM_OPTIONS_LIB}
+molecuilder_LDADD = UIElements/libMolecuilderUI.a Actions/libMolecuilderActions.a libmolecuilder.a libparser.a libgslwrapper.a $(BOOST_LIB) ${BOOST_THREAD_LIB} ${BOOST_PROGRAM_OPTIONS_LIB}
 
 #Stuff for building the GUI using QT
@@ -289 +293 @@
 molecuildergui_CXXFLAGS = ${QT_CXXFLAGS} ${GLU_CXXFLAGS} -DUSE_GUI_QT
 molecuildergui_LDFLAGS = $(BOOST_LIB) ${QT_LDFLAGS} ${GLU_LDFLAGS}
-molecuildergui_LDADD = UIElements/libMolecuilderUI.a Actions/libMolecuilderActions.a libmolecuilder.a libgslwrapper.a $(BOOST_LIB) ${BOOST_THREAD_LIB} ${BOOST_PROGRAM_OPTIONS_LIB}  ${GUI_LIBS}
+molecuildergui_LDADD = UIElements/libMolecuilderUI.a Actions/libMolecuilderActions.a libmolecuilder.a libparser.a libgslwrapper.a $(BOOST_LIB) ${BOOST_THREAD_LIB} ${BOOST_PROGRAM_OPTIONS_LIB}  ${GUI_LIBS}
 
 joiner_SOURCES = joiner.cpp datacreator.cpp parser.cpp datacreator.hpp helpers.hpp parser.hpp periodentafel.hpp

src/Parser/ChangeTracker.cpp

@@ -7 +7 @@
 
 #include "Helpers/MemDebug.hpp"
+#include "Parser/ChangeTracker.hpp"
+#include "Patterns/Singleton_impl.hpp"
 
-#include "ChangeTracker.hpp"
-
-ChangeTracker* ChangeTracker::instance = NULL;
 
 /**
@@ -27 +26 @@
 ChangeTracker::~ChangeTracker() {
   World::getInstance().signOff(this);
-}
-
-/**
- * Returns the change tracker instance.
- *
- * \return this
- */
-ChangeTracker* ChangeTracker::get() {
-  if (instance == NULL) {
-    instance = new ChangeTracker();
-  }
-
-  return instance;
-}
-
-/**
- * Destroys the change tracker instance. Be careful, the change tracker is a
- * singleton and destruction might lead to a loss of consistency.
- */
-void ChangeTracker::destroy() {
-  delete instance;
-  instance = NULL;
 }
 
@@ -76 +53 @@
   }
 }
+
+CONSTRUCT_SINGLETON(ChangeTracker)

src/Parser/ChangeTracker.hpp

@@ -18 +18 @@
  * changes to it.
  */
-class ChangeTracker : public Observable {
+class ChangeTracker : public Singleton<ChangeTracker>, public Observable {
+  friend class Singleton<ChangeTracker>;
 public:
   void saveStatus();
@@ -31 +32 @@
   bool isConsistent;
   static ChangeTracker* instance;
+
+  // private constructor and destructor due to singleton
   ChangeTracker();
   ~ChangeTracker();

src/Parser/FormatParser.cpp

@@ -19 +19 @@
   Observer("FormatParser")
 {
-  ChangeTracker::get()->signOn(this);
+  ChangeTracker::getInstance().signOn(this);
   saveStream = NULL;
 }
@@ -27 +27 @@
  */
 FormatParser::~FormatParser() {
-  ChangeTracker::get()->signOff(this);
+  ChangeTracker::getInstance().signOff(this);
 }
 

src/Parser/TremoloParser.cpp

@@ -49 +49 @@
   knownKeys["GrpTypeNo"] = TremoloKey::GrpTypeNo;
   knownKeys["torsion"] = TremoloKey::torsion;
+
+  // default behavior: use all possible keys on output
+  for (std::map<std::string, TremoloKey::atomDataKey>::iterator iter = knownKeys.begin(); iter != knownKeys.end(); ++iter)
+    usedFields.push_back(iter->first);
 }
 
@@ -193 +197 @@
 
   lineStream << line.substr(offset);
+  usedFields.clear();
   while (lineStream.good()) {
     lineStream >> keyword;

src/Parser/TremoloParser.hpp

@@ -10 +10 @@
 
 #include <string>
-#include "FormatParser.hpp"
+#include "Parser/FormatParser.hpp"
 
 /**

src/Parser/XyzParser.cpp

@@ -27 +27 @@
  */
 XyzParser::~XyzParser() {
-  delete(&comment);
 }
 
@@ -58 +57 @@
  */
 void XyzParser::save(ostream* file) {
+  if (comment == "") {
+    time_t now = time((time_t *)NULL);   // Get the system time and put it into 'now' as 'calender time'
+    comment = "\tCreated by molecuilder on ";
+    // ctime ends in \n\0, we have to cut away the newline
+    std::string time(ctime(&now));
+    size_t pos = time.find('\n');
+    if (pos != 0)
+      comment += time.substr(0,pos);
+    else
+      comment += time;
+  }
   *file << World::getInstance().numAtoms() << endl << comment << endl;
 
   vector<atom*> atoms = World::getInstance().getAllAtoms();
   for(vector<atom*>::iterator it = atoms.begin(); it != atoms.end(); it++) {
-    *file << fixed << (*it)->getType()->symbol << "\t" << (*it)->x[0] << "\t" << (*it)->x[1] << "\t" << (*it)->x[2] << endl;
+    *file << noshowpoint << (*it)->getType()->symbol << "\t" << (*it)->x[0] << "\t" << (*it)->x[1] << "\t" << (*it)->x[2] << endl;
   }
 }

src/Parser/XyzParser.hpp

@@ -10 +10 @@
 
 #include <string>
-#include "FormatParser.hpp"
+#include "Parser/FormatParser.hpp"
 
 /**

src/Patterns/Singleton.hpp

@@ -9 +9 @@
 #define SINGLETON_HPP_
 
-#include <cassert>
 #include <boost/thread.hpp>
 
+#include "Helpers/Assert.hpp"
 #include "defs.hpp"
 
@@ -181 +181 @@
 
     inline static void set(creator_T*&,creator_T*){
-      assert(0 && "Cannot set the Instance for a singleton of this type");
+      ASSERT(0, "Cannot set the Instance for a singleton of this type");
     }
   };
@@ -191 +191 @@
   struct creator_t<creator_T,false>{
     inline static creator_T* make(){
-      assert(0 && "Cannot create a singleton of this type directly");
+      ASSERT(0, "Cannot create a singleton of this type directly");
+      return 0;
     }
     inline static void set(ptr_t& dest,creator_T* src){

src/Patterns/Singleton_impl.hpp

@@ -72 +72 @@
 template <class T,bool _may_create>
 void Singleton<T,_may_create>::setInstance(T* newInstance){
-  assert(!theInstance.get() && "Trying to set the instance of an already created singleton");
+  ASSERT(!theInstance.get(), "Trying to set the instance of an already created singleton");
   boost::recursive_mutex::scoped_lock guard(instanceLock);
   theInstance.reset(newInstance);
@@ -83 +83 @@
 template <class T, bool _may_create>
 Singleton<T,_may_create>::Singleton(const Singleton<T,_may_create>&){
-  assert(0 && "Copy constructor of singleton template called");
+  ASSERT(0, "Copy constructor of singleton template called");
 }
 

src/UIElements/CommandLineUI/CommandLineWindow.cpp

@@ -47 +47 @@
 #include "Actions/WorldAction/CenterOnEdgeAction.hpp"
 #include "Actions/WorldAction/ChangeBoxAction.hpp"
+#include "Actions/WorldAction/InputAction.hpp"
+#include "Actions/WorldAction/OutputAction.hpp"
 #include "Actions/WorldAction/RemoveSphereOfAtomsAction.hpp"
 #include "Actions/WorldAction/RepeatBoxAction.hpp"
@@ -81 +83 @@
 
 void CommandLineWindow::display() {
+  //cout << ActionRegistry::getInstance();
+
   // go through all possible actions
   for (std::list<std::string>::iterator CommandRunner = CommandLineParser::getInstance().SequenceOfActions.begin(); CommandRunner != CommandLineParser::getInstance().SequenceOfActions.end(); ++CommandRunner) {
-    cout << "Checking presence of " << *CommandRunner << endl;
-    if (ActionRegistry::getInstance().isActionByNamePresent(*CommandRunner))
+    cout << "Checking presence of " << *CommandRunner << ": ";
+    if (ActionRegistry::getInstance().isActionByNamePresent(*CommandRunner)) {
+      cout << "calling " << *CommandRunner << endl;
       ActionRegistry::getInstance().getActionByName(*CommandRunner)->call();
+    } else {
+      cout << "absent." << endl;
+    }
   }
 }
@@ -152 +160 @@
   new WorldCenterOnEdgeAction();
   new WorldChangeBoxAction();
+  new WorldInputAction();
+  new WorldOutputAction();
   new WorldRemoveSphereOfAtomsAction();
   new WorldRepeatBoxAction();

src/UIElements/Dialog.cpp

@@ -7 +7 @@
 
 #include "Helpers/MemDebug.hpp"
-
-#include <cassert>
 
 #include "Dialog.hpp"

src/World.cpp

@@ -14 +14 @@
 #include "molecule.hpp"
 #include "periodentafel.hpp"
+#include "ThermoStatContainer.hpp"
 #include "Descriptors/AtomDescriptor.hpp"
 #include "Descriptors/AtomDescriptor_impl.hpp"
@@ -20 +21 @@
 #include "Descriptors/SelectiveIterator_impl.hpp"
 #include "Actions/ManipulateAtomsProcess.hpp"
+#include "Helpers/Assert.hpp"
 
 #include "Patterns/Singleton_impl.hpp"
@@ -90 +92 @@
 };
 
+class ThermoStatContainer * World::getThermostats()
+{
+  return Thermostats;
+}
+
+
 int World::getExitFlag() {
   return ExitFlag;
@@ -105 +113 @@
   molecule *mol = NULL;
   mol = NewMolecule();
-  assert(!molecules.count(currMoleculeId));
+  ASSERT(!molecules.count(currMoleculeId),"currMoleculeId did not specify an unused ID");
   mol->setId(currMoleculeId++);
   // store the molecule by ID
@@ -121 +129 @@
   OBSERVE;
   molecule *mol = molecules[id];
-  assert(mol);
+  ASSERT(mol,"Molecule id that was meant to be destroyed did not exist");
   DeleteMolecule(mol);
   molecules.erase(id);
@@ -157 +165 @@
   OBSERVE;
   atom *atom = atoms[id];
-  assert(atom);
+  ASSERT(atom,"Atom ID that was meant to be destroyed did not exist");
   DeleteAtom(atom);
   atoms.erase(id);
@@ -169 +177 @@
   if(!target){
     target = atoms[oldId];
-    assert(target && "Atom with that ID not found");
+    ASSERT(target,"Atom with that ID not found");
     return target->changeId(newId);
   }
@@ -281 +289 @@
     periode(new periodentafel),
     configuration(new config),
+    Thermostats(new ThermoStatContainer),
     ExitFlag(0),
     atoms(),
@@ -306 +315 @@
   delete periode;
   delete configuration;
+  delete Thermostats;
   MoleculeSet::iterator molIter;
   for(molIter=molecules.begin();molIter!=molecules.end();++molIter){

src/World.hpp

@@ -30 +30 @@
 
 // forward declarations
-class config;
-class periodentafel;
-class MoleculeListClass;
 class atom;
-class molecule;
 class AtomDescriptor;
 class AtomDescriptor_impl;
+template<typename T> class AtomsCalculation;
+class config;
+class ManipulateAtomsProcess;
+class molecule;
 class MoleculeDescriptor;
 class MoleculeDescriptor_impl;
-class ManipulateAtomsProcess;
-template<typename T>
-class AtomsCalculation;
+class MoleculeListClass;
+class periodentafel;
+class ThermoStatContainer;
 
 /****************************************** forward declarations *****************************/
@@ -142 +142 @@
   void setDefaultName(std::string name);
 
+  /**
+   * get pointer to World's ThermoStatContainer
+   */
+  ThermoStatContainer * getThermostats();
+
   /*
    * get the ExitFlag
@@ -254 +259 @@
   static double *cell_size;
   std::string defaultName;
+  class ThermoStatContainer *Thermostats;
   int ExitFlag;
 public:

src/atom.cpp

@@ -159 +159 @@
   * \return true - \a *out present, false - \a *out is NULL
  */
-bool atom::OutputArrayIndexed(ofstream * const out, const int *ElementNo, int *AtomNo, const char *comment) const
+bool atom::OutputArrayIndexed(ostream * const out, const int *ElementNo, int *AtomNo, const char *comment) const
 {
   AtomNo[type->Z]++;  // increment number
@@ -236 +236 @@
  * \param *AtomNo pointer to atom counter that is increased by one
  */
-void atom::OutputMPQCLine(ofstream * const out, const Vector *center, int *AtomNo = NULL) const
+void atom::OutputMPQCLine(ostream * const out, const Vector *center, int *AtomNo = NULL) const
 {
   *out << "\t\t" << type->symbol << " [ " << x[0]-center->at(0) << "\t" << x[1]-center->at(1) << "\t" << x[2]-center->at(2) << " ]" << endl;

src/atom.hpp

@@ -51 +51 @@
 
   bool OutputIndexed(ofstream * const out, const int ElementNo, const int AtomNo, const char *comment = NULL) const;
-  bool OutputArrayIndexed(ofstream * const out, const int *ElementNo, int *AtomNo, const char *comment = NULL) const;
+  bool OutputArrayIndexed(ostream * const out, const int *ElementNo, int *AtomNo, const char *comment = NULL) const;
   bool OutputXYZLine(ofstream *out) const;
   bool OutputTrajectory(ofstream * const out, const int *ElementNo, int *AtomNo, const int step) const;
   bool OutputTrajectoryXYZ(ofstream * const out, const int step) const;
-  void OutputMPQCLine(ofstream * const out, const Vector *center, int *AtomNo) const;
+  void OutputMPQCLine(ostream * const out, const Vector *center, int *AtomNo) const;
 
   void InitComponentNr();

src/atom_trajectoryparticle.cpp

@@ -15 +15 @@
 #include "log.hpp"
 #include "parser.hpp"
+#include "ThermoStatContainer.hpp"
 #include "verbose.hpp"
 
@@ -197 +198 @@
 void TrajectoryParticle::Thermostat_Langevin(int Step, gsl_rng * r, double *ekin, config *configuration)
 {
-  double sigma  = sqrt(configuration->TargetTemp/type->mass); // sigma = (k_b T)/m (Hartree/atomicmass = atomiclength/atomictime)
+  double sigma  = sqrt(configuration->Thermostats->TargetTemp/type->mass); // sigma = (k_b T)/m (Hartree/atomicmass = atomiclength/atomictime)
   Vector &U = Trajectory.U.at(Step);
   if (FixedIon == 0) { // even FixedIon moves, only not by other's forces
     // throw a dice to determine whether it gets hit by a heat bath particle
-    if (((((rand()/(double)RAND_MAX))*configuration->TempFrequency) < 1.)) {
+    if (((((rand()/(double)RAND_MAX))*configuration->Thermostats->TempFrequency) < 1.)) {
       DoLog(3) && (Log() << Verbose(3) << "Particle " << *this << " was hit (sigma " << sigma << "): " << sqrt(U[0]*U[0]+U[1]*U[1]+U[2]*U[2]) << " -> ");
       // pick three random numbers from a Boltzmann distribution around the desired temperature T for each momenta axis
@@ -225 +226 @@
   if (FixedIon == 0) { // even FixedIon moves, only not by other's forces
     for (int d=0; d<NDIM; d++) {
-      U[d] *= sqrt(1+(configuration->Deltat/configuration->TempFrequency)*(ScaleTempFactor-1));
+      U[d] *= sqrt(1+(configuration->Deltat/configuration->Thermostats->TempFrequency)*(ScaleTempFactor-1));
       *ekin += 0.5*type->mass * U[d]*U[d];
     }
@@ -255 +256 @@
   if (FixedIon == 0) { // even FixedIon moves, only not by other's forces
     for (int d=0; d<NDIM; d++) {
-        U[d] += configuration->Deltat/type->mass * (configuration->alpha * (U[d] * type->mass));
+        U[d] += configuration->Deltat/type->mass * (configuration->Thermostats->alpha * (U[d] * type->mass));
         *ekin += (0.5*type->mass) * U[d]*U[d];
       }

src/builder.cpp

@@ -1 +1 @@
 /** \file builder.cpp
  *
- * By stating absolute positions or binding angles and distances atomic positions of a molecule can be constructed.
- * The output is the complete configuration file for PCP for direct use.
- * Features:
- * -# Atomic data is retrieved from a file, if not found requested and stored there for later re-use
- * -# step-by-step construction of the molecule beginning either at a centre of with a certain atom
+ *  date: Jan 1, 2007
+ *  author: heber
  *
  */
 
-/*! \mainpage Molecuilder - a molecular set builder
- *
- * This introductory shall briefly make aquainted with the program, helping in installing and a first run.
+/*! \mainpage MoleCuilder - a molecular set builder
+ *
+ * This introductory shall briefly make acquainted with the program, helping in installing and a first run.
  *
  * \section about About the Program
  *
- *  Molecuilder is a short program, written in C++, that enables the construction of a coordinate set for the
- *  atoms making up an molecule by the successive statement of binding angles and distances and referencing to
- *  already constructed atoms.
- *
- *  A configuration file may be written that is compatible to the format used by PCP - a parallel Car-Parrinello
- *  molecular dynamics implementation.
+ *  MoleCuilder is a program, written entirely in C++, that enables the construction of a coordinate set for the
+ *  atoms making up an molecule. It allows for both building of simple molecules by adding atom-wise giving bond
+ *  angles and distances or absolute coordinates, but also using them as templates. Regions can be specified and
+ *  ordered to be filled with a molecule in a certain manner. Greater conglomerations of molecules can be tesselated
+ *  and recognized as a region themselves to be subsequently surrounded by other (surface solvated) molecules.
+ *  In the end, MoleCuilder allows the construction of arbitrary nano structures, whether they be crystalline or
+ *  amorphic in nature.
+ *
  *
  * \section install Installation
@@ -32 +31 @@
  *  -# make clean uninstall: deletes .o-files and removes executable from the given binary directory\n
  *  -# make doxygen-doc: Creates these html pages out of the documented source
+ *  -# make check: Runs an extensive set of unit tests and a testsuite which also gives a good overview on the set of
+ *                 functions.
  *
  * \section run Running
@@ -37 +38 @@
  *  The program can be executed by running: ./molecuilder
  *
- *  Note, that it uses a database, called "elements.db", in the executable's directory. If the file is not found,
- *  it is created and any given data on elements of the periodic table will be stored therein and re-used on
- *  later re-execution.
- *
- * \section ref References
- *
- *  For the special configuration file format, see the documentation of pcp.
+ *  MoleCuilder has three interfaces at your disposal:
+ *  -# Textmenu: A simple interactive console-based menu, where awaits your choices and inputs in order to set atoms
+ *               as you like
+ *  -# CommandLineUI: Every command can also be chained up as a sequence of actions on the command line to be executed
+ *               with any user interaction.
+ *  -# GraphicalUI: A graphical user interface that also display the molecular structure being built and lots of other
+ *               informations to ease the construction of bigger geometries.
+ *
+ *  The supported output formats right now are:
+ *  -# mpqc: Configuration files of the Massively Parallel Quantum Chemistry package (Sandia labs)
+ *  -# pcp: Configuration files of the Parallel Car-Parrinello program (Institute for Numerical Simulation)
+ *  -# tremolo: Configuration files of TREMOLO (Institute for Numerical Simulation)
+ *  -# xyz: the most basic format for the 3d arrangement of atoms consisting of a list of element and 3 coordinates.
  *
  */
@@ -49 +56 @@
 #include "Helpers/MemDebug.hpp"
 
-#include <boost/bind.hpp>
-
-using namespace std;
-
-#include <cstring>
-#include <cstdlib>
-
-#include "analysis_bonds.hpp"
-#include "analysis_correlation.hpp"
-#include "atom.hpp"
-#include "bond.hpp"
 #include "bondgraph.hpp"
-#include "boundary.hpp"
 #include "CommandLineParser.hpp"
 #include "config.hpp"
-#include "element.hpp"
-#include "ellipsoid.hpp"
-#include "helpers.hpp"
-#include "leastsquaremin.hpp"
-#include "linkedcell.hpp"
 #include "log.hpp"
 #include "memoryusageobserver.hpp"
@@ -82 +72 @@
 #include "UIElements/Dialog.hpp"
 #include "Menu/ActionMenuItem.hpp"
+#include "verbose.hpp"
+#include "World.hpp"
+
 #include "Actions/ActionRegistry.hpp"
 #include "Actions/ActionHistory.hpp"
 #include "Actions/MapOfActions.hpp"
-#include "Actions/MethodAction.hpp"
-#include "Actions/MoleculeAction/ChangeNameAction.hpp"
-#include "World.hpp"
+
+#include "Parser/ChangeTracker.hpp"
+#include "Parser/FormatParserStorage.hpp"
+
+#include "UIElements/UIFactory.hpp"
+#include "UIElements/TextUI/TextUIFactory.hpp"
+#include "UIElements/CommandLineUI/CommandLineUIFactory.hpp"
+#include "UIElements/MainWindow.hpp"
+#include "UIElements/Dialog.hpp"
+
 #include "version.h"
-#include "World.hpp"
-
-
-/********************************************* Subsubmenu routine ************************************/
-#if 0
-/** Submenu for adding atoms to the molecule.
- * \param *periode periodentafel
- * \param *molecule molecules with atoms
- */
-static void AddAtoms(periodentafel *periode, molecule *mol)
-{
-  atom *first, *second, *third, *fourth;
-  Vector **atoms;
-  Vector x,y,z,n;  // coordinates for absolute point in cell volume
-  double a,b,c;
-  char choice;  // menu choice char
-  bool valid;
-
-  cout << Verbose(0) << "===========ADD ATOM============================" << endl;
-  cout << Verbose(0) << " a - state absolute coordinates of atom" << endl;
-  cout << Verbose(0) << " b - state relative coordinates of atom wrt to reference point" << endl;
-  cout << Verbose(0) << " c - state relative coordinates of atom wrt to already placed atom" << endl;
-  cout << Verbose(0) << " d - state two atoms, two angles and a distance" << endl;
-  cout << Verbose(0) << " e - least square distance position to a set of atoms" << endl;
-  cout << Verbose(0) << "all else - go back" << endl;
-  cout << Verbose(0) << "===============================================" << endl;
-  cout << Verbose(0) << "Note: Specifiy angles in degrees not multiples of Pi!" << endl;
-  cout << Verbose(0) << "INPUT: ";
-  cin >> choice;
-
-  switch (choice) {
-    default:
-      DoeLog(2) && (eLog()<< Verbose(2) << "Not a valid choice." << endl);
-      break;
-      case 'a': // absolute coordinates of atom
-        cout << Verbose(0) << "Enter absolute coordinates." << endl;
-        first = new atom;
-        first->x.AskPosition(World::getInstance().getDomain(), false);
-        first->type = periode->AskElement();  // give type
-        mol->AddAtom(first);  // add to molecule
-        break;
-
-      case 'b': // relative coordinates of atom wrt to reference point
-        first = new atom;
-        valid = true;
-        do {
-          if (!valid) DoeLog(2) && (eLog()<< Verbose(2) << "Resulting position out of cell." << endl);
-          cout << Verbose(0) << "Enter reference coordinates." << endl;
-          x.AskPosition(World::getInstance().getDomain(), true);
-          cout << Verbose(0) << "Enter relative coordinates." << endl;
-          first->x.AskPosition(World::getInstance().getDomain(), false);
-          first->x.AddVector((const Vector *)&x);
-          cout << Verbose(0) << "\n";
-        } while (!(valid = mol->CheckBounds((const Vector *)&first->x)));
-        first->type = periode->AskElement();  // give type
-        mol->AddAtom(first);  // add to molecule
-        break;
-
-      case 'c': // relative coordinates of atom wrt to already placed atom
-        first = new atom;
-        valid = true;
-        do {
-          if (!valid) DoeLog(2) && (eLog()<< Verbose(2) << "Resulting position out of cell." << endl);
-          second = mol->AskAtom("Enter atom number: ");
-          DoLog(0) && (Log() << Verbose(0) << "Enter relative coordinates." << endl);
-          first->x.AskPosition(World::getInstance().getDomain(), false);
-          for (int i=NDIM;i--;) {
-            first->x.x[i] += second->x.x[i];
-          }
-        } while (!(valid = mol->CheckBounds((const Vector *)&first->x)));
-        first->type = periode->AskElement();  // give type
-        mol->AddAtom(first);  // add to molecule
-        break;
-
-    case 'd': // two atoms, two angles and a distance
-        first = new atom;
-        valid = true;
-        do {
-          if (!valid) {
-            DoeLog(2) && (eLog()<< Verbose(2) << "Resulting coordinates out of cell - " << first->x << endl);
-          }
-          cout << Verbose(0) << "First, we need two atoms, the first atom is the central, while the second is the outer one." << endl;
-          second = mol->AskAtom("Enter central atom: ");
-          third = mol->AskAtom("Enter second atom (specifying the axis for first angle): ");
-          fourth = mol->AskAtom("Enter third atom (specifying a plane for second angle): ");
-          a = ask_value("Enter distance between central (first) and new atom: ");
-          b = ask_value("Enter angle between new, first and second atom (degrees): ");
-          b *= M_PI/180.;
-          bound(&b, 0., 2.*M_PI);
-          c = ask_value("Enter second angle between new and normal vector of plane defined by first, second and third atom (degrees): ");
-          c *= M_PI/180.;
-          bound(&c, -M_PI, M_PI);
-          cout << Verbose(0) << "radius: " << a << "\t phi: " << b*180./M_PI << "\t theta: " << c*180./M_PI << endl;
-/*
-          second->Output(1,1,(ofstream *)&cout);
-          third->Output(1,2,(ofstream *)&cout);
-          fourth->Output(1,3,(ofstream *)&cout);
-          n.MakeNormalvector((const vector *)&second->x, (const vector *)&third->x, (const vector *)&fourth->x);
-          x.Copyvector(&second->x);
-          x.SubtractVector(&third->x);
-          x.Copyvector(&fourth->x);
-          x.SubtractVector(&third->x);
-
-          if (!z.SolveSystem(&x,&y,&n, b, c, a)) {
-         coutg() << Verbose(0) << "Failure solving self-dependent linear system!" << endl;
-            continue;
-          }
-          DoLog(0) && (Log() << Verbose(0) << "resulting relative coordinates: ");
-          z.Output();
-          DoLog(0) && (Log() << Verbose(0) << endl);
-          */
-          // calc axis vector
-          x.CopyVector(&second->x);
-          x.SubtractVector(&third->x);
-          x.Normalize();
-          Log() << Verbose(0) << "x: ",
-          x.Output();
-          DoLog(0) && (Log() << Verbose(0) << endl);
-          z.MakeNormalVector(&second->x,&third->x,&fourth->x);
-          Log() << Verbose(0) << "z: ",
-          z.Output();
-          DoLog(0) && (Log() << Verbose(0) << endl);
-          y.MakeNormalVector(&x,&z);
-          Log() << Verbose(0) << "y: ",
-          y.Output();
-          DoLog(0) && (Log() << Verbose(0) << endl);
-
-          // rotate vector around first angle
-          first->x.CopyVector(&x);
-          first->x.RotateVector(&z,b - M_PI);
-          Log() << Verbose(0) << "Rotated vector: ",
-          first->x.Output();
-          DoLog(0) && (Log() << Verbose(0) << endl);
-          // remove the projection onto the rotation plane of the second angle
-          n.CopyVector(&y);
-          n.Scale(first->x.ScalarProduct(&y));
-          Log() << Verbose(0) << "N1: ",
-          n.Output();
-          DoLog(0) && (Log() << Verbose(0) << endl);
-          first->x.SubtractVector(&n);
-          Log() << Verbose(0) << "Subtracted vector: ",
-          first->x.Output();
-          DoLog(0) && (Log() << Verbose(0) << endl);
-          n.CopyVector(&z);
-          n.Scale(first->x.ScalarProduct(&z));
-          Log() << Verbose(0) << "N2: ",
-          n.Output();
-          DoLog(0) && (Log() << Verbose(0) << endl);
-          first->x.SubtractVector(&n);
-          Log() << Verbose(0) << "2nd subtracted vector: ",
-          first->x.Output();
-          DoLog(0) && (Log() << Verbose(0) << endl);
-
-          // rotate another vector around second angle
-          n.CopyVector(&y);
-          n.RotateVector(&x,c - M_PI);
-          Log() << Verbose(0) << "2nd Rotated vector: ",
-          n.Output();
-          DoLog(0) && (Log() << Verbose(0) << endl);
-
-          // add the two linear independent vectors
-          first->x.AddVector(&n);
-          first->x.Normalize();
-          first->x.Scale(a);
-          first->x.AddVector(&second->x);
-
-          DoLog(0) && (Log() << Verbose(0) << "resulting coordinates: ");
-          first->x.Output();
-          DoLog(0) && (Log() << Verbose(0) << endl);
-        } while (!(valid = mol->CheckBounds((const Vector *)&first->x)));
-        first->type = periode->AskElement();  // give type
-        mol->AddAtom(first);  // add to molecule
-        break;
-
-      case 'e': // least square distance position to a set of atoms
-        first = new atom;
-        atoms = new (Vector*[128]);
-        valid = true;
-        for(int i=128;i--;)
-          atoms[i] = NULL;
-        int i=0, j=0;
-        cout << Verbose(0) << "Now we need at least three molecules.\n";
-        do {
-          cout << Verbose(0) << "Enter " << i+1 << "th atom: ";
-          cin >> j;
-          if (j != -1) {
-            second = mol->FindAtom(j);
-            atoms[i++] = &(second->x);
-          }
-        } while ((j != -1) && (i<128));
-        if (i >= 2) {
-          first->x.LSQdistance((const Vector **)atoms, i);
-          first->x.Output();
-          first->type = periode->AskElement();  // give type
-          mol->AddAtom(first);  // add to molecule
-        } else {
-          delete first;
-          cout << Verbose(0) << "Please enter at least two vectors!\n";
-        }
-        break;
-  };
-};
-
-/** Submenu for centering the atoms in the molecule.
- * \param *mol molecule with all the atoms
- */
-static void CenterAtoms(molecule *mol)
-{
-  Vector x, y, helper;
-  char choice;  // menu choice char
-
-  cout << Verbose(0) << "===========CENTER ATOMS=========================" << endl;
-  cout << Verbose(0) << " a - on origin" << endl;
-  cout << Verbose(0) << " b - on center of gravity" << endl;
-  cout << Verbose(0) << " c - within box with additional boundary" << endl;
-  cout << Verbose(0) << " d - within given simulation box" << endl;
-  cout << Verbose(0) << "all else - go back" << endl;
-  cout << Verbose(0) << "===============================================" << endl;
-  cout << Verbose(0) << "INPUT: ";
-  cin >> choice;
-
-  switch (choice) {
-    default:
-      cout << Verbose(0) << "Not a valid choice." << endl;
-      break;
-    case 'a':
-      cout << Verbose(0) << "Centering atoms in config file on origin." << endl;
-      mol->CenterOrigin();
-      break;
-    case 'b':
-      cout << Verbose(0) << "Centering atoms in config file on center of gravity." << endl;
-      mol->CenterPeriodic();
-      break;
-    case 'c':
-      cout << Verbose(0) << "Centering atoms in config file within given additional boundary." << endl;
-      for (int i=0;i<NDIM;i++) {
-        cout << Verbose(0) << "Enter axis " << i << " boundary: ";
-        cin >> y.x[i];
-      }
-      mol->CenterEdge(&x);  // make every coordinate positive
-      mol->Center.AddVector(&y); // translate by boundary
-      helper.CopyVector(&y);
-      helper.Scale(2.);
-      helper.AddVector(&x);
-      mol->SetBoxDimension(&helper);  // update Box of atoms by boundary
-      break;
-    case 'd':
-      cout << Verbose(1) << "Centering atoms in config file within given simulation box." << endl;
-      for (int i=0;i<NDIM;i++) {
-        cout << Verbose(0) << "Enter axis " << i << " boundary: ";
-        cin >> x.x[i];
-      }
-      // update Box of atoms by boundary
-      mol->SetBoxDimension(&x);
-      // center
-      mol->CenterInBox();
-      break;
-  }
-};
-
-/** Submenu for aligning the atoms in the molecule.
- * \param *periode periodentafel
- * \param *mol molecule with all the atoms
- */
-static void AlignAtoms(periodentafel *periode, molecule *mol)
-{
-  atom *first, *second, *third;
-  Vector x,n;
-  char choice;  // menu choice char
-
-  cout << Verbose(0) << "===========ALIGN ATOMS=========================" << endl;
-  cout << Verbose(0) << " a - state three atoms defining align plane" << endl;
-  cout << Verbose(0) << " b - state alignment vector" << endl;
-  cout << Verbose(0) << " c - state two atoms in alignment direction" << endl;
-  cout << Verbose(0) << " d - align automatically by least square fit" << endl;
-  cout << Verbose(0) << "all else - go back" << endl;
-  cout << Verbose(0) << "===============================================" << endl;
-  cout << Verbose(0) << "INPUT: ";
-  cin >> choice;
-
-  switch (choice) {
-    default:
-    case 'a': // three atoms defining mirror plane
-      first = mol->AskAtom("Enter first atom: ");
-      second = mol->AskAtom("Enter second atom: ");
-      third = mol->AskAtom("Enter third atom: ");
-
-      n.MakeNormalVector((const Vector *)&first->x,(const Vector *)&second->x,(const Vector *)&third->x);
-      break;
-    case 'b': // normal vector of mirror plane
-      cout << Verbose(0) << "Enter normal vector of mirror plane." << endl;
-      n.AskPosition(World::getInstance().getDomain(),0);
-      n.Normalize();
-      break;
-    case 'c': // three atoms defining mirror plane
-      first = mol->AskAtom("Enter first atom: ");
-      second = mol->AskAtom("Enter second atom: ");
-
-      n.CopyVector((const Vector *)&first->x);
-      n.SubtractVector((const Vector *)&second->x);
-      n.Normalize();
-      break;
-    case 'd':
-      char shorthand[4];
-      Vector a;
-      struct lsq_params param;
-      do {
-        fprintf(stdout, "Enter the element of atoms to be chosen: ");
-        fscanf(stdin, "%3s", shorthand);
-      } while ((param.type = periode->FindElement(shorthand)) == NULL);
-      cout << Verbose(0) << "Element is " << param.type->name << endl;
-      mol->GetAlignvector(&param);
-      for (int i=NDIM;i--;) {
-        x.x[i] = gsl_vector_get(param.x,i);
-        n.x[i] = gsl_vector_get(param.x,i+NDIM);
-      }
-      gsl_vector_free(param.x);
-      cout << Verbose(0) << "Offset vector: ";
-      x.Output();
-      DoLog(0) && (Log() << Verbose(0) << endl);
-      n.Normalize();
-      break;
-  };
-  DoLog(0) && (Log() << Verbose(0) << "Alignment vector: ");
-  n.Output();
-  DoLog(0) && (Log() << Verbose(0) << endl);
-  mol->Align(&n);
-};
-
-/** Submenu for mirroring the atoms in the molecule.
- * \param *mol molecule with all the atoms
- */
-static void MirrorAtoms(molecule *mol)
-{
-  atom *first, *second, *third;
-  Vector n;
-  char choice;  // menu choice char
-
-  DoLog(0) && (Log() << Verbose(0) << "===========MIRROR ATOMS=========================" << endl);
-  DoLog(0) && (Log() << Verbose(0) << " a - state three atoms defining mirror plane" << endl);
-  DoLog(0) && (Log() << Verbose(0) << " b - state normal vector of mirror plane" << endl);
-  DoLog(0) && (Log() << Verbose(0) << " c - state two atoms in normal direction" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "all else - go back" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "===============================================" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "INPUT: ");
-  cin >> choice;
-
-  switch (choice) {
-    default:
-    case 'a': // three atoms defining mirror plane
-      first = mol->AskAtom("Enter first atom: ");
-      second = mol->AskAtom("Enter second atom: ");
-      third = mol->AskAtom("Enter third atom: ");
-
-      n.MakeNormalVector((const Vector *)&first->x,(const Vector *)&second->x,(const Vector *)&third->x);
-      break;
-    case 'b': // normal vector of mirror plane
-      DoLog(0) && (Log() << Verbose(0) << "Enter normal vector of mirror plane." << endl);
-      n.AskPosition(World::getInstance().getDomain(),0);
-      n.Normalize();
-      break;
-    case 'c': // three atoms defining mirror plane
-      first = mol->AskAtom("Enter first atom: ");
-      second = mol->AskAtom("Enter second atom: ");
-
-      n.CopyVector((const Vector *)&first->x);
-      n.SubtractVector((const Vector *)&second->x);
-      n.Normalize();
-      break;
-  };
-  DoLog(0) && (Log() << Verbose(0) << "Normal vector: ");
-  n.Output();
-  DoLog(0) && (Log() << Verbose(0) << endl);
-  mol->Mirror((const Vector *)&n);
-};
-
-/** Submenu for removing the atoms from the molecule.
- * \param *mol molecule with all the atoms
- */
-static void RemoveAtoms(molecule *mol)
-{
-  atom *first, *second;
-  int axis;
-  double tmp1, tmp2;
-  char choice;  // menu choice char
-
-  DoLog(0) && (Log() << Verbose(0) << "===========REMOVE ATOMS=========================" << endl);
-  DoLog(0) && (Log() << Verbose(0) << " a - state atom for removal by number" << endl);
-  DoLog(0) && (Log() << Verbose(0) << " b - keep only in radius around atom" << endl);
-  DoLog(0) && (Log() << Verbose(0) << " c - remove this with one axis greater value" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "all else - go back" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "===============================================" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "INPUT: ");
-  cin >> choice;
-
-  switch (choice) {
-    default:
-    case 'a':
-      if (mol->RemoveAtom(mol->AskAtom("Enter number of atom within molecule: ")))
-        DoLog(1) && (Log() << Verbose(1) << "Atom removed." << endl);
-      else
-        DoLog(1) && (Log() << Verbose(1) << "Atom not found." << endl);
-      break;
-    case 'b':
-      second = mol->AskAtom("Enter number of atom as reference point: ");
-      DoLog(0) && (Log() << Verbose(0) << "Enter radius: ");
-      cin >> tmp1;
-      first = mol->start;
-      second = first->next;
-      while(second != mol->end) {
-        first = second;
-        second = first->next;
-        if (first->x.DistanceSquared((const Vector *)&second->x) > tmp1*tmp1) // distance to first above radius ...
-          mol->RemoveAtom(first);
-      }
-      break;
-    case 'c':
-      DoLog(0) && (Log() << Verbose(0) << "Which axis is it: ");
-      cin >> axis;
-      DoLog(0) && (Log() << Verbose(0) << "Lower boundary: ");
-      cin >> tmp1;
-      DoLog(0) && (Log() << Verbose(0) << "Upper boundary: ");
-      cin >> tmp2;
-      first = mol->start;
-      second = first->next;
-      while(second != mol->end) {
-        first = second;
-        second = first->next;
-        if ((first->x.x[axis] < tmp1) || (first->x.x[axis] > tmp2)) {// out of boundary ...
-          //Log() << Verbose(0) << "Atom " << *first << " with " << first->x.x[axis] << " on axis " << axis << " is out of bounds [" << tmp1 << "," << tmp2 << "]." << endl;
-          mol->RemoveAtom(first);
-        }
-      }
-      break;
-  };
-  //mol->Output();
-  choice = 'r';
-};
-
-/** Submenu for measuring out the atoms in the molecule.
- * \param *periode periodentafel
- * \param *mol molecule with all the atoms
- */
-static void MeasureAtoms(periodentafel *periode, molecule *mol, config *configuration)
-{
-  atom *first, *second, *third;
-  Vector x,y;
-  double min[256], tmp1, tmp2, tmp3;
-  int Z;
-  char choice;  // menu choice char
-
-  DoLog(0) && (Log() << Verbose(0) << "===========MEASURE ATOMS=========================" << endl);
-  DoLog(0) && (Log() << Verbose(0) << " a - calculate bond length between one atom and all others" << endl);
-  DoLog(0) && (Log() << Verbose(0) << " b - calculate bond length between two atoms" << endl);
-  DoLog(0) && (Log() << Verbose(0) << " c - calculate bond angle" << endl);
-  DoLog(0) && (Log() << Verbose(0) << " d - calculate principal axis of the system" << endl);
-  DoLog(0) && (Log() << Verbose(0) << " e - calculate volume of the convex envelope" << endl);
-  DoLog(0) && (Log() << Verbose(0) << " f - calculate temperature from current velocity" << endl);
-  DoLog(0) && (Log() << Verbose(0) << " g - output all temperatures per step from velocities" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "all else - go back" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "===============================================" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "INPUT: ");
-  cin >> choice;
-
-  switch(choice) {
-    default:
-      DoLog(1) && (Log() << Verbose(1) << "Not a valid choice." << endl);
-      break;
-    case 'a':
-      first = mol->AskAtom("Enter first atom: ");
-      for (int i=MAX_ELEMENTS;i--;)
-        min[i] = 0.;
-
-      second = mol->start;
-      while ((second->next != mol->end)) {
-        second = second->next; // advance
-        Z = second->type->Z;
-        tmp1 = 0.;
-        if (first != second) {
-          x.CopyVector((const Vector *)&first->x);
-          x.SubtractVector((const Vector *)&second->x);
-          tmp1 = x.Norm();
-        }
-        if ((tmp1 != 0.) && ((min[Z] == 0.) || (tmp1 < min[Z]))) min[Z] = tmp1;
-        //Log() << Verbose(0) << "Bond length between Atom " << first->nr << " and " << second->nr << ": " << tmp1 << " a.u." << endl;
-      }
-      for (int i=MAX_ELEMENTS;i--;)
-        if (min[i] != 0.) Log() << Verbose(0) << "Minimum Bond length between " << first->type->name << " Atom " << first->nr << " and next Ion of type " << (periode->FindElement(i))->name << ": " << min[i] << " a.u." << endl;
-      break;
-
-    case 'b':
-      first = mol->AskAtom("Enter first atom: ");
-      second = mol->AskAtom("Enter second atom: ");
-      for (int i=NDIM;i--;)
-        min[i] = 0.;
-      x.CopyVector((const Vector *)&first->x);
-      x.SubtractVector((const Vector *)&second->x);
-      tmp1 = x.Norm();
-      DoLog(1) && (Log() << Verbose(1) << "Distance vector is ");
-      x.Output();
-      DoLog(0) && (Log() << Verbose(0) << "." << endl << "Norm of distance is " << tmp1 << "." << endl);
-      break;
-
-    case 'c':
-      DoLog(0) && (Log() << Verbose(0) << "Evaluating bond angle between three - first, central, last - atoms." << endl);
-      first = mol->AskAtom("Enter first atom: ");
-      second = mol->AskAtom("Enter central atom: ");
-      third  = mol->AskAtom("Enter last atom: ");
-      tmp1 = tmp2 = tmp3 = 0.;
-      x.CopyVector((const Vector *)&first->x);
-      x.SubtractVector((const Vector *)&second->x);
-      y.CopyVector((const Vector *)&third->x);
-      y.SubtractVector((const Vector *)&second->x);
-      DoLog(0) && (Log() << Verbose(0) << "Bond angle between first atom Nr." << first->nr << ", central atom Nr." << second->nr << " and last atom Nr." << third->nr << ": ");
-      DoLog(0) && (Log() << Verbose(0) << (acos(x.ScalarProduct((const Vector *)&y)/(y.Norm()*x.Norm()))/M_PI*180.) << " degrees" << endl);
-      break;
-    case 'd':
-      DoLog(0) && (Log() << Verbose(0) << "Evaluating prinicipal axis." << endl);
-      DoLog(0) && (Log() << Verbose(0) << "Shall we rotate? [0/1]: ");
-      cin >> Z;
-      if ((Z >=0) && (Z <=1))
-        mol->PrincipalAxisSystem((bool)Z);
-      else
-        mol->PrincipalAxisSystem(false);
-      break;
-    case 'e':
-      {
-        DoLog(0) && (Log() << Verbose(0) << "Evaluating volume of the convex envelope.");
-        class Tesselation *TesselStruct = NULL;
-        const LinkedCell *LCList = NULL;
-        LCList = new LinkedCell(mol, 10.);
-        FindConvexBorder(mol, TesselStruct, LCList, NULL);
-        double clustervolume = VolumeOfConvexEnvelope(TesselStruct, configuration);
-        DoLog(0) && (Log() << Verbose(0) << "The tesselated surface area is " << clustervolume << "." << endl);\
-        delete(LCList);
-        delete(TesselStruct);
-      }
-      break;
-    case 'f':
-      mol->OutputTemperatureFromTrajectories((ofstream *)&cout, mol->MDSteps-1, mol->MDSteps);
-      break;
-    case 'g':
-      {
-        char filename[255];
-        DoLog(0) && (Log() << Verbose(0) << "Please enter filename: " << endl);
-        cin >> filename;
-        DoLog(1) && (Log() << Verbose(1) << "Storing temperatures in " << filename << "." << endl);
-        ofstream *output = new ofstream(filename, ios::trunc);
-        if (!mol->OutputTemperatureFromTrajectories(output, 0, mol->MDSteps))
-          DoLog(2) && (Log() << Verbose(2) << "File could not be written." << endl);
-        else
-          DoLog(2) && (Log() << Verbose(2) << "File stored." << endl);
-        output->close();
-        delete(output);
-      }
-      break;
-  }
-};
-
-/** Submenu for measuring out the atoms in the molecule.
- * \param *mol molecule with all the atoms
- * \param *configuration configuration structure for the to be written config files of all fragments
- */
-static void FragmentAtoms(molecule *mol, config *configuration)
-{
-  int Order1;
-  clock_t start, end;
-
-  DoLog(0) && (Log() << Verbose(0) << "Fragmenting molecule with current connection matrix ..." << endl);
-  DoLog(0) && (Log() << Verbose(0) << "What's the desired bond order: ");
-  cin >> Order1;
-  if (mol->first->next != mol->last) {  // there are bonds
-    start = clock();
-    mol->FragmentMolecule(Order1, configuration);
-    end = clock();
-    DoLog(0) && (Log() << Verbose(0) << "Clocks for this operation: " << (end-start) << ", time: " << ((double)(end-start)/CLOCKS_PER_SEC) << "s." << endl);
-  } else
-    DoLog(0) && (Log() << Verbose(0) << "Connection matrix has not yet been generated!" << endl);
-};
-
-/********************************************** Submenu routine **************************************/
-
-/** Submenu for manipulating atoms.
- * \param *periode periodentafel
- * \param *molecules list of molecules whose atoms are to be manipulated
- */
-static void ManipulateAtoms(periodentafel *periode, MoleculeListClass *molecules, config *configuration)
-{
-  atom *first, *second, *third;
-  molecule *mol = NULL;
-  Vector x,y,z,n; // coordinates for absolute point in cell volume
-  double *factor; // unit factor if desired
-  double bond, minBond;
-  char choice;  // menu choice char
-  bool valid;
-
-  DoLog(0) && (Log() << Verbose(0) << "=========MANIPULATE ATOMS======================" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "a - add an atom" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "r - remove an atom" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "b - scale a bond between atoms" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "t - turn an atom round another bond" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "u - change an atoms element" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "l - measure lengths, angles, ... for an atom" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "all else - go back" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "===============================================" << endl);
-  if (molecules->NumberOfActiveMolecules() > 1)
-    DoeLog(2) && (eLog()<< Verbose(2) << "There is more than one molecule active! Atoms will be added to each." << endl);
-  DoLog(0) && (Log() << Verbose(0) << "INPUT: ");
-  cin >> choice;
-
-  switch (choice) {
-    default:
-      DoLog(0) && (Log() << Verbose(0) << "Not a valid choice." << endl);
-      break;
-
-    case 'a': // add atom
-      for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++)
-        if ((*ListRunner)->ActiveFlag) {
-        mol = *ListRunner;
-        DoLog(0) && (Log() << Verbose(0) << "Current molecule is: " << mol->IndexNr << "\t" << mol->name << endl);
-        AddAtoms(periode, mol);
-      }
-      break;
-
-    case 'b': // scale a bond
-      for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++)
-        if ((*ListRunner)->ActiveFlag) {
-        mol = *ListRunner;
-        DoLog(0) && (Log() << Verbose(0) << "Current molecule is: " << mol->IndexNr << "\t" << mol->name << endl);
-        DoLog(0) && (Log() << Verbose(0) << "Scaling bond length between two atoms." << endl);
-        first = mol->AskAtom("Enter first (fixed) atom: ");
-        second = mol->AskAtom("Enter second (shifting) atom: ");
-        minBond = 0.;
-        for (int i=NDIM;i--;)
-          minBond += (first->x.x[i]-second->x.x[i])*(first->x.x[i] - second->x.x[i]);
-        minBond = sqrt(minBond);
-        DoLog(0) && (Log() << Verbose(0) << "Current Bond length between " << first->type->name << " Atom " << first->nr << " and " << second->type->name << " Atom " << second->nr << ": " << minBond << " a.u." << endl);
-        DoLog(0) && (Log() << Verbose(0) << "Enter new bond length [a.u.]: ");
-        cin >> bond;
-        for (int i=NDIM;i--;) {
-          second->x.x[i] -= (second->x.x[i]-first->x.x[i])/minBond*(minBond-bond);
-        }
-        //Log() << Verbose(0) << "New coordinates of Atom " << second->nr << " are: ";
-        //second->Output(second->type->No, 1);
-      }
-      break;
-
-    case 'c': // unit scaling of the metric
-      for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++)
-        if ((*ListRunner)->ActiveFlag) {
-        mol = *ListRunner;
-        DoLog(0) && (Log() << Verbose(0) << "Current molecule is: " << mol->IndexNr << "\t" << mol->name << endl);
-       DoLog(0) && (Log() << Verbose(0) << "Angstroem -> Bohrradius: 1.8897261\t\tBohrradius -> Angstroem: 0.52917721" << endl);
-       DoLog(0) && (Log() << Verbose(0) << "Enter three factors: ");
-       factor = new double[NDIM];
-       cin >> factor[0];
-       cin >> factor[1];
-       cin >> factor[2];
-       valid = true;
-       mol->Scale((const double ** const)&factor);
-       delete[](factor);
-      }
-     break;
-
-    case 'l': // measure distances or angles
-      for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++)
-        if ((*ListRunner)->ActiveFlag) {
-        mol = *ListRunner;
-        DoLog(0) && (Log() << Verbose(0) << "Current molecule is: " << mol->IndexNr << "\t" << mol->name << endl);
-        MeasureAtoms(periode, mol, configuration);
-      }
-      break;
-
-    case 'r': // remove atom
-      for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++)
-        if ((*ListRunner)->ActiveFlag) {
-        mol = *ListRunner;
-        DoLog(0) && (Log() << Verbose(0) << "Current molecule is: " << mol->IndexNr << "\t" << mol->name << endl);
-        RemoveAtoms(mol);
-      }
-      break;
-
-    case 't': // turn/rotate atom
-      for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++)
-        if ((*ListRunner)->ActiveFlag) {
-          mol = *ListRunner;
-          DoLog(0) && (Log() << Verbose(0) << "Turning atom around another bond - first is atom to turn, second (central) and third specify bond" << endl);
-          first = mol->AskAtom("Enter turning atom: ");
-          second = mol->AskAtom("Enter central atom: ");
-          third  = mol->AskAtom("Enter bond atom: ");
-          cout << Verbose(0) << "Enter new angle in degrees: ";
-          double tmp = 0.;
-          cin >> tmp;
-          // calculate old angle
-          x.CopyVector((const Vector *)&first->x);
-          x.SubtractVector((const Vector *)&second->x);
-          y.CopyVector((const Vector *)&third->x);
-          y.SubtractVector((const Vector *)&second->x);
-          double alpha = (acos(x.ScalarProduct((const Vector *)&y)/(y.Norm()*x.Norm()))/M_PI*180.);
-          cout << Verbose(0) << "Bond angle between first atom Nr." << first->nr << ", central atom Nr." << second->nr << " and last atom Nr." << third->nr << ": ";
-          cout << Verbose(0) << alpha << " degrees" << endl;
-          // rotate
-          z.MakeNormalVector(&x,&y);
-          x.RotateVector(&z,(alpha-tmp)*M_PI/180.);
-          x.AddVector(&second->x);
-          first->x.CopyVector(&x);
-          // check new angle
-          x.CopyVector((const Vector *)&first->x);
-          x.SubtractVector((const Vector *)&second->x);
-          alpha = (acos(x.ScalarProduct((const Vector *)&y)/(y.Norm()*x.Norm()))/M_PI*180.);
-          cout << Verbose(0) << "new Bond angle between first atom Nr." << first->nr << ", central atom Nr." << second->nr << " and last atom Nr." << third->nr << ": ";
-          cout << Verbose(0) << alpha << " degrees" << endl;
-        }
-      break;
-
-    case 'u': // change an atom's element
-      for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++)
-        if ((*ListRunner)->ActiveFlag) {
-        int Z;
-        mol = *ListRunner;
-        DoLog(0) && (Log() << Verbose(0) << "Current molecule is: " << mol->IndexNr << "\t" << mol->name << endl);
-        first = NULL;
-        do {
-          DoLog(0) && (Log() << Verbose(0) << "Change the element of which atom: ");
-          cin >> Z;
-        } while ((first = mol->FindAtom(Z)) == NULL);
-        DoLog(0) && (Log() << Verbose(0) << "New element by atomic number Z: ");
-        cin >> Z;
-        first->type = periode->FindElement(Z);
-        DoLog(0) && (Log() << Verbose(0) << "Atom " << first->nr << "'s element is " << first->type->name << "." << endl);
-      }
-      break;
-  }
-};
-
-/** Submenu for manipulating molecules.
- * \param *periode periodentafel
- * \param *molecules list of molecule to manipulate
- */
-static void ManipulateMolecules(periodentafel *periode, MoleculeListClass *molecules, config *configuration)
-{
-  atom *first = NULL;
-  Vector x,y,z,n; // coordinates for absolute point in cell volume
-  int j, axis, count, faktor;
-  char choice;  // menu choice char
-  molecule *mol = NULL;
-  element **Elements;
-  Vector **vectors;
-  MoleculeLeafClass *Subgraphs = NULL;
-
-  DoLog(0) && (Log() << Verbose(0) << "=========MANIPULATE GLOBALLY===================" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "c - scale by unit transformation" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "d - duplicate molecule/periodic cell" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "f - fragment molecule many-body bond order style" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "g - center atoms in box" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "i - realign molecule" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "m - mirror all molecules" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "o - create connection matrix" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "t - translate molecule by vector" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "all else - go back" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "===============================================" << endl);
-  if (molecules->NumberOfActiveMolecules() > 1)
-    DoeLog(2) && (eLog()<< Verbose(2) << "There is more than one molecule active! Atoms will be added to each." << endl);
-  DoLog(0) && (Log() << Verbose(0) << "INPUT: ");
-  cin >> choice;
-
-  switch (choice) {
-    default:
-      DoLog(0) && (Log() << Verbose(0) << "Not a valid choice." << endl);
-      break;
-
-    case 'd': // duplicate the periodic cell along a given axis, given times
-      for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++) 
-        if ((*ListRunner)->ActiveFlag) {
-        mol = *ListRunner;
-        DoLog(0) && (Log() << Verbose(0) << "Current molecule is: " << mol->IndexNr << "\t" << mol->name << endl);
-        DoLog(0) && (Log() << Verbose(0) << "State the axis [(+-)123]: ");
-        cin >> axis;
-        DoLog(0) && (Log() << Verbose(0) << "State the factor: ");
-        cin >> faktor;
-
-        mol->CountAtoms(); // recount atoms
-        if (mol->getAtomCount() != 0) {  // if there is more than none
-          count = mol->getAtomCount();  // is changed becausing of adding, thus has to be stored away beforehand
-          Elements = new element *[count];
-          vectors = new Vector *[count];
-          j = 0;
-          first = mol->start;
-          while (first->next != mol->end) { // make a list of all atoms with coordinates and element
-            first = first->next;
-            Elements[j] = first->type;
-            vectors[j] = &first->x;
-            j++;
-          }
-          if (count != j)
-            DoeLog(1) && (eLog()<< Verbose(1) << "AtomCount " << count << " is not equal to number of atoms in molecule " << j << "!" << endl);
-          x.Zero();
-          y.Zero();
-          y.x[abs(axis)-1] = World::getInstance().getDomain()[(abs(axis) == 2) ? 2 : ((abs(axis) == 3) ? 5 : 0)] * abs(axis)/axis; // last term is for sign, first is for magnitude
-          for (int i=1;i<faktor;i++) {  // then add this list with respective translation factor times
-            x.AddVector(&y); // per factor one cell width further
-            for (int k=count;k--;) { // go through every atom of the original cell
-              first = new atom(); // create a new body
-              first->x.CopyVector(vectors[k]);  // use coordinate of original atom
-              first->x.AddVector(&x);     // translate the coordinates
-              first->type = Elements[k];  // insert original element
-              mol->AddAtom(first);        // and add to the molecule (which increments ElementsInMolecule, AtomCount, ...)
-            }
-          }
-          if (mol->first->next != mol->last) // if connect matrix is present already, redo it
-            mol->CreateAdjacencyList(mol->BondDistance, configuration->GetIsAngstroem(), &BondGraph::CovalentMinMaxDistance, NULL);
-          // free memory
-          delete[](Elements);
-          delete[](vectors);
-          // correct cell size
-          if (axis < 0) { // if sign was negative, we have to translate everything
-            x.Zero();
-            x.AddVector(&y);
-            x.Scale(-(faktor-1));
-            mol->Translate(&x);
-          }
-          World::getInstance().getDomain()[(abs(axis) == 2) ? 2 : ((abs(axis) == 3) ? 5 : 0)] *= faktor;
-        }
-      }
-      break;
-
-    case 'f':
-      FragmentAtoms(mol, configuration);
-      break;
-
-    case 'g': // center the atoms
-      for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++)
-        if ((*ListRunner)->ActiveFlag) {
-        mol = *ListRunner;
-        DoLog(0) && (Log() << Verbose(0) << "Current molecule is: " << mol->IndexNr << "\t" << mol->name << endl);
-        CenterAtoms(mol);
-      }
-      break;
-
-    case 'i': // align all atoms
-      for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++)
-        if ((*ListRunner)->ActiveFlag) {
-        mol = *ListRunner;
-        DoLog(0) && (Log() << Verbose(0) << "Current molecule is: " << mol->IndexNr << "\t" << mol->name << endl);
-        AlignAtoms(periode, mol);
-      }
-      break;
-
-    case 'm': // mirror atoms along a given axis
-      for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++)
-        if ((*ListRunner)->ActiveFlag) {
-        mol = *ListRunner;
-        DoLog(0) && (Log() << Verbose(0) << "Current molecule is: " << mol->IndexNr << "\t" << mol->name << endl);
-        MirrorAtoms(mol);
-      }
-      break;
-
-    case 'o': // create the connection matrix
-      for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++)
-        if ((*ListRunner)->ActiveFlag) {
-          mol = *ListRunner;
-          double bonddistance;
-          clock_t start,end;
-          DoLog(0) && (Log() << Verbose(0) << "What's the maximum bond distance: ");
-          cin >> bonddistance;
-          start = clock();
-          mol->CreateAdjacencyList(bonddistance, configuration->GetIsAngstroem(), &BondGraph::CovalentMinMaxDistance, NULL);
-          end = clock();
-          DoLog(0) && (Log() << Verbose(0) << "Clocks for this operation: " << (end-start) << ", time: " << ((double)(end-start)/CLOCKS_PER_SEC) << "s." << endl);
-        }
-      break;
-
-    case 't': // translate all atoms
-      for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++)
-        if ((*ListRunner)->ActiveFlag) {
-        mol = *ListRunner;
-        DoLog(0) && (Log() << Verbose(0) << "Current molecule is: " << mol->IndexNr << "\t" << mol->name << endl);
-        DoLog(0) && (Log() << Verbose(0) << "Enter translation vector." << endl);
-        x.AskPosition(World::getInstance().getDomain(),0);
-        mol->Center.AddVector((const Vector *)&x);
-     }
-     break;
-  }
-  // Free all
-  if (Subgraphs != NULL) {  // free disconnected subgraph list of DFS analysis was performed
-    while (Subgraphs->next != NULL) {
-      Subgraphs = Subgraphs->next;
-      delete(Subgraphs->previous);
-    }
-    delete(Subgraphs);
-  }
-};
-
-
-/** Submenu for creating new molecules.
- * \param *periode periodentafel
- * \param *molecules list of molecules to add to
- */
-static void EditMolecules(periodentafel *periode, MoleculeListClass *molecules)
-{
-  char choice;  // menu choice char
-  Vector center;
-  int nr, count;
-  molecule *mol = NULL;
-
-  DoLog(0) && (Log() << Verbose(0) << "==========EDIT MOLECULES=====================" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "c - create new molecule" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "l - load molecule from xyz file" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "n - change molecule's name" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "N - give molecules filename" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "p - parse atoms in xyz file into molecule" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "r - remove a molecule" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "all else - go back" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "===============================================" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "INPUT: ");
-  cin >> choice;
-
-  switch (choice) {
-    default:
-      DoLog(0) && (Log() << Verbose(0) << "Not a valid choice." << endl);
-      break;
-    case 'c':
-      mol = World::getInstance().createMolecule();
-      molecules->insert(mol);
-      break;
-
-    case 'l': // load from XYZ file
-      {
-        char filename[MAXSTRINGSIZE];
-        DoLog(0) && (Log() << Verbose(0) << "Format should be XYZ with: ShorthandOfElement\tX\tY\tZ" << endl);
-        mol = World::getInstance().createMolecule();
-        do {
-          DoLog(0) && (Log() << Verbose(0) << "Enter file name: ");
-          cin >> filename;
-        } while (!mol->AddXYZFile(filename));
-        mol->SetNameFromFilename(filename);
-        // center at set box dimensions
-        mol->CenterEdge(&center);
-        double * const cell_size = World::getInstance().getDomain();
-        cell_size[0] = center.x[0];
-        cell_size[1] = 0;
-        cell_size[2] = center.x[1];
-        cell_size[3] = 0;
-        cell_size[4] = 0;
-        cell_size[5] = center.x[2];
-        molecules->insert(mol);
-      }
-      break;
-
-    case 'n':
-      {
-        char filename[MAXSTRINGSIZE];
-        do {
-          DoLog(0) && (Log() << Verbose(0) << "Enter index of molecule: ");
-          cin >> nr;
-          mol = molecules->ReturnIndex(nr);
-        } while (mol == NULL);
-        DoLog(0) && (Log() << Verbose(0) << "Enter name: ");
-        cin >> filename;
-        strcpy(mol->name, filename);
-      }
-      break;
-
-    case 'N':
-      {
-        char filename[MAXSTRINGSIZE];
-        do {
-          DoLog(0) && (Log() << Verbose(0) << "Enter index of molecule: ");
-          cin >> nr;
-          mol = molecules->ReturnIndex(nr);
-        } while (mol == NULL);
-        DoLog(0) && (Log() << Verbose(0) << "Enter name: ");
-        cin >> filename;
-        mol->SetNameFromFilename(filename);
-      }
-      break;
-
-    case 'p': // parse XYZ file
-      {
-        char filename[MAXSTRINGSIZE];
-        mol = NULL;
-        do {
-          DoLog(0) && (Log() << Verbose(0) << "Enter index of molecule: ");
-          cin >> nr;
-          mol = molecules->ReturnIndex(nr);
-        } while (mol == NULL);
-        DoLog(0) && (Log() << Verbose(0) << "Format should be XYZ with: ShorthandOfElement\tX\tY\tZ" << endl);
-        do {
-          DoLog(0) && (Log() << Verbose(0) << "Enter file name: ");
-          cin >> filename;
-        } while (!mol->AddXYZFile(filename));
-        mol->SetNameFromFilename(filename);
-      }
-      break;
-
-    case 'r':
-      DoLog(0) && (Log() << Verbose(0) << "Enter index of molecule: ");
-      cin >> nr;
-      count = 1;
-      for(MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++)
-        if (nr == (*ListRunner)->IndexNr) {
-          mol = *ListRunner;
-          molecules->ListOfMolecules.erase(ListRunner);
-          delete(mol);
-          break;
-        }
-      break;
-  }
-};
-
-
-/** Submenu for merging molecules.
- * \param *periode periodentafel
- * \param *molecules list of molecules to add to
- */
-static void MergeMolecules(periodentafel *periode, MoleculeListClass *molecules)
-{
-  char choice;  // menu choice char
-
-  DoLog(0) && (Log() << Verbose(0) << "===========MERGE MOLECULES=====================" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "a - simple add of one molecule to another" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "b - count the number of bonds of two elements" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "B - count the number of bonds of three elements " << endl);
-  DoLog(0) && (Log() << Verbose(0) << "e - embedding merge of two molecules" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "h - count the number of hydrogen bonds" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "b - count the number of hydrogen bonds" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "m - multi-merge of all molecules" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "s - scatter merge of two molecules" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "t - simple merge of two molecules" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "all else - go back" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "===============================================" << endl);
-  DoLog(0) && (Log() << Verbose(0) << "INPUT: ");
-  cin >> choice;
-
-  switch (choice) {
-    default:
-      DoLog(0) && (Log() << Verbose(0) << "Not a valid choice." << endl);
-      break;
-
-    case 'a':
-      {
-        int src, dest;
-        molecule *srcmol = NULL, *destmol = NULL;
-        {
-          do {
-            DoLog(0) && (Log() << Verbose(0) << "Enter index of destination molecule: ");
-            cin >> dest;
-            destmol = molecules->ReturnIndex(dest);
-          } while ((destmol == NULL) && (dest != -1));
-          do {
-            DoLog(0) && (Log() << Verbose(0) << "Enter index of source molecule to add from: ");
-            cin >> src;
-            srcmol = molecules->ReturnIndex(src);
-          } while ((srcmol == NULL) && (src != -1));
-          if ((src != -1) && (dest != -1))
-            molecules->SimpleAdd(srcmol, destmol);
-        }
-      }
-      break;
-
-    case 'b':
-      {
-        const int nr = 2;
-        char *names[nr] = {"first", "second"};
-        int Z[nr];
-        element *elements[nr];
-        for (int i=0;i<nr;i++) {
-          Z[i] = 0;
-          do {
-            cout << "Enter " << names[i] << " element: ";
-            cin >> Z[i];
-          } while ((Z[i] <= 0) && (Z[i] > MAX_ELEMENTS));
-          elements[i] = periode->FindElement(Z[i]);
-        }
-        const int count = CountBondsOfTwo(molecules, elements[0], elements[1]);
-        cout << endl << "There are " << count << " ";
-        for (int i=0;i<nr;i++) {
-          if (i==0)
-            cout << elements[i]->symbol;
-          else
-            cout << "-" << elements[i]->symbol;
-        }
-        cout << " bonds." << endl;
-      }
-    break;
-
-    case 'B':
-      {
-        const int nr = 3;
-        char *names[nr] = {"first", "second", "third"};
-        int Z[nr];
-        element *elements[nr];
-        for (int i=0;i<nr;i++) {
-          Z[i] = 0;
-          do {
-            cout << "Enter " << names[i] << " element: ";
-            cin >> Z[i];
-          } while ((Z[i] <= 0) && (Z[i] > MAX_ELEMENTS));
-          elements[i] = periode->FindElement(Z[i]);
-        }
-        const int count = CountBondsOfThree(molecules, elements[0], elements[1], elements[2]);
-        cout << endl << "There are " << count << " ";
-        for (int i=0;i<nr;i++) {
-          if (i==0)
-            cout << elements[i]->symbol;
-          else
-            cout << "-" << elements[i]->symbol;
-        }
-        cout << " bonds." << endl;
-      }
-    break;
-
-    case 'e':
-      {
-        int src, dest;
-        molecule *srcmol = NULL, *destmol = NULL;
-        do {
-          DoLog(0) && (Log() << Verbose(0) << "Enter index of matrix molecule (the variable one): ");
-          cin >> src;
-          srcmol = molecules->ReturnIndex(src);
-        } while ((srcmol == NULL) && (src != -1));
-        do {
-          DoLog(0) && (Log() << Verbose(0) << "Enter index of molecule to merge into (the fixed one): ");
-          cin >> dest;
-          destmol = molecules->ReturnIndex(dest);
-        } while ((destmol == NULL) && (dest != -1));
-        if ((src != -1) && (dest != -1))
-          molecules->EmbedMerge(destmol, srcmol);
-      }
-      break;
-
-    case 'h':
-      {
-        int Z;
-        cout << "Please enter interface element: ";
-        cin >> Z;
-        element * const InterfaceElement = periode->FindElement(Z);
-        cout << endl << "There are " << CountHydrogenBridgeBonds(molecules, InterfaceElement) << " hydrogen bridges with connections to " << (InterfaceElement != 0 ? InterfaceElement->name : "None") << "." << endl;
-      }
-      break;
-
-    case 'm':
-      {
-        int nr;
-        molecule *mol = NULL;
-        do {
-          DoLog(0) && (Log() << Verbose(0) << "Enter index of molecule to merge into: ");
-          cin >> nr;
-          mol = molecules->ReturnIndex(nr);
-        } while ((mol == NULL) && (nr != -1));
-        if (nr != -1) {
-          int N = molecules->ListOfMolecules.size()-1;
-          int *src = new int(N);
-          for(MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++)
-            if ((*ListRunner)->IndexNr != nr)
-              src[N++] = (*ListRunner)->IndexNr;        
-          molecules->SimpleMultiMerge(mol, src, N);
-          delete[](src);
-        }
-      }
-      break;
-
-    case 's':
-      DoLog(0) && (Log() << Verbose(0) << "Not implemented yet." << endl);
-      break;
-
-    case 't':
-      {
-        int src, dest;
-        molecule *srcmol = NULL, *destmol = NULL;
-        {
-          do {
-            DoLog(0) && (Log() << Verbose(0) << "Enter index of destination molecule: ");
-            cin >> dest;
-            destmol = molecules->ReturnIndex(dest);
-          } while ((destmol == NULL) && (dest != -1));
-          do {
-            DoLog(0) && (Log() << Verbose(0) << "Enter index of source molecule to merge into: ");
-            cin >> src;
-            srcmol = molecules->ReturnIndex(src);
-          } while ((srcmol == NULL) && (src != -1));
-          if ((src != -1) && (dest != -1))
-            molecules->SimpleMerge(srcmol, destmol);
-        }
-      }
-      break;
-  }
-};
-
-/********************************************** Test routine **************************************/
-
-/** Is called always as option 'T' in the menu.
- * \param *molecules list of molecules
- */
-static void testroutine(MoleculeListClass *molecules)
-{
-  // the current test routine checks the functionality of the KeySet&Graph concept:
-  // We want to have a multiindex (the KeySet) describing a unique subgraph
-  int i, comp, counter=0;
-
-  // create a clone
-  molecule *mol = NULL;
-  if (molecules->ListOfMolecules.size() != 0) // clone
-    mol = (molecules->ListOfMolecules.front())->CopyMolecule();
-  else {
-    DoeLog(0) && (eLog()<< Verbose(0) << "I don't have anything to test on ... ");
-    performCriticalExit();
-    return;
-  }
-  atom *Walker = mol->start;
-
-  // generate some KeySets
-  DoLog(0) && (Log() << Verbose(0) << "Generating KeySets." << endl);
-  KeySet TestSets[mol->getAtomCount()+1];
-  i=1;
-  while (Walker->next != mol->end) {
-    Walker = Walker->next;
-    for (int j=0;j<i;j++) {
-      TestSets[j].insert(Walker->nr);
-    }
-    i++;
-  }
-  DoLog(0) && (Log() << Verbose(0) << "Testing insertion of already present item in KeySets." << endl);
-  KeySetTestPair test;
-  test = TestSets[mol->getAtomCount()-1].insert(Walker->nr);
-  if (test.second) {
-    DoLog(1) && (Log() << Verbose(1) << "Insertion worked?!" << endl);
-  } else {
-    DoLog(1) && (Log() << Verbose(1) << "Insertion rejected: Present object is " << (*test.first) << "." << endl);
-  }
-  TestSets[mol->getAtomCount()].insert(mol->end->previous->nr);
-  TestSets[mol->getAtomCount()].insert(mol->end->previous->previous->previous->nr);
-
-  // constructing Graph structure
-  DoLog(0) && (Log() << Verbose(0) << "Generating Subgraph class." << endl);
-  Graph Subgraphs;
-
-  // insert KeySets into Subgraphs
-  DoLog(0) && (Log() << Verbose(0) << "Inserting KeySets into Subgraph class." << endl);
-  for (int j=0;j<mol->getAtomCount();j++) {
-    Subgraphs.insert(GraphPair (TestSets[j],pair<int, double>(counter++, 1.)));
-  }
-  DoLog(0) && (Log() << Verbose(0) << "Testing insertion of already present item in Subgraph." << endl);
-  GraphTestPair test2;
-  test2 = Subgraphs.insert(GraphPair (TestSets[mol->getAtomCount()],pair<int, double>(counter++, 1.)));
-  if (test2.second) {
-    DoLog(1) && (Log() << Verbose(1) << "Insertion worked?!" << endl);
-  } else {
-    DoLog(1) && (Log() << Verbose(1) << "Insertion rejected: Present object is " << (*(test2.first)).second.first << "." << endl);
-  }
-
-  // show graphs
-  DoLog(0) && (Log() << Verbose(0) << "Showing Subgraph's contents, checking that it's sorted." << endl);
-  Graph::iterator A = Subgraphs.begin();
-  while (A !=  Subgraphs.end()) {
-    DoLog(0) && (Log() << Verbose(0) << (*A).second.first << ": ");
-    KeySet::iterator key = (*A).first.begin();
-    comp = -1;
-    while (key != (*A).first.end()) {
-      if ((*key) > comp)
-        DoLog(0) && (Log() << Verbose(0) << (*key) << " ");
-      else
-        DoLog(0) && (Log() << Verbose(0) << (*key) << "! ");
-      comp = (*key);
-      key++;
-    }
-    DoLog(0) && (Log() << Verbose(0) << endl);
-    A++;
-  }
-  delete(mol);
-};
-
-
-/** Tries given filename or standard on saving the config file.
- * \param *ConfigFileName name of file
- * \param *configuration pointer to configuration structure with all the values
- * \param *periode pointer to periodentafel structure with all the elements
- * \param *molecules list of molecules structure with all the atoms and coordinates
- */
-static void SaveConfig(char *ConfigFileName, config *configuration, periodentafel *periode, MoleculeListClass *molecules)
-{
-  char filename[MAXSTRINGSIZE];
-  ofstream output;
-  molecule *mol = World::getInstance().createMolecule();
-  mol->SetNameFromFilename(ConfigFileName);
-
-  if (!strcmp(configuration->configpath, configuration->GetDefaultPath())) {
-    DoeLog(2) && (eLog()<< Verbose(2) << "config is found under different path then stated in config file::defaultpath!" << endl);
-  }
-
-
-  // first save as PDB data
-  if (ConfigFileName != NULL)
-    strcpy(filename, ConfigFileName);
-  if (output == NULL)
-    strcpy(filename,"main_pcp_linux");
-  DoLog(0) && (Log() << Verbose(0) << "Saving as pdb input ");
-  if (configuration->SavePDB(filename, molecules))
-    DoLog(0) && (Log() << Verbose(0) << "done." << endl);
-  else
-    DoLog(0) && (Log() << Verbose(0) << "failed." << endl);
-
-  // then save as tremolo data file
-  if (ConfigFileName != NULL)
-    strcpy(filename, ConfigFileName);
-  if (output == NULL)
-    strcpy(filename,"main_pcp_linux");
-  DoLog(0) && (Log() << Verbose(0) << "Saving as tremolo data input ");
-  if (configuration->SaveTREMOLO(filename, molecules))
-    DoLog(0) && (Log() << Verbose(0) << "done." << endl);
-  else
-    DoLog(0) && (Log() << Verbose(0) << "failed." << endl);
-
-  // translate each to its center and merge all molecules in MoleculeListClass into this molecule
-  int N = molecules->ListOfMolecules.size();
-  int *src = new int[N];
-  N=0;
-  for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++) {
-    src[N++] = (*ListRunner)->IndexNr;
-    (*ListRunner)->Translate(&(*ListRunner)->Center);
-  }
-  molecules->SimpleMultiAdd(mol, src, N);
-  delete[](src);
-
-  // ... and translate back
-  for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++) {
-    (*ListRunner)->Center.Scale(-1.);
-    (*ListRunner)->Translate(&(*ListRunner)->Center);
-    (*ListRunner)->Center.Scale(-1.);
-  }
-
-  DoLog(0) && (Log() << Verbose(0) << "Storing configuration ... " << endl);
-  // get correct valence orbitals
-  mol->CalculateOrbitals(*configuration);
-  configuration->InitMaxMinStopStep = configuration->MaxMinStopStep = configuration->MaxPsiDouble;
-  if (ConfigFileName != NULL) { // test the file name
-    strcpy(filename, ConfigFileName);
-    output.open(filename, ios::trunc);
-  } else if (strlen(configuration->configname) != 0) {
-    strcpy(filename, configuration->configname);
-    output.open(configuration->configname, ios::trunc);
-    } else {
-      strcpy(filename, DEFAULTCONFIG);
-      output.open(DEFAULTCONFIG, ios::trunc);
-    }
-  output.close();
-  output.clear();
-  DoLog(0) && (Log() << Verbose(0) << "Saving of config file ");
-  if (configuration->Save(filename, periode, mol))
-    DoLog(0) && (Log() << Verbose(0) << "successful." << endl);
-  else
-    DoLog(0) && (Log() << Verbose(0) << "failed." << endl);
-
-  // and save to xyz file
-  if (ConfigFileName != NULL) {
-    strcpy(filename, ConfigFileName);
-    strcat(filename, ".xyz");
-    output.open(filename, ios::trunc);
-  }
-  if (output == NULL) {
-    strcpy(filename,"main_pcp_linux");
-    strcat(filename, ".xyz");
-    output.open(filename, ios::trunc);
-  }
-  DoLog(0) && (Log() << Verbose(0) << "Saving of XYZ file ");
-  if (mol->MDSteps <= 1) {
-    if (mol->OutputXYZ(&output))
-      DoLog(0) && (Log() << Verbose(0) << "successful." << endl);
-    else
-      DoLog(0) && (Log() << Verbose(0) << "failed." << endl);
-  } else {
-    if (mol->OutputTrajectoriesXYZ(&output))
-      DoLog(0) && (Log() << Verbose(0) << "successful." << endl);
-    else
-      DoLog(0) && (Log() << Verbose(0) << "failed." << endl);
-  }
-  output.close();
-  output.clear();
-
-  // and save as MPQC configuration
-  if (ConfigFileName != NULL)
-    strcpy(filename, ConfigFileName);
-  if (output == NULL)
-    strcpy(filename,"main_pcp_linux");
-  DoLog(0) && (Log() << Verbose(0) << "Saving as mpqc input ");
-  if (configuration->SaveMPQC(filename, mol))
-    DoLog(0) && (Log() << Verbose(0) << "done." << endl);
-  else
-    DoLog(0) && (Log() << Verbose(0) << "failed." << endl);
-
-  if (!strcmp(configuration->configpath, configuration->GetDefaultPath())) {
-    DoeLog(2) && (eLog()<< Verbose(2) << "config is found under different path then stated in config file::defaultpath!" << endl);
-  }
-
-  World::getInstance().destroyMolecule(mol);
-};
-
-#endif
-
-/** Parses the command line options.
- * Note that this function is from now on transitional. All commands that are not passed
- * here are handled by CommandLineParser and the actions of CommandLineUIFactory.
- * \param argc argument count
- * \param **argv arguments array
- * \param *molecules list of molecules structure
- * \param *periode elements structure
- * \param configuration config file structure
- * \param *ConfigFileName pointer to config file name in **argv
- * \param *PathToDatabases pointer to db's path in **argv
- * \param &ArgcList list of arguments that we do not parse here
- * \return exit code (0 - successful, all else - something's wrong)
- */
-static int ParseCommandLineOptions(int argc, char **argv, MoleculeListClass *&molecules, periodentafel *&periode,
-                                   config& configuration, char **ConfigFileName, set<int> &ArgcList)
-{
-  Vector x,y,z,n;  // coordinates for absolute point in cell volume
-  ifstream test;
-  ofstream output;
-  string line;
-  bool SaveFlag = false;
-  int ExitFlag = 0;
-  int j;
-  double volume = 0.;
-  enum ConfigStatus configPresent = absent;
-  int argptr;
-  molecule *mol = NULL;
-  string BondGraphFileName("\n");
-
-  if (argc > 1) { // config file specified as option
-    // 1. : Parse options that just set variables or print help
-    argptr = 1;
-    do {
-      if (argv[argptr][0] == '-') {
-        DoLog(0) && (Log() << Verbose(0) << "Recognized command line argument: " << argv[argptr][1] << ".\n");
-        argptr++;
-        switch(argv[argptr-1][1]) {
-          case 'h':
-          case 'H':
-          case '?':
-            ArgcList.insert(argptr-1);
-            return(1);
-            break;
-          case 'v':
-            if ((argptr >= argc) || (argv[argptr][0] == '-')) {
-              DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments for specifying verbosity: -v <level>" << endl);
-              performCriticalExit();
-            } else {
-              setVerbosity(atoi(argv[argptr]));
-              ArgcList.insert(argptr-1);
-              ArgcList.insert(argptr);
-              argptr++;
-            }
-            break;
-          case 'V':
-            ArgcList.insert(argptr-1);
-            return(1);
-            break;
-          case 'B':
-            if (ExitFlag == 0) ExitFlag = 1;
-            if ((argptr+5 >= argc)) {
-              DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments for setting Box: -B <xx> <<xy> <<xz> <yy> <yz> <zz>" << endl);
-              performCriticalExit();
-            } else {
-              ArgcList.insert(argptr-1);
-              ArgcList.insert(argptr);
-              ArgcList.insert(argptr+1);
-              ArgcList.insert(argptr+2);
-              ArgcList.insert(argptr+3);
-              ArgcList.insert(argptr+4);
-              ArgcList.insert(argptr+5);
-              argptr+=6;
-            }
-            break;
-          case 'e':
-            if ((argptr >= argc) || (argv[argptr][0] == '-')) {
-              DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments for specifying element db: -e <db file>" << endl);
-              performCriticalExit();
-            } else {
-              ArgcList.insert(argptr-1);
-              ArgcList.insert(argptr);
-              argptr+=1;
-            }
-            break;
-          case 'g':
-            if ((argptr >= argc) || (argv[argptr][0] == '-')) {
-              DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments for specifying bond length table: -g <table file>" << endl);
-              performCriticalExit();
-            } else {
-              ArgcList.insert(argptr-1);
-              ArgcList.insert(argptr);
-              argptr+=1;
-            }
-            break;
-          case 'M':
-            if ((argptr >= argc) || (argv[argptr][0] == '-')) {
-              ExitFlag = 255;
-              DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for setting MPQC basis: -M <basis name>" << endl);
-              performCriticalExit();
-            } else {
-              ArgcList.insert(argptr-1);
-              ArgcList.insert(argptr);
-              argptr+=1;
-            }
-            break;
-          case 'n':
-            if ((argptr >= argc) || (argv[argptr][0] == '-')) {
-              ExitFlag = 255;
-              DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for setting fast-parsing: -n <0/1>" << endl);
-              performCriticalExit();
-            } else {
-              ArgcList.insert(argptr-1);
-              ArgcList.insert(argptr);
-              argptr+=1;
-            }
-            break;
-          case 'X':
-            if ((argptr >= argc) || (argv[argptr][0] == '-')) {
-              ExitFlag = 255;
-              DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for setting default molecule name: -X <name>" << endl);
-              performCriticalExit();
-            } else {
-              ArgcList.insert(argptr-1);
-              ArgcList.insert(argptr);
-              argptr+=1;
-            }
-            break;
-          default:   // no match? Step on
-            argptr++;
-            break;
-        }
-      } else
-        argptr++;
-    } while (argptr < argc);
-
-    // 3b. Find config file name and parse if possible, also BondGraphFileName
-    if (argv[1][0] != '-') {
-      // simply create a new molecule, wherein the config file is loaded and the manipulation takes place
-      DoLog(0) && (Log() << Verbose(0) << "Config file given." << endl);
-      test.open(argv[1], ios::in);
-      if (test == NULL) {
-        //return (1);
-        output.open(argv[1], ios::out);
-        if (output == NULL) {
-          DoLog(1) && (Log() << Verbose(1) << "Specified config file " << argv[1] << " not found." << endl);
-          configPresent = absent;
-        } else {
-          DoLog(0) && (Log() << Verbose(0) << "Empty configuration file." << endl);
-          strcpy(*ConfigFileName, argv[1]);
-          configPresent = empty;
-          output.close();
-        }
-      } else {
-        test.close();
-        strcpy(*ConfigFileName, argv[1]);
-        DoLog(1) && (Log() << Verbose(1) << "Specified config file found, parsing ... ");
-        switch (configuration.TestSyntax(*ConfigFileName, periode)) {
-          case 1:
-            DoLog(0) && (Log() << Verbose(0) << "new syntax." << endl);
-            configuration.Load(*ConfigFileName, BondGraphFileName, periode, molecules);
-            configPresent = present;
-            break;
-          case 0:
-            DoLog(0) && (Log() << Verbose(0) << "old syntax." << endl);
-            configuration.LoadOld(*ConfigFileName, BondGraphFileName, periode, molecules);
-            configPresent = present;
-            break;
-          default:
-            DoLog(0) && (Log() << Verbose(0) << "Unknown syntax or empty, yet present file." << endl);
-            configPresent = empty;
-       }
-      }
-    } else
-      configPresent = absent;
-     // set mol to first active molecule
-     if (molecules->ListOfMolecules.size() != 0) {
-       for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++)
-         if ((*ListRunner)->ActiveFlag) {
-           mol = *ListRunner;
-           break;
-         }
-     }
-     if (mol == NULL) {
-       mol = World::getInstance().createMolecule();
-       mol->ActiveFlag = true;
-       if (*ConfigFileName != NULL)
-         mol->SetNameFromFilename(*ConfigFileName);
-       molecules->insert(mol);
-     }
-
-    // 4. parse again through options, now for those depending on elements db and config presence
-    argptr = 1;
-    do {
-      DoLog(0) && (Log() << Verbose(0) << "Current Command line argument: " << argv[argptr] << "." << endl);
-      if (argv[argptr][0] == '-') {
-        argptr++;
-        if ((configPresent == present) || (configPresent == empty)) {
-          switch(argv[argptr-1][1]) {
-            case 'p':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr >= argc) || (argv[argptr][0] == '-')) {
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough arguments for parsing: -p <xyz file>" << endl);
-                performCriticalExit();
-              } else {
-                SaveFlag = true;
-                DoLog(1) && (Log() << Verbose(1) << "Parsing xyz file for new atoms." << endl);
-                if (!mol->AddXYZFile(argv[argptr]))
-                  DoLog(2) && (Log() << Verbose(2) << "File not found." << endl);
-                else {
-                  DoLog(2) && (Log() << Verbose(2) << "File found and parsed." << endl);
-                  configPresent = present;
-                }
-              }
-              break;
-            case 'a':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr+4 >= argc) || (argv[argptr][0] == '-')) {
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough arguments for adding atom: -a <Z> --position <x> <y> <z>" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                ArgcList.insert(argptr+1);
-                ArgcList.insert(argptr+2);
-                ArgcList.insert(argptr+3);
-                ArgcList.insert(argptr+4);
-                argptr+=5;
-              }
-              break;
-            default:   // no match? Don't step on (this is done in next switch's default)
-              break;
-          }
-        }
-        if (configPresent == present) {
-          switch(argv[argptr-1][1]) {
-            case 'D':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr >= argc) || (argv[argptr][0] == '-')) {
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for depth-first-search analysis: -D <max. bond distance>" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                argptr+=1;
-              }
-              break;
-            case 'I':
-              DoLog(1) && (Log() << Verbose(1) << "Dissecting molecular system into a set of disconnected subgraphs ... " << endl);
-              ArgcList.insert(argptr-1);
-              argptr+=0;
-              break;
-            case 'C':
-              {
-                if (ExitFlag == 0) ExitFlag = 1;
-                if ((argptr+11 >= argc)) {
-                  ExitFlag = 255;
-                  DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for pair correlation analysis: -C[p] <type: E/P/S> [more params] <output> <bin output> <BinStart> <BinEnd>" << endl);
-                  performCriticalExit();
-                } else {
-                  switch(argv[argptr][0]) {
-                    case 'E':
-                      ArgcList.insert(argptr-1);
-                      ArgcList.insert(argptr);
-                      ArgcList.insert(argptr+1);
-                      ArgcList.insert(argptr+2);
-                      ArgcList.insert(argptr+3);
-                      ArgcList.insert(argptr+4);
-                      ArgcList.insert(argptr+5);
-                      ArgcList.insert(argptr+6);
-                      ArgcList.insert(argptr+7);
-                      ArgcList.insert(argptr+8);
-                      ArgcList.insert(argptr+9);
-                      ArgcList.insert(argptr+10);
-                      ArgcList.insert(argptr+11);
-                      argptr+=12;
-                      break;
-
-                    case 'P':
-                      ArgcList.insert(argptr-1);
-                      ArgcList.insert(argptr);
-                      ArgcList.insert(argptr+1);
-                      ArgcList.insert(argptr+2);
-                      ArgcList.insert(argptr+3);
-                      ArgcList.insert(argptr+4);
-                      ArgcList.insert(argptr+5);
-                      ArgcList.insert(argptr+6);
-                      ArgcList.insert(argptr+7);
-                      ArgcList.insert(argptr+8);
-                      ArgcList.insert(argptr+9);
-                      ArgcList.insert(argptr+10);
-                      ArgcList.insert(argptr+11);
-                      ArgcList.insert(argptr+12);
-                      ArgcList.insert(argptr+13);
-                      ArgcList.insert(argptr+14);
-                      argptr+=15;
-                      break;
-
-                    case 'S':
-                      ArgcList.insert(argptr-1);
-                      ArgcList.insert(argptr);
-                      ArgcList.insert(argptr+1);
-                      ArgcList.insert(argptr+2);
-                      ArgcList.insert(argptr+3);
-                      ArgcList.insert(argptr+4);
-                      ArgcList.insert(argptr+5);
-                      ArgcList.insert(argptr+6);
-                      ArgcList.insert(argptr+7);
-                      ArgcList.insert(argptr+8);
-                      ArgcList.insert(argptr+9);
-                      ArgcList.insert(argptr+10);
-                      ArgcList.insert(argptr+11);
-                      ArgcList.insert(argptr+12);
-                      ArgcList.insert(argptr+13);
-                      ArgcList.insert(argptr+14);
-                      argptr+=15;
-                      break;
-
-                    default:
-                      ExitFlag = 255;
-                      DoeLog(0) && (eLog()<< Verbose(0) << "Invalid type given for pair correlation analysis: -C <type: E/P/S> [more params] <output> <bin output>" << endl);
-                      performCriticalExit();
-                      break;
-                  }
-                }
-                break;
-              }
-            case 'E':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr+2 >= argc) || (!IsValidNumber(argv[argptr]))) {
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for changing element: -E <atom nr.> --element <Z>" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                ArgcList.insert(argptr+1);
-                ArgcList.insert(argptr+2);
-                argptr+=3;
-              }
-              break;
-            case 'F':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr+12 >= argc) || (argv[argptr][0] == '-')) {
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for filling with molecule: -F <filler xyz file> --MaxDistance <distance or -1> --distances <x> <y> <z>  --lengths <surface> <randatm> <randmol> --DoRotate <0/1>" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                ArgcList.insert(argptr+1);
-                ArgcList.insert(argptr+2);
-                ArgcList.insert(argptr+3);
-                ArgcList.insert(argptr+4);
-                ArgcList.insert(argptr+5);
-                ArgcList.insert(argptr+6);
-                ArgcList.insert(argptr+7);
-                ArgcList.insert(argptr+8);
-                ArgcList.insert(argptr+9);
-                ArgcList.insert(argptr+10);
-                ArgcList.insert(argptr+11);
-                ArgcList.insert(argptr+12);
-                argptr+=13;
-              }
-              break;
-            case 'A':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr+2 >= argc) || (argv[argptr][0] == '-')) {
-                ExitFlag =255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Missing source file for bonds in molecule: -A <bond sourcefile> --molecule-by-id <molecule_id>" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                ArgcList.insert(argptr+1);
-                ArgcList.insert(argptr+2);
-                argptr+=3;
-              }
-              break;
-
-            case 'J':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr+2 >= argc) || (argv[argptr][0] == '-')) {
-                ExitFlag =255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Missing path of adjacency file: -J <path> --molecule-by-id <molecule_id>" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                ArgcList.insert(argptr+1);
-                ArgcList.insert(argptr+2);
-                argptr+=3;
-              }
-              break;
-
-            case 'j':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr >= argc) || (argv[argptr][0] == '-')) {
-                ExitFlag =255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Missing path of bonds file: -j <path> --molecule-by-id <molecule_id>" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                ArgcList.insert(argptr+1);
-                ArgcList.insert(argptr+2);
-                argptr+=3;
-              }
-              break;
-
-            case 'N':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr+4 >= argc) || (argv[argptr][0] == '-')){
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for non-convex envelope: -N <molecule_id> --sphere-radius <radius> --nonconvex-file <output prefix>" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                ArgcList.insert(argptr+1);
-                ArgcList.insert(argptr+2);
-                ArgcList.insert(argptr+3);
-                ArgcList.insert(argptr+4);
-                argptr+=5;
-              }
-              break;
-            case 'S':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr+2 >= argc) || (argv[argptr][0] == '-')) {
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for storing tempature: -S <temperature file> --molecule-by-id 0" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                ArgcList.insert(argptr+1);
-                ArgcList.insert(argptr+2);
-                argptr+=3;
-              }
-              break;
-            case 'L':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr+8 >= argc) || (argv[argptr][0] == '-')) {
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for linear interpolation: -L <prefix> --start-step <step0> --end-step <step1> --molecule-by-id 0 --id-mapping <0/1>" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                ArgcList.insert(argptr+1);
-                ArgcList.insert(argptr+2);
-                ArgcList.insert(argptr+3);
-                ArgcList.insert(argptr+4);
-                ArgcList.insert(argptr+5);
-                ArgcList.insert(argptr+6);
-                ArgcList.insert(argptr+7);
-                ArgcList.insert(argptr+8);
-                argptr+=9;
-              }
-              break;
-            case 'P':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr+2 >= argc) || (argv[argptr][0] == '-')) {
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for parsing and integrating forces: -P <forces file> --molecule-by-id <molecule_id>" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                ArgcList.insert(argptr+1);
-                ArgcList.insert(argptr+2);
-                argptr+=3;
-              }
-              break;
-            case 'R':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr+4 >= argc) || (argv[argptr][0] == '-'))  {
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for removing atoms: -R <distance> --position <x> <y> <z>" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                ArgcList.insert(argptr+1);
-                ArgcList.insert(argptr+2);
-                ArgcList.insert(argptr+3);
-                ArgcList.insert(argptr+4);
-                argptr+=5;
-              }
-              break;
-            case 't':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr+4 >= argc) || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) ) {
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for translation: -t <x> <y> <z> --molecule-by-id <molecule_id> --periodic <0/1>" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                ArgcList.insert(argptr+1);
-                ArgcList.insert(argptr+2);
-                ArgcList.insert(argptr+3);
-                ArgcList.insert(argptr+4);
-                ArgcList.insert(argptr+5);
-                ArgcList.insert(argptr+6);
-                argptr+=7;
-              }
-              break;
-            case 's':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr+2 >= argc) || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) ) {
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for scaling: -s <factor_x> [factor_y] [factor_z]" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                ArgcList.insert(argptr+1);
-                ArgcList.insert(argptr+2);
-                argptr+=3;
-              }
-              break;
-            case 'b':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr+5 >= argc) || (argv[argptr][0] == '-') || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) || (!IsValidNumber(argv[argptr+3])) || (!IsValidNumber(argv[argptr+4])) || (!IsValidNumber(argv[argptr+5])) ) {
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for centering in box: -b <xx> <xy> <xz> <yy> <yz> <zz>" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                ArgcList.insert(argptr+1);
-                ArgcList.insert(argptr+2);
-                ArgcList.insert(argptr+3);
-                ArgcList.insert(argptr+4);
-                ArgcList.insert(argptr+5);
-                argptr+=6;
-              }
-              break;
-            case 'B':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr+5 >= argc) || (argv[argptr][0] == '-') || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) || (!IsValidNumber(argv[argptr+3])) || (!IsValidNumber(argv[argptr+4])) || (!IsValidNumber(argv[argptr+5])) ) {
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for bounding in box: -B <xx> <xy> <xz> <yy> <yz> <zz>" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                ArgcList.insert(argptr+1);
-                ArgcList.insert(argptr+2);
-                ArgcList.insert(argptr+3);
-                ArgcList.insert(argptr+4);
-                ArgcList.insert(argptr+5);
-                argptr+=6;
-              }
-              break;
-            case 'c':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr+2 >= argc) || (argv[argptr][0] == '-') || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) ) {
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for centering with boundary: -c <boundary_x> <boundary_y> <boundary_z>" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                ArgcList.insert(argptr+1);
-                ArgcList.insert(argptr+2);
-                argptr+=3;
-              }
-              break;
-            case 'O':
-              if (ExitFlag == 0) ExitFlag = 1;
-              ArgcList.insert(argptr-1);
-              argptr+=0;
-              break;
-            case 'r':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr >= argc) || (argv[argptr][0] == '-') || (!IsValidNumber(argv[argptr])))  {
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for removing atoms: -r <id>" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                argptr+=1;
-              }
-              break;
-            case 'f':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr+1 >= argc) || (argv[argptr][0] == '-')) {
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments for fragmentation: -f <max. bond distance> <bond order>" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                ArgcList.insert(argptr+1);
-                ArgcList.insert(argptr+2);
-                ArgcList.insert(argptr+3);
-                ArgcList.insert(argptr+4);
-                argptr+=5;
-              }
-              break;
-            case 'm':
-              if (ExitFlag == 0) ExitFlag = 1;
-              j = atoi(argv[argptr++]);
-              if ((j<0) || (j>1)) {
-                DoeLog(1) && (eLog()<< Verbose(1) << "Argument of '-m' should be either 0 for no-rotate or 1 for rotate." << endl);
-                j = 0;
-              }
-              if (j) {
-                SaveFlag = true;
-                DoLog(0) && (Log() << Verbose(0) << "Converting to prinicipal axis system." << endl);
-                mol->PrincipalAxisSystem((bool)j);
-              } else
-                ArgcList.insert(argptr-1);
-                argptr+=0;
-              break;
-            case 'o':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr+4 >= argc) || (argv[argptr][0] == '-')){
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for convex envelope: -o <molecule_id> --output-file <output file> --output-file <binned output file>" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                ArgcList.insert(argptr+1);
-                ArgcList.insert(argptr+2);
-                ArgcList.insert(argptr+3);
-                ArgcList.insert(argptr+4);
-                argptr+=5;
-              }
-              break;
-            case 'U':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr+1 >= argc) || (argv[argptr][0] == '-') || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) ) {
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for suspension with specified volume: -U <volume> <density>" << endl);
-                performCriticalExit();
-              } else {
-                volume = atof(argv[argptr++]);
-                DoLog(0) && (Log() << Verbose(0) << "Using " << volume << " angstrom^3 as the volume instead of convex envelope one's." << endl);
-              }
-            case 'u':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr >= argc) || (argv[argptr][0] == '-') || (!IsValidNumber(argv[argptr])) ) {
-                if (volume != -1)
-                  ExitFlag = 255;
-                  DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for suspension: -u <density>" << endl);
-                  performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                argptr+=1;
-              }
-              break;
-            case 'd':
-              if (ExitFlag == 0) ExitFlag = 1;
-              if ((argptr+2 >= argc) || (argv[argptr][0] == '-') || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) ) {
-                ExitFlag = 255;
-                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for repeating cells: -d <repeat_x> <repeat_y> <repeat_z>" << endl);
-                performCriticalExit();
-              } else {
-                ArgcList.insert(argptr-1);
-                ArgcList.insert(argptr);
-                ArgcList.insert(argptr+1);
-                ArgcList.insert(argptr+2);
-                argptr+=3;
-              }
-              break;
-            default:   // no match? Step on
-              if ((argptr < argc) && (argv[argptr][0] != '-')) // if it started with a '-' we've already made a step!
-                argptr++;
-              break;
-          }
-        }
-      } else argptr++;
-    } while (argptr < argc);
-    if (SaveFlag)
-      configuration.SaveAll(*ConfigFileName, periode, molecules);
-  } else {  // no arguments, hence scan the elements db
-    if (periode->LoadPeriodentafel(configuration.databasepath))
-      DoLog(0) && (Log() << Verbose(0) << "Element list loaded successfully." << endl);
-    else
-      DoLog(0) && (Log() << Verbose(0) << "Element list loading failed." << endl);
-    configuration.RetrieveConfigPathAndName("main_pcp_linux");
-  }
-  return(ExitFlag);
-};
+
 
 /********************************************** Main routine **************************************/
 
 void cleanUp(){
+  FormatParserStorage::purgeInstance();
+  ChangeTracker::purgeInstance();
   World::purgeInstance();
   logger::purgeInstance();
@@ -2200 +112 @@
 int main(int argc, char **argv)
 {
-    // while we are non interactive, we want to abort from asserts
-    //ASSERT_DO(Assert::Abort);
-    Vector x, y, z, n;
-    ifstream test;
-    ofstream output;
-    string line;
-    char **Arguments = NULL;
-    int ArgcSize = 0;
-    int ExitFlag = 0;
-    bool ArgumentsCopied = false;
-    char *ConfigFileName = new char[MAXSTRINGSIZE];
-
-    // print version check whether arguments are present at all
-    cout << ESPACKVersion << endl;
-    if (argc < 2) {
-      cout << "Obtain help with " << argv[0] << " -h." << endl;
-      cleanUp();
-      Memory::getState();
-      return(1);
-    }
-
-
-    setVerbosity(0);
-    // need to init the history before any action is created
-    ActionHistory::init();
-
-    // In the interactive mode, we can leave the user the choice in case of error
-    ASSERT_DO(Assert::Ask);
-
-    // from this moment on, we need to be sure to deeinitialize in the correct order
-    // this is handled by the cleanup function
-    atexit(cleanUp);
-
-    // Parse command line options and if present create respective UI
-    {
-      set<int> ArgcList;
-      ArgcList.insert(0); // push back program!
-      ArgcList.insert(1); // push back config file name
-      // handle arguments by ParseCommandLineOptions()
-      ExitFlag = ParseCommandLineOptions(argc,argv,World::getInstance().getMolecules(),World::getInstance().getPeriode(),*World::getInstance().getConfig(), &ConfigFileName, ArgcList);
-      World::getInstance().setExitFlag(ExitFlag);
-      // copy all remaining arguments to a new argv
-      Arguments = new char *[ArgcList.size()];
-      cout << "The following arguments are handled by CommandLineParser: ";
-      for (set<int>::iterator ArgcRunner = ArgcList.begin(); ArgcRunner != ArgcList.end(); ++ArgcRunner) {
-        Arguments[ArgcSize] = new char[strlen(argv[*ArgcRunner])+2];
-        strcpy(Arguments[ArgcSize], argv[*ArgcRunner]);
-        cout << " " << argv[*ArgcRunner];
-        ArgcSize++;
-      }
-      cout << endl;
-      ArgumentsCopied = true;
-      // handle remaining arguments by CommandLineParser
-      MapOfActions::getInstance().AddOptionsToParser();
-      map <std::string, std::string> ShortFormToActionMap = MapOfActions::getInstance().getShortFormToActionMap();
-      CommandLineParser::getInstance().Run(ArgcSize,Arguments, ShortFormToActionMap);
-      if (!CommandLineParser::getInstance().isEmpty()) {
-        DoLog(0) && (Log() << Verbose(0) << "Setting UI to CommandLine." << endl);
-        UIFactory::registerFactory(new CommandLineUIFactory::description());
-        UIFactory::makeUserInterface("CommandLine");
+  // while we are non interactive, we want to abort from asserts
+  //ASSERT_DO(Assert::Abort);
+  string line;
+  char **Arguments = NULL;
+  int ArgcSize = 0;
+  int ExitFlag = 0;
+  bool ArgumentsCopied = false;
+  std::string BondGraphFileName("\n");
+  FormatParserStorage::getInstance().addMpqc();
+  FormatParserStorage::getInstance().addPcp();
+  FormatParserStorage::getInstance().addXyz();
+
+  // print version check whether arguments are present at all
+  cout << ESPACKVersion << endl;
+  if (argc < 2) {
+    cout << "Obtain help with " << argv[0] << " -h." << endl;
+    cleanUp();
+    Memory::getState();
+    return(1);
+  }
+
+
+  setVerbosity(0);
+  // need to init the history before any action is created
+  ActionHistory::init();
+
+  // In the interactive mode, we can leave the user the choice in case of error
+  ASSERT_DO(Assert::Ask);
+
+  // from this moment on, we need to be sure to deeinitialize in the correct order
+  // this is handled by the cleanup function
+  atexit(cleanUp);
+
+  // Parse command line options and if present create respective UI
+  {
+    // construct bond graph
+    if (World::getInstance().getConfig()->BG == NULL) {
+      World::getInstance().getConfig()->BG = new BondGraph(World::getInstance().getConfig()->GetIsAngstroem());
+      if (World::getInstance().getConfig()->BG->LoadBondLengthTable(BondGraphFileName)) {
+        DoLog(0) && (Log() << Verbose(0) << "Bond length table loaded successfully." << endl);
       } else {
-          #ifdef USE_GUI_QT
-            DoLog(0) && (Log() << Verbose(0) << "Setting UI to QT4." << endl);
-            UIFactory::registerFactory(new QTUIFactory::description());
-            UIFactory::makeUserInterface("QT4");
-          #else
-            DoLog(0) && (Log() << Verbose(0) << "Setting UI to Text." << endl);
-            cout << ESPACKVersion << endl;
-            UIFactory::registerFactory(new TextUIFactory::description());
-            UIFactory::makeUserInterface("Text");
-          #endif
+        DoeLog(1) && (eLog()<< Verbose(1) << "Bond length table loading failed." << endl);
       }
     }
-
-    {
-      MainWindow *mainWindow = UIFactory::getInstance().makeMainWindow();
-      mainWindow->display();
-      delete mainWindow;
+    // handle remaining arguments by CommandLineParser
+    MapOfActions::getInstance().AddOptionsToParser();
+    map <std::string, std::string> ShortFormToActionMap = MapOfActions::getInstance().getShortFormToActionMap();
+    CommandLineParser::getInstance().Run(argc,argv, ShortFormToActionMap);
+    if (!CommandLineParser::getInstance().isEmpty()) {
+      DoLog(0) && (Log() << Verbose(0) << "Setting UI to CommandLine." << endl);
+      UIFactory::registerFactory(new CommandLineUIFactory::description());
+      UIFactory::makeUserInterface("CommandLine");
+    } else {
+      #ifdef USE_GUI_QT
+        DoLog(0) && (Log() << Verbose(0) << "Setting UI to QT4." << endl);
+        UIFactory::registerFactory(new QTUIFactory::description());
+        UIFactory::makeUserInterface("QT4");
+      #else
+        DoLog(0) && (Log() << Verbose(0) << "Setting UI to Text." << endl);
+        cout << ESPACKVersion << endl;
+        UIFactory::registerFactory(new TextUIFactory::description());
+        UIFactory::makeUserInterface("Text");
+      #endif
     }
-
-    Log() << Verbose(0) << "Saving to " << ConfigFileName << "." << endl;
-    World::getInstance().getConfig()->SaveAll(ConfigFileName, World::getInstance().getPeriode(), World::getInstance().getMolecules());
+  }
+
+  {
+    MainWindow *mainWindow = UIFactory::getInstance().makeMainWindow();
+    mainWindow->display();
+    delete mainWindow;
+  }
+
+  FormatParserStorage::getInstance().SaveAll();
+  ChangeTracker::getInstance().saveStatus();
 
   // free the new argv
@@ -2289 +193 @@
     delete[](Arguments);
   }
-  delete[](ConfigFileName);
+  //delete[](ConfigFileName);
 
   ExitFlag = World::getInstance().getExitFlag();

src/config.cpp

@@ -10 +10 @@
 #include <cstring>
 
-#include "World.hpp"
 #include "atom.hpp"
 #include "bond.hpp"
+#include "bondgraph.hpp"
 #include "config.hpp"
+#include "ConfigFileBuffer.hpp"
 #include "element.hpp"
 #include "helpers.hpp"
@@ -23 +24 @@
 #include "molecule.hpp"
 #include "periodentafel.hpp"
+#include "ThermoStatContainer.hpp"
 #include "World.hpp"
-
-/******************************** Functions for class ConfigFileBuffer **********************/
-
-/** Structure containing compare function for Ion_Type sorting.
- */
-struct IonTypeCompare {
-  bool operator()(const char* s1, const char *s2) const {
-    char number1[8];
-    char number2[8];
-    const char *dummy1, *dummy2;
-    //Log() << Verbose(0) << s1 << "  " << s2 << endl;
-    dummy1 = strchr(s1, '_')+sizeof(char)*5;  // go just after "Ion_Type"
-    dummy2 = strchr(dummy1, '_');
-    strncpy(number1, dummy1, dummy2-dummy1); // copy the number
-    number1[dummy2-dummy1]='\0';
-    dummy1 = strchr(s2, '_')+sizeof(char)*5;  // go just after "Ion_Type"
-    dummy2 = strchr(dummy1, '_');
-    strncpy(number2, dummy1, dummy2-dummy1); // copy the number
-    number2[dummy2-dummy1]='\0';
-    if (atoi(number1) != atoi(number2))
-      return (atoi(number1) < atoi(number2));
-    else {
-      dummy1 = strchr(s1, '_')+sizeof(char);
-      dummy1 = strchr(dummy1, '_')+sizeof(char);
-      dummy2 = strchr(dummy1, ' ') < strchr(dummy1, '\t') ? strchr(dummy1, ' ') : strchr(dummy1, '\t');
-      strncpy(number1, dummy1, dummy2-dummy1); // copy the number
-      number1[dummy2-dummy1]='\0';
-      dummy1 = strchr(s2, '_')+sizeof(char);
-      dummy1 = strchr(dummy1, '_')+sizeof(char);
-      dummy2 = strchr(dummy1, ' ') < strchr(dummy1, '\t') ? strchr(dummy1, ' ') : strchr(dummy1, '\t');
-      strncpy(number2, dummy1, dummy2-dummy1); // copy the number
-      number2[dummy2-dummy1]='\0';
-      return (atoi(number1) < atoi(number2));
-    }
-  }
-};
-
-/** Constructor for ConfigFileBuffer class.
- */
-ConfigFileBuffer::ConfigFileBuffer() : buffer(NULL), LineMapping(NULL), CurrentLine(0), NoLines(0)
-{
-};
-
-/** Constructor for ConfigFileBuffer class with filename to be parsed.
- * \param *filename file name
- */
-ConfigFileBuffer::ConfigFileBuffer(const char * const filename) : buffer(NULL), LineMapping(NULL), CurrentLine(0), NoLines(0)
-{
-  ifstream *file = NULL;
-  char line[MAXSTRINGSIZE];
-
-  // prescan number of lines
-  file= new ifstream(filename);
-  if (file == NULL) {
-    DoeLog(1) && (eLog()<< Verbose(1) << "config file " << filename << " missing!" << endl);
-    return;
-  }
-  NoLines = 0; // we're overcounting by one
-  long file_position = file->tellg(); // mark current position
-  do {
-    file->getline(line, 256);
-    NoLines++;
-  } while (!file->eof());
-  file->clear();
-  file->seekg(file_position, ios::beg);
-  DoLog(1) && (Log() << Verbose(1) << NoLines-1 << " lines were recognized." << endl);
-
-  // allocate buffer's 1st dimension
-  if (buffer != NULL) {
-    DoeLog(1) && (eLog()<< Verbose(1) << "FileBuffer->buffer is not NULL!" << endl);
-    return;
-  } else
-    buffer = new char *[NoLines];
-
-  // scan each line and put into buffer
-  int lines=0;
-  int i;
-  do {
-    buffer[lines] = new char[MAXSTRINGSIZE];
-    file->getline(buffer[lines], MAXSTRINGSIZE-1);
-    i = strlen(buffer[lines]);
-    buffer[lines][i] = '\n';
-    buffer[lines][i+1] = '\0';
-    lines++;
-  } while((!file->eof()) && (lines < NoLines));
-  DoLog(1) && (Log() << Verbose(1) << lines-1 << " lines were read into the buffer." << endl);
-
-  // close and exit
-  file->close();
-  file->clear();
-  delete(file);
-}
-
-/** Destructor for ConfigFileBuffer class.
- */
-ConfigFileBuffer::~ConfigFileBuffer()
-{
-  for(int i=0;i<NoLines;++i)
-    delete[](buffer[i]);
-  delete[](buffer);
-  delete[](LineMapping);
-}
-
-
-/** Create trivial mapping.
- */
-void ConfigFileBuffer::InitMapping()
-{
-  LineMapping = new int[NoLines];
-  for (int i=0;i<NoLines;i++)
-    LineMapping[i] = i;
-}
-
-/** Creates a mapping for the \a *FileBuffer's lines containing the Ion_Type keyword such that they are sorted.
- * \a *map on return contains a list of NoAtom entries such that going through the list, yields indices to the
- * lines in \a *FileBuffer in a sorted manner of the Ion_Type?_? keywords. We assume that ConfigFileBuffer::CurrentLine
- * points to first Ion_Type entry.
- * \param *FileBuffer pointer to buffer structure
- * \param NoAtoms of subsequent lines to look at
- */
-void ConfigFileBuffer::MapIonTypesInBuffer(const int NoAtoms)
-{
-  map<const char *, int, IonTypeCompare> IonTypeLineMap;
-  if (LineMapping == NULL) {
-    DoeLog(0) && (eLog()<< Verbose(0) << "map pointer is NULL: " << LineMapping << endl);
-    performCriticalExit();
-    return;
-  }
-
-  // put all into hashed map
-  for (int i=0; i<NoAtoms; ++i) {
-    IonTypeLineMap.insert(pair<const char *, int> (buffer[CurrentLine+i], CurrentLine+i));
-  }
-
-  // fill map
-  int nr=0;
-  for (map<const char *, int, IonTypeCompare>::iterator runner = IonTypeLineMap.begin(); runner != IonTypeLineMap.end(); ++runner) {
-    if (CurrentLine+nr < NoLines)
-      LineMapping[CurrentLine+(nr++)] = runner->second;
-    else {
-      DoeLog(0) && (eLog()<< Verbose(0) << "config::MapIonTypesInBuffer - NoAtoms is wrong: We are past the end of the file!" << endl);
-      performCriticalExit();
-    }
-  }
-}
 
 /************************************* Functions for class config ***************************/
@@ -174 +31 @@
 /** Constructor for config file class.
  */
-config::config() : BG(NULL), PsiType(0), MaxPsiDouble(0), PsiMaxNoUp(0), PsiMaxNoDown(0), MaxMinStopStep(1), InitMaxMinStopStep(1), ProcPEGamma(8), ProcPEPsi(1), configpath(NULL),
-    configname(NULL), FastParsing(false), Deltat(0.01), basis(""), databasepath(NULL), DoConstrainedMD(0), MaxOuterStep(0), Thermostat(4), ThermostatImplemented(NULL),
-    ThermostatNames(NULL), TempFrequency(2.5), alpha(0.), HooverMass(0.), TargetTemp(0.00095004455), ScaleTempStep(25),  mainname(NULL), defaultpath(NULL), pseudopotpath(NULL),
+config::config() : BG(NULL), Thermostats(0), PsiType(0), MaxPsiDouble(0), PsiMaxNoUp(0), PsiMaxNoDown(0), MaxMinStopStep(1), InitMaxMinStopStep(1), ProcPEGamma(8), ProcPEPsi(1),
+    configname(NULL), FastParsing(false), Deltat(0.01), basis(""), databasepath(NULL), DoConstrainedMD(0), MaxOuterStep(0), mainname(NULL), defaultpath(NULL), pseudopotpath(NULL),
     DoOutVis(0), DoOutMes(1), DoOutNICS(0), DoOutOrbitals(0), DoOutCurrent(0), DoFullCurrent(0), DoPerturbation(0), DoWannier(0), CommonWannier(0), SawtoothStart(0.01),
     VectorPlane(0), VectorCut(0.), UseAddGramSch(1), Seed(1), OutVisStep(10), OutSrcStep(5), MaxPsiStep(0), EpsWannier(1e-7), MaxMinStep(100), RelEpsTotalEnergy(1e-7),
@@ -186 +42 @@
   pseudopotpath = new char[MAXSTRINGSIZE];
   databasepath = new char[MAXSTRINGSIZE];
-  configpath = new char[MAXSTRINGSIZE];
   configname = new char[MAXSTRINGSIZE];
+  Thermostats = new ThermoStatContainer();
   strcpy(mainname,"pcp");
   strcpy(defaultpath,"not specified");
   strcpy(pseudopotpath,"not specified");
-  configpath[0]='\0';
   configname[0]='\0';
   basis = "3-21G";
-
-  InitThermostats();
 };
 
@@ -206 +59 @@
   delete[](pseudopotpath);
   delete[](databasepath);
-  delete[](configpath);
   delete[](configname);
-  delete[](ThermostatImplemented);
-  for (int j=0;j<MaxThermostats;j++)
-    delete[](ThermostatNames[j]);
-  delete[](ThermostatNames);
+  if (Thermostats != NULL)
+    delete(Thermostats);
 
   if (BG != NULL)
     delete(BG);
 };
-
-/** Initialises variables in class config for Thermostats.
- */
-void config::InitThermostats()
-{
-  ThermostatImplemented = new int[MaxThermostats];
-  ThermostatNames = new char *[MaxThermostats];
-  for (int j=0;j<MaxThermostats;j++)
-    ThermostatNames[j] = new char[12];
-
-  strcpy(ThermostatNames[0],"None");
-  ThermostatImplemented[0] = 1;
-  strcpy(ThermostatNames[1],"Woodcock");
-  ThermostatImplemented[1] = 1;
-  strcpy(ThermostatNames[2],"Gaussian");
-  ThermostatImplemented[2] = 1;
-  strcpy(ThermostatNames[3],"Langevin");
-  ThermostatImplemented[3] = 1;
-  strcpy(ThermostatNames[4],"Berendsen");
-  ThermostatImplemented[4] = 1;
-  strcpy(ThermostatNames[5],"NoseHoover");
-  ThermostatImplemented[5] = 1;
-};
-
-/** Readin of Thermostat related values from parameter file.
- * \param *fb file buffer containing the config file
- */
-void config::ParseThermostats(class ConfigFileBuffer * const fb)
-{
-  char * const thermo = new char[12];
-  const int verbose = 0;
-
-  // read desired Thermostat from file along with needed additional parameters
-  if (ParseForParameter(verbose,fb,"Thermostat", 0, 1, 1, string_type, thermo, 1, optional)) {
-    if (strcmp(thermo, ThermostatNames[0]) == 0) { // None
-      if (ThermostatImplemented[0] == 1) {
-        Thermostat = None;
-      } else {
-        DoLog(1) && (Log() << Verbose(1) << "Warning: " << ThermostatNames[0] << " thermostat not implemented, falling back to None." << endl);
-        Thermostat = None;
-      }
-    } else if (strcmp(thermo, ThermostatNames[1]) == 0) { // Woodcock
-      if (ThermostatImplemented[1] == 1) {
-        Thermostat = Woodcock;
-        ParseForParameter(verbose,fb,"Thermostat", 0, 2, 1, int_type, &ScaleTempStep, 1, critical); // read scaling frequency
-      } else {
-        DoLog(1) && (Log() << Verbose(1) << "Warning: " << ThermostatNames[0] << " thermostat not implemented, falling back to None." << endl);
-        Thermostat = None;
-      }
-    } else if (strcmp(thermo, ThermostatNames[2]) == 0) { // Gaussian
-      if (ThermostatImplemented[2] == 1) {
-        Thermostat = Gaussian;
-        ParseForParameter(verbose,fb,"Thermostat", 0, 2, 1, int_type, &ScaleTempStep, 1, critical); // read collision rate
-      } else {
-        DoLog(1) && (Log() << Verbose(1) << "Warning: " << ThermostatNames[0] << " thermostat not implemented, falling back to None." << endl);
-        Thermostat = None;
-      }
-    } else if (strcmp(thermo, ThermostatNames[3]) == 0) { // Langevin
-      if (ThermostatImplemented[3] == 1) {
-        Thermostat = Langevin;
-        ParseForParameter(verbose,fb,"Thermostat", 0, 2, 1, double_type, &TempFrequency, 1, critical); // read gamma
-        if (ParseForParameter(verbose,fb,"Thermostat", 0, 3, 1, double_type, &alpha, 1, optional)) {
-          DoLog(2) && (Log() << Verbose(2) << "Extended Stochastic Thermostat detected with interpolation coefficient " << alpha << "." << endl);
-        } else {
-          alpha = 1.;
-        }
-      } else {
-        DoLog(1) && (Log() << Verbose(1) << "Warning: " << ThermostatNames[0] << " thermostat not implemented, falling back to None." << endl);
-        Thermostat = None;
-      }
-    } else if (strcmp(thermo, ThermostatNames[4]) == 0) { // Berendsen
-      if (ThermostatImplemented[4] == 1) {
-        Thermostat = Berendsen;
-        ParseForParameter(verbose,fb,"Thermostat", 0, 2, 1, double_type, &TempFrequency, 1, critical); // read \tau_T
-      } else {
-        DoLog(1) && (Log() << Verbose(1) << "Warning: " << ThermostatNames[0] << " thermostat not implemented, falling back to None." << endl);
-        Thermostat = None;
-      }
-    } else if (strcmp(thermo, ThermostatNames[5]) == 0) { // Nose-Hoover
-      if (ThermostatImplemented[5] == 1) {
-        Thermostat = NoseHoover;
-        ParseForParameter(verbose,fb,"Thermostat", 0, 2, 1, double_type, &HooverMass, 1, critical); // read Hoovermass
-        alpha = 0.;
-      } else {
-        DoLog(1) && (Log() << Verbose(1) << "Warning: " << ThermostatNames[0] << " thermostat not implemented, falling back to None." << endl);
-        Thermostat = None;
-      }
-    } else {
-      DoLog(1) && (Log() << Verbose(1) << " Warning: thermostat name was not understood!" << endl);
-      Thermostat = None;
-    }
-  } else {
-    if ((MaxOuterStep > 0) && (TargetTemp != 0)) 
-      DoLog(2) && (Log() << Verbose(2) <<  "No thermostat chosen despite finite temperature MD, falling back to None." << endl);
-    Thermostat = None;
-  }
-  delete[](thermo);
-};
-
 
 /** Displays menu for editing each entry of the config file.
@@ -629 +380 @@
 };
 
-/** Retrieves the path in the given config file name.
- * \param filename config file string
- */
-void config::RetrieveConfigPathAndName(const string filename)
-{
-  char *ptr = NULL;
-  char *buffer = new char[MAXSTRINGSIZE];
-  strncpy(buffer, filename.c_str(), MAXSTRINGSIZE);
-  int last = -1;
-  for(last=MAXSTRINGSIZE;last--;) {
-    if (buffer[last] == '/')
-      break;
-  }
-  if (last == -1) { // no path in front, set to local directory.
-    strcpy(configpath, "./");
-    ptr = buffer;
-  } else {
-    strncpy(configpath, buffer, last+1);
-    ptr = &buffer[last+1];
-    if (last < 254)
-      configpath[last+1]='\0';
-  }
-  strcpy(configname, ptr);
-  DoLog(0) && (Log() << Verbose(0) << "Found configpath: " << configpath << ", dir slash was found at " << last << ", config name is " << configname << "." << endl);
-  delete[](buffer);
-};
-
-/** Initializes ConfigFileBuffer from a file.
- * \param *file input file stream being the opened config file
- * \param *FileBuffer pointer to FileBuffer on return, should point to NULL
- */
-void PrepareFileBuffer(const char * const filename, struct ConfigFileBuffer *&FileBuffer)
-{
-  if (FileBuffer != NULL) {
-    DoeLog(2) && (eLog()<< Verbose(2) << "deleting present FileBuffer in PrepareFileBuffer()." << endl);
-    delete(FileBuffer);
-  }
-  FileBuffer = new ConfigFileBuffer(filename);
-
-  FileBuffer->InitMapping();
-};
-
 /** Loads a molecule from a ConfigFileBuffer.
  * \param *mol molecule to load
@@ -864 +573 @@
   file->close();
   delete(file);
-  RetrieveConfigPathAndName(filename);
 
   // ParseParameterFile
-  struct ConfigFileBuffer *FileBuffer = NULL;
-  PrepareFileBuffer(filename,FileBuffer);
+  class ConfigFileBuffer *FileBuffer = new ConfigFileBuffer(filename);
 
   /* Oeffne Hauptparameterdatei */
@@ -877 +584 @@
   int verbose = 0;
   
+  //TODO: This is actually sensible?: if (MaxOuterStep > 0)
   ParseThermostats(FileBuffer);
   
@@ -941 +649 @@
   ParseForParameter(verbose,FileBuffer,"OutVisStep", 0, 1, 1, int_type, &(config::OutVisStep), 1, optional);
   ParseForParameter(verbose,FileBuffer,"OutSrcStep", 0, 1, 1, int_type, &(config::OutSrcStep), 1, optional);
-  ParseForParameter(verbose,FileBuffer,"TargetTemp", 0, 1, 1, double_type, &(config::TargetTemp), 1, optional);
+  ParseForParameter(verbose,FileBuffer,"TargetTemp", 0, 1, 1, double_type, &(Thermostats->TargetTemp), 1, optional);
   //ParseForParameter(verbose,FileBuffer,"Thermostat", 0, 1, 1, int_type, &(config::ScaleTempStep), 1, optional);
   if (!ParseForParameter(verbose,FileBuffer,"EpsWannier", 0, 1, 1, double_type, &(config::EpsWannier), 1, optional))
@@ -1101 +809 @@
     return;
   }
-  RetrieveConfigPathAndName(filename);
   // ParseParameters
 
@@ -1150 +857 @@
   ParseForParameter(verbose,file,"VisOuterStep", 0, 1, 1, int_type, &(config::OutVisStep), 1, optional);
   ParseForParameter(verbose,file,"VisSrcOuterStep", 0, 1, 1, int_type, &(config::OutSrcStep), 1, optional);
-  ParseForParameter(verbose,file,"TargetTemp", 0, 1, 1, double_type, &(config::TargetTemp), 1, optional);
-  ParseForParameter(verbose,file,"ScaleTempStep", 0, 1, 1, int_type, &(config::ScaleTempStep), 1, optional);
+  ParseForParameter(verbose,file,"TargetTemp", 0, 1, 1, double_type, &(Thermostats->TargetTemp), 1, optional);
+  ParseForParameter(verbose,file,"ScaleTempStep", 0, 1, 1, int_type, &(Thermostats->ScaleTempStep), 1, optional);
   config::EpsWannier = 1e-8;
 
@@ -1339 +1046 @@
     *output << "DoFullCurrent\t" << config::DoFullCurrent << "\t# Do full perturbation" << endl;
     *output << "DoConstrainedMD\t" << config::DoConstrainedMD << "\t# Do perform a constrained (>0, relating to current MD step) instead of unconstrained (0) MD" << endl;
-    *output << "Thermostat\t" << ThermostatNames[Thermostat] << "\t";
-    switch(Thermostat) {
+    *output << "Thermostat\t" << Thermostats->ThermostatNames[Thermostats->Thermostat] << "\t";
+    switch(Thermostats->Thermostat) {
       default:
       case None:
         break;
       case Woodcock:
-        *output << ScaleTempStep;
+        *output << Thermostats->ScaleTempStep;
         break;
       case Gaussian:
-        *output << ScaleTempStep;
+        *output << Thermostats->ScaleTempStep;
         break;
       case Langevin:
-        *output << TempFrequency << "\t" << alpha;
+        *output << Thermostats->TempFrequency << "\t" << Thermostats->alpha;
         break;
       case Berendsen:
-        *output << TempFrequency;
+        *output << Thermostats->TempFrequency;
         break;
       case NoseHoover:
-        *output << HooverMass;
+        *output << Thermostats->HooverMass;
         break;
     };
@@ -1372 +1079 @@
     *output << "OutVisStep\t" << config::OutVisStep << "\t# Output visual data every ...th step" << endl;
     *output << "OutSrcStep\t" << config::OutSrcStep << "\t# Output \"restart\" data every ..th step" << endl;
-    *output << "TargetTemp\t" << config::TargetTemp << "\t# Target temperature" << endl;
+    *output << "TargetTemp\t" << Thermostats->TargetTemp << "\t# Target temperature" << endl;
     *output << "MaxPsiStep\t" << config::MaxPsiStep << "\t# number of Minimisation steps per state (0 - default)" << endl;
     *output << "EpsWannier\t" << config::EpsWannier << "\t# tolerance value for spread minimisation of orbitals" << endl;
@@ -1483 +1190 @@
     // output of atoms
     AtomNo = 0;
-    mol->ActOnAllAtoms( &atom::OutputMPQCLine, output, (const Vector *)center, &AtomNo );
+    mol->ActOnAllAtoms( &atom::OutputMPQCLine, (ostream * const) output, (const Vector *)center, &AtomNo );
     delete(center);
     *output << "\t}" << endl;
@@ -1525 +1232 @@
     // output of atoms
     AtomNo = 0;
-    mol->ActOnAllAtoms( &atom::OutputMPQCLine, output, (const Vector *)center, &AtomNo );
+    mol->ActOnAllAtoms( &atom::OutputMPQCLine, (ostream * const) output, (const Vector *)center, &AtomNo );
     delete(center);
     *output << "\t}" << endl;
@@ -1786 +1493 @@
   molecule *mol = NULL;
 
-  if (!strcmp(configpath, GetDefaultPath())) {
-    eLog() << Verbose(2) << "config is found under different path then stated in config file::defaultpath!" << endl;
-  }
-
-
   // first save as PDB data
   if (ConfigFileName != NULL)
@@ -1827 +1529 @@
     mol->doCountAtoms();
     mol->CountElements();
-    mol->CalculateOrbitals(*this);
+    //mol->CalculateOrbitals(*this);
     delete[](src);
   } else {
@@ -1833 +1535 @@
       mol = *(molecules->ListOfMolecules.begin());
       mol->doCountAtoms();
-      mol->CalculateOrbitals(*this);
+      //mol->CalculateOrbitals(*this);
     } else {
       DoeLog(1) && (eLog() << Verbose(1) << "There are no molecules to save!" << endl);
@@ -1895 +1597 @@
   else
     Log() << Verbose(0) << "\t... failed." << endl;
-
-  if (!strcmp(configpath, GetDefaultPath())) {
-    eLog() << Verbose(2) << "config is found under different path then stated in config file::defaultpath!" << endl;
-  }
 
   // don't destroy molecule as it contains all our atoms
@@ -2346 +2044 @@
   return (found); // true if found, false if not
 }
+
+/** Reading of Thermostat related values from parameter file.
+ * \param *fb file buffer containing the config file
+ */
+void config::ParseThermostats(class ConfigFileBuffer * const fb)
+{
+  char * const thermo = new char[12];
+  const int verbose = 0;
+
+  // read desired Thermostat from file along with needed additional parameters
+  if (ParseForParameter(verbose,fb,"Thermostat", 0, 1, 1, string_type, thermo, 1, optional)) {
+    if (strcmp(thermo, Thermostats->ThermostatNames[0]) == 0) { // None
+      if (Thermostats->ThermostatImplemented[0] == 1) {
+        Thermostats->Thermostat = None;
+      } else {
+        DoLog(1) && (Log() << Verbose(1) << "Warning: " << Thermostats->ThermostatNames[0] << " thermostat not implemented, falling back to None." << endl);
+        Thermostats->Thermostat = None;
+      }
+    } else if (strcmp(thermo, Thermostats->ThermostatNames[1]) == 0) { // Woodcock
+      if (Thermostats->ThermostatImplemented[1] == 1) {
+        Thermostats->Thermostat = Woodcock;
+        ParseForParameter(verbose,fb,"Thermostat", 0, 2, 1, int_type, &Thermostats->ScaleTempStep, 1, critical); // read scaling frequency
+      } else {
+        DoLog(1) && (Log() << Verbose(1) << "Warning: " << Thermostats->ThermostatNames[0] << " thermostat not implemented, falling back to None." << endl);
+        Thermostats->Thermostat = None;
+      }
+    } else if (strcmp(thermo, Thermostats->ThermostatNames[2]) == 0) { // Gaussian
+      if (Thermostats->ThermostatImplemented[2] == 1) {
+        Thermostats->Thermostat = Gaussian;
+        ParseForParameter(verbose,fb,"Thermostat", 0, 2, 1, int_type, &Thermostats->ScaleTempStep, 1, critical); // read collision rate
+      } else {
+        DoLog(1) && (Log() << Verbose(1) << "Warning: " << Thermostats->ThermostatNames[0] << " thermostat not implemented, falling back to None." << endl);
+        Thermostats->Thermostat = None;
+      }
+    } else if (strcmp(thermo, Thermostats->ThermostatNames[3]) == 0) { // Langevin
+      if (Thermostats->ThermostatImplemented[3] == 1) {
+        Thermostats->Thermostat = Langevin;
+        ParseForParameter(verbose,fb,"Thermostat", 0, 2, 1, double_type, &Thermostats->TempFrequency, 1, critical); // read gamma
+        if (ParseForParameter(verbose,fb,"Thermostat", 0, 3, 1, double_type, &Thermostats->alpha, 1, optional)) {
+          DoLog(2) && (Log() << Verbose(2) << "Extended Stochastic Thermostat detected with interpolation coefficient " << Thermostats->alpha << "." << endl);
+        } else {
+          Thermostats->alpha = 1.;
+        }
+      } else {
+        DoLog(1) && (Log() << Verbose(1) << "Warning: " << Thermostats->ThermostatNames[0] << " thermostat not implemented, falling back to None." << endl);
+        Thermostats->Thermostat = None;
+      }
+    } else if (strcmp(thermo, Thermostats->ThermostatNames[4]) == 0) { // Berendsen
+      if (Thermostats->ThermostatImplemented[4] == 1) {
+        Thermostats->Thermostat = Berendsen;
+        ParseForParameter(verbose,fb,"Thermostat", 0, 2, 1, double_type, &Thermostats->TempFrequency, 1, critical); // read \tau_T
+      } else {
+        DoLog(1) && (Log() << Verbose(1) << "Warning: " << Thermostats->ThermostatNames[0] << " thermostat not implemented, falling back to None." << endl);
+        Thermostats->Thermostat = None;
+      }
+    } else if (strcmp(thermo, Thermostats->ThermostatNames[5]) == 0) { // Nose-Hoover
+      if (Thermostats->ThermostatImplemented[5] == 1) {
+        Thermostats->Thermostat = NoseHoover;
+        ParseForParameter(verbose,fb,"Thermostat", 0, 2, 1, double_type, &Thermostats->HooverMass, 1, critical); // read Hoovermass
+        Thermostats->alpha = 0.;
+      } else {
+        DoLog(1) && (Log() << Verbose(1) << "Warning: " << Thermostats->ThermostatNames[0] << " thermostat not implemented, falling back to None." << endl);
+        Thermostats->Thermostat = None;
+      }
+    } else {
+      DoLog(1) && (Log() << Verbose(1) << " Warning: thermostat name was not understood!" << endl);
+      Thermostats->Thermostat = None;
+    }
+  } else {
+    if ((Thermostats->TargetTemp != 0))
+      DoLog(2) && (Log() << Verbose(2) <<  "No thermostat chosen despite finite temperature MD, falling back to None." << endl);
+    Thermostats->Thermostat = None;
+  }
+  delete[](thermo);
+};
+

src/config.hpp

@@ -20 +20 @@
 #include <string>
 
-#include "bondgraph.hpp"
-
 /****************************************** forward declarations *****************************/
 
+class BondGraph;
+class ConfigFileBuffer;
 class molecule;
 class MoleculeListClass;
 class periodentafel;
+class ThermoStatContainer;
 
 /********************************************** declarations *******************************/
-
-class ConfigFileBuffer {
-  public:
-    char **buffer;
-    int *LineMapping;
-    int CurrentLine;
-    int NoLines;
-
-    ConfigFileBuffer();
-    ConfigFileBuffer(const char * const filename);
-    ~ConfigFileBuffer();
-
-    void InitMapping();
-    void MapIonTypesInBuffer(const int NoAtoms);
-};
 
 /** The config file.
@@ -51 +37 @@
   public:
     class BondGraph *BG;
+    class ThermoStatContainer *Thermostats;
 
     int PsiType;
@@ -60 +47 @@
     int ProcPEGamma;
     int ProcPEPsi;
-    char *configpath;
     char *configname;
     bool FastParsing;
@@ -70 +56 @@
     int DoConstrainedMD;
     int MaxOuterStep;
-    int Thermostat;
-    int *ThermostatImplemented;
-    char **ThermostatNames;
-    double TempFrequency;
-    double alpha;
-    double HooverMass;
-    double TargetTemp;
-    int ScaleTempStep;
 
   private:
@@ -138 +116 @@
   void Load(const char * const filename, const string &BondGraphFileName, const periodentafel * const periode, MoleculeListClass * const &MolList);
   void LoadOld(const char * const filename, const string &BondGraphFileName, const periodentafel * const periode, MoleculeListClass * const &MolList);
-  void RetrieveConfigPathAndName(const string filename);
   bool Save(const char * const filename, const periodentafel * const periode, molecule * const mol) const;
   bool SaveMPQC(const char * const filename, const molecule * const mol) const;
@@ -152 +129 @@
   char *GetDefaultPath() const;
   void SetDefaultPath(const char * const path);
-  void InitThermostats();
   void ParseThermostats(class ConfigFileBuffer * const fb);
 };

src/molecule.cpp

@@ -79 +79 @@
 void molecule::setName(const std::string _name){
   OBSERVE;
+  cout << "Set name of molecule " << getId() << " to " << _name << endl;
   strncpy(name,_name.c_str(),MAXSTRINGSIZE);
 }
@@ -154 +155 @@
 molecule::const_iterator molecule::erase( atom * key )
 {
-  cout << "trying to erase atom" << endl;
   molecule::const_iterator iter = find(key);
   if (iter != end()){
@@ -740 +740 @@
   else
     length = strlen(molname) - strlen(endname);
+  cout << "Set name of molecule " << getId() << " to " << molname << endl;
   strncpy(name, molname, length);
   name[length]='\0';
@@ -881 +882 @@
         ElementNo[i] = current++;
     }
-    ActOnAllAtoms( &atom::OutputArrayIndexed, output, (const int *)ElementNo, (int *)AtomNo, (const char *) NULL );
+    ActOnAllAtoms( &atom::OutputArrayIndexed, (ostream * const) output, (const int *)ElementNo, (int *)AtomNo, (const char *) NULL );
     return true;
   }
@@ -1004 +1005 @@
   for(int i=MAX_ELEMENTS;i--;)
     ElementCount += (ElementsInMolecule[i] != 0 ? 1 : 0);
-};
-
-
-/** Counts necessary number of valence electrons and returns number and SpinType.
- * \param configuration containing everything
- */
-void molecule::CalculateOrbitals(class config &configuration)
-{
-  configuration.MaxPsiDouble = configuration.PsiMaxNoDown = configuration.PsiMaxNoUp = configuration.PsiType = 0;
-  for(int i=MAX_ELEMENTS;i--;) {
-    if (ElementsInMolecule[i] != 0) {
-      //Log() << Verbose(0) << "CalculateOrbitals: " << elemente->FindElement(i)->name << " has a valence of " << (int)elemente->FindElement(i)->Valence << " and there are " << ElementsInMolecule[i] << " of it." << endl;
-      configuration.MaxPsiDouble += ElementsInMolecule[i]*((int)elemente->FindElement(i)->Valence);
-    }
-  }
-  configuration.PsiMaxNoDown = configuration.MaxPsiDouble/2 + (configuration.MaxPsiDouble % 2);
-  configuration.PsiMaxNoUp = configuration.MaxPsiDouble/2;
-  configuration.MaxPsiDouble /= 2;
-  configuration.PsiType = (configuration.PsiMaxNoDown == configuration.PsiMaxNoUp) ? 0 : 1;
-  if ((configuration.PsiType == 1) && (configuration.ProcPEPsi < 2) && ((configuration.PsiMaxNoDown != 1) || (configuration.PsiMaxNoUp != 0))) {
-    configuration.ProcPEGamma /= 2;
-    configuration.ProcPEPsi *= 2;
-  } else {
-    configuration.ProcPEGamma *= configuration.ProcPEPsi;
-    configuration.ProcPEPsi = 1;
-  }
-  cout << configuration.PsiMaxNoDown << ">" << configuration.PsiMaxNoUp << endl;
-  if (configuration.PsiMaxNoDown > configuration.PsiMaxNoUp) {
-    configuration.InitMaxMinStopStep = configuration.MaxMinStopStep = configuration.PsiMaxNoDown;
-    cout << configuration.PsiMaxNoDown << " " << configuration.InitMaxMinStopStep << endl;
-  } else {
-    configuration.InitMaxMinStopStep = configuration.MaxMinStopStep = configuration.PsiMaxNoUp;
-    cout << configuration.PsiMaxNoUp << " " << configuration.InitMaxMinStopStep << endl;
-  }
 };
 

src/molecule.hpp

@@ -80 +80 @@
   double *PenaltyConstants;   //!<  penalty constant in front of each term
 };
-
-#define MaxThermostats 6      //!< maximum number of thermostat entries in Ions#ThermostatNames and Ions#ThermostatImplemented 
-enum thermostats { None, Woodcock, Gaussian, Langevin, Berendsen, NoseHoover };   //!< Thermostat names for output
-
 
 /** The complete molecule.
@@ -265 +261 @@
   /// Count and change present atoms' coordination.
   void CountElements();
-  void CalculateOrbitals(class config &configuration);
   bool CenterInBox();
   bool BoundInBox();
@@ -292 +287 @@
   double MinimiseConstrainedPotential(atom **&permutation, int startstep, int endstep, bool IsAngstroem);
   void EvaluateConstrainedForces(int startstep, int endstep, atom **PermutationMap, ForceMatrix *Force);
-  bool LinearInterpolationBetweenConfiguration(int startstep, int endstep, const char *prefix, config &configuration, bool MapByIdentity);
+  bool LinearInterpolationBetweenConfiguration(int startstep, int endstep, std::string &prefix, config &configuration, bool MapByIdentity);
         
   bool CheckBounds(const Vector *x) const;
@@ -324 +319 @@
 
   /// Fragment molecule by two different approaches:
-  int FragmentMolecule(int Order, config *configuration);
-  bool CheckOrderAtSite(bool *AtomMask, Graph *GlobalKeySetList, int Order, int *MinimumRingSize, char *path = NULL);
-  bool StoreBondsToFile(char *path, char *filename);
-  bool StoreAdjacencyToFile(char *path, char *filename);
-  bool CheckAdjacencyFileAgainstMolecule(char *path, atom **ListOfAtoms);
-  bool ParseOrderAtSiteFromFile(char *path);
-  bool StoreOrderAtSiteFile(char *path);
-  bool StoreForcesFile(MoleculeListClass *BondFragments, char *path, int *SortIndex);
+  int FragmentMolecule(int Order, std::string &prefix);
+  bool CheckOrderAtSite(bool *AtomMask, Graph *GlobalKeySetList, int Order, int *MinimumRingSize, std::string path = "");
+  bool StoreBondsToFile(std::string &filename, std::string path = "");
+  bool StoreAdjacencyToFile(std::string &filename, std::string path = "");
+  bool CheckAdjacencyFileAgainstMolecule(std::string &path, atom **ListOfAtoms);
+  bool ParseOrderAtSiteFromFile(std::string &path);
+  bool StoreOrderAtSiteFile(std::string &path);
+  bool StoreForcesFile(MoleculeListClass *BondFragments, std::string &path, int *SortIndex);
   bool CreateMappingLabelsToConfigSequence(int *&SortIndex);
   bool CreateFatherLookupTable(atom **&LookupTable, int count = 0);
@@ -380 +375 @@
   ~MoleculeListClass();
 
-  bool AddHydrogenCorrection(char *path);
-  bool StoreForcesFile(char *path, int *SortIndex);
+  bool AddHydrogenCorrection(std::string &path);
+  bool StoreForcesFile(std::string &path, int *SortIndex);
   void insert(molecule *mol);
   void erase(molecule *mol);
   molecule * ReturnIndex(int index);
-  bool OutputConfigForListOfFragments(config *configuration, int *SortIndex);
+  bool OutputConfigForListOfFragments(std::string &prefix, int *SortIndex);
   int NumberOfActiveMolecules();
   void Enumerate(ostream *out);

src/molecule_dynamics.cpp

@@ -18 +18 @@
 #include "parser.hpp"
 #include "Plane.hpp"
+#include "ThermoStatContainer.hpp"
 
 /************************************* Functions for class molecule *********************************/
@@ -472 +473 @@
  * \param startstep stating initial configuration in molecule::Trajectories
  * \param endstep stating final configuration in molecule::Trajectories
+ * \param &prefix path and prefix
  * \param &config configuration structure
  * \param MapByIdentity if true we just use the identity to map atoms in start config to end config, if not we find mapping by \sa MinimiseConstrainedPotential()
  * \return true - success in writing step files, false - error writing files or only one step in molecule::Trajectories
  */
-bool molecule::LinearInterpolationBetweenConfiguration(int startstep, int endstep, const char *prefix, config &configuration, bool MapByIdentity)
+bool molecule::LinearInterpolationBetweenConfiguration(int startstep, int endstep, std::string &prefix, config &configuration, bool MapByIdentity)
 {
   molecule *mol = NULL;
@@ -524 +526 @@
   for (int i=getAtomCount(); i--; )
     SortIndex[i] = i;
-  status = MoleculePerStep->OutputConfigForListOfFragments(&configuration, SortIndex);
+
+  status = MoleculePerStep->OutputConfigForListOfFragments(prefix, SortIndex);
   delete[](SortIndex);
 
@@ -643 +646 @@
 
   // calculate scale configuration
-  ScaleTempFactor = configuration.TargetTemp/ActualTemp;
+  ScaleTempFactor = configuration.Thermostats->TargetTemp/ActualTemp;
 
   // differentating between the various thermostats
@@ -651 +654 @@
       break;
      case Woodcock:
-      if ((configuration.ScaleTempStep > 0) && ((MDSteps-1) % configuration.ScaleTempStep == 0)) {
+      if ((configuration.Thermostats->ScaleTempStep > 0) && ((MDSteps-1) % configuration.Thermostats->ScaleTempStep == 0)) {
         DoLog(2) && (Log() << Verbose(2) <<  "Applying Woodcock thermostat..." << endl);
         ActOnAllAtoms( &atom::Thermostat_Woodcock, sqrt(ScaleTempFactor), MDSteps, &ekin );
@@ -684 +687 @@
       delta_alpha = 0.;
       ActOnAllAtoms( &atom::Thermostat_NoseHoover_init, MDSteps, &delta_alpha );
-      delta_alpha = (delta_alpha - (3.*getAtomCount()+1.) * configuration.TargetTemp)/(configuration.HooverMass*Units2Electronmass);
-      configuration.alpha += delta_alpha*configuration.Deltat;
-      DoLog(3) && (Log() << Verbose(3) << "alpha = " << delta_alpha << " * " << configuration.Deltat << " = " << configuration.alpha << "." << endl);
+      delta_alpha = (delta_alpha - (3.*getAtomCount()+1.) * configuration.Thermostats->TargetTemp)/(configuration.Thermostats->HooverMass*Units2Electronmass);
+      configuration.Thermostats->alpha += delta_alpha*configuration.Deltat;
+      DoLog(3) && (Log() << Verbose(3) << "alpha = " << delta_alpha << " * " << configuration.Deltat << " = " << configuration.Thermostats->alpha << "." << endl);
       // apply updated alpha as additional force
       ActOnAllAtoms( &atom::Thermostat_NoseHoover_scale, MDSteps, &ekin, &configuration );

src/molecule_fragmentation.cpp

@@ -82 +82 @@
  * -# Scans TEFactors file and sets the TEFactor of each key set in the temporary graph accordingly
  * Finally, the temporary graph is inserted into the given \a FragmentList for return.
- * \param *out output stream for debugging
- * \param *path path to file
+ * \param &path path to file
  * \param *FragmentList empty, filled on return
  * \return true - parsing successfully, false - failure on parsing (FragmentList will be NULL)
  */
-bool ParseKeySetFile(char *path, Graph *&FragmentList)
+bool ParseKeySetFile(std::string &path, Graph *&FragmentList)
 {
   bool status = true;
@@ -94 +93 @@
   GraphTestPair testGraphInsert;
   int NumberOfFragments = 0;
-  char filename[MAXSTRINGSIZE];
+  string filename;
 
   if (FragmentList == NULL) { // check list pointer
@@ -102 +101 @@
   // 1st pass: open file and read
   DoLog(1) && (Log() << Verbose(1) << "Parsing the KeySet file ... " << endl);
-  sprintf(filename, "%s/%s%s", path, FRAGMENTPREFIX, KEYSETFILE);
-  InputFile.open(filename);
-  if (InputFile != NULL) {
+  filename = path + KEYSETFILE;
+  InputFile.open(filename.c_str());
+  if (InputFile.good()) {
     // each line represents a new fragment
     char buffer[MAXSTRINGSIZE];
@@ -181 +180 @@
 
 /** Stores key sets to file.
- * \param *out output stream for debugging
  * \param KeySetList Graph with Keysets
- * \param *path path to file
+ * \param &path path to file
  * \return true - file written successfully, false - writing failed
  */
-bool StoreKeySetFile(Graph &KeySetList, char *path)
-{
-  ofstream output;
+bool StoreKeySetFile(Graph &KeySetList, std::string &path)
+{
   bool status =  true;
-  string line;
+  string line = path + KEYSETFILE;
+  ofstream output(line.c_str());
 
   // open KeySet file
-  line = path;
-  line.append("/");
-  line += FRAGMENTPREFIX;
-  line += KEYSETFILE;
-  output.open(line.c_str(), ios::out);
   DoLog(1) && (Log() << Verbose(1) << "Saving key sets of the total graph ... ");
-  if(output != NULL) {
+  if(output.good()) {
     for(Graph::iterator runner = KeySetList.begin(); runner != KeySetList.end(); runner++) {
       for (KeySet::iterator sprinter = (*runner).first.begin();sprinter != (*runner).first.end(); sprinter++) {
@@ -302 +295 @@
 
 /** Scans the adaptive order file and insert (index, value) into map.
- * \param *out output stream for debugging
- * \param *path path to ENERGYPERFRAGMENT file (may be NULL if Order is non-negative)
+ * \param &path path to ENERGYPERFRAGMENT file (may be NULL if Order is non-negative)
  * \param &IndexedKeySetList list to find key set for a given index \a No
  * \return adaptive criteria list from file
  */
-map<int, pair<double,int> > * ScanAdaptiveFileIntoMap(char *path, map<int,KeySet> &IndexKeySetList)
+map<int, pair<double,int> > * ScanAdaptiveFileIntoMap(std::string &path, map<int,KeySet> &IndexKeySetList)
 {
   map<int, pair<double,int> > *AdaptiveCriteriaList = new map<int, pair<double,int> >;
@@ -313 +305 @@
   double Value = 0.;
   char buffer[MAXSTRINGSIZE];
-  sprintf(buffer, "%s/%s%s.dat", path, FRAGMENTPREFIX, ENERGYPERFRAGMENT);
-  ifstream InputFile(buffer, ios::in);
+  string filename = path + ENERGYPERFRAGMENT;
+  ifstream InputFile(filename.c_str());
+
+  if (InputFile.fail()) {
+    DoeLog(1) && (eLog() << Verbose(1) << "Cannot find file " << filename << "." << endl);
+    return AdaptiveCriteriaList;
+  }
 
   if (CountLinesinFile(InputFile) > 0) {
@@ -419 +416 @@
 
 /** Checks whether the OrderAtSite is still below \a Order at some site.
- * \param *out output stream for debugging
  * \param *AtomMask defines true/false per global Atom::nr to mask in/out each nuclear site, used to activate given number of site to increment order adaptively
  * \param *GlobalKeySetList list of keysets with global ids (valid in "this" molecule) needed for adaptive increase
  * \param Order desired Order if positive, desired exponent in threshold criteria if negative (0 is single-step)
  * \param *MinimumRingSize array of max. possible order to avoid loops
- * \param *path path to ENERGYPERFRAGMENT file (may be NULL if Order is non-negative)
+ * \param path path to ENERGYPERFRAGMENT file (may be NULL if Order is non-negative)
  * \return true - needs further fragmentation, false - does not need fragmentation
  */
-bool molecule::CheckOrderAtSite(bool *AtomMask, Graph *GlobalKeySetList, int Order, int *MinimumRingSize, char *path)
+bool molecule::CheckOrderAtSite(bool *AtomMask, Graph *GlobalKeySetList, int Order, int *MinimumRingSize, std::string path)
 {
   bool status = false;
@@ -585 +581 @@
  * of vertex indices: Global always means the index in "this" molecule, whereas local refers to the molecule or
  * subgraph in the MoleculeListClass.
- * \param *out output stream for debugging
  * \param Order up to how many neighbouring bonds a fragment contains in BondOrderScheme::BottumUp scheme
- * \param *configuration configuration for writing config files for each fragment
+ * \param &prefix path and prefix of the bond order configs to be written
  * \return 1 - continue, 2 - stop (no fragmentation occured)
  */
-int molecule::FragmentMolecule(int Order, config *configuration)
+int molecule::FragmentMolecule(int Order, std::string &prefix)
 {
   MoleculeListClass *BondFragments = NULL;
@@ -624 +619 @@
 
   // === compare it with adjacency file ===
-  FragmentationToDo = FragmentationToDo && CheckAdjacencyFileAgainstMolecule(configuration->configpath, ListOfAtoms);
+  FragmentationToDo = FragmentationToDo && CheckAdjacencyFileAgainstMolecule(prefix, ListOfAtoms);
   delete[](ListOfAtoms);
 
@@ -658 +653 @@
 
   // ===== 3. if structure still valid, parse key set file and others =====
-  FragmentationToDo = FragmentationToDo && ParseKeySetFile(configuration->configpath, ParsedFragmentList);
+  FragmentationToDo = FragmentationToDo && ParseKeySetFile(prefix, ParsedFragmentList);
 
   // ===== 4. check globally whether there's something to do actually (first adaptivity check)
-  FragmentationToDo = FragmentationToDo && ParseOrderAtSiteFromFile(configuration->configpath);
+  FragmentationToDo = FragmentationToDo && ParseOrderAtSiteFromFile(prefix);
 
   // =================================== Begin of FRAGMENTATION ===============================
@@ -672 +667 @@
   AtomMask[getAtomCount()] = false;
   FragmentationToDo = false;  // if CheckOrderAtSite just ones recommends fragmentation, we will save fragments afterwards
-  while ((CheckOrder = CheckOrderAtSite(AtomMask, ParsedFragmentList, Order, MinimumRingSize, configuration->configpath))) {
+  while ((CheckOrder = CheckOrderAtSite(AtomMask, ParsedFragmentList, Order, MinimumRingSize, prefix))) {
     FragmentationToDo = FragmentationToDo || CheckOrder;
     AtomMask[getAtomCount()] = true;   // last plus one entry is used as marker that we have been through this loop once already in CheckOrderAtSite()
@@ -727 +722 @@
     KeySet test = (*runner).first;
     DoLog(0) && (Log() << Verbose(0) << "Fragment No." << (*runner).second.first << " with TEFactor " << (*runner).second.second << "." << endl);
-    BondFragments->insert(StoreFragmentFromKeySet(test, configuration));
+    BondFragments->insert(StoreFragmentFromKeySet(test, World::getInstance().getConfig()));
     k++;
   }
@@ -739 +734 @@
 
     DoLog(1) && (Log() << Verbose(1) << "Writing " << BondFragments->ListOfMolecules.size() << " possible bond fragmentation configs" << endl);
-    if (BondFragments->OutputConfigForListOfFragments(configuration, SortIndex))
+    if (BondFragments->OutputConfigForListOfFragments(prefix, SortIndex))
       DoLog(1) && (Log() << Verbose(1) << "All configs written." << endl);
     else
@@ -745 +740 @@
 
     // store force index reference file
-    BondFragments->StoreForcesFile(configuration->configpath, SortIndex);
+    BondFragments->StoreForcesFile(prefix, SortIndex);
 
     // store keysets file
-    StoreKeySetFile(TotalGraph, configuration->configpath);
+    StoreKeySetFile(TotalGraph, prefix);
 
     {
       // store Adjacency file
-      char filename[MAXSTRINGSIZE];
-      strcpy(filename, FRAGMENTPREFIX);
-      strcat(filename, ADJACENCYFILE);
-      StoreAdjacencyToFile(configuration->configpath, filename);
+      std::string filename = prefix + ADJACENCYFILE;
+      StoreAdjacencyToFile(filename);
     }
 
     // store Hydrogen saturation correction file
-    BondFragments->AddHydrogenCorrection(configuration->configpath);
+    BondFragments->AddHydrogenCorrection(prefix);
 
     // store adaptive orders into file
-    StoreOrderAtSiteFile(configuration->configpath);
+    StoreOrderAtSiteFile(prefix);
 
     // restore orbital and Stop values
-    CalculateOrbitals(*configuration);
+    //CalculateOrbitals(*configuration);
 
     // free memory for bond part
@@ -782 +775 @@
 /** Stores pairs (Atom::nr, Atom::AdaptiveOrder) into file.
  * Atoms not present in the file get "-1".
- * \param *out output stream for debugging
- * \param *path path to file ORDERATSITEFILE
+ * \param &path path to file ORDERATSITEFILE
  * \return true - file writable, false - not writable
  */
-bool molecule::StoreOrderAtSiteFile(char *path)
-{
-  stringstream line;
+bool molecule::StoreOrderAtSiteFile(std::string &path)
+{
+  string line;
   ofstream file;
 
-  line << path << "/" << FRAGMENTPREFIX << ORDERATSITEFILE;
-  file.open(line.str().c_str());
+  line = path + ORDERATSITEFILE;
+  file.open(line.c_str());
   DoLog(1) && (Log() << Verbose(1) << "Writing OrderAtSite " << ORDERATSITEFILE << " ... " << endl);
-  if (file != NULL) {
+  if (file.good()) {
     ActOnAllAtoms( &atom::OutputOrder, &file );
     file.close();
@@ -800 +792 @@
     return true;
   } else {
-    DoLog(1) && (Log() << Verbose(1) << "failed to open file " << line.str() << "." << endl);
+    DoLog(1) && (Log() << Verbose(1) << "failed to open file " << line << "." << endl);
     return false;
   }
@@ -807 +799 @@
 /** Parses pairs(Atom::nr, Atom::AdaptiveOrder) from file and stores in molecule's Atom's.
  * Atoms not present in the file get "0".
- * \param *out output stream for debugging
- * \param *path path to file ORDERATSITEFILEe
+ * \param &path path to file ORDERATSITEFILEe
  * \return true - file found and scanned, false - file not found
  * \sa ParseKeySetFile() and CheckAdjacencyFileAgainstMolecule() as this is meant to be used in conjunction with the two
  */
-bool molecule::ParseOrderAtSiteFromFile(char *path)
+bool molecule::ParseOrderAtSiteFromFile(std::string &path)
 {
   unsigned char *OrderArray = new unsigned char[getAtomCount()];
@@ -818 +809 @@
   bool status;
   int AtomNr, value;
-  stringstream line;
+  string line;
   ifstream file;
 
@@ -827 +818 @@
 
   DoLog(1) && (Log() << Verbose(1) << "Begin of ParseOrderAtSiteFromFile" << endl);
-  line << path << "/" << FRAGMENTPREFIX << ORDERATSITEFILE;
-  file.open(line.str().c_str());
-  if (file != NULL) {
+  line = path + ORDERATSITEFILE;
+  file.open(line.c_str());
+  if (file.good()) {
     while (!file.eof()) { // parse from file
       AtomNr = -1;
@@ -850 +841 @@
     status = true;
   } else {
-    DoLog(1) && (Log() << Verbose(1) << "\t ... failed to open file " << line.str() << "." << endl);
+    DoLog(1) && (Log() << Verbose(1) << "\t ... failed to open file " << line << "." << endl);
     status = false;
   }

src/molecule_graph.cpp

@@ -12 +12 @@
 #include "bondgraph.hpp"
 #include "config.hpp"
+#include "defs.hpp"
 #include "element.hpp"
 #include "helpers.hpp"
@@ -1024 +1025 @@
 /** Storing the bond structure of a molecule to file.
  * Simply stores Atom::nr and then the Atom::nr of all bond partners per line.
- * \param *path path to file
- * \param *filename name of file
+ * \param &filename name of file
+ * \param path path to file, defaults to empty
  * \return true - file written successfully, false - writing failed
  */
-bool molecule::StoreAdjacencyToFile(char *path, char *filename)
+bool molecule::StoreAdjacencyToFile(std::string &filename, std::string path)
 {
   ofstream AdjacencyFile;
-  stringstream line;
+  string line;
   bool status = true;
 
-  if (path != NULL)
-    line << path << "/" << filename;
+  if (path != "")
+    line = path + "/" + filename;
   else
-    line << filename;
-  AdjacencyFile.open(line.str().c_str(), ios::out);
+    line = filename;
+  AdjacencyFile.open(line.c_str(), ios::out);
   DoLog(1) && (Log() << Verbose(1) << "Saving adjacency list ... " << endl);
-  if (AdjacencyFile != NULL) {
+  if (AdjacencyFile.good()) {
     AdjacencyFile << "m\tn" << endl;
     ActOnAllAtoms(&atom::OutputAdjacency, &AdjacencyFile);
@@ -1046 +1047 @@
     DoLog(1) && (Log() << Verbose(1) << "\t... done." << endl);
   } else {
-    DoLog(1) && (Log() << Verbose(1) << "\t... failed to open file " << line.str() << "." << endl);
+    DoLog(1) && (Log() << Verbose(1) << "\t... failed to open file " << line << "." << endl);
     status = false;
   }
@@ -1056 +1057 @@
 /** Storing the bond structure of a molecule to file.
  * Simply stores Atom::nr and then the Atom::nr of all bond partners, one per line.
- * \param *path path to file
- * \param *filename name of file
+ * \param &filename name of file
+ * \param path path to file, defaults to empty
  * \return true - file written successfully, false - writing failed
  */
-bool molecule::StoreBondsToFile(char *path, char *filename)
+bool molecule::StoreBondsToFile(std::string &filename, std::string path)
 {
   ofstream BondFile;
-  stringstream line;
+  string line;
   bool status = true;
 
-  if (path != NULL)
-    line << path << "/" << filename;
+  if (path != "")
+    line = path + "/" + filename;
   else
-    line << filename;
-  BondFile.open(line.str().c_str(), ios::out);
+    line = filename;
+  BondFile.open(line.c_str(), ios::out);
   DoLog(1) && (Log() << Verbose(1) << "Saving adjacency list ... " << endl);
-  if (BondFile != NULL) {
+  if (BondFile.good()) {
     BondFile << "m\tn" << endl;
     ActOnAllAtoms(&atom::OutputBonds, &BondFile);
@@ -1078 +1079 @@
     DoLog(1) && (Log() << Verbose(1) << "\t... done." << endl);
   } else {
-    DoLog(1) && (Log() << Verbose(1) << "\t... failed to open file " << line.str() << "." << endl);
+    DoLog(1) && (Log() << Verbose(1) << "\t... failed to open file " << line << "." << endl);
     status = false;
   }
@@ -1086 +1087 @@
 ;
 
-bool CheckAdjacencyFileAgainstMolecule_Init(char *path, ifstream &File, int *&CurrentBonds)
-{
-  stringstream filename;
-  filename << path << "/" << FRAGMENTPREFIX << ADJACENCYFILE;
-  File.open(filename.str().c_str(), ios::out);
+bool CheckAdjacencyFileAgainstMolecule_Init(std::string &path, ifstream &File, int *&CurrentBonds)
+{
+  string filename;
+  filename = path + ADJACENCYFILE;
+  File.open(filename.c_str(), ios::out);
   DoLog(1) && (Log() << Verbose(1) << "Looking at bond structure stored in adjacency file and comparing to present one ... " << endl);
-  if (File == NULL)
+  if (File.fail())
     return false;
 
@@ -1146 +1147 @@
  * \return true - structure is equal, false - not equivalence
  */
-bool molecule::CheckAdjacencyFileAgainstMolecule(char *path, atom **ListOfAtoms)
+bool molecule::CheckAdjacencyFileAgainstMolecule(std::string &path, atom **ListOfAtoms)
 {
   ifstream File;

src/moleculelist.cpp

@@ -12 +12 @@
 #include "atom.hpp"
 #include "bond.hpp"
+#include "bondgraph.hpp"
 #include "boundary.hpp"
 #include "config.hpp"
@@ -379 +380 @@
  * bonded to the same atom, then we add for this pair a correction term constructed from a Morse
  * potential function fit to QM calculations with respecting to the interatomic hydrogen distance.
- * \param *out output stream for debugging
- * \param *path path to file
- */
-bool MoleculeListClass::AddHydrogenCorrection(char *path)
+ * \param &path path to file
+ */
+bool MoleculeListClass::AddHydrogenCorrection(std::string &path)
 {
   bond *Binder = NULL;
@@ -400 +400 @@
   // 0a. find dimension of matrices with constants
   line = path;
-  line.append("/");
-  line += FRAGMENTPREFIX;
   line += "1";
   line += FITCONSTANTSUFFIX;
   input.open(line.c_str());
-  if (input == NULL) {
+  if (input.fail()) {
     DoLog(1) && (Log() << Verbose(1) << endl << "Unable to open " << line << ", is the directory correct?" << endl);
     return false;
@@ -569 +567 @@
 
 /** Store force indices, i.e. the connection between the nuclear index in the total molecule config and the respective atom in fragment config.
- * \param *out output stream for debugging
- * \param *path path to file
+ * \param &path path to file
  * \param *SortIndex Index to map from the BFS labeling to the sequence how of Ion_Type in the config
  * \return true - file written successfully, false - writing failed
  */
-bool MoleculeListClass::StoreForcesFile(char *path,
-    int *SortIndex)
+bool MoleculeListClass::StoreForcesFile(std::string &path, int *SortIndex)
 {
   bool status = true;
-  ofstream ForcesFile;
-  stringstream line;
+  string filename(path);
+  filename += FORCESFILE;
+  ofstream ForcesFile(filename.c_str());
   periodentafel *periode=World::getInstance().getPeriode();
 
   // open file for the force factors
   DoLog(1) && (Log() << Verbose(1) << "Saving  force factors ... ");
-  line << path << "/" << FRAGMENTPREFIX << FORCESFILE;
-  ForcesFile.open(line.str().c_str(), ios::out);
-  if (ForcesFile != NULL) {
+  if (!ForcesFile.fail()) {
     //Log() << Verbose(1) << "Final AtomicForcesList: ";
     //output << prefix << "Forces" << endl;
@@ -611 +606 @@
   } else {
     status = false;
-    DoLog(1) && (Log() << Verbose(1) << "failed to open file " << line.str() << "." << endl);
+    DoLog(1) && (Log() << Verbose(1) << "failed to open file " << filename << "." << endl);
   }
   ForcesFile.close();
@@ -620 +615 @@
 /** Writes a config file for each molecule in the given \a **FragmentList.
  * \param *out output stream for debugging
- * \param *configuration standard configuration to attach atoms in fragment molecule to.
+ * \param &prefix path and prefix to the fragment config files
  * \param *SortIndex Index to map from the BFS labeling to the sequence how of Ion_Type in the config
  * \return true - success (each file was written), false - something went wrong.
  */
-bool MoleculeListClass::OutputConfigForListOfFragments(config *configuration, int *SortIndex)
+bool MoleculeListClass::OutputConfigForListOfFragments(std::string &prefix, int *SortIndex)
 {
   ofstream outputFragment;
-  char FragmentName[MAXSTRINGSIZE];
+  std::string FragmentName;
   char PathBackup[MAXSTRINGSIZE];
   bool result = true;
@@ -645 +640 @@
   for (MoleculeList::iterator ListRunner = ListOfMolecules.begin(); ListRunner != ListOfMolecules.end(); ListRunner++) {
     // save default path as it is changed for each fragment
-    path = configuration->GetDefaultPath();
+    path = World::getInstance().getConfig()->GetDefaultPath();
     if (path != NULL)
       strcpy(PathBackup, path);
@@ -658 +653 @@
     // output xyz file
     FragmentNumber = FixedDigitNumber(ListOfMolecules.size(), FragmentCounter++);
-    sprintf(FragmentName, "%s/%s%s.conf.xyz", configuration->configpath, FRAGMENTPREFIX, FragmentNumber);
-    outputFragment.open(FragmentName, ios::out);
+    FragmentName = prefix + FragmentNumber + ".conf.xyz";
+    outputFragment.open(FragmentName.c_str(), ios::out);
     DoLog(2) && (Log() << Verbose(2) << "Saving bond fragment No. " << FragmentNumber << "/" << FragmentCounter - 1 << " as XYZ ...");
     if ((intermediateResult = (*ListRunner)->OutputXYZ(&outputFragment)))
@@ -682 +677 @@
     for (int k = 0; k < NDIM; k++) {
       j += k + 1;
-      BoxDimension[k] = 2.5 * (configuration->GetIsAngstroem() ? 1. : 1. / AtomicLengthToAngstroem);
+      BoxDimension[k] = 2.5 * (World::getInstance().getConfig()->GetIsAngstroem() ? 1. : 1. / AtomicLengthToAngstroem);
       cell_size[j] = BoxDimension[k] * 2.;
     }
@@ -689 +684 @@
     // also calculate necessary orbitals
     (*ListRunner)->CountElements(); // this is a bugfix, atoms should shoulds actually be added correctly to this fragment
-    (*ListRunner)->CalculateOrbitals(*configuration);
+    //(*ListRunner)->CalculateOrbitals(*World::getInstance().getConfig);
 
     // change path in config
-    //strcpy(PathBackup, configuration->configpath);
-    sprintf(FragmentName, "%s/%s%s/", PathBackup, FRAGMENTPREFIX, FragmentNumber);
-    configuration->SetDefaultPath(FragmentName);
+    FragmentName = PathBackup;
+    FragmentName += "/";
+    FragmentName += FRAGMENTPREFIX;
+    FragmentName += FragmentNumber;
+    FragmentName += "/";
+    World::getInstance().getConfig()->SetDefaultPath(FragmentName.c_str());
 
     // and save as config
-    sprintf(FragmentName, "%s/%s%s.conf", configuration->configpath, FRAGMENTPREFIX, FragmentNumber);
+    FragmentName = prefix + FragmentNumber + ".conf";
     DoLog(2) && (Log() << Verbose(2) << "Saving bond fragment No. " << FragmentNumber << "/" << FragmentCounter - 1 << " as config ...");
-    if ((intermediateResult = configuration->Save(FragmentName, (*ListRunner)->elemente, (*ListRunner))))
+    if ((intermediateResult = World::getInstance().getConfig()->Save(FragmentName.c_str(), (*ListRunner)->elemente, (*ListRunner))))
       DoLog(0) && (Log() << Verbose(0) << " done." << endl);
     else
@@ -706 +704 @@
 
     // restore old config
-    configuration->SetDefaultPath(PathBackup);
+    World::getInstance().getConfig()->SetDefaultPath(PathBackup);
 
     // and save as mpqc input file
-    sprintf(FragmentName, "%s/%s%s.conf", configuration->configpath, FRAGMENTPREFIX, FragmentNumber);
+    FragmentName = prefix + FragmentNumber + ".conf";
     DoLog(2) && (Log() << Verbose(2) << "Saving bond fragment No. " << FragmentNumber << "/" << FragmentCounter - 1 << " as mpqc input ...");
-    if ((intermediateResult = configuration->SaveMPQC(FragmentName, (*ListRunner))))
+    if ((intermediateResult = World::getInstance().getConfig()->SaveMPQC(FragmentName.c_str(), (*ListRunner))))
       DoLog(2) && (Log() << Verbose(2) << " done." << endl);
     else
@@ -767 +765 @@
 
   // 1. dissect the molecule into connected subgraphs
-  if (!configuration->BG->ConstructBondGraph(mol)) {
-    World::getInstance().destroyMolecule(mol);
-    DoeLog(1) && (eLog()<< Verbose(1) << "There are no bonds." << endl);
+  if (configuration->BG != NULL) {
+    if (!configuration->BG->ConstructBondGraph(mol)) {
+      World::getInstance().destroyMolecule(mol);
+      DoeLog(1) && (eLog()<< Verbose(1) << "There are no bonds." << endl);
+      return;
+    }
+  } else {
+    DoeLog(1) && (eLog()<< Verbose(1) << "There is no BondGraph class present to create bonds." << endl);
     return;
   }

src/periodentafel.cpp

@@ -32 +32 @@
 periodentafel::periodentafel()
 {
-  bool status = true;
-  status = LoadElementsDatabase(new stringstream(elementsDB,ios_base::in));
-  ASSERT(status,  "General element initialization failed");
-  status = LoadValenceDatabase(new stringstream(valenceDB,ios_base::in));
-  ASSERT(status, "Valence entry of element initialization failed");
-  status = LoadOrbitalsDatabase(new stringstream(orbitalsDB,ios_base::in));
-  ASSERT(status, "Orbitals entry of element initialization failed");
-  status = LoadHBondAngleDatabase(new stringstream(HbondangleDB,ios_base::in));
-  ASSERT(status, "HBond angle entry of element initialization failed");
-  status = LoadHBondLengthsDatabase(new stringstream(HbonddistanceDB,ios_base::in));
-  ASSERT(status, "HBond distance entry of element initialization failed");
+  {
+    stringstream input(elementsDB,ios_base::in);
+    bool status = LoadElementsDatabase(&input);
+    ASSERT(status,  "General element initialization failed");
+  }
+  {
+    stringstream input(valenceDB,ios_base::in);
+    bool status = LoadValenceDatabase(&input);
+    ASSERT(status, "Valence entry of element initialization failed");
+  }
+  {
+    stringstream input(orbitalsDB,ios_base::in);
+    bool status = LoadOrbitalsDatabase(&input);
+    ASSERT(status, "Orbitals entry of element initialization failed");
+  }
+  {
+    stringstream input(HbondangleDB,ios_base::in);
+    bool status = LoadHBondAngleDatabase(&input);
+    ASSERT(status, "HBond angle entry of element initialization failed");
+  }
+  {
+    stringstream input(HbonddistanceDB,ios_base::in);
+    bool status = LoadHBondLengthsDatabase(&input);
+    ASSERT(status, "HBond distance entry of element initialization failed");
+  }
 };
 
@@ -333 +347 @@
           ASSERT(Elemental != NULL, "element should be present but is not??");
           *Elemental = *neues;
+          delete(neues);
+          neues = Elemental;
         } else {
           InserterTest = elements.insert(pair <atomicNumber_t,element*> (neues->getNumber(), neues));

src/tesselation.cpp

@@ -9 +9 @@
 
 #include <fstream>
-#include <assert.h>
 
 #include "helpers.hpp"
@@ -24 +23 @@
 #include "Plane.hpp"
 #include "Exceptions/LinearDependenceException.hpp"
-#include "Helpers/Assert.hpp"
-
 #include "Helpers/Assert.hpp"
 
@@ -2527 +2524 @@
     baseline = Runner->second;
     if (baseline->pointlist.empty()) {
-      assert((baseline->BaseLine->triangles.size() == 1) && ("Open line without exactly one attached triangle"));
+      ASSERT((baseline->BaseLine->triangles.size() == 1),"Open line without exactly one attached triangle");
       T = (((baseline->BaseLine->triangles.begin()))->second);
       DoLog(1) && (Log() << Verbose(1) << "Finding best candidate for open line " << *baseline->BaseLine << " of triangle " << *T << endl);

src/unittests/LinkedCellUnitTest.cpp

@@ -264 +264 @@
   Vector tester;
   LinkedCell::LinkedNodes *ListOfPoints = NULL;
-  atom *Walker = NULL;
   size_t size = 0;
 
@@ -326 +325 @@
   Vector tester;
   LinkedCell::LinkedNodes *ListOfPoints = NULL;
-  atom *Walker = NULL;
   size_t size = 0;
 

src/unittests/Makefile.am

@@ -52 +52 @@
 GSLLIBS = ../libgslwrapper.a $(BOOST_LIB) ${BOOST_THREAD_LIB}
 ALLLIBS = ../libmolecuilder.a ${GSLLIBS} 
+PARSERLIBS = ../libparser.a ${ALLLIBS} 
 
 TESTSOURCES = \
@@ -201 +202 @@
 
 ParserUnitTest_SOURCES = UnitTestMain.cpp ParserUnitTest.cpp ParserUnitTest.hpp
-ParserUnitTest_LDADD = ${ALLLIBS}
+ParserUnitTest_LDADD = ${PARSERLIBS}
 
 periodentafelTest_SOURCES = UnitTestMain.cpp periodentafelTest.cpp periodentafelTest.hpp
@@ -210 +211 @@
 
 SingletonTest_SOURCES = UnitTestMain.cpp SingletonTest.cpp SingletonTest.hpp
-SingletonTest_LDADD = $(BOOST_LIB) ${BOOST_THREAD_LIB}
+SingletonTest_LDADD = ${ALLLIBS} $(BOOST_LIB) ${BOOST_THREAD_LIB}
 
 StackClassUnitTest_SOURCES = UnitTestMain.cpp stackclassunittest.cpp stackclassunittest.hpp

src/unittests/ObserverTest.cpp

@@ -382 +382 @@
   // make this Observable its own subject. NEVER DO THIS IN ACTUAL CODE
   simpleObservable1->signOn(simpleObservable1);
+#ifndef NDEBUG
   CPPUNIT_ASSERT_THROW(simpleObservable1->changeMethod(),Assert::AssertionFailure);
+#else
+  simpleObservable1->changeMethod();
+#endif
 
   // more complex test
@@ -388 +392 @@
   simpleObservable1->signOn(simpleObservable2);
   simpleObservable2->signOn(simpleObservable1);
+#ifndef NDEBUG
   CPPUNIT_ASSERT_THROW(simpleObservable1->changeMethod(),Assert::AssertionFailure);
+#else
+  simpleObservable1->changeMethod();
+#endif
+
+
   simpleObservable1->signOff(simpleObservable2);
   simpleObservable2->signOff(simpleObservable1);

src/unittests/ParserUnitTest.cpp

@@ -12 +12 @@
 #include <cppunit/ui/text/TestRunner.h>
 
+#include "Parser/MpqcParser.hpp"
+#include "Parser/PcpParser.hpp"
+#include "Parser/TremoloParser.hpp"
 #include "Parser/XyzParser.hpp"
-#include "Parser/TremoloParser.hpp"
 #include "World.hpp"
 #include "atom.hpp"
@@ -29 +31 @@
 CPPUNIT_TEST_SUITE_REGISTRATION( ParserUnitTest );
 
+static string waterPcp = "# ParallelCarParinello - main configuration file - created with molecuilder\n\
+\n\
+mainname\tpcp\t# programm name (for runtime files)\n\
+defaultpath\not specified\t# where to put files during runtime\n\
+pseudopotpath\not specified\t# where to find pseudopotentials\n\
+\n\
+ProcPEGamma\t8\t# for parallel computing: share constants\n\
+ProcPEPsi\t1\t# for parallel computing: share wave functions\n\
+DoOutVis\t0\t# Output data for OpenDX\n\
+DoOutMes\t1\t# Output data for measurements\n\
+DoOutOrbitals\t0\t# Output all Orbitals\n\
+DoOutCurr\t0\t# Ouput current density for OpenDx\n\
+DoOutNICS\t0\t# Output Nucleus independent current shieldings\n\
+DoPerturbation\t0\t# Do perturbation calculate and determine susceptibility and shielding\n\
+DoFullCurrent\t0\t# Do full perturbation\n\
+DoConstrainedMD\t0\t# Do perform a constrained (>0, relating to current MD step) instead of unconstrained (0) MD\n\
+Thermostat\tBerendsen\t2.5\t# Which Thermostat and its parameters to use in MD case.\n\
+CommonWannier\t0\t# Put virtual centers at indivual orbits, all common, merged by variance, to grid point, to cell center\n\
+SawtoothStart\t0.01\t# Absolute value for smooth transition at cell border \n\
+VectorPlane\t0\t# Cut plane axis (x, y or z: 0,1,2) for two-dim current vector plot\n\
+VectorCut\t0\t# Cut plane axis value\n\
+AddGramSch\t1\t# Additional GramSchmidtOrtogonalization to be safe\n\
+Seed\t1\t# initial value for random seed for Psi coefficients\n\
+\n\
+MaxOuterStep\t0\t# number of MolecularDynamics/Structure optimization steps\n\
+Deltat\t0.01\t# time per MD step\n\
+OutVisStep\t10\t# Output visual data every ...th step\n\
+OutSrcStep\t5\t# Output \"restart\" data every ..th step\n\
+TargetTemp\t0.000950045\t# Target temperature\n\
+MaxPsiStep\t3\t# number of Minimisation steps per state (0 - default)\n\
+EpsWannier\t1e-07\t# tolerance value for spread minimisation of orbitals\n\
+# Values specifying when to stop\n\
+MaxMinStep\t100\t# Maximum number of steps\n\
+RelEpsTotalE\t1e-07\t# relative change in total energy\n\
+RelEpsKineticE\t1e-05\t# relative change in kinetic energy\n\
+MaxMinStopStep\t2\t# check every ..th steps\n\
+MaxMinGapStopStep\t1\t# check every ..th steps\n\
+\n\
+# Values specifying when to stop for INIT, otherwise same as above\n\
+MaxInitMinStep\t100\t# Maximum number of steps\n\
+InitRelEpsTotalE\t1e-05\t# relative change in total energy\n\
+InitRelEpsKineticE\t0.0001\t# relative change in kinetic energy\n\
+InitMaxMinStopStep\t2\t# check every ..th steps\n\
+InitMaxMinGapStopStep\t1\t# check every ..th steps\n\
+\n\
+BoxLength\t# (Length of a unit cell)\n\
+20\n\
+0\t20\n\
+0\t0\t20\n\
+\n\
+ECut\t128\t# energy cutoff for discretization in Hartrees\n\
+MaxLevel\t5\t# number of different levels in the code, >=2\n\
+Level0Factor\t2\t# factor by which node number increases from S to 0 level\n\
+RiemannTensor\t0\t# (Use metric)\n\
+PsiType\t0\t# 0 - doubly occupied, 1 - SpinUp,SpinDown\n\
+MaxPsiDouble\t2\t# here: specifying both maximum number of SpinUp- and -Down-states\n\
+PsiMaxNoUp\t2\t# here: specifying maximum number of SpinUp-states\n\
+PsiMaxNoDown\t2\t# here: specifying maximum number of SpinDown-states\n\
+AddPsis\t0\t# Additional unoccupied Psis for bandgap determination\n\
+\n\
+RCut\t20\t# R-cut for the ewald summation\n\
+StructOpt\t0\t# Do structure optimization beforehand\n\
+IsAngstroem\t1\t# 0 - Bohr, 1 - Angstroem\n\
+RelativeCoord\t0\t# whether ion coordinates are relative (1) or absolute (0)\n\
+MaxTypes\t2\t# maximum number of different ion types\n\
+\n\
+# Ion type data (PP = PseudoPotential, Z = atomic number)\n\
+#Ion_TypeNr.\tAmount\tZ\tRGauss\tL_Max(PP)L_Loc(PP)IonMass\t# chemical name, symbol\n\
+Ion_Type1\t2\t1\t1.0\t3\t3\t1.008\tHydrogen\tH\n\
+Ion_Type2\t1\t8\t1.0\t3\t3\t15.999\tOxygen\tO\n\
+#Ion_TypeNr._Nr.R[0]\tR[1]\tR[2]\tMoveType (0 MoveIon, 1 FixedIon)\n\
+Ion_Type2_1\t0.000000000\t0.000000000\t0.000000000\t0 # molecule nr 0\n\
+Ion_Type1_1\t0.758602\t0.000000000\t0.504284\t0 # molecule nr 1\n\
+Ion_Type1_2\t0.758602\t0.000000000\t-0.504284\t0 # molecule nr 2\n";
+static string waterMpqc ="% Created by MoleCuilder\n\
+mpqc: (\n\
+\tsavestate = no\n\
+\tdo_gradient = yes\n\
+\tmole<MBPT2>: (\n\
+\t\tmaxiter = 200\n\
+\t\tbasis = $:basis\n\
+\t\tmolecule = $:molecule\n\
+\t\treference<CLHF>: (\n\
+\t\t\tbasis = $:basis\n\
+\t\t\tmolecule = $:molecule\n\
+\t\t)\n\
+\t)\n\
+)\n\
+molecule<Molecule>: (\n\
+\tunit = angstrom\n\
+\t{ atoms geometry } = {\n\
+\t\tO [ -0.505735\t0\t0 ]\n\
+\t\tH [ 0.252867\t0\t0.504284 ]\n\
+\t\tH [ 0.252867\t0\t-0.504284 ]\n\
+\t}\n\
+)\n\
+basis<GaussianBasisSet>: (\n\
+\tname = \"3-21G\"\n\
+\tmolecule = $:molecule\n\
+)\n";
+static string waterXyz = "3\nH2O: water molecule\nO\t0\t0\t0\nH\t0.758602\t0\t0.504284\nH\t0.758602\t0\t-0.504284\n";
+static string Tremolo_Atomdata1 = "# ATOMDATA\tId\tname\tType\tx=3\n";
+static string Tremolo_Atomdata2 = "#\n#ATOMDATA Id name Type x=3\n1 hydrogen H 3.0 4.5 0.1\n\n";
+static string Tremolo_invalidkey = "#\n#ATOMDATA Id name foo Type x=3\n\n\n";
+static string Tremolo_velocity = "#\n#ATOMDATA Id name Type u=3\n1 hydrogen H 3.0 4.5 0.1\n\n";
+static string Tremolo_neighbours = "#\n#ATOMDATA Id Type neighbors=2\n1 H 3 0\n2 H 3 0\n3 O 1 2\n";
+static string Tremolo_improper = "#\n#ATOMDATA Id Type imprData\n8 H 9-10\n9 H 10-8,8-10\n10 O -\n";
+static string Tremolo_torsion = "#\n#ATOMDATA Id Type torsion\n8 H 9-10\n9 H 10-8,8-10\n10 O -\n";
+static string Tremolo_full = "# ATOMDATA\tx=3\tu=3\tF\tstress\tId\tneighbors=5\timprData\tGroupMeasureTypeNo\tType\textType\tname\tresName\tchainID\tresSeq\toccupancy\ttempFactor\tsegID\tCharge\tcharge\tGrpTypeNo\ttorsion\n0\t0\t0\t0\t0\t0\t0\t0\t0\t0\t0\t0\t0\t0\t-\t0\tH\t-\t-\t-\t0\t0\t0\t0\t0\t0\t0\t0\t-\t\n";
 
 void ParserUnitTest::setUp() {
   World::getInstance();
+
+  // we need hydrogens and oxygens in the following tests
+  CPPUNIT_ASSERT(World::getInstance().getPeriode()->FindElement(1) != NULL);
+  CPPUNIT_ASSERT(World::getInstance().getPeriode()->FindElement(8) != NULL);
 }
 
 void ParserUnitTest::tearDown() {
+  ChangeTracker::purgeInstance();
   World::purgeInstance();
 }
@@ -43 +159 @@
   cout << "Testing the XYZ parser." << endl;
   XyzParser* testParser = new XyzParser();
-  string waterXyz = "3\nH2O: water molecule\nO\t0.000000\t0.000000\t0.000000\nH\t0.758602\t0.000000\t0.504284\nH\t0.758602\t0.000000\t-0.504284\n";
   stringstream input;
   input << waterXyz;
@@ -62 +177 @@
   TremoloParser* testParser = new TremoloParser();
   stringstream input, output;
-  string waterTremolo;
 
   // Atomdata beginning with "# ATOMDATA"
-  waterTremolo = "# ATOMDATA\tId\tname\tType\tx=3\n";
-  input << waterTremolo;
-  testParser->load(&input);
-  testParser->save(&output);
-  CPPUNIT_ASSERT(waterTremolo == output.str());
+  input << Tremolo_Atomdata1;
+  testParser->load(&input);
+  testParser->save(&output);
+  CPPUNIT_ASSERT(Tremolo_Atomdata1 == output.str());
   input.clear();
   output.clear();
 
   // Atomdata beginning with "#ATOMDATA"
-  waterTremolo = "#\n#ATOMDATA Id name Type x=3\n1 hydrogen H 3.0 4.5 0.1\n\n";
-  input << waterTremolo;
+  input << Tremolo_Atomdata2;
   testParser->load(&input);
   testParser->save(&output);
@@ -83 +195 @@
 
   // Invalid key in Atomdata line
-  waterTremolo = "#\n#ATOMDATA Id name foo Type x=3\n\n\n";
-  input << waterTremolo;
+  input << Tremolo_invalidkey;
   testParser->load(&input);
   //TODO: proove invalidity
@@ -93 +204 @@
   TremoloParser* testParser = new TremoloParser();
   stringstream input;
-  string waterTremolo;
 
   // One simple data line
-  waterTremolo = "#\n#ATOMDATA Id name Type x=3\n1 hydrogen H 3.0 4.5 0.1\n\n";
-  input << waterTremolo;
+  input << Tremolo_Atomdata2;
   testParser->load(&input);
   CPPUNIT_ASSERT(World::getInstance().getAtom(AtomByType(1))->x[0] == 3.0);
@@ -106 +215 @@
   TremoloParser* testParser = new TremoloParser();
   stringstream input;
-  string waterTremolo;
 
   // One simple data line
-  waterTremolo = "#\n#ATOMDATA Id name Type u=3\n1 hydrogen H 3.0 4.5 0.1\n\n";
-  input << waterTremolo;
+  input << Tremolo_velocity;
   testParser->load(&input);
   CPPUNIT_ASSERT(World::getInstance().getAtom(AtomByType(1))->v[0] == 3.0);
@@ -119 +226 @@
   TremoloParser* testParser = new TremoloParser();
   stringstream input;
-  string waterTremolo;
 
   // Neighbor data
-  waterTremolo = "#\n#ATOMDATA Id Type neighbors=2\n1 H 3 0\n2 H 3 0\n3 O 1 2\n";
-  input << waterTremolo;
+  input << Tremolo_neighbours;
   testParser->load(&input);
 
@@ -135 +240 @@
   TremoloParser* testParser = new TremoloParser();
   stringstream input, output;
-  string waterTremolo;
 
   // Neighbor data
-  waterTremolo = "#\n#ATOMDATA Id Type imprData\n8 H 9-10\n9 H 10-8,8-10\n10 O -\n";
-  input << waterTremolo;
+  input << Tremolo_improper;
   testParser->load(&input);
   testParser->save(&output);
@@ -151 +254 @@
   TremoloParser* testParser = new TremoloParser();
   stringstream input, output;
-  string waterTremolo;
 
   // Neighbor data
-  waterTremolo = "#\n#ATOMDATA Id Type torsion\n8 H 9-10\n9 H 10-8,8-10\n10 O -\n";
-  input << waterTremolo;
+  input << Tremolo_torsion;
   testParser->load(&input);
   testParser->save(&output);
@@ -173 +274 @@
   testParser->setFieldsForSave("x=3 u=3 F stress Id neighbors=5 imprData GroupMeasureTypeNo Type extType name resName chainID resSeq occupancy tempFactor segID Charge charge GrpTypeNo torsion");
   testParser->save(&output);
-  CPPUNIT_ASSERT(output.str() == "# ATOMDATA\tx=3\tu=3\tF\tstress\tId\tneighbors=5\timprData\tGroupMeasureTypeNo\tType\textType\tname\tresName\tchainID\tresSeq\toccupancy\ttempFactor\tsegID\tCharge\tcharge\tGrpTypeNo\ttorsion\n0\t0\t0\t0\t0\t0\t0\t0\t0\t0\t0\t0\t0\t0\t-\t0\tH\t-\t-\t-\t0\t0\t0\t0\t0\t0\t0\t0\t-\t\n");
+  CPPUNIT_ASSERT(output.str() == Tremolo_full);
 
   cout << "testing the tremolo parser is done" << endl;
 }
+
+void ParserUnitTest::readwritePcpTest() {
+  stringstream input(waterPcp);
+  PcpParser* testParser = new PcpParser();
+  testParser->load(&input);
+  input.clear();
+
+  CPPUNIT_ASSERT_EQUAL(3, World::getInstance().numAtoms());
+
+  string newWaterPcp = "";
+  stringstream output;
+  testParser->save(&output);
+
+  input << output;
+  PcpParser* testParser2 = new PcpParser();
+  testParser2->load(&input);
+
+  CPPUNIT_ASSERT_EQUAL(6, World::getInstance().numAtoms());
+
+  CPPUNIT_ASSERT(*testParser == *testParser2);
+}
+
+void ParserUnitTest::writeMpqcTest() {
+  // build up water molecule
+  atom *Walker = NULL;
+  Walker = World::getInstance().createAtom();
+  Walker->type = World::getInstance().getPeriode()->FindElement(8);
+  Walker->x = Vector(0,0,0);
+  Walker = World::getInstance().createAtom();
+  Walker->type = World::getInstance().getPeriode()->FindElement(1);
+  Walker->x = Vector(0.758602,0,0.504284);
+  Walker = World::getInstance().createAtom();
+  Walker->type = World::getInstance().getPeriode()->FindElement(1);
+  Walker->x = Vector(0.758602,0,-0.504284);
+  CPPUNIT_ASSERT_EQUAL(3, World::getInstance().numAtoms());
+
+  // create two stringstreams, one stored, one created
+  stringstream input(waterMpqc);
+  MpqcParser* testParser = new MpqcParser();
+  stringstream output;
+  testParser->save(&output);
+
+  // compare both configs
+  string first = input.str();
+  string second = output.str();
+  CPPUNIT_ASSERT(first == second);
+}

src/unittests/ParserUnitTest.hpp

@@ -22 +22 @@
   CPPUNIT_TEST ( readAndWriteTremoloTorsionInformationTest );
   CPPUNIT_TEST ( writeTremoloTest );
+  CPPUNIT_TEST ( readwritePcpTest );
+  CPPUNIT_TEST ( writeMpqcTest );
   CPPUNIT_TEST_SUITE_END();
 
@@ -36 +38 @@
   void readAndWriteTremoloTorsionInformationTest();
   void writeTremoloTest();
+  void readwritePcpTest();
+  void writeMpqcTest();
 };
 

src/unittests/manipulateAtomsTest.cpp

@@ -87 +87 @@
 static void operation(atom* _atom){
   AtomStub *atom = dynamic_cast<AtomStub*>(_atom);
-  assert(atom);
+  CPPUNIT_ASSERT(atom);
   atom->doSomething();
 }
@@ -100 +100 @@
     AtomStub *atom;
     atom = dynamic_cast<AtomStub*>(*iter);
-    assert(atom);
+    CPPUNIT_ASSERT(atom);
     CPPUNIT_ASSERT(atom->manipulated);
   }
@@ -114 +114 @@
     AtomStub *atom;
     atom = dynamic_cast<AtomStub*>(*iter);
-    assert(atom);
+    CPPUNIT_ASSERT(atom);
     if(atom->getId()!=(int)ATOM_COUNT/2)
       CPPUNIT_ASSERT(atom->manipulated);

tests/regression/Domain/6/post/test-x.conf.xyz

@@ -1 +1 @@
 22
         Created by molecuilder on Mon May 31 19:07:16 2010
-H       9.78209 2.64589 2.64589 
-C       27.2836 3.27519 3.53589 
-C       28.5328 4.15859 3.53589 
-C       29.7821 3.27519 3.53589 
-H       27.2836 2.64589 4.42589 
-H       9.78209 2.64589 4.42589 
-H       10.672  3.90454 3.53589 
-H       8.53279 4.78789 2.64589 
-H       27.2836 2.64589 2.64589 
-H       26.3936 3.90454 3.53589 
-H       28.5328 4.78789 4.42589 
-H       28.5328 4.78789 2.64589 
-H       30.672  3.90454 3.53589 
-H       29.7821 2.64589 4.42589 
-H       29.7821 2.64589 2.64589 
-H       8.53279 4.78789 4.42589 
-H       6.39363 3.90454 3.53589 
-H       7.28359 2.64589 2.64589 
-H       7.28359 2.64589 4.42589 
-C       9.78209 3.27519 3.53589 
-C       8.53279 4.15859 3.53589 
-C       7.28359 3.27519 3.53589 
+H       9.78209 2.64589 2.64589
+C       27.2836 3.27519 3.53589
+C       28.5328 4.15859 3.53589
+C       29.7821 3.27519 3.53589
+H       27.2836 2.64589 4.42589
+H       9.78209 2.64589 4.42589
+H       10.672  3.90454 3.53589
+H       8.53279 4.78789 2.64589
+H       27.2836 2.64589 2.64589
+H       26.3936 3.90454 3.53589
+H       28.5328 4.78789 4.42589
+H       28.5328 4.78789 2.64589
+H       30.672  3.90454 3.53589
+H       29.7821 2.64589 4.42589
+H       29.7821 2.64589 2.64589
+H       8.53279 4.78789 4.42589
+H       6.39363 3.90454 3.53589
+H       7.28359 2.64589 2.64589
+H       7.28359 2.64589 4.42589
+C       9.78209 3.27519 3.53589
+C       8.53279 4.15859 3.53589
+C       7.28359 3.27519 3.53589

tests/regression/Domain/6/post/test-y.conf.xyz

@@ -1 +1 @@
 22
         Created by molecuilder on Mon May 31 19:07:16 2010
-H       9.78209 2.64589 2.64589 
-C       7.28359 23.2752 3.53589 
-C       8.53279 24.1586 3.53589 
-C       9.78209 23.2752 3.53589 
-H       7.28359 22.6459 4.42589 
-H       9.78209 2.64589 4.42589 
-H       10.672  3.90454 3.53589 
-H       8.53279 4.78789 2.64589 
-H       7.28359 22.6459 2.64589 
-H       6.39363 23.9045 3.53589 
-H       8.53279 24.7879 4.42589 
-H       8.53279 24.7879 2.64589 
-H       10.672  23.9045 3.53589 
-H       9.78209 22.6459 4.42589 
-H       9.78209 22.6459 2.64589 
-H       8.53279 4.78789 4.42589 
-H       6.39363 3.90454 3.53589 
-H       7.28359 2.64589 2.64589 
-H       7.28359 2.64589 4.42589 
-C       9.78209 3.27519 3.53589 
-C       8.53279 4.15859 3.53589 
-C       7.28359 3.27519 3.53589 
+H       9.78209 2.64589 2.64589
+C       7.28359 23.2752 3.53589
+C       8.53279 24.1586 3.53589
+C       9.78209 23.2752 3.53589
+H       7.28359 22.6459 4.42589
+H       9.78209 2.64589 4.42589
+H       10.672  3.90454 3.53589
+H       8.53279 4.78789 2.64589
+H       7.28359 22.6459 2.64589
+H       6.39363 23.9045 3.53589
+H       8.53279 24.7879 4.42589
+H       8.53279 24.7879 2.64589
+H       10.672  23.9045 3.53589
+H       9.78209 22.6459 4.42589
+H       9.78209 22.6459 2.64589
+H       8.53279 4.78789 4.42589
+H       6.39363 3.90454 3.53589
+H       7.28359 2.64589 2.64589
+H       7.28359 2.64589 4.42589
+C       9.78209 3.27519 3.53589
+C       8.53279 4.15859 3.53589
+C       7.28359 3.27519 3.53589

tests/regression/Domain/6/post/test-z.conf.xyz

@@ -1 +1 @@
 22
         Created by molecuilder on Mon May 31 19:07:16 2010
-H       9.78209 2.64589 2.64589 
-C       7.28359 3.27519 23.5359 
-C       8.53279 4.15859 23.5359 
-C       9.78209 3.27519 23.5359 
-H       7.28359 2.64589 24.4259 
-H       9.78209 2.64589 4.42589 
-H       10.672  3.90454 3.53589 
-H       8.53279 4.78789 2.64589 
-H       7.28359 2.64589 22.6459 
-H       6.39363 3.90454 23.5359 
-H       8.53279 4.78789 24.4259 
-H       8.53279 4.78789 22.6459 
-H       10.672  3.90454 23.5359 
-H       9.78209 2.64589 24.4259 
-H       9.78209 2.64589 22.6459 
-H       8.53279 4.78789 4.42589 
-H       6.39363 3.90454 3.53589 
-H       7.28359 2.64589 2.64589 
-H       7.28359 2.64589 4.42589 
-C       9.78209 3.27519 3.53589 
-C       8.53279 4.15859 3.53589 
-C       7.28359 3.27519 3.53589 
+H       9.78209 2.64589 2.64589
+C       7.28359 3.27519 23.5359
+C       8.53279 4.15859 23.5359
+C       9.78209 3.27519 23.5359
+H       7.28359 2.64589 24.4259
+H       9.78209 2.64589 4.42589
+H       10.672  3.90454 3.53589
+H       8.53279 4.78789 2.64589
+H       7.28359 2.64589 22.6459
+H       6.39363 3.90454 23.5359
+H       8.53279 4.78789 24.4259
+H       8.53279 4.78789 22.6459
+H       10.672  3.90454 23.5359
+H       9.78209 2.64589 24.4259
+H       9.78209 2.64589 22.6459
+H       8.53279 4.78789 4.42589
+H       6.39363 3.90454 3.53589
+H       7.28359 2.64589 2.64589
+H       7.28359 2.64589 4.42589
+C       9.78209 3.27519 3.53589
+C       8.53279 4.15859 3.53589
+C       7.28359 3.27519 3.53589

tests/regression/Domain/6/post/test.conf.xyz

@@ -1 +1 @@
 11
         Created by molecuilder on Mon May 31 19:21:12 2010
-H       9.78209 2.64589 2.64589 
-H       9.78209 2.64589 4.42589 
-H       10.672  3.90454 3.53589 
-H       8.53279 4.78789 2.64589 
-H       8.53279 4.78789 4.42589 
-H       6.39363 3.90454 3.53589 
-H       7.28359 2.64589 2.64589 
-H       7.28359 2.64589 4.42589 
-C       9.78209 3.27519 3.53589 
-C       8.53279 4.15859 3.53589 
-C       7.28359 3.27519 3.53589 
+H       9.78209 2.64589 2.64589
+H       9.78209 2.64589 4.42589
+H       10.672  3.90454 3.53589
+H       8.53279 4.78789 2.64589
+H       8.53279 4.78789 4.42589
+H       6.39363 3.90454 3.53589
+H       7.28359 2.64589 2.64589
+H       7.28359 2.64589 4.42589
+C       9.78209 3.27519 3.53589
+C       8.53279 4.15859 3.53589
+C       7.28359 3.27519 3.53589

tests/regression/Filling/1/post/test.conf

@@ -35 +35 @@
 RelEpsTotalE    1e-07   # relative change in total energy
 RelEpsKineticE  1e-05   # relative change in kinetic energy
-MaxMinStopStep  1350    # check every ..th steps
+MaxMinStopStep  680     # check every ..th steps
 MaxMinGapStopStep       1       # check every ..th steps
 
@@ -42 +42 @@
 InitRelEpsTotalE        1e-05   # relative change in total energy
 InitRelEpsKineticE      0.0001  # relative change in kinetic energy
-InitMaxMinStopStep      1350    # check every ..th steps
+InitMaxMinStopStep      680     # check every ..th steps
 InitMaxMinGapStopStep   1       # check every ..th steps
 
@@ -55 +55 @@
 RiemannTensor   0       # (Use metric)
 PsiType         0       # 0 - doubly occupied, 1 - SpinUp,SpinDown
-MaxPsiDouble    1350    # here: specifying both maximum number of SpinUp- and -Down-states
-PsiMaxNoUp      1350    # here: specifying maximum number of SpinUp-states
-PsiMaxNoDown    1350    # here: specifying maximum number of SpinDown-states
+MaxPsiDouble    680     # here: specifying both maximum number of SpinUp- and -Down-states
+PsiMaxNoUp      680     # here: specifying maximum number of SpinUp-states
+PsiMaxNoDown    680     # here: specifying maximum number of SpinDown-states
 AddPsis         0       # Additional unoccupied Psis for bandgap determination
 

tests/regression/Simple_configuration/2/post/test.conf.xyz

@@ -1 +1 @@
 1
         Created by molecuilder on Tue Oct  6 18:34:23 2009
-H       10      10      10      
+H       10      10      10

tests/regression/Simple_configuration/3/post/test.conf.xyz

@@ -1 +1 @@
 1
         Created by molecuilder on Tue Oct  6 18:34:23 2009
-H       10      10      10      
+H       10      10      10

tests/regression/Simple_configuration/4/post/test.conf

@@ -35 +35 @@
 RelEpsTotalE    1e-07   # relative change in total energy
 RelEpsKineticE  1e-05   # relative change in kinetic energy
-MaxMinStopStep  1       # check every ..th steps
+MaxMinStopStep  2       # check every ..th steps
 MaxMinGapStopStep       1       # check every ..th steps
 
@@ -42 +42 @@
 InitRelEpsTotalE        1e-05   # relative change in total energy
 InitRelEpsKineticE      0.0001  # relative change in kinetic energy
-InitMaxMinStopStep      1       # check every ..th steps
+InitMaxMinStopStep      2       # check every ..th steps
 InitMaxMinGapStopStep   1       # check every ..th steps
 
@@ -54 +54 @@
 Level0Factor    2       # factor by which node number increases from S to 0 level
 RiemannTensor   0       # (Use metric)
-PsiType         1       # 0 - doubly occupied, 1 - SpinUp,SpinDown
-MaxPsiDouble    0       # here: specifying both maximum number of SpinUp- and -Down-states
-PsiMaxNoUp      0       # here: specifying maximum number of SpinUp-states
-PsiMaxNoDown    1       # here: specifying maximum number of SpinDown-states
+PsiType         0       # 0 - doubly occupied, 1 - SpinUp,SpinDown
+MaxPsiDouble    2       # here: specifying both maximum number of SpinUp- and -Down-states
+PsiMaxNoUp      2       # here: specifying maximum number of SpinUp-states
+PsiMaxNoDown    2       # here: specifying maximum number of SpinDown-states
 AddPsis         0       # Additional unoccupied Psis for bandgap determination
 

tests/regression/Simple_configuration/4/post/test.conf.xyz

@@ -1 +1 @@
 1
         Created by molecuilder on Tue Oct  6 18:31:23 2009
-C       10      10      10      
+C       10      10      10

tests/regression/Simple_configuration/5/post/test.conf

@@ -35 +35 @@
 RelEpsTotalE    1e-07   # relative change in total energy
 RelEpsKineticE  1e-05   # relative change in kinetic energy
-MaxMinStopStep  1       # check every ..th steps
+MaxMinStopStep  0       # check every ..th steps
 MaxMinGapStopStep       1       # check every ..th steps
 
@@ -42 +42 @@
 InitRelEpsTotalE        1e-05   # relative change in total energy
 InitRelEpsKineticE      0.0001  # relative change in kinetic energy
-InitMaxMinStopStep      1       # check every ..th steps
+InitMaxMinStopStep      0       # check every ..th steps
 InitMaxMinGapStopStep   1       # check every ..th steps
 
@@ -54 +54 @@
 Level0Factor    2       # factor by which node number increases from S to 0 level
 RiemannTensor   0       # (Use metric)
-PsiType         1       # 0 - doubly occupied, 1 - SpinUp,SpinDown
+PsiType         0       # 0 - doubly occupied, 1 - SpinUp,SpinDown
 MaxPsiDouble    0       # here: specifying both maximum number of SpinUp- and -Down-states
 PsiMaxNoUp      0       # here: specifying maximum number of SpinUp-states
-PsiMaxNoDown    1       # here: specifying maximum number of SpinDown-states
+PsiMaxNoDown    0       # here: specifying maximum number of SpinDown-states
 AddPsis         0       # Additional unoccupied Psis for bandgap determination
 

tests/regression/Simple_configuration/8/post/test.conf.xyz

@@ -1 +1 @@
 144
         Created by molecuilder on Mon May 31 18:50:20 2010
-H       9.78209 2.64589 2.64589 
-H       9.78209 2.64589 4.42589 
-H       10.672  3.90454 3.53589 
-H       8.53279 4.78789 2.64589 
-H       8.53279 4.78789 4.42589 
-H       6.39363 3.90454 3.53589 
-H       7.28359 2.64589 2.64589 
-H       7.28359 2.64589 4.42589 
-H       0.758602        2.85714 3.86571 
-H       0.758602        2.85714 2.85714 
-H       0.758602        2.85714 5.71429 
-H       0.758602        5.71429 3.86571 
-H       0.758602        5.71429 2.85714 
-H       3.61574 0       1.00857 
-H       3.61574 0       0       
-H       3.61574 0       3.86571 
-H       3.61574 0       2.85714 
-H       3.61574 0       6.72285 
-H       3.61574 0       5.71429 
-H       3.61574 2.85714 1.00857 
-H       3.61574 2.85714 0       
-H       3.61574 2.85714 3.86571 
-H       3.61574 2.85714 2.85714 
-H       3.61574 2.85714 6.72285 
-H       3.61574 2.85714 5.71429 
-H       3.61574 2.85714 8.57143 
-H       3.61574 5.71429 1.00857 
-H       3.61574 5.71429 0       
-H       3.61574 5.71429 3.86571 
-H       3.61574 5.71429 2.85714 
-H       3.61574 5.71429 6.72285 
-H       3.61574 5.71429 5.71429 
-H       3.61574 5.71429 8.57143 
-H       3.61574 8.57143 1.00857 
-H       3.61574 8.57143 3.86571 
-H       3.61574 8.57143 2.85714 
-H       3.61574 8.57143 5.71429 
-H       6.47289 0       1.00857 
-H       6.47289 0       0       
-H       6.47289 0       3.86571 
-H       6.47289 0       2.85714 
-H       6.47289 0       6.72285 
-H       6.47289 0       5.71429 
-H       6.47289 0       9.58    
-H       6.47289 0       8.57143 
-H       6.47289 2.85714 0       
-H       6.47289 2.85714 6.72285 
-H       6.47289 2.85714 9.58    
-H       6.47289 2.85714 8.57143 
-H       6.47289 5.71429 1.00857 
-H       6.47289 5.71429 0       
-H       6.47289 5.71429 6.72285 
-H       6.47289 5.71429 5.71429 
-H       6.47289 5.71429 9.58    
-H       6.47289 5.71429 8.57143 
-H       6.47289 8.57143 1.00857 
-H       6.47289 8.57143 0       
-H       6.47289 8.57143 3.86571 
-H       6.47289 8.57143 2.85714 
-H       6.47289 8.57143 6.72285 
-H       6.47289 8.57143 5.71429 
-H       9.33003 0       1.00857 
-H       9.33003 0       0       
-H       9.33003 0       3.86571 
-H       9.33003 0       2.85714 
-H       9.33003 0       6.72285 
-H       9.33003 0       5.71429 
-H       9.33003 0       8.57143 
-H       9.33003 2.85714 0       
-H       9.33003 2.85714 6.72285 
-H       9.33003 2.85714 9.58    
-H       9.33003 2.85714 8.57143 
-H       9.33003 5.71429 0       
-H       9.33003 5.71429 6.72285 
-H       9.33003 5.71429 9.58    
-H       9.33003 5.71429 8.57143 
-H       9.33003 8.57143 1.00857 
-H       9.33003 8.57143 0       
-H       9.33003 8.57143 3.86571 
-H       9.33003 8.57143 2.85714 
-H       9.33003 8.57143 6.72285 
-H       9.33003 8.57143 5.71429 
-H       12.1872 0       1.00857 
-H       12.1872 0       0       
-H       12.1872 0       3.86571 
-H       12.1872 0       2.85714 
-H       12.1872 0       6.72285 
-H       12.1872 0       5.71429 
-H       12.1872 2.85714 1.00857 
-H       12.1872 2.85714 0       
-H       12.1872 2.85714 6.72285 
-H       12.1872 2.85714 5.71429 
-H       12.1872 5.71429 1.00857 
-H       12.1872 5.71429 0       
-H       12.1872 5.71429 3.86571 
-H       12.1872 5.71429 2.85714 
-H       12.1872 5.71429 6.72285 
-H       12.1872 5.71429 5.71429 
-C       9.78209 3.27519 3.53589 
-C       8.53279 4.15859 3.53589 
-C       7.28359 3.27519 3.53589 
-O       2.85714 0       0.504284        
-O       2.85714 0       3.36143 
-O       2.85714 0       6.21857 
-O       2.85714 2.85714 0.504284        
-O       2.85714 2.85714 3.36143 
-O       2.85714 2.85714 6.21857 
-O       2.85714 5.71429 0.504284        
-O       2.85714 5.71429 3.36143 
-O       2.85714 5.71429 6.21857 
-O       2.85714 8.57143 3.36143 
-O       5.71429 0       0.504284        
-O       5.71429 0       3.36143 
-O       5.71429 0       6.21857 
-O       5.71429 0       9.07571 
-O       5.71429 2.85714 0.504284        
-O       5.71429 2.85714 6.21857 
-O       5.71429 2.85714 9.07571 
-O       5.71429 5.71429 0.504284        
-O       5.71429 5.71429 6.21857 
-O       5.71429 5.71429 9.07571 
-O       5.71429 8.57143 0.504284        
-O       5.71429 8.57143 3.36143 
-O       5.71429 8.57143 6.21857 
-O       8.57143 0       0.504284        
-O       8.57143 0       3.36143 
-O       8.57143 0       6.21857 
-O       8.57143 0       9.07571 
-O       8.57143 2.85714 0.504284        
-O       8.57143 2.85714 9.07571 
-O       8.57143 5.71429 0.504284        
-O       8.57143 5.71429 9.07571 
-O       8.57143 8.57143 0.504284        
-O       8.57143 8.57143 3.36143 
-O       8.57143 8.57143 6.21857 
-O       11.4286 0       0.504284        
-O       11.4286 0       3.36143 
-O       11.4286 0       6.21857 
-O       11.4286 2.85714 0.504284        
-O       11.4286 2.85714 6.21857 
-O       11.4286 2.85714 9.07571 
-O       11.4286 5.71429 0.504284        
-O       11.4286 5.71429 6.21857 
-O       11.4286 8.57143 3.36143 
+H       9.78209 2.64589 2.64589
+H       9.78209 2.64589 4.42589
+H       10.672  3.90454 3.53589
+H       8.53279 4.78789 2.64589
+H       8.53279 4.78789 4.42589
+H       6.39363 3.90454 3.53589
+H       7.28359 2.64589 2.64589
+H       7.28359 2.64589 4.42589
+H       0.758602        2.85714 3.86571
+H       0.758602        2.85714 2.85714
+H       0.758602        2.85714 5.71429
+H       0.758602        5.71429 3.86571
+H       0.758602        5.71429 2.85714
+H       3.61574 0       1.00857
+H       3.61574 0       0
+H       3.61574 0       3.86571
+H       3.61574 0       2.85714
+H       3.61574 0       6.72285
+H       3.61574 0       5.71429
+H       3.61574 2.85714 1.00857
+H       3.61574 2.85714 0
+H       3.61574 2.85714 3.86571
+H       3.61574 2.85714 2.85714
+H       3.61574 2.85714 6.72285
+H       3.61574 2.85714 5.71429
+H       3.61574 2.85714 8.57143
+H       3.61574 5.71429 1.00857
+H       3.61574 5.71429 0
+H       3.61574 5.71429 3.86571
+H       3.61574 5.71429 2.85714
+H       3.61574 5.71429 6.72285
+H       3.61574 5.71429 5.71429
+H       3.61574 5.71429 8.57143
+H       3.61574 8.57143 1.00857
+H       3.61574 8.57143 3.86571
+H       3.61574 8.57143 2.85714
+H       3.61574 8.57143 5.71429
+H       6.47289 0       1.00857
+H       6.47289 0       0
+H       6.47289 0       3.86571
+H       6.47289 0       2.85714
+H       6.47289 0       6.72285
+H       6.47289 0       5.71429
+H       6.47289 0       9.58
+H       6.47289 0       8.57143
+H       6.47289 2.85714 0
+H       6.47289 2.85714 6.72285
+H       6.47289 2.85714 9.58
+H       6.47289 2.85714 8.57143
+H       6.47289 5.71429 1.00857
+H       6.47289 5.71429 0
+H       6.47289 5.71429 6.72285
+H       6.47289 5.71429 5.71429
+H       6.47289 5.71429 9.58
+H       6.47289 5.71429 8.57143
+H       6.47289 8.57143 1.00857
+H       6.47289 8.57143 0
+H       6.47289 8.57143 3.86571
+H       6.47289 8.57143 2.85714
+H       6.47289 8.57143 6.72285
+H       6.47289 8.57143 5.71429
+H       9.33003 0       1.00857
+H       9.33003 0       0
+H       9.33003 0       3.86571
+H       9.33003 0       2.85714
+H       9.33003 0       6.72285
+H       9.33003 0       5.71429
+H       9.33003 0       8.57143
+H       9.33003 2.85714 0
+H       9.33003 2.85714 6.72285
+H       9.33003 2.85714 9.58
+H       9.33003 2.85714 8.57143
+H       9.33003 5.71429 0
+H       9.33003 5.71429 6.72285
+H       9.33003 5.71429 9.58
+H       9.33003 5.71429 8.57143
+H       9.33003 8.57143 1.00857
+H       9.33003 8.57143 0
+H       9.33003 8.57143 3.86571
+H       9.33003 8.57143 2.85714
+H       9.33003 8.57143 6.72285
+H       9.33003 8.57143 5.71429
+H       12.1872 0       1.00857
+H       12.1872 0       0
+H       12.1872 0       3.86571
+H       12.1872 0       2.85714
+H       12.1872 0       6.72285
+H       12.1872 0       5.71429
+H       12.1872 2.85714 1.00857
+H       12.1872 2.85714 0
+H       12.1872 2.85714 6.72285
+H       12.1872 2.85714 5.71429
+H       12.1872 5.71429 1.00857
+H       12.1872 5.71429 0
+H       12.1872 5.71429 3.86571
+H       12.1872 5.71429 2.85714
+H       12.1872 5.71429 6.72285
+H       12.1872 5.71429 5.71429
+C       9.78209 3.27519 3.53589
+C       8.53279 4.15859 3.53589
+C       7.28359 3.27519 3.53589
+O       2.85714 0       0.504284
+O       2.85714 0       3.36143
+O       2.85714 0       6.21857
+O       2.85714 2.85714 0.504284
+O       2.85714 2.85714 3.36143
+O       2.85714 2.85714 6.21857
+O       2.85714 5.71429 0.504284
+O       2.85714 5.71429 3.36143
+O       2.85714 5.71429 6.21857
+O       2.85714 8.57143 3.36143
+O       5.71429 0       0.504284
+O       5.71429 0       3.36143
+O       5.71429 0       6.21857
+O       5.71429 0       9.07571
+O       5.71429 2.85714 0.504284
+O       5.71429 2.85714 6.21857
+O       5.71429 2.85714 9.07571
+O       5.71429 5.71429 0.504284
+O       5.71429 5.71429 6.21857
+O       5.71429 5.71429 9.07571
+O       5.71429 8.57143 0.504284
+O       5.71429 8.57143 3.36143
+O       5.71429 8.57143 6.21857
+O       8.57143 0       0.504284
+O       8.57143 0       3.36143
+O       8.57143 0       6.21857
+O       8.57143 0       9.07571
+O       8.57143 2.85714 0.504284
+O       8.57143 2.85714 9.07571
+O       8.57143 5.71429 0.504284
+O       8.57143 5.71429 9.07571
+O       8.57143 8.57143 0.504284
+O       8.57143 8.57143 3.36143
+O       8.57143 8.57143 6.21857
+O       11.4286 0       0.504284
+O       11.4286 0       3.36143
+O       11.4286 0       6.21857
+O       11.4286 2.85714 0.504284
+O       11.4286 2.85714 6.21857
+O       11.4286 2.85714 9.07571
+O       11.4286 5.71429 0.504284
+O       11.4286 5.71429 6.21857
+O       11.4286 8.57143 3.36143

tests/regression/testsuite-fragmentation.at

@@ -12 +12 @@
 AT_SETUP([Fragmentation - Fragmentation])
 AT_CHECK([/bin/cp -f ${abs_top_srcdir}/${AUTOTEST_PATH}/Fragmentation/2/pre/test.conf .], 0)
-AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -v 1 -f 0 --distance 1.55 --order 2], 0, [ignore], [ignore])
+AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -v 1 -f ./BondFragment --molecule-by-id 0 --distance 1.55 --order 2], 0, [ignore], [ignore])
 AT_CHECK([ls -l BondFragment*.conf | wc -l], 0, [5
 ], [ignore])
@@ -21 +21 @@
 AT_SETUP([Fragmentation - BROKEN: Fragmentation is at MaxOrder])
 AT_CHECK([/bin/cp -f ${abs_top_srcdir}/${AUTOTEST_PATH}/Fragmentation/3/pre/test.conf .], 0)
-AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -v 1 -f 0 --distance 1.55 --order 2], 0, [ignore], [ignore])
-AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -v 1 -f 0 --distance 1.55 --order 2], [ignore], [ignore], [ignore])
+AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -v 1 -f ./BondFragment --molecule-by-id 0 --distance 1.55 --order 2], 0, [ignore], [ignore])
+AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -v 1 -f ./BondFragment --molecule-by-id 0 --distance 1.55 --order 2], [ignore], [ignore], [ignore])
 AT_CLEANUP

tests/regression/testsuite-simple_configuration.at

@@ -52 +52 @@
 AT_SETUP([Simple configuration - invalid commands on empty configs])
 AT_KEYWORDS([configuration])
-AT_CHECK([../../molecuilder empty.conf -e ${abs_top_srcdir}/src/ -t -s -b -E -c -b -a -U -T -u], 255, [ignore], [stderr])
-AT_CHECK([fgrep -c "Not enough" stderr], 0, [1
-], [ignore])
+AT_CHECK([../../molecuilder empty.conf -e ${abs_top_srcdir}/src/ -t -s -b -E -c -b -a -U -T -u], 134, [ignore], [stderr])
 AT_CLEANUP
 
@@ -61 +59 @@
 AT_KEYWORDS([configuration])
 AT_CHECK([/bin/cp -f ${abs_top_srcdir}/${AUTOTEST_PATH}/Simple_configuration/7/pre/test.conf .], 0)
-AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -t], 255, [ignore], [stderr])
-AT_CHECK([fgrep -c "CRITICAL: Not enough" stderr], 0, [ignore], [ignore])
-AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -s -b -E -c -b -a -U -T -u], 255, [ignore], [stderr])
-AT_CHECK([fgrep -c "CRITICAL: Not enough" stderr], 0, [ignore], [ignore])
-AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -b -E -c -b -a -U -T -u], 255, [ignore], [stderr])
-AT_CHECK([fgrep -c "CRITICAL: Not enough" stderr], 0, [ignore], [ignore])
-AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -E -c -b -a -U -T -u], 255, [ignore], [stderr])
-AT_CHECK([fgrep -c "CRITICAL: Not enough" stderr], 0, [ignore], [ignore])
-AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -c -b -a -U -T -u], 255, [ignore], [stderr])
-AT_CHECK([fgrep -c "CRITICAL: Not enough" stderr], 0, [ignore], [ignore])
-AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -b -a -U -T -u], 255, [ignore], [stderr])
-AT_CHECK([fgrep -c "CRITICAL: Not enough" stderr], 0, [ignore], [ignore])
-AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -a -U -T -u], 255, [ignore], [stderr])
-AT_CHECK([fgrep -c "CRITICAL: Not enough" stderr], 0, [ignore], [ignore])
-AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -U -T -u], 255, [ignore], [stderr])
-AT_CHECK([fgrep -c "CRITICAL: Not enough" stderr], 0, [ignore], [ignore])
-AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -T -u], 255, [ignore], [stderr])
-AT_CHECK([fgrep -c "CRITICAL: Not enough" stderr], 0, [ignore], [ignore])
-AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -u], 255, [ignore], [stderr])
-AT_CHECK([fgrep -c "CRITICAL: Not enough" stderr], 0, [ignore], [ignore])
+AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -t], 134, [ignore], [stderr])
+AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -s -b -E -c -b -a -U -T -u], 134, [ignore], [stderr])
+AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -b -E -c -b -a -U -T -u], 134, [ignore], [stderr])
+AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -E -c -b -a -U -T -u], 134, [ignore], [stderr])
+AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -c -b -a -U -T -u], 134, [ignore], [stderr])
+AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -b -a -U -T -u], 134, [ignore], [stderr])
+AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -a -U -T -u], 134, [ignore], [stderr])
+AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -U -T -u], 134, [ignore], [stderr])
+AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -T -u], 134, [ignore], [stderr])
+AT_CHECK([../../molecuilder test.conf -e ${abs_top_srcdir}/src/ -u], 134, [ignore], [stderr])
 AT_CLEANUP
 

tests/regression/testsuite-standard_options.at

@@ -18 +18 @@
 AT_SETUP([Standard Options - no element database])
 AT_KEYWORDS([options])
-AT_CHECK([../../molecuilder -e], 255, [ignore], [stderr])
-AT_CHECK([fgrep "Not enough or invalid arguments" stderr], 0, [ignore], [ignore])
+AT_CHECK([../../molecuilder -e], 134, [ignore], [stderr])
 AT_CLEANUP