Changeset 6ac7ee for src/parser.cpp
- Timestamp:
- Feb 9, 2009, 5:24:10 PM (17 years ago)
- 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:
- d8b94a
- Parents:
- 124df1
- git-author:
- Frederik Heber <heber@…> (02/09/09 15:55:37)
- git-committer:
- Frederik Heber <heber@…> (02/09/09 17:24:10)
- File:
-
- 1 edited
-
src/parser.cpp (modified) (26 diffs, 1 prop)
Legend:
- Unmodified
- Added
- Removed
-
src/parser.cpp
-
Property mode
changed from
100644to100755
r124df1 r6ac7ee 2 2 * 3 3 * Declarations of various class functions for the parsing of value files. 4 * 4 * 5 5 */ 6 6 7 7 // ======================================= INCLUDES ========================================== 8 8 9 #include "helpers.hpp" 9 #include "helpers.hpp" 10 10 #include "parser.hpp" 11 11 … … 24 24 bool FilePresent(const char *filename, bool test) 25 25 { 26 27 28 29 30 31 32 33 34 35 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; 36 36 }; 37 37 … … 43 43 bool TestParams(int argc, char **argv) 44 44 { 45 46 47 48 49 45 ifstream input; 46 stringstream line; 47 48 line << argv[1] << FRAGMENTPREFIX << KEYSETFILE; 49 return FilePresent(line.str().c_str(), false); 50 50 }; 51 51 … … 55 55 */ 56 56 MatrixContainer::MatrixContainer() { 57 58 59 60 61 62 63 64 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; 65 65 }; 66 66 … … 68 68 */ 69 69 MatrixContainer::~MatrixContainer() { 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 Free((void **)&Indices, "MatrixContainer::~MatrixContainer: **Indices"); 90 91 92 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"); 93 93 }; 94 94 95 95 96 96 /** Parsing a number of matrices. 97 * 98 * 99 * 100 * 101 * 102 * 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 103 103 * \param *name directory with files 104 104 * \param skiplines number of inital lines to skip 105 105 * \param skiplines number of inital columns to skip 106 106 * \param MatrixNr index number in Matrix array to parse into 107 * \return parsing successful 107 * \return parsing successful 108 108 */ 109 109 bool MatrixContainer::ParseMatrix(const char *name, int skiplines, int skipcolumns, int MatrixNr) 110 110 { 111 112 113 114 115 116 117 //cout << "Opening " << name << " ... "<< input << endl;118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 RowCounter[MatrixNr]=-1;// counts one line too much144 while (!input.eof()) { 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 input.getline(Header, 1023);// skip header165 166 for(int k=skipcolumns;k--;)// skip columns in header too167 168 strncpy(Header, line.str().c_str(), 1023); 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 cerr << "ERROR: Matrix nr. " << MatrixNr << " has column and row count of (" << ColumnCounter << "," << RowCounter[MatrixNr] << "), could not allocate nor parse!" << endl; 187 188 189 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; 190 190 }; 191 191 192 192 /** Parsing a number of matrices. 193 193 * -# First, count the number of matrices by counting lines in KEYSETFILE 194 * -# Then, 195 * 196 * 197 * 198 * 199 * 200 * 201 * 202 * -# Finally, allocate one additional matrix (\a MatrixCounter) containing combined or temporary values 194 * -# 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 202 * -# Finally, allocate one additional matrix (\a MatrixCounter) containing combined or temporary values 203 203 * \param *name directory with files 204 204 * \param *prefix prefix of each matrix file … … 206 206 * \param skiplines number of inital lines to skip 207 207 * \param skiplines number of inital columns to skip 208 * \return parsing successful 208 * \return parsing successful 209 209 */ 210 210 bool MatrixContainer::ParseFragmentMatrix(char *name, char *prefix, string suffix, int skiplines, int skipcolumns) 211 211 { 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 input.close(); 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 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; 251 251 }; 252 252 … … 256 256 * \param *RCounter number of rows for each matrix 257 257 * \param CCounter number of columns for all matrices 258 * \return Allocation successful 258 * \return Allocation successful 259 259 */ 260 260 bool MatrixContainer::AllocateMatrix(char *GivenHeader, int MCounter, int *RCounter, int CCounter) 261 261 { 262 263 264 265 266 267 268 269 270 271 Matrix[i] = (double **) Malloc(sizeof(double *)*(RowCounter[i]+1), "MatrixContainer::ParseFragmentMatrix: **Matrix[]"); 272 273 274 275 276 277 278 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; 279 279 }; 280 280 … … 284 284 bool MatrixContainer::ResetMatrix() 285 285 { 286 287 288 289 290 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; 291 291 }; 292 292 … … 296 296 double MatrixContainer::FindMaxValue() 297 297 { 298 299 300 301 302 303 304 305 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; 306 306 return max; 307 307 }; … … 312 312 double MatrixContainer::FindMinValue() 313 313 { 314 315 316 317 318 319 320 321 322 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; 323 323 }; 324 324 … … 330 330 bool MatrixContainer::SetLastMatrix(double value, int skipcolumns) 331 331 { 332 333 334 335 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; 336 336 }; 337 337 … … 343 343 bool MatrixContainer::SetLastMatrix(double **values, int skipcolumns) 344 344 { 345 346 347 348 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; 349 349 }; 350 350 … … 358 358 bool MatrixContainer::SumSubManyBodyTerms(class MatrixContainer &MatrixValues, class KeySetsContainer &KeySet, int Order) 359 359 { 360 361 362 363 364 365 366 367 368 369 370 371 372 373 374 375 376 377 break; 378 379 380 381 382 cerr << "In fragment No. " << KeySet.OrderSet[Order][CurrentFragment]<< " current force index " << m << " is greater than " << RowCounter[ KeySet.OrderSet[Order][CurrentFragment] ] << "!" << endl;383 384 385 386 for(int l=ColumnCounter;l--;) // then adds/subtract each column 387 388 389 390 391 392 393 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 397 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 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; 405 405 }; 406 406 … … 412 412 bool MatrixContainer::WriteTotalFragments(const char *name, const char *prefix) 413 413 { 414 415 416 417 418 419 420 421 422 423 424 425 426 427 428 429 430 431 432 433 434 435 436 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; 437 437 }; 438 438 … … 445 445 bool MatrixContainer::WriteLastMatrix(const char *name, const char *prefix, const char *suffix) 446 446 { 447 448 449 450 451 452 453 454 455 456 457 458 459 460 461 462 463 464 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; 465 465 }; 466 466 … … 471 471 * \return creation sucessful 472 472 */ 473 bool EnergyMatrix::ParseIndices() 474 { 475 476 477 478 479 480 481 } 482 return true; 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 483 }; 484 484 … … 494 494 bool EnergyMatrix::SumSubEnergy(class EnergyMatrix &Fragments, class EnergyMatrix *CorrectionFragments, class KeySetsContainer &KeySet, int Order, double sign) 495 495 { 496 497 498 499 500 501 502 503 504 505 506 507 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; 508 508 }; 509 509 … … 515 515 * \param skiplines number of inital columns to skip 516 516 * \return parsing successful 517 */ 517 */ 518 518 bool EnergyMatrix::ParseFragmentMatrix(char *name, char *prefix, string suffix, int skiplines, int skipcolumns) 519 519 { 520 521 522 523 524 525 526 527 528 529 530 531 532 533 534 535 536 537 538 539 540 541 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; 542 542 }; 543 543 … … 546 546 /** Parsing force Indices of each fragment 547 547 * \param *name directory with \a ForcesFile 548 * \return parsing successful 549 */ 550 bool ForceMatrix::ParseIndices(char *name) 551 { 552 553 554 555 556 557 558 559 560 561 //cout << "Opening " << line.str() << " ... "<< input << endl;562 563 564 565 566 567 568 569 570 571 572 573 574 575 576 577 578 579 580 581 582 583 584 585 586 548 * \return parsing successful 549 */ 550 bool ForceMatrix::ParseIndices(char *name) 551 { 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; 587 587 }; 588 588 … … 592 592 * \param KeySet KeySetContainer with bond Order and association mapping of each fragment to an order 593 593 * \param Order bond order 594 * 594 * \param sign +1 or -1 595 595 * \return true if summing was successful 596 596 */ 597 597 bool ForceMatrix::SumSubForces(class ForceMatrix &Fragments, class KeySetsContainer &KeySet, int Order, double sign) 598 598 { 599 600 601 602 603 604 605 606 607 608 609 610 611 612 613 614 615 616 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; 617 617 }; 618 618 … … 625 625 * \param skiplines number of inital columns to skip 626 626 * \return parsing successful 627 */ 627 */ 628 628 bool ForceMatrix::ParseFragmentMatrix(char *name, char *prefix, string suffix, int skiplines, int skipcolumns) 629 629 { 630 631 632 633 634 635 636 637 638 639 640 641 642 643 644 645 646 647 648 649 650 651 652 653 654 655 RowCounter[MatrixCounter]++;// nr start at 0, count starts at 1656 657 658 659 660 661 662 663 664 665 666 667 668 669 670 671 672 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; 673 673 }; 674 674 … … 678 678 */ 679 679 KeySetsContainer::KeySetsContainer() { 680 681 682 683 684 685 680 KeySets = NULL; 681 AtomCounter = NULL; 682 FragmentCounter = 0; 683 Order = 0; 684 FragmentsPerOrder = 0; 685 OrderSet = NULL; 686 686 }; 687 687 … … 689 689 */ 690 690 KeySetsContainer::~KeySetsContainer() { 691 692 693 694 695 696 697 698 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"); 699 699 }; 700 700 … … 706 706 */ 707 707 bool KeySetsContainer::ParseKeySets(const char *name, const int *ACounter, const int FCounter) { 708 709 710 711 712 713 714 715 716 717 718 719 720 721 722 723 724 725 726 727 728 729 730 731 732 733 734 735 736 737 738 739 740 741 742 743 744 745 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; 746 746 }; 747 747 … … 751 751 bool KeySetsContainer::ParseManyBodyTerms() 752 752 { 753 754 755 756 757 758 759 760 761 762 763 764 765 766 767 768 769 770 771 772 773 774 775 776 777 778 779 780 781 782 783 784 785 786 787 788 789 790 791 792 793 794 795 796 797 798 799 800 801 802 803 804 805 806 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; 807 807 }; 808 808 … … 814 814 bool KeySetsContainer::Contains(const int GreaterSet, const int SmallerSet) 815 815 { 816 817 818 819 820 821 822 823 824 825 826 827 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; 828 828 }; 829 829 -
Property mode
changed from
Note:
See TracChangeset
for help on using the changeset viewer.
