ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/trunk/src/math/Quaternion.hpp
(Generate patch)

Comparing trunk/src/math/Quaternion.hpp (file contents):
Revision 93 by tim, Sun Oct 17 01:19:11 2004 UTC vs.
Revision 246 by gezelter, Wed Jan 12 22:41:40 2005 UTC

# Line 1 | Line 1
1 < /*
2 < * Copyright (C) 2000-2004  Object Oriented Parallel Simulation Engine (OOPSE) project
3 < *
4 < * Contact: oopse@oopse.org
5 < *
6 < * This program is free software; you can redistribute it and/or
7 < * modify it under the terms of the GNU Lesser General Public License
8 < * as published by the Free Software Foundation; either version 2.1
9 < * of the License, or (at your option) any later version.
10 < * All we ask is that proper credit is given for our work, which includes
11 < * - but is not limited to - adding the above copyright notice to the beginning
12 < * of your source code files, and to any copyright notice that you may distribute
13 < * with programs based on this work.
14 < *
15 < * This program is distributed in the hope that it will be useful,
16 < * but WITHOUT ANY WARRANTY; without even the implied warranty of
17 < * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
18 < * GNU Lesser General Public License for more details.
19 < *
20 < * You should have received a copy of the GNU Lesser General Public License
21 < * along with this program; if not, write to the Free Software
22 < * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
1 > /*
2 > * Copyright (c) 2005 The University of Notre Dame. All Rights Reserved.
3   *
4 + * The University of Notre Dame grants you ("Licensee") a
5 + * non-exclusive, royalty free, license to use, modify and
6 + * redistribute this software in source and binary code form, provided
7 + * that the following conditions are met:
8 + *
9 + * 1. Acknowledgement of the program authors must be made in any
10 + *    publication of scientific results based in part on use of the
11 + *    program.  An acceptable form of acknowledgement is citation of
12 + *    the article in which the program was described (Matthew
13 + *    A. Meineke, Charles F. Vardeman II, Teng Lin, Christopher
14 + *    J. Fennell and J. Daniel Gezelter, "OOPSE: An Object-Oriented
15 + *    Parallel Simulation Engine for Molecular Dynamics,"
16 + *    J. Comput. Chem. 26, pp. 252-271 (2005))
17 + *
18 + * 2. Redistributions of source code must retain the above copyright
19 + *    notice, this list of conditions and the following disclaimer.
20 + *
21 + * 3. Redistributions in binary form must reproduce the above copyright
22 + *    notice, this list of conditions and the following disclaimer in the
23 + *    documentation and/or other materials provided with the
24 + *    distribution.
25 + *
26 + * This software is provided "AS IS," without a warranty of any
27 + * kind. All express or implied conditions, representations and
28 + * warranties, including any implied warranty of merchantability,
29 + * fitness for a particular purpose or non-infringement, are hereby
30 + * excluded.  The University of Notre Dame and its licensors shall not
31 + * be liable for any damages suffered by licensee as a result of
32 + * using, modifying or distributing the software or its
33 + * derivatives. In no event will the University of Notre Dame or its
34 + * licensors be liable for any lost revenue, profit or data, or for
35 + * direct, indirect, special, consequential, incidental or punitive
36 + * damages, however caused and regardless of the theory of liability,
37 + * arising out of the use of or inability to use software, even if the
38 + * University of Notre Dame has been advised of the possibility of
39 + * such damages.
40   */
41 <
41 >
42   /**
43   * @file Quaternion.hpp
44   * @author Teng Lin
# Line 34 | Line 50
50   #define MATH_QUATERNION_HPP
51  
52   #include "math/Vector.hpp"
53 + #include "math/SquareMatrix.hpp"
54  
55   namespace oopse{
56  
# Line 48 | Line 65 | namespace oopse{
65      template<typename Real>
66      class Quaternion : public Vector<Real, 4> {
67          public:
68 <            Quaternion();
68 >            Quaternion() : Vector<Real, 4>() {}
69  
70              /** Constructs and initializes a Quaternion from w, x, y, z values */    
71              Quaternion(Real w, Real x, Real y, Real z) {
# Line 58 | Line 75 | namespace oopse{
75                  data_[3] = z;                
76              }
77              
78 <            /**
62 <             *  
63 <             */
78 >            /** Constructs and initializes a Quaternion from a  Vector<Real,4> */    
79              Quaternion(const Vector<Real,4>& v)
80                  : Vector<Real, 4>(v){
81              }
82  
83 <            /** */
83 >            /** copy assignment */
84              Quaternion& operator =(const Vector<Real, 4>& v){
85                  if (this == & v)
86                      return *this;
# Line 139 | Line 154 | namespace oopse{
154              }
155  
156              /**
157 +             * Tests if this quaternion is equal to other quaternion
158 +             * @return true if equal, otherwise return false
159 +             * @param q quaternion to be compared
160 +             */
161 +             inline bool operator ==(const Quaternion<Real>& q) {
162 +
163 +                for (unsigned int i = 0; i < 4; i ++) {
164 +                    if (!equal(data_[i], q[i])) {
165 +                        return false;
166 +                    }
167 +                }
168 +                
169 +                return true;
170 +            }
171 +            
172 +            /**
173               * Returns the inverse of this quaternion
174               * @return inverse
175               * @note since quaternion is a complex number, the inverse of quaternion
176               * q = w + xi + yj+ zk is inv_q = (w -xi - yj - zk)/(|q|^2)
177               */
178 <            Quaternion<Real> inverse(){
178 >            Quaternion<Real> inverse() {
179                  Quaternion<Real> q;
180 <                Real d = this->lengthSquared();
180 >                Real d = this->lengthSquare();
181                  
182                  q.w() = w() / d;
183                  q.x() = -x() / d;
# Line 161 | Line 192 | namespace oopse{
192               * @param q the other quaternion
193               */
194              void mul(const Quaternion<Real>& q) {
195 +                Quaternion<Real> tmp(*this);
196  
197 <                Real a0( (z() - y()) * (q.y() - q.z()) );
198 <                Real a1( (w() + x()) * (q.w() + q.x()) );
199 <                Real a2( (w() - x()) * (q.y() + q.z()) );
200 <                Real a3( (y() + z()) * (q.w() - q.x()) );
201 <                Real b0( -(x() - z()) * (q.x() - q.y()) );
170 <                Real b1( -(x() + z()) * (q.x() + q.y()) );
171 <                Real b2( (w() + y()) * (q.w() - q.z()) );
172 <                Real b3( (w() - y()) * (q.w() + q.z()) );
197 >                data_[0] = (tmp[0]*q[0]) -(tmp[1]*q[1]) - (tmp[2]*q[2]) - (tmp[3]*q[3]);
198 >                data_[1] = (tmp[0]*q[1]) + (tmp[1]*q[0]) + (tmp[2]*q[3]) - (tmp[3]*q[2]);
199 >                data_[2] = (tmp[0]*q[2]) + (tmp[2]*q[0]) + (tmp[3]*q[1]) - (tmp[1]*q[3]);
200 >                data_[3] = (tmp[0]*q[3]) + (tmp[3]*q[0]) + (tmp[1]*q[2]) - (tmp[2]*q[1]);                
201 >            }
202  
203 <                data_[0] = a0 + 0.5*(b0 + b1 + b2 + b3),;
204 <                data_[1] = a1 + 0.5*(b0 + b1 - b2 - b3);
205 <                data_[2] = a2 + 0.5*(b0 - b1 + b2 - b3),
206 <                data_[3] = a3 + 0.5*(b0 - b1 - b2 + b3) );
203 >            void mul(const Real& s) {
204 >                data_[0] *= s;
205 >                data_[1] *= s;
206 >                data_[2] *= s;
207 >                data_[3] *= s;
208              }
209  
180 –
210              /** Set the value of this quaternion to the division of itself by another quaternion */
211 <            void div(const Quaternion<Real>& q) {
211 >            void div(Quaternion<Real>& q) {
212                  mul(q.inverse());
213              }
214 +
215 +            void div(const Real& s) {
216 +                data_[0] /= s;
217 +                data_[1] /= s;
218 +                data_[2] /= s;
219 +                data_[3] /= s;
220 +            }
221              
222              Quaternion<Real>& operator *=(const Quaternion<Real>& q) {
223                  mul(q);
224                  return *this;
225              }
226 <                        
227 <            Quaternion<Real>& operator /=(const Quaternion<Real>& q) {
228 <                mul(q.inverse());
226 >
227 >            Quaternion<Real>& operator *=(const Real& s) {
228 >                mul(s);
229                  return *this;
230              }
231              
232 +            Quaternion<Real>& operator /=(Quaternion<Real>& q) {                
233 +                *this *= q.inverse();
234 +                return *this;
235 +            }
236 +
237 +            Quaternion<Real>& operator /=(const Real& s) {
238 +                div(s);
239 +                return *this;
240 +            }            
241              /**
242               * Returns the conjugate quaternion of this quaternion
243               * @return the conjugate quaternion of this quaternion
# Line 205 | Line 250 | namespace oopse{
250               * Returns the corresponding rotation matrix (3x3)
251               * @return a 3x3 rotation matrix
252               */
253 <            SquareMatrix<Real, 3, 3> toRotationMatrix3() {
254 <                SquareMatrix<Real, 3, 3> rotMat3;
253 >            SquareMatrix<Real, 3> toRotationMatrix3() {
254 >                SquareMatrix<Real, 3> rotMat3;
255  
256                  Real w2;
257                  Real x2;
# Line 232 | Line 277 | namespace oopse{
277                  rotMat3(2, 0) = 2.0 * ( x() * z() + w() * y() );
278                  rotMat3(2, 1) = 2.0 * ( y() * z() - w() * x() );
279                  rotMat3(2, 2) = w2 - x2 -y2 +z2;
280 +
281 +                return rotMat3;
282              }
283  
284      };//end Quaternion
285  
286 +
287      /**
288 +     * Returns the vaule of scalar multiplication of this quaterion q (q * s).
289 +     * @return  the vaule of scalar multiplication of this vector
290 +     * @param q the source quaternion
291 +     * @param s the scalar value
292 +     */
293 +    template<typename Real, unsigned int Dim>                
294 +    Quaternion<Real> operator * ( const Quaternion<Real>& q, Real s) {      
295 +        Quaternion<Real> result(q);
296 +        result.mul(s);
297 +        return result;          
298 +    }
299 +    
300 +    /**
301 +     * Returns the vaule of scalar multiplication of this quaterion q (q * s).
302 +     * @return  the vaule of scalar multiplication of this vector
303 +     * @param s the scalar value
304 +     * @param q the source quaternion
305 +     */  
306 +    template<typename Real, unsigned int Dim>
307 +    Quaternion<Real> operator * ( const Real& s, const Quaternion<Real>& q ) {
308 +        Quaternion<Real> result(q);
309 +        result.mul(s);
310 +        return result;          
311 +    }    
312 +
313 +    /**
314       * Returns the multiplication of two quaternion
315       * @return the multiplication of two quaternion
316       * @param q1 the first quaternion
# Line 256 | Line 330 | namespace oopse{
330       */
331  
332      template<typename Real>
333 <    inline Quaternion<Real> operator /(const Quaternion<Real>& q1, const Quaternion<Real>& q2) {
333 >    inline Quaternion<Real> operator /( Quaternion<Real>& q1,  Quaternion<Real>& q2) {
334          return q1 * q2.inverse();
335      }
336  
# Line 268 | Line 342 | namespace oopse{
342       * @note for a quaternion q, 1/q = q.inverse()
343       */
344      template<typename Real>
345 <    Quaternion<Real> operator /(const Quaternion<Real>& s, const Quaternion<Real>& q) {
345 >    Quaternion<Real> operator /(const Real& s,  Quaternion<Real>& q) {
346  
347 <        Quaternion<Real> x = q.inv();
348 <        return x * s;
347 >        Quaternion<Real> x;
348 >        x = q.inverse();
349 >        x *= s;
350 >        return x;
351      }
352 <
352 >    
353 >    template <class T>
354 >    inline bool operator==(const Quaternion<T>& lhs, const Quaternion<T>& rhs) {
355 >        return equal(lhs[0] ,rhs[0]) && equal(lhs[1] , rhs[1]) && equal(lhs[2], rhs[2]) && equal(lhs[3], rhs[3]);
356 >    }
357 >    
358      typedef Quaternion<double> Quat4d;
359   }
360   #endif //MATH_QUATERNION_HPP

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines