 |
Reference documentation for deal.II version 9.2.0
|
\(\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\}}\)
Go to the documentation of this file.
16 #ifndef dealii_point_h
17 #define dealii_point_h
26 #include <boost/geometry/core/cs.hpp>
27 #include <boost/geometry/geometries/point.hpp>
110 template <
int dim,
typename Number =
double>
138 Point(
const Number x);
150 Point(
const Number x,
const Number y);
162 Point(
const Number x,
const Number y,
const Number z);
167 template <std::size_t dummy_dim,
168 typename std::enable_if<(dim == dummy_dim) && (dummy_dim != 0),
170 Point(
const boost::geometry::model::
171 point<Number, dummy_dim, boost::geometry::cs::cartesian> &boost_pt);
204 template <
typename OtherNumber>
268 template <
typename OtherNumber>
280 template <
typename OtherNumber>
336 template <
class Archive>
338 serialize(Archive &ar,
const unsigned int version);
347 template <
int dim,
typename Number>
354 template <
int dim,
typename Number>
357 :
Tensor<1, dim, Number>(t)
362 template <
int dim,
typename Number>
366 # ifndef __CUDA_ARCH__
369 "You can only initialize Point<1> objects using the constructor "
370 "that takes only one argument. Point<dim> objects with dim!=1 "
371 "require initialization with the constructor that takes 'dim' "
390 template <
int dim,
typename Number>
394 # ifndef __CUDA_ARCH__
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' "
406 constexpr
unsigned int y_index = (dim < 2) ? 0 : 1;
408 this->values[y_index] = y;
413 template <
int dim,
typename Number>
417 # ifndef __CUDA_ARCH__
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' "
429 constexpr
unsigned int y_index = (dim < 2) ? 0 : 1;
430 constexpr
unsigned int z_index = (dim < 3) ? 0 : 2;
432 this->values[y_index] = y;
433 this->values[z_index] = z;
438 template <
int dim,
typename Number>
440 std::size_t dummy_dim,
441 typename std::enable_if<(dim == dummy_dim) && (dummy_dim != 0),
int>::type>
443 const boost::geometry::model::
444 point<Number, dummy_dim, boost::geometry::cs::cartesian> &boost_pt)
447 this->values[0] = boost::geometry::get<0>(boost_pt);
448 constexpr
unsigned int y_index = (dim < 2) ? 0 : 1;
449 constexpr
unsigned int z_index = (dim < 3) ? 0 : 2;
452 this->values[y_index] = boost::geometry::get<y_index>(boost_pt);
455 this->values[z_index] = boost::geometry::get<z_index>(boost_pt);
460 template <
int dim,
typename Number>
470 template <
int dim,
typename Number>
474 # ifndef __CUDA_ARCH__
477 return this->values[index];
482 template <
int dim,
typename Number>
486 # ifndef __CUDA_ARCH__
489 return this->values[index];
494 template <
int dim,
typename Number>
495 template <
typename OtherNumber>
505 template <
int dim,
typename Number>
516 template <
int dim,
typename Number>
525 template <
int dim,
typename Number>
536 template <
int dim,
typename Number>
541 for (
unsigned int i = 0; i < dim; ++i)
542 result.
values[i] = -this->values[i];
548 template <
int dim,
typename Number>
549 template <
typename OtherNumber>
557 for (
unsigned int i = 0; i < dim; ++i)
558 tmp[i] = this->
operator[](i) * factor;
564 template <
int dim,
typename Number>
565 template <
typename OtherNumber>
577 ::operator/(base_object, factor));
582 template <
int dim,
typename Number>
586 Number res = Number();
587 for (
unsigned int i = 0; i < dim; ++i)
588 res += this->
operator[](i) * p[i];
593 template <
int dim,
typename Number>
597 return this->norm_square();
602 template <
int dim,
typename Number>
611 template <
int dim,
typename Number>
616 for (
unsigned int i = 0; i < dim; ++i)
618 const Number diff =
static_cast<Number
>(this->values[i]) - p(i);
627 template <
int dim,
typename Number>
628 template <
class Archive>
651 template <
int dim,
typename Number,
typename OtherNumber>
668 template <
int dim,
typename Number>
669 inline std::ostream &
672 for (
unsigned int i = 0; i < dim - 1; ++i)
685 template <
int dim,
typename Number>
686 inline std::istream &
689 for (
unsigned int i = 0; i < dim; ++i)
703 template <
typename Number>
704 inline std::ostream &
numbers::NumberTraits< Number >::real_type distance_square(const Point< dim, Number > &p) const
static ::ExceptionBase & ExcNotImplemented()
static Point< dim, Number > unit_vector(const unsigned int i)
static constexpr const T & value(const T &t)
Point< dim, Number > & operator=(const Tensor< 1, dim, OtherNumber > &p)
#define AssertIndexRange(index, range)
static constexpr std::enable_if< std::is_same< Dummy, number >::value &&is_cuda_compatible< Dummy >::value, real_type >::type abs_square(const number &x)
Point< dim, Number > operator+(const Tensor< 1, dim, Number > &) const
#define DEAL_II_ALWAYS_INLINE
void serialize(Archive &ar, const unsigned int version)
#define DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
Point< dim, typename ProductType< Number, typename EnableIfScalar< OtherNumber >::type >::type > operator*(const OtherNumber) const
Number operator()(const unsigned int index) const
static ::ExceptionBase & ExcMessage(std::string arg1)
Point< dim, Number > operator-() const
#define DEAL_II_NAMESPACE_OPEN
SymmetricTensor< rank, dim, Number > sum(const SymmetricTensor< rank, dim, Number > &local, const MPI_Comm &mpi_communicator)
numbers::NumberTraits< Number >::real_type distance(const Point< dim, Number > &p) const
Point< dim, typename ProductType< Number, typename EnableIfScalar< OtherNumber >::type >::type > operator/(const OtherNumber) const
inline ::VectorizedArray< Number, width > sqrt(const ::VectorizedArray< Number, width > &x)
#define Assert(cond, exc)
OutputOperator< VectorType > & operator<<(OutputOperator< VectorType > &out, unsigned int step)
constexpr Tensor & operator=(const Tensor< rank_, dim, OtherNumber > &rhs)
numbers::NumberTraits< Number >::real_type square() const
#define DEAL_II_CUDA_HOST_DEV
#define DEAL_II_NAMESPACE_CLOSE
std::istream & operator>>(std::istream &in, Point< dim, Number > &p)
Tensor< rank_ - 1, dim, double > values[(dim !=0) ? dim :1]
#define DEAL_II_DISABLE_EXTRA_DIAGNOSTICS