Changeset 1f91f4 for src/Actions


Ignore:
Timestamp:
Feb 3, 2011, 9:51:19 AM (14 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, 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:
3f9eba
Parents:
4782599
git-author:
Frederik Heber <heber@…> (01/10/11 22:31:40)
git-committer:
Frederik Heber <heber@…> (02/03/11 09:51:19)
Message:

factored functionality in PrincipalAxisSystemAction() and RotateToPrincipalAxisSystemAction() into class molecule:

Location:
src/Actions
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • src/Actions/AnalysisAction/PrincipalAxisSystemAction.cpp

    r4782599 r1f91f4  
    4040/** =========== define the function ====================== */
    4141Action::state_ptr AnalysisPrincipalAxisSystemAction::performCall() {
    42   RealSpaceMatrix InertiaTensor;
    43 
    4442  DoLog(0) && (Log() << Verbose(0) << "Evaluating prinicipal axis." << endl);
    4543  for (World::MoleculeSelectionIterator iter = World::getInstance().beginMoleculeSelection(); iter != World::getInstance().endMoleculeSelection(); ++iter) {
    4644    molecule *mol = iter->second;
    47     Vector *CenterOfGravity = mol->DetermineCenterOfGravity();
    4845
    49     // reset inertia tensor
    50     InertiaTensor.setZero();
     46    // get inertia tensor
     47    RealSpaceMatrix InertiaTensor = mol->getInertiaTensor();
    5148
    52     // sum up inertia tensor
    53     for (molecule::const_iterator iter = mol->begin(); iter != mol->end(); ++iter) {
    54       Vector x = (*iter)->getPosition();
    55       x -= *CenterOfGravity;
    56       const double mass = (*iter)->getType()->getMass();
    57       InertiaTensor.at(0,0) += mass*(x[1]*x[1] + x[2]*x[2]);
    58       InertiaTensor.at(0,1) += mass*(-x[0]*x[1]);
    59       InertiaTensor.at(0,2) += mass*(-x[0]*x[2]);
    60       InertiaTensor.at(1,0) += mass*(-x[1]*x[0]);
    61       InertiaTensor.at(1,1) += mass*(x[0]*x[0] + x[2]*x[2]);
    62       InertiaTensor.at(1,2) += mass*(-x[1]*x[2]);
    63       InertiaTensor.at(2,0) += mass*(-x[2]*x[0]);
    64       InertiaTensor.at(2,1) += mass*(-x[2]*x[1]);
    65       InertiaTensor.at(2,2) += mass*(x[0]*x[0] + x[1]*x[1]);
    66     }
    67     // print InertiaTensor
    68     DoLog(0) && (Log() << Verbose(0) << "The inertia tensor of molecule "
    69         << mol->getName() <<  " is:"
    70         << InertiaTensor << endl);
    71 
    72     // diagonalize matrix
     49                // diagonalize matrix
    7350    Vector EigenValues = InertiaTensor.transformToEigenbasis();
    7451

  • src/Actions/MoleculeAction/RotateToPrincipalAxisSystemAction.cpp

    r4782599 r1f91f4  
    5050    mol = iter->second;
    5151    DoLog(0) && (Log() << Verbose(0) << "Converting to prinicipal axis system." << endl);
    52     RealSpaceMatrix InertiaTensor;
    53     Vector *CenterOfGravity = mol->DetermineCenterOfGravity();
    5452
    55     // reset inertia tensor
    56     InertiaTensor.setZero();
     53    RealSpaceMatrix InertiaTensor = mol->getInertiaTensor();
    5754
    58     // sum up inertia tensor
    59     for (molecule::const_iterator iter = mol->begin(); iter != mol->end(); ++iter) {
    60       Vector x = (*iter)->getPosition();
    61       x -= *CenterOfGravity;
    62       const double mass = (*iter)->getType()->getMass();
    63       InertiaTensor.at(0,0) += mass*(x[1]*x[1] + x[2]*x[2]);
    64       InertiaTensor.at(0,1) += mass*(-x[0]*x[1]);
    65       InertiaTensor.at(0,2) += mass*(-x[0]*x[2]);
    66       InertiaTensor.at(1,0) += mass*(-x[1]*x[0]);
    67       InertiaTensor.at(1,1) += mass*(x[0]*x[0] + x[2]*x[2]);
    68       InertiaTensor.at(1,2) += mass*(-x[1]*x[2]);
    69       InertiaTensor.at(2,0) += mass*(-x[2]*x[0]);
    70       InertiaTensor.at(2,1) += mass*(-x[2]*x[1]);
    71       InertiaTensor.at(2,2) += mass*(x[0]*x[0] + x[1]*x[1]);
    72     }
    73     // print InertiaTensor for debugging
    74     DoLog(0) && (Log() << Verbose(0) << "The inertia tensor is:" << InertiaTensor << endl);
    75 
    76     // diagonalize to determine principal axis system
    77     Vector Eigenvalues = InertiaTensor.transformToEigenbasis();
    78 
    79     for(int i=0;i<NDIM;i++)
    80       DoLog(0) && (Log() << Verbose(0) << "eigenvalue = " << Eigenvalues[i] << ", eigenvector = " << InertiaTensor.column(i) << endl);
    81 
    82     // check whether we rotate or not
    83     DoLog(0) && (Log() << Verbose(0) << "Transforming molecule into PAS ... ");
    84 
    85     // obtain first column, eigenvector to biggest eigenvalue
    86     Vector BiggestEigenvector(InertiaTensor.column(Eigenvalues.SmallestComponent()));
    87     Vector DesiredAxis(params.Axis);
    88 
    89     // Creation Line that is the rotation axis
    90     DesiredAxis.VectorProduct(BiggestEigenvector);
    91     Line RotationAxis(Vector(0.,0.,0.), DesiredAxis);
    92 
    93     // determine angle
    94     const double alpha = BiggestEigenvector.Angle(params.Axis);
    95 
    96     DoLog(0) && (Log() << Verbose(0) << alpha << endl);
    97 
    98     for (molecule::iterator iter = mol->begin(); iter != mol->end(); ++iter) {
    99       *(*iter) -= *CenterOfGravity;
    100       (*iter)->setPosition(RotationAxis.rotateVector((*iter)->getPosition(), alpha));
    101       *(*iter) += *CenterOfGravity;
    102     }
    103     DoLog(0) && (Log() << Verbose(0) << "done." << endl);
     55    mol->RotateToPrincipalAxisSystem(params.Axis);
    10456
    10557    // summing anew for debugging (resulting matrix has to be diagonal!)
    106     // reset inertia tensor
    107     InertiaTensor.setZero();
    108 
    109     // sum up inertia tensor
    110     for (molecule::const_iterator iter = mol->begin(); iter != mol->end(); ++iter) {
    111       Vector x = (*iter)->getPosition();
    112       x -= *CenterOfGravity;
    113       const double mass = (*iter)->getType()->getMass();
    114       InertiaTensor.at(0,0) += mass*(x[1]*x[1] + x[2]*x[2]);
    115       InertiaTensor.at(0,1) += mass*(-x[0]*x[1]);
    116       InertiaTensor.at(0,2) += mass*(-x[0]*x[2]);
    117       InertiaTensor.at(1,0) += mass*(-x[1]*x[0]);
    118       InertiaTensor.at(1,1) += mass*(x[0]*x[0] + x[2]*x[2]);
    119       InertiaTensor.at(1,2) += mass*(-x[1]*x[2]);
    120       InertiaTensor.at(2,0) += mass*(-x[2]*x[0]);
    121       InertiaTensor.at(2,1) += mass*(-x[2]*x[1]);
    122       InertiaTensor.at(2,2) += mass*(x[0]*x[0] + x[1]*x[1]);
    123       // print InertiaTensor for debugging
    124       DoLog(0) && (Log() << Verbose(0) << "The inertia tensor is:" << InertiaTensor << endl);
    125     }
    126 
    127     // free everything
    128     delete(CenterOfGravity);
     58    InertiaTensor = mol->getInertiaTensor();
    12959  }
    13060  return Action::success;
Note: See TracChangeset for help on using the changeset viewer.