Changeset 042f82 for src/joiner.cpp
- Timestamp:
- Jul 23, 2009, 2:21:57 PM (16 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, 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:
- 36ec71
- Parents:
- 205ccd
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
src/joiner.cpp
r205ccd r042f82 17 17 int main(int argc, char **argv) 18 18 { 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 19 periodentafel *periode = NULL; // and a period table of all elements 20 EnergyMatrix Energy; 21 EnergyMatrix Hcorrection; 22 ForceMatrix Force; 23 EnergyMatrix EnergyFragments; 24 EnergyMatrix HcorrectionFragments; 25 ForceMatrix ForceFragments; 26 ForceMatrix Shielding; 27 ForceMatrix ShieldingPAS; 28 ForceMatrix ShieldingFragments; 29 ForceMatrix ShieldingPASFragments; 30 KeySetsContainer KeySet; 31 stringstream prefix; 32 char *dir = NULL; 33 bool Hcorrected = true; 34 34 35 36 35 cout << "Joiner" << endl; 36 cout << "======" << endl; 37 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 38 // Get the command line options 39 if (argc < 3) { 40 cout << "Usage: " << argv[0] << " <inputdir> <prefix> [elementsdb]" << endl; 41 cout << "<inputdir>\ttherein the output of a molecuilder fragmentation is expected, each fragment with a subdir containing an energy.all and a forces.all file." << endl; 42 cout << "<prefix>\tprefix of energy and forces file." << endl; 43 cout << "[elementsdb]\tpath to elements database, needed for shieldings." << endl; 44 return 1; 45 } else { 46 dir = (char *) Malloc(sizeof(char)*(strlen(argv[2])+2), "main: *dir"); 47 strcpy(dir, "/"); 48 strcat(dir, argv[2]); 49 } 50 if (argc > 3) { 51 periode = new periodentafel; 52 periode->LoadPeriodentafel(argv[3]); 53 } 54 54 55 56 57 55 // Test the given directory 56 if (!TestParams(argc, argv)) 57 return 1; 58 58 59 59 // +++++++++++++++++ PARSING +++++++++++++++++++++++++++++++ 60 60 61 62 63 64 65 66 67 68 61 // ------------- Parse through all Fragment subdirs -------- 62 if (!Energy.ParseFragmentMatrix(argv[1], dir, EnergySuffix, 0,0)) return 1; 63 Hcorrected = Hcorrection.ParseFragmentMatrix(argv[1], "", HCORRECTIONSUFFIX, 0,0); 64 if (!Force.ParseFragmentMatrix(argv[1], dir, ForcesSuffix, 0,0)) return 1; 65 if (periode != NULL) { // also look for PAS values 66 if (!Shielding.ParseFragmentMatrix(argv[1], dir, ShieldingSuffix, 1, 0)) return 1; 67 if (!ShieldingPAS.ParseFragmentMatrix(argv[1], dir, ShieldingPASSuffix, 1, 0)) return 1; 68 } 69 69 70 71 72 70 // ---------- Parse the TE Factors into an array ----------------- 71 if (!Energy.ParseIndices()) return 1; 72 if (Hcorrected) Hcorrection.ParseIndices(); 73 73 74 75 74 // ---------- Parse the Force indices into an array --------------- 75 if (!Force.ParseIndices(argv[1])) return 1; 76 76 77 78 79 80 81 77 // ---------- Parse the shielding indices into an array --------------- 78 if (periode != NULL) { // also look for PAS values 79 if(!Shielding.ParseIndices(argv[1])) return 1; 80 if(!ShieldingPAS.ParseIndices(argv[1])) return 1; 81 } 82 82 83 84 83 // ---------- Parse the KeySets into an array --------------- 84 if (!KeySet.ParseKeySets(argv[1], Force.RowCounter, Force.MatrixCounter)) return 1; 85 85 86 87 88 if (Hcorrected)HcorrectionFragments.AllocateMatrix(Hcorrection.Header, Hcorrection.MatrixCounter, Hcorrection.RowCounter, Hcorrection.ColumnCounter);89 90 91 92 93 86 if (!KeySet.ParseManyBodyTerms()) return 1; 87 if (!EnergyFragments.AllocateMatrix(Energy.Header, Energy.MatrixCounter, Energy.RowCounter, Energy.ColumnCounter)) return 1; 88 if (Hcorrected) HcorrectionFragments.AllocateMatrix(Hcorrection.Header, Hcorrection.MatrixCounter, Hcorrection.RowCounter, Hcorrection.ColumnCounter); 89 if (!ForceFragments.AllocateMatrix(Force.Header, Force.MatrixCounter, Force.RowCounter, Force.ColumnCounter)) return 1; 90 if (periode != NULL) { // also look for PAS values 91 if (!ShieldingFragments.AllocateMatrix(Shielding.Header, Shielding.MatrixCounter, Shielding.RowCounter, Shielding.ColumnCounter)) return 1; 92 if (!ShieldingPASFragments.AllocateMatrix(ShieldingPAS.Header, ShieldingPAS.MatrixCounter, ShieldingPAS.RowCounter, ShieldingPAS.ColumnCounter)) return 1; 93 } 94 94 95 96 97 98 99 100 101 95 // ----------- Resetting last matrices (where full QM values are stored right now) 96 if(!Energy.SetLastMatrix(0., 0)) return 1; 97 if(!Force.SetLastMatrix(0., 2)) return 1; 98 if (periode != NULL) { // also look for PAS values 99 if(!Shielding.SetLastMatrix(0., 2)) return 1; 100 if(!ShieldingPAS.SetLastMatrix(0., 2)) return 1; 101 } 102 102 103 103 // +++++++++++++++++ SUMMING +++++++++++++++++++++++++++++++ 104 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 105 // --------- sum up and write for each order---------------- 106 for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) { 107 // --------- sum up energy -------------------- 108 cout << "Summing energy of order " << BondOrder+1 << " ..." << endl; 109 if (!EnergyFragments.SumSubManyBodyTerms(Energy, KeySet, BondOrder)) return 1; 110 if (Hcorrected) { 111 HcorrectionFragments.SumSubManyBodyTerms(Hcorrection, KeySet, BondOrder); 112 if (!Energy.SumSubEnergy(EnergyFragments, &HcorrectionFragments, KeySet, BondOrder, 1.)) return 1; 113 if (Hcorrected) Hcorrection.SumSubEnergy(HcorrectionFragments, NULL, KeySet, BondOrder, 1.); 114 } else 115 if (!Energy.SumSubEnergy(EnergyFragments, NULL, KeySet, BondOrder, 1.)) return 1; 116 // --------- sum up Forces -------------------- 117 cout << "Summing forces of order " << BondOrder+1 << " ..." << endl; 118 if (!ForceFragments.SumSubManyBodyTerms(Force, KeySet, BondOrder)) return 1; 119 if (!Force.SumSubForces(ForceFragments, KeySet, BondOrder, 1.)) return 1; 120 if (periode != NULL) { // also look for PAS values 121 cout << "Summing shieldings of order " << BondOrder+1 << " ..." << endl; 122 if (!ShieldingFragments.SumSubManyBodyTerms(Shielding, KeySet, BondOrder)) return 1; 123 if (!Shielding.SumSubForces(ShieldingFragments, KeySet, BondOrder, 1.)) return 1; 124 if (!ShieldingPASFragments.SumSubManyBodyTerms(ShieldingPAS, KeySet, BondOrder)) return 1; 125 if (!ShieldingPAS.SumSubForces(ShieldingPASFragments, KeySet, BondOrder, 1.)) return 1; 126 } 127 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 128 // --------- write the energy and forces file -------------------- 129 prefix.str(" "); 130 prefix << dir << OrderSuffix << (BondOrder+1); 131 cout << "Writing files " << argv[1] << prefix.str() << ". ..." << endl; 132 // energy 133 if (!Energy.WriteLastMatrix(argv[1], (prefix.str()).c_str(), EnergySuffix)) return 1; 134 // forces 135 if (!Force.WriteLastMatrix(argv[1], (prefix.str()).c_str(), ForcesSuffix)) return 1; 136 // shieldings 137 if (periode != NULL) { // also look for PAS values 138 if (!Shielding.WriteLastMatrix(argv[1], (prefix.str()).c_str(), ShieldingSuffix)) return 1; 139 if (!ShieldingPAS.WriteLastMatrix(argv[1], (prefix.str()).c_str(), ShieldingPASSuffix)) return 1; 140 } 141 } 142 // fragments 143 prefix.str(" "); 144 prefix << dir << EnergyFragmentSuffix; 145 if (!EnergyFragments.WriteTotalFragments(argv[1], (prefix.str()).c_str())) return 1; 146 if (Hcorrected) { 147 prefix.str(" "); 148 prefix << dir << HcorrectionFragmentSuffix; 149 if (!HcorrectionFragments.WriteTotalFragments(argv[1], (prefix.str()).c_str())) return 1; 150 } 151 prefix.str(" "); 152 prefix << dir << ForceFragmentSuffix; 153 if (!ForceFragments.WriteTotalFragments(argv[1], (prefix.str()).c_str())) return 1; 154 if (!CreateDataFragment(EnergyFragments, KeySet, argv[1], FRAGMENTPREFIX ENERGYPERFRAGMENT, "fragment energy versus the Fragment No", "today", CreateEnergy)) return 1; 155 if (periode != NULL) { // also look for PAS values 156 prefix.str(" "); 157 prefix << dir << ShieldingFragmentSuffix; 158 if (!ShieldingFragments.WriteTotalFragments(argv[1], (prefix.str()).c_str())) return 1; 159 prefix.str(" "); 160 prefix << dir << ShieldingPASFragmentSuffix; 161 if (!ShieldingPASFragments.WriteTotalFragments(argv[1], (prefix.str()).c_str())) return 1; 162 } 163 163 164 165 166 167 168 169 170 171 164 // write last matrices as fragments into central dir (not subdir as above), for analyzer to know index bounds 165 if (!Energy.WriteLastMatrix(argv[1], dir, EnergyFragmentSuffix)) return 1; 166 if (Hcorrected) Hcorrection.WriteLastMatrix(argv[1], dir, HcorrectionFragmentSuffix); 167 if (!Force.WriteLastMatrix(argv[1], dir, ForceFragmentSuffix)) return 1; 168 if (periode != NULL) { // also look for PAS values 169 if (!Shielding.WriteLastMatrix(argv[1], dir, ShieldingFragmentSuffix)) return 1; 170 if (!ShieldingPAS.WriteLastMatrix(argv[1], dir, ShieldingPASFragmentSuffix)) return 1; 171 } 172 172 173 174 175 176 177 173 // exit 174 delete(periode); 175 Free((void **)&dir, "main: *dir"); 176 cout << "done." << endl; 177 return 0; 178 178 }; 179 179
Note:
See TracChangeset
for help on using the changeset viewer.