1#ifndef POINT_HPP
2#define POINT_HPP
3
4#include <boost/fusion/sequence/intrinsic/at_c.hpp>
5#include <boost/fusion/include/at_c.hpp>
6#include <boost/fusion/container/vector.hpp>
7#include <boost/fusion/include/vector.hpp>
8#include <boost/fusion/container/vector/vector_fwd.hpp>
9#include <boost/fusion/include/vector_fwd.hpp>
10#include "boost/multi_array.hpp"
11#include "memory_ly/Encap.hpp"
12#include "Point_operators.hpp"
13
14
15/*! \brief This class implement the point shape in an N-dimensional space
16 *
17 * \param T type of the space
18 * \param dim dimensionality
19 *
20 */
21
22template<unsigned int dim ,typename T> class Point
23{
24 public:
25
26 //! Indicate that this object is vtk writable
27 typedef int is_vtk_writable;
28
29 //! coordinate type
30 typedef T coord_type;
31
32 //! boost fusion that store the point
33 typedef boost::fusion::vector<T[dim]> type;
34
35 //! Indicate that is a Point
36 typedef int yes_is_point;
37
38 //! structure that store the data of the point
39 type data;
40
41 //! Property id of the point
42 static const unsigned int x = 0;
43
44
45 /*! \brief Evaluate the expression and save the result on the point
46 *
47 * \param p_exp point expression to evaluate
48 *
49 */
50 template<typename orig, typename exp1, typename exp2, unsigned int op> __device__ __host__ Point(const point_expression_op<orig,exp1,exp2,op> & p_exp)
51 {
52 this->operator=(p_exp);
53 }
54
55 /*! \brief Point constructor from point
56 *
57 * \param p the point
58 *
59 */
60 __device__ __host__ inline Point(const Point<dim,T> && p)
61 {
62 for(size_t i = 0; i < dim ; i++)
63 {get(i) = p.get(i);}
64 }
65
66 /*! \brief Point constructor from point
67 *
68 * \param p the point
69 *
70 */
71 __device__ __host__ inline Point(const Point<dim,T> & p)
72 {
73 for(size_t i = 0; i < dim ; i++)
74 {get(i) = p.get(i);}
75 }
76
77 /*! \brief Constructor from an array
78 *
79 * \param p array with the coordinate of the point
80 *
81 */
82 __device__ __host__ inline Point(const T (&p)[dim])
83 {
84 for(size_t i = 0; i < dim ; i++)
85 {get(i) = p[i];}
86 }
87
88 /*! \brief Constructor from scalar
89 *
90 * \param d scalar
91 *
92 */
93 __device__ __host__ inline Point(T d)
94 {
95 this->operator=(d);
96 }
97
98 /*! \brief Point constructor
99 *
100 * \param p Point
101 *
102 */
103 template <typename S> __device__ __host__ inline Point(const Point<dim,S> & p)
104 {
105 for (size_t i = 0 ; i < dim ; i++)
106 {get(i) = static_cast<S>(p.get(i));}
107 }
108
109 /*! \brief Point constructor
110 *
111 * \param p encapc Point
112 *
113 */
114 template <unsigned int d, typename M> inline Point(const encapc<d,Point<dim,T>,M> & p)
115 {
116 for (size_t i = 0 ; i < dim ; i++)
117 get(i) = p.template get<0>()[i];
118 }
119
120 /*! \brief Point constructor from multi array
121 *
122 *
123 *
124 */
125 template <typename vmpl> inline __device__ __host__ Point(const openfpm::detail::multi_array::sub_array_openfpm<T,1,vmpl> & mar)
126 {
127 for (unsigned int i = 0 ; i < dim ; i++)
128 {get(i) = mar[i];}
129 }
130
131 /*! \brief Point constructor from multi array
132 *
133 *
134 *
135 */
136 template <typename vmpl> inline __device__ __host__ Point(const openfpm::detail::multi_array::const_sub_array_openfpm<T,1,vmpl> & mar)
137 {
138 for (unsigned int i = 0 ; i < dim ; i++)
139 {get(i) = mar[i];}
140 }
141
142 /*! \brief Constructor from a list
143 *
144 * [Example] Point<3,float> p({0.0,0.0,1.0})
145 *
146 * \param p1 initializer list
147 *
148 */
149 __device__ __host__ inline Point(std::initializer_list<T> p1)
150 {
151 size_t i = 0;
152 for(T x : p1)
153 {get(i) = x;i++;}
154 }
155
156 //! Default contructor
157 __device__ __host__ inline Point()
158 {}
159
160 /*! \brief Get coordinate
161 *
162 * \param i dimension
163 * \return the i-coordinate of the point
164 *
165 */
166 __device__ __host__ inline const T & get(int i) const
167 {
168 return boost::fusion::at_c<x>(data)[i];
169 }
170
171 /*! \brief Get coordinate
172 *
173 * \param i dimension
174 * \return the i-coordinate of the point
175 *
176 */
177
178 inline T get_vtk(size_t i) const
179 {
180 return boost::fusion::at_c<x>(data)[i];
181 }
182
183 /*! \brief Get coordinate
184 *
185 * \param i dimension
186 * \return the i-coordinate of the point
187 *
188 */
189 __device__ __host__ inline T& get(int i)
190 {
191 return boost::fusion::at_c<x>(data)[i];
192 }
193
194 /*! \brief Get the component i
195 *
196 * \param i component
197 *
198 * \return the i-component
199 *
200 */
201
202 __device__ __host__ inline T& operator[](size_t i)
203 {
204 return get(i);
205 }
206
207 /*! \brief Get the component i
208 *
209 * \param i component
210 *
211 * \return the i-component
212 *
213 */
214
215 __device__ __host__ inline const T& operator[](size_t i) const
216 {
217 return get(i);
218 }
219
220 /*! \brief norm of the vector
221 *
222 * \return the norm of the vector
223 *
224 */
225 __device__ __host__ T norm()
226 {
227 T n = 0.0;
228
229 for (size_t i = 0 ; i < dim ; i++)
230 n+=get(i) * get(i);
231
232 return sqrt(n);
233 }
234
235 /*! \brief It calculate the distance between 2 points
236 *
237 * The distance between itself (p) and the other point (q)
238 *
239 * \param q target point
240 *
241 * \return the distance
242 *
243 */
244 __device__ __host__ T distance(const Point<dim,T> & q) const
245 {
246 T tot = 0.0;
247
248 for (size_t i = 0 ; i < dim ; i++)
249 tot += (this->get(i) - q.get(i)) * (this->get(i) - q.get(i));
250
251 return sqrt(tot);
252 }
253
254 /*! \brief It calculate the square distance between 2 points
255 *
256 * The distance between itself (p) and the other point (q)
257 *
258 * \param q target point
259 *
260 * \return the square of the distance
261 *
262 */
263 T distance2(const Point<dim,T> & q) const
264 {
265 T tot = 0.0;
266
267 for (size_t i = 0 ; i < dim ; i++)
268 tot += (this->get(i) - q.get(i)) * (this->get(i) - q.get(i));
269
270 return tot;
271 }
272
273
274 /*! \brief Set to zero the point coordinate
275 *
276 *
277 */
278 __device__ __host__ inline void zero()
279 {
280 for (size_t i = 0 ; i < dim ; i++)
281 {
282 get(i) = 0;
283 }
284 }
285
286 /*! \brief Set to one the point coordinate
287 *
288 *
289 */
290 inline void one()
291 {
292 for (size_t i = 0 ; i < dim ; i++)
293 {
294 get(i) = 1;
295 }
296 }
297
298 /*! \brief Create a point set to zero
299 *
300 * \return a point with all coorfinate set to 0
301 *
302 */
303 inline static Point<dim,T> zero_p()
304 {
305 Point<dim,T> p;
306
307 for (size_t i = 0 ; i < dim ; i++)
308 {
309 p.get(i) = 0;
310 }
311
312 return p;
313 }
314
315 /*! \brief Convert the point into a string
316 *
317 * \return the string
318 *
319 */
320 std::string toPointString() const
321 {
322 std::stringstream ps;
323
324 for (size_t i = 0 ; i < dim ; i++)
325 ps << "x[" << i << "]=" << get(i) << " ";
326
327 ps << "\n";
328
329 return ps.str();
330 }
331
332 /*! \brief exchange the data of two points
333 *
334 * \param p Point to swap with
335 *
336 */
337 void swap(Point<dim,T> & p)
338 {
339 for (size_t i = 0 ; i < dim ; i++)
340 {
341 T tmp = get(i);
342 get(i) = p.get(i);
343 p.get(i) = tmp;
344 }
345 }
346
347 /*! \brief Check if two points match
348 *
349 * \param p point to compare with
350 *
351 * \return true if two points match
352 *
353 */
354 __device__ __host__ inline bool operator==(const Point<dim,T> & p) const
355 {
356 for (size_t i = 0 ; i < dim ; i++)
357 {
358 if (p.get(i) != get(i))
359 {return false;}
360 }
361
362 return true;
363 }
364
365 /*! \brief Check if two points match
366 *
367 * \param p point to compare with
368 *
369 * \return true if two points match
370 *
371 */
372 __device__ __host__ inline bool operator!=(const Point<dim,T> & p) const
373 {
374 return !this->operator==(p);
375 }
376
377 /*! \brief Return the string with the point coordinate
378 *
379 * \return the string
380 *
381 */
382 std::string to_string() const
383 {
384 return toString();
385 }
386
387 /*! \brief Return the string with the point coordinate
388 *
389 * \return the string
390 *
391 */
392 std::string toString() const
393 {
394 std::string str;
395
396 for (size_t i = 0 ; i < dim - 1 ; i++)
397 {
398 // coverty[dead_error_line]
399 str += std::to_string(static_cast<double>(get(i))) + " ";
400 }
401 str += std::to_string(static_cast<double>(get(dim-1)));
402
403 return str;
404 }
405
406 /*! \brief Return the reference to the value at coordinate i
407 *
408 * \param i coordinate to return
409 *
410 * \return the reference
411 *
412 */
413 T & value(size_t i)
414 {
415 return get(i);
416 }
417
418 /*! \brief Return the value at coordinate i
419 *
420 * \param i coordinate to return
421 *
422 * \return the value
423 *
424 */
425 __device__ __host__ inline T value(size_t i) const
426 {
427 return get(i);
428 }
429
430 /*! \brief Return the coordinated of the point as reference array
431 *
432 * \return the reference array
433 *
434 */
435 T (&asArray())[dim]
436 {
437 return boost::fusion::at_c<x>(data);
438 }
439
440 /*! Convert the point from Point<dim,T> to Point<dim,A>
441 *
442 * \return the converted point
443 *
444 */
445 template<typename A> Point<dim,A> convertPoint() const
446 {
447 Point<dim,A> p;
448
449 for (size_t i = 0; i < dim ; i++)
450 p.get(i) = static_cast<A>(get(i));
451
452 return p;
453 }
454
455
456 //! This structure has no internal pointers
457 static bool noPointers()
458 {
459 return true;
460 }
461
462 ////////////////////////////////////////////////////////////////
463 ////////////////////// ARITMETIC OPERATORS /////////////////////
464 ////////////////////////////////////////////////////////////////
465
466 /*! \brief Fill the vector with the evaluated expression
467 *
468 * \param p_exp expression to evaluate
469 *
470 * \return itself
471 *
472 */
473 template<typename orig, typename exp1, typename exp2, unsigned int op> __device__ __host__ Point<dim,T> & operator=(const point_expression_op<orig,exp1,exp2,op> & p_exp)
474 {
475 p_exp.init();
476
477 for (size_t i = 0; i < dim ; i++)
478 {get(i) = p_exp.value(i);}
479
480 return *this;
481 }
482
483 /*! \brief Fill the vector property with the evaluated expression
484 *
485 * \param p_exp expression to evaluate
486 *
487 * \return itself
488 *
489 */
490 __device__ __host__ Point<dim,T> & operator=(const point_expression<T[dim]> & p_exp)
491 {
492 p_exp.init();
493
494 for (size_t i = 0; i < dim ; i++)
495 get(i) = p_exp.value(i);
496
497 return *this;
498 }
499
500 /*! \brief Fill the point with the value specified in the array
501 *
502 * \param p array
503 *
504 * \return itself
505 *
506 */
507 __device__ __host__ Point<dim,T> & operator=(const T (& p)[dim])
508 {
509 for (size_t i = 0; i < dim ; i++)
510 get(i) = p[i];
511
512 return *this;
513 }
514
515 /*! \brief Fill the vector property with the evaluated expression
516 *
517 * \tparam check disable this method if T is a constant
518 *
519 * \param p_exp expression to evaluate
520 *
521 * \return itself
522 *
523 */
524 template<typename T1, typename check = typename std::enable_if<std::is_const<T1>::value == false>::type> __device__ __host__ Point<dim,T> & operator=(const point_expression<const T1[dim]> & p_exp)
525 {
526 p_exp.init();
527
528 for (size_t i = 0; i < dim ; i++)
529 get(i) = p_exp.value(i);
530
531 return *this;
532 }
533
534 /*! \brief divide each component by an array
535 *
536 * \param ar Component wise division
537 *
538 * \return itself
539 *
540 */
541 template<typename aT> __device__ __host__ inline Point<dim,T> operator/(const aT (&ar)[dim])
542 {
543 Point<dim,T> result;
544
545 for (size_t i = 0 ; i < dim ; i++)
546 result.get(i) = get(i) / ar[i];
547
548 return result;
549 }
550
551 /*! \brief divide each component by a constant
552 *
553 * \param c constant
554 *
555 * \return itself
556 *
557 */
558 template<typename aT> __device__ __host__ inline Point<dim,T> operator/=(const aT c)
559 {
560 for (size_t i = 0 ; i < dim ; i++)
561 get(i) = get(i) / c;
562
563 return *this;
564 }
565
566
567 /*! \brief Fill the vector property with some value
568 *
569 * \param d value to fill
570 *
571 * \return itself
572 *
573 */
574 __device__ __host__ Point<dim,T> & operator=(T d)
575 {
576 for (size_t i = 0; i < dim ; i++)
577 get(i) = d;
578
579 return *this;
580 }
581
582 /*! \brief operator= between points
583 *
584 * \param p Point
585 *
586 * \return itself
587 *
588 */
589 __device__ __host__ inline Point<dim,T> & operator=(const Point<dim,T> & p)
590 {
591 for (size_t i = 0 ; i < dim ; i++)
592 get(i) = p.get(i);
593
594 return *this;
595 }
596
597 /*! \brief Subtract two points
598 *
599 * \param p point to subtract
600 *
601 * \return itself
602 *
603 */
604 __device__ __host__ inline Point<dim,T> & operator-=(const Point<dim,T> & p)
605 {
606 for (size_t i = 0 ; i < dim ; i++)
607 get(i) -= p.get(i);
608
609 return *this;
610 }
611
612 /*! \brief Sum two points
613 *
614 * \param p point to sum
615 *
616 * \return itself
617 *
618 */
619 __device__ __host__ inline Point<dim,T> & operator+=(const Point<dim,T> & p)
620 {
621 for (size_t i = 0 ; i < dim ; i++)
622 {get(i) += p.get(i);}
623
624 return *this;
625 }
626
627 /*! \brief Do nothing stub operation
628 *
629 * Required to make the code compilable
630 *
631 */
632 __device__ __host__ inline void init() const
633 {}
634
635 ///////////////////////////////////////////////////////////////////////
636 ///////////////////////////////////////////////////////////////////////
637
638 //! The point has one property
639 static const unsigned int max_prop = 1;
640 static const unsigned int max_prop_real = 1;
641
642 //! expose the dimension
643 static const unsigned int dims = dim;
644
645 //! expose the dimension with a different name
646 static const unsigned int nvals = dim;
647};
648
649/*! \brief Convert an array of point coordinate into string
650 *
651 * \param p coordinate on each dimension
652 *
653 * \return the string
654 *
655 */
656template <unsigned int N, typename T> std::string toPointString(const T (&p)[N] )
657{
658 std::stringstream ps;
659
660 for (size_t i = 0 ; i < N ; i++)
661 ps << "x[" << i << "]=" << p[i] << " ";
662
663 ps << "\n";
664
665 return ps.str();
666}
667
668/*! \brief Convert an encapsulated point into string
669 *
670 * \param p coordinate on each dimension
671 *
672 * \return the string
673 *
674 */
675template <unsigned int N, typename T, typename Mem> std::string toPointString(const encapc<1,Point<N,T>,Mem> & p )
676{
677 std::stringstream ps;
678
679 for (size_t i = 0 ; i < N ; i++)
680 ps << "x[" << i << "]=" << p.template get<Point<N,T>::x>()[i] << " ";
681
682 ps << "\n";
683
684 return ps.str();
685}
686
687
688//! A point is a vector on a computer (But do not say this to a Mathematician)
689
690template<unsigned int dim, typename T> using VectorS = Point<dim,T>;
691
692template<typename T, typename Sfinae = void>
693struct is_Point: std::false_type {};
694
695
696/*! \brief Check if a type T is an aggregate
697 *
698 * return true if T is an aggregate
699 *
700 */
701template<typename T>
702struct is_Point<T, typename Void< typename T::yes_is_point>::type> : std::true_type
703{};
704
705#endif
706