Changeset 4e855e for src

Timestamp:

Feb 4, 2011, 6:56:38 PM (14 years ago)

Author:

Frederik Heber <heber@…>

Branches:

Action_Thermostats, Add_AtomRandomPerturbation, Add_FitFragmentPartialChargesAction, Add_RotateAroundBondAction, Add_SelectAtomByNameAction, Added_ParseSaveFragmentResults, AddingActions_SaveParseParticleParameters, Adding_Graph_to_ChangeBondActions, Adding_MD_integration_tests, Adding_ParticleName_to_Atom, Adding_StructOpt_integration_tests, AtomFragments, Automaking_mpqc_open, AutomationFragmentation_failures, Candidate_v1.5.4, Candidate_v1.6.0, Candidate_v1.6.1, ChangeBugEmailaddress, ChangingTestPorts, ChemicalSpaceEvaluator, CombiningParticlePotentialParsing, Combining_Subpackages, Debian_Package_split, Debian_package_split_molecuildergui_only, Disabling_MemDebug, Docu_Python_wait, EmpiricalPotential_contain_HomologyGraph, EmpiricalPotential_contain_HomologyGraph_documentation, Enable_parallel_make_install, Enhance_userguide, Enhanced_StructuralOptimization, Enhanced_StructuralOptimization_continued, Example_ManyWaysToTranslateAtom, Exclude_Hydrogens_annealWithBondGraph, FitPartialCharges_GlobalError, Fix_BoundInBox_CenterInBox_MoleculeActions, Fix_ChargeSampling_PBC, Fix_ChronosMutex, Fix_FitPartialCharges, Fix_FitPotential_needs_atomicnumbers, Fix_ForceAnnealing, Fix_IndependentFragmentGrids, Fix_ParseParticles, Fix_ParseParticles_split_forward_backward_Actions, Fix_PopActions, Fix_QtFragmentList_sorted_selection, Fix_Restrictedkeyset_FragmentMolecule, Fix_StatusMsg, Fix_StepWorldTime_single_argument, Fix_Verbose_Codepatterns, Fix_fitting_potentials, Fixes, ForceAnnealing_goodresults, ForceAnnealing_oldresults, ForceAnnealing_tocheck, ForceAnnealing_with_BondGraph, ForceAnnealing_with_BondGraph_continued, ForceAnnealing_with_BondGraph_continued_betteresults, ForceAnnealing_with_BondGraph_contraction-expansion, FragmentAction_writes_AtomFragments, FragmentMolecule_checks_bonddegrees, GeometryObjects, Gui_Fixes, Gui_displays_atomic_force_velocity, ImplicitCharges, IndependentFragmentGrids, IndependentFragmentGrids_IndividualZeroInstances, IndependentFragmentGrids_IntegrationTest, IndependentFragmentGrids_Sole_NN_Calculation, JobMarket_RobustOnKillsSegFaults, JobMarket_StableWorkerPool, JobMarket_unresolvable_hostname_fix, MoreRobust_FragmentAutomation, ODR_violation_mpqc_open, PartialCharges_OrthogonalSummation, PdbParser_setsAtomName, PythonUI_with_named_parameters, QtGui_reactivate_TimeChanged_changes, Recreated_GuiChecks, Rewrite_FitPartialCharges, RotateToPrincipalAxisSystem_UndoRedo, SaturateAtoms_findBestMatching, SaturateAtoms_singleDegree, StoppableMakroAction, Subpackage_CodePatterns, Subpackage_JobMarket, Subpackage_LinearAlgebra, Subpackage_levmar, Subpackage_mpqc_open, Subpackage_vmg, Switchable_LogView, ThirdParty_MPQC_rebuilt_buildsystem, TrajectoryDependenant_MaxOrder, TremoloParser_IncreasedPrecision, TremoloParser_MultipleTimesteps, TremoloParser_setsAtomName, Ubuntu_1604_changes, stable

Children:

febef3

Parents:

fd3788

git-author:

Frederik Heber <heber@…> (02/04/11 18:54:12)

git-committer:

Frederik Heber <heber@…> (02/04/11 18:56:38)

Message:

BondGraph::LoadBondLengthTable() now accepts istream instead of const char *.

• This has been quite a lot of work because the matrix parsing is involved.
• Hence, the fix is not really clean, this should get fixed as soon as we deal with the BondGraph implementation.
• This fixes ticket #130.

