Changeset f731ae


Ignore:
Timestamp:
Oct 18, 2008, 3:05:54 PM (16 years ago)
Author:
Frederik Heber <heber@…>
Branches:
Action_Thermostats, Add_AtomRandomPerturbation, Add_FitFragmentPartialChargesAction, Add_RotateAroundBondAction, Add_SelectAtomByNameAction, Added_ParseSaveFragmentResults, AddingActions_SaveParseParticleParameters, Adding_Graph_to_ChangeBondActions, Adding_MD_integration_tests, Adding_ParticleName_to_Atom, Adding_StructOpt_integration_tests, AtomFragments, Automaking_mpqc_open, AutomationFragmentation_failures, Candidate_v1.5.4, Candidate_v1.6.0, Candidate_v1.6.1, ChangeBugEmailaddress, ChangingTestPorts, ChemicalSpaceEvaluator, CombiningParticlePotentialParsing, Combining_Subpackages, Debian_Package_split, Debian_package_split_molecuildergui_only, Disabling_MemDebug, Docu_Python_wait, EmpiricalPotential_contain_HomologyGraph, EmpiricalPotential_contain_HomologyGraph_documentation, Enable_parallel_make_install, Enhance_userguide, Enhanced_StructuralOptimization, Enhanced_StructuralOptimization_continued, Example_ManyWaysToTranslateAtom, Exclude_Hydrogens_annealWithBondGraph, FitPartialCharges_GlobalError, Fix_BoundInBox_CenterInBox_MoleculeActions, Fix_ChargeSampling_PBC, Fix_ChronosMutex, Fix_FitPartialCharges, Fix_FitPotential_needs_atomicnumbers, Fix_ForceAnnealing, Fix_IndependentFragmentGrids, Fix_ParseParticles, Fix_ParseParticles_split_forward_backward_Actions, Fix_PopActions, Fix_QtFragmentList_sorted_selection, Fix_Restrictedkeyset_FragmentMolecule, Fix_StatusMsg, Fix_StepWorldTime_single_argument, Fix_Verbose_Codepatterns, Fix_fitting_potentials, Fixes, ForceAnnealing_goodresults, ForceAnnealing_oldresults, ForceAnnealing_tocheck, ForceAnnealing_with_BondGraph, ForceAnnealing_with_BondGraph_continued, ForceAnnealing_with_BondGraph_continued_betteresults, ForceAnnealing_with_BondGraph_contraction-expansion, FragmentAction_writes_AtomFragments, FragmentMolecule_checks_bonddegrees, GeometryObjects, Gui_Fixes, Gui_displays_atomic_force_velocity, ImplicitCharges, IndependentFragmentGrids, IndependentFragmentGrids_IndividualZeroInstances, IndependentFragmentGrids_IntegrationTest, IndependentFragmentGrids_Sole_NN_Calculation, JobMarket_RobustOnKillsSegFaults, JobMarket_StableWorkerPool, JobMarket_unresolvable_hostname_fix, MoreRobust_FragmentAutomation, ODR_violation_mpqc_open, PartialCharges_OrthogonalSummation, PdbParser_setsAtomName, PythonUI_with_named_parameters, QtGui_reactivate_TimeChanged_changes, Recreated_GuiChecks, Rewrite_FitPartialCharges, RotateToPrincipalAxisSystem_UndoRedo, SaturateAtoms_findBestMatching, SaturateAtoms_singleDegree, StoppableMakroAction, Subpackage_CodePatterns, Subpackage_JobMarket, Subpackage_LinearAlgebra, Subpackage_levmar, Subpackage_mpqc_open, Subpackage_vmg, Switchable_LogView, ThirdParty_MPQC_rebuilt_buildsystem, TrajectoryDependenant_MaxOrder, TremoloParser_IncreasedPrecision, TremoloParser_MultipleTimesteps, TremoloParser_setsAtomName, Ubuntu_1604_changes, stable
Children:
b12a35, eeec8f
Parents:
85d278
Message:

ColumnCounter -> *ColumnCounter

columns of the MatrixContainer class may now be variable over each matrix. This was changed globally to allow for varying column numbers of the new class HessianMatrix.

Location:
src
Files:
5 edited

Legend:

