Rheolef  7.2
an efficient C++ finite element environment
 
Loading...
Searching...
No Matches
point.h
Go to the documentation of this file.
1# ifndef _RHEO_BASIC_POINT_H
2# define _RHEO_BASIC_POINT_H
3//
4// This file is part of Rheolef.
5//
6// Copyright (C) 2000-2009 Pierre Saramito <Pierre.Saramito@imag.fr>
7//
8// Rheolef is free software; you can redistribute it and/or modify
9// it under the terms of the GNU General Public License as published by
10// the Free Software Foundation; either version 2 of the License, or
11// (at your option) any later version.
12//
13// Rheolef is distributed in the hope that it will be useful,
14// but WITHOUT ANY WARRANTY; without even the implied warranty of
15// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16// GNU General Public License for more details.
17//
18// You should have received a copy of the GNU General Public License
19// along with Rheolef; if not, write to the Free Software
20// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
21//
22// =========================================================================
23// author: Pierre.Saramito@imag.fr
24
25namespace rheolef {
77} // namespace rheolef
78
79#include "rheolef/compiler_mpi.h"
80#include <sstream>
81
82namespace rheolef {
83
84// [verbatim_point_basic]
85template <class T>
87public:
88
89// typedefs:
90
91 typedef size_t size_type;
92 typedef T element_type;
93 typedef T scalar_type;
94 typedef T float_type;
95
96// allocators:
97
98 explicit point_basic();
99 explicit point_basic (const T& x0, const T& x1 = 0, const T& x2 = 0);
100
101 template <class T1>
103
104 template <class T1>
106
107 point_basic (const std::initializer_list<T>& il);
108
109// accessors:
110
111 T& operator[](int i_coord) { return _x[i_coord%3]; }
112 T& operator()(int i_coord) { return _x[i_coord%3]; }
113 const T& operator[](int i_coord) const { return _x[i_coord%3]; }
114 const T& operator()(int i_coord) const { return _x[i_coord%3]; }
115
116// algebra:
117
118 bool operator== (const point_basic<T>& v) const;
119 bool operator!= (const point_basic<T>& v) const;
127
128 template <class U>
129 typename
130 std::enable_if<
133 >::type
134 operator* (const U& a) const;
135 point_basic<T> operator/ (const T& a) const;
137
138// i/o:
139
140 std::istream& get (std::istream& s, int d = 3);
141 std::ostream& put (std::ostream& s, int d = 3) const;
142// [verbatim_point_basic]
143
144// interface for CGAL library inter-operability:
145
146 const T& x() const { return _x[0]; }
147 const T& y() const { return _x[1]; }
148 const T& z() const { return _x[2]; }
149 T& x(){ return _x[0]; }
150 T& y(){ return _x[1]; }
151 T& z(){ return _x[2]; }
152
153// data:
154// protected:
155 T _x[3];
156// internal:
157 static T _my_abs(const T& x) { return (x > T(0)) ? x : -x; }
158// [verbatim_point_basic_cont]
159};
160// [verbatim_point_basic_cont]
161
162// [verbatim_point]
164// [verbatim_point]
165
166// [verbatim_point_basic_cont2]
167template<class T>
168std::istream& operator >> (std::istream& s, point_basic<T>& p);
169
170template<class T>
171std::ostream& operator << (std::ostream& s, const point_basic<T>& p);
172
173template <class T, class U>
174typename
175std::enable_if<
176 details::is_rheolef_arithmetic<U>::value
178>::type
179operator* (const U& a, const point_basic<T>& u);
180
181template<class T>
183vect (const point_basic<T>& v, const point_basic<T>& w);
184
185// metrics:
186template<class T>
187T dot (const point_basic<T>& x, const point_basic<T>& y);
188
189template<class T>
190T norm2 (const point_basic<T>& x);
191
192template<class T>
193T norm (const point_basic<T>& x);
194
195template<class T>
196T dist2 (const point_basic<T>& x, const point_basic<T>& y);
197
198template<class T>
199T dist (const point_basic<T>& x, const point_basic<T>& y);
200
201template<class T>
202T dist_infty (const point_basic<T>& x, const point_basic<T>& y);
203
204template <class T>
205T vect2d (const point_basic<T>& v, const point_basic<T>& w);
206
207template <class T>
208T mixt (const point_basic<T>& u, const point_basic<T>& v, const point_basic<T>& w);
209
210// robust(exact) floating point predicates: return the sign of the value as (0, > 0, < 0)
211// formally: orient2d(a,b,x) = vect2d(a-x,b-x)
212template <class T>
213int sign_orient2d (
214 const point_basic<T>& a,
215 const point_basic<T>& b,
216 const point_basic<T>& c);
217
218template <class T>
219int sign_orient3d (
220 const point_basic<T>& a,
221 const point_basic<T>& b,
222 const point_basic<T>& c,
223 const point_basic<T>& d);
224
225// compute also the value:
226template <class T>
227T orient2d(
228 const point_basic<T>& a,
229 const point_basic<T>& b,
230 const point_basic<T>& c);
231
232// formally: orient3d(a,b,c,x) = mixt3d(a-x,b-x,c-x)
233template <class T>
234T orient3d(
235 const point_basic<T>& a,
236 const point_basic<T>& b,
237 const point_basic<T>& c,
238 const point_basic<T>& d);
239
240template <class T>
241std::string ptos (const point_basic<T>& x, int d = 3);
242
243// ccomparators: lexicographic order
244template<class T, size_t d>
246// [verbatim_point_basic_cont2]
247// -------------------------------------------------------------------------------------
248// inline'd
249// -------------------------------------------------------------------------------------
250template <class T, class U>
251inline
252typename
253std::enable_if<
254 details::is_rheolef_arithmetic<U>::value
256>::type
257operator* (const U& a, const point_basic<T>& u)
258{
259 return point_basic<T> (a*u[0], a*u[1], a*u[2]);
260}
261template<class T>
262inline
263point_basic<T>
265{
266 return point_basic<T> (
267 v[1]*w[2]-v[2]*w[1],
268 v[2]*w[0]-v[0]*w[2],
269 v[0]*w[1]-v[1]*w[0]);
270}
271// metrics:
272template<class T>
273inline
274T dot (const point_basic<T>& x, const point_basic<T>& y)
275{
276 return x[0]*y[0]+x[1]*y[1]+x[2]*y[2];
277}
278template<class T>
279inline
281{
282 return dot(x,x);
283}
284template<class T>
285inline
287{
288 return sqrt(norm2(x));
289}
290template<class T>
291inline
293{
294 return norm2(x-y);
295}
296template<class T>
297inline
299{
300 return norm(x-y);
301}
302template<class T>
303inline
305{
306 return max(point_basic<T>::_my_abs(x[0]-y[0]),
307 max(point_basic<T>::_my_abs(x[1]-y[1]),
308 point_basic<T>::_my_abs(x[2]-y[2])));
309}
310// ccomparators: lexicographic order
311template<class T, size_t d>
312inline
313bool
315{
316 for (typename point_basic<T>::size_type i = 0; i < d; i++) {
317 if (a[i] < b[i]) return true;
318 if (a[i] > b[i]) return false;
319 }
320 return false; // equality
321}
323template<class T> struct scalar_traits { typedef T type; };
324template<class T> struct scalar_traits<point_basic<T> > { typedef T type; };
325template<class T> struct float_traits<point_basic<T> > { typedef typename float_traits<T>::type type; };
326
327template<class T>
329 _x[0] = T();
330 _x[1] = T();
331 _x[2] = T();
332}
333template<class T>
335 const T& x0,
336 const T& x1,
337 const T& x2)
338{
339 _x[0] = x0;
340 _x[1] = x1;
341 _x[2] = x2;
342}
343template<class T>
344template<class T1>
346{
347 _x[0] = p._x[0];
348 _x[1] = p._x[1];
349 _x[2] = p._x[2];
350}
351template<class T>
352template <class T1>
355{
356 _x[0] = p._x[0];
357 _x[1] = p._x[1];
358 _x[2] = p._x[2];
359 return *this;
360}
361template<class T>
362point_basic<T>::point_basic (const std::initializer_list<T>& il) : _x() {
363 typedef typename std::initializer_list<T>::const_iterator const_iterator;
364 check_macro (il.size() <= 3, "unexpected initializer list size=" << il.size() << " > 3");
365 size_type i = 0;
366 for (const_iterator iter = il.begin(); iter != il.end(); ++iter, ++i) {
367 _x[i] = *iter;
368 }
369 for (i = il.size(); i < 3; ++i) {
370 _x[i] = T();
371 }
372}
373// input/output
374template<class T>
375std::istream&
376point_basic<T>::get (std::istream& s, int d)
377{
378 switch (d) {
379 case 0 : _x[0] = _x[1] = _x[2] = T(0); return s;
380 case 1 : _x[1] = _x[2] = T(0); return s >> _x[0];
381 case 2 : _x[2] = T(0); return s >> _x[0] >> _x[1];
382 default: return s >> _x[0] >> _x[1] >> _x[2];
383 }
384}
385template<class T>
386inline
387std::ostream&
388point_basic<T>::put (std::ostream& s, int d) const
389{
390 switch (d) {
391 case 0 : return s;
392 case 1 : return s << _x[0];
393 case 2 : return s << _x[0] << " " << _x[1];
394 default: return s << _x[0] << " " << _x[1] << " " << _x[2];
395 }
396}
397template<class T>
398inline
399std::istream&
400operator >> (std::istream& s, point_basic<T>& p)
401{
402 return s >> p[0] >> p[1] >> p[2];
403}
404template<class T>
405inline
406std::ostream&
407operator << (std::ostream& s, const point_basic<T>& p)
408{
409 return s << p[0] << " " << p[1] << " " << p[2];
410}
411template<class T>
412std::string
413ptos (const point_basic<T>& x, int d)
414{
415 std::ostringstream ostrstr;
416 x.put (ostrstr, d);
417 return ostrstr.str();
418}
419// ----------------------------------------------------------
420// point extra: inlined
421// ----------------------------------------------------------
422#define def_point_function2(f,F) \
423template<class T> \
424inline \
425point_basic<T> \
426f (const point_basic<T>& x) \
427{ \
428 point_basic<T> y; \
429 for (size_t i = 0; i < 3; i++) \
430 y[i] = F(x[i]); \
431 return y; \
432}
433
434#define def_point_function(f) def_point_function2(f,f)
435
441#undef def_point_function2
442#undef def_point_function
443
444template <class T>
445bool
447{
448 return _x[0] == v[0] && _x[1] == v[1] && _x[2] == v[2];
449}
450template <class T>
451bool
452point_basic<T>::operator!= (const point_basic<T>& v) const
453{
454 return !operator==(v);
455}
456template <class T>
457point_basic<T>&
458point_basic<T>::operator+= (const point_basic<T>& v)
459{
460 _x[0] += v[0];
461 _x[1] += v[1];
462 _x[2] += v[2];
463 return *this;
464}
465template <class T>
466point_basic<T>&
468{
469 _x[0] -= v[0];
470 _x[1] -= v[1];
471 _x[2] -= v[2];
472 return *this;
473}
474template <class T>
477{
478 _x[0] *= a;
479 _x[1] *= a;
480 _x[2] *= a;
481 return *this;
482}
483template <class T>
486{
487 _x[0] /= a;
488 _x[1] /= a;
489 _x[2] /= a;
490 return *this;
491}
492template <class T>
495{
496 return point_basic<T> (_x[0]+v[0],
497 _x[1]+v[1],
498 _x[2]+v[2]);
499}
500template <class T>
503{
504 return point_basic<T> (-_x[0],
505 -_x[1],
506 -_x[2]);
507}
508template <class T>
511{
512 return point_basic<T> (_x[0]-v[0],
513 _x[1]-v[1],
514 _x[2]-v[2]);
515}
516template<class T1, class T2>
517inline
519operator/ (const T2& a, const point_basic<T1>& x)
520{
522 for (size_t i = 0; i < 3; i++)
523 if (x[i] != 0) y[i] = a/x[i];
524 return y;
525}
526template <class T>
527template <class U>
528typename
529std::enable_if<
530 details::is_rheolef_arithmetic<U>::value
531 ,point_basic<T>
532>::type
534{
535 return point_basic<T> (_x[0]*a,
536 _x[1]*a,
537 _x[2]*a);
538}
539template <class T>
542{
543 return operator* (T(1)/T(a));
544}
545template <class T>
548{
549 return point_basic<T> (_x[0]/v[0],
550 _x[1]/v[1],
551 _x[2]/v[2]);
552}
553// vect2d() and mixt() are deduced from:
554template <class T>
555inline
556T
558{
559 return orient2d (point_basic<T>(), v, w);
560}
561template <class T>
562inline
563T
565{
566 return orient3d (point_basic<T>(), u, v, w);
567}
568
569}// namespace rheolef
570// -------------------------------------------------------------
571// serialization
572// -------------------------------------------------------------
573#ifdef _RHEOLEF_HAVE_MPI
574#include <boost/serialization/serialization.hpp>
575
576namespace boost {
577 namespace serialization {
578 template<class Archive, class T>
579 void serialize (Archive& ar, class rheolef::point_basic<T>& x, const unsigned int version) {
580 ar & x[0];
581 ar & x[1];
582 ar & x[2];
583 }
584 } // namespace serialization
585} // namespace boost
586
587// Some serializable types, like geo_element, have a fixed amount of data stored at fixed field positions.
588// When this is the case, boost::mpi can optimize their serialization and transmission to avoid extraneous copy operations.
589// To enable this optimization, we specialize the type trait is_mpi_datatype, e.g.:
590namespace boost {
591 namespace mpi {
592 // TODO: when point_basic<T> is not a simple type, such as T=bigfloat or T=gmp, etc
593 template <>
594 struct is_mpi_datatype<rheolef::point_basic<double> > : mpl::true_ { };
595 } // namespace mpi
596} // namespace boost
597#endif // _RHEOLEF_HAVE_MPI
598
599#endif /* _RHEO_BASIC_POINT_H */
point_basic< T > operator-() const
Definition point.h:502
const T & operator()(int i_coord) const
Definition point.h:114
std::istream & get(std::istream &s, int d=3)
Definition point.h:376
bool operator==(const point_basic< T > &v) const
T & operator()(int i_coord)
Definition point.h:112
point_basic< T > & operator+=(const point_basic< T > &v)
point_basic< T > operator+(const point_basic< T > &v) const
Definition point.h:494
std::ostream & put(std::ostream &s, int d=3) const
Definition point.h:388
bool operator!=(const point_basic< T > &v) const
point_basic< T > operator/(const T &a) const
Definition point.h:541
point_basic< T > & operator=(const point_basic< T1 > &p)
Definition point.h:354
const T & z() const
Definition point.h:148
std::enable_if< details::is_rheolef_arithmetic< U >::value, point_basic< T > >::type operator*(const U &a) const
Definition point.h:533
static T _my_abs(const T &x)
Definition point.h:157
point_basic< T > & operator*=(const T &a)
Definition point.h:476
point_basic< T > & operator/=(const T &a)
Definition point.h:485
const T & operator[](int i_coord) const
Definition point.h:113
point_basic(const T &x0, const T &x1=0, const T &x2=0)
Definition point.h:334
const T & x() const
Definition point.h:146
T & operator[](int i_coord)
Definition point.h:111
point_basic(const std::initializer_list< T > &il)
Definition point.h:362
point_basic< T > & operator-=(const point_basic< T > &v)
Definition point.h:467
const T & y() const
Definition point.h:147
point_basic< Float > point
Definition point.h:163
Expr1::float_type T
Definition field_expr.h:230
check_macro(expr1.have_homogeneous_space(Xh1), "dual(expr1,expr2); expr1 should have homogeneous space. HINT: use dual(interpolate(Xh, expr1),expr2)")
void serialize(Archive &ar, class rheolef::point_basic< T > &x, const unsigned int version)
Definition point.h:579
This file is part of Rheolef.
bool lexicographically_less(const point_basic< T > &a, const point_basic< T > &b)
Definition point.h:314
T orient2d(const point_basic< T > &a, const point_basic< T > &b, const point_basic< T > &c)
bool operator==(const heap_allocator< T1 > &lhs, const heap_allocator< T1 > &rhs)
std::ostream & operator<<(std::ostream &os, const catchmark &m)
Definition catchmark.h:99
point_basic< T > vect(const point_basic< T > &v, const point_basic< T > &w)
Definition point.h:264
dia< T, M > operator/(const T &lambda, const dia< T, M > &d)
Definition dia.h:145
T norm2(const vec< T, M > &x)
norm2(x): see the expression page for the full documentation
Definition vec.h:379
T dist_infty(const point_basic< T > &x, const point_basic< T > &y)
Definition point.h:304
T orient3d(const point_basic< T > &a, const point_basic< T > &b, const point_basic< T > &c, const point_basic< T > &d)
int sign_orient2d(const point_basic< T > &a, const point_basic< T > &b, const point_basic< T > &c)
int sign_orient3d(const point_basic< T > &a, const point_basic< T > &b, const point_basic< T > &c, const point_basic< T > &d)
T mixt(const point_basic< T > &u, const point_basic< T > &v, const point_basic< T > &w)
Definition point.h:564
rheolef::std enable_if ::type dot const Expr1 expr1, const Expr2 expr2 dot(const Expr1 &expr1, const Expr2 &expr2)
T dist2(const point_basic< T > &x, const point_basic< T > &y)
Definition point.h:292
std::string ptos(const point_basic< T > &x, int d=3)
Definition point.h:413
std::istream & operator>>(std::istream &is, const catchmark &m)
Definition catchmark.h:88
T vect2d(const point_basic< T > &v, const point_basic< T > &w)
Definition point.h:557
T dist(const point_basic< T > &x, const point_basic< T > &y)
Definition point.h:298
csr< T, sequential > operator*(const T &lambda, const csr< T, sequential > &a)
Definition csr.h:437
T norm(const vec< T, M > &x)
norm(x): see the expression page for the full documentation
Definition vec.h:387
#define def_point_function(f)
Definition point.h:434
Definition sphere.icc:25
float_traits< T >::type type
Definition point.h:325
helper for std::complex<T>: get basic T type
Definition Float.h:93
helper for point_basic<T> & tensor_basic<T>: get basic T type
Definition point.h:323
Definition leveque.h:25