Changeset b38b64 for src/joiner.cpp

Timestamp:

Jul 23, 2009, 9:14:18 AM (16 years ago)

Author:

Frederik Heber <heber@…>

Branches:

Action_Thermostats, Add_AtomRandomPerturbation, Add_FitFragmentPartialChargesAction, Add_RotateAroundBondAction, Add_SelectAtomByNameAction, Added_ParseSaveFragmentResults, AddingActions_SaveParseParticleParameters, Adding_Graph_to_ChangeBondActions, Adding_MD_integration_tests, Adding_ParticleName_to_Atom, Adding_StructOpt_integration_tests, AtomFragments, Automaking_mpqc_open, AutomationFragmentation_failures, Candidate_v1.5.4, Candidate_v1.6.0, Candidate_v1.6.1, Candidate_v1.7.0, ChangeBugEmailaddress, ChangingTestPorts, ChemicalSpaceEvaluator, CombiningParticlePotentialParsing, Combining_Subpackages, Debian_Package_split, Debian_package_split_molecuildergui_only, Disabling_MemDebug, Docu_Python_wait, EmpiricalPotential_contain_HomologyGraph, EmpiricalPotential_contain_HomologyGraph_documentation, Enable_parallel_make_install, Enhance_userguide, Enhanced_StructuralOptimization, Enhanced_StructuralOptimization_continued, Example_ManyWaysToTranslateAtom, Exclude_Hydrogens_annealWithBondGraph, FitPartialCharges_GlobalError, Fix_BoundInBox_CenterInBox_MoleculeActions, Fix_ChargeSampling_PBC, Fix_ChronosMutex, Fix_FitPartialCharges, Fix_FitPotential_needs_atomicnumbers, Fix_ForceAnnealing, Fix_IndependentFragmentGrids, Fix_ParseParticles, Fix_ParseParticles_split_forward_backward_Actions, Fix_PopActions, Fix_QtFragmentList_sorted_selection, Fix_Restrictedkeyset_FragmentMolecule, Fix_StatusMsg, Fix_StepWorldTime_single_argument, Fix_Verbose_Codepatterns, Fix_fitting_potentials, Fixes, ForceAnnealing_goodresults, ForceAnnealing_oldresults, ForceAnnealing_tocheck, ForceAnnealing_with_BondGraph, ForceAnnealing_with_BondGraph_continued, ForceAnnealing_with_BondGraph_continued_betteresults, ForceAnnealing_with_BondGraph_contraction-expansion, FragmentAction_writes_AtomFragments, FragmentMolecule_checks_bonddegrees, GeometryObjects, Gui_Fixes, Gui_displays_atomic_force_velocity, ImplicitCharges, IndependentFragmentGrids, IndependentFragmentGrids_IndividualZeroInstances, IndependentFragmentGrids_IntegrationTest, IndependentFragmentGrids_Sole_NN_Calculation, JobMarket_RobustOnKillsSegFaults, JobMarket_StableWorkerPool, JobMarket_unresolvable_hostname_fix, MoreRobust_FragmentAutomation, ODR_violation_mpqc_open, PartialCharges_OrthogonalSummation, PdbParser_setsAtomName, PythonUI_with_named_parameters, QtGui_reactivate_TimeChanged_changes, Recreated_GuiChecks, Rewrite_FitPartialCharges, RotateToPrincipalAxisSystem_UndoRedo, SaturateAtoms_findBestMatching, SaturateAtoms_singleDegree, StoppableMakroAction, Subpackage_CodePatterns, Subpackage_JobMarket, Subpackage_LinearAlgebra, Subpackage_levmar, Subpackage_mpqc_open, Subpackage_vmg, Switchable_LogView, ThirdParty_MPQC_rebuilt_buildsystem, TrajectoryDependenant_MaxOrder, TremoloParser_IncreasedPrecision, TremoloParser_MultipleTimesteps, TremoloParser_setsAtomName, Ubuntu_1604_changes, stable

Children:

631dcb

Parents:

b77416 (diff), 72744a (diff)
Note: this is a merge changeset, the changes displayed below correspond to the merge itself.
Use the (diff) links above to see all the changes relative to each parent.

Message:

Merge branch 'HessianMatrix'

File:

• src/joiner.cpp (modified) (12 diffs)

