1 // Copyright (C) 2007-2016 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, or (at your option) any later version.
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__
31 namespace INTERP_KERNEL
33 /// Precision used for tests of 3D part of INTERP_KERNEL
34 const double VOL_PREC = 1.0e-6;
36 /// Default relative tolerance in epsilonEqualRelative
37 const double DEFAULT_REL_TOL = 1.0e-6;
39 /// Default absolute tolerance in epsilonEqual and epsilonEqualRelative
40 const double DEFAULT_ABS_TOL = 5.0e-12;
43 * @param a first point. Should point on a array of size at least equal to SPACEDIM.
44 * @param b second point. Should point on a array of size at least equal to SPACEDIM.
46 template<int SPACEDIM>
47 inline double getDistanceBtw2Pts(const double *a, const double *b)
50 for(int i=0;i<SPACEDIM;i++)
51 ret2+=(a[i]-b[i])*(a[i]-b[i]);
55 // -------------------------------------------------------------------
56 // Math operations for vectors represented by double[3] - arrays
57 // -------------------------------------------------------------------
60 * Copies a double[3] vector from src to dest
62 * @param src source vector
63 * @param dest destination vector
66 inline void copyVector3(const double* src, double* dest)
68 for(int i = 0 ; i < 3 ; ++i)
73 * Creates a string representation of a double[3] vector
75 * @param pt a 3-vector
76 * @return a string of the form [x, y, z]
78 inline const std::string vToStr(const double* pt)
80 std::stringstream ss(std::ios::out);
81 ss << "[" << pt[0] << ", " << pt[1] << ", " << pt[2] << "]";
86 * Adds a double[3] - vector to another one.
89 * @param res vector in which to store the result res + v.
91 inline void add(const double* v, double* res)
99 * Calculates the cross product of two double[3] - vectors.
101 * @param v1 vector v1
102 * @param v2 vector v2
103 * @param res vector in which to store the result v1 x v2. It should not be one of v1 and v2.
105 inline void cross(const double* v1, const double* v2,double* res)
107 res[0] = v1[1]*v2[2] - v1[2]*v2[1];
108 res[1] = v1[2]*v2[0] - v1[0]*v2[2];
109 res[2] = v1[0]*v2[1] - v1[1]*v2[0];
113 * Calculates the dot product of two double[3] - vectors
115 * @param v1 vector v1
116 * @param v2 vector v2
117 * @return dot (scalar) product v1.v2
119 inline double dot(const double* v1, const double* v2)
121 return v1[0]*v2[0] + v1[1]*v2[1] + v1[2]*v2[2];
125 * Calculates norm of a double[3] vector
127 * @param v a vector v
128 * @return euclidean norm of v
130 inline double norm(const double* v)
132 return sqrt(dot(v,v));
136 * Compares doubles using an absolute tolerance
137 * This is suitable mainly for comparisons with 0.0
139 * @param x first value
140 * @param y second value
141 * @param errTol maximum allowed absolute difference that is to be treated as equality
142 * @return true if |x - y| < errTol, false otherwise
144 inline bool epsilonEqual(const double x, const double y, const double errTol = DEFAULT_ABS_TOL)
146 return y < x ? x - y < errTol : y - x < errTol;
147 // return std::fabs(x - y) < errTol;
152 * Test whether two 3D vectors are colinear. The two vectors are expected to be of unit norm (not checked)
153 * Implemented by checking that the norm of the cross product is null.
155 inline bool isColinear3D(const double *v1, const double *v2, const double eps = DEFAULT_ABS_TOL)
159 return epsilonEqual(dot(cros, cros), 0.0, eps);
164 * Compares doubles using a relative tolerance
165 * This is suitable mainly for comparing larger values to each other. Before performing the relative test,
166 * an absolute test is performed to guard from problems when comparing to 0.0
168 * @param x first value
169 * @param y second value
170 * @param relTol maximum allowed relative difference that is to be treated as equality
171 * @param absTol maximum allowed absolute difference that is to be treated as equality
172 * @return true if |x - y| <= absTol or |x - y|/max(|x|,|y|) <= relTol, false otherwise
174 inline bool epsilonEqualRelative(const double x, const double y, const double relTol = DEFAULT_REL_TOL, const double absTol = DEFAULT_ABS_TOL)
176 // necessary for comparing values close to zero
177 // in order to avoid division by very small numbers
178 if(std::fabs(x - y) < absTol)
183 const double relError = std::fabs((x - y) / std::max(std::fabs(x), std::fabs(y)));
185 return relError < relTol;