FreeFOAM The Cross-Platform CFD Toolkit
linearAxialAngularSpring.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
29 #include <OpenFOAM/transform.H>
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace sixDoFRigidBodyMotionRestraints
36 {
39  (
43  );
44 };
45 };
46 
47 
48 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
49 
52 (
53  const dictionary& sDoFRBMRDict
54 )
55 :
56  sixDoFRigidBodyMotionRestraint(sDoFRBMRDict),
57  refQ_(),
58  axis_(),
59  stiffness_(),
60  damping_()
61 {
62  read(sDoFRBMRDict);
63 }
64 
65 
66 // * * * * * * * * * * * * * * * * Destructors * * * * * * * * * * * * * * * //
67 
70 {}
71 
72 
73 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
74 
75 void
77 (
78  const sixDoFRigidBodyMotion& motion,
79  vector& restraintPosition,
80  vector& restraintForce,
81  vector& restraintMoment
82 ) const
83 {
84  vector refDir = rotationTensor(vector(1, 0 ,0), axis_) & vector(0, 1, 0);
85 
86  vector oldDir = refQ_ & refDir;
87 
88  vector newDir = motion.orientation() & refDir;
89 
90  if (mag(oldDir & axis_) > 0.95 || mag(newDir & axis_) > 0.95)
91  {
92  // Directions getting close to the axis, change reference
93 
94  refDir = rotationTensor(vector(1, 0 ,0), axis_) & vector(0, 0, 1);
95 
96  vector oldDir = refQ_ & refDir;
97 
98  vector newDir = motion.orientation() & refDir;
99  }
100 
101  // Removing any axis component from oldDir and newDir and normalising
102  oldDir -= (axis_ & oldDir)*axis_;
103  oldDir /= (mag(oldDir) + VSMALL);
104 
105  newDir -= (axis_ & newDir)*axis_;
106  newDir /= (mag(newDir) + VSMALL);
107 
108  scalar theta = mag(acos(min(oldDir & newDir, 1.0)));
109 
110  // Temporary axis with sign information.
111  vector a = (oldDir ^ newDir);
112 
113  // Remove any component that is not along axis that may creep in
114  a = (a & axis_)*axis_;
115 
116  scalar magA = mag(a);
117 
118  if (magA > VSMALL)
119  {
120  a /= magA;
121  }
122  else
123  {
124  a = vector::zero;
125  }
126 
127  // Damping of along axis angular velocity only
128  restraintMoment = -stiffness_*theta*a - damping_*(motion.omega() & a)*a;
129 
130  restraintForce = vector::zero;
131 
132  // Not needed to be altered as restraintForce is zero, but set to
133  // centreOfMass to be sure of no spurious moment
134  restraintPosition = motion.centreOfMass();
135 
136  if (motion.report())
137  {
138  Info<< " angle " << theta*sign(a & axis_)
139  << " force " << restraintForce
140  << " moment " << restraintMoment
141  << endl;
142  }
143 }
144 
145 
147 (
148  const dictionary& sDoFRBMRDict
149 )
150 {
152 
153  refQ_ = sDoFRBMRCoeffs_.lookupOrDefault<tensor>("referenceOrientation", I);
154 
155  if (mag(mag(refQ_) - sqrt(3.0)) > 1e-9)
156  {
158  (
159  "Foam::sixDoFRigidBodyMotionRestraints::linearAxialAngularSpring::"
160  "read"
161  "("
162  "const dictionary& sDoFRBMRDict"
163  ")"
164  )
165  << "referenceOrientation " << refQ_ << " is not a rotation tensor. "
166  << "mag(referenceOrientation) - sqrt(3) = "
167  << mag(refQ_) - sqrt(3.0) << nl
168  << exit(FatalError);
169  }
170 
171  axis_ = sDoFRBMRCoeffs_.lookup("axis");
172 
173  scalar magAxis(mag(axis_));
174 
175  if (magAxis > VSMALL)
176  {
177  axis_ /= magAxis;
178  }
179  else
180  {
182  (
183  "Foam::sixDoFRigidBodyMotionRestraints::linearAxialAngularSpring::"
184  "read"
185  "("
186  "const dictionary& sDoFRBMCDict"
187  ")"
188  )
189  << "axis has zero length"
190  << abort(FatalError);
191  }
192 
193  sDoFRBMRCoeffs_.lookup("stiffness") >> stiffness_;
194 
195  sDoFRBMRCoeffs_.lookup("damping") >> damping_;
196 
197  return true;
198 }
199 
200 
202 (
203  Ostream& os
204 ) const
205 {
206  os.writeKeyword("referenceOrientation")
207  << refQ_ << token::END_STATEMENT << nl;
208 
209  os.writeKeyword("axis")
210  << axis_ << token::END_STATEMENT << nl;
211 
212  os.writeKeyword("stiffness")
213  << stiffness_ << token::END_STATEMENT << nl;
214 
215  os.writeKeyword("damping")
216  << damping_ << token::END_STATEMENT << nl;
217 }
218 
219 // ************************ vim: set sw=4 sts=4 et: ************************ //