src/joiner.cpp

@@ -19 +19 @@
   periodentafel *periode = NULL; // and a period table of all elements
   EnergyMatrix Energy;
+  EnergyMatrix EnergyFragments;
+  
   EnergyMatrix Hcorrection;
+  EnergyMatrix HcorrectionFragments;
+  
   ForceMatrix Force;
-  EnergyMatrix EnergyFragments;
-  EnergyMatrix HcorrectionFragments;
   ForceMatrix ForceFragments;
+
+  HessianMatrix Hessian;
+  HessianMatrix HessianFragments;
+  
   ForceMatrix Shielding;
   ForceMatrix ShieldingPAS;
@@ -35 +41 @@
   stringstream prefix;
   char *dir = NULL;
-  bool Hcorrected = true;
+  bool NoHCorrection = false;
+  bool NoHessian = false;
 
   cout << "Joiner" << endl;
@@ -65 +72 @@
   // ------------- Parse through all Fragment subdirs --------
   if (!Energy.ParseFragmentMatrix(argv[1], dir, EnergySuffix, 0,0)) return 1;
-  Hcorrected = Hcorrection.ParseFragmentMatrix(argv[1], "", HCORRECTIONSUFFIX, 0,0);
+  if (!Hcorrection.ParseFragmentMatrix(argv[1], "", HCORRECTIONSUFFIX, 0,0)) {
+    NoHCorrection = true;
+    cout << "No HCorrection matrices found, skipping these." << endl;
+  }
   if (!Force.ParseFragmentMatrix(argv[1], dir, ForcesSuffix, 0,0)) return 1;