Location:

src

Files:

• Actions/CommandAction/BondLengthTableAction.cpp (modified) (2 diffs)
• Actions/CommandAction/BondLengthTableAction.def (modified) (2 diffs)
• Actions/FragmentationAction/DepthFirstSearchAction.cpp (modified) (1 diff)
• bondgraph.cpp (modified) (3 diffs)
• bondgraph.hpp (modified) (1 diff)
• builder_init.cpp (modified) (2 diffs)
• molecule_dynamics.cpp (modified) (2 diffs)
• molecule_graph.cpp (modified) (2 diffs)
• parser.cpp (modified) (14 diffs)
• parser.hpp (modified) (1 diff)
• unittests/AnalysisBondsUnitTest.cpp (modified) (3 diffs)
• unittests/BondGraphUnitTest.cpp (modified) (6 diffs)
• unittests/BondGraphUnitTest.hpp (modified) (1 diff)

src/Actions/CommandAction/BondLengthTableAction.cpp

@@ -45 +45 @@
   DoLog(0) && (Log() << Verbose(0) << "Using " << params.BondGraphFileName << " as bond length table." << endl);
   config *configuration = World::getInstance().getConfig();
-  if (configuration->BG == NULL) {
-    configuration->BG = new BondGraph(configuration->GetIsAngstroem());
-    if ((!params.BondGraphFileName.empty())
-        && boost::filesystem::exists(params.BondGraphFileName)
-        && (configuration->BG->LoadBondLengthTable(params.BondGraphFileName.string()))) {
-      DoLog(0) && (Log() << Verbose(0) << "Bond length table loaded successfully." << endl);
-      return Action::success;
+  BondGraph *OldBG;
+  if (configuration->BG != NULL) {
+    OldBG = configuration->BG;
+    DoLog(0) && (Log() << Verbose(0) << "There is a bond length table already present." << endl);
+  }
+  configuration->BG = new BondGraph(configuration->GetIsAngstroem());
+  if ((!params.BondGraphFileName.empty())
+      && boost::filesystem::exists(params.BondGraphFileName)) {
+    std::ifstream input(params.BondGraphFileName.string().c_str());
+    if ((input.good()) && (configuration->BG->LoadBondLengthTable(input))) {
+      DoLog(0) && (Log() << Verbose(0) << "Bond length table parsed successfully." << endl);
+      input.close();
+      return Action::state_ptr(new CommandBondLengthTableState(*OldBG, params));
     } else {
-      DoeLog(1) && (eLog()<< Verbose(1) << "Bond length table loading failed." << endl);
+      DoeLog(1) && (eLog()<< Verbose(1) << "Bond length table parsing failed." << endl);
+      input.close();
+      configuration->BG = new BondGraph(*OldBG);
       return Action::failure;
     }
   } else {
-    DoLog(0) && (Log() << Verbose(0) << "Bond length table already present." << endl);
+    DoeLog(1) && (eLog()<< Verbose(1) << "Bond length table loading failed." << endl);
+    configuration->BG = new BondGraph(*OldBG);
     return Action::failure;
   }
@@ -63 +72 @@
 
 Action::state_ptr CommandBondLengthTableAction::performUndo(Action::state_ptr _state) {
-//  ParserLoadXyzState *state = assert_cast<ParserLoadXyzState*>(_state.get());
+  CommandBondLengthTableState *state = assert_cast<CommandBondLengthTableState*>(_state.get());
 
-  return Action::failure;
-//  string newName = state->mol->getName();
-//  state->mol->setName(state->lastName);
-//
-//  return Action::state_ptr(new ParserLoadXyzState(state->mol,newName));
+  config *configuration = World::getInstance().getConfig();
+  BondGraph *OldBG = new BondGraph(*configuration->BG);
+  configuration->BG = new BondGraph(state->OldBG);
+
+  return Action::state_ptr(new CommandBondLengthTableState(*OldBG,params));
 }
 
 Action::state_ptr CommandBondLengthTableAction::performRedo(Action::state_ptr _state){
-  return Action::failure;
+  CommandBondLengthTableState *state = assert_cast<CommandBondLengthTableState*>(_state.get());
+
+  config *configuration = World::getInstance().getConfig();
+  BondGraph *OldBG = new BondGraph(*configuration->BG);
+  configuration->BG = new BondGraph(state->OldBG);
+
+  return Action::state_ptr(new CommandBondLengthTableState(*OldBG,params));
 }
 

src/Actions/CommandAction/BondLengthTableAction.def

@@ -8 +8 @@
 // all includes and forward declarations necessary for non-integral types below
 
+#include "bondgraph.hpp"
 
 // i.e. there is an integer with variable name Z that can be found in
@@ -18 +19 @@
 #define paramreferences (BondGraphFileName)
 
-#undef statetypes
-#undef statereferences
+#define statetypes (BondGraph)
+#define statereferences (OldBG)
 
 // some defines for all the names, you may use ACTION, STATE and PARAMS

src/Actions/FragmentationAction/DepthFirstSearchAction.cpp

@@ -52 +52 @@
   std::deque<bond *> *BackEdgeStack = NULL;
   std::deque<bond *> *LocalBackEdgeStack = NULL;
-  mol->CreateAdjacencyList(params.distance, World::getInstance().getConfig()->GetIsAngstroem(), &BondGraph::CovalentMinMaxDistance, NULL);
+  BondGraph *BG = World::getInstance().getConfig()->BG;
+  if (BG != NULL)
+    mol->CreateAdjacencyList(params.distance, World::getInstance().getConfig()->GetIsAngstroem(), &BondGraph::BondLengthMatrixMinMaxDistance, BG);
+  else
+    mol->CreateAdjacencyList(params.distance, World::getInstance().getConfig()->GetIsAngstroem(), &BondGraph::CovalentMinMaxDistance, BG);
   Subgraphs = mol->DepthFirstSearchAnalysis(BackEdgeStack);
   if (Subgraphs != NULL) {

src/bondgraph.cpp

@@ -57 +57 @@
  * Allocates \a MatrixContainer for BondGraph::BondLengthMatrix, using MatrixContainer::ParseMatrix(),
  * but only if parsing is successful. Otherwise variable is left as NULL.
- * \param *out output stream for debugging
- * \param filename file with bond lengths to parse
+ * \param &input input stream to parse table from
  * \return true - success in parsing file, false - failed to parse the file
  */
-bool BondGraph::LoadBondLengthTable(const string &filename)
+bool BondGraph::LoadBondLengthTable(std::istream &input)
 {
   Info FunctionInfo(__func__);
@@ -75 +74 @@
 
   // parse in matrix
-  if ((status = TempContainer->ParseMatrix(filename.c_str(), 0, 1, 0))) {
+  if ((input.good()) && (status = TempContainer->ParseMatrix(input, 0, 1, 0))) {
     DoLog(1) && (Log() << Verbose(1) << "Parsing bond length matrix successful." << endl);
   } else {
     DoeLog(1) && (eLog()<< Verbose(1) << "Parsing bond length matrix failed." << endl);
+    status = false;
   }
 
@@ -129 +129 @@
 double BondGraph::GetBondLength(int firstZ, int secondZ)
 {
-  if (BondLengthMatrix == NULL)
+  std::cout << "Request for length between " << firstZ << " and " << secondZ << ": ";
+  if (BondLengthMatrix == NULL) {
+    std::cout << "-1." << std::endl;
     return( -1. );
-  else
+  } else {
+    std::cout << BondLengthMatrix->Matrix[0][firstZ][secondZ] << std::endl;
     return (BondLengthMatrix->Matrix[0][firstZ][secondZ]);
+  }
 };
 

src/bondgraph.hpp

@@ -36 +36 @@
   BondGraph(bool IsA);
   ~BondGraph();
-  bool LoadBondLengthTable(const string &filename);
+  bool LoadBondLengthTable(std::istream &input);
   bool ConstructBondGraph(molecule * const mol);
   double GetBondLength(int firstelement, int secondelement);

src/builder_init.cpp

@@ -19 +19 @@
 
 #include "CodePatterns/MemDebug.hpp"
+
+#include <iostream>
 
 #include "bondgraph.hpp"
@@ -107 +109 @@
     World::getInstance().getConfig()->BG = new BondGraph(World::getInstance().getConfig()->GetIsAngstroem());
     if (boost::filesystem::exists(BondGraphFileName)) {
-      if (World::getInstance().getConfig()->BG->LoadBondLengthTable(BondGraphFileName)) {
+      std::ifstream input(BondGraphFileName.c_str());
+      if ((input.good()) && (World::getInstance().getConfig()->BG->LoadBondLengthTable(input))) {
         DoLog(0) && (Log() << Verbose(0) << "Bond length table loaded successfully." << endl);
       } else {
         DoeLog(1) && (eLog()<< Verbose(1) << "Bond length table loading failed." << endl);
       }
+      input.close();
     }
   }

src/molecule_dynamics.cpp

@@ -574 +574 @@
 {
   Info FunctionInfo(__func__);
-  ifstream input(file);
   string token;
   stringstream item;
@@ -582 +581 @@
 
   const int AtomCount = getAtomCount();
-  // check file
-  if (input == NULL) {
+  // parse file into ForceMatrix
+  std::ifstream input(file);
+  if ((input.good()) && (!Force.ParseMatrix(input, 0,0,0))) {
+    DoeLog(0) && (eLog()<< Verbose(0) << "Could not parse Force Matrix file " << file << "." << endl);
+    performCriticalExit();
     return false;
-  } else {
-    // parse file into ForceMatrix
-    if (!Force.ParseMatrix(file, 0,0,0)) {
-      DoeLog(0) && (eLog()<< Verbose(0) << "Could not parse Force Matrix file " << file << "." << endl);
-      performCriticalExit();
-      return false;
-    }
-    if (Force.RowCounter[0] != AtomCount) {
-      DoeLog(0) && (eLog()<< Verbose(0) << "Mismatch between number of atoms in file " << Force.RowCounter[0] << " and in molecule " << getAtomCount() << "." << endl);
-      performCriticalExit();
-      return false;
-    }
-    // correct Forces
-    Velocity.Zero();
-    for(int i=0;i<AtomCount;i++)
-      for(int d=0;d<NDIM;d++) {
-        Velocity[d] += Force.Matrix[0][i][d+offset];
-      }
-    for(int i=0;i<AtomCount;i++)
-      for(int d=0;d<NDIM;d++) {
-        Force.Matrix[0][i][d+offset] -= Velocity[d]/static_cast<double>(AtomCount);
-      }
-    // solve a constrained potential if we are meant to
-    if (configuration.DoConstrainedMD) {
-      // calculate forces and potential
-      atom **PermutationMap = NULL;
-      ConstrainedPotentialEnergy = MinimiseConstrainedPotential(PermutationMap,configuration.DoConstrainedMD, 0, configuration.GetIsAngstroem());
-      EvaluateConstrainedForces(configuration.DoConstrainedMD, 0, PermutationMap, &Force);
-      delete[](PermutationMap);
-    }
-
-    // and perform Verlet integration for each atom with position, velocity and force vector
-    // check size of vectors
-    for_each(atoms.begin(),
-             atoms.end(),
-             boost::bind(&atom::VelocityVerletUpdate,_1,MDSteps+1, &configuration, &Force, (const size_t) 0));
-  }
+  }
+  input.close();
+  if (Force.RowCounter[0] != AtomCount) {
+    DoeLog(0) && (eLog()<< Verbose(0) << "Mismatch between number of atoms in file " << Force.RowCounter[0] << " and in molecule " << getAtomCount() << "." << endl);
+    performCriticalExit();
+    return false;
+  }
+  // correct Forces
+  Velocity.Zero();
+  for(int i=0;i<AtomCount;i++)
+    for(int d=0;d<NDIM;d++) {
+      Velocity[d] += Force.Matrix[0][i][d+offset];
+    }
+  for(int i=0;i<AtomCount;i++)
+    for(int d=0;d<NDIM;d++) {
+      Force.Matrix[0][i][d+offset] -= Velocity[d]/static_cast<double>(AtomCount);
+    }
+  // solve a constrained potential if we are meant to
+  if (configuration.DoConstrainedMD) {
+    // calculate forces and potential
+    atom **PermutationMap = NULL;
+    ConstrainedPotentialEnergy = MinimiseConstrainedPotential(PermutationMap,configuration.DoConstrainedMD, 0, configuration.GetIsAngstroem());
+    EvaluateConstrainedForces(configuration.DoConstrainedMD, 0, PermutationMap, &Force);
+    delete[](PermutationMap);
+  }
+
+  // and perform Verlet integration for each atom with position, velocity and force vector
+  // check size of vectors
+  for_each(atoms.begin(),
+           atoms.end(),
+           boost::bind(&atom::VelocityVerletUpdate,_1,MDSteps+1, &configuration, &Force, (const size_t) 0));
+
   // correct velocities (rather momenta) so that center of mass remains motionless
   Velocity = atoms.totalMomentumAtStep(MDSteps+1);

src/molecule_graph.cpp

@@ -174 +174 @@
         for (LC->n[2] = 0; LC->n[2] < LC->N[2]; LC->n[2]++) {
           const LinkedCell::LinkedNodes *List = LC->GetCurrentCell();
-//          Log() << Verbose(2) << "Current cell is " << LC->n[0] << ", " << LC->n[1] << ", " << LC->n[2] << " with No. " << LC->index << " containing " << List->size() << " points." << endl;
+          Log() << Verbose(2) << "Current cell is " << LC->n[0] << ", " << LC->n[1] << ", " << LC->n[2] << " with No. " << LC->index << " containing " << List->size() << " points." << endl;
           if (List != NULL) {
             for (LinkedCell::LinkedNodes::const_iterator Runner = List->begin(); Runner != List->end(); Runner++) {
               Walker = dynamic_cast<atom*>(*Runner);
               ASSERT(Walker,"Tesselpoint that was not an atom retrieved from LinkedNode");
-              //Log() << Verbose(0) << "Current Atom is " << *Walker << "." << endl;
+              Log() << Verbose(0) << "Current Atom is " << *Walker << "." << endl;
               // 3c. check for possible bond between each atom in this and every one in the 27 cells
               for (n[0] = -1; n[0] <= 1; n[0]++)
@@ -185 +185 @@
                   for (n[2] = -1; n[2] <= 1; n[2]++) {
                     const LinkedCell::LinkedNodes *OtherList = LC->GetRelativeToCurrentCell(n);
-//                    Log() << Verbose(2) << "Current relative cell is " << LC->n[0] << ", " << LC->n[1] << ", " << LC->n[2] << " with No. " << LC->index << " containing " << List->size() << " points." << endl;
                     if (OtherList != NULL) {
+                      Log() << Verbose(2) << "Current relative cell is " << LC->n[0] << ", " << LC->n[1] << ", " << LC->n[2] << " with No. " << LC->index << " containing " << List->size() << " points." << endl;
                       for (LinkedCell::LinkedNodes::const_iterator OtherRunner = OtherList->begin(); OtherRunner != OtherList->end(); OtherRunner++) {
                         if ((*OtherRunner)->nr > Walker->nr) {
                           OtherWalker = dynamic_cast<atom*>(*OtherRunner);
                           ASSERT(OtherWalker,"TesselPoint that was not an atom retrieved from LinkedNode");
-                          //Log() << Verbose(1) << "Checking distance " << OtherWalker->x.PeriodicDistanceSquared(&(Walker->x), cell_size) << " against typical bond length of " << bonddistance*bonddistance << "." << endl;
                           (BG->*minmaxdistance)(Walker, OtherWalker, MinDistance, MaxDistance, IsAngstroem);
                           const double distance = domain.periodicDistanceSquared(OtherWalker->getPosition(),Walker->getPosition());
+                          Log() << Verbose(1) << "Checking distance " << distance << " against typical bond length of " << bonddistance*bonddistance << "." << endl;
                           const bool status = (distance <= MaxDistance * MaxDistance) && (distance >= MinDistance * MinDistance);
-//                          Log() << Verbose(1) << "MinDistance is " << MinDistance << " and MaxDistance is " << MaxDistance << "." << endl;
+                          Log() << Verbose(1) << "MinDistance is " << MinDistance << " and MaxDistance is " << MaxDistance << "." << endl;
                           if (OtherWalker->father->nr > Walker->father->nr) {
                             if (status) { // create bond if distance is smaller
-//                              Log() << Verbose(1) << "Adding Bond between " << *Walker << " and " << *OtherWalker << " in distance " << sqrt(distance) << "." << endl;
+                              Log() << Verbose(1) << "Adding Bond between " << *Walker << " and " << *OtherWalker << " in distance " << sqrt(distance) << "." << endl;
                               AddBond(Walker->father, OtherWalker->father, 1); // also increases molecule::BondCount
                             } else {
-//                              Log() << Verbose(1) << "Not Adding: distance too great." << endl;
+                              Log() << Verbose(1) << "Not Adding: distance too great." << endl;
                             }
                           } else {
-//                            Log() << Verbose(1) << "Not Adding: Wrong order of labels." << endl;
+                            Log() << Verbose(1) << "Not Adding: Wrong order of labels." << endl;
                           }
                         }

src/parser.cpp

@@ -22 +22 @@
 
 #include <cstring>
+#include <iomanip>
 
 #include "Helpers/helpers.hpp"
@@ -158 +159 @@
  *    -# allocate matrix
  *    -# loop over found column and row counts and parse in each entry
- * \param *name directory with files
+ * \param &input input stream
  * \param skiplines number of inital lines to skip
  * \param skiplines number of inital columns to skip
@@ -164 +165 @@
  * \return parsing successful
  */
-bool MatrixContainer::ParseMatrix(const char *name, int skiplines, int skipcolumns, int MatrixNr)
-{
-  ifstream input;
+bool MatrixContainer::ParseMatrix(std::istream &input, int skiplines, int skipcolumns, int MatrixNr)
+{
   stringstream line;
   string token;
   char filename[1023];
 
-  input.open(name, ios::in);
   //Log() << Verbose(1) << "Opening " << name << " ... "  << input << endl;
-  if (input == NULL) {
-    DoeLog(1) && (eLog()<< Verbose(1) << endl << "MatrixContainer::ParseMatrix: Unable to open " << name << ", is the directory correct?" << endl);
+  if (input.bad()) {
+    DoeLog(1) && (eLog()<< Verbose(1) << endl << "MatrixContainer::ParseMatrix: Unable to parse istream." << endl);
     //performCriticalExit();
     return false;
@@ -190 +189 @@
   for(int k=skipcolumns;k--;)
     line >> Header[MatrixNr];
-  //Log() << Verbose(0) << line.str() << endl;
+  Log() << Verbose(0) << line.str() << endl;
   ColumnCounter[MatrixNr]=0;
   while ( getline(line,token, '\t') ) {
@@ -196 +195 @@
       ColumnCounter[MatrixNr]++;
   }
-  //Log() << Verbose(0) << line.str() << endl;
-  //Log() << Verbose(1) << "ColumnCounter[" << MatrixNr << "]: " << ColumnCounter[MatrixNr] << "." << endl;
+  Log() << Verbose(0) << line.str() << endl;
+  Log() << Verbose(1) << "ColumnCounter[" << MatrixNr << "]: " << ColumnCounter[MatrixNr] << "." << endl;
   if (ColumnCounter[MatrixNr] == 0) {
-    DoeLog(0) && (eLog()<< Verbose(0) << "ColumnCounter[" << MatrixNr << "]: " << ColumnCounter[MatrixNr] << " from file " << name << ", this is probably an error!" << endl);
+    DoeLog(0) && (eLog()<< Verbose(0) << "ColumnCounter[" << MatrixNr << "]: " << ColumnCounter[MatrixNr] << " from ostream." << endl);
     performCriticalExit();
   }
@@ -207 +206 @@
   while (!input.eof()) {
     input.getline(filename, 1023);
-    //Log() << Verbose(0) << "Comparing: " << strncmp(filename,"MeanForce",9) << endl;
+    Log() << Verbose(0) << "Comparing: " << strncmp(filename,"MeanForce",9) << endl;
     RowCounter[MatrixNr]++; // then line was not last MeanForce
     if (strncmp(filename,"MeanForce", 9) == 0) {
@@ -213 +212 @@
     }
   }
-  //Log() << Verbose(1) << "RowCounter[" << MatrixNr << "]: " << RowCounter[MatrixNr] << " from file " << name << "." << endl;
+  Log() << Verbose(1) << "RowCounter[" << MatrixNr << "]: " << RowCounter[MatrixNr] << " from input stream." << endl;
   if (RowCounter[MatrixNr] == 0) {
-    DoeLog(0) && (eLog()<< Verbose(0) << "RowCounter[" << MatrixNr << "]: " << RowCounter[MatrixNr] << " from file " << name << ", this is probably an error!" << endl);
+    DoeLog(0) && (eLog()<< Verbose(0) << "RowCounter[" << MatrixNr << "]: " << RowCounter[MatrixNr] << " from input stream." << endl);
     performCriticalExit();
   }
@@ -233 +232 @@
       input.getline(Header[MatrixNr], 1023);    // skip header
     line.str(Header[MatrixNr]);
+    Log() << Verbose(0) << "Header: " << line.str() << endl;
     for(int k=skipcolumns;k--;)  // skip columns in header too
       line >> filename;
@@ -245 +245 @@
       input.getline(filename, 1023);
       stringstream lines(filename);
-      //Log() << Verbose(2) << "Matrix at level " << j << ":";// << filename << endl;
+      Log() << Verbose(2) << "Matrix at level " << j << ":";// << filename << endl;
       for(int k=skipcolumns;k--;)
         lines >> filename;
       for(int k=0;(k<ColumnCounter[MatrixNr]) && (!lines.eof());k++) {
         lines >> Matrix[MatrixNr][j][k];
-        //Log() << Verbose(1) << " " << setprecision(2) << Matrix[MatrixNr][j][k] << endl;
+        Log() << Verbose(1) << " " << setprecision(2) << Matrix[MatrixNr][j][k] << endl;
       }
       if (Matrix[MatrixNr][ RowCounter[MatrixNr] ] != NULL)
@@ -261 +261 @@
     DoeLog(1) && (eLog()<< Verbose(1) << "Matrix nr. " << MatrixNr << " has column and row count of (" << ColumnCounter[MatrixNr] << "," << RowCounter[MatrixNr] << "), could not allocate nor parse!" << endl); 
   }
-  input.close();
   return true;
 };
@@ -327 +326 @@
     file.str(" ");
     file << name << FRAGMENTPREFIX << FragmentNumber << prefix << suffix;
-    if (!ParseMatrix(file.str().c_str(), skiplines, skipcolumns, i))
+    std::ifstream input(file.str().c_str());
+    if (!ParseMatrix(input, skiplines, skipcolumns, i)) {
+      input.close();
       return false;
+    }
+    input.close();
     delete[](FragmentNumber);
   }
@@ -628 +631 @@
     strncat(filename, prefix, 1023-strlen(filename));
     strncat(filename, suffix.c_str(), 1023-strlen(filename));
-    ParseMatrix(filename, skiplines, skipcolumns, MatrixCounter);
+    std::ifstream input(filename);
+    ParseMatrix(input, skiplines, skipcolumns, MatrixCounter);
+    input.close();
   }
   return status;
@@ -764 +769 @@
     strncat(filename, prefix, 1023-strlen(filename));
     strncat(filename, suffix.c_str(), 1023-strlen(filename));
-    ParseMatrix(filename, skiplines, skipcolumns, MatrixCounter);
+    std::ifstream input(filename);
+    ParseMatrix(input, skiplines, skipcolumns, MatrixCounter);
+    input.close();
   }
  
@@ -992 +999 @@
     strncat(filename, prefix, 1023-strlen(filename));
     strncat(filename, suffix.c_str(), 1023-strlen(filename));
-    ParseMatrix(filename, skiplines, skipcolumns, MatrixCounter);
+    std::ifstream input(filename);
+    ParseMatrix(input, skiplines, skipcolumns, MatrixCounter);
+    input.close();
   }
 

src/parser.hpp

@@ -86 +86 @@
   
   bool InitialiseIndices(class MatrixContainer *Matrix = NULL);
-  bool ParseMatrix(const char *name, int skiplines, int skipcolumns, int MatrixNr);
+  bool ParseMatrix(std::istream &input, int skiplines, int skipcolumns, int MatrixNr);
   virtual bool ParseFragmentMatrix(const char *name, const char *prefix, string suffix, int skiplines, int skipcolumns);
   bool AllocateMatrix(char **GivenHeader, int MCounter, int *RCounter, int *CCounter);

src/unittests/AnalysisBondsUnitTest.cpp

@@ -92 +92 @@
   CPPUNIT_ASSERT_EQUAL( TestMolecule->getAtomCount(), 5 );
 
-  // create a small file with table
-  filename = new string("test.dat");
-  CPPUNIT_ASSERT(filename != NULL && "could not create string");
-  ofstream test(filename->c_str());
+  // create stream with table
+  std::stringstream test;
   test << ".\tH\tHe\tLi\tBe\tB\tC\n";
   test << "H\t1.\t1.\t1.\t1.\t1.\t1.2\n";
@@ -103 +101 @@
   test << "B\t1.\t1.\t1.\t1.\t1.\t1.\n";
   test << "C\t1.2\t1.\t1.\t1.\t1.\t1.5\n";
-  test.close();
   BG = new BondGraph(true);
   CPPUNIT_ASSERT(BG != NULL && "could not create BondGraph");
 
-  CPPUNIT_ASSERT_EQUAL( true , BG->LoadBondLengthTable(*filename) );
+  CPPUNIT_ASSERT_EQUAL( true , BG->LoadBondLengthTable(test) );
   CPPUNIT_ASSERT_EQUAL( 1., BG->GetBondLength(0,0) );
   CPPUNIT_ASSERT_EQUAL( 1.2, BG->GetBondLength(0,5) );
@@ -119 +116 @@
 {
   // remove the file
-  remove(filename->c_str());
-  delete(filename);
   delete(BG);
 

src/unittests/BondGraphUnitTest.cpp

@@ -92 +92 @@
   CPPUNIT_ASSERT_EQUAL( TestMolecule->getAtomCount(), 4 );
 
-  // create a small file with table
-  dummyname = new string("dummy.dat");
-  CPPUNIT_ASSERT(dummyname != NULL && "could not create string");
-  filename = new string("test.dat");
-  CPPUNIT_ASSERT(filename != NULL && "could not create string");
-  ofstream test(filename->c_str());
+  // create stream with table
   test << ".\tH\tHe\tLi\tBe\tB\tC\n";
   test << "H\t1.\t1.\t1.\t1.\t1.\t1.2\n";
@@ -105 +100 @@
   test << "B\t1.\t1.\t1.\t1.\t1.\t1.\n";
   test << "C\t1.2\t1.\t1.\t1.\t1.\t1.5\n";
-  test.close();
+  // created bad stream (i.e. non-present file)
+  dummy.setstate(ios::eofbit);
+  CPPUNIT_ASSERT(dummy.eof());
   BG = new BondGraph(true);
   CPPUNIT_ASSERT(BG != NULL && "could not create BondGraph");
@@ -114 +111 @@
 {
   // remove the file
-  remove(filename->c_str());
-  delete(filename);
-  delete(dummyname);
   delete(BG);
 
@@ -139 +133 @@
 void BondGraphTest::LoadTableTest()
 {
-  CPPUNIT_ASSERT_EQUAL( true , BG->LoadBondLengthTable(*filename) );
+  CPPUNIT_ASSERT_EQUAL( true , BG->LoadBondLengthTable(test) );
   CPPUNIT_ASSERT_EQUAL( 1., BG->GetBondLength(0,0) );
   CPPUNIT_ASSERT_EQUAL( 1.2, BG->GetBondLength(0,5) );
@@ -152 +146 @@
   molecule::iterator Runner = TestMolecule->begin();
   Runner++;
-  CPPUNIT_ASSERT_EQUAL( true , BG->LoadBondLengthTable(*filename) );
+  CPPUNIT_ASSERT_EQUAL( true , BG->LoadBondLengthTable(test) );
   CPPUNIT_ASSERT_EQUAL( true , BG->ConstructBondGraph(TestMolecule) );
   CPPUNIT_ASSERT_EQUAL( true , (*Walker)->IsBondedTo((*Runner)) );
@@ -165 +159 @@
   //atom *Runner = TestMolecule->end->previous;
   //CPPUNIT_ASSERT( TestMolecule->end != Walker );
-  CPPUNIT_ASSERT_EQUAL( false , BG->LoadBondLengthTable(*dummyname) );
+  CPPUNIT_ASSERT_EQUAL( false , BG->LoadBondLengthTable(dummy) );
   CPPUNIT_ASSERT_EQUAL( true , BG->ConstructBondGraph(TestMolecule) );
 

src/unittests/BondGraphUnitTest.hpp

@@ -49 +49 @@
 
       BondGraph *BG;
-      string *filename;
-      string *dummyname;
+      std::stringstream dummy;
+      std::stringstream test;
 };