1 // Copyright (C) 2007-2013 CEA/DEN, EDF R&D
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License.
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 // Lesser General Public License for more details.
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
20 #ifndef __VECTORUTILS_HXX__
21 #define __VECTORUTILS_HXX__
29 namespace INTERP_KERNEL
31 /// Precision used for tests of 3D part of INTERP_KERNEL
32 const double VOL_PREC = 1.0e-6;
34 /// Default relative tolerance in epsilonEqualRelative
35 const double DEFAULT_REL_TOL = 1.0e-6;
37 /// Default absolute tolerance in epsilonEqual and epsilonEqualRelative
38 const double DEFAULT_ABS_TOL = 5.0e-12;
41 * @param a first point. Should point on a array of size at least equal to SPACEDIM.
42 * @param b second point. Should point on a array of size at least equal to SPACEDIM.
44 template<int SPACEDIM>
45 inline double getDistanceBtw2Pts(const double *a, const double *b)
48 for(int i=0;i<SPACEDIM;i++)
49 ret2+=(a[i]-b[i])*(a[i]-b[i]);
53 // -------------------------------------------------------------------
54 // Math operations for vectors represented by double[3] - arrays
55 // -------------------------------------------------------------------
58 * Copies a double[3] vector from src to dest
60 * @param src source vector
61 * @param dest destination vector
64 inline void copyVector3(const double* src, double* dest)
66 for(int i = 0 ; i < 3 ; ++i)
71 * Creates a string representation of a double[3] vector
73 * @param pt a 3-vector
74 * @return a string of the form [x, y, z]
76 inline const std::string vToStr(const double* pt)
78 std::stringstream ss(std::ios::out);
79 ss << "[" << pt[0] << ", " << pt[1] << ", " << pt[2] << "]";
84 * Adds a double[3] - vector to another one.
87 * @param res vector in which to store the result res + v.
89 inline void add(const double* v, double* res)
97 * Calculates the cross product of two double[3] - vectors.
100 * @param v2 vector v2
101 * @param res vector in which to store the result v1 x v2. It should not be one of v1 and v2.
103 inline void cross(const double* v1, const double* v2,double* res)
105 res[0] = v1[1]*v2[2] - v1[2]*v2[1];
106 res[1] = v1[2]*v2[0] - v1[0]*v2[2];
107 res[2] = v1[0]*v2[1] - v1[1]*v2[0];
111 * Calculates the dot product of two double[3] - vectors
113 * @param v1 vector v1
114 * @param v2 vector v2
115 * @return dot (scalar) product v1.v2
117 inline double dot(const double* v1, const double* v2)
119 return v1[0]*v2[0] + v1[1]*v2[1] + v1[2]*v2[2];
123 * Calculates norm of a double[3] vector
125 * @param v a vector v
126 * @return euclidean norm of v
128 inline double norm(const double* v)
130 return sqrt(dot(v,v));
134 * Compares doubles using an absolute tolerance
135 * This is suitable mainly for comparisons with 0.0
137 * @param x first value
138 * @param y second value
139 * @param errTol maximum allowed absolute difference that is to be treated as equality
140 * @return true if |x - y| < errTol, false otherwise
142 inline bool epsilonEqual(const double x, const double y, const double errTol = DEFAULT_ABS_TOL)
144 return y < x ? x - y < errTol : y - x < errTol;
145 // return std::fabs(x - y) < errTol;
149 * Compares doubles using a relative tolerance
150 * This is suitable mainly for comparing larger values to each other. Before performing the relative test,
151 * an absolute test is performed to guard from problems when comparing to 0.0
153 * @param x first value
154 * @param y second value
155 * @param relTol maximum allowed relative difference that is to be treated as equality
156 * @param absTol maximum allowed absolute difference that is to be treated as equality
157 * @return true if |x - y| <= absTol or |x - y|/max(|x|,|y|) <= relTol, false otherwise
159 inline bool epsilonEqualRelative(const double x, const double y, const double relTol = DEFAULT_REL_TOL, const double absTol = DEFAULT_ABS_TOL)
161 // necessary for comparing values close to zero
162 // in order to avoid division by very small numbers
163 if(std::fabs(x - y) < absTol)
168 const double relError = std::fabs((x - y) / std::max(std::fabs(x), std::fabs(y)));
170 return relError < relTol;