source: src/Actions/MoleculeAction/RotateToPrincipalAxisSystemAction.cpp@ 9c1c89

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
Last change on this file since 9c1c89 was 83f176, checked in by Frederik Heber <heber@…>, 15 years ago

Made all member variables of class element private, added accessor functions and periodentafel is friend.

  • Property mode set to 100644
File size: 6.5 KB
Line 
1/*
2 * Project: MoleCuilder
3 * Description: creates and alters molecular systems
4 * Copyright (C) 2010 University of Bonn. All rights reserved.
5 * Please see the LICENSE file or "Copyright notice" in builder.cpp for details.
6 */
7
8/*
9 * RotateToPrincipalAxisSystemAction.cpp
10 *
11 * Created on: May 10, 2010
12 * Author: heber
13 */
14
15// include config.h
16#ifdef HAVE_CONFIG_H
17#include <config.h>
18#endif
19
20#include "Helpers/MemDebug.hpp"
21
22#include "Actions/MoleculeAction/RotateToPrincipalAxisSystemAction.hpp"
23#include "Actions/ActionRegistry.hpp"
24#include "Helpers/Log.hpp"
25#include "Helpers/Verbose.hpp"
26#include "LinearAlgebra/Line.hpp"
27#include "LinearAlgebra/Matrix.hpp"
28#include "LinearAlgebra/Vector.hpp"
29#include "element.hpp"
30#include "molecule.hpp"
31
32
33#include <iostream>
34#include <fstream>
35#include <string>
36
37using namespace std;
38
39#include "UIElements/UIFactory.hpp"
40#include "UIElements/Dialog.hpp"
41#include "Actions/ValueStorage.hpp"
42
43/****** MoleculeRotateToPrincipalAxisSystemAction *****/
44
45// memento to remember the state when undoing
46
47//class MoleculeRotateToPrincipalAxisSystemState : public ActionState {
48//public:
49// MoleculeRotateToPrincipalAxisSystemState(molecule* _mol,std::string _lastName) :
50// mol(_mol),
51// lastName(_lastName)
52// {}
53// molecule* mol;
54// std::string lastName;
55//};
56
57const char MoleculeRotateToPrincipalAxisSystemAction::NAME[] = "rotate-to-pas";
58
59MoleculeRotateToPrincipalAxisSystemAction::MoleculeRotateToPrincipalAxisSystemAction() :
60 Action(NAME)
61{}
62
63MoleculeRotateToPrincipalAxisSystemAction::~MoleculeRotateToPrincipalAxisSystemAction()
64{}
65
66void MoleculeRotateToPrincipalAxisSystem(Vector &Axis) {
67 ValueStorage::getInstance().setCurrentValue(MoleculeRotateToPrincipalAxisSystemAction::NAME, Axis);
68 ActionRegistry::getInstance().getActionByName(MoleculeRotateToPrincipalAxisSystemAction::NAME)->call(Action::NonInteractive);
69};
70
71Dialog* MoleculeRotateToPrincipalAxisSystemAction::fillDialog(Dialog *dialog) {
72 ASSERT(dialog,"No Dialog given when filling action dialog");
73
74 dialog->queryVector(NAME, false, ValueStorage::getInstance().getDescription(NAME));
75
76 return dialog;
77}
78
79Action::state_ptr MoleculeRotateToPrincipalAxisSystemAction::performCall() {
80 molecule *mol = NULL;
81 Vector Axis;
82
83 // obtain axis to rotate to
84 ValueStorage::getInstance().queryCurrentValue(NAME, Axis);
85
86 for (World::MoleculeSelectionIterator iter = World::getInstance().beginMoleculeSelection(); iter != World::getInstance().endMoleculeSelection(); ++iter) {
87 mol = iter->second;
88 DoLog(0) && (Log() << Verbose(0) << "Converting to prinicipal axis system." << endl);
89 Matrix InertiaTensor;
90 Vector *CenterOfGravity = mol->DetermineCenterOfGravity();
91
92 // reset inertia tensor
93 InertiaTensor.zero();
94
95 // sum up inertia tensor
96 for (molecule::const_iterator iter = mol->begin(); iter != mol->end(); ++iter) {
97 Vector x = (*iter)->getPosition();
98 x -= *CenterOfGravity;
99 const double mass = (*iter)->getType()->getMass();
100 InertiaTensor.at(0,0) += mass*(x[1]*x[1] + x[2]*x[2]);
101 InertiaTensor.at(0,1) += mass*(-x[0]*x[1]);
102 InertiaTensor.at(0,2) += mass*(-x[0]*x[2]);
103 InertiaTensor.at(1,0) += mass*(-x[1]*x[0]);
104 InertiaTensor.at(1,1) += mass*(x[0]*x[0] + x[2]*x[2]);
105 InertiaTensor.at(1,2) += mass*(-x[1]*x[2]);
106 InertiaTensor.at(2,0) += mass*(-x[2]*x[0]);
107 InertiaTensor.at(2,1) += mass*(-x[2]*x[1]);
108 InertiaTensor.at(2,2) += mass*(x[0]*x[0] + x[1]*x[1]);
109 }
110 // print InertiaTensor for debugging
111 DoLog(0) && (Log() << Verbose(0) << "The inertia tensor is:" << InertiaTensor << endl);
112
113 // diagonalize to determine principal axis system
114 Vector Eigenvalues = InertiaTensor.transformToEigenbasis();
115
116 for(int i=0;i<NDIM;i++)
117 DoLog(0) && (Log() << Verbose(0) << "eigenvalue = " << Eigenvalues[i] << ", eigenvector = " << InertiaTensor.column(i) << endl);
118
119 // check whether we rotate or not
120 DoLog(0) && (Log() << Verbose(0) << "Transforming molecule into PAS ... ");
121
122 // obtain first column, eigenvector to biggest eigenvalue
123 Vector BiggestEigenvector(InertiaTensor.column(Eigenvalues.SmallestComponent()));
124 Vector DesiredAxis(Axis);
125
126 // Creation Line that is the rotation axis
127 DesiredAxis.VectorProduct(BiggestEigenvector);
128 Line RotationAxis(Vector(0.,0.,0.), DesiredAxis);
129
130 // determine angle
131 const double alpha = BiggestEigenvector.Angle(Axis);
132
133 DoLog(0) && (Log() << Verbose(0) << alpha << endl);
134
135 for (molecule::iterator iter = mol->begin(); iter != mol->end(); ++iter) {
136 *(*iter) -= *CenterOfGravity;
137 (*iter)->setPosition(RotationAxis.rotateVector((*iter)->getPosition(), alpha));
138 *(*iter) += *CenterOfGravity;
139 }
140 DoLog(0) && (Log() << Verbose(0) << "done." << endl);
141
142 // summing anew for debugging (resulting matrix has to be diagonal!)
143 // reset inertia tensor
144 InertiaTensor.zero();
145
146 // sum up inertia tensor
147 for (molecule::const_iterator iter = mol->begin(); iter != mol->end(); ++iter) {
148 Vector x = (*iter)->getPosition();
149 x -= *CenterOfGravity;
150 const double mass = (*iter)->getType()->getMass();
151 InertiaTensor.at(0,0) += mass*(x[1]*x[1] + x[2]*x[2]);
152 InertiaTensor.at(0,1) += mass*(-x[0]*x[1]);
153 InertiaTensor.at(0,2) += mass*(-x[0]*x[2]);
154 InertiaTensor.at(1,0) += mass*(-x[1]*x[0]);
155 InertiaTensor.at(1,1) += mass*(x[0]*x[0] + x[2]*x[2]);
156 InertiaTensor.at(1,2) += mass*(-x[1]*x[2]);
157 InertiaTensor.at(2,0) += mass*(-x[2]*x[0]);
158 InertiaTensor.at(2,1) += mass*(-x[2]*x[1]);
159 InertiaTensor.at(2,2) += mass*(x[0]*x[0] + x[1]*x[1]);
160 // print InertiaTensor for debugging
161 DoLog(0) && (Log() << Verbose(0) << "The inertia tensor is:" << InertiaTensor << endl);
162 }
163
164 // free everything
165 delete(CenterOfGravity);
166 }
167 return Action::success;
168}
169
170Action::state_ptr MoleculeRotateToPrincipalAxisSystemAction::performUndo(Action::state_ptr _state) {
171// MoleculeRotateToPrincipalAxisSystemState *state = assert_cast<MoleculeRotateToPrincipalAxisSystemState*>(_state.get());
172
173// string newName = state->mol->getName();
174// state->mol->setName(state->lastName);
175
176 return Action::failure;
177}
178
179Action::state_ptr MoleculeRotateToPrincipalAxisSystemAction::performRedo(Action::state_ptr _state){
180 // Undo and redo have to do the same for this action
181 return performUndo(_state);
182}
183
184bool MoleculeRotateToPrincipalAxisSystemAction::canUndo() {
185 return false;
186}
187
188bool MoleculeRotateToPrincipalAxisSystemAction::shouldUndo() {
189 return false;
190}
191
192const string MoleculeRotateToPrincipalAxisSystemAction::getName() {
193 return NAME;
194}
Note: See TracBrowser for help on using the repository browser.