Reference documentation for deal.II version GIT c00406db16 2023-09-28 18:00:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
point.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2023 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_point_h
17 #define dealii_point_h
18 
19 
20 #include <deal.II/base/config.h>
21 
23 #include <deal.II/base/tensor.h>
24 
25 #include <boost/geometry/core/cs.hpp>
26 #include <boost/geometry/geometries/point.hpp>
27 
28 #include <cmath>
29 
31 
109 template <int dim, typename Number = double>
110 DEAL_II_CXX20_REQUIRES(dim >= 0)
111 class Point : public Tensor<1, dim, Number>
112 {
113 public:
121  Point();
122 
126  explicit DEAL_II_HOST_DEVICE
128 
137  explicit DEAL_II_HOST_DEVICE
138  Point(const Number x);
139 
150  Point(const Number x, const Number y);
151 
162  Point(const Number x, const Number y, const Number z);
163 
167  template <std::size_t dummy_dim,
168  std::enable_if_t<(dim == dummy_dim) && (dummy_dim != 0), int> = 0>
169  Point(const boost::geometry::model::
170  point<Number, dummy_dim, boost::geometry::cs::cartesian> &boost_pt);
171 
180  unit_vector(const unsigned int i);
181 
187  DEAL_II_HOST_DEVICE Number
188  operator()(const unsigned int index) const;
189 
195  DEAL_II_HOST_DEVICE Number &
196  operator()(const unsigned int index);
197 
203  template <typename OtherNumber>
206 
219 
231 
242 
249  operator-() const;
250 
267  template <typename OtherNumber>
269  dim,
270  typename ProductType<Number,
271  typename EnableIfScalar<OtherNumber>::type>::type>
272  operator*(const OtherNumber) const;
273 
279  template <typename OtherNumber>
281  dim,
282  typename ProductType<Number,
283  typename EnableIfScalar<OtherNumber>::type>::type>
284  operator/(const OtherNumber) const;
285 
291  DEAL_II_HOST_DEVICE Number
293 
307  square() const;
308 
317  distance(const Point<dim, Number> &p) const;
318 
327 
337  template <class Archive>
338  void
339  serialize(Archive &ar, const unsigned int version);
340 };
341 
342 /*--------------------------- Inline functions: Point -----------------------*/
343 
344 #ifndef DOXYGEN
345 
346 // At least clang-3.7 requires us to have a user-defined constructor
347 // and we can't use 'Point<dim,Number>::Point () = default' here.
348 template <int dim, typename Number>
349 DEAL_II_CXX20_REQUIRES(dim >= 0)
350 inline DEAL_II_HOST_DEVICE Point<dim, Number>::Point() // NOLINT
351 {}
352 
353 
354 
355 template <int dim, typename Number>
356 DEAL_II_CXX20_REQUIRES(dim >= 0)
357 inline DEAL_II_HOST_DEVICE Point<dim, Number>::Point(
358  const Tensor<1, dim, Number> &t)
359  : Tensor<1, dim, Number>(t)
360 {}
361 
362 
363 
364 template <int dim, typename Number>
365 DEAL_II_CXX20_REQUIRES(dim >= 0)
366 inline DEAL_II_HOST_DEVICE Point<dim, Number>::Point(const Number x)
367 {
368  Assert(dim == 1,
369  ExcMessage(
370  "You can only initialize Point<1> objects using the constructor "
371  "that takes only one argument. Point<dim> objects with dim!=1 "
372  "require initialization with the constructor that takes 'dim' "
373  "arguments."));
374 
375  // we can only get here if we pass the assertion. use the switch anyway so
376  // as to avoid compiler warnings about uninitialized elements or writing
377  // beyond the end of the 'values' array
378  switch (dim)
379  {
380  case 1:
381  this->values[0] = x;
382  break;
383 
384  default:;
385  }
386 }
387 
388 
389 
390 template <int dim, typename Number>
391 DEAL_II_CXX20_REQUIRES(dim >= 0)
392 inline DEAL_II_HOST_DEVICE Point<dim, Number>::Point(const Number x,
393  const Number y)
394 {
395  Assert(dim == 2,
396  ExcMessage(
397  "You can only initialize Point<2> objects using the constructor "
398  "that takes two arguments. Point<dim> objects with dim!=2 "
399  "require initialization with the constructor that takes 'dim' "
400  "arguments."));
401 
402  // we can only get here if we pass the assertion. use the indirection anyway
403  // so as to avoid compiler warnings about uninitialized elements or writing
404  // beyond the end of the 'values' array
405  constexpr unsigned int y_index = (dim < 2) ? 0 : 1;
406  this->values[0] = x;
407  this->values[y_index] = y;
408 }
409 
410 
411 
412 template <int dim, typename Number>
413 DEAL_II_CXX20_REQUIRES(dim >= 0)
414 inline DEAL_II_HOST_DEVICE Point<dim, Number>::Point(const Number x,
415  const Number y,
416  const Number z)
417 {
418  Assert(dim == 3,
419  ExcMessage(
420  "You can only initialize Point<3> objects using the constructor "
421  "that takes three arguments. Point<dim> objects with dim!=3 "
422  "require initialization with the constructor that takes 'dim' "
423  "arguments."));
424 
425  // we can only get here if we pass the assertion. use the indirection anyway
426  // so as to avoid compiler warnings about uninitialized elements or writing
427  // beyond the end of the 'values' array
428  constexpr unsigned int y_index = (dim < 2) ? 0 : 1;
429  constexpr unsigned int z_index = (dim < 3) ? 0 : 2;
430  this->values[0] = x;
431  this->values[y_index] = y;
432  this->values[z_index] = z;
433 }
434 
435 
436 
437 template <int dim, typename Number>
438 DEAL_II_CXX20_REQUIRES(dim >= 0)
439 template <std::size_t dummy_dim,
440  std::enable_if_t<(dim == dummy_dim) && (dummy_dim != 0), int>>
441 inline Point<dim, Number>::Point(
442  const boost::geometry::model::
443  point<Number, dummy_dim, boost::geometry::cs::cartesian> &boost_pt)
444 {
445  Assert(dim <= 3, ExcNotImplemented());
446  this->values[0] = boost::geometry::get<0>(boost_pt);
447  constexpr unsigned int y_index = (dim < 2) ? 0 : 1;
448  constexpr unsigned int z_index = (dim < 3) ? 0 : 2;
449 
450  if (dim >= 2)
451  this->values[y_index] = boost::geometry::get<y_index>(boost_pt);
452 
453  if (dim >= 3)
454  this->values[z_index] = boost::geometry::get<z_index>(boost_pt);
455 }
456 
457 
458 
459 template <int dim, typename Number>
460 DEAL_II_CXX20_REQUIRES(dim >= 0)
461 inline DEAL_II_HOST_DEVICE Point<dim, Number> Point<dim, Number>::unit_vector(
462  unsigned int i)
463 {
465  p[i] = 1.;
466  return p;
467 }
468 
469 
470 template <int dim, typename Number>
471 DEAL_II_CXX20_REQUIRES(dim >= 0)
472 inline DEAL_II_HOST_DEVICE Number Point<dim, Number>::operator()(
473  const unsigned int index) const
474 {
475  AssertIndexRange(static_cast<int>(index), dim);
476  return this->values[index];
477 }
478 
479 
480 
481 template <int dim, typename Number>
482 DEAL_II_CXX20_REQUIRES(dim >= 0)
483 inline DEAL_II_HOST_DEVICE Number &Point<dim, Number>::operator()(
484  const unsigned int index)
485 {
486  AssertIndexRange(static_cast<int>(index), dim);
487  return this->values[index];
488 }
489 
490 
491 
492 template <int dim, typename Number>
493 DEAL_II_CXX20_REQUIRES(dim >= 0)
494 template <typename OtherNumber>
495 inline DEAL_II_ALWAYS_INLINE Point<dim, Number> &Point<dim, Number>::operator=(
496  const Tensor<1, dim, OtherNumber> &p)
497 {
499  return *this;
500 }
501 
502 
503 
504 template <int dim, typename Number>
505 DEAL_II_CXX20_REQUIRES(dim >= 0)
506 inline DEAL_II_HOST_DEVICE Point<dim, Number> Point<dim, Number>::operator+(
507  const Tensor<1, dim, Number> &p) const
508 {
509  Point<dim, Number> tmp = *this;
510  tmp += p;
511  return tmp;
512 }
513 
514 
515 
516 template <int dim, typename Number>
517 DEAL_II_CXX20_REQUIRES(dim >= 0)
518 inline DEAL_II_HOST_DEVICE Tensor<1, dim, Number> Point<dim, Number>::operator-(
519  const Point<dim, Number> &p) const
520 {
521  return (Tensor<1, dim, Number>(*this) -= p);
522 }
523 
524 
525 
526 template <int dim, typename Number>
527 DEAL_II_CXX20_REQUIRES(dim >= 0)
528 inline DEAL_II_HOST_DEVICE Point<dim, Number> Point<dim, Number>::operator-(
529  const Tensor<1, dim, Number> &p) const
530 {
531  Point<dim, Number> tmp = *this;
532  tmp -= p;
533  return tmp;
534 }
535 
536 
537 
538 template <int dim, typename Number>
539 DEAL_II_CXX20_REQUIRES(dim >= 0)
540 inline DEAL_II_HOST_DEVICE Point<dim, Number> Point<dim, Number>::operator-()
541  const
542 {
543  Point<dim, Number> result;
544  for (unsigned int i = 0; i < dim; ++i)
545  result.values[i] = -this->values[i];
546  return result;
547 }
548 
549 
550 
551 template <int dim, typename Number>
552 DEAL_II_CXX20_REQUIRES(dim >= 0)
553 template <typename OtherNumber>
555  dim,
556  typename ProductType<Number, typename EnableIfScalar<OtherNumber>::type>::
557  type> Point<dim, Number>::operator*(const OtherNumber factor) const
558 {
560  for (unsigned int i = 0; i < dim; ++i)
561  tmp[i] = this->operator[](i) * factor;
562  return tmp;
563 }
564 
565 
566 
567 template <int dim, typename Number>
568 DEAL_II_CXX20_REQUIRES(dim >= 0)
569 template <typename OtherNumber>
571  dim,
572  typename ProductType<Number, typename EnableIfScalar<OtherNumber>::type>::
573  type> Point<dim, Number>::operator/(const OtherNumber factor) const
574 {
575  const Tensor<1, dim, Number> &base_object = *this;
576  return Point<
577  dim,
578  typename ProductType<Number,
579  typename EnableIfScalar<OtherNumber>::type>::type>(
580  ::operator/(base_object, factor));
581 }
582 
583 
584 
585 template <int dim, typename Number>
586 DEAL_II_CXX20_REQUIRES(dim >= 0)
587 inline DEAL_II_HOST_DEVICE Number Point<dim, Number>::operator*(
588  const Tensor<1, dim, Number> &p) const
589 {
590  Number res = Number();
591  for (unsigned int i = 0; i < dim; ++i)
592  res += this->operator[](i) * p[i];
593  return res;
594 }
595 
596 
597 template <int dim, typename Number>
598 DEAL_II_CXX20_REQUIRES(dim >= 0)
599 inline DEAL_II_HOST_DEVICE typename numbers::NumberTraits<Number>::real_type
600  Point<dim, Number>::square() const
601 {
602  return this->norm_square();
603 }
604 
605 
606 
607 template <int dim, typename Number>
608 DEAL_II_CXX20_REQUIRES(dim >= 0)
609 inline DEAL_II_HOST_DEVICE typename numbers::NumberTraits<Number>::real_type
610  Point<dim, Number>::distance(const Point<dim, Number> &p) const
611 {
612  return std::sqrt(distance_square(p));
613 }
614 
615 
616 
617 template <int dim, typename Number>
618 DEAL_II_CXX20_REQUIRES(dim >= 0)
619 inline DEAL_II_HOST_DEVICE typename numbers::NumberTraits<Number>::real_type
620  Point<dim, Number>::distance_square(const Point<dim, Number> &p) const
621 {
623  for (unsigned int i = 0; i < dim; ++i)
624  {
625  const Number diff = static_cast<Number>(this->values[i]) - p(i);
627  }
628 
629  return sum;
630 }
631 
632 
633 
634 template <int dim, typename Number>
635 DEAL_II_CXX20_REQUIRES(dim >= 0)
636 template <class Archive>
637 inline void Point<dim, Number>::serialize(Archive &ar, const unsigned int)
638 {
639  // forward to serialization
640  // function in the base class
641  ar &static_cast<Tensor<1, dim, Number> &>(*this);
642 }
643 
644 #endif // DOXYGEN
645 
646 
647 /*--------------------------- Global functions: Point -----------------------*/
648 
649 
657 template <int dim, typename Number, typename OtherNumber>
658 inline DEAL_II_HOST_DEVICE
659  Point<dim,
660  typename ProductType<Number,
661  typename EnableIfScalar<OtherNumber>::type>::type>
662  operator*(const OtherNumber factor, const Point<dim, Number> &p)
663 {
664  return p * factor;
665 }
666 
667 
668 
676 template <int dim, typename Number>
677 DEAL_II_CXX20_REQUIRES(dim >= 0)
678 inline std::ostream &
679 operator<<(std::ostream &out, const Point<dim, Number> &p)
680 {
681  for (unsigned int i = 0; i < dim - 1; ++i)
682  out << p[i] << ' ';
683  out << p[dim - 1];
684 
685  return out;
686 }
687 
688 
689 
696 template <int dim, typename Number>
697 DEAL_II_CXX20_REQUIRES(dim >= 0)
698 inline std::istream &
699 operator>>(std::istream &in, Point<dim, Number> &p)
700 {
701  for (unsigned int i = 0; i < dim; ++i)
702  in >> p[i];
703 
704  return in;
705 }
706 
707 
708 #ifndef DOXYGEN
709 
715 template <typename Number>
716 inline std::ostream &
717 operator<<(std::ostream &out, const Point<1, Number> &p)
718 {
719  out << p[0];
720 
721  return out;
722 }
723 
724 #endif // DOXYGEN
726 
727 #endif
OutputOperator< VectorType > & operator<<(OutputOperator< VectorType > &out, unsigned int step)
Definition: operator.h:165
Definition: point.h:112
numbers::NumberTraits< Number >::real_type distance_square(const Point< dim, Number > &p) const
Point< dim, typename ProductType< Number, typename EnableIfScalar< OtherNumber >::type >::type > operator/(const OtherNumber) const
Tensor< 1, dim, Number > operator-(const Point< dim, Number > &) const
Point< dim, typename ProductType< Number, typename EnableIfScalar< OtherNumber >::type >::type > operator*(const OtherNumber) const
Number operator*(const Tensor< 1, dim, Number > &p) const
Number & operator()(const unsigned int index)
numbers::NumberTraits< Number >::real_type distance(const Point< dim, Number > &p) const
Point(const boost::geometry::model::point< Number, dummy_dim, boost::geometry::cs::cartesian > &boost_pt)
Point(const Tensor< 1, dim, Number > &)
Point(const Number x, const Number y, const Number z)
Point< dim, Number > & operator=(const Tensor< 1, dim, OtherNumber > &p)
Point< dim, Number > operator-(const Tensor< 1, dim, Number > &) const
Point(const Number x)
void serialize(Archive &ar, const unsigned int version)
numbers::NumberTraits< Number >::real_type square() const
Point< dim, Number > operator-() const
Point< dim, Number > operator+(const Tensor< 1, dim, Number > &) const
Point(const Number x, const Number y)
static Point< dim, Number > unit_vector(const unsigned int i)
Number operator()(const unsigned int index) const
Point< dim, typename ProductType< Number, typename EnableIfScalar< OtherNumber >::type >::type > operator*(const OtherNumber factor, const Point< dim, Number > &p)
Definition: point.h:662
SymmetricTensor< rank, dim, Number > sum(const SymmetricTensor< rank, dim, Number > &local, const MPI_Comm mpi_communicator)
Definition: tensor.h:516
constexpr Tensor & operator=(const Tensor< rank_, dim, OtherNumber > &rhs)
Tensor< rank_ - 1, dim, Number > values[(dim !=0) ? dim :1]
Definition: tensor.h:884
#define DEAL_II_ALWAYS_INLINE
Definition: config.h:110
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_CXX20_REQUIRES(condition)
Definition: config.h:166
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
#define Assert(cond, exc)
Definition: exceptions.h:1616
static ::ExceptionBase & ExcNotImplemented()
#define AssertIndexRange(index, range)
Definition: exceptions.h:1857
static ::ExceptionBase & ExcMessage(std::string arg1)
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
Definition: utilities.cc:192
#define DEAL_II_HOST_DEVICE
Definition: numbers.h:35
::VectorizedArray< Number, width > sqrt(const ::VectorizedArray< Number, width > &)
static constexpr DEAL_II_HOST_DEVICE_ALWAYS_INLINE const T & value(const T &t)
Definition: numbers.h:703
static constexpr real_type abs_square(const number &x)
Definition: numbers.h:585