Reference documentation for deal.II version Git 3f1f337db3 2021-10-23 13:19:02 -0600
\(\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 - 2021 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 
26 #include <boost/geometry/core/cs.hpp>
27 #include <boost/geometry/geometries/point.hpp>
29 
30 #include <cmath>
31 
33 
109 template <int dim, typename Number = double>
110 class Point : public Tensor<1, dim, Number>
111 {
112 public:
120  Point();
121 
125  explicit DEAL_II_CUDA_HOST_DEV
126  Point(const Tensor<1, dim, Number> &);
127 
136  explicit DEAL_II_CUDA_HOST_DEV
137  Point(const Number x);
138 
149  Point(const Number x, const Number y);
150 
161  Point(const Number x, const Number y, const Number z);
162 
166  template <std::size_t dummy_dim,
167  typename std::enable_if<(dim == dummy_dim) && (dummy_dim != 0),
168  int>::type = 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_CUDA_HOST_DEV Number
188  operator()(const unsigned int index) const;
189 
195  DEAL_II_CUDA_HOST_DEV Number &
196  operator()(const unsigned int index);
197 
203  template <typename OtherNumber>
206 
218  operator+(const Tensor<1, dim, Number> &) const;
219 
230  operator-(const Point<dim, Number> &) const;
231 
241  operator-(const Tensor<1, dim, Number> &) const;
242 
249  operator-() const;
250 
267  template <typename OtherNumber>
268  DEAL_II_CUDA_HOST_DEV Point<
269  dim,
270  typename ProductType<Number,
271  typename EnableIfScalar<OtherNumber>::type>::type>
272  operator*(const OtherNumber) const;
273 
279  template <typename OtherNumber>
280  DEAL_II_CUDA_HOST_DEV Point<
281  dim,
282  typename ProductType<Number,
283  typename EnableIfScalar<OtherNumber>::type>::type>
284  operator/(const OtherNumber) const;
285 
291  DEAL_II_CUDA_HOST_DEV Number
292  operator*(const Tensor<1, dim, Number> &p) const;
293 
307  square() const;
308 
317  distance(const Point<dim, Number> &p) const;
318 
326  distance_square(const Point<dim, Number> &p) const;
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>
350 Point<dim, Number>::Point() // NOLINT
351 {}
352 
353 
354 
355 template <int dim, typename Number>
359 {}
360 
361 
362 
363 template <int dim, typename Number>
365 Point<dim, Number>::Point(const Number x)
366 {
367 # ifndef __CUDA_ARCH__
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 # endif
375 
376  // we can only get here if we pass the assertion. use the switch anyway so
377  // as to avoid compiler warnings about uninitialized elements or writing
378  // beyond the end of the 'values' array
379  switch (dim)
380  {
381  case 1:
382  this->values[0] = x;
383  break;
384 
385  default:;
386  }
387 }
388 
389 
390 
391 template <int dim, typename Number>
393 Point<dim, Number>::Point(const Number x, const Number y)
394 {
395 # ifndef __CUDA_ARCH__
396  Assert(dim == 2,
397  ExcMessage(
398  "You can only initialize Point<2> objects using the constructor "
399  "that takes two arguments. Point<dim> objects with dim!=2 "
400  "require initialization with the constructor that takes 'dim' "
401  "arguments."));
402 # endif
403 
404  // we can only get here if we pass the assertion. use the indirection anyway
405  // so as to avoid compiler warnings about uninitialized elements or writing
406  // beyond the end of the 'values' array
407  constexpr unsigned int y_index = (dim < 2) ? 0 : 1;
408  this->values[0] = x;
409  this->values[y_index] = y;
410 }
411 
412 
413 
414 template <int dim, typename Number>
416 Point<dim, Number>::Point(const Number x, const Number y, const Number z)
417 {
418 # ifndef __CUDA_ARCH__
419  Assert(dim == 3,
420  ExcMessage(
421  "You can only initialize Point<3> objects using the constructor "
422  "that takes three arguments. Point<dim> objects with dim!=3 "
423  "require initialization with the constructor that takes 'dim' "
424  "arguments."));
425 # endif
426 
427  // we can only get here if we pass the assertion. use the indirection anyway
428  // so as to avoid compiler warnings about uninitialized elements or writing
429  // beyond the end of the 'values' array
430  constexpr unsigned int y_index = (dim < 2) ? 0 : 1;
431  constexpr unsigned int z_index = (dim < 3) ? 0 : 2;
432  this->values[0] = x;
433  this->values[y_index] = y;
434  this->values[z_index] = z;
435 }
436 
437 
438 
439 template <int dim, typename Number>
440 template <
441  std::size_t dummy_dim,
442  typename std::enable_if<(dim == dummy_dim) && (dummy_dim != 0), int>::type>
444  const boost::geometry::model::
445  point<Number, dummy_dim, boost::geometry::cs::cartesian> &boost_pt)
446 {
447  Assert(dim <= 3, ExcNotImplemented());
448  this->values[0] = boost::geometry::get<0>(boost_pt);
449  constexpr unsigned int y_index = (dim < 2) ? 0 : 1;
450  constexpr unsigned int z_index = (dim < 3) ? 0 : 2;
451 
452  if (dim >= 2)
453  this->values[y_index] = boost::geometry::get<y_index>(boost_pt);
454 
455  if (dim >= 3)
456  this->values[z_index] = boost::geometry::get<z_index>(boost_pt);
457 }
458 
459 
460 
461 template <int dim, typename Number>
463  Point<dim, Number>::unit_vector(unsigned int i)
464 {
466  p[i] = 1.;
467  return p;
468 }
469 
470 
471 template <int dim, typename Number>
472 inline DEAL_II_CUDA_HOST_DEV Number
473 Point<dim, Number>::operator()(const unsigned int index) const
474 {
475 # ifndef __CUDA_ARCH__
476  AssertIndexRange(static_cast<int>(index), dim);
477 # endif
478  return this->values[index];
479 }
480 
481 
482 
483 template <int dim, typename Number>
484 inline DEAL_II_CUDA_HOST_DEV Number &
485 Point<dim, Number>::operator()(const unsigned int index)
486 {
487 # ifndef __CUDA_ARCH__
488  AssertIndexRange(static_cast<int>(index), dim);
489 # endif
490  return this->values[index];
491 }
492 
493 
494 
495 template <int dim, typename Number>
496 template <typename OtherNumber>
499 {
501  return *this;
502 }
503 
504 
505 
506 template <int dim, typename Number>
509 {
510  Point<dim, Number> tmp = *this;
511  tmp += p;
512  return tmp;
513 }
514 
515 
516 
517 template <int dim, typename Number>
520 {
521  return (Tensor<1, dim, Number>(*this) -= p);
522 }
523 
524 
525 
526 template <int dim, typename Number>
529 {
530  Point<dim, Number> tmp = *this;
531  tmp -= p;
532  return tmp;
533 }
534 
535 
536 
537 template <int dim, typename Number>
540 {
541  Point<dim, Number> result;
542  for (unsigned int i = 0; i < dim; ++i)
543  result.values[i] = -this->values[i];
544  return result;
545 }
546 
547 
548 
549 template <int dim, typename Number>
550 template <typename OtherNumber>
552  Point<dim,
553  typename ProductType<Number,
554  typename EnableIfScalar<OtherNumber>::type>::type>
555  Point<dim, Number>::operator*(const OtherNumber factor) const
556 {
558  for (unsigned int i = 0; i < dim; ++i)
559  tmp[i] = this->operator[](i) * factor;
560  return tmp;
561 }
562 
563 
564 
565 template <int dim, typename Number>
566 template <typename OtherNumber>
568  Point<dim,
569  typename ProductType<Number,
570  typename EnableIfScalar<OtherNumber>::type>::type>
571  Point<dim, Number>::operator/(const OtherNumber factor) const
572 {
573  const Tensor<1, dim, Number> &base_object = *this;
574  return Point<
575  dim,
576  typename ProductType<Number,
577  typename EnableIfScalar<OtherNumber>::type>::type>(
578  ::operator/(base_object, factor));
579 }
580 
581 
582 
583 template <int dim, typename Number>
584 inline DEAL_II_CUDA_HOST_DEV Number
586 {
587  Number res = Number();
588  for (unsigned int i = 0; i < dim; ++i)
589  res += this->operator[](i) * p[i];
590  return res;
591 }
592 
593 
594 template <int dim, typename Number>
597 {
598  return this->norm_square();
599 }
600 
601 
602 
603 template <int dim, typename Number>
606 {
607  return std::sqrt(distance_square(p));
608 }
609 
610 
611 
612 template <int dim, typename Number>
615 {
617  for (unsigned int i = 0; i < dim; ++i)
618  {
619  const Number diff = static_cast<Number>(this->values[i]) - p(i);
621  }
622 
623  return sum;
624 }
625 
626 
627 
628 template <int dim, typename Number>
629 template <class Archive>
630 inline void
631 Point<dim, Number>::serialize(Archive &ar, const unsigned int)
632 {
633  // forward to serialization
634  // function in the base class
635  ar &static_cast<Tensor<1, dim, Number> &>(*this);
636 }
637 
638 #endif // DOXYGEN
639 
640 
641 /*--------------------------- Global functions: Point -----------------------*/
642 
643 
652 template <int dim, typename Number, typename OtherNumber>
654  Point<dim,
655  typename ProductType<Number,
656  typename EnableIfScalar<OtherNumber>::type>::type>
657  operator*(const OtherNumber factor, const Point<dim, Number> &p)
658 {
659  return p * factor;
660 }
661 
662 
663 
669 template <int dim, typename Number>
670 inline std::ostream &
671 operator<<(std::ostream &out, const Point<dim, Number> &p)
672 {
673  for (unsigned int i = 0; i < dim - 1; ++i)
674  out << p[i] << ' ';
675  out << p[dim - 1];
676 
677  return out;
678 }
679 
680 
681 
686 template <int dim, typename Number>
687 inline std::istream &
688 operator>>(std::istream &in, Point<dim, Number> &p)
689 {
690  for (unsigned int i = 0; i < dim; ++i)
691  in >> p[i];
692 
693  return in;
694 }
695 
696 
697 #ifndef DOXYGEN
698 
704 template <typename Number>
705 inline std::ostream &
706 operator<<(std::ostream &out, const Point<1, Number> &p)
707 {
708  out << p[0];
709 
710  return out;
711 }
712 
713 #endif // DOXYGEN
715 
716 #endif
Tensor< rank, dim, Number > sum(const Tensor< rank, dim, Number > &local, const MPI_Comm &mpi_communicator)
Point< dim, Number > operator+(const Tensor< 1, dim, Number > &) const
static constexpr const T & value(const T &t)
Definition: numbers.h:703
Point< dim, Number > & operator=(const Tensor< 1, dim, OtherNumber > &p)
static constexpr std::enable_if< std::is_same< Dummy, number >::value &&is_cuda_compatible< Dummy >::value, real_type >::type abs_square(const number &x)
Definition: numbers.h:587
#define AssertIndexRange(index, range)
Definition: exceptions.h:1720
static Point< dim, Number > unit_vector(const unsigned int i)
Tensor< rank_ - 1, dim, Number > values[(dim !=0) ? dim :1]
Definition: tensor.h:869
Number operator()(const unsigned int index) const
Definition: point.h:110
constexpr Tensor & operator=(const Tensor< rank_, dim, OtherNumber > &rhs)
#define DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
Definition: config.h:414
static ::ExceptionBase & ExcMessage(std::string arg1)
Point< dim, Number > operator-() const
#define Assert(cond, exc)
Definition: exceptions.h:1461
Point< dim, typename ProductType< Number, typename EnableIfScalar< OtherNumber >::type >::type > operator*(const OtherNumber) const
void serialize(Archive &ar, const unsigned int version)
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
#define DEAL_II_ALWAYS_INLINE
Definition: config.h:98
constexpr numbers::NumberTraits< Number >::real_type norm_square() const
numbers::NumberTraits< Number >::real_type distance(const Point< dim, Number > &p) const
Definition: tensor.h:506
#define DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
Definition: config.h:452
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:400
numbers::NumberTraits< Number >::real_type square() const
Point< dim, typename ProductType< Number, typename EnableIfScalar< OtherNumber >::type >::type > operator/(const OtherNumber) const
#define DEAL_II_CUDA_HOST_DEV
Definition: numbers.h:34
static ::ExceptionBase & ExcNotImplemented()
numbers::NumberTraits< Number >::real_type distance_square(const Point< dim, Number > &p) const
std::istream & operator>>(std::istream &in, Point< dim, Number > &p)
Definition: point.h:688