Changeset 36ec71 for src/datacreator.cpp

Timestamp:

Jul 24, 2009, 10:38:32 AM (16 years ago)

Author:

Frederik Heber <heber@…>

Branches:

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

Children:

d30402

Parents:

042f82 (diff), 51c910 (diff)
Note: this is a merge changeset, the changes displayed below correspond to the merge itself.
Use the (diff) links above to see all the changes relative to each parent.

git-author:

Frederik Heber <heber@…> (07/23/09 14:23:32)

git-committer:

Frederik Heber <heber@…> (07/24/09 10:38:32)

Message:

Merge branch 'master' into ConcaveHull

Conflicts:

molecuilder/src/analyzer.cpp
molecuilder/src/bond.cpp
molecuilder/src/boundary.cpp
molecuilder/src/boundary.hpp
molecuilder/src/builder.cpp
molecuilder/src/config.cpp
molecuilder/src/datacreator.cpp
molecuilder/src/datacreator.hpp
molecuilder/src/defs.hpp
molecuilder/src/ellipsoid.cpp
molecuilder/src/joiner.cpp
molecuilder/src/molecules.cpp
molecuilder/src/molecules.hpp
molecuilder/src/parser.cpp
molecuilder/src/parser.hpp

merges:

• analyzer.cpp: all from master, due to Hessian; only Chi(PAS)Fragments from HEAD
• bond.cpp: bond::GetOtherAtom gives cerr output if neither atom matches from master
• boundary.cpp: all from HEAD, is tesselation stuff
• boundary.hpp: all from HEAD, is tesselation stuff
• builder.cpp: Center...(), RemoveAtoms(), SaveConfig() from master, is MultipleMolecules, rest from HEAD (Translation and ...)
• config.cpp: FileBuffer stuff from HEAD, only some taken over from HEAD
• datacreator.cpp: all from master, is Hessian
• datacreator.hpp: all from master, is Hessian
• defs.hpp: all from master, is Thermostat
• ellipsoid.cpp: all from HEAD, created during tesselation
• joiner.cpp: all from master, is Hessian
• molecules.cpp: molecule::molecule(), ConstrainedPotential() stuff, Thermostats(), RemoveAtom(),CreateAdjacencyList(), FragmentMolecule() from master, CenterInBox(), CenterEdge(), CenterOrigin(), CenterPeriodic() (instead of CenterAtVector()), TranslatePeriodically(), DeterminePeriodicCenter(), Output...(), ReturnFullForSymmnetricMatrix() from HEAD
• molecules.hpp: see above, ConfigfileBuffer struct from HEAD, class config joined
• parser.cpp: all from master, is Hessian
• parser.hpp: all from master, is Hessian

compilation fixes:

• molecules.cpp: PrincipalAxisSystem() - CenterGravity() -> CenterOrigin(), DetermineCenter() -> DeterminePeriodicCenter
• datacreator.cpp: all char * to const char *
• builder.cpp: ParseCommandLineOptions() - case 'c': had argptr++ and argptr+=3; case 'O' - x was not zero'd

File:

• src/datacreator.cpp (modified) (37 diffs)

src/datacreator.cpp

@@ -63 +63 @@
   cout << msg << endl;
   output << "# " << msg << ", created on " << datum;
-  output << "#Order\tFrag.No.\t" << Fragments.Header << endl;
+  output << "#Order\tFrag.No.\t" << Fragments.Header[Fragments.MatrixCounter] << endl;
   for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
     for(int i=KeySet.FragmentsPerOrder[BondOrder];i--;) {
       for(int j=Fragments.RowCounter[ KeySet.OrderSet[BondOrder][i] ];j--;)
-        for(int k=Fragments.ColumnCounter;k--;)
+        for(int k=Fragments.ColumnCounter[ KeySet.OrderSet[BondOrder][i] ];k--;)
           Fragments.Matrix[Fragments.MatrixCounter][j][k] += Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
     }
     output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
