Changeset 437922 for src/parser.cpp


Ignore:
Timestamp:
Jul 23, 2009, 12:14:13 PM (16 years ago)
Author:
Frederik Heber <heber@…>
Branches:
Action_Thermostats, Add_AtomRandomPerturbation, Add_FitFragmentPartialChargesAction, Add_RotateAroundBondAction, Add_SelectAtomByNameAction, Added_ParseSaveFragmentResults, AddingActions_SaveParseParticleParameters, Adding_Graph_to_ChangeBondActions, Adding_MD_integration_tests, Adding_ParticleName_to_Atom, Adding_StructOpt_integration_tests, AtomFragments, Automaking_mpqc_open, AutomationFragmentation_failures, Candidate_v1.5.4, Candidate_v1.6.0, Candidate_v1.6.1, Candidate_v1.7.0, 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:
d067d45
Parents:
178f92
Message:

Fix indentation from tab to two spaces.

The trouble was caused at the merge e08f45e4539ffcc30e039dec5606cf06b45ab6be. Seemingly, I thought eclipse had pulled some shit which i didn't

File:
1 edited

Legend:

Unmodified
Added
Removed

  • src/parser.cpp

    r178f92 r437922  
    2424bool FilePresent(const char *filename, bool test)
    2525{
    26         ifstream input;
    27 
    28         input.open(filename, ios::in);
    29         if (input == NULL) {
    30                 if (!test)
    31                         cout << endl << "Unable to open " << filename << ", is the directory correct?" << endl;
    32                 return false;
    33         }
    34         input.close();
    35         return true;
     26  ifstream input;
     27
     28  input.open(filename, ios::in);
     29  if (input == NULL) {
     30    if (!test)
     31      cout << endl << "Unable to open " << filename << ", is the directory correct?" << endl;
     32    return false;
     33  }
     34  input.close();
     35  return true;
    3636};
    3737
    … …  
    4343bool TestParams(int argc, char **argv)
    4444{
    45         ifstream input;
    46         stringstream line;
    47 
    48         line << argv[1] << FRAGMENTPREFIX << KEYSETFILE;
    49         return FilePresent(line.str().c_str(), false);
     45  ifstream input;
     46  stringstream line;
     47
     48  line << argv[1] << FRAGMENTPREFIX << KEYSETFILE;
     49  return FilePresent(line.str().c_str(), false);
    5050};
    5151
    … …  
    5555 */
    5656MatrixContainer::MatrixContainer() {
    57         Indices = NULL;
    58         Header = (char *) Malloc(sizeof(char)*1023, "MatrixContainer::MatrixContainer: *Header");
    59         Matrix = (double ***) Malloc(sizeof(double **)*(1), "MatrixContainer::MatrixContainer: ***Matrix"); // one more each for the total molecule
    60         RowCounter = (int *) Malloc(sizeof(int)*(1), "MatrixContainer::MatrixContainer: *RowCounter");
    61         Matrix[0] = NULL;
    62         RowCounter[0] = -1;
    63         MatrixCounter = 0;
    64         ColumnCounter = 0;
     57  Indices = NULL;
     58  Header = (char *) Malloc(sizeof(char)*1023, "MatrixContainer::MatrixContainer: *Header");
     59  Matrix = (double ***) Malloc(sizeof(double **)*(1), "MatrixContainer::MatrixContainer: ***Matrix"); // one more each for the total molecule
     60  RowCounter = (int *) Malloc(sizeof(int)*(1), "MatrixContainer::MatrixContainer: *RowCounter");
     61  Matrix[0] = NULL;
     62  RowCounter[0] = -1;
     63  MatrixCounter = 0;
     64  ColumnCounter = 0;
    6565};
    6666
    … …  
    6868 */
    6969MatrixContainer::~MatrixContainer() {
    70         if (Matrix != NULL) {
    71                 for(int i=MatrixCounter;i--;) {
    72                         if (RowCounter != NULL) {
    73                                         for(int j=RowCounter[i]+1;j--;)
    74                                                 Free((void **)&Matrix[i][j], "MatrixContainer::~MatrixContainer: *Matrix[][]");
    75                                 Free((void **)&Matrix[i], "MatrixContainer::~MatrixContainer: **Matrix[]");
    76                         }
    77                 }
    78                 if ((RowCounter != NULL) && (Matrix[MatrixCounter] != NULL))
    79                         for(int j=RowCounter[MatrixCounter]+1;j--;)
    80                                 Free((void **)&Matrix[MatrixCounter][j], "MatrixContainer::~MatrixContainer: *Matrix[MatrixCounter][]");
    81                 if (MatrixCounter != 0)
    82                         Free((void **)&Matrix[MatrixCounter], "MatrixContainer::~MatrixContainer: **Matrix[MatrixCounter]");
    83                 Free((void **)&Matrix, "MatrixContainer::~MatrixContainer: ***Matrix");
    84         }
    85         if (Indices != NULL)
    86                 for(int i=MatrixCounter+1;i--;) {
    87                         Free((void **)&Indices[i], "MatrixContainer::~MatrixContainer: *Indices[]");
    88                 }
    89         Free((void **)&Indices, "MatrixContainer::~MatrixContainer: **Indices");
    90 
    91         Free((void **)&Header, "MatrixContainer::~MatrixContainer: *Header");
    92         Free((void **)&RowCounter, "MatrixContainer::~MatrixContainer: *RowCounter");
     70  if (Matrix != NULL) {
     71    for(int i=MatrixCounter;i--;) {
     72      if (RowCounter != NULL) {
     73          for(int j=RowCounter[i]+1;j--;)
     74            Free((void **)&Matrix[i][j], "MatrixContainer::~MatrixContainer: *Matrix[][]");
     75        Free((void **)&Matrix[i], "MatrixContainer::~MatrixContainer: **Matrix[]");
     76      }
     77    }
     78    if ((RowCounter != NULL) && (Matrix[MatrixCounter] != NULL))
     79      for(int j=RowCounter[MatrixCounter]+1;j--;)
     80        Free((void **)&Matrix[MatrixCounter][j], "MatrixContainer::~MatrixContainer: *Matrix[MatrixCounter][]");
     81    if (MatrixCounter != 0)
     82      Free((void **)&Matrix[MatrixCounter], "MatrixContainer::~MatrixContainer: **Matrix[MatrixCounter]");
     83    Free((void **)&Matrix, "MatrixContainer::~MatrixContainer: ***Matrix");
     84  }
     85  if (Indices != NULL)
     86    for(int i=MatrixCounter+1;i--;) {
     87      Free((void **)&Indices[i], "MatrixContainer::~MatrixContainer: *Indices[]");
     88    }
     89  Free((void **)&Indices, "MatrixContainer::~MatrixContainer: **Indices");
     90
     91  Free((void **)&Header, "MatrixContainer::~MatrixContainer: *Header");
     92  Free((void **)&RowCounter, "MatrixContainer::~MatrixContainer: *RowCounter");
    9393};
    9494
    9595
    9696/** Parsing a number of matrices.
    97  *              -# open the matrix file
    98  *              -# skip some lines (\a skiplines)
    99  *              -# scan header lines for number of columns
    100  *              -# scan lines for number of rows
    101  *              -# allocate matrix
    102  *              -# loop over found column and row counts and parse in each entry
     97 *    -# open the matrix file
     98 *    -# skip some lines (\a skiplines)
     99 *    -# scan header lines for number of columns
     100 *    -# scan lines for number of rows
     101 *    -# allocate matrix
     102 *    -# loop over found column and row counts and parse in each entry
    103103 * \param *name directory with files
    104104 * \param skiplines number of inital lines to skip
    … …  
    109109bool MatrixContainer::ParseMatrix(const char *name, int skiplines, int skipcolumns, int MatrixNr)
    110110{
    111         ifstream input;
    112         stringstream line;
    113         string token;
    114         char filename[1023];
    115 
    116         input.open(name, ios::in);
    117         //cout << "Opening " << name << " ... " << input << endl;
    118         if (input == NULL) {
    119                 cerr << endl << "Unable to open " << name << ", is the directory correct?" << endl;
    120                 return false;
    121         }
    122 
    123         // skip some initial lines
    124         for (int m=skiplines+1;m--;)
    125                 input.getline(Header, 1023);
    126 
    127         // scan header for number of columns
    128         line.str(Header);
    129         for(int k=skipcolumns;k--;)
    130                 line >> Header;
    131         //cout << line.str() << endl;
    132         ColumnCounter=0;
    133         while ( getline(line,token, '\t') ) {
    134                 if (token.length() > 0)
    135                         ColumnCounter++;
    136         }
    137         //cout << line.str() << endl;
    138         //cout << "ColumnCounter: " << ColumnCounter << "." << endl;
    139         if (ColumnCounter == 0)
    140                 cerr << "ColumnCounter: " << ColumnCounter << " from file " << name << ", this is probably an error!" << endl;
    141 
    142         // scan rest for number of rows/lines
    143         RowCounter[MatrixNr]=-1;                // counts one line too much
    144         while (!input.eof()) {
    145                 input.getline(filename, 1023);
    146                 //cout << "Comparing: " << strncmp(filename,"MeanForce",9) << endl;
    147                 RowCounter[MatrixNr]++; // then line was not last MeanForce
    148                 if (strncmp(filename,"MeanForce", 9) == 0) {
    149                         break;
    150                 }
    151         }
    152         //cout << "RowCounter[" << MatrixNr << "]: " << RowCounter[MatrixNr] << " from file " << name << "." << endl;
    153         if (RowCounter[MatrixNr] == 0)
    154                 cerr << "RowCounter[" << MatrixNr << "]: " << RowCounter[MatrixNr] << " from file " << name << ", this is probably an error!" << endl;
    155 
    156         // allocate matrix if it's not zero dimension in one direction
    157         if ((ColumnCounter > 0) && (RowCounter[MatrixNr] > -1)) {
    158                 Matrix[MatrixNr] = (double **) Malloc(sizeof(double *)*(RowCounter[MatrixNr]+1), "MatrixContainer::ParseFragmentMatrix: **Matrix[]");
    159 
    160                 // parse in each entry for this matrix
    161                 input.clear();
    162                 input.seekg(ios::beg);
    163                 for (int m=skiplines+1;m--;)
    164                         input.getline(Header, 1023);            // skip header
    165                 line.str(Header);
    166                 for(int k=skipcolumns;k--;)     // skip columns in header too
    167                         line >> filename;
    168                 strncpy(Header, line.str().c_str(), 1023);
    169                 for(int j=0;j<RowCounter[MatrixNr];j++) {
    170                         Matrix[MatrixNr][j] = (double *) Malloc(sizeof(double)*ColumnCounter, "MatrixContainer::ParseFragmentMatrix: *Matrix[][]");
    171                         input.getline(filename, 1023);
    172                         stringstream lines(filename);
    173                         //cout << "Matrix at level " << j << ":";// << filename << endl;
    174                         for(int k=skipcolumns;k--;)
    175                                 lines >> filename;
    176                         for(int k=0;(k<ColumnCounter) && (!lines.eof());k++) {
    177                                 lines >> Matrix[MatrixNr][j][k];
    178                                 //cout << " " << setprecision(2) << Matrix[MatrixNr][j][k];;
    179                         }
    180                         //cout << endl;
    181                         Matrix[MatrixNr][ RowCounter[MatrixNr] ] = (double *) Malloc(sizeof(double)*ColumnCounter, "MatrixContainer::ParseFragmentMatrix: *Matrix[RowCounter[MatrixNr]][]");
    182                         for(int j=ColumnCounter;j--;)
    183                                 Matrix[MatrixNr][ RowCounter[MatrixNr] ][j] = 0.;
    184                 }
    185         } else {
    186                 cerr << "ERROR: Matrix nr. " << MatrixNr << " has column and row count of (" << ColumnCounter << "," << RowCounter[MatrixNr] << "), could not allocate nor parse!" << endl;
    187         }
    188         input.close();
    189         return true;
     111  ifstream input;
     112  stringstream line;
     113  string token;
     114  char filename[1023];
     115
     116  input.open(name, ios::in);
     117  //cout << "Opening " << name << " ... "  << input << endl;
     118  if (input == NULL) {
     119    cerr << endl << "Unable to open " << name << ", is the directory correct?" << endl;
     120    return false;
     121  }
     122
     123  // skip some initial lines
     124  for (int m=skiplines+1;m--;)
     125    input.getline(Header, 1023);
     126
     127  // scan header for number of columns
     128  line.str(Header);
     129  for(int k=skipcolumns;k--;)
     130    line >> Header;
     131  //cout << line.str() << endl;
     132  ColumnCounter=0;
     133  while ( getline(line,token, '\t') ) {
     134    if (token.length() > 0)
     135      ColumnCounter++;
     136  }
     137  //cout << line.str() << endl;
     138  //cout << "ColumnCounter: " << ColumnCounter << "." << endl;
     139  if (ColumnCounter == 0)
     140    cerr << "ColumnCounter: " << ColumnCounter << " from file " << name << ", this is probably an error!" << endl;
     141
     142  // scan rest for number of rows/lines
     143  RowCounter[MatrixNr]=-1;    // counts one line too much
     144  while (!input.eof()) {
     145    input.getline(filename, 1023);
     146    //cout << "Comparing: " << strncmp(filename,"MeanForce",9) << endl;
     147    RowCounter[MatrixNr]++; // then line was not last MeanForce
     148    if (strncmp(filename,"MeanForce", 9) == 0) {
     149      break;
     150    }
     151  }
     152  //cout << "RowCounter[" << MatrixNr << "]: " << RowCounter[MatrixNr] << " from file " << name << "." << endl;
     153  if (RowCounter[MatrixNr] == 0)
     154    cerr << "RowCounter[" << MatrixNr << "]: " << RowCounter[MatrixNr] << " from file " << name << ", this is probably an error!" << endl;
     155
     156  // allocate matrix if it's not zero dimension in one direction
     157  if ((ColumnCounter > 0) && (RowCounter[MatrixNr] > -1)) {
     158    Matrix[MatrixNr] = (double **) Malloc(sizeof(double *)*(RowCounter[MatrixNr]+1), "MatrixContainer::ParseFragmentMatrix: **Matrix[]");
     159
     160    // parse in each entry for this matrix
     161    input.clear();
     162    input.seekg(ios::beg);
     163    for (int m=skiplines+1;m--;)
     164      input.getline(Header, 1023);    // skip header
     165    line.str(Header);
     166    for(int k=skipcolumns;k--;)  // skip columns in header too
     167      line >> filename;
     168    strncpy(Header, line.str().c_str(), 1023);
     169    for(int j=0;j<RowCounter[MatrixNr];j++) {
     170      Matrix[MatrixNr][j] = (double *) Malloc(sizeof(double)*ColumnCounter, "MatrixContainer::ParseFragmentMatrix: *Matrix[][]");
     171      input.getline(filename, 1023);
     172      stringstream lines(filename);
     173      //cout << "Matrix at level " << j << ":";// << filename << endl;
     174      for(int k=skipcolumns;k--;)
     175        lines >> filename;
     176      for(int k=0;(k<ColumnCounter) && (!lines.eof());k++) {
     177        lines >> Matrix[MatrixNr][j][k];
     178        //cout << " " << setprecision(2) << Matrix[MatrixNr][j][k];;
     179      }
     180      //cout << endl;
     181      Matrix[MatrixNr][ RowCounter[MatrixNr] ] = (double *) Malloc(sizeof(double)*ColumnCounter, "MatrixContainer::ParseFragmentMatrix: *Matrix[RowCounter[MatrixNr]][]");
     182      for(int j=ColumnCounter;j--;)
     183        Matrix[MatrixNr][ RowCounter[MatrixNr] ][j] = 0.;
     184    }
     185  } else {
     186    cerr << "ERROR: Matrix nr. " << MatrixNr << " has column and row count of (" << ColumnCounter << "," << RowCounter[MatrixNr] << "), could not allocate nor parse!" << endl;
     187  }
     188  input.close();
     189  return true;
    190190};
    191191
    … …  
    193193 * -# First, count the number of matrices by counting lines in KEYSETFILE
    194194 * -# Then,
    195  *              -# construct the fragment number
    196  *              -# open the matrix file
    197  *              -# skip some lines (\a skiplines)
    198  *              -# scan header lines for number of columns
    199  *              -# scan lines for number of rows
    200  *              -# allocate matrix
    201  *              -# loop over found column and row counts and parse in each entry
     195 *    -# construct the fragment number
     196 *    -# open the matrix file
     197 *    -# skip some lines (\a skiplines)
     198 *    -# scan header lines for number of columns
     199 *    -# scan lines for number of rows
     200 *    -# allocate matrix
     201 *    -# loop over found column and row counts and parse in each entry
    202202 * -# Finally, allocate one additional matrix (\a MatrixCounter) containing combined or temporary values
    203203 * \param *name directory with files
    … …  
    210210bool MatrixContainer::ParseFragmentMatrix(char *name, char *prefix, string suffix, int skiplines, int skipcolumns)
    211211{
    212         char filename[1023];
    213         ifstream input;
    214         char *FragmentNumber = NULL;
    215         stringstream file;
    216         string token;
    217 
    218         // count the number of matrices
    219         MatrixCounter = -1; // we count one too much
    220         file << name << FRAGMENTPREFIX << KEYSETFILE;
    221         input.open(file.str().c_str(), ios::in);
    222         if (input == NULL) {
    223                 cout << endl << "Unable to open " << file.str() << ", is the directory correct?" << endl;
    224                 return false;
    225         }
    226         while (!input.eof()) {
    227                 input.getline(filename, 1023);
    228                 stringstream zeile(filename);
    229                 MatrixCounter++;
    230         }
    231         input.close();
    232         cout << "Determined " << MatrixCounter << " fragments." << endl;
    233 
    234         cout << "Parsing through each fragment and retrieving " << prefix << suffix << "." << endl;
    235         Matrix = (double ***) ReAlloc(Matrix, sizeof(double **)*(MatrixCounter+1), "MatrixContainer::ParseFragmentMatrix: ***Matrix"); // one more each for the total molecule
    236         RowCounter = (int *) ReAlloc(RowCounter, sizeof(int)*(MatrixCounter+1), "MatrixContainer::ParseFragmentMatrix: *RowCounter");
    237         for(int i=MatrixCounter+1;i--;) {
    238                 Matrix[i] = NULL;
    239                 RowCounter[i] = -1;
    240         }
    241         for(int i=0; i < MatrixCounter;i++) {
    242                 // open matrix file
    243                 FragmentNumber = FixedDigitNumber(MatrixCounter, i);
    244                 file.str(" ");
    245                 file << name << FRAGMENTPREFIX << FragmentNumber << prefix << suffix;
    246                 if (!ParseMatrix(file.str().c_str(), skiplines, skipcolumns, i))
    247                         return false;
    248                 Free((void **)&FragmentNumber, "MatrixContainer::ParseFragmentMatrix: *FragmentNumber");
    249         }
    250         return true;
     212  char filename[1023];
     213  ifstream input;
     214  char *FragmentNumber = NULL;
     215  stringstream file;
     216  string token;
     217
     218  // count the number of matrices
     219  MatrixCounter = -1; // we count one too much
     220  file << name << FRAGMENTPREFIX << KEYSETFILE;
     221  input.open(file.str().c_str(), ios::in);
     222  if (input == NULL) {
     223    cout << endl << "Unable to open " << file.str() << ", is the directory correct?" << endl;
     224    return false;
     225  }
     226  while (!input.eof()) {
     227    input.getline(filename, 1023);
     228    stringstream zeile(filename);
     229    MatrixCounter++;
     230  }
     231  input.close();
     232  cout << "Determined " << MatrixCounter << " fragments." << endl;
     233
     234  cout << "Parsing through each fragment and retrieving " << prefix << suffix << "." << endl;
     235  Matrix = (double ***) ReAlloc(Matrix, sizeof(double **)*(MatrixCounter+1), "MatrixContainer::ParseFragmentMatrix: ***Matrix"); // one more each for the total molecule
     236  RowCounter = (int *) ReAlloc(RowCounter, sizeof(int)*(MatrixCounter+1), "MatrixContainer::ParseFragmentMatrix: *RowCounter");
     237  for(int i=MatrixCounter+1;i--;) {
     238    Matrix[i] = NULL;
     239    RowCounter[i] = -1;
     240  }
     241  for(int i=0; i < MatrixCounter;i++) {
     242    // open matrix file
     243    FragmentNumber = FixedDigitNumber(MatrixCounter, i);
     244    file.str(" ");
     245    file << name << FRAGMENTPREFIX << FragmentNumber << prefix << suffix;
     246    if (!ParseMatrix(file.str().c_str(), skiplines, skipcolumns, i))
     247      return false;
     248    Free((void **)&FragmentNumber, "MatrixContainer::ParseFragmentMatrix: *FragmentNumber");
     249  }
     250  return true;
    251251};
    252252
    … …  
    260260bool MatrixContainer::AllocateMatrix(char *GivenHeader, int MCounter, int *RCounter, int CCounter)
    261261{
    262         Header = (char *) Malloc(sizeof(char)*1024, "MatrixContainer::ParseFragmentMatrix: *EnergyHeader");
    263         strncpy(Header, GivenHeader, 1023);
    264 
    265         MatrixCounter = MCounter;
    266         ColumnCounter = CCounter;
    267         Matrix = (double ***) Malloc(sizeof(double **)*(MatrixCounter+1), "MatrixContainer::ParseFragmentMatrix: ***Matrix"); // one more each for the total molecule
    268         RowCounter = (int *) Malloc(sizeof(int)*(MatrixCounter+1), "MatrixContainer::ParseFragmentMatrix: *RowCounter");
    269         for(int i=MatrixCounter+1;i--;) {
    270                 RowCounter[i] = RCounter[i];
    271                 Matrix[i] = (double **) Malloc(sizeof(double *)*(RowCounter[i]+1), "MatrixContainer::ParseFragmentMatrix: **Matrix[]");
    272                 for(int j=RowCounter[i]+1;j--;) {
    273                         Matrix[i][j] = (double *) Malloc(sizeof(double)*ColumnCounter, "MatrixContainer::ParseFragmentMatrix: *Matrix[][]");
    274                         for(int k=ColumnCounter;k--;)
    275                                 Matrix[i][j][k] = 0.;
    276                 }
    277         }
    278         return true;
     262  Header = (char *) Malloc(sizeof(char)*1024, "MatrixContainer::ParseFragmentMatrix: *EnergyHeader");
     263  strncpy(Header, GivenHeader, 1023);
     264
     265  MatrixCounter = MCounter;
     266  ColumnCounter = CCounter;
     267  Matrix = (double ***) Malloc(sizeof(double **)*(MatrixCounter+1), "MatrixContainer::ParseFragmentMatrix: ***Matrix"); // one more each for the total molecule
     268  RowCounter = (int *) Malloc(sizeof(int)*(MatrixCounter+1), "MatrixContainer::ParseFragmentMatrix: *RowCounter");
     269  for(int i=MatrixCounter+1;i--;) {
     270    RowCounter[i] = RCounter[i];
     271    Matrix[i] = (double **) Malloc(sizeof(double *)*(RowCounter[i]+1), "MatrixContainer::ParseFragmentMatrix: **Matrix[]");
     272    for(int j=RowCounter[i]+1;j--;) {
     273      Matrix[i][j] = (double *) Malloc(sizeof(double)*ColumnCounter, "MatrixContainer::ParseFragmentMatrix: *Matrix[][]");
     274      for(int k=ColumnCounter;k--;)
     275        Matrix[i][j][k] = 0.;
     276    }
     277  }
     278  return true;
    279279};
    280280
    … …  
    284284bool MatrixContainer::ResetMatrix()
    285285{
    286         for(int i=MatrixCounter+1;i--;)
    287                 for(int j=RowCounter[i]+1;j--;)
    288                         for(int k=ColumnCounter;k--;)
    289                                 Matrix[i][j][k] = 0.;
    290         return true;
     286  for(int i=MatrixCounter+1;i--;)
     287    for(int j=RowCounter[i]+1;j--;)
     288      for(int k=ColumnCounter;k--;)
     289        Matrix[i][j][k] = 0.;
     290  return true;
    291291};
    292292
    … …  
    296296double MatrixContainer::FindMaxValue()
    297297{
    298         double max = Matrix[0][0][0];
    299         for(int i=MatrixCounter+1;i--;)
    300                 for(int j=RowCounter[i]+1;j--;)
    301                         for(int k=ColumnCounter;k--;)
    302                                 if (fabs(Matrix[i][j][k]) > max)
    303                                         max = fabs(Matrix[i][j][k]);
    304         if (fabs(max) < MYEPSILON)
    305                 max += MYEPSILON;
     298  double max = Matrix[0][0][0];
     299  for(int i=MatrixCounter+1;i--;)
     300    for(int j=RowCounter[i]+1;j--;)
     301      for(int k=ColumnCounter;k--;)
     302        if (fabs(Matrix[i][j][k]) > max)
     303          max = fabs(Matrix[i][j][k]);
     304  if (fabs(max) < MYEPSILON)
     305    max += MYEPSILON;
    306306 return max;
    307307};
    … …  
    312312double MatrixContainer::FindMinValue()
    313313{
    314         double min = Matrix[0][0][0];
    315         for(int i=MatrixCounter+1;i--;)
    316                 for(int j=RowCounter[i]+1;j--;)
    317                         for(int k=ColumnCounter;k--;)
    318                                 if (fabs(Matrix[i][j][k]) < min)
    319                                         min = fabs(Matrix[i][j][k]);
    320         if (fabs(min) < MYEPSILON)
    321                 min += MYEPSILON;
    322         return min;
     314  double min = Matrix[0][0][0];
     315  for(int i=MatrixCounter+1;i--;)
     316    for(int j=RowCounter[i]+1;j--;)
     317      for(int k=ColumnCounter;k--;)
     318        if (fabs(Matrix[i][j][k]) < min)
     319          min = fabs(Matrix[i][j][k]);
     320  if (fabs(min) < MYEPSILON)
     321    min += MYEPSILON;
     322  return min;
    323323};
    324324
    … …  
    330330bool MatrixContainer::SetLastMatrix(double value, int skipcolumns)
    331331{
    332         for(int j=RowCounter[MatrixCounter]+1;j--;)
    333                 for(int k=skipcolumns;k<ColumnCounter;k++)
    334                         Matrix[MatrixCounter][j][k] = value;
    335         return true;
     332  for(int j=RowCounter[MatrixCounter]+1;j--;)
     333    for(int k=skipcolumns;k<ColumnCounter;k++)
     334      Matrix[MatrixCounter][j][k] = value;
     335  return true;
    336336};
    337337
    … …  
    343343bool MatrixContainer::SetLastMatrix(double **values, int skipcolumns)
    344344{
    345         for(int j=RowCounter[MatrixCounter]+1;j--;)
    346                 for(int k=skipcolumns;k<ColumnCounter;k++)
    347                         Matrix[MatrixCounter][j][k] = values[j][k];
    348         return true;
     345  for(int j=RowCounter[MatrixCounter]+1;j--;)
     346    for(int k=skipcolumns;k<ColumnCounter;k++)
     347      Matrix[MatrixCounter][j][k] = values[j][k];
     348  return true;
    349349};
    350350
    … …  
    358358bool MatrixContainer::SumSubManyBodyTerms(class MatrixContainer &MatrixValues, class KeySetsContainer &KeySet, int Order)
    359359{
    360         // go through each order
    361         for (int CurrentFragment=0;CurrentFragment<KeySet.FragmentsPerOrder[Order];CurrentFragment++) {
    362                 //cout << "Current Fragment is " << CurrentFragment << "/" << KeySet.OrderSet[Order][CurrentFragment] << "." << endl;
    363                 // then go per order through each suborder and pick together all the terms that contain this fragment
    364                 for(int SubOrder=0;SubOrder<=Order;SubOrder++) { // go through all suborders up to the desired order
    365                         for (int j=0;j<KeySet.FragmentsPerOrder[SubOrder];j++) { // go through all possible fragments of size suborder
    366                                 if (KeySet.Contains(KeySet.OrderSet[Order][CurrentFragment], KeySet.OrderSet[SubOrder][j])) {
    367                                         //cout << "Current other fragment is " << j << "/" << KeySet.OrderSet[SubOrder][j] << "." << endl;
    368                                         // if the fragment's indices are all in the current fragment
    369                                         for(int k=0;k<RowCounter[ KeySet.OrderSet[SubOrder][j] ];k++) { // go through all atoms in this fragment
    370                                                 int m = MatrixValues.Indices[ KeySet.OrderSet[SubOrder][j] ][k];
    371                                                 //cout << "Current index is " << k << "/" << m << "." << endl;
    372                                                 if (m != -1) { // if it's not an added hydrogen
    373                                                         for (int l=0;l<RowCounter[ KeySet.OrderSet[Order][CurrentFragment] ];l++) { // look for the corresponding index in the current fragment
    374                                                                 //cout << "Comparing " << m << " with " << MatrixValues.Indices[ KeySet.OrderSet[Order][CurrentFragment] ][l] << "." << endl;
    375                                                                 if (m == MatrixValues.Indices[ KeySet.OrderSet[Order][CurrentFragment] ][l]) {
    376                                                                         m = l;
    377                                                                         break;
    378                                                                 }
    379                                                         }
    380                                                         //cout << "Corresponding index in CurrentFragment is " << m << "." << endl;
    381                                                         if (m > RowCounter[ KeySet.OrderSet[Order][CurrentFragment] ]) {
    382                                                                 cerr << "In fragment No. " << KeySet.OrderSet[Order][CurrentFragment]    << " current force index " << m << " is greater than " << RowCounter[ KeySet.OrderSet[Order][CurrentFragment] ] << "!" << endl;
    383                                                                 return false;
    384                                                         }
    385                                                         if (Order == SubOrder) { // equal order is always copy from Energies
    386                                                                 for(int l=ColumnCounter;l--;) // then adds/subtract each column
    387                                                                         Matrix[ KeySet.OrderSet[Order][CurrentFragment] ][m][l] += MatrixValues.Matrix[ KeySet.OrderSet[SubOrder][j] ][k][l];
    388                                                         } else {
    389                                                                 for(int l=ColumnCounter;l--;)
    390                                                                         Matrix[ KeySet.OrderSet[Order][CurrentFragment] ][m][l] -= Matrix[ KeySet.OrderSet[SubOrder][j] ][k][l];
    391                                                         }
    392                                                 }
    393                                                 //if ((ColumnCounter>1) && (RowCounter[0]-1 >= 1))
    394                                                  //cout << "Fragments[ KeySet.OrderSet[" << Order << "][" << CurrentFragment << "]=" << KeySet.OrderSet[Order][CurrentFragment] << " ][" << RowCounter[0]-1 << "][" << 1 << "] = " <<   Matrix[ KeySet.OrderSet[Order][CurrentFragment] ][RowCounter[0]-1][1] << endl;
    395                                         }
    396                                 } else {
    397                                         //cout << "Fragment " << KeySet.OrderSet[SubOrder][j] << " is not contained in fragment " << KeySet.OrderSet[Order][CurrentFragment] << "." << endl;
    398                                 }
    399                         }
    400                 }
    401          //cout << "Final Fragments[ KeySet.OrderSet[" << Order << "][" << CurrentFragment << "]=" << KeySet.OrderSet[Order][CurrentFragment] << " ][" << KeySet.AtomCounter[0]-1 << "][" << 1 << "] = " <<     Matrix[ KeySet.OrderSet[Order][CurrentFragment] ][KeySet.AtomCounter[0]-1][1] << endl;
    402         }
    403 
    404         return true;
     360  // go through each order
     361  for (int CurrentFragment=0;CurrentFragment<KeySet.FragmentsPerOrder[Order];CurrentFragment++) {
     362    //cout << "Current Fragment is " << CurrentFragment << "/" << KeySet.OrderSet[Order][CurrentFragment] << "." << endl;
     363    // then go per order through each suborder and pick together all the terms that contain this fragment
     364    for(int SubOrder=0;SubOrder<=Order;SubOrder++) { // go through all suborders up to the desired order
     365      for (int j=0;j<KeySet.FragmentsPerOrder[SubOrder];j++) { // go through all possible fragments of size suborder
     366        if (KeySet.Contains(KeySet.OrderSet[Order][CurrentFragment], KeySet.OrderSet[SubOrder][j])) {
     367          //cout << "Current other fragment is " << j << "/" << KeySet.OrderSet[SubOrder][j] << "." << endl;
     368          // if the fragment's indices are all in the current fragment
     369          for(int k=0;k<RowCounter[ KeySet.OrderSet[SubOrder][j] ];k++) { // go through all atoms in this fragment
     370            int m = MatrixValues.Indices[ KeySet.OrderSet[SubOrder][j] ][k];
     371            //cout << "Current index is " << k << "/" << m << "." << endl;
     372            if (m != -1) { // if it's not an added hydrogen
     373              for (int l=0;l<RowCounter[ KeySet.OrderSet[Order][CurrentFragment] ];l++) { // look for the corresponding index in the current fragment
     374                //cout << "Comparing " << m << " with " << MatrixValues.Indices[ KeySet.OrderSet[Order][CurrentFragment] ][l] << "." << endl;
     375                if (m == MatrixValues.Indices[ KeySet.OrderSet[Order][CurrentFragment] ][l]) {
     376                  m = l;
     377                  break;
     378                }
     379              }
     380              //cout << "Corresponding index in CurrentFragment is " << m << "." << endl;
     381              if (m > RowCounter[ KeySet.OrderSet[Order][CurrentFragment] ]) {
     382                cerr << "In fragment No. " << KeySet.OrderSet[Order][CurrentFragment]  << " current force index " << m << " is greater than " << RowCounter[ KeySet.OrderSet[Order][CurrentFragment] ] << "!" << endl;
     383                return false;
     384              }
     385              if (Order == SubOrder) { // equal order is always copy from Energies
     386                for(int l=ColumnCounter;l--;) // then adds/subtract each column
     387                  Matrix[ KeySet.OrderSet[Order][CurrentFragment] ][m][l] += MatrixValues.Matrix[ KeySet.OrderSet[SubOrder][j] ][k][l];
     388              } else {
     389                for(int l=ColumnCounter;l--;)
     390                  Matrix[ KeySet.OrderSet[Order][CurrentFragment] ][m][l] -= Matrix[ KeySet.OrderSet[SubOrder][j] ][k][l];
     391              }
     392            }
     393            //if ((ColumnCounter>1) && (RowCounter[0]-1 >= 1))
     394             //cout << "Fragments[ KeySet.OrderSet[" << Order << "][" << CurrentFragment << "]=" << KeySet.OrderSet[Order][CurrentFragment] << " ][" << RowCounter[0]-1 << "][" << 1 << "] = " <<  Matrix[ KeySet.OrderSet[Order][CurrentFragment] ][RowCounter[0]-1][1] << endl;
     395          }
     396        } else {
     397          //cout << "Fragment " << KeySet.OrderSet[SubOrder][j] << " is not contained in fragment " << KeySet.OrderSet[Order][CurrentFragment] << "." << endl;
     398        }
     399      }
     400    }
     401   //cout << "Final Fragments[ KeySet.OrderSet[" << Order << "][" << CurrentFragment << "]=" << KeySet.OrderSet[Order][CurrentFragment] << " ][" << KeySet.AtomCounter[0]-1 << "][" << 1 << "] = " <<  Matrix[ KeySet.OrderSet[Order][CurrentFragment] ][KeySet.AtomCounter[0]-1][1] << endl;
     402  }
     403
     404  return true;
    405405};
    406406
    … …  
    412412bool MatrixContainer::WriteTotalFragments(const char *name, const char *prefix)
    413413{
    414         ofstream output;
    415         char *FragmentNumber = NULL;
    416 
    417         cout << "Writing fragment files." << endl;
    418         for(int i=0;i<MatrixCounter;i++) {
    419                 stringstream line;
    420                 FragmentNumber = FixedDigitNumber(MatrixCounter, i);
    421                 line << name << FRAGMENTPREFIX << FragmentNumber << "/" << prefix;
    422                 Free((void **)&FragmentNumber, "*FragmentNumber");
    423                 output.open(line.str().c_str(), ios::out);
    424                 if (output == NULL) {
    425                         cerr << "Unable to open output energy file " << line.str() << "!" << endl;
    426                         return false;
    427                 }
    428                 output << Header << endl;
    429                 for(int j=0;j<RowCounter[i];j++) {
    430                         for(int k=0;k<ColumnCounter;k++)
    431                                 output << scientific << Matrix[i][j][k] << "\t";
    432                         output << endl;
    433                 }
    434                 output.close();
    435         }
    436         return true;
     414  ofstream output;
     415  char *FragmentNumber = NULL;
     416
     417  cout << "Writing fragment files." << endl;
     418  for(int i=0;i<MatrixCounter;i++) {
     419    stringstream line;
     420    FragmentNumber = FixedDigitNumber(MatrixCounter, i);
     421    line << name << FRAGMENTPREFIX << FragmentNumber << "/" << prefix;
     422    Free((void **)&FragmentNumber, "*FragmentNumber");
     423    output.open(line.str().c_str(), ios::out);
     424    if (output == NULL) {
     425      cerr << "Unable to open output energy file " << line.str() << "!" << endl;
     426      return false;
     427    }
     428    output << Header << endl;
     429    for(int j=0;j<RowCounter[i];j++) {
     430      for(int k=0;k<ColumnCounter;k++)
     431        output << scientific << Matrix[i][j][k] << "\t";
     432      output << endl;
     433    }
     434    output.close();
     435  }
     436  return true;
    437437};
    438438
    … …  
    445445bool MatrixContainer::WriteLastMatrix(const char *name, const char *prefix, const char *suffix)
    446446{
    447         ofstream output;
    448         stringstream line;
    449 
    450         cout << "Writing matrix values of " << suffix << "." << endl;
    451         line << name << prefix << suffix;
    452         output.open(line.str().c_str(), ios::out);
    453         if (output == NULL) {
    454                 cerr << "Unable to open output matrix file " << line.str() << "!" << endl;
    455                 return false;
    456         }
    457         output << Header << endl;
    458         for(int j=0;j<RowCounter[MatrixCounter];j++) {
    459                 for(int k=0;k<ColumnCounter;k++)
    460                         output << scientific << Matrix[MatrixCounter][j][k] << "\t";
    461                 output << endl;
    462         }
    463         output.close();
    464         return true;
     447  ofstream output;
     448  stringstream line;
     449
     450  cout << "Writing matrix values of " << suffix << "." << endl;
     451  line << name << prefix << suffix;
     452  output.open(line.str().c_str(), ios::out);
     453  if (output == NULL) {
     454    cerr << "Unable to open output matrix file " << line.str() << "!" << endl;
     455    return false;
     456  }
     457  output << Header << endl;
     458  for(int j=0;j<RowCounter[MatrixCounter];j++) {
     459    for(int k=0;k<ColumnCounter;k++)
     460      output << scientific << Matrix[MatrixCounter][j][k] << "\t";
     461    output << endl;
     462  }
     463  output.close();
     464  return true;
    465465};
    466466
    … …  
    473473bool EnergyMatrix::ParseIndices()
    474474{
    475         cout << "Parsing energy indices." << endl;
    476         Indices = (int **) Malloc(sizeof(int *)*(MatrixCounter+1), "EnergyMatrix::ParseIndices: **Indices");
    477         for(int i=MatrixCounter+1;i--;) {
    478                 Indices[i] = (int *) Malloc(sizeof(int)*RowCounter[i], "EnergyMatrix::ParseIndices: *Indices[]");
    479                 for(int j=RowCounter[i];j--;)
    480                         Indices[i][j] = j;
    481         }
    482         return true;
     475  cout << "Parsing energy indices." << endl;
     476  Indices = (int **) Malloc(sizeof(int *)*(MatrixCounter+1), "EnergyMatrix::ParseIndices: **Indices");
     477  for(int i=MatrixCounter+1;i--;) {
     478    Indices[i] = (int *) Malloc(sizeof(int)*RowCounter[i], "EnergyMatrix::ParseIndices: *Indices[]");
     479    for(int j=RowCounter[i];j--;)
     480      Indices[i][j] = j;
     481  }
     482  return true;
    483483};
    484484
    … …  
    494494bool EnergyMatrix::SumSubEnergy(class EnergyMatrix &Fragments, class EnergyMatrix *CorrectionFragments, class KeySetsContainer &KeySet, int Order, double sign)
    495495{
    496         // sum energy
    497         if (CorrectionFragments == NULL)
    498                 for(int i=KeySet.FragmentsPerOrder[Order];i--;)
    499                         for(int j=RowCounter[ KeySet.OrderSet[Order][i] ];j--;)
    500                                 for(int k=ColumnCounter;k--;)
    501                                         Matrix[MatrixCounter][j][k] += sign*Fragments.Matrix[ KeySet.OrderSet[Order][i] ][j][k];
    502         else
    503                 for(int i=KeySet.FragmentsPerOrder[Order];i--;)
    504                         for(int j=RowCounter[ KeySet.OrderSet[Order][i] ];j--;)
    505                                 for(int k=ColumnCounter;k--;)
    506                                         Matrix[MatrixCounter][j][k] += sign*(Fragments.Matrix[ KeySet.OrderSet[Order][i] ][j][k] + CorrectionFragments->Matrix[ KeySet.OrderSet[Order][i] ][j][k]);
    507         return true;
     496  // sum energy
     497  if (CorrectionFragments == NULL)
     498    for(int i=KeySet.FragmentsPerOrder[Order];i--;)
     499      for(int j=RowCounter[ KeySet.OrderSet[Order][i] ];j--;)
     500        for(int k=ColumnCounter;k--;)
     501          Matrix[MatrixCounter][j][k] += sign*Fragments.Matrix[ KeySet.OrderSet[Order][i] ][j][k];
     502  else
     503    for(int i=KeySet.FragmentsPerOrder[Order];i--;)
     504      for(int j=RowCounter[ KeySet.OrderSet[Order][i] ];j--;)
     505        for(int k=ColumnCounter;k--;)
     506          Matrix[MatrixCounter][j][k] += sign*(Fragments.Matrix[ KeySet.OrderSet[Order][i] ][j][k] + CorrectionFragments->Matrix[ KeySet.OrderSet[Order][i] ][j][k]);
     507  return true;
    508508};
    509509
    … …  
    518518bool EnergyMatrix::ParseFragmentMatrix(char *name, char *prefix, string suffix, int skiplines, int skipcolumns)
    519519{
    520         char filename[1024];
    521         bool status = MatrixContainer::ParseFragmentMatrix(name, prefix, suffix, skiplines, skipcolumns);
    522 
    523         if (status) {
    524                 // count maximum of columns
    525                 RowCounter[MatrixCounter] = 0;
    526                 for(int j=0; j < MatrixCounter;j++) // (energy matrix might be bigger than number of atoms in terms of rows)
    527                         if (RowCounter[j] > RowCounter[MatrixCounter])
    528                                 RowCounter[MatrixCounter] = RowCounter[j];
    529                 // allocate last plus one matrix
    530                 cout << "Allocating last plus one matrix with " << (RowCounter[MatrixCounter]+1) << " rows and " << ColumnCounter << " columns." << endl;
    531                 Matrix[MatrixCounter] = (double **) Malloc(sizeof(double *)*(RowCounter[MatrixCounter]+1), "MatrixContainer::ParseFragmentMatrix: **Matrix[]");
    532                 for(int j=0;j<=RowCounter[MatrixCounter];j++)
    533                         Matrix[MatrixCounter][j] = (double *) Malloc(sizeof(double)*ColumnCounter, "MatrixContainer::ParseFragmentMatrix: *Matrix[][]");
    534 
    535                 // try independently to parse global energysuffix file if present
    536                 strncpy(filename, name, 1023);
    537                 strncat(filename, prefix, 1023-strlen(filename));
    538                 strncat(filename, suffix.c_str(), 1023-strlen(filename));
    539                 ParseMatrix(filename, skiplines, skipcolumns, MatrixCounter);
    540         }
    541         return status;
     520  char filename[1024];
     521  bool status = MatrixContainer::ParseFragmentMatrix(name, prefix, suffix, skiplines, skipcolumns);
     522
     523  if (status) {
     524    // count maximum of columns
     525    RowCounter[MatrixCounter] = 0;
     526    for(int j=0; j < MatrixCounter;j++) // (energy matrix might be bigger than number of atoms in terms of rows)
     527      if (RowCounter[j] > RowCounter[MatrixCounter])
     528        RowCounter[MatrixCounter] = RowCounter[j];
     529    // allocate last plus one matrix
     530    cout << "Allocating last plus one matrix with " << (RowCounter[MatrixCounter]+1) << " rows and " << ColumnCounter << " columns." << endl;
     531    Matrix[MatrixCounter] = (double **) Malloc(sizeof(double *)*(RowCounter[MatrixCounter]+1), "MatrixContainer::ParseFragmentMatrix: **Matrix[]");
     532    for(int j=0;j<=RowCounter[MatrixCounter];j++)
     533      Matrix[MatrixCounter][j] = (double *) Malloc(sizeof(double)*ColumnCounter, "MatrixContainer::ParseFragmentMatrix: *Matrix[][]");
     534
     535    // try independently to parse global energysuffix file if present
     536    strncpy(filename, name, 1023);
     537    strncat(filename, prefix, 1023-strlen(filename));
     538    strncat(filename, suffix.c_str(), 1023-strlen(filename));
     539    ParseMatrix(filename, skiplines, skipcolumns, MatrixCounter);
     540  }
     541  return status;
    542542};
    543543
    … …  
    550550bool ForceMatrix::ParseIndices(char *name)
    551551{
    552         ifstream input;
    553         char *FragmentNumber = NULL;
    554         char filename[1023];
    555         stringstream line;
    556 
    557         cout << "Parsing force indices for " << MatrixCounter << " matrices." << endl;
    558         Indices = (int **) Malloc(sizeof(int *)*(MatrixCounter+1), "ForceMatrix::ParseIndices: **Indices");
    559         line << name << FRAGMENTPREFIX << FORCESFILE;
    560         input.open(line.str().c_str(), ios::in);
    561         //cout << "Opening " << line.str() << " ... "   << input << endl;
    562         if (input == NULL) {
    563                 cout << endl << "Unable to open " << line.str() << ", is the directory correct?" << endl;
    564                 return false;
    565         }
    566         for (int i=0;(i<MatrixCounter) && (!input.eof());i++) {
    567                 // get the number of atoms for this fragment
    568                 input.getline(filename, 1023);
    569                 line.str(filename);
    570                 // parse the values
    571                 Indices[i] = (int *) Malloc(sizeof(int)*RowCounter[i], "ForceMatrix::ParseIndices: *Indices[]");
    572                 FragmentNumber = FixedDigitNumber(MatrixCounter, i);
    573                 //cout << FRAGMENTPREFIX << FragmentNumber << "[" << RowCounter[i] << "]:";
    574                 Free((void **)&FragmentNumber, "ForceMatrix::ParseIndices: *FragmentNumber");
    575                 for(int j=0;(j<RowCounter[i]) && (!line.eof());j++) {
    576                         line >> Indices[i][j];
    577                         //cout << " " << Indices[i][j];
    578                 }
    579                 //cout << endl;
    580         }
    581         Indices[MatrixCounter] = (int *) Malloc(sizeof(int)*RowCounter[MatrixCounter], "ForceMatrix::ParseIndices: *Indices[]");
    582         for(int j=RowCounter[MatrixCounter];j--;) {
    583                 Indices[MatrixCounter][j] = j;
    584         }
    585         input.close();
    586         return true;
     552  ifstream input;
     553  char *FragmentNumber = NULL;
     554  char filename[1023];
     555  stringstream line;
     556
     557  cout << "Parsing force indices for " << MatrixCounter << " matrices." << endl;
     558  Indices = (int **) Malloc(sizeof(int *)*(MatrixCounter+1), "ForceMatrix::ParseIndices: **Indices");
     559  line << name << FRAGMENTPREFIX << FORCESFILE;
     560  input.open(line.str().c_str(), ios::in);
     561  //cout << "Opening " << line.str() << " ... "  << input << endl;
     562  if (input == NULL) {
     563    cout << endl << "Unable to open " << line.str() << ", is the directory correct?" << endl;
     564    return false;
     565  }
     566  for (int i=0;(i<MatrixCounter) && (!input.eof());i++) {
     567    // get the number of atoms for this fragment
     568    input.getline(filename, 1023);
     569    line.str(filename);
     570    // parse the values
     571    Indices[i] = (int *) Malloc(sizeof(int)*RowCounter[i], "ForceMatrix::ParseIndices: *Indices[]");
     572    FragmentNumber = FixedDigitNumber(MatrixCounter, i);
     573    //cout << FRAGMENTPREFIX << FragmentNumber << "[" << RowCounter[i] << "]:";
     574    Free((void **)&FragmentNumber, "ForceMatrix::ParseIndices: *FragmentNumber");
     575    for(int j=0;(j<RowCounter[i]) && (!line.eof());j++) {
     576      line >> Indices[i][j];
     577      //cout << " " << Indices[i][j];
     578    }
     579    //cout << endl;
     580  }
     581  Indices[MatrixCounter] = (int *) Malloc(sizeof(int)*RowCounter[MatrixCounter], "ForceMatrix::ParseIndices: *Indices[]");
     582  for(int j=RowCounter[MatrixCounter];j--;) {
     583    Indices[MatrixCounter][j] = j;
     584  }
     585  input.close();
     586  return true;
    587587};
    588588
    … …  
    592592 * \param KeySet KeySetContainer with bond Order and association mapping of each fragment to an order
    593593 * \param Order bond order
    594  *      \param sign +1 or -1
     594 *  \param sign +1 or -1
    595595 * \return true if summing was successful
    596596 */
    597597bool ForceMatrix::SumSubForces(class ForceMatrix &Fragments, class KeySetsContainer &KeySet, int Order, double sign)
    598598{
    599         int FragmentNr;
    600         // sum forces
    601         for(int i=0;i<KeySet.FragmentsPerOrder[Order];i++) {
    602                 FragmentNr = KeySet.OrderSet[Order][i];
    603                 for(int l=0;l<RowCounter[ FragmentNr ];l++) {
    604                         int j = Indices[ FragmentNr ][l];
    605                         if (j > RowCounter[MatrixCounter]) {
    606                                 cerr << "Current force index " << j << " is greater than " << RowCounter[MatrixCounter] << "!" << endl;
    607                                 return false;
    608                         }
    609                         if (j != -1) {
    610                                 //if (j == 0) cout << "Summing onto ion 0, type 0 from fragment " << FragmentNr << ", ion " << l << "." << endl;
    611                                 for(int k=2;k<ColumnCounter;k++)
    612                                         Matrix[MatrixCounter][j][k] += sign*Fragments.Matrix[ FragmentNr ][l][k];
    613                         }
    614                 }
    615         }
    616         return true;
     599  int FragmentNr;
     600  // sum forces
     601  for(int i=0;i<KeySet.FragmentsPerOrder[Order];i++) {
     602    FragmentNr = KeySet.OrderSet[Order][i];
     603    for(int l=0;l<RowCounter[ FragmentNr ];l++) {
     604      int j = Indices[ FragmentNr ][l];
     605      if (j > RowCounter[MatrixCounter]) {
     606        cerr << "Current force index " << j << " is greater than " << RowCounter[MatrixCounter] << "!" << endl;
     607        return false;
     608      }
     609      if (j != -1) {
     610        //if (j == 0) cout << "Summing onto ion 0, type 0 from fragment " << FragmentNr << ", ion " << l << "." << endl;
     611        for(int k=2;k<ColumnCounter;k++)
     612          Matrix[MatrixCounter][j][k] += sign*Fragments.Matrix[ FragmentNr ][l][k];
     613      }
     614    }
     615  }
     616  return true;
    617617};
    618618
    … …  
    628628bool ForceMatrix::ParseFragmentMatrix(char *name, char *prefix, string suffix, int skiplines, int skipcolumns)
    629629{
    630         char filename[1023];
    631         ifstream input;
    632         stringstream file;
    633         int nr;
    634         bool status = MatrixContainer::ParseFragmentMatrix(name, prefix, suffix, skiplines, skipcolumns);
    635 
    636         if (status) {
    637                 // count number of atoms for last plus one matrix
    638                 file << name << FRAGMENTPREFIX << KEYSETFILE;
    639                 input.open(file.str().c_str(), ios::in);
    640                 if (input == NULL) {
    641                         cout << endl << "Unable to open " << file.str() << ", is the directory correct?" << endl;
    642                         return false;
    643                 }
    644                 RowCounter[MatrixCounter] = 0;
    645                 while (!input.eof()) {
    646                         input.getline(filename, 1023);
    647                         stringstream zeile(filename);
    648                         while (!zeile.eof()) {
    649                                 zeile >> nr;
    650                                 //cout << "Current index: " << nr << "." << endl;
    651                                 if (nr > RowCounter[MatrixCounter])
    652                                         RowCounter[MatrixCounter] = nr;
    653                         }
    654                 }
    655                 RowCounter[MatrixCounter]++;            // nr start at 0, count starts at 1
    656                 input.close();
    657 
    658                 // allocate last plus one matrix
    659                 cout << "Allocating last plus one matrix with " << (RowCounter[MatrixCounter]+1) << " rows and " << ColumnCounter << " columns." << endl;
    660                 Matrix[MatrixCounter] = (double **) Malloc(sizeof(double *)*(RowCounter[MatrixCounter]+1), "MatrixContainer::ParseFragmentMatrix: **Matrix[]");
    661                 for(int j=0;j<=RowCounter[MatrixCounter];j++)
    662                         Matrix[MatrixCounter][j] = (double *) Malloc(sizeof(double)*ColumnCounter, "MatrixContainer::ParseFragmentMatrix: *Matrix[][]");
    663 
    664                 // try independently to parse global forcesuffix file if present
    665                 strncpy(filename, name, 1023);
    666                 strncat(filename, prefix, 1023-strlen(filename));
    667                 strncat(filename, suffix.c_str(), 1023-strlen(filename));
    668                 ParseMatrix(filename, skiplines, skipcolumns, MatrixCounter);
    669         }
    670 
    671 
    672         return status;
     630  char filename[1023];
     631  ifstream input;
     632  stringstream file;
     633  int nr;
     634  bool status = MatrixContainer::ParseFragmentMatrix(name, prefix, suffix, skiplines, skipcolumns);
     635
     636  if (status) {
     637    // count number of atoms for last plus one matrix
     638    file << name << FRAGMENTPREFIX << KEYSETFILE;
     639    input.open(file.str().c_str(), ios::in);
     640    if (input == NULL) {
     641      cout << endl << "Unable to open " << file.str() << ", is the directory correct?" << endl;
     642      return false;
     643    }
     644    RowCounter[MatrixCounter] = 0;
     645    while (!input.eof()) {
     646      input.getline(filename, 1023);
     647      stringstream zeile(filename);
     648      while (!zeile.eof()) {
     649        zeile >> nr;
     650        //cout << "Current index: " << nr << "." << endl;
     651        if (nr > RowCounter[MatrixCounter])
     652          RowCounter[MatrixCounter] = nr;
     653      }
     654    }
     655    RowCounter[MatrixCounter]++;    // nr start at 0, count starts at 1
     656    input.close();
     657
     658    // allocate last plus one matrix
     659    cout << "Allocating last plus one matrix with " << (RowCounter[MatrixCounter]+1) << " rows and " << ColumnCounter << " columns." << endl;
     660    Matrix[MatrixCounter] = (double **) Malloc(sizeof(double *)*(RowCounter[MatrixCounter]+1), "MatrixContainer::ParseFragmentMatrix: **Matrix[]");
     661    for(int j=0;j<=RowCounter[MatrixCounter];j++)
     662      Matrix[MatrixCounter][j] = (double *) Malloc(sizeof(double)*ColumnCounter, "MatrixContainer::ParseFragmentMatrix: *Matrix[][]");
     663
     664    // try independently to parse global forcesuffix file if present
     665    strncpy(filename, name, 1023);
     666    strncat(filename, prefix, 1023-strlen(filename));
     667    strncat(filename, suffix.c_str(), 1023-strlen(filename));
     668    ParseMatrix(filename, skiplines, skipcolumns, MatrixCounter);
     669  }
     670
     671
     672  return status;
    673673};
    674674
    … …  
    678678 */
    679679KeySetsContainer::KeySetsContainer() {
    680         KeySets = NULL;
    681         AtomCounter = NULL;
    682         FragmentCounter = 0;
    683         Order = 0;
    684         FragmentsPerOrder = 0;
    685         OrderSet = NULL;
     680  KeySets = NULL;
     681  AtomCounter = NULL;
     682  FragmentCounter = 0;
     683  Order = 0;
     684  FragmentsPerOrder = 0;
     685  OrderSet = NULL;
    686686};
    687687
    … …  
    689689 */
    690690KeySetsContainer::~KeySetsContainer() {
    691         for(int i=FragmentCounter;i--;)
    692                 Free((void **)&KeySets[i], "KeySetsContainer::~KeySetsContainer: *KeySets[]");
    693         for(int i=Order;i--;)
    694                 Free((void **)&OrderSet[i], "KeySetsContainer::~KeySetsContainer: *OrderSet[]");
    695         Free((void **)&KeySets, "KeySetsContainer::~KeySetsContainer: **KeySets");
    696         Free((void **)&OrderSet, "KeySetsContainer::~KeySetsContainer: **OrderSet");
    697         Free((void **)&AtomCounter, "KeySetsContainer::~KeySetsContainer: *AtomCounter");
    698         Free((void **)&FragmentsPerOrder, "KeySetsContainer::~KeySetsContainer: *FragmentsPerOrder");
     691  for(int i=FragmentCounter;i--;)
     692    Free((void **)&KeySets[i], "KeySetsContainer::~KeySetsContainer: *KeySets[]");
     693  for(int i=Order;i--;)
     694    Free((void **)&OrderSet[i], "KeySetsContainer::~KeySetsContainer: *OrderSet[]");
     695  Free((void **)&KeySets, "KeySetsContainer::~KeySetsContainer: **KeySets");
     696  Free((void **)&OrderSet, "KeySetsContainer::~KeySetsContainer: **OrderSet");
     697  Free((void **)&AtomCounter, "KeySetsContainer::~KeySetsContainer: *AtomCounter");
     698  Free((void **)&FragmentsPerOrder, "KeySetsContainer::~KeySetsContainer: *FragmentsPerOrder");
    699699};
    700700
    … …  
    706706 */
    707707bool KeySetsContainer::ParseKeySets(const char *name, const int *ACounter, const int FCounter) {
    708         ifstream input;
    709         char *FragmentNumber = NULL;
    710         stringstream file;
    711         char filename[1023];
    712 
    713         FragmentCounter = FCounter;
    714         cout << "Parsing key sets." << endl;
    715         KeySets = (int **) Malloc(sizeof(int *)*FragmentCounter, "KeySetsContainer::ParseKeySets: **KeySets");
    716         for(int i=FragmentCounter;i--;)
    717                 KeySets[i] = NULL;
    718         file << name << FRAGMENTPREFIX << KEYSETFILE;
    719         input.open(file.str().c_str(), ios::in);
    720         if (input == NULL) {
    721                 cout << endl << "Unable to open " << file.str() << ", is the directory correct?" << endl;
    722                 return false;
    723         }
    724 
    725         AtomCounter = (int *) Malloc(sizeof(int)*FragmentCounter, "KeySetsContainer::ParseKeySets: *RowCounter");
    726         for(int i=0;(i<FragmentCounter) && (!input.eof());i++) {
    727                 stringstream line;
    728                 AtomCounter[i] = ACounter[i];
    729                 // parse the values
    730                 KeySets[i] = (int *) Malloc(sizeof(int)*AtomCounter[i], "KeySetsContainer::ParseKeySets: *KeySets[]");
    731                 for(int j=AtomCounter[i];j--;)
    732                         KeySets[i][j] = -1;
    733                 FragmentNumber = FixedDigitNumber(FragmentCounter, i);
    734                 //cout << FRAGMENTPREFIX << FragmentNumber << "[" << AtomCounter[i] << "]:";
    735                 Free((void **)&FragmentNumber, "KeySetsContainer::ParseKeySets: *FragmentNumber");
    736                 input.getline(filename, 1023);
    737                 line.str(filename);
    738                 for(int j=0;(j<AtomCounter[i]) && (!line.eof());j++) {
    739                         line >> KeySets[i][j];
    740                         //cout << " " << KeySets[i][j];
    741                 }
    742                 //cout << endl;
    743         }
    744         input.close();
    745         return true;
     708  ifstream input;
     709  char *FragmentNumber = NULL;
     710  stringstream file;
     711  char filename[1023];
     712
     713  FragmentCounter = FCounter;
     714  cout << "Parsing key sets." << endl;
     715  KeySets = (int **) Malloc(sizeof(int *)*FragmentCounter, "KeySetsContainer::ParseKeySets: **KeySets");
     716  for(int i=FragmentCounter;i--;)
     717    KeySets[i] = NULL;
     718  file << name << FRAGMENTPREFIX << KEYSETFILE;
     719  input.open(file.str().c_str(), ios::in);
     720  if (input == NULL) {
     721    cout << endl << "Unable to open " << file.str() << ", is the directory correct?" << endl;
     722    return false;
     723  }
     724
     725  AtomCounter = (int *) Malloc(sizeof(int)*FragmentCounter, "KeySetsContainer::ParseKeySets: *RowCounter");
     726  for(int i=0;(i<FragmentCounter) && (!input.eof());i++) {
     727    stringstream line;
     728    AtomCounter[i] = ACounter[i];
     729    // parse the values
     730    KeySets[i] = (int *) Malloc(sizeof(int)*AtomCounter[i], "KeySetsContainer::ParseKeySets: *KeySets[]");
     731    for(int j=AtomCounter[i];j--;)
     732      KeySets[i][j] = -1;
     733    FragmentNumber = FixedDigitNumber(FragmentCounter, i);
     734    //cout << FRAGMENTPREFIX << FragmentNumber << "[" << AtomCounter[i] << "]:";
     735    Free((void **)&FragmentNumber, "KeySetsContainer::ParseKeySets: *FragmentNumber");
     736    input.getline(filename, 1023);
     737    line.str(filename);
     738    for(int j=0;(j<AtomCounter[i]) && (!line.eof());j++) {
     739      line >> KeySets[i][j];
     740      //cout << " " << KeySets[i][j];
     741    }
     742    //cout << endl;
     743  }
     744  input.close();
     745  return true;
    746746};
    747747
    … …  
    751751bool KeySetsContainer::ParseManyBodyTerms()
    752752{
    753         int Counter;
    754 
    755         cout << "Creating Fragment terms." << endl;
    756         // scan through all to determine maximum order
    757         Order=0;
    758         for(int i=FragmentCounter;i--;) {
    759                 Counter=0;
    760                 for(int j=AtomCounter[i];j--;)
    761                         if (KeySets[i][j] != -1)
    762                                 Counter++;
    763                 if (Counter > Order)
    764                         Order = Counter;
    765         }
    766         cout << "Found Order is " << Order << "." << endl;
    767 
    768         // scan through all to determine fragments per order
    769         FragmentsPerOrder = (int *) Malloc(sizeof(int)*Order, "KeySetsContainer::ParseManyBodyTerms: *FragmentsPerOrder");
    770         for(int i=Order;i--;)
    771                 FragmentsPerOrder[i] = 0;
    772         for(int i=FragmentCounter;i--;) {
    773                 Counter=0;
    774                 for(int j=AtomCounter[i];j--;)
    775                         if (KeySets[i][j] != -1)
    776                                 Counter++;
    777                 FragmentsPerOrder[Counter-1]++;
    778         }
    779         for(int i=0;i<Order;i++)
    780                 cout << "Found No. of Fragments of Order " << i+1 << " is " << FragmentsPerOrder[i] << "." << endl;
    781 
    782         // scan through all to gather indices to each order set
    783         OrderSet = (int **) Malloc(sizeof(int *)*Order, "KeySetsContainer::ParseManyBodyTerms: **OrderSet");
    784         for(int i=Order;i--;)
    785                 OrderSet[i] = (int *) Malloc(sizeof(int)*FragmentsPerOrder[i], "KeySetsContainer::ParseManyBodyTermsKeySetsContainer::ParseManyBodyTerms: *OrderSet[]");
    786         for(int i=Order;i--;)
    787                 FragmentsPerOrder[i] = 0;
    788         for(int i=FragmentCounter;i--;) {
    789                 Counter=0;
    790                 for(int j=AtomCounter[i];j--;)
    791                         if (KeySets[i][j] != -1)
    792                                 Counter++;
    793                 OrderSet[Counter-1][FragmentsPerOrder[Counter-1]] = i;
    794                 FragmentsPerOrder[Counter-1]++;
    795         }
    796         cout << "Printing OrderSet." << endl;
    797         for(int i=0;i<Order;i++) {
    798                 for (int j=0;j<FragmentsPerOrder[i];j++) {
    799                         cout << " " << OrderSet[i][j];
    800                 }
    801                 cout << endl;
    802         }
    803         cout << endl;
    804 
    805 
    806         return true;
     753  int Counter;
     754
     755  cout << "Creating Fragment terms." << endl;
     756  // scan through all to determine maximum order
     757  Order=0;
     758  for(int i=FragmentCounter;i--;) {
     759    Counter=0;
     760    for(int j=AtomCounter[i];j--;)
     761      if (KeySets[i][j] != -1)
     762        Counter++;
     763    if (Counter > Order)
     764      Order = Counter;
     765  }
     766  cout << "Found Order is " << Order << "." << endl;
     767
     768  // scan through all to determine fragments per order
     769  FragmentsPerOrder = (int *) Malloc(sizeof(int)*Order, "KeySetsContainer::ParseManyBodyTerms: *FragmentsPerOrder");
     770  for(int i=Order;i--;)
     771    FragmentsPerOrder[i] = 0;
     772  for(int i=FragmentCounter;i--;) {
     773    Counter=0;
     774    for(int j=AtomCounter[i];j--;)
     775      if (KeySets[i][j] != -1)
     776        Counter++;
     777    FragmentsPerOrder[Counter-1]++;
     778  }
     779  for(int i=0;i<Order;i++)
     780    cout << "Found No. of Fragments of Order " << i+1 << " is " << FragmentsPerOrder[i] << "." << endl;
     781
     782  // scan through all to gather indices to each order set
     783  OrderSet = (int **) Malloc(sizeof(int *)*Order, "KeySetsContainer::ParseManyBodyTerms: **OrderSet");
     784  for(int i=Order;i--;)
     785    OrderSet[i] = (int *) Malloc(sizeof(int)*FragmentsPerOrder[i], "KeySetsContainer::ParseManyBodyTermsKeySetsContainer::ParseManyBodyTerms: *OrderSet[]");
     786  for(int i=Order;i--;)
     787    FragmentsPerOrder[i] = 0;
     788  for(int i=FragmentCounter;i--;) {
     789    Counter=0;
     790    for(int j=AtomCounter[i];j--;)
     791      if (KeySets[i][j] != -1)
     792        Counter++;
     793    OrderSet[Counter-1][FragmentsPerOrder[Counter-1]] = i;
     794    FragmentsPerOrder[Counter-1]++;
     795  }
     796  cout << "Printing OrderSet." << endl;
     797  for(int i=0;i<Order;i++) {
     798    for (int j=0;j<FragmentsPerOrder[i];j++) {
     799      cout << " " << OrderSet[i][j];
     800    }
     801    cout << endl;
     802  }
     803  cout << endl;
     804
     805
     806  return true;
    807807};
    808808
    … …  
    814814bool KeySetsContainer::Contains(const int GreaterSet, const int SmallerSet)
    815815{
    816         bool result = true;
    817         bool intermediate;
    818         if ((GreaterSet < 0) || (SmallerSet < 0) || (GreaterSet > FragmentCounter) || (SmallerSet > FragmentCounter)) // index out of bounds
    819                 return false;
    820         for(int i=AtomCounter[SmallerSet];i--;) {
    821                 intermediate = false;
    822                 for (int j=AtomCounter[GreaterSet];j--;)
    823                         intermediate = (intermediate || ((KeySets[SmallerSet][i] == KeySets[GreaterSet][j]) || (KeySets[SmallerSet][i] == -1)));
    824                 result = result && intermediate;
    825         }
    826 
    827         return result;
     816  bool result = true;
     817  bool intermediate;
     818  if ((GreaterSet < 0) || (SmallerSet < 0) || (GreaterSet > FragmentCounter) || (SmallerSet > FragmentCounter)) // index out of bounds
     819    return false;
     820  for(int i=AtomCounter[SmallerSet];i--;) {
     821    intermediate = false;
     822    for (int j=AtomCounter[GreaterSet];j--;)
     823      intermediate = (intermediate || ((KeySets[SmallerSet][i] == KeySets[GreaterSet][j]) || (KeySets[SmallerSet][i] == -1)));
     824    result = result && intermediate;
     825  }
     826
     827  return result;
    828828};
    829829
Note: See TracChangeset for help on using the changeset viewer.