point(2rheolef) | rheolef-7.0 | point(2rheolef) |
point - vertex of a mesh
Defines geometrical vertex as an array of coordinates. This array is also used as a vector of the three dimensional physical space.
template <class T> class point_basic {
public: // typedefs:
typedef size_t size_type;
typedef T element_type;
typedef T scalar_type;
typedef T float_type; // allocators:
explicit point_basic () { _x[0] = T(); _x[1] = T(); _x[2] = T(); }
explicit point_basic (
const T& x0,
const T& x1 = 0,
const T& x2 = 0)
{ _x[0] = x0; _x[1] = x1; _x[2] = x2; }
template <class T1>
point_basic<T>(const point_basic<T1>& p)
{ _x[0] = p._x[0]; _x[1] = p._x[1]; _x[2] = p._x[2]; }
template <class T1>
point_basic<T>& operator = (const point_basic<T1>& p)
{ _x[0] = p._x[0]; _x[1] = p._x[1]; _x[2] = p._x[2]; return *this; } #ifdef _RHEOLEF_HAVE_STD_INITIALIZER_LIST
point_basic (const std::initializer_list<T>& il); #endif // _RHEOLEF_HAVE_STD_INITIALIZER_LIST // accessors:
T& operator[](int i_coord) { return _x[i_coord%3]; }
const T& operator[](int i_coord) const { return _x[i_coord%3]; }
T& operator()(int i_coord) { return _x[i_coord%3]; }
const T& operator()(int i_coord) const { return _x[i_coord%3]; }
// interface for CGAL library inter-operability:
const T& x() const { return _x[0]; }
const T& y() const { return _x[1]; }
const T& z() const { return _x[2]; }
T& x(){ return _x[0]; }
T& y(){ return _x[1]; }
T& z(){ return _x[2]; } // inputs/outputs:
std::istream& get (std::istream& s, int d = 3)
{
switch (d) {
case 0 : _x[0] = _x[1] = _x[2] = T(0); return s;
case 1 : _x[1] = _x[2] = T(0); return s >> _x[0];
case 2 : _x[2] = T(0); return s >> _x[0] >> _x[1];
default: return s >> _x[0] >> _x[1] >> _x[2];
}
}
// output
std::ostream& put (std::ostream& s, int d = 3) const; // algebra:
bool operator== (const point_basic<T>& v) const
{ return _x[0] == v[0] && _x[1] == v[1] && _x[2] == v[2]; }
bool operator!= (const point_basic<T>& v) const
{ return !operator==(v); }
point_basic<T>& operator+= (const point_basic<T>& v)
{ _x[0] += v[0]; _x[1] += v[1]; _x[2] += v[2]; return *this; }
point_basic<T>& operator-= (const point_basic<T>& v)
{ _x[0] -= v[0]; _x[1] -= v[1]; _x[2] -= v[2]; return *this; }
point_basic<T>& operator*= (const T& a)
{ _x[0] *= a; _x[1] *= a; _x[2] *= a; return *this; }
point_basic<T>& operator/= (const T& a)
{ _x[0] /= a; _x[1] /= a; _x[2] /= a; return *this; }
point_basic<T> operator+ (const point_basic<T>& v) const
{ return point_basic<T> (_x[0]+v[0], _x[1]+v[1], _x[2]+v[2]); }
point_basic<T> operator- () const
{ return point_basic<T> (-_x[0], -_x[1], -_x[2]); }
point_basic<T> operator- (const point_basic<T>& v) const
{ return point_basic<T> (_x[0]-v[0], _x[1]-v[1], _x[2]-v[2]); }
template <class U>
typename
std::enable_if<
details::is_rheolef_arithmetic<U>::value
,point_basic<T>
>::type
operator* (const U& a) const
{ return point_basic<T> (_x[0]*a, _x[1]*a, _x[2]*a); }
point_basic<T> operator/ (const T& a) const
{ return operator* (T(1)/T(a)); }
point_basic<T> operator/ (point_basic<T> v) const
{ return point_basic<T> (_x[0]/v[0], _x[1]/v[1], _x[2]/v[2]); } // data: // protected:
T _x[3]; // internal:
static T _my_abs(const T& x) { return (x > T(0)) ? x : -x; } }; typedef point_basic<Float> point; // algebra: template <class T, class U> inline typename std::enable_if<
details::is_rheolef_arithmetic<U>::value
,point_basic<T> >::type operator* (const U& a, const point_basic<T>& u) {
return point_basic<T> (a*u[0], a*u[1], a*u[2]); } template<class T> inline point_basic<T> vect (const point_basic<T>& v, const point_basic<T>& w) {
return point_basic<T> (
v[1]*w[2]-v[2]*w[1],
v[2]*w[0]-v[0]*w[2],
v[0]*w[1]-v[1]*w[0]); } // metrics: template<class T> inline T dot (const point_basic<T>& x, const point_basic<T>& y) {
return x[0]*y[0]+x[1]*y[1]+x[2]*y[2]; } template<class T> inline T norm2 (const point_basic<T>& x) {
return dot(x,x); } template<class T> inline T norm (const point_basic<T>& x) {
return sqrt(norm2(x)); } template<class T> inline T dist2 (const point_basic<T>& x, const point_basic<T>& y) {
return norm2(x-y); } template<class T> inline T dist (const point_basic<T>& x, const point_basic<T>& y) {
return norm(x-y); } template<class T> inline T dist_infty (const point_basic<T>& x, const point_basic<T>& y) {
return max(point_basic<T>::_my_abs(x[0]-y[0]),
max(point_basic<T>::_my_abs(x[1]-y[1]),
point_basic<T>::_my_abs(x[2]-y[2]))); } template <class T> T vect2d (const point_basic<T>& v, const point_basic<T>& w); template <class T> T mixt (const point_basic<T>& u, const point_basic<T>& v, const point_basic<T>& w); // robust(exact) floating point predicates: return the sign of the value as (0, > 0, < 0) // formally: orient2d(a,b,x) = vect2d(a-x,b-x) template <class T> int sign_orient2d (
const point_basic<T>& a,
const point_basic<T>& b,
const point_basic<T>& c); template <class T> int sign_orient3d (
const point_basic<T>& a,
const point_basic<T>& b,
const point_basic<T>& c,
const point_basic<T>& d); // compute also the value: template <class T> T orient2d(
const point_basic<T>& a,
const point_basic<T>& b,
const point_basic<T>& c); // formally: orient3d(a,b,c,x) = mixt3d(a-x,b-x,c-x) template <class T> T orient3d(
const point_basic<T>& a,
const point_basic<T>& b,
const point_basic<T>& c,
const point_basic<T>& d); template <class T> std::string ptos (const point_basic<T>& x, int d = 3); // ccomparators: lexicographic order template<class T, size_t d> bool lexicographically_less (const point_basic<T>& a, const point_basic<T>& b) {
for (typename point_basic<T>::size_type i = 0; i < d; i++) {
if (a[i] < b[i]) return true;
if (a[i] > b[i]) return false;
}
return false; // equality }
Copyright (C) 2000-2018 Pierre Saramito <Pierre.Saramito@imag.fr> GPLv3+: GNU GPL version 3 or later <http://gnu.org/licenses/gpl.html>. This is free software: you are free to change and redistribute it. There is NO WARRANTY, to the extent permitted by law.
rheolef-7.0 | rheolef-7.0 |