+  if (!Hessian.ParseFragmentMatrix(argv[1], dir, HessianSuffix, 0,0)) {
+    NoHessian = true;
+    cout << "No hessian matrices found, skipping these." << endl;
+  }
   if (periode != NULL) { // also look for PAS values
     if (!Shielding.ParseFragmentMatrix(argv[1], dir, ShieldingSuffix, 1, 0)) return 1;
@@ -75 +89 @@
 
   // ---------- Parse the TE Factors into an array -----------------
-  if (!Energy.ParseIndices()) return 1;
-  if (Hcorrected) Hcorrection.ParseIndices();
+  if (!Energy.InitialiseIndices()) return 1;
+  if (!NoHCorrection) 
+    Hcorrection.InitialiseIndices();
   
   // ---------- Parse the Force indices into an array ---------------
   if (!Force.ParseIndices(argv[1])) return 1;
+
+  // ---------- Parse the Hessian (=force) indices into an array ---------------
+  if (!NoHessian)
+    if (!Hessian.InitialiseIndices((class MatrixContainer *)&Force)) return 1;
 
   // ---------- Parse the shielding indices into an array ---------------
@@ -94 +113 @@
 
   if (!EnergyFragments.AllocateMatrix(Energy.Header, Energy.MatrixCounter, Energy.RowCounter, Energy.ColumnCounter)) return 1;
-  if (Hcorrected)  HcorrectionFragments.AllocateMatrix(Hcorrection.Header, Hcorrection.MatrixCounter, Hcorrection.RowCounter, Hcorrection.ColumnCounter);
+  if (!NoHCorrection)  
+    HcorrectionFragments.AllocateMatrix(Hcorrection.Header, Hcorrection.MatrixCounter, Hcorrection.RowCounter, Hcorrection.ColumnCounter);
   if (!ForceFragments.AllocateMatrix(Force.Header, Force.MatrixCounter, Force.RowCounter, Force.ColumnCounter)) return 1;
+  if (!NoHessian)
+    if (!HessianFragments.AllocateMatrix(Hessian.Header, Hessian.MatrixCounter, Hessian.RowCounter, Hessian.ColumnCounter)) return 1;
   if (periode != NULL) { // also look for PAS values
     if (!ShieldingFragments.AllocateMatrix(Shielding.Header, Shielding.MatrixCounter, Shielding.RowCounter, Shielding.ColumnCounter)) return 1;
@@ -106 +128 @@
   if(!Energy.SetLastMatrix(0., 0)) return 1;
   if(!Force.SetLastMatrix(0., 2)) return 1;
+  if (!NoHessian)
+    if (!Hessian.SetLastMatrix(0., 0)) return 1;
   if (periode != NULL) { // also look for PAS values
     if(!Shielding.SetLastMatrix(0., 2)) return 1;
@@ -120 +144 @@
     cout << "Summing energy of order " << BondOrder+1 << " ..." << endl;
     if (!EnergyFragments.SumSubManyBodyTerms(Energy, KeySet, BondOrder)) return 1;
-    if (Hcorrected) { 
+    if (!NoHCorrection) { 
       HcorrectionFragments.SumSubManyBodyTerms(Hcorrection, KeySet, BondOrder);
       if (!Energy.SumSubEnergy(EnergyFragments, &HcorrectionFragments, KeySet, BondOrder, 1.)) return 1;
-      if (Hcorrected) Hcorrection.SumSubEnergy(HcorrectionFragments, NULL, KeySet, BondOrder, 1.);
+      Hcorrection.SumSubEnergy(HcorrectionFragments, NULL, KeySet, BondOrder, 1.);
     } else 
       if (!Energy.SumSubEnergy(EnergyFragments, NULL, KeySet, BondOrder, 1.)) return 1;
@@ -130 +154 @@
     if (!ForceFragments.SumSubManyBodyTerms(Force, KeySet, BondOrder)) return 1;
     if (!Force.SumSubForces(ForceFragments, KeySet, BondOrder, 1.)) return 1;
+    // --------- sum up Hessian --------------------
+    if (!NoHessian) {
+      cout << "Summing Hessian of order " << BondOrder+1 << " ..." << endl;
+      if (!HessianFragments.SumSubManyBodyTerms(Hessian, KeySet, BondOrder)) return 1;
+      if (!Hessian.SumSubHessians(HessianFragments, KeySet, BondOrder, 1.)) return 1;
+    }
     if (periode != NULL) { // also look for PAS values
       cout << "Summing shieldings and susceptibilities of order " << BondOrder+1 << " ..." << endl;
@@ -150 +180 @@
     // forces
     if (!Force.WriteLastMatrix(argv[1], (prefix.str()).c_str(), ForcesSuffix)) return 1;
+    // hessian
+    if (!NoHessian)
+      if (!Hessian.WriteLastMatrix(argv[1], (prefix.str()).c_str(), HessianSuffix)) return 1;
     // shieldings
     if (periode != NULL) { // also look for PAS values
@@ -162 +195 @@
   prefix << dir << EnergyFragmentSuffix;
   if (!EnergyFragments.WriteTotalFragments(argv[1], (prefix.str()).c_str())) return 1;
-  if (Hcorrected) {
+  if (!NoHCorrection) {
     prefix.str(" ");
     prefix << dir << HcorrectionFragmentSuffix;
@@ -171 +204 @@
   if (!ForceFragments.WriteTotalFragments(argv[1], (prefix.str()).c_str())) return 1;
   if (!CreateDataFragment(EnergyFragments, KeySet, argv[1], FRAGMENTPREFIX ENERGYPERFRAGMENT, "fragment energy versus the Fragment No", "today", CreateEnergy)) return 1;
+  if (!NoHessian) {
+    prefix.str(" ");
+    prefix << dir << HessianFragmentSuffix;
+    if (!HessianFragments.WriteTotalFragments(argv[1], (prefix.str()).c_str())) return 1;
+  }
   if (periode != NULL) { // also look for PAS values
     prefix.str(" ");
@@ -188 +226 @@
   // write last matrices as fragments into central dir (not subdir as above), for analyzer to know index bounds
   if (!Energy.WriteLastMatrix(argv[1], dir, EnergyFragmentSuffix)) return 1;
-  if (Hcorrected) Hcorrection.WriteLastMatrix(argv[1], dir, HcorrectionFragmentSuffix);
+  if (!NoHCorrection) Hcorrection.WriteLastMatrix(argv[1], dir, HcorrectionFragmentSuffix);
   if (!Force.WriteLastMatrix(argv[1], dir, ForceFragmentSuffix)) return 1;
+  if (!NoHessian)
+    if (!Hessian.WriteLastMatrix(argv[1], dir, HessianFragmentSuffix)) return 1;
   if (periode != NULL) { // also look for PAS values
     if (!Shielding.WriteLastMatrix(argv[1], dir, ShieldingFragmentSuffix)) return 1;