Changes in src/analyzer.cpp [5bc4d0:f05407]
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
src/analyzer.cpp
-
Property mode
changed from
100755
to100644
r5bc4d0 rf05407 29 29 ForceMatrix Shielding; 30 30 ForceMatrix ShieldingPAS; 31 ForceMatrix Chi;32 ForceMatrix ChiPAS;33 31 EnergyMatrix Time; 34 32 EnergyMatrix EnergyFragments; … … 37 35 ForceMatrix ShieldingFragments; 38 36 ForceMatrix ShieldingPASFragments; 39 ForceMatrix ChiFragments;40 ForceMatrix ChiPASFragments;41 37 KeySetsContainer KeySet; 42 38 ofstream output; … … 94 90 if (!Shielding.ParseFragmentMatrix(argv[1], dir, ShieldingSuffix, 1, 0)) return 1; 95 91 if (!ShieldingPAS.ParseFragmentMatrix(argv[1], dir, ShieldingPASSuffix, 1, 0)) return 1; 96 if (!Chi.ParseFragmentMatrix(argv[1], dir, ChiSuffix, 1, 0)) return 1;97 if (!ChiPAS.ParseFragmentMatrix(argv[1], dir, ChiPASSuffix, 1, 0)) return 1;98 92 } 99 93 … … 111 105 if(!Shielding.ParseIndices(argv[1])) return 1; 112 106 if(!ShieldingPAS.ParseIndices(argv[1])) return 1; 113 if (!ShieldingFragments.AllocateMatrix(Shielding.Header, Shielding.MatrixCounter, Shielding.RowCounter, Shielding.ColumnCounter)) return 1;114 if (!ShieldingPASFragments.AllocateMatrix(ShieldingPAS.Header, ShieldingPAS.MatrixCounter, ShieldingPAS.RowCounter, ShieldingPAS.ColumnCounter)) return 1;115 if(!ShieldingFragments.ParseIndices(argv[1])) return 1;116 if(!ShieldingPASFragments.ParseIndices(argv[1])) return 1;117 if(!Chi.ParseIndices(argv[1])) return 1;118 if(!ChiPAS.ParseIndices(argv[1])) return 1;119 if (!ChiFragments.AllocateMatrix(Chi.Header, Chi.MatrixCounter, Chi.RowCounter, Chi.ColumnCounter)) return 1;120 if (!ChiPASFragments.AllocateMatrix(ChiPAS.Header, ChiPAS.MatrixCounter, ChiPAS.RowCounter, ChiPAS.ColumnCounter)) return 1;121 if(!ChiFragments.ParseIndices(argv[1])) return 1;122 if(!ChiPASFragments.ParseIndices(argv[1])) return 1;123 107 } 124 108 … … 132 116 if (!ForceFragments.ParseFragmentMatrix(argv[1], dir, ForceFragmentSuffix,0,0)) return 1; 133 117 if (periode != NULL) { // also look for PAS values 134 if (!ShieldingFragments.ParseFragmentMatrix(argv[1], dir, ShieldingFragmentSuffix, 1, 0)) return 1; 135 if (!ShieldingPASFragments.ParseFragmentMatrix(argv[1], dir, ShieldingPASFragmentSuffix, 1, 0)) return 1; 136 if (!ChiFragments.ParseFragmentMatrix(argv[1], dir, ChiFragmentSuffix, 1, 0)) return 1; 137 if (!ChiPASFragments.ParseFragmentMatrix(argv[1], dir, ChiPASFragmentSuffix, 1, 0)) return 1; 118 if (!ShieldingFragments.ParseFragmentMatrix(argv[1], dir, ShieldingSuffix, 1, 0)) return 1; 119 if (!ShieldingPASFragments.ParseFragmentMatrix(argv[1], dir, ShieldingPASSuffix, 1, 0)) return 1; 138 120 } 139 121 … … 173 155 for(int k=0;k<ShieldingPAS.ColumnCounter;k++) 174 156 output << scientific << ShieldingPAS.Matrix[ ShieldingPAS.MatrixCounter ][j][k] << "\t"; 175 output << endl;176 }177 output << endl;178 179 output << endl << "Total Chis" << endl << "===============" << endl << Chi.Header[Chi.MatrixCounter] << endl;180 for(int j=0;j<Chi.RowCounter[Chi.MatrixCounter];j++) {181 for(int k=0;k<Chi.ColumnCounter;k++)182 output << scientific << Chi.Matrix[ Chi.MatrixCounter ][j][k] << "\t";183 output << endl;184 }185 output << endl;186 187 output << endl << "Total Chis PAS" << endl << "===============" << endl << ChiPAS.Header[ChiPAS.MatrixCounter] << endl;188 for(int j=0;j<ChiPAS.RowCounter[ChiPAS.MatrixCounter];j++) {189 for(int k=0;k<ChiPAS.ColumnCounter;k++)190 output << scientific << ChiPAS.Matrix[ ChiPAS.MatrixCounter ][j][k] << "\t";191 157 output << endl; 192 158 } … … 245 211 output2.close(); 246 212 247 // +++++++++++++++++++++++++++++++++++++++ Plotting shieldings 248 249 if (periode != NULL) { // also look for PAS values 250 if (!CreateDataDeltaForcesOrderPerAtom(ShieldingPAS, ShieldingPASFragments, KeySet, argv[3], "DeltaShieldingsPAS-Order", "Plot of error between approximated shieldings and full shieldings versus the Bond Order", datum)) return 1; 251 if (!CreateDataForcesOrderPerAtom(ShieldingPASFragments, KeySet, argv[3], "ShieldingsPAS-Order", "Plot of approximated shieldings versus the Bond Order", datum)) return 1; 252 if (!AppendOutputFile(output, argv[3], "ShieldingsPAS-Order.dat" )) return false; 253 output << endl << "# Full" << endl; 254 for(int j=0;j<ShieldingPAS.RowCounter[ShieldingPAS.MatrixCounter];j++) { 255 output << j << "\t"; 256 for(int k=0;k<ShieldingPAS.ColumnCounter;k++) 257 output << scientific << ShieldingPAS.Matrix[ ShieldingPAS.MatrixCounter ][j][k] << "\t"; //*(((k>1) && (k<6))? 1.e6 : 1.) << "\t"; 258 output << endl; 259 } 260 output.close(); 261 if (!CreateDataDeltaForcesOrderPerAtom(ChiPAS, ChiPASFragments, KeySet, argv[3], "DeltaChisPAS-Order", "Plot of error between approximated Chis and full Chis versus the Bond Order", datum)) return 1; 262 if (!CreateDataForcesOrderPerAtom(ChiPASFragments, KeySet, argv[3], "ChisPAS-Order", "Plot of approximated Chis versus the Bond Order", datum)) return 1; 263 if (!AppendOutputFile(output, argv[3], "ChisPAS-Order.dat" )) return false; 264 output << endl << "# Full" << endl; 265 for(int j=0;j<ChiPAS.RowCounter[ChiPAS.MatrixCounter];j++) { 266 output << j << "\t"; 267 for(int k=0;k<ChiPAS.ColumnCounter;k++) 268 output << scientific << ChiPAS.Matrix[ ChiPAS.MatrixCounter ][j][k] << "\t"; //*(((k>1) && (k<6))? 1.e6 : 1.) << "\t"; 269 output << endl; 270 } 271 output.close(); 272 } 273 274 213 275 214 // +++++++++++++++++++++++++++++++++++++++ Plotting deviation in energy to full QM 276 215 if (!CreateDataDeltaEnergyOrder(Energy, EnergyFragments, KeySet, argv[3], "DeltaEnergies-Order", "Plot of error between approximated and full energies energies versus the Bond Order", datum)) return 1; … … 293 232 // ++++++++++++++++++++++++++++++++++++++Plotting Forces vs. Order 294 233 if (!CreateDataForcesOrderPerAtom(ForceFragments, KeySet, argv[3], "Forces-Order", "Plot of approximated forces versus the Bond Order", datum)) return 1; 295 if (!AppendOutputFile(output, argv[3], "Forces-Order.dat" )) return false; 296 output << endl << "# Full" << endl; 297 for(int j=0;j<Force.RowCounter[Force.MatrixCounter];j++) { 298 output << j << "\t"; 299 for(int k=0;k<Force.ColumnCounter;k++) 300 output << scientific << Force.Matrix[ Force.MatrixCounter ][j][k] << "\t"; 301 output << endl; 302 } 303 output.close(); 234 304 235 // min force 305 236 if (!CreateDataForcesOrder(ForceFragments, KeySet, argv[3], "MinForces-Order", "Plot of min approximated forces versus the Bond Order", datum, CreateMinimumForce)) return 1; … … 428 359 if (!CreatePlotOrder(ForceFragments, KeySet, argv[3], "MaxForces-FragmentOrder", 5, "below", "y", "set boxwidth 0.2", 1, 1, "bond order", "maximum of approximated forces [Ht/a.u.]", Orderxrange.str().c_str(), yrange.str().c_str(), "1" , "with boxes fillcolor", BoxesForcePlotLine)) return 1; 429 360 430 // +++++++++++++++++++++++++++++++=Ploting approximated and true shielding for each atom431 if (periode != NULL) { // also look for PAS values432 if(!OpenOutputFile(output, argv[3], "ShieldingsPAS-Order.pyx")) return 1;433 if(!OpenOutputFile(output2, argv[3], "DeltaShieldingsPAS-Order.pyx")) return 1;434 CreatePlotHeader(output, "ShieldingsPAS-Order", 1, "top right", NULL, NULL, 1, 5, "nuclei index", "iso chemical shielding value [ppm]");435 CreatePlotHeader(output2, "DeltaShieldingsPAS-Order", 1, "top right", NULL, NULL, 1, 5, "nuclei index", "iso chemical shielding value [ppm]");436 double step=0.8/KeySet.Order;437 output << "set boxwidth " << step << endl;438 output << "plot [0:" << ShieldingPAS.RowCounter[ShieldingPAS.MatrixCounter]+10 << "]\\" << endl;439 output2 << "plot [0:" << ShieldingPAS.RowCounter[ShieldingPAS.MatrixCounter]+10 << "]\\" << endl;440 for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {441 output << "'ShieldingsPAS-Order.dat' index " << BondOrder << " title 'Order " << BondOrder+1 << "' using ($1+" << step*(double)BondOrder << "):7 with boxes, \\" << endl;442 output2 << "'DeltaShieldingsPAS-Order.dat' index " << BondOrder << " title 'Order " << BondOrder+1 << "' using ($1):7 with linespoints";443 if (BondOrder-1 != KeySet.Order)444 output2 << ", \\" << endl;445 }446 output << "'ShieldingsPAS-Order.dat' index " << KeySet.Order << " title 'Full' using ($1+" << step*(double)KeySet.Order << "):7 with boxes" << endl;447 output.close();448 output2.close();449 450 if(!OpenOutputFile(output, argv[3], "ChisPAS-Order.pyx")) return 1;451 if(!OpenOutputFile(output2, argv[3], "DeltaChisPAS-Order.pyx")) return 1;452 CreatePlotHeader(output, "ChisPAS-Order", 1, "top right", NULL, NULL, 1, 5, "nuclei index", "iso chemical Chi value [ppm]");453 CreatePlotHeader(output2, "DeltaChisPAS-Order", 1, "top right", NULL, NULL, 1, 5, "nuclei index", "iso chemical Chi value [ppm]");454 output << "set boxwidth " << step << endl;455 output << "plot [0:" << ChiPAS.RowCounter[ChiPAS.MatrixCounter]+10 << "]\\" << endl;456 output2 << "plot [0:" << ChiPAS.RowCounter[ChiPAS.MatrixCounter]+10 << "]\\" << endl;457 for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {458 output << "'ChisPAS-Order.dat' index " << BondOrder << " title 'Order " << BondOrder+1 << "' using ($1+" << step*(double)BondOrder << "):7 with boxes, \\" << endl;459 output2 << "'DeltaChisPAS-Order.dat' index " << BondOrder << " title 'Order " << BondOrder+1 << "' using ($1):7 with linespoints";460 if (BondOrder-1 != KeySet.Order)461 output2 << ", \\" << endl;462 }463 output << "'ChisPAS-Order.dat' index " << KeySet.Order << " title 'Full' using ($1+" << step*(double)KeySet.Order << "):7 with boxes" << endl;464 output.close();465 output2.close();466 }467 468 361 // create Makefile 469 362 if(!OpenOutputFile(output, argv[3], "Makefile")) return 1; -
Property mode
changed from
Note:
See TracChangeset
for help on using the changeset viewer.