-    for (int l=0;l<Fragments.ColumnCounter;l++)
+    for (int l=0;l<Fragments.ColumnCounter[Fragments.MatrixCounter];l++)
       output << scientific << "\t" << Fragments.Matrix[Fragments.MatrixCounter][ Fragments.RowCounter[Fragments.MatrixCounter]-1 ][l];
     output << endl;
@@ -96 +96 @@
   cout << msg << endl;
   output << "# " << msg << ", created on " << datum;
-  output << "#Order\tFrag.No.\t" << Fragments.Header << endl;
+  output << "#Order\tFrag.No.\t" << Fragments.Header[Fragments.MatrixCounter] << endl;
   Fragments.SetLastMatrix(Energy.Matrix[Energy.MatrixCounter],0);
   for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
     for(int i=KeySet.FragmentsPerOrder[BondOrder];i--;) {
       for(int j=Fragments.RowCounter[ KeySet.OrderSet[BondOrder][i] ];j--;)
-        for(int k=Fragments.ColumnCounter;k--;)
+        for(int k=Fragments.ColumnCounter[ KeySet.OrderSet[BondOrder][i] ];k--;)
           Fragments.Matrix[Fragments.MatrixCounter][j][k] -= Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
     }
     output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
-    for (int l=0;l<Fragments.ColumnCounter;l++)
+    for (int l=0;l<Fragments.ColumnCounter[Energy.MatrixCounter];l++)
       if (fabs(Energy.Matrix[Energy.MatrixCounter][ Energy.RowCounter[Energy.MatrixCounter]-1 ][l]) < MYEPSILON)
         output << scientific << "\t" << Fragments.Matrix[Fragments.MatrixCounter][ Fragments.RowCounter[Fragments.MatrixCounter]-1 ][l];
@@ -133 +133 @@
   cout << msg << endl;
   output << "# " << msg << ", created on " << datum;
-  output << "# Order\tFrag.No.\t" << Fragments.Header << endl;
+  output << "# Order\tFrag.No.\t" << Fragments.Header[Fragments.MatrixCounter] << endl;
   Fragments.SetLastMatrix(0.,0);
   for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
@@ -139 +139 @@
     output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
     CreateForce(Fragments, Fragments.MatrixCounter);
-    for (int l=0;l<Fragments.ColumnCounter;l++)
+    for (int l=0;l<Fragments.ColumnCounter[Fragments.MatrixCounter];l++)
        output << scientific << "\t" << Fragments.Matrix[Fragments.MatrixCounter][ Fragments.RowCounter[Fragments.MatrixCounter] ][l];
     output << endl;
@@ -165 +165 @@
   cout << msg << endl;
   output << "# " << msg << ", created on " << datum;
-  output << "# Order\tFrag.No.\t" << Fragments.Header << endl;
+  output << "# Order\tFrag.No.\t" << Fragments.Header[Fragments.MatrixCounter] << endl;
   Fragments.SetLastMatrix(Force.Matrix[Force.MatrixCounter],0);
   for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
@@ -171 +171 @@
     output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
     CreateForce(Fragments, Fragments.MatrixCounter);
-    for (int l=0;l<Fragments.ColumnCounter;l++)
+    for (int l=0;l<Fragments.ColumnCounter[Fragments.MatrixCounter];l++)
        output << scientific << "\t" << Fragments.Matrix[Fragments.MatrixCounter][ Fragments.RowCounter[Fragments.MatrixCounter] ][l];
     output << endl;
@@ -198 +198 @@
   cout << msg << endl;
   output << "# " << msg << ", created on " << datum;
-  output << "# AtomNo\t" << Fragments.Header << endl;
+  output << "# AtomNo\t" << Fragments.Header[Fragments.MatrixCounter] << endl;
   Fragments.SetLastMatrix(Force.Matrix[Force.MatrixCounter], 0);
   for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
