Changeset 8cbb97 for src/molecule.cpp


Ignore:
Timestamp:
Apr 29, 2010, 1:55:21 PM (15 years ago)
Author:
Tillmann Crueger <crueger@…>
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:
d79639
Parents:
632bc3 (diff), 753f02 (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 'VectorRefactoring' into StructureRefactoring

Conflicts:

molecuilder/src/Legacy/oldmenu.cpp
molecuilder/src/Makefile.am
molecuilder/src/analysis_correlation.cpp
molecuilder/src/boundary.cpp
molecuilder/src/builder.cpp
molecuilder/src/config.cpp
molecuilder/src/ellipsoid.cpp
molecuilder/src/linkedcell.cpp
molecuilder/src/molecule.cpp
molecuilder/src/molecule_fragmentation.cpp
molecuilder/src/molecule_geometry.cpp
molecuilder/src/molecule_graph.cpp
molecuilder/src/moleculelist.cpp
molecuilder/src/tesselation.cpp
molecuilder/src/tesselationhelpers.cpp
molecuilder/src/unittests/AnalysisCorrelationToSurfaceUnitTest.cpp
molecuilder/src/unittests/bondgraphunittest.cpp
molecuilder/src/vector.cpp
molecuilder/src/vector.hpp

File:
1 edited

Legend:

Unmodified
Added
Removed

  • src/molecule.cpp

    r632bc3 r8cbb97  
    2626#include "vector.hpp"
    2727#include "World.hpp"
     28#include "Plane.hpp"
     29#include "Exceptions/LinearDependenceException.hpp"
     30
    2831
    2932/************************************* Functions for class molecule *********************************/
    … …  
    218221//  Log() << Verbose(3) << "Begin of AddHydrogenReplacementAtom." << endl;
    219222  // create vector in direction of bond
    220   InBondvector.CopyVector(&TopReplacement->x);
    221   InBondvector.SubtractVector(&TopOrigin->x);
     223  InBondvector = TopReplacement->x - TopOrigin->x;
    222224  bondlength = InBondvector.Norm();
    223225
    … …  
    231233    Orthovector1.Zero();
    232234    for (int i=NDIM;i--;) {
    233       l = TopReplacement->x.x[i] - TopOrigin->x.x[i];
     235      l = TopReplacement->x[i] - TopOrigin->x[i];
    234236      if (fabs(l) > BondDistance) { // is component greater than bond distance
    235         Orthovector1.x[i] = (l < 0) ? -1. : +1.;
     237        Orthovector1[i] = (l < 0) ? -1. : +1.;
    236238      } // (signs are correct, was tested!)
    237239    }
    238240    matrix = ReturnFullMatrixforSymmetric(cell_size);
    239241    Orthovector1.MatrixMultiplication(matrix);
    240     InBondvector.SubtractVector(&Orthovector1); // subtract just the additional translation
     242    InBondvector -= Orthovector1; // subtract just the additional translation
    241243    Free(&matrix);
    242244    bondlength = InBondvector.Norm();
    … …  
    263265      FirstOtherAtom = World::getInstance().createAtom();    // new atom
    264266      FirstOtherAtom->type = elemente->FindElement(1);  // element is Hydrogen
    265       FirstOtherAtom->v.CopyVector(&TopReplacement->v); // copy velocity
     267      FirstOtherAtom->v = TopReplacement->v; // copy velocity
    266268      FirstOtherAtom->FixedIon = TopReplacement->FixedIon;
    267269      if (TopReplacement->type->Z == 1) { // neither rescale nor replace if it's already hydrogen
    … …  
    271273        FirstOtherAtom->father = NULL;  // if we replace hydrogen, we mark it as our father, otherwise we are just an added hydrogen with no father
    272274      }
    273       InBondvector.Scale(&BondRescale);   // rescale the distance vector to Hydrogen bond length
    274       FirstOtherAtom->x.CopyVector(&TopOrigin->x); // set coordination to origin ...
    275       FirstOtherAtom->x.AddVector(&InBondvector);  // ... and add distance vector to replacement atom
     275      InBondvector *= BondRescale;   // rescale the distance vector to Hydrogen bond length
     276      FirstOtherAtom->x = TopOrigin->x; // set coordination to origin ...
     277      FirstOtherAtom->x = InBondvector;  // ... and add distance vector to replacement atom
    276278      AllWentWell = AllWentWell && AddAtom(FirstOtherAtom);
    277279//      Log() << Verbose(4) << "Added " << *FirstOtherAtom << " at: ";
    … …  
    305307
    306308        // determine the plane of these two with the *origin
    307         AllWentWell = AllWentWell && Orthovector1.MakeNormalVector(&TopOrigin->x, &FirstOtherAtom->x, &SecondOtherAtom->x);
     309        try {
     310          Orthovector1 =Plane(TopOrigin->x, FirstOtherAtom->x, SecondOtherAtom->x).getNormal();
     311        }
     312        catch(LinearDependenceException &excp){
     313          Log() << Verbose(0) << excp;
     314          // TODO: figure out what to do with the Orthovector in this case
     315          AllWentWell = false;
     316        }
    308317      } else {
    309         Orthovector1.GetOneNormalVector(&InBondvector);
     318        Orthovector1.GetOneNormalVector(InBondvector);
    310319      }
    311320      //Log() << Verbose(3)<< "Orthovector1: ";
    … …  
    313322      //Log() << Verbose(0) << endl;
    314323      // orthogonal vector and bond vector between origin and replacement form the new plane
    315       Orthovector1.MakeNormalVector(&InBondvector);
     324      Orthovector1.MakeNormalTo(InBondvector);
    316325      Orthovector1.Normalize();
    317326      //Log() << Verbose(3) << "ReScaleCheck: " << Orthovector1.Norm() << " and " << InBondvector.Norm() << "." << endl;
    … …  
    322331      FirstOtherAtom->type = elemente->FindElement(1);
    323332      SecondOtherAtom->type = elemente->FindElement(1);
    324       FirstOtherAtom->v.CopyVector(&TopReplacement->v); // copy velocity
     333      FirstOtherAtom->v = TopReplacement->v; // copy velocity
    325334      FirstOtherAtom->FixedIon = TopReplacement->FixedIon;
    326       SecondOtherAtom->v.CopyVector(&TopReplacement->v); // copy velocity
     335      SecondOtherAtom->v = TopReplacement->v; // copy velocity
    327336      SecondOtherAtom->FixedIon = TopReplacement->FixedIon;
    328337      FirstOtherAtom->father = NULL;  // we are just an added hydrogen with no father
    … …  
    345354      SecondOtherAtom->x.Zero();
    346355      for(int i=NDIM;i--;) { // rotate by half the bond angle in both directions (InBondvector is bondangle = 0 direction)
    347         FirstOtherAtom->x.x[i] = InBondvector.x[i] * cos(bondangle) + Orthovector1.x[i] * (sin(bondangle));
    348         SecondOtherAtom->x.x[i] = InBondvector.x[i] * cos(bondangle) + Orthovector1.x[i] * (-sin(bondangle));
    349       }
    350       FirstOtherAtom->x.Scale(&BondRescale);  // rescale by correct BondDistance
    351       SecondOtherAtom->x.Scale(&BondRescale);
     356        FirstOtherAtom->x[i] = InBondvector[i] * cos(bondangle) + Orthovector1[i] * (sin(bondangle));
     357        SecondOtherAtom->x[i] = InBondvector[i] * cos(bondangle) + Orthovector1[i] * (-sin(bondangle));
     358      }
     359      FirstOtherAtom->x *= BondRescale;  // rescale by correct BondDistance
     360      SecondOtherAtom->x *= BondRescale;
    352361      //Log() << Verbose(3) << "ReScaleCheck: " << FirstOtherAtom->x.Norm() << " and " << SecondOtherAtom->x.Norm() << "." << endl;
    353362      for(int i=NDIM;i--;) { // and make relative to origin atom
    354         FirstOtherAtom->x.x[i] += TopOrigin->x.x[i];
    355         SecondOtherAtom->x.x[i] += TopOrigin->x.x[i];
     363        FirstOtherAtom->x[i] += TopOrigin->x[i];
     364        SecondOtherAtom->x[i] += TopOrigin->x[i];
    356365      }
    357366      // ... and add to molecule
    … …  
    379388      SecondOtherAtom->type = elemente->FindElement(1);
    380389      ThirdOtherAtom->type = elemente->FindElement(1);
    381       FirstOtherAtom->v.CopyVector(&TopReplacement->v); // copy velocity
     390      FirstOtherAtom->v = TopReplacement->v; // copy velocity
    382391      FirstOtherAtom->FixedIon = TopReplacement->FixedIon;
    383       SecondOtherAtom->v.CopyVector(&TopReplacement->v); // copy velocity
     392      SecondOtherAtom->v = TopReplacement->v; // copy velocity
    384393      SecondOtherAtom->FixedIon = TopReplacement->FixedIon;
    385       ThirdOtherAtom->v.CopyVector(&TopReplacement->v); // copy velocity
     394      ThirdOtherAtom->v = TopReplacement->v; // copy velocity
    386395      ThirdOtherAtom->FixedIon = TopReplacement->FixedIon;
    387396      FirstOtherAtom->father = NULL;  //  we are just an added hydrogen with no father
    … …  
    390399
    391400      // we need to vectors orthonormal the InBondvector
    392       AllWentWell = AllWentWell && Orthovector1.GetOneNormalVector(&InBondvector);
     401      AllWentWell = AllWentWell && Orthovector1.GetOneNormalVector(InBondvector);
    393402//      Log() << Verbose(3) << "Orthovector1: ";
    394403//      Orthovector1.Output(out);
    395404//      Log() << Verbose(0) << endl;
    396       AllWentWell = AllWentWell && Orthovector2.MakeNormalVector(&InBondvector, &Orthovector1);
     405      try{
     406        Orthovector2 = Plane(InBondvector, Orthovector1,0).getNormal();
     407      }
     408      catch(LinearDependenceException &excp) {
     409        Log() << Verbose(0) << excp;
     410        AllWentWell = false;
     411      }
    397412//      Log() << Verbose(3) << "Orthovector2: ";
    398413//      Orthovector2.Output(out);
    … …  
    411426      factors[1] = f;
    412427      factors[2] = 0.;
    413       FirstOtherAtom->x.LinearCombinationOfVectors(&InBondvector, &Orthovector1, &Orthovector2, factors);
     428      FirstOtherAtom->x.LinearCombinationOfVectors(InBondvector, Orthovector1, Orthovector2, factors);
    414429      factors[1] = -0.5*f;
    415430      factors[2] = g;
    416       SecondOtherAtom->x.LinearCombinationOfVectors(&InBondvector, &Orthovector1, &Orthovector2, factors);
     431      SecondOtherAtom->x.LinearCombinationOfVectors(InBondvector, Orthovector1, Orthovector2, factors);
    417432      factors[2] = -g;
    418       ThirdOtherAtom->x.LinearCombinationOfVectors(&InBondvector, &Orthovector1, &Orthovector2, factors);
     433      ThirdOtherAtom->x.LinearCombinationOfVectors(InBondvector, Orthovector1, Orthovector2, factors);
    419434
    420435      // rescale each to correct BondDistance
    … …  
    424439
    425440      // and relative to *origin atom
    426       FirstOtherAtom->x.AddVector(&TopOrigin->x);
    427       SecondOtherAtom->x.AddVector(&TopOrigin->x);
    428       ThirdOtherAtom->x.AddVector(&TopOrigin->x);
     441      FirstOtherAtom->x += TopOrigin->x;
     442      SecondOtherAtom->x += TopOrigin->x;
     443      ThirdOtherAtom->x += TopOrigin->x;
    429444
    430445      // ... and add to molecule
    … …  
    514529    }
    515530    for(j=NDIM;j--;) {
    516       Walker->x.x[j] = x[j];
    517       Walker->Trajectory.R.at(MDSteps-1).x[j] = x[j];
    518       Walker->Trajectory.U.at(MDSteps-1).x[j] = 0;
    519       Walker->Trajectory.F.at(MDSteps-1).x[j] = 0;
     531      Walker->x[j] = x[j];
     532      Walker->Trajectory.R.at(MDSteps-1)[j] = x[j];
     533      Walker->Trajectory.U.at(MDSteps-1)[j] = 0;
     534      Walker->Trajectory.F.at(MDSteps-1)[j] = 0;
    520535    }
    521536    AddAtom(Walker);  // add to molecule
    … …  
    662677{
    663678  double * const cell_size = World::getInstance().getDomain();
    664   cell_size[0] = dim->x[0];
     679  cell_size[0] = dim->at(0);
    665680  cell_size[1] = 0.;
    666   cell_size[2] = dim->x[1];
     681  cell_size[2] = dim->at(1);
    667682  cell_size[3] = 0.;
    668683  cell_size[4] = 0.;
    669   cell_size[5] = dim->x[2];
     684  cell_size[5] = dim->at(2);
    670685};
    671686
    … …  
    757772  for (int i=0;i<NDIM;i++) {
    758773    j += i+1;
    759     result = result && ((x->x[i] >= 0) && (x->x[i] < cell_size[j]));
     774    result = result && ((x->at(i) >= 0) && (x->at(i) < cell_size[j]));
    760775  }
    761776  //return result;
    … …  
    10071022    DeterminePeriodicCenter(CenterOfGravity);
    10081023    OtherMolecule->DeterminePeriodicCenter(OtherCenterOfGravity);
    1009     DoLog(5) && (Log() << Verbose(5) << "Center of Gravity: ");
    1010     CenterOfGravity.Output();
    1011     DoLog(0) && (Log() << Verbose(0) << endl << Verbose(5) << "Other Center of Gravity: ");
    1012     OtherCenterOfGravity.Output();
    1013     DoLog(0) && (Log() << Verbose(0) << endl);
    1014     if (CenterOfGravity.DistanceSquared(&OtherCenterOfGravity) > threshold*threshold) {
     1024    DoLog(5) && (Log() << Verbose(5) << "Center of Gravity: " << CenterOfGravity << endl);
     1025    DoLog(5) && (Log() << Verbose(5) << "Other Center of Gravity: " << OtherCenterOfGravity << endl);
     1026    if (CenterOfGravity.DistanceSquared(OtherCenterOfGravity) > threshold*threshold) {
    10151027      DoLog(4) && (Log() << Verbose(4) << "Centers of gravity don't match." << endl);
    10161028      result = false;
Note: See TracChangeset for help on using the changeset viewer.