Unmodified
Added
Removed
  • src/analyzer.cpp

    r85d278 rf731ae  
    9393
    9494  // ---------- Parse the TE Factors into an array -----------------
    95   if (!Energy.ParseIndices()) return 1;
    96   if (Hcorrected) Hcorrection.ParseIndices();
     95  if (!Energy.InitialiseIndices()) return 1;
     96  if (Hcorrected) Hcorrection.InitialiseIndices();
    9797 
    9898  // ---------- Parse the Force indices into an array ---------------
     
    132132  output << endl << "Total Energy" << endl << "==============" << endl << Energy.Header << endl;
    133133  for(int j=0;j<Energy.RowCounter[Energy.MatrixCounter];j++) {
    134     for(int k=0;k<Energy.ColumnCounter;k++)
     134    for(int k=0;k<Energy.ColumnCounter[Energy.MatrixCounter];k++)
    135135      output << scientific << Energy.Matrix[ Energy.MatrixCounter ][j][k] << "\t";
    136136    output << endl;
     
    140140  output << endl << "Total Forces" << endl << "===============" << endl << Force.Header << endl;
    141141  for(int j=0;j<Force.RowCounter[Force.MatrixCounter];j++) {
    142     for(int k=0;k<Force.ColumnCounter;k++)
     142    for(int k=0;k<Force.ColumnCounter[Force.MatrixCounter];k++)
    143143      output << scientific << Force.Matrix[ Force.MatrixCounter ][j][k] << "\t";
    144144    output << endl;
     
    149149    output << endl << "Total Shieldings" << endl << "===============" << endl << Shielding.Header << endl;
    150150    for(int j=0;j<Shielding.RowCounter[Shielding.MatrixCounter];j++) {
    151       for(int k=0;k<Shielding.ColumnCounter;k++)
     151      for(int k=0;k<Shielding.ColumnCounter[Shielding.MatrixCounter];k++)
    152152        output << scientific << Shielding.Matrix[ Shielding.MatrixCounter ][j][k] << "\t";
    153153      output << endl;
     
    157157    output << endl << "Total Shieldings PAS" << endl << "===============" << endl << ShieldingPAS.Header << endl;
    158158    for(int j=0;j<ShieldingPAS.RowCounter[ShieldingPAS.MatrixCounter];j++) {
    159       for(int k=0;k<ShieldingPAS.ColumnCounter;k++)
     159      for(int k=0;k<ShieldingPAS.ColumnCounter[ShieldingPAS.MatrixCounter];k++)
    160160        output << scientific << ShieldingPAS.Matrix[ ShieldingPAS.MatrixCounter ][j][k] << "\t";
    161161      output << endl;
     
    166166  output << endl << "Total Times" << endl << "===============" << endl << Time.Header << endl;
    167167  for(int j=0;j<Time.RowCounter[Time.MatrixCounter];j++) {
    168     for(int k=0;k<Time.ColumnCounter;k++) {
     168    for(int k=0;k<Time.ColumnCounter[Time.MatrixCounter];k++) {
    169169      output << scientific << Time.Matrix[ Time.MatrixCounter ][j][k] << "\t";
    170170    }
     
    173173  output << endl;
    174174  output.close();
    175   for(int k=0;k<Time.ColumnCounter;k++)
     175  for(int k=0;k<Time.ColumnCounter[Time.MatrixCounter];k++)
    176176    Time.Matrix[ Time.MatrixCounter ][ Time.RowCounter[Time.MatrixCounter] ][k] = Time.Matrix[ Time.MatrixCounter ][Time.RowCounter[Time.MatrixCounter]-1][k];
    177177
     
    187187  if (!OpenOutputFile(output2, argv[3], "DeltaSimTime-Order.dat" )) return false;
    188188  for(int j=Time.RowCounter[Time.MatrixCounter];j--;)
    189     for(int k=Time.ColumnCounter;k--;) {
     189    for(int k=Time.ColumnCounter[Time.MatrixCounter];k--;) {
    190190      Time.Matrix[ Time.MatrixCounter ][j][k] = 0.;
    191191    }
     
    196196    for(int i=KeySet.FragmentsPerOrder[BondOrder];i--;)
    197197      for(int j=Time.RowCounter[Time.MatrixCounter];j--;)
    198         for(int k=Time.ColumnCounter;k--;) {
     198        for(int k=Time.ColumnCounter[Time.MatrixCounter];k--;) {
    199199          Time.Matrix[ Time.MatrixCounter ][j][k] += Time.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
    200200        }
     
    202202    output << BondOrder+1 << "\t" << counter;
    203203    output2 << BondOrder+1 << "\t" << counter;
    204     for(int k=0;k<Time.ColumnCounter;k++) {
     204    for(int k=0;k<Time.ColumnCounter[Time.MatrixCounter];k++) {
    205205      output << "\t" << scientific << Time.Matrix[ Time.MatrixCounter ][ Time.RowCounter[Time.MatrixCounter]-1 ][k];
    206206      if (fabs(Time.Matrix[ Time.MatrixCounter ][ Time.RowCounter[Time.MatrixCounter] ][k]) > MYEPSILON)
     
    224224    for(int j=0;j<ShieldingPAS.RowCounter[ShieldingPAS.MatrixCounter];j++) {
    225225      output << j << "\t";
    226       for(int k=0;k<ShieldingPAS.ColumnCounter;k++)
     226      for(int k=0;k<ShieldingPAS.ColumnCounter[ShieldingPAS.MatrixCounter];k++)
    227227        output << scientific <<  ShieldingPAS.Matrix[ ShieldingPAS.MatrixCounter ][j][k] << "\t"; //*(((k>1) && (k<6))? 1.e6 : 1.) << "\t";
    228228      output << endl;
     
    256256  for(int j=0;j<Force.RowCounter[Force.MatrixCounter];j++) {
    257257    output << j << "\t";
    258     for(int k=0;k<Force.ColumnCounter;k++)
     258    for(int k=0;k<Force.ColumnCounter[Force.MatrixCounter];k++)
    259259      output << scientific <<  Force.Matrix[ Force.MatrixCounter ][j][k] << "\t";
    260260    output << endl;
  • src/datacreator.cpp

    r85d278 rf731ae  
    6767    for(int i=KeySet.FragmentsPerOrder[BondOrder];i--;) {
    6868      for(int j=Fragments.RowCounter[ KeySet.OrderSet[BondOrder][i] ];j--;)
    69         for(int k=Fragments.ColumnCounter;k--;)
     69        for(int k=Fragments.ColumnCounter[ KeySet.OrderSet[BondOrder][i] ];k--;)
    7070          Fragments.Matrix[Fragments.MatrixCounter][j][k] += Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
    7171    }
    7272    output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
    73     for (int l=0;l<Fragments.ColumnCounter;l++)
     73    for (int l=0;l<Fragments.ColumnCounter[Fragments.MatrixCounter];l++)
    7474      output << scientific << "\t" << Fragments.Matrix[Fragments.MatrixCounter][ Fragments.RowCounter[Fragments.MatrixCounter]-1 ][l];
    7575    output << endl;
     
    101101    for(int i=KeySet.FragmentsPerOrder[BondOrder];i--;) {
    102102      for(int j=Fragments.RowCounter[ KeySet.OrderSet[BondOrder][i] ];j--;)
    103         for(int k=Fragments.ColumnCounter;k--;)
     103        for(int k=Fragments.ColumnCounter[ KeySet.OrderSet[BondOrder][i] ];k--;)
    104104          Fragments.Matrix[Fragments.MatrixCounter][j][k] -= Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
    105105    }
    106106    output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
    107     for (int l=0;l<Fragments.ColumnCounter;l++)
     107    for (int l=0;l<Fragments.ColumnCounter[Energy.MatrixCounter];l++)
    108108      if (fabs(Energy.Matrix[Energy.MatrixCounter][ Energy.RowCounter[Energy.MatrixCounter]-1 ][l]) < MYEPSILON)
    109109        output << scientific << "\t" << Fragments.Matrix[Fragments.MatrixCounter][ Fragments.RowCounter[Fragments.MatrixCounter]-1 ][l];
     
    139139    output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
    140140    CreateForce(Fragments, Fragments.MatrixCounter);
    141     for (int l=0;l<Fragments.ColumnCounter;l++)
     141    for (int l=0;l<Fragments.ColumnCounter[Fragments.MatrixCounter];l++)
    142142       output << scientific << "\t" << Fragments.Matrix[Fragments.MatrixCounter][ Fragments.RowCounter[Fragments.MatrixCounter] ][l];
    143143    output << endl;
     
    171171    output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
    172172    CreateForce(Fragments, Fragments.MatrixCounter);
    173     for (int l=0;l<Fragments.ColumnCounter;l++)
     173    for (int l=0;l<Fragments.ColumnCounter[Fragments.MatrixCounter];l++)
    174174       output << scientific << "\t" << Fragments.Matrix[Fragments.MatrixCounter][ Fragments.RowCounter[Fragments.MatrixCounter] ][l];
    175175    output << endl;
     
    207207    for(int j=0;j<Fragments.RowCounter[ Fragments.MatrixCounter ];j++) {
    208208      output << Fragments.Indices[Fragments.MatrixCounter][j] << "\t";
    209       for (int l=0;l<Fragments.ColumnCounter;l++) {
     209      for (int l=0;l<Fragments.ColumnCounter[ Fragments.MatrixCounter ];l++) {
    210210        if (((l+1) % 3) == 0) {
    211211          norm = 0.;
     
    252252    for(int j=0;j<Fragments.RowCounter[ Fragments.MatrixCounter ];j++) {
    253253      output << Fragments.Indices[Fragments.MatrixCounter][j] << "\t";
    254       for (int l=0;l<Fragments.ColumnCounter;l++)
     254      for (int l=0;l<Fragments.ColumnCounter[ Fragments.MatrixCounter ];l++)
    255255        output << scientific << Fragments.Matrix[Fragments.MatrixCounter][ j ][l] << "\t";
    256256      output << endl;
     
    278278      output << BondOrder+1 << "\t" << KeySet.OrderSet[BondOrder][i]+1;
    279279      CreateFragment(Fragment, KeySet.OrderSet[BondOrder][i]);
    280       for (int l=0;l<Fragment.ColumnCounter;l++)
     280      for (int l=0;l<Fragment.ColumnCounter[ KeySet.OrderSet[BondOrder][i] ];l++)
    281281        output << scientific << "\t" << Fragment.Matrix[ KeySet.OrderSet[BondOrder][i] ][ Fragment.RowCounter[ KeySet.OrderSet[BondOrder][i] ] ][l];
    282282      output << endl;
     
    297297    for(int i=KeySet.FragmentsPerOrder[BondOrder];i--;) {
    298298      if (fabs(Fragments.Matrix[ Fragments.MatrixCounter ][j][1]) < fabs(Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][1])) {
    299         for (int k=Fragments.ColumnCounter;k--;)
     299        for (int k=Fragments.ColumnCounter[ Fragments.MatrixCounter ];k--;)
    300300          Fragments.Matrix[ Fragments.MatrixCounter ][j][k] = Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
    301301      }
     
    314314    int i=0;
    315315    do {  // first get a minimum value unequal to 0
    316       for (int k=Fragments.ColumnCounter;k--;)
     316      for (int k=Fragments.ColumnCounter[ Fragments.MatrixCounter ];k--;)
    317317        Fragments.Matrix[ Fragments.MatrixCounter ][j][k] = Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
    318318      i++;
     
    320320    for(;i<KeySet.FragmentsPerOrder[BondOrder];i++) { // then find lowest
    321321      if (fabs(Fragments.Matrix[ Fragments.MatrixCounter ][j][1]) > fabs(Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][1])) {
    322         for (int k=Fragments.ColumnCounter;k--;)
     322        for (int k=Fragments.ColumnCounter[ Fragments.MatrixCounter ];k--;)
    323323          Fragments.Matrix[ Fragments.MatrixCounter ][j][k] = Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
    324324      }
     
    344344    CreateFragmentOrder(Fragment, KeySet, BondOrder);
    345345    output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
    346     for (int l=0;l<Fragment.ColumnCounter;l++)
     346    for (int l=0;l<Fragment.ColumnCounter[ Fragment.MatrixCounter ];l++)
    347347      output << scientific << "\t" << Fragment.Matrix[ Fragment.MatrixCounter ][ Fragment.RowCounter[ Fragment.MatrixCounter ]-1 ][l];
    348348    output << endl;
     
    358358void CreateEnergy(class MatrixContainer &Energy, int MatrixNumber)
    359359{
    360   for(int k=0;k<Energy.ColumnCounter;k++)
     360  for(int k=0;k<Energy.ColumnCounter[MatrixNumber];k++)
    361361    Energy.Matrix[MatrixNumber][ Energy.RowCounter[MatrixNumber] ] [k] =  Energy.Matrix[MatrixNumber][ Energy.RowCounter[MatrixNumber]-1 ] [k];
    362362};
     
    369369void CreateMinimumForce(class MatrixContainer &Force, int MatrixNumber)
    370370{
    371   for (int l=Force.ColumnCounter;l--;)
     371  for (int l=Force.ColumnCounter[MatrixNumber];l--;)
    372372    Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l] = 0.;
    373   for (int l=5;l<Force.ColumnCounter;l+=3) {
     373  for (int l=5;l<Force.ColumnCounter[MatrixNumber];l+=3) {
    374374    double stored = 0;
    375375    int k=0;
     
    404404{
    405405  int divisor = 0;
    406   for (int l=Force.ColumnCounter;l--;)
     406  for (int l=Force.ColumnCounter[MatrixNumber];l--;)
    407407    Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l] = 0.;
    408   for (int l=5;l<Force.ColumnCounter;l+=3) {
     408  for (int l=5;l<Force.ColumnCounter[MatrixNumber];l+=3) {
    409409    double tmp = 0;
    410410    for (int k=Force.RowCounter[MatrixNumber];k--;) {
     
    428428void CreateMaximumForce(class MatrixContainer &Force, int MatrixNumber)
    429429{
    430   for (int l=5;l<Force.ColumnCounter;l+=3) {
     430  for (int l=5;l<Force.ColumnCounter[MatrixNumber];l+=3) {
    431431    double stored = 0;
    432432    for (int k=Force.RowCounter[MatrixNumber];k--;) {
     
    460460void CreateVectorSumForce(class MatrixContainer &Force, int MatrixNumber)
    461461{
    462   for (int l=Force.ColumnCounter;l--;)
     462  for (int l=Force.ColumnCounter[MatrixNumber];l--;)
    463463    Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l] = 0.;
    464   for (int l=0;l<Force.ColumnCounter;l++) {
     464  for (int l=0;l<Force.ColumnCounter[MatrixNumber];l++) {
    465465    for (int k=Force.RowCounter[MatrixNumber];k--;)
    466466      Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l] += Force.Matrix[MatrixNumber][k][l];
     
    542542
    543543  getline(line, token, '\t');
    544   for (int i=2; i<= Energy.ColumnCounter;i++) {
     544  for (int i=2; i<= Energy.ColumnCounter[Energy.MatrixCounter];i++) {
    545545    getline(line, token, '\t');
    546546    while (token[0] == ' ') // remove leading white spaces
    547547      token.erase(0,1);
    548548    output << "'" << prefix << ".dat' title '" << token << "' using " << xargument << ":(abs($" << i+2 << ")) " << uses;
    549     if (i != (Energy.ColumnCounter))
     549    if (i != (Energy.ColumnCounter[Energy.MatrixCounter]))
    550550      output << ", \\";
    551551    output << endl;
     
    566566
    567567  getline(line, token, '\t');
    568   for (int i=1; i<= Energy.ColumnCounter;i++) {
     568  for (int i=1; i<= Energy.ColumnCounter[Energy.MatrixCounter];i++) {
    569569    getline(line, token, '\t');
    570570    while (token[0] == ' ') // remove leading white spaces
    571571      token.erase(0,1);
    572572    output << "'" << prefix << ".dat' title '" << token << "' using " << xargument << ":" << i+2 << " " << uses;
    573     if (i != (Energy.ColumnCounter))
     573    if (i != (Energy.ColumnCounter[Energy.MatrixCounter]))
    574574      output << ", \\";
    575575    output << endl;
     
    594594  getline(line, token, '\t');
    595595  getline(line, token, '\t');
    596   for (int i=7; i< Force.ColumnCounter;i+=NDIM) {
     596  for (int i=7; i< Force.ColumnCounter[Force.MatrixCounter];i+=NDIM) {
    597597    getline(line, token, '\t');
    598598    while (token[0] == ' ') // remove leading white spaces
     
    600600    token.erase(token.length(), 1);  // kill residual index char (the '0')
    601601    output << "'" << prefix << ".dat' title '" << token << "' using " << xargument << ":(sqrt($" << i+1 << "*$" << i+1 << "+$" << i+2 << "*$" << i+2 << "+$" << i+3 << "*$" << i+3 << ")) " << uses;
    602     if (i != (Force.ColumnCounter-1))
     602    if (i != (Force.ColumnCounter[Force.MatrixCounter]-1))
    603603      output << ", \\";
    604604    output << endl;
     
    625625  getline(line, token, '\t');
    626626  getline(line, token, '\t');
    627   for (int i=7; i< Force.ColumnCounter;i+=NDIM) {
     627  for (int i=7; i< Force.ColumnCounter[Force.MatrixCounter];i+=NDIM) {
    628628    getline(line, token, '\t');
    629629    while (token[0] == ' ') // remove leading white spaces
     
    631631    token.erase(token.length(), 1);  // kill residual index char (the '0')
    632632    output << "'" << prefix << ".dat' title '" << token << "' using " << xargument << ":(abs($" << i+1 << ")) " << uses;
    633     if (i != (Force.ColumnCounter-1))
     633    if (i != (Force.ColumnCounter[Force.MatrixCounter]-1))
    634634      output << ", \\";
    635635    output << endl;
     
    657657  getline(line, token, '\t');
    658658  getline(line, token, '\t');
    659   for (int i=7; i< Force.ColumnCounter;i+=NDIM) {
     659  for (int i=7; i< Force.ColumnCounter[Force.MatrixCounter];i+=NDIM) {
    660660    getline(line, token, '\t');
    661661    while (token[0] == ' ') // remove leading white spaces
     
    663663    token.erase(token.length(), 1);  // kill residual index char (the '0')
    664664    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];
    665     if (i != (Force.ColumnCounter-1))
     665    if (i != (Force.ColumnCounter[Force.MatrixCounter]-1))
    666666      output << ", \\";
    667667    output << endl;
     
    689689  getline(line, token, '\t');
    690690  getline(line, token, '\t');
    691   for (int i=7; i< Force.ColumnCounter;i+=NDIM) {
     691  for (int i=7; i< Force.ColumnCounter[Force.MatrixCounter];i+=NDIM) {
    692692    getline(line, token, '\t');
    693693    while (token[0] == ' ') // remove leading white spaces
     
    695695    token.erase(token.length(), 1);  // kill residual index char (the '0')
    696696    output << "'" << prefix << ".dat' title '" << token << "' using ($" << xargument << "+" << fixed << setprecision(1) << (double)((i-7)/3)*0.2 << "):" << i+1 << " " << uses << " " << fillcolor[(i-7)/3];
    697     if (i != (Force.ColumnCounter-1))
     697    if (i != (Force.ColumnCounter[Force.MatrixCounter]-1))
    698698      output << ", \\";
    699699    output << endl;
  • src/joiner.cpp

    r85d278 rf731ae  
    8080
    8181  // ---------- Parse the TE Factors into an array -----------------
    82   if (!Energy.ParseIndices()) return 1;
    83   if (Hcorrected) Hcorrection.ParseIndices();
     82  if (!Energy.InitialiseIndices()) return 1;
     83  if (Hcorrected) Hcorrection.InitialiseIndices();
    8484 
    8585  // ---------- Parse the Force indices into an array ---------------
     
    8787
    8888  // ---------- Parse the Hessian (=force) indices into an array ---------------
    89   if (!Hessian.ParseIndices(argv[1])) return 1;
     89  if (!Hessian.InitialiseIndices((class MatrixContainer *)&Force)) return 1;
    9090
    9191  // ---------- Parse the shielding indices into an array ---------------
  • src/parser.cpp

    r85d278 rf731ae  
    5959  Matrix = (double ***) Malloc(sizeof(double **)*(1), "MatrixContainer::MatrixContainer: ***Matrix"); // one more each for the total molecule
    6060  RowCounter = (int *) Malloc(sizeof(int)*(1), "MatrixContainer::MatrixContainer: *RowCounter");
     61  ColumnCounter = (int *) Malloc(sizeof(int)*(1), "MatrixContainer::MatrixContainer: *ColumnCounter");
    6162  Matrix[0] = NULL;
    6263  RowCounter[0] = -1;
    6364  MatrixCounter = 0;
    64   ColumnCounter = 0;
     65  ColumnCounter[0] = -1;
    6566};
    6667
     
    7071  if (Matrix != NULL) {
    7172    for(int i=MatrixCounter;i--;) {
    72       if (RowCounter != NULL) {
     73      if ((ColumnCounter != NULL) && (RowCounter != NULL)) {
    7374          for(int j=RowCounter[i]+1;j--;)
    7475            Free((void **)&Matrix[i][j], "MatrixContainer::~MatrixContainer: *Matrix[][]");
     
    7677      }
    7778    }
    78     if ((RowCounter != NULL) && (Matrix[MatrixCounter] != NULL))
     79    if ((ColumnCounter != NULL) && (RowCounter != NULL) && (Matrix[MatrixCounter] != NULL))
    7980      for(int j=RowCounter[MatrixCounter]+1;j--;)
    8081        Free((void **)&Matrix[MatrixCounter][j], "MatrixContainer::~MatrixContainer: *Matrix[MatrixCounter][]");
     
    9192  Free((void **)&Header, "MatrixContainer::~MatrixContainer: *Header");
    9293  Free((void **)&RowCounter, "MatrixContainer::~MatrixContainer: *RowCounter");
    93 };
    94 
     94  Free((void **)&ColumnCounter, "MatrixContainer::~MatrixContainer: *RowCounter");
     95};
     96
     97/** Either copies index matrix from another MatrixContainer or initialises with trivial mapping if NULL.
     98 * This either copies the index matrix or just maps 1 to 1, 2 to 2 and so on for all fragments.
     99 * \param *Matrix pointer to other MatrixContainer
     100 * \return true - copy/initialisation sucessful, false - dimension false for copying
     101 */
     102bool MatrixContainer::InitialiseIndices(class MatrixContainer *Matrix)
     103{
     104  cout << "Initialising indices";
     105  if (Matrix == NULL) {
     106    cout << " with trivial mapping." << endl;
     107    Indices = (int **) Malloc(sizeof(int *)*(MatrixCounter+1), "MatrixContainer::InitialiseIndices: **Indices");
     108    for(int i=MatrixCounter+1;i--;) {
     109      Indices[i] = (int *) Malloc(sizeof(int)*RowCounter[i], "MatrixContainer::InitialiseIndices: *Indices[]");
     110      for(int j=RowCounter[i];j--;)
     111        Indices[i][j] = j;
     112    }
     113  } else {
     114    cout << " from other MatrixContainer." << endl;
     115    if (MatrixCounter != Matrix->MatrixCounter)
     116      return false;
     117    Indices = (int **) Malloc(sizeof(int *)*(MatrixCounter+1), "MatrixContainer::InitialiseIndices: **Indices");
     118    for(int i=MatrixCounter+1;i--;) {
     119      if (RowCounter[i] != Matrix->RowCounter[i])
     120        return false;
     121      Indices[i] = (int *) Malloc(sizeof(int)*Matrix->RowCounter[i], "MatrixContainer::InitialiseIndices: *Indices[]");
     122      for(int j=Matrix->RowCounter[i];j--;)
     123        Indices[i][j] = Matrix->Indices[i][j];
     124    }
     125  }
     126  return true;
     127};
    95128
    96129/** Parsing a number of matrices.
     
    130163    line >> Header;
    131164  //cout << line.str() << endl;
    132   ColumnCounter=0;
     165  ColumnCounter[MatrixNr]=0;
    133166  while ( getline(line,token, '\t') ) {
    134167    if (token.length() > 0)
    135       ColumnCounter++;
     168      ColumnCounter[MatrixNr]++;
    136169  }
    137170  //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;
     171  cout << "ColumnCounter[" << MatrixNr << "]: " << ColumnCounter[MatrixNr] << "." << endl;
     172  if (ColumnCounter[MatrixNr] == 0)
     173    cerr << "ColumnCounter[" << MatrixNr << "]: " << ColumnCounter[MatrixNr] << " from file " << name << ", this is probably an error!" << endl;
    141174 
    142175  // scan rest for number of rows/lines
     
    150183    }
    151184  }
    152   //cout << "RowCounter[" << MatrixNr << "]: " << RowCounter[MatrixNr] << " from file " << name << "." << endl;
     185  cout << "RowCounter[" << MatrixNr << "]: " << RowCounter[MatrixNr] << " from file " << name << "." << endl;
    153186  if (RowCounter[MatrixNr] == 0)
    154187    cerr << "RowCounter[" << MatrixNr << "]: " << RowCounter[MatrixNr] << " from file " << name << ", this is probably an error!" << endl;
    155188 
    156189  // allocate matrix if it's not zero dimension in one direction
    157   if ((ColumnCounter > 0) && (RowCounter[MatrixNr] > -1)) {
     190  if ((ColumnCounter[MatrixNr] > 0) && (RowCounter[MatrixNr] > -1)) {
    158191    Matrix[MatrixNr] = (double **) Malloc(sizeof(double *)*(RowCounter[MatrixNr]+1), "MatrixContainer::ParseFragmentMatrix: **Matrix[]");
    159192 
     
    168201    strncpy(Header, line.str().c_str(), 1023); 
    169202    for(int j=0;j<RowCounter[MatrixNr];j++) {
    170       Matrix[MatrixNr][j] = (double *) Malloc(sizeof(double)*ColumnCounter, "MatrixContainer::ParseFragmentMatrix: *Matrix[][]");
     203      Matrix[MatrixNr][j] = (double *) Malloc(sizeof(double)*ColumnCounter[MatrixNr], "MatrixContainer::ParseFragmentMatrix: *Matrix[][]");
    171204      input.getline(filename, 1023);
    172205      stringstream lines(filename);
     
    174207      for(int k=skipcolumns;k--;)
    175208        lines >> filename;
    176       for(int k=0;(k<ColumnCounter) && (!lines.eof());k++) {
     209      for(int k=0;(k<ColumnCounter[MatrixNr]) && (!lines.eof());k++) {
    177210        lines >> Matrix[MatrixNr][j][k];
    178211        //cout << " " << setprecision(2) << Matrix[MatrixNr][j][k];;
    179212      }
    180213      //cout << endl;
    181       Matrix[MatrixNr][ RowCounter[MatrixNr] ] = (double *) Malloc(sizeof(double)*ColumnCounter, "MatrixContainer::ParseFragmentMatrix: *Matrix[RowCounter[MatrixNr]][]");
    182       for(int j=ColumnCounter;j--;)
     214      Matrix[MatrixNr][ RowCounter[MatrixNr] ] = (double *) Malloc(sizeof(double)*ColumnCounter[MatrixNr], "MatrixContainer::ParseFragmentMatrix: *Matrix[RowCounter[MatrixNr]][]");
     215      for(int j=ColumnCounter[MatrixNr];j--;)
    183216        Matrix[MatrixNr][ RowCounter[MatrixNr] ][j] = 0.;
    184217    }
    185218  } else {
    186     cerr << "ERROR: Matrix nr. " << MatrixNr << " has column and row count of (" << ColumnCounter << "," << RowCounter[MatrixNr] << "), could not allocate nor parse!" << endl;
     219    cerr << "ERROR: Matrix nr. " << MatrixNr << " has column and row count of (" << ColumnCounter[MatrixNr] << "," << RowCounter[MatrixNr] << "), could not allocate nor parse!" << endl;
    187220  }
    188221  input.close();
     
    235268  Matrix = (double ***) ReAlloc(Matrix, sizeof(double **)*(MatrixCounter+1), "MatrixContainer::ParseFragmentMatrix: ***Matrix"); // one more each for the total molecule
    236269  RowCounter = (int *) ReAlloc(RowCounter, sizeof(int)*(MatrixCounter+1), "MatrixContainer::ParseFragmentMatrix: *RowCounter");
     270  ColumnCounter = (int *) ReAlloc(ColumnCounter, sizeof(int)*(MatrixCounter+1), "MatrixContainer::ParseFragmentMatrix: *ColumnCounter");
    237271  for(int i=MatrixCounter+1;i--;) {
    238272    Matrix[i] = NULL;
    239273    RowCounter[i] = -1;
     274    ColumnCounter[i] = -1;
    240275  }
    241276  for(int i=0; i < MatrixCounter;i++) {
     
    255290 * \param MCounter number of matrices
    256291 * \param *RCounter number of rows for each matrix
    257  * \param CCounter number of columns for all matrices
     292 * \param *CCounter number of columns for each matrices
    258293 * \return Allocation successful
    259294 */
    260 bool MatrixContainer::AllocateMatrix(char *GivenHeader, int MCounter, int *RCounter, int CCounter)
     295bool MatrixContainer::AllocateMatrix(char *GivenHeader, int MCounter, int *RCounter, int *CCounter)
    261296{
    262297  Header = (char *) Malloc(sizeof(char)*1024, "MatrixContainer::ParseFragmentMatrix: *EnergyHeader");
     
    264299
    265300  MatrixCounter = MCounter;
    266   ColumnCounter = CCounter;
    267301  Matrix = (double ***) Malloc(sizeof(double **)*(MatrixCounter+1), "MatrixContainer::ParseFragmentMatrix: ***Matrix"); // one more each for the total molecule
    268302  RowCounter = (int *) Malloc(sizeof(int)*(MatrixCounter+1), "MatrixContainer::ParseFragmentMatrix: *RowCounter");
     303  ColumnCounter = (int *) Malloc(sizeof(int)*(MatrixCounter+1), "MatrixContainer::ParseFragmentMatrix: *ColumnCounter");
    269304  for(int i=MatrixCounter+1;i--;) {
    270305    RowCounter[i] = RCounter[i];
     306    ColumnCounter[i] = CCounter[i];
    271307    Matrix[i] = (double **) Malloc(sizeof(double *)*(RowCounter[i]+1), "MatrixContainer::ParseFragmentMatrix: **Matrix[]"); 
    272308    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--;)
     309      Matrix[i][j] = (double *) Malloc(sizeof(double)*ColumnCounter[i], "MatrixContainer::ParseFragmentMatrix: *Matrix[][]");
     310      for(int k=ColumnCounter[i];k--;)
    275311        Matrix[i][j][k] = 0.;
    276312    }
     
    286322  for(int i=MatrixCounter+1;i--;)
    287323    for(int j=RowCounter[i]+1;j--;)
    288       for(int k=ColumnCounter;k--;)
     324      for(int k=ColumnCounter[i];k--;)
    289325        Matrix[i][j][k] = 0.;
    290326   return true;
     
    299335  for(int i=MatrixCounter+1;i--;)
    300336    for(int j=RowCounter[i]+1;j--;)
    301       for(int k=ColumnCounter;k--;)
     337      for(int k=ColumnCounter[i];k--;)
    302338        if (fabs(Matrix[i][j][k]) > max)
    303339          max = fabs(Matrix[i][j][k]);
     
    315351  for(int i=MatrixCounter+1;i--;)
    316352    for(int j=RowCounter[i]+1;j--;)
    317       for(int k=ColumnCounter;k--;)
     353      for(int k=ColumnCounter[i];k--;)
    318354        if (fabs(Matrix[i][j][k]) < min)
    319355          min = fabs(Matrix[i][j][k]);
     
    331367{
    332368  for(int j=RowCounter[MatrixCounter]+1;j--;)
    333     for(int k=skipcolumns;k<ColumnCounter;k++)
     369    for(int k=skipcolumns;k<ColumnCounter[MatrixCounter];k++)
    334370      Matrix[MatrixCounter][j][k] = value;
    335371   return true;
     
    344380{
    345381  for(int j=RowCounter[MatrixCounter]+1;j--;)
    346     for(int k=skipcolumns;k<ColumnCounter;k++)
     382    for(int k=skipcolumns;k<ColumnCounter[MatrixCounter];k++)
    347383      Matrix[MatrixCounter][j][k] = values[j][k];
    348384   return true;
     
    351387/** Sums the energy with each factor and put into last element of \a ***Energies.
    352388 * Sums over "E"-terms to create the "F"-terms
    353  * \param Matrix MatrixContainer with matrices (LevelCounter by ColumnCounter) with all the energies.
     389 * \param Matrix MatrixContainer with matrices (LevelCounter by *ColumnCounter) with all the energies.
    354390 * \param KeySet KeySetContainer with bond Order and association mapping of each fragment to an order
    355391 * \param Order bond order
     
    384420              }
    385421              if (Order == SubOrder) { // equal order is always copy from Energies
    386                 for(int l=ColumnCounter;l--;) // then adds/subtract each column
     422                for(int l=ColumnCounter[ KeySet.OrderSet[SubOrder][j] ];l--;) // then adds/subtract each column
    387423                  Matrix[ KeySet.OrderSet[Order][CurrentFragment] ][m][l] += MatrixValues.Matrix[ KeySet.OrderSet[SubOrder][j] ][k][l];
    388424              } else {
    389                 for(int l=ColumnCounter;l--;)
     425                for(int l=ColumnCounter[ KeySet.OrderSet[SubOrder][j] ];l--;)
    390426                  Matrix[ KeySet.OrderSet[Order][CurrentFragment] ][m][l] -= Matrix[ KeySet.OrderSet[SubOrder][j] ][k][l];
    391427              }
    392428            }
    393             //if ((ColumnCounter>1) && (RowCounter[0]-1 >= 1))
     429            //if ((ColumnCounter[ KeySet.OrderSet[SubOrder][j] ]>1) && (RowCounter[0]-1 >= 1))
    394430             //cout << "Fragments[ KeySet.OrderSet[" << Order << "][" << CurrentFragment << "]=" << KeySet.OrderSet[Order][CurrentFragment] << " ][" << RowCounter[0]-1 << "][" << 1 << "] = " <<  Matrix[ KeySet.OrderSet[Order][CurrentFragment] ][RowCounter[0]-1][1] << endl;
    395431          }
     
    428464    output << Header << endl;
    429465    for(int j=0;j<RowCounter[i];j++) {
    430       for(int k=0;k<ColumnCounter;k++)
     466      for(int k=0;k<ColumnCounter[i];k++)
    431467        output << scientific << Matrix[i][j][k] << "\t";
    432468      output << endl;
     
    457493  output << Header << endl;
    458494  for(int j=0;j<RowCounter[MatrixCounter];j++) {
    459     for(int k=0;k<ColumnCounter;k++)
     495    for(int k=0;k<ColumnCounter[MatrixCounter];k++)
    460496      output << scientific << Matrix[MatrixCounter][j][k] << "\t";
    461497    output << endl;
     
    467503// ======================================= CLASS EnergyMatrix =============================
    468504
    469 /** Create a trivial energy index mapping.
    470  * This just maps 1 to 1, 2 to 2 and so on for all fragments.
    471  * \return creation sucessful
    472  */
    473 bool EnergyMatrix::ParseIndices()
    474 {
    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;
    483 };
    484 
    485505/** Sums the energy with each factor and put into last element of \a EnergyMatrix::Matrix.
    486506 * Sums over the "F"-terms in ANOVA decomposition.
    487  * \param Matrix MatrixContainer with matrices (LevelCounter by ColumnCounter) with all the energies.
     507 * \param Matrix MatrixContainer with matrices (LevelCounter by *ColumnCounter) with all the energies.
    488508 * \param CorrectionFragments MatrixContainer with hydrogen saturation correction per fragments
    489509 * \param KeySet KeySetContainer with bond Order and association mapping of each fragment to an order
     
    498518    for(int i=KeySet.FragmentsPerOrder[Order];i--;)
    499519      for(int j=RowCounter[ KeySet.OrderSet[Order][i] ];j--;)
    500         for(int k=ColumnCounter;k--;)
     520        for(int k=ColumnCounter[ KeySet.OrderSet[Order][i] ];k--;)
    501521          Matrix[MatrixCounter][j][k] += sign*Fragments.Matrix[ KeySet.OrderSet[Order][i] ][j][k];
    502522  else
    503523    for(int i=KeySet.FragmentsPerOrder[Order];i--;)
    504524      for(int j=RowCounter[ KeySet.OrderSet[Order][i] ];j--;)
    505         for(int k=ColumnCounter;k--;)
     525        for(int k=ColumnCounter[ KeySet.OrderSet[Order][i] ];k--;)
    506526          Matrix[MatrixCounter][j][k] += sign*(Fragments.Matrix[ KeySet.OrderSet[Order][i] ][j][k] + CorrectionFragments->Matrix[ KeySet.OrderSet[Order][i] ][j][k]);
    507527  return true;
     
    524544    // count maximum of columns
    525545    RowCounter[MatrixCounter] = 0;
    526     for(int j=0; j < MatrixCounter;j++) // (energy matrix might be bigger than number of atoms in terms of rows)
     546    ColumnCounter[MatrixCounter] = 0;
     547    for(int j=0; j < MatrixCounter;j++) { // (energy matrix might be bigger than number of atoms in terms of rows)
    527548      if (RowCounter[j] > RowCounter[MatrixCounter])
    528549        RowCounter[MatrixCounter] = RowCounter[j];
     550      if (ColumnCounter[j] > ColumnCounter[MatrixCounter])  // take maximum of all for last matrix
     551        ColumnCounter[MatrixCounter] = ColumnCounter[j];
     552    }
    529553    // allocate last plus one matrix
    530     cout << "Allocating last plus one matrix with " << (RowCounter[MatrixCounter]+1) << " rows and " << ColumnCounter << " columns." << endl;
     554    cout << "Allocating last plus one matrix with " << (RowCounter[MatrixCounter]+1) << " rows and " << ColumnCounter[MatrixCounter] << " columns." << endl;
    531555    Matrix[MatrixCounter] = (double **) Malloc(sizeof(double *)*(RowCounter[MatrixCounter]+1), "MatrixContainer::ParseFragmentMatrix: **Matrix[]");
    532556    for(int j=0;j<=RowCounter[MatrixCounter];j++)
    533       Matrix[MatrixCounter][j] = (double *) Malloc(sizeof(double)*ColumnCounter, "MatrixContainer::ParseFragmentMatrix: *Matrix[][]");
     557      Matrix[MatrixCounter][j] = (double *) Malloc(sizeof(double)*ColumnCounter[MatrixCounter], "MatrixContainer::ParseFragmentMatrix: *Matrix[][]");
    534558   
    535559    // try independently to parse global energysuffix file if present
     
    589613
    590614/** Sums the forces and puts into last element of \a ForceMatrix::Matrix.
    591  * \param Matrix MatrixContainer with matrices (LevelCounter by ColumnCounter) with all the energies.
     615 * \param Matrix MatrixContainer with matrices (LevelCounter by *ColumnCounter) with all the energies.
    592616 * \param KeySet KeySetContainer with bond Order and association mapping of each fragment to an order
    593617 * \param Order bond order
     
    609633      if (j != -1) {
    610634        //if (j == 0) cout << "Summing onto ion 0, type 0 from fragment " << FragmentNr << ", ion " << l << "." << endl;
    611         for(int k=2;k<ColumnCounter;k++)
     635        for(int k=2;k<ColumnCounter[MatrixCounter];k++)
    612636          Matrix[MatrixCounter][j][k] += sign*Fragments.Matrix[ FragmentNr ][l][k];
    613637      }
     
    655679    RowCounter[MatrixCounter]++;    // nr start at 0, count starts at 1
    656680    input.close();
     681
     682    ColumnCounter[MatrixCounter] = 0;
     683    for(int j=0; j < MatrixCounter;j++) { // (energy matrix might be bigger than number of atoms in terms of rows)
     684      if (ColumnCounter[j] > ColumnCounter[MatrixCounter])  // take maximum of all for last matrix
     685        ColumnCounter[MatrixCounter] = ColumnCounter[j];
     686    }
    657687 
    658688    // allocate last plus one matrix
    659     cout << "Allocating last plus one matrix with " << (RowCounter[MatrixCounter]+1) << " rows and " << ColumnCounter << " columns." << endl;
     689    cout << "Allocating last plus one matrix with " << (RowCounter[MatrixCounter]+1) << " rows and " << ColumnCounter[MatrixCounter] << " columns." << endl;
    660690    Matrix[MatrixCounter] = (double **) Malloc(sizeof(double *)*(RowCounter[MatrixCounter]+1), "MatrixContainer::ParseFragmentMatrix: **Matrix[]");
    661691    for(int j=0;j<=RowCounter[MatrixCounter];j++)
    662       Matrix[MatrixCounter][j] = (double *) Malloc(sizeof(double)*ColumnCounter, "MatrixContainer::ParseFragmentMatrix: *Matrix[][]");
     692      Matrix[MatrixCounter][j] = (double *) Malloc(sizeof(double)*ColumnCounter[MatrixCounter], "MatrixContainer::ParseFragmentMatrix: *Matrix[][]");
    663693
    664694    // try independently to parse global forcesuffix file if present
     
    686716  stringstream line;
    687717 
    688   cout << "Parsing force indices for " << MatrixCounter << " matrices." << endl;
     718  cout << "Parsing hessian indices for " << MatrixCounter << " matrices." << endl;
    689719  Indices = (int **) Malloc(sizeof(int *)*(MatrixCounter+1), "HessianMatrix::ParseIndices: **Indices");
    690720  line << name << FRAGMENTPREFIX << FORCESFILE;
     
    720750
    721751/** Sums the hessian entries and puts into last element of \a HessianMatrix::Matrix.
    722  * \param Matrix MatrixContainer with matrices (LevelCounter by ColumnCounter) with all the energies.
     752 * \param Matrix MatrixContainer with matrices (LevelCounter by *ColumnCounter) with all the energies.
    723753 * \param KeySet KeySetContainer with bond Order and association mapping of each fragment to an order
    724754 * \param Order bond order
     
    739769      }
    740770      if (j != -1) {
    741         //for(int m=0;m<ColumnCounter[ FragmentNr ];m++) {
    742         for(int m=0;m<ColumnCounter;m++) {
     771        for(int m=0;m<ColumnCounter[ FragmentNr ];m++) {
    743772          int k = Indices[ FragmentNr ][m];
    744 //          if (k > ColumnCounter[MatrixCounter]) {
    745 //            cerr << "Current hessian index " << k << " is greater than " << ColumnCounter[MatrixCounter] << "!" << endl;
    746           if (k > ColumnCounter) {
    747             cerr << "Current hessian index " << k << " is greater than " << ColumnCounter << "!" << endl;
     773          if (k > ColumnCounter[MatrixCounter]) {
     774            cerr << "Current hessian index " << k << " is greater than " << ColumnCounter[MatrixCounter] << "!" << endl;
    748775            return false;
    749776          }
     
    784811    }
    785812    RowCounter[MatrixCounter] = 0;
     813    ColumnCounter[MatrixCounter] = 0;
    786814    while (!input.eof()) {
    787815      input.getline(filename, 1023);
     
    790818        zeile >> nr;
    791819        //cout << "Current index: " << nr << "." << endl;
    792         if (nr > RowCounter[MatrixCounter])
     820        if (nr > RowCounter[MatrixCounter]) {
    793821          RowCounter[MatrixCounter] = nr;
     822          ColumnCounter[MatrixCounter] = nr;
     823        }
    794824      }
    795825    }
    796826    RowCounter[MatrixCounter]++;    // nr start at 0, count starts at 1
     827    ColumnCounter[MatrixCounter]++;    // nr start at 0, count starts at 1
    797828    input.close();
    798829 
    799830    // allocate last plus one matrix
    800     cout << "Allocating last plus one matrix with " << (RowCounter[MatrixCounter]+1) << " rows and " << ColumnCounter << " columns." << endl;
     831    cout << "Allocating last plus one matrix with " << (RowCounter[MatrixCounter]+1) << " rows and " << ColumnCounter[MatrixCounter] << " columns." << endl;
    801832    Matrix[MatrixCounter] = (double **) Malloc(sizeof(double *)*(RowCounter[MatrixCounter]+1), "MatrixContainer::ParseFragmentMatrix: **Matrix[]");
    802833    for(int j=0;j<=RowCounter[MatrixCounter];j++)
    803       Matrix[MatrixCounter][j] = (double *) Malloc(sizeof(double)*ColumnCounter, "MatrixContainer::ParseFragmentMatrix: *Matrix[][]");
     834      Matrix[MatrixCounter][j] = (double *) Malloc(sizeof(double)*ColumnCounter[MatrixCounter], "MatrixContainer::ParseFragmentMatrix: *Matrix[][]");
    804835
    805836    // try independently to parse global forcesuffix file if present
  • src/parser.hpp

    r85d278 rf731ae  
    4848    int MatrixCounter;
    4949    int *RowCounter;
    50     int ColumnCounter;
     50    int *ColumnCounter;
    5151 
    5252  MatrixContainer();
    5353  ~MatrixContainer();
    5454 
     55  bool InitialiseIndices(class MatrixContainer *Matrix = NULL);
    5556  bool ParseMatrix(const char *name, int skiplines, int skipcolumns, int MatrixNr);
    5657  virtual bool ParseFragmentMatrix(char *name, char *prefix, string suffix, int skiplines, int skipcolumns);
    57   bool AllocateMatrix(char *GivenHeader, int MCounter, int *RCounter, int CCounter);
     58  bool AllocateMatrix(char *GivenHeader, int MCounter, int *RCounter, int *CCounter);
    5859  bool ResetMatrix();
    5960  double FindMinValue();
     
    7273class EnergyMatrix : public MatrixContainer {
    7374  public:
    74     bool ParseIndices();
    7575    bool SumSubEnergy(class EnergyMatrix &Fragments, class EnergyMatrix *CorrectionFragments, class KeySetsContainer &KeySet, int Order, double sign);
    7676    bool ParseFragmentMatrix(char *name, char *prefix, string suffix, int skiplines, int skipcolumns);
Note: See TracChangeset for help on using the changeset viewer.