@@ -207 +207 @@
     for(int j=0;j<Fragments.RowCounter[ Fragments.MatrixCounter ];j++) {
       output << Fragments.Indices[Fragments.MatrixCounter][j] << "\t";
-      for (int l=0;l<Fragments.ColumnCounter;l++) {
+      for (int l=0;l<Fragments.ColumnCounter[ Fragments.MatrixCounter ];l++) {
         if (((l+1) % 3) == 0) {
           norm = 0.;
@@ -213 +213 @@
             norm += Force.Matrix[Force.MatrixCounter][ j ][l+m]*Force.Matrix[Force.MatrixCounter][ j ][l+m];
           norm = sqrt(norm);
-        }
+        }                                                                                                    
 //        if (norm < MYEPSILON)
           output << scientific << Fragments.Matrix[Fragments.MatrixCounter][ j ][l] << "\t";
@@ -244 +244 @@
   cout << msg << endl;
   output << "# " << msg << ", created on " << datum;
-  output << "# AtomNo\t" << Fragments.Header << endl;
+  output << "# AtomNo\t" << Fragments.Header[Fragments.MatrixCounter] << endl;
   for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
     //cout << "Current order is " << BondOrder << "." << endl;
@@ -252 +252 @@
     for(int j=0;j<Fragments.RowCounter[ Fragments.MatrixCounter ];j++) {
       output << Fragments.Indices[Fragments.MatrixCounter][j] << "\t";
-      for (int l=0;l<Fragments.ColumnCounter;l++)
+      for (int l=0;l<Fragments.ColumnCounter[ Fragments.MatrixCounter ];l++)
         output << scientific << Fragments.Matrix[Fragments.MatrixCounter][ j ][l] << "\t";
       output << endl;
@@ -262 +262 @@
 };
 
+
+/** Plot hessian error vs. vs atom vs. order.
+ * \param &Hessian HessianMatrix containing reference values (in MatrixCounter matrix)
+ * \param &Fragments HessianMatrix class containing matrix values
+ * \param KeySet KeySetContainer class holding bond KeySetContainer::Order
+ * \param *prefix prefix in filename (without ending)
+ * \param *msg message to be place in first line as a comment
+ * \param *datum current date and time
+ * \return true if file was written successfully
+ */
+bool CreateDataDeltaHessianOrderPerAtom(class HessianMatrix &Hessian, class HessianMatrix &Fragments, class KeySetsContainer &KeySet, const char *dir, const char *prefix, const char *msg, const char *datum) 
+{
+  stringstream filename;
+  ofstream output;
+
+  filename << prefix << ".dat";
+  if (!OpenOutputFile(output, dir, filename.str().c_str())) return false; 
+  cout << msg << endl;
+  output << "# " << msg << ", created on " << datum;
+  output << "# AtomNo\t" << Fragments.Header[Fragments.MatrixCounter] << endl;
+  Fragments.SetLastMatrix(Hessian.Matrix[Hessian.MatrixCounter], 0);
+  for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
+    //cout << "Current order is " << BondOrder << "." << endl;
+    Fragments.SumSubHessians(Fragments, KeySet, BondOrder, -1.);
+    // errors per atom
+    output << endl << "#Order\t" << BondOrder+1 << endl;
+    for(int j=0;j<Fragments.RowCounter[ Fragments.MatrixCounter ];j++) { 
+      output << Fragments.Indices[Fragments.MatrixCounter][j] << "\t";
+      for (int l=0;l<Fragments.ColumnCounter[ Fragments.MatrixCounter ];l++) {
+        output << scientific << Fragments.Matrix[Fragments.MatrixCounter][ j ][l] << "\t";
+      }
+      output << endl;
+    }
+    output << endl;
+  }
+  output.close();
+  return true;
+};
+
+/** Plot hessian error vs. vs atom vs. order in the frobenius norm.
+ * \param &Hessian HessianMatrix containing reference values (in MatrixCounter matrix)
+ * \param &Fragments HessianMatrix class containing matrix values
+ * \param KeySet KeySetContainer class holding bond KeySetContainer::Order
+ * \param *prefix prefix in filename (without ending)
+ * \param *msg message to be place in first line as a comment
+ * \param *datum current date and time
+ * \return true if file was written successfully
+ */
+bool CreateDataDeltaFrobeniusOrderPerAtom(class HessianMatrix &Hessian, class HessianMatrix &Fragments, class KeySetsContainer &KeySet, const char *dir, const char *prefix, const char *msg, const char *datum) 
+{
+  stringstream filename;
+  ofstream output;
+  double norm = 0;
+  double tmp;
+
+  filename << prefix << ".dat";
+  if (!OpenOutputFile(output, dir, filename.str().c_str())) return false; 
+  cout << msg << endl;
+  output << "# " << msg << ", created on " << datum;
+  output << "# AtomNo\t";
+  Fragments.SetLastMatrix(Hessian.Matrix[Hessian.MatrixCounter], 0);
+  for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
+    output << "Order" << BondOrder+1 << "\t";
+  }
+  output << endl;
+  output << Fragments.RowCounter[ Fragments.MatrixCounter ] << "\t";
+  for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
+    //cout << "Current order is " << BondOrder << "." << endl;
+    Fragments.SumSubHessians(Fragments, KeySet, BondOrder, -1.);
+    // frobenius norm of errors per atom
+    norm = 0.;
+    for(int j=0;j<Fragments.RowCounter[ Fragments.MatrixCounter ];j++) { 
+      for (int l=0;l<Fragments.ColumnCounter[ Fragments.MatrixCounter ];l++) {
+        tmp = Fragments.Matrix[Fragments.MatrixCounter][ j ][l];
+        norm += tmp*tmp;
+      }
+    }
+    output << scientific << sqrt(norm)/(Fragments.RowCounter[ Fragments.MatrixCounter ]*Fragments.ColumnCounter[ Fragments.MatrixCounter] ) << "\t";
+  }
+  output << endl;
+  output.close();
+  return true;
+};
+
+/** Plot hessian error vs. vs atom vs. order.
+ * \param &Fragments HessianMatrix class containing matrix values
+ * \param KeySet KeySetContainer class holding bond KeySetContainer::Order
+ * \param *prefix prefix in filename (without ending)
+ * \param *msg message to be place in first line as a comment
+ * \param *datum current date and time
+ * \return true if file was written successfully
+ */
+bool CreateDataHessianOrderPerAtom(class HessianMatrix &Fragments, class KeySetsContainer &KeySet, const char *dir, const char *prefix, const char *msg, const char *datum) 
+{
+  stringstream filename;
+  ofstream output;
+
+  filename << prefix << ".dat";
+  if (!OpenOutputFile(output, dir, filename.str().c_str())) return false; 
+  cout << msg << endl;
+  output << "# " << msg << ", created on " << datum;
+  output << "# AtomNo\t" << Fragments.Header[ Fragments.MatrixCounter ] << endl;
+  Fragments.SetLastMatrix(0., 0);
+  for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
+    //cout << "Current order is " << BondOrder << "." << endl;
+    Fragments.SumSubHessians(Fragments, KeySet, BondOrder, 1.);
+    // errors per atom
+    output << endl << "#Order\t" << BondOrder+1 << endl;
+    for(int j=0;j<Fragments.RowCounter[ Fragments.MatrixCounter ];j++) { 
+      output << Fragments.Indices[Fragments.MatrixCounter][j] << "\t";
+      for (int l=0;l<Fragments.ColumnCounter[ Fragments.MatrixCounter ];l++)
+        output << scientific << Fragments.Matrix[Fragments.MatrixCounter][ j ][l] << "\t";
+      output << endl;
+    }
+    output << endl;
+  }
+  output.close();
+  return true;
+};
+
 /** Plot matrix vs. fragment.
  */
