OpenStructure
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
vecmat3_op.hh
Go to the documentation of this file.
1 //------------------------------------------------------------------------------
2 // This file is part of the OpenStructure project <www.openstructure.org>
3 //
4 // Copyright (C) 2008-2011 by the OpenStructure authors
5 //
6 // This library is free software; you can redistribute it and/or modify it under
7 // the terms of the GNU Lesser General Public License as published by the Free
8 // Software Foundation; either version 3.0 of the License, or (at your option)
9 // any later version.
10 // This library is distributed in the hope that it will be useful, but WITHOUT
11 // ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
12 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
13 // details.
14 //
15 // You should have received a copy of the GNU Lesser General Public License
16 // along with this library; if not, write to the Free Software Foundation, Inc.,
17 // 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
18 //------------------------------------------------------------------------------
19 #ifndef GEOM_VECMAT3_OP_HH
20 #define GEOM_VECMAT3_OP_HH
21 
22 #include <ostream>
23 #include "constants.hh"
24 
26 #include <ost/geom/vec3.hh>
27 #include <ost/geom/mat3.hh>
28 
29 namespace geom {
30 
31 
33 inline Real Length2(const Vec3& v)
34 {
35  return v[0]*v[0]+v[1]*v[1]+v[2]*v[2];
36 }
37 
39 inline Real Length(const Vec3& v)
40 {
41  return std::sqrt(Length2(v));
42 }
43 
45 inline bool Equal(const Vec3& v1, const Vec3& v2,
46  Real ephilon=EPSILON)
47 {
48  return std::fabs(v1[0]-v2[0])<ephilon &&
49  std::fabs(v1[1]-v2[1])<ephilon &&
50  std::fabs(v1[2]-v2[2])<ephilon;
51 }
52 
54 inline bool Equal(const Mat3& m1, const Mat3& m2,
55  Real ephilon=EPSILON)
56 {
57  return std::fabs(m1(0,0)-m2(0,0))<ephilon &&
58  std::fabs(m1(0,1)-m2(0,1))<ephilon &&
59  std::fabs(m1(0,2)-m2(0,2))<ephilon &&
60  std::fabs(m1(1,0)-m2(1,0))<ephilon &&
61  std::fabs(m1(1,1)-m2(1,1))<ephilon &&
62  std::fabs(m1(1,2)-m2(1,2))<ephilon &&
63  std::fabs(m1(2,0)-m2(2,0))<ephilon &&
64  std::fabs(m1(2,1)-m2(2,1))<ephilon &&
65  std::fabs(m1(2,2)-m2(2,2))<ephilon;
66 }
67 
69 inline Real Dot(const Vec3& v1, const Vec3& v2)
70 {
71  return v1[0]*v2[0]+v1[1]*v2[1]+v1[2]*v2[2];
72 }
73 
75 inline Vec3 Normalize(const Vec3& v)
76 {
77  Real l=Length(v);
78  if(l==0.0) {
79  return v;
80  }
81  return v/l;
82 }
83 
85 inline Vec3 Cross(const Vec3& v1, const Vec3& v2)
86 {
87  Vec3 nrvo(v1[1]*v2[2]-v2[1]*v1[2],
88  v1[2]*v2[0]-v2[2]*v1[0],
89  v1[0]*v2[1]-v2[0]*v1[1]);
90  return nrvo;
91 }
92 
94 inline Vec3 CompMultiply(const Vec3& v1, const Vec3& v2)
95 {
96  Vec3 nrvo(v1[0]*v2[0],v1[1]*v2[1],v1[2]*v2[2]);
97  return nrvo;
98 }
99 
101 inline Vec3 CompDivide(const Vec3& v1, const Vec3& v2)
102 {
103  Vec3 nrvo(v1[0]/v2[0],v1[1]/v2[1],v1[2]/v2[2]);
104  return nrvo;
105 }
106 
108 inline Vec3 operator*(const Vec3& v,const Mat3& m)
109 {
110  Vec3 nrvo(v[0]*m(0,0)+v[1]*m(1,0)+v[2]*m(2,0),
111  v[0]*m(0,1)+v[1]*m(1,1)+v[2]*m(2,1),
112  v[0]*m(0,2)+v[1]*m(1,2)+v[2]*m(2,2));
113  return nrvo;
114 }
115 
117 inline Vec3 operator*(const Mat3& m, const Vec3& v)
118 {
119  Vec3 nrvo(v[0]*m(0,0)+v[1]*m(0,1)+v[2]*m(0,2),
120  v[0]*m(1,0)+v[1]*m(1,1)+v[2]*m(1,2),
121  v[0]*m(2,0)+v[1]*m(2,1)+v[2]*m(2,2));
122  return nrvo;
123 }
124 
126 inline Mat3 operator*(const Mat3& m1, const Mat3& m2)
127 {
128  Mat3 nrvo(m1(0,0)*m2(0,0)+m1(0,1)*m2(1,0)+m1(0,2)*m2(2,0),
129  m1(0,0)*m2(0,1)+m1(0,1)*m2(1,1)+m1(0,2)*m2(2,1),
130  m1(0,0)*m2(0,2)+m1(0,1)*m2(1,2)+m1(0,2)*m2(2,2),
131 
132  m1(1,0)*m2(0,0)+m1(1,1)*m2(1,0)+m1(1,2)*m2(2,0),
133  m1(1,0)*m2(0,1)+m1(1,1)*m2(1,1)+m1(1,2)*m2(2,1),
134  m1(1,0)*m2(0,2)+m1(1,1)*m2(1,2)+m1(1,2)*m2(2,2),
135 
136  m1(2,0)*m2(0,0)+m1(2,1)*m2(1,0)+m1(2,2)*m2(2,0),
137  m1(2,0)*m2(0,1)+m1(2,1)*m2(1,1)+m1(2,2)*m2(2,1),
138  m1(2,0)*m2(0,2)+m1(2,1)*m2(1,2)+m1(2,2)*m2(2,2));
139  return nrvo;
140 }
141 
142 Mat3 DLLEXPORT_OST_GEOM Invert(const Mat3& m);
143 Mat3 DLLEXPORT_OST_GEOM Transpose(const Mat3& m);
144 
145 // algebraic complement
146 Real DLLEXPORT_OST_GEOM Comp(const Mat3& m, unsigned int i, unsigned int j);
147 
148 // minor
149 Real DLLEXPORT_OST_GEOM Minor(const Mat3& m, unsigned int i, unsigned int j);
150 
151 // determinant
152 Real DLLEXPORT_OST_GEOM Det(const Mat3& m);
153 
154 // angle between two vectors
155 Real DLLEXPORT_OST_GEOM Angle(const Vec3& v1, const Vec3& v2);
156 
157 // signed angle between two vectors, based on a reference normal
158 Real DLLEXPORT_OST_GEOM SignedAngle(const Vec3& v1, const Vec3& v2, const Vec3& ref);
159 
161 
162 Mat3 DLLEXPORT_OST_GEOM AxisRotation(const Vec3& axis, Real angle);
163 
167 Vec3 DLLEXPORT_OST_GEOM OrthogonalVector(const Vec3& axis);
168 
169 
170 Real DLLEXPORT_OST_GEOM DihedralAngle(const Vec3& p1, const Vec3& p2,
171  const Vec3& p3, const Vec3&p4);
172 
174 inline Vec3 Min(const Vec3& v1, const Vec3& v2)
175 {
176  Vec3 nrvo(std::min(v1[0],v2[0]),
177  std::min(v1[1],v2[1]),
178  std::min(v1[2],v2[2]));
179  return nrvo;
180 }
181 
183 inline Vec3 Max(const Vec3& v1, const Vec3& v2)
184 {
185  Vec3 nrvo(std::max(v1[0],v2[0]),
186  std::max(v1[1],v2[1]),
187  std::max(v1[2],v2[2]));
188  return nrvo;
189 }
190 
192 inline Real Distance(const Vec3& p1, const Vec3& p2)
193 {
194  return Length(p1-p2);
195 }
196 
197 
199 inline Real Distance2WithPBC(const Vec3& v1, const Vec3& v2, const Vec3& basis_vec){
200  Vec3 v;
201  v=v1-v2;
202  for (int i=0; i<3; i++) {
203  if (std::fabs(v[i])>basis_vec[i]/2.){
204  v[i]=std::fabs(v[i])-basis_vec[i]*int(std::fabs(v[i])/basis_vec[i]+0.5);
205  }
206  }
207  return Length2(v);
208 }
210 inline Real DistanceWithPBC(const Vec3& v1, const Vec3& v2, const Vec3& basis_vec){
211  return sqrt(Distance2WithPBC(v1, v2, basis_vec));
212 }
214 Real DLLEXPORT_OST_GEOM MinDistance(const Vec3List& l1, const Vec3List& l2);
216 // with periodic boundaries in x,y,z given in basis_vec
217 Real DLLEXPORT_OST_GEOM MinDistanceWithPBC(const Vec3List& l1, const Vec3List& l2, Vec3& basis_vec);
220 std::vector<unsigned int> DLLEXPORT_OST_GEOM MinDistanceIndices(const Vec3List& l1, const Vec3List& l2);
221 
223 Vec3 DLLEXPORT_OST_GEOM WrapVec3(const Vec3& v1,const Vec3& box_center,const Vec3& basis_vec);
225 Vec3List DLLEXPORT_OST_GEOM WrapVec3List(const Vec3List& vl,const Vec3& box_center,const Vec3& basis_vec);
226 
227 
228 } // ns
229 
230 #endif