@@ -273 +393 @@
   cout << msg << endl;
   output << "# " << msg << ", created on " << datum << endl;
-  output << "#Order\tFrag.No.\t" << Fragment.Header << endl;
+  output << "#Order\tFrag.No.\t" << Fragment.Header[ Fragment.MatrixCounter ] << endl;
   for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
     for(int i=0;i<KeySet.FragmentsPerOrder[BondOrder];i++) {
       output << BondOrder+1 << "\t" << KeySet.OrderSet[BondOrder][i]+1;
       CreateFragment(Fragment, KeySet.OrderSet[BondOrder][i]);
-      for (int l=0;l<Fragment.ColumnCounter;l++)
+      for (int l=0;l<Fragment.ColumnCounter[ KeySet.OrderSet[BondOrder][i] ];l++)
         output << scientific << "\t" << Fragment.Matrix[ KeySet.OrderSet[BondOrder][i] ][ Fragment.RowCounter[ KeySet.OrderSet[BondOrder][i] ] ][l];
       output << endl;
@@ -297 +417 @@
     for(int i=KeySet.FragmentsPerOrder[BondOrder];i--;) {
       if (fabs(Fragments.Matrix[ Fragments.MatrixCounter ][j][1]) < fabs(Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][1])) {
-        for (int k=Fragments.ColumnCounter;k--;)
+        for (int k=Fragments.ColumnCounter[ Fragments.MatrixCounter ];k--;)
           Fragments.Matrix[ Fragments.MatrixCounter ][j][k] = Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
       }
@@ -314 +434 @@
     int i=0;
     do {  // first get a minimum value unequal to 0
-      for (int k=Fragments.ColumnCounter;k--;)
+      for (int k=Fragments.ColumnCounter[ Fragments.MatrixCounter ];k--;)
         Fragments.Matrix[ Fragments.MatrixCounter ][j][k] = Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
       i++;
@@ -320 +440 @@
     for(;i<KeySet.FragmentsPerOrder[BondOrder];i++) { // then find lowest
       if (fabs(Fragments.Matrix[ Fragments.MatrixCounter ][j][1]) > fabs(Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][1])) {
-        for (int k=Fragments.ColumnCounter;k--;)
+        for (int k=Fragments.ColumnCounter[ Fragments.MatrixCounter ];k--;)
           Fragments.Matrix[ Fragments.MatrixCounter ][j][k] = Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
       }
@@ -338 +458 @@
   cout << msg << endl;
   output << "# " << msg << ", created on " << datum;
-  output << "#Order\tFrag.No.\t" << Fragment.Header << endl;
+  output << "#Order\tFrag.No.\t" << Fragment.Header[ Fragment.MatrixCounter ] << endl;
   // max
   for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
@@ -344 +464 @@
     CreateFragmentOrder(Fragment, KeySet, BondOrder);
     output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
-    for (int l=0;l<Fragment.ColumnCounter;l++)
+    for (int l=0;l<Fragment.ColumnCounter[ Fragment.MatrixCounter ];l++)
       output << scientific << "\t" << Fragment.Matrix[ Fragment.MatrixCounter ][ Fragment.RowCounter[ Fragment.MatrixCounter ]-1 ][l];
     output << endl;
@@ -358 +478 @@
 void CreateEnergy(class MatrixContainer &Energy, int MatrixNumber)
 {
-  for(int k=0;k<Energy.ColumnCounter;k++)
+  for(int k=0;k<Energy.ColumnCounter[MatrixNumber];k++)
     Energy.Matrix[MatrixNumber][ Energy.RowCounter[MatrixNumber] ] [k] =  Energy.Matrix[MatrixNumber][ Energy.RowCounter[MatrixNumber]-1 ] [k];
 };
@@ -369 +489 @@
 void CreateMinimumForce(class MatrixContainer &Force, int MatrixNumber)
 {
-  for (int l=Force.ColumnCounter;l--;)
+  for (int l=Force.ColumnCounter[MatrixNumber];l--;)
     Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l] = 0.;
-  for (int l=5;l<Force.ColumnCounter;l+=3) {
+  for (int l=5;l<Force.ColumnCounter[MatrixNumber];l+=3) {
     double stored = 0;
     int k=0;
@@ -404 +524 @@
 {
   int divisor = 0;
-  for (int l=Force.ColumnCounter;l--;)
+  for (int l=Force.ColumnCounter[MatrixNumber];l--;)
     Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l] = 0.;
-  for (int l=5;l<Force.ColumnCounter;l+=3) {
+  for (int l=5;l<Force.ColumnCounter[MatrixNumber];l+=3) {
     double tmp = 0;
     for (int k=Force.RowCounter[MatrixNumber];k--;) {
@@ -428 +548 @@
 void CreateMaximumForce(class MatrixContainer &Force, int MatrixNumber)
 {
-  for (int l=5;l<Force.ColumnCounter;l+=3) {
+  for (int l=5;l<Force.ColumnCounter[MatrixNumber];l+=3) {
     double stored = 0;
     for (int k=Force.RowCounter[MatrixNumber];k--;) {
@@ -460 +580 @@
 void CreateVectorSumForce(class MatrixContainer &Force, int MatrixNumber)
 {
-  for (int l=Force.ColumnCounter;l--;)
+  for (int l=Force.ColumnCounter[MatrixNumber];l--;)
     Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l] = 0.;
-  for (int l=0;l<Force.ColumnCounter;l++) {
+  for (int l=0;l<Force.ColumnCounter[MatrixNumber];l++) {
     for (int k=Force.RowCounter[MatrixNumber];k--;)
       Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l] += Force.Matrix[MatrixNumber][k][l];
@@ -538 +658 @@
 void AbsEnergyPlotLine(ofstream &output, class MatrixContainer &Energy, const char *prefix, const char *xargument, const char *uses)
 {
-  stringstream line(Energy.Header);
+  stringstream line(Energy.Header[ Energy.MatrixCounter ]);
   string token;
 
   getline(line, token, '\t');
-  for (int i=2; i<= Energy.ColumnCounter;i++) {
+  for (int i=2; i<= Energy.ColumnCounter[Energy.MatrixCounter];i++) {
     getline(line, token, '\t');
     while (token[0] == ' ') // remove leading white spaces
       token.erase(0,1);
     output << "'" << prefix << ".dat' title '" << token << "' using " << xargument << ":(abs($" << i+2 << ")) " << uses;
-    if (i != (Energy.ColumnCounter))
+    if (i != (Energy.ColumnCounter[Energy.MatrixCounter]))
       output << ", \\";
     output << endl;
@@ -562 +682 @@
 void EnergyPlotLine(ofstream &output, class MatrixContainer &Energy, const char *prefix, const char *xargument, const char *uses)
 {
-  stringstream line(Energy.Header);
+  stringstream line(Energy.Header[Energy.MatrixCounter]);
   string token;
 
   getline(line, token, '\t');
-  for (int i=1; i<= Energy.ColumnCounter;i++) {
+  for (int i=1; i<= Energy.ColumnCounter[Energy.MatrixCounter];i++) {
     getline(line, token, '\t');
     while (token[0] == ' ') // remove leading white spaces
       token.erase(0,1);
     output << "'" << prefix << ".dat' title '" << token << "' using " << xargument << ":" << i+2 << " " << uses;
-    if (i != (Energy.ColumnCounter))
+    if (i != (Energy.ColumnCounter[Energy.MatrixCounter]))
       output << ", \\";
     output << endl;
@@ -586 +706 @@
 void ForceMagnitudePlotLine(ofstream &output, class MatrixContainer &Force, const char *prefix, const char *xargument, const char *uses)
 {
-  stringstream line(Force.Header);
+  stringstream line(Force.Header[Force.MatrixCounter]);
   string token;
 
@@ -594 +714 @@
   getline(line, token, '\t');
   getline(line, token, '\t');
-  for (int i=7; i< Force.ColumnCounter;i+=NDIM) {
+  for (int i=7; i< Force.ColumnCounter[Force.MatrixCounter];i+=NDIM) {
     getline(line, token, '\t');
     while (token[0] == ' ') // remove leading white spaces
@@ -600 +720 @@
     token.erase(token.length(), 1);  // kill residual index char (the '0')
     output << "'" << prefix << ".dat' title '" << token << "' using " << xargument << ":(sqrt($" << i+1 << "*$" << i+1 << "+$" << i+2 << "*$" << i+2 << "+$" << i+3 << "*$" << i+3 << ")) " << uses;
-    if (i != (Force.ColumnCounter-1))
+    if (i != (Force.ColumnCounter[Force.MatrixCounter]-1))
       output << ", \\";
     output << endl;
@@ -617 +737 @@
 void AbsFirstForceValuePlotLine(ofstream &output, class MatrixContainer &Force, const char *prefix, const char *xargument, const char *uses)
 {
-  stringstream line(Force.Header);
+  stringstream line(Force.Header[Force.MatrixCounter]);
   string token;
 
@@ -625 +745 @@
   getline(line, token, '\t');
   getline(line, token, '\t');
-  for (int i=7; i< Force.ColumnCounter;i+=NDIM) {
+  for (int i=7; i< Force.ColumnCounter[Force.MatrixCounter];i+=NDIM) {
     getline(line, token, '\t');
     while (token[0] == ' ') // remove leading white spaces
@@ -631 +751 @@
     token.erase(token.length(), 1);  // kill residual index char (the '0')
     output << "'" << prefix << ".dat' title '" << token << "' using " << xargument << ":(abs($" << i+1 << ")) " << uses;
-    if (i != (Force.ColumnCounter-1))
+    if (i != (Force.ColumnCounter[Force.MatrixCounter]-1))
       output << ", \\";
     output << endl;
@@ -648 +768 @@
 void BoxesForcePlotLine(ofstream &output, class MatrixContainer &Force, const char *prefix, const char *xargument, const char *uses)
 {
-  stringstream line(Force.Header);
-  const char *fillcolor[5] = {"black", "red", "blue", "green", "cyan"};
+  stringstream line(Force.Header[Force.MatrixCounter]);
+  char *fillcolor[5] = {"black", "red", "blue", "green", "cyan"};
   string token;
 
@@ -657 +777 @@
   getline(line, token, '\t');
   getline(line, token, '\t');
-  for (int i=7; i< Force.ColumnCounter;i+=NDIM) {
+  for (int i=7; i< Force.ColumnCounter[Force.MatrixCounter];i+=NDIM) {
     getline(line, token, '\t');
     while (token[0] == ' ') // remove leading white spaces
@@ -663 +783 @@
     token.erase(token.length(), 1);  // kill residual index char (the '0')
     output << "'" << prefix << ".dat' title '" << token << "' using ($" << xargument << "+" << fixed << setprecision(1) << (double)((i-7)/3)*0.2 << "):(sqrt($" << i+1 << "*$" << i+1 << "+$" << i+2 << "*$" << i+2 << "+$" << i+3 << "*$" << i+3 << ")) " << uses << " " << fillcolor[(i-7)/3];
-    if (i != (Force.ColumnCounter-1))
+    if (i != (Force.ColumnCounter[Force.MatrixCounter]-1))
       output << ", \\";
     output << endl;
@@ -680 +800 @@
 void BoxesFirstForceValuePlotLine(ofstream &output, class MatrixContainer &Force, const char *prefix, const char *xargument, const char *uses)
 {
-  stringstream line(Force.Header);
-  const char *fillcolor[5] = {"black", "red", "blue", "green", "cyan"};
+  stringstream line(Force.Header[Force.MatrixCounter]);
+  char *fillcolor[5] = {"black", "red", "blue", "green", "cyan"};
   string token;
 
@@ -689 +809 @@
   getline(line, token, '\t');
   getline(line, token, '\t');
-  for (int i=7; i< Force.ColumnCounter;i+=NDIM) {
+  for (int i=7; i< Force.ColumnCounter[Force.MatrixCounter];i+=NDIM) {
     getline(line, token, '\t');
     while (token[0] == ' ') // remove leading white spaces
@@ -695 +815 @@
     token.erase(token.length(), 1);  // kill residual index char (the '0')
     output << "'" << prefix << ".dat' title '" << token << "' using ($" << xargument << "+" << fixed << setprecision(1) << (double)((i-7)/3)*0.2 << "):" << i+1 << " " << uses << " " << fillcolor[(i-7)/3];
-    if (i != (Force.ColumnCounter-1))
+    if (i != (Force.ColumnCounter[Force.MatrixCounter]-1))
       output << ", \\";
     output << endl;