v0.8.16
DataStructures.hpp
Go to the documentation of this file.
1 /** \file DataStructures.hpp
2 
3 \brief Data structures for accessing information about finite element and its
4 degrees of freedom.
5 
6 */
7 
8 /* This file is part of MoFEM.
9  * MoFEM is free software: you can redistribute it and/or modify it under
10  * the terms of the GNU Lesser General Public License as published by the
11  * Free Software Foundation, either version 3 of the License, or (at your
12  * option) any later version.
13  *
14  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
15  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
16  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
17  * License for more details.
18  *
19  * You should have received a copy of the GNU Lesser General Public
20  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
21 
22 #ifndef __DATASTRUCTURES_HPP
23 #define __DATASTRUCTURES_HPP
24 
25 using namespace boost::numeric;
26 
27 namespace MoFEM {
28 
29 typedef ublas::unbounded_array<
30  boost::shared_ptr<const FEDofEntity>,
31  std::allocator<boost::shared_ptr<const FEDofEntity> > >
33 typedef ublas::vector<boost::shared_ptr<const FEDofEntity>, DofsAllocator>
35 
36 /**
37 * \brief Get tensor rank 0 (scalar) form data vector
38 * \ingroup mofem_forces_and_sources_user_data_operators
39 
40 Example how to use it.
41 \code
42 VectorDouble vec;
43 vec.resize(nb_gauss_pts,false);
44 vec.clear();
45 auto t0 = getFTensor0FromData(data);
46 for(int gg = 0;gg!=nb_gauss_pts;gg++) {
47 
48  ++t0;
49 }
50 \endcode
51 
52 */
53 template <class T, class A>
55 getFTensor0FromVec(ublas::vector<T, A> &data) {
56  static_assert(1, "not implemented");
57 }
58 
59 /**
60  * @deprecated Name change to getFTensor0FromVec
61  */
62 template <class T, class A>
64 getTensor0FormData(ublas::vector<T, A> &data) {
65  return getFTensor0FromVec(data);
66 }
67 
68 template <>
71  ublas::vector<double, DoubleAllocator> &data);
72 
73 /**
74  * \brief Get tensor rank 1 (vector) form data matrix
75  * \ingroup mofem_forces_and_sources_user_data_operators
76  */
77 template <int Tensor_Dim, class T, class L, class A>
79 getFTensor1FromMat(ublas::matrix<T, L, A> &data) {
80  static_assert(1, "not implemented");
81 }
82 
83 /**
84  * \brief Get tensor rank 1 (vector) form data matrix (specialization)
85  * \ingroup mofem_forces_and_sources_user_data_operators
86  */
87 template <int Tensor_Dim>
90  return getFTensor1FromMat<Tensor_Dim, double, ublas::row_major,
91  DoubleAllocator>(data);
92 }
93 
94 template <>
96 getFTensor1FromMat<3, double, ublas::row_major, DoubleAllocator>(
97  MatrixDouble &data);
98 
99 template <>
101 getFTensor1FromMat<2, double, ublas::row_major, DoubleAllocator>(
102  MatrixDouble &data);
103 
104 /**
105  * @deprecated Name change to getFTensor1FromMat
106  */
107 template <int Tensor_Dim>
110  return getFTensor1FromMat<Tensor_Dim>(data);
111 }
112 
113 /**
114  * \brief Get tensor rank 2 (matrix) form data matrix
115  * \ingroup mofem_forces_and_sources_user_data_operators
116  */
117 template <int Tensor_Dim0, int Tensor_Dim1, class T, class L, class A>
118 FTensor::Tensor2<FTensor::PackPtr<T *, 1>, Tensor_Dim0, Tensor_Dim1>
119 getFTensor2FromMat(ublas::matrix<T, L, A> &data) {
120  static_assert(1, "not implemented");
121 }
122 
123 /**
124  * Template specialization for getFTensor2FromMat
125  * \ingroup mofem_forces_and_sources_user_data_operators
126  *
127  */
128 template <>
131 
132 /**
133  * Template specialization for getFTensor2FromMat
134  */
135 template <>
138 
139 /**
140  * \brief Get tensor rank 2 (matrix) form data matrix (specialization)
141  * \ingroup mofem_forces_and_sources_user_data_operators
142  */
143 template <int Tensor_Dim0, int Tensor_Dim1>
144 FTensor::Tensor2<FTensor::PackPtr<double *, 1>, Tensor_Dim0, Tensor_Dim1>
146  return getFTensor2FromMat<Tensor_Dim0, Tensor_Dim1, double, ublas::row_major,
147  DoubleAllocator>(data);
148 }
149 
150 /**
151  * @deprecated Name change to getFTensor1FromMat
152  */
153 template <int Tensor_Dim0, int Tensor_Dim1>
156  return getFTensor2FromMat<Tensor_Dim0, Tensor_Dim1>(data);
157 }
158 
159 /**
160  * \brief Get symmetric tensor rank 2 (matrix) form data matrix
161  * \ingroup mofem_forces_and_sources_user_data_operators
162  */
163 template <int Tensor_Dim, class T, class L, class A>
165 getFTensor2SymmetricFromMat(ublas::matrix<T, L, A> &data) {
166  static_assert(1, "not implemented");
167 }
168 
169 /**
170  * @brief Get symmetric tensor rank 2 form matrix of for dimension 3
171  *
172  * Specialisation for symmetric tensor 2
173  *
174  * @tparam
175  * @param data
176  * @return FTensor::Tensor2_symmetric<FTensor::PackPtr<double *, 1>, 3>
177  */
178 template <>
181 
182 /**
183  * @brief Get symmetric tensor rank 2 form matrix
184  *
185  * Specialisation for symmetric tensor 2
186  *
187  * @tparam Tensor_Dim
188  * @param data
189  * @return FTensor::Tensor2_symmetric<FTensor::PackPtr<double *, 1>, Tensor_Dim>
190  */
191 template <int Tensor_Dim>
195  DoubleAllocator>(data);
196 }
197 
198 /** \brief data structure for finite element entity
199  * \ingroup mofem_forces_and_sources_user_data_operators
200  *
201  * It keeps that about indices of degrees of freedom, dofs data, base functions
202  * functions, entity side number, type of entities, approximation order, etc.
203  *
204  */
206 
207  /** \brief Data on single entity (This is passed as argument to
208  * DataOperator::doWork) \ingroup mofem_forces_and_sources_user_data_operators
209  * \nosubgrouping
210  *
211  * \todo Hdiv and Hcurl functions should be accessed through common interface.
212  */
213  struct EntData {
214 
215  /** \name Constructor and destructor */
216 
217  /**@{*/
218 
219  EntData();
220  virtual ~EntData() {}
221 
222  /**@}*/
223 
224  /** \name Sense, order and indices */
225 
226  /**@{*/
227 
228  /// \brief get entity sense, need to calculate base functions with
229  /// conforming approximation fields
230  virtual int getSense() const { return sEnse; }
231 
232  /// \brief get approximation order
233  inline ApproximationOrder getOrder() const { return oRder; }
234 
235  /// \brief Get global indices of dofs on entity
236  inline const VectorInt &getIndices() const { return iNdices; }
237 
238  /// \brief get global indices of dofs on entity up to given order
239  inline const VectorIntAdaptor getIndicesUpToOrder(int order) {
240  unsigned int size = 0;
241  if (iNdices.size()) {
242  size = dOfs[0]->getOrderNbDofs(order) * dOfs[0]->getNbOfCoeffs();
243  size = size < iNdices.size() ? size : iNdices.size();
244  }
245  int *data = &*iNdices.data().begin();
246  return VectorIntAdaptor(size,
247  ublas::shallow_array_adaptor<int>(size, data));
248  }
249 
250  /// \brief get local indices of dofs on entity
251  inline const VectorInt &getLocalIndices() const { return localIndices; }
252 
253  /// \brief get local indices of dofs on entity up to given order
254  inline const VectorIntAdaptor getLocalIndicesUpToOrder(int order) {
255  unsigned int size = 0;
256  if (localIndices.size()) {
257  size = dOfs[0]->getOrderNbDofs(order) * dOfs[0]->getNbOfCoeffs();
258  size = size < localIndices.size() ? size : localIndices.size();
259  }
260  int *data = &*localIndices.data().begin();
261  return VectorIntAdaptor(size,
262  ublas::shallow_array_adaptor<int>(size, data));
263  }
264 
265  inline int &getSense() { return sEnse; }
266  inline ApproximationOrder &getDataOrder() { return oRder; }
267  inline VectorInt &getIndices() { return iNdices; }
268  inline VectorInt &getLocalIndices() { return localIndices; }
269 
270  /**@}*/
271 
272  /** \name Data on entity */
273 
274  /**@{*/
275 
276  /// \brief get dofs values
277  inline const VectorDouble &getFieldData() const { return fieldData; }
278 
279  /// \brief get dofs values up to given order
280  inline const VectorAdaptor getFieldDataUpToOrder(int order) {
281  unsigned int size = 0;
282  if (fieldData.size()) {
283  size = dOfs[0]->getOrderNbDofs(order) * dOfs[0]->getNbOfCoeffs();
284  size = size < fieldData.size() ? size : fieldData.size();
285  }
286  double *data = &*fieldData.data().begin();
287  return getVectorAdaptor(data, size);
288  }
289 
290  /// \brief get dofs data stature FEDofEntity
291  inline const VectorDofs &getFieldDofs() const { return dOfs; }
292 
293  inline VectorDouble &getFieldData() { return fieldData; }
294 
295  /**
296  * @brief Return FTensor of rank 1, i.e. vector from filed data coeffinects
297  *
298  * \code
299  * auto t_vec = data.getFTensor1FieldData<3>();
300  * \endcode
301  *
302  * @tparam Tensor_Dim size of vector
303  * @return FTensor::Tensor1<FTensor::PackPtr<double *, Tensor_Dim>,
304  * Tensor_Dim>
305  */
306  template <int Tensor_Dim>
309  std::stringstream s;
310  s << "Not implemented for this dimension dim = " << Tensor_Dim;
311  THROW_MESSAGE(s.str());
312  }
313 
314  /**
315  * @brief Return FTensor rank 2, i.e. matrix from filed data coeffinects
316  *
317  * \code
318  * auto t_mat = data.getFTensor2FieldData<3,3>();
319  * \endcode
320  *
321  * @tparam Tensor_Dim0
322  * @tparam Tensor_Dim1
323  * @return FTensor::Tensor2<FTensor::PackPtr<double *, Tensor_Dim0 * Tensor_Dim1>, Tensor_Dim0, Tensor_Dim1>
324  */
325  template <int Tensor_Dim0, int Tensor_Dim1>
327  Tensor_Dim0, Tensor_Dim1>
329  std::stringstream s;
330  s << "Not implemented for this dimension dim0 = " << Tensor_Dim0;
331  s << " and dim1 " << Tensor_Dim1;
332  THROW_MESSAGE(s.str());
333  }
334 
335  /**
336  * @brief Return symmetric FTensor rank 2, i.e. matrix from filed data coeffinects
337  *
338  * \code
339  * auto t_mat = data.getFTensor2SymmetricFieldData<3>();
340  * \endcode
341  *
342  * @tparam Tensor_Dim dimension of the tensor
343  * @return FTensor::Tensor2_symmetric<FTensor::PackPtr<double *, (Tensor_Dim * (Tensor_Dim + 1)) / 2>, Tensor_Dim>
344  */
345  template <int Tensor_Dim>
347  FTensor::PackPtr<double *, (Tensor_Dim * (Tensor_Dim + 1)) / 2>,
348  Tensor_Dim>
350  std::stringstream s;
351  s << "Not implemented for this dimension dim = " << Tensor_Dim;
352  THROW_MESSAGE(s.str());
353  }
354 
355  /**
356  * @brief Resturn scalar files as a FTensor of rank 0
357  *
358  * @return FTensor::Tensor0<FTensor::PackPtr<double *,1> >
359  */
360  FTensor::Tensor0<FTensor::PackPtr<double *,1> > getFTensor0FieldData();
361 
362  inline VectorDofs &getFieldDofs() { return dOfs; }
363 
364  /**@}*/
365 
366  /** \name Base and space */
367 
368  /**@{*/
369 
370  /**
371  * \brief Get approximation base
372  * @return Approximation base
373  */
374  inline FieldApproximationBase &getBase() { return bAse; }
375 
376  /**
377  * \brief Get field space
378  * @return Field space
379  */
380  inline FieldSpace &getSpace() { return sPace; }
381 
382  /**
383  * Get shared pointer to base base functions
384  */
385  virtual boost::shared_ptr<MatrixDouble> &
387  return N[base];
388  }
389 
390  /**
391  * Get shared pointer to base base functions
392  */
393  virtual const boost::shared_ptr<MatrixDouble> &
395  return N[base];
396  }
397 
398  /**
399  * Get shared pointer to derivatives of base base functions
400  */
401  virtual boost::shared_ptr<MatrixDouble> &
403  return diffN[base];
404  }
405 
406  /**
407  * Get shared pointer to derivatives of base base functions
408  */
409  virtual const boost::shared_ptr<MatrixDouble> &
411  return diffN[base];
412  }
413 
414  /**@}*/
415 
416  /** \name Get base functions for H1/L2 */
417 
418  /**@{*/
419 
420  /** \brief get base functions
421  * this return matrix (nb. of rows is equal to nb. of Gauss pts, nb. of
422  * columns is equal to number of base functions on this entity.
423  *
424  * Note that for vectorial base, like Hdiv or Hcurl, in columns are
425  * vectorial base functions. For tonsorial would be tonsorial base
426  * functions. Interpretation depends on type of base, scalar, vectorial or
427  * tonsorial and dimension fo problem.
428  *
429  */
430  virtual const MatrixDouble &getN(const FieldApproximationBase base) const {
431  return *(getNSharedPtr(base));
432  }
433 
434  inline const MatrixDouble &getN() const { return getN(bAse); }
435 
436  /** \brief get derivatives of base functions
437  *
438  * Matrix at rows has nb. of Gauss pts, at columns it has derivative of
439  * base functions. Columns are structured as follows, [ dN1/dx, dN1/dy,
440  * dN1/dz, dN2/dx, dN2/dy, dN2/dz, ... ]
441  *
442  * Scalar base functions:
443  * Note that base functions are calculated in file H1.c
444  * Above description not apply for derivatives of nodal functions, since
445  * derivative of nodal functions in case of simplexes, EDGES, TRIANGLES and
446  * TETS are constant. So that matrix rows represents nb. of base
447  * functions, columns are derivatives. Nb. of columns depend on element
448  * dimension, for EDGES is one, for TRIS is 2 and TETS is 3.
449  *
450  * Note that for node element this function make no sense.
451  *
452  * Tonsorial base functions:
453  * Note: In rows ale integration pts, columns are formatted that that
454  * components of vectors and then derivatives, for example row for given
455  * integration points is formatted in array
456  * \f[
457  * t_{0,0}, t_{1,0}, t_{1,0}, t_{0,1}, t_{1,1}, t_{1,1}, t_{0,2}, t_{1,2},
458  * t_{1,2} \f] where comma express derivative, i.e. \f$t_{2,1} =
459  * \frac{\partial t_2}{\partial \xi_1}\f$
460  *
461  */
462  virtual const MatrixDouble &
463  getDiffN(const FieldApproximationBase base) const {
464  return *(getDiffNSharedPtr(base));
465  }
466 
467  inline const MatrixDouble &getDiffN() const { return getDiffN(bAse); }
468 
469  /**
470  * \brief Get base functions
471  * @param base Approximation base
472  * @return Error code
473  */
475  return *(getNSharedPtr(base));
476  }
477 
478  /**
479  * \brief Get base functions
480  *
481  * It assumed that approximation base for given field is known and stored in
482  * this data structure
483  *
484  * @return Error code
485  */
486  inline MatrixDouble &getN() { return getN(bAse); }
487 
488  /**
489  * \brief Get derivatives of base functions
490  * @param base Approximation base
491  * @return Error code
492  */
494  return *(getDiffNSharedPtr(base));
495  }
496 
497  /**
498  * \brief Get derivatives of base functions
499  *
500  * It assumed that approximation base for given field is known and stored in
501  * this data structure
502  *
503  * @return Error code
504  */
505  inline MatrixDouble &getDiffN() { return getDiffN(bAse); }
506 
507  /// \brief get base functions at Gauss pts
508  inline const VectorAdaptor getN(const FieldApproximationBase base,
509  const int gg) {
510  int size = getN(base).size2();
511  double *data = &getN(base)(gg, 0);
512  return VectorAdaptor(size,
513  ublas::shallow_array_adaptor<double>(size, data));
514  }
515 
516  /// \brief get base functions at Gauss pts
517  inline const VectorAdaptor getN(const int gg) { return getN(bAse, gg); }
518 
519  /** \brief get derivative of base functions at Gauss pts
520 
521  * returned matrix on rows has base functions, in column its derivatives.
522  *
523  * \param base Approximation base
524  * \param gg Nb. of Gauss pts.
525  *
526  */
528  const int gg) {
529  // FIXME: That is bug, it will not work if number of integration pts is
530  // equal to number of nodes on entity. User who not implementing low level
531  // DataOperator will not experience this.
532  if (getN(base).size1() == getDiffN(base).size1()) {
533  int size = getN(base).size2();
534  int dim = getDiffN(base).size2() / size;
535  double *data = &getDiffN(base)(gg, 0);
536  return MatrixAdaptor(
537  getN(base).size2(), dim,
538  ublas::shallow_array_adaptor<double>(getDiffN(base).size2(), data));
539  } else {
540  // in some cases, f.e. for derivatives of nodal base functions at only
541  // one gauss point is needed
542  return MatrixAdaptor(
543  getN(base).size1(), getN(base).size2(),
544  ublas::shallow_array_adaptor<double>(getDiffN(base).data().size(),
545  &getDiffN(base).data()[0]));
546  }
547  }
548 
549  /** \brief get derivative of base functions at Gauss pts
550 
551  * returned matrix on rows has base functions, in column its derivatives.
552  *
553  * \param gg nb. of Gauss pts.
554  *
555  */
556  inline const MatrixAdaptor getDiffN(const int gg) {
557  return getDiffN(bAse, gg);
558  }
559 
560  /** \brief get base functions at Gauss pts
561 
562  * Note that multi field element, two different field can have different
563  * approximation orders. Since we use hierarchical approximation basis,
564  * base functions are calculated once for element, using maximal
565  * approximation order on given entity.
566  *
567  * \param base Approximation base
568  * \param gg number of Gauss point
569  * \param nb_base_functions number of of base functions returned
570 
571  */
572  inline const VectorAdaptor getN(const FieldApproximationBase base,
573  const int gg, const int nb_base_functions) {
574  (void)getN()(gg, nb_base_functions -
575  1); // throw error if nb_base_functions is to big
576  double *data = &getN(base)(gg, 0);
577  return VectorAdaptor(
578  nb_base_functions,
579  ublas::shallow_array_adaptor<double>(nb_base_functions, data));
580  }
581 
582  /** \brief get base functions at Gauss pts
583 
584  * Note that multi field element, two different field can have different
585  * approximation orders. Since we use hierarchical approximation basis,
586  * base functions are calculated once for element, using maximal
587  * approximation order on given entity.
588  *
589  * \param gg number of Gauss point
590  * \param nb_base_functions number of of base functions returned
591 
592  */
593  inline const VectorAdaptor getN(const int gg, const int nb_base_functions) {
594  return getN(bAse, gg, nb_base_functions);
595  }
596 
597  /** \brief get derivatives of base functions at Gauss pts
598  *
599  * Note that multi field element, two different field can have different
600  * approximation orders. Since we use hierarchical approximation basis,
601  * base functions are calculated once for element, using maximal
602  * approximation order on given entity.
603  *
604  * \param base Approximation base
605  * \param gg nb. of Gauss point
606  * \param nb_base_functions number of of base functions
607  *
608  */
610  const int gg,
611  const int nb_base_functions) {
612  // FIXME: That is bug, it will not work if number of integration pts is
613  // equal to number of nodes on entity. User who not implementing low level
614  // DataOperator will not experience this.
615  if (getN(base).size1() == getDiffN(base).size1()) {
616  (void)getN(base)(gg,
617  nb_base_functions -
618  1); // throw error if nb_base_functions is to big
619  int dim = getDiffN(base).size2() / getN(base).size2();
620  double *data = &getDiffN(base)(gg, 0);
621  return MatrixAdaptor(nb_base_functions, dim,
622  ublas::shallow_array_adaptor<double>(
623  dim * nb_base_functions, data));
624  } else {
625  // in some cases, f.e. for derivatives of nodal base functions only one
626  // gauss point is needed
627  return MatrixAdaptor(
628  getN(base).size1(), getN(base).size2(),
629  ublas::shallow_array_adaptor<double>(getDiffN(base).data().size(),
630  &getDiffN(base).data()[0]));
631  }
632  }
633 
634  /** \brief get derivatives of base functions at Gauss pts
635  *
636  * Note that multi field element, two different field can have different
637  * approximation orders. Since we use hierarchical approximation basis,
638  * base functions are calculated once for element, using maximal
639  * approximation order on given entity.
640  *
641  * \param gg nb. of Gauss point
642  * \param nb_base_functions number of of base functions
643  *
644  */
645  inline const MatrixAdaptor getDiffN(const int gg,
646  const int nb_base_functions) {
647  return getDiffN(bAse, gg, nb_base_functions);
648  }
649 
650  /**@}*/
651 
652  /** \name Get base functions for vectorial approximation basese, i.e.
653  * Hdiv/Hcurl */
654 
655  /**@{*/
656 
657  /** \brief get Hdiv of base functions at Gauss pts
658  *
659  * \param base Approximation base
660  * \param gg nb. of Gauss point
661  *
662  */
663  template <int DIM>
665  const int gg) {
666  if (PetscUnlikely(getN(base).size2() % DIM)) {
667  THROW_MESSAGE("Wrong dimension");
668  }
669 
670  const int nb_base_functions = getN(base).size2() / DIM;
671  double *data = &getN(base)(gg, 0);
672  return MatrixAdaptor(
673  nb_base_functions, DIM,
674  ublas::shallow_array_adaptor<double>(DIM * nb_base_functions, data));
675  }
676 
677  /** \brief get Hdiv of base functions at Gauss pts
678  *
679  * \param gg nb. of Gauss point
680  * \param number of of base functions
681  *
682  */
683  template<int DIM>
684  inline const MatrixAdaptor getVectorN(const int gg) {
685  return getVectorN<DIM>(bAse, gg);
686  }
687 
688  /** \brief get DiffHdiv of base functions at Gauss pts
689  *
690  * \param base Approximation base
691  * \param gg nb. of Gauss point
692  * \param number of of base functions
693  *
694  */
695  template <int DIM0, int DIM1>
697  const int gg) {
698  if (PetscUnlikely(getDiffN(base).size2() % (DIM0*DIM1))) {
699  THROW_MESSAGE("Wrong dimension");
700  }
701 
702  const int nb_base_functions = getN(base).size2() / (DIM0 * DIM1);
703  double *data = &getN(base)(gg, 0);
704  return MatrixAdaptor(nb_base_functions, DIM0 * DIM1,
705  ublas::shallow_array_adaptor<double>(
706  DIM0 * DIM1 * nb_base_functions, data));
707  }
708 
709  /** \brief get DiffHdiv of base functions at Gauss pts
710  *
711  * \param gg nb. of Gauss point
712  * \param number of of base functions
713  *
714  */
715  template <int DIM0, int DIM1>
716  inline const MatrixAdaptor getVectorDiffN(const int gg) {
717  return getVectorDiffN<DIM0,DIM1>(bAse, gg);
718  }
719 
720  /** \brief get DiffHdiv of base functions at Gauss pts
721  *
722  * \param base Approximation base
723  * \param gg nb. of Gauss point
724  * \param number of of base functions
725  *
726  */
727  template <int DIM0, int DIM1>
729  const int dof, const int gg) {
730  double *data = &getDiffN(base)(gg, DIM0 * DIM1 * dof);
731  return MatrixAdaptor(
732  DIM0, DIM1, ublas::shallow_array_adaptor<double>(DIM0 * DIM1, data));
733  }
734 
735  /** \brief get DiffHdiv of base functions at Gauss pts
736  *
737  * \param gg nb. of Gauss point
738  * \param number of of base functions
739  *
740  */
741  template <int DIM0, int DIM1>
742  inline const MatrixAdaptor getVectorDiffN(const int dof, const int gg) {
743  return getVectorDiffN<DIM0,DIM1>(bAse, dof, gg);
744  }
745 
746  /**@}*/
747 
748  /** \name Get base functions with FTensor */
749 
750  /**@{*/
751 
752  /**
753  * \brief Get base function as Tensor0
754  *
755  * \param base
756  * \return Tensor0
757  *
758  */
761  double *ptr = &*getN(base).data().begin();
763  };
764 
765  /**
766  * \brief Get base function as Tensor0
767  *
768  * Return base functions for field base
769  *
770  * \return Tensor0
771  *
772  */
774  return getFTensor0N(bAse);
775  };
776 
777  /**
778  * \brief Get base function as Tensor0 (Loop by integration points)
779  *
780  * \param base
781  * \param bb base function
782  * \return Tensor0
783 
784  Note that:
785  \code
786  t0 = data.getFTensor0N(base,bb);
787  ++t0
788  \endcode
789  Increment in above code will move pointer to base function in next
790  integration point.
791 
792  *
793  */
795  getFTensor0N(const FieldApproximationBase base, const int bb) {
796  double *ptr = &getN(base)(0, bb);
797  return FTensor::Tensor0<double *>(ptr, getN(base).size2());
798  };
799 
800  /**
801  * \brief Get base function as Tensor0 (Loop by integration points)
802  *
803  * Return base functions for field base
804  *
805  * \param bb base function
806  * \return Tensor0
807  *
808  *
809  */
811  return getFTensor0N(bAse, bb);
812  };
813 
814  /**
815  * \brief Get base function as Tensor0 (Loop by integration points)
816  *
817  * \param base
818  * \param gg integration points
819  * \param bb base function
820  * \return Tensor0
821 
822  Note that:
823  \code
824  t0 = data.getFTensor0N(base,bb);
825  ++t0
826  \endcode
827  Increment in above code will move pointer to base function in next
828  integration point.
829 
830  *
831  */
833  getFTensor0N(const FieldApproximationBase base, const int gg,
834  const int bb) {
835  double *ptr = &getN(base)(gg, bb);
837  };
838 
839  /**
840  * \brief Get base function as Tensor0 (Loop by integration points)
841  *
842  * Return base functions for field base
843  *
844  * \param bb base function
845  * \return Tensor0
846  *
847  */
849  getFTensor0N(const int gg, const int bb) {
850  return getFTensor0N(bAse, gg, bb);
851  };
852 
853  /**
854  * \brief Get derivatives of base functions
855  *
856  * For volume element like tetrahedral or prism,
857  * \code
858  * Tensor1<double*,3> diff_base = data.getFTensor1DiffN<3>();
859  * \endcode
860  *
861  * For face element like triangle or quad
862  * \code
863  * Tensor1<double*,2> diff_base = data.getFTensor1DiffN<2>();
864  * \endcode
865  *
866  * \param base functions
867  * \return Tensor rank 1 (vector)
868  *
869  */
870  template <int Tensor_Dim>
872  getFTensor1DiffN(const FieldApproximationBase base);
873 
874  /**
875  * \brief Get derivatives of base functions
876  *
877  * For volume element like tetrahedral or prism,
878  * \code
879  * Tensor1<double*,3> diff_base = data.getFTensor1DiffN<3>();
880  * \endcode
881  *
882  * For face element like triangle or quad
883  * \code
884  * Tensor1<double*,2> diff_base = data.getFTensor1DiffN<2>();
885  * \endcode
886  *
887  * \return Tensor rank 1 (vector)
888  *
889  */
890  template <int Tensor_Dim>
891  FTensor::Tensor1<double *, Tensor_Dim> getFTensor1DiffN();
892 
893  /**
894  * \brief Get derivatives of base functions (Loop by integration points)
895  *
896  * For volume element like tetrahedral or prism,
897  * \code
898  * Tensor1<double*,3> diff_base = data.getFTensor1DiffN<3>(base,bb);
899  * \endcode
900  * where bb is base function. Operator ++diff_base will move tensor pointer
901  * to next integration point.
902  *
903  * For face element like triangle or quad
904  * \code
905  * Tensor1<double*,2> diff_base = data.getFTensor1DiffN<2>(base,bb);
906  * \endcode
907  *
908  * \param base functions
909  * \return Tensor rank 1 (vector)
910  *
911  */
912  template <int Tensor_Dim>
914  getFTensor1DiffN(const FieldApproximationBase base, const int bb);
915 
916  /**
917  * \brief Get derivatives of base functions (Loop by integration points)
918  *
919  * For volume element like tetrahedral or prism,
920  * \code
921  * Tensor1<double*,3> diff_base = data.getFTensor1DiffN<3>(bb);
922  * \endcode
923  * where bb is base function. Operator ++diff_base will move tensor pointer
924  * to next integration point.
925  *
926  * For face element like triangle or quad
927  * \code
928  * Tensor1<double*,2> diff_base = data.getFTensor1DiffN<2>(bb);
929  * \endcode
930  *
931  * \return Tensor rank 1 (vector)
932  *
933  */
934  template <int Tensor_Dim>
935  FTensor::Tensor1<double *, Tensor_Dim> getFTensor1DiffN(const int bb);
936 
937  /**
938  * \brief Get derivatives of base functions (Loop by integration points)
939  *
940  * For volume element like tetrahedral or prism,
941  * \code
942  * Tensor1<double*,3> diff_base = data.getFTensor1DiffN<3>(base,gg,bb);
943  * \endcode
944  * where bb is base function and gg is integration pt. Operator ++diff_base
945  * will move tensor pointer to next integration point.
946  *
947  * For face element like triangle or quad
948  * \code
949  * Tensor1<double*,2> diff_base = data.getFTensor1DiffN<2>(base,gg,bb);
950  * \endcode
951  *
952  * \return Tensor rank 1 (vector)
953  *
954  */
955  template <int Tensor_Dim>
957  getFTensor1DiffN(const FieldApproximationBase base, const int gg,
958  const int bb);
959 
960  /**
961  * \brief Get derivatives of base functions (Loop by integration points)
962  *
963  * For volume element like tetrahedral or prism,
964  * \code
965  * Tensor1<double*,3> diff_base = data.getFTensor1DiffN<3>(gg,bb);
966  * \endcode
967  * where bb is base function and gg is integration pt. Operator ++diff_base
968  * will move tensor pointer to next integration point.
969  *
970  * For face element like triangle or quad
971  * \code
972  * Tensor1<double*,2> diff_base = data.getFTensor1DiffN<2>(gg,bb);
973  * \endcode
974  *
975  * \return Tensor rank 1 (vector)
976  *
977  */
978  template <int Tensor_Dim>
979  FTensor::Tensor1<double *, Tensor_Dim> getFTensor1DiffN(const int gg,
980  const int bb);
981 
982  /** \brief Get base functions for Hdiv/Hcurl spaces
983 
984  \note You probably like to use getFTensor1N(), in typical use base is
985  set automatically based on base set to field.
986 
987  * @param base Approximation base
988 
989  Example:
990  \code
991  FTensor::Index<'i',3> i;
992  int nb_dofs = data.getFieldData().size();
993  auto t_n_hdiv = data.getFTensor1N<3>();
994  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
995  int ll = 0;
996  for(;ll!=nb_dofs;ll++) {
997  double dot_product = t_n_hdiv(i)*t_n_hdiv(i);
998  ++t_n_hdiv;
999  }
1000  for(;ll!=data.getVectorN().size2()/3;ll++) {
1001  ++t_n_hdiv;
1002  }
1003  }
1004  \endcode
1005 
1006  */
1007  template <int Tensor_Dim>
1009  getFTensor1N(FieldApproximationBase base);
1010 
1011  /** \brief Get base functions for Hdiv space
1012 
1013  Example:
1014  \code
1015  FTensor::Index<'i',3> i;
1016  int nb_dofs = data.getFieldData().size();
1017  auto t_n_hdiv = data.getFTensor1N<3>();
1018  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
1019  int ll = 0;
1020  for(;ll!=nb_dofs;ll++) {
1021  double dot_product = t_n_hdiv(i)*t_n_hdiv(i);
1022  ++t_n_hdiv;
1023  }
1024  for(;ll!=data.getVectorN().size2()/3;ll++) {
1025  ++t_n_hdiv;
1026  }
1027  }
1028  \endcode
1029 
1030  */
1031  template <int Tensor_Dim> auto getFTensor1N() {
1032  return getFTensor1N<Tensor_Dim>(bAse);
1033  }
1034 
1035  /** \brief Get derivatives of base functions for Hdiv space
1036  */
1037  template <int Tensor_Dim0, int Tensor_Dim1>
1039  Tensor_Dim0, Tensor_Dim1>
1040  getFTensor2DiffN(FieldApproximationBase base);
1041 
1042  /** \brief Get derivatives of base functions for Hdiv space at integration
1043  * pts
1044  */
1045  template <int Tensor_Dim0, int Tensor_Dim1>
1047  Tensor_Dim0, Tensor_Dim1>
1048  getFTensor2DiffN(FieldApproximationBase base, const int gg,
1049  const int bb);
1050 
1051  /** \brief Get derivatives of base functions for Hdiv space
1052  */
1053  template <int Tensor_Dim0, int Tensor_Dim1>
1054  inline FTensor::Tensor2<
1056  Tensor_Dim1>
1058  return getFTensor2DiffN<Tensor_Dim0, Tensor_Dim1>(bAse);
1059  }
1060 
1061  /** \brief Get derivatives of base functions for Hdiv space at integration
1062  * pts
1063  */
1064  template <int Tensor_Dim0, int Tensor_Dim1>
1065  inline FTensor::Tensor2<
1067  Tensor_Dim1>
1068  getFTensor2DiffN(const int gg, const int bb) {
1069  return getFTensor2DiffN<Tensor_Dim0, Tensor_Dim1>(bAse, gg, bb);
1070  }
1071 
1072  /**
1073  * \brief Get Hdiv base functions at integration point
1074 
1075  \code
1076  FTensor::Index<'i',3> i;
1077  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
1078  auto t_base = data.getFTensor1N(base,gg,bb);
1079  for(int bb = 0;bb!=nb_base_functions;bb++) {
1080  auto dot = t_base(i)*t_base(i);
1081  }
1082  }
1083  \endcode
1084 
1085  */
1086  template <int Tensor_Dim>
1088  getFTensor1N(FieldApproximationBase base, const int gg, const int bb);
1089 
1090  /**
1091  * \brief Get Hdiv base functions at integration point
1092 
1093  \code
1094  FTensor::Index<'i',3> i;
1095  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
1096  auto t_base = data.getFTensor1N(gg,0);
1097  for(int bb = 0;bb!=nb_base_functions;bb++) {
1098  double dot = t_base(i)*t_base(i);
1099  }
1100  }
1101  \endcode
1102 
1103  */
1104  template <int Tensor_Dim>
1105  inline auto getFTensor1N(const int gg, const int bb) {
1106  return getFTensor1N<Tensor_Dim>(bAse, gg, bb);
1107  }
1108 
1109  /** \brief Get base functions for Hdiv/Hcurl spaces
1110 
1111  \note You probably like to use getFTensor1N(), in typical use base is
1112  set automatically based on base set to field.
1113 
1114  * @param base Approximation base
1115 
1116  Example:
1117  \code
1118  FTensor::Index<'i',3> i;
1119  FTensor::Index<'i',3> j;
1120  int nb_dofs = data.getFieldData().size();
1121  auto t_n_hdiv = data.getFTensor2N<3,3>();
1122  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
1123  int ll = 0;
1124  for(;ll!=nb_dofs;ll++) {
1125  double dot_product = t_n_hdiv(i,j)*t_n_hdiv(i,j);
1126  ++t_n_hdiv;
1127  }
1128  for(;ll!=data.getVectorN().size2()/3;ll++) {
1129  ++t_n_hdiv;
1130  }
1131  }
1132  \endcode
1133 
1134  */
1135  template <int Tensor_Dim0, int Tensor_Dim1>
1137  Tensor_Dim0, Tensor_Dim1>
1138  getFTensor2N(FieldApproximationBase base);
1139 
1140  /** \brief Get base functions for Hdiv space
1141 
1142  Example:
1143  \code
1144  FTensor::Index<'i',3> i;
1145  FTensor::Index<'j',3> j;
1146 
1147  int nb_dofs = data.getFieldData().size();
1148  auto t_n_hdiv = data.getFTensor2N<3,3>();
1149  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
1150  int ll = 0;
1151  for(;ll!=nb_dofs;ll++) {
1152  double dot_product = t_n_hdiv(i,j)*t_n_hdiv(i,j);
1153  ++t_n_hdiv;
1154  }
1155  for(;ll!=data.getVectorN().size2()/3;ll++) {
1156  ++t_n_hdiv;
1157  }
1158  }
1159  \endcode
1160 
1161  */
1162  template <int Tensor_Dim0, int Tensor_Dim1> auto getFTensor2N() {
1163  return getFTensor2N<Tensor_Dim0, Tensor_Dim1>(bAse);
1164  }
1165 
1166  /** \brief Get base functions for tensor Hdiv/Hcurl spaces
1167 
1168  \note You probably like to use getFTensor2N(), in typical use base is
1169  set automatically based on base set to field.
1170 
1171  @param base Approximation base
1172 
1173  Example:
1174  \code
1175  FTensor::Index<'i',3> i;
1176  FTensor::Index<'j',3> i;
1177  int nb_dofs = data.getFieldData().size();
1178  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
1179  auto t_n_hdiv = data.getFTensor2N<3>(base,gg,bb);
1180  int ll = 0;
1181  for(;ll!=nb_dofs;ll++) {
1182  double dot_product = t_n_hdiv(i,j)*t_n_hdiv(i,j);
1183  ++t_n_hdiv;
1184  }
1185  for(;ll!=data.getVectorN().size2()/3;ll++) {
1186  ++t_n_hdiv;
1187  }
1188  }
1189  \endcode
1190 
1191  */
1192  template <int Tensor_Dim0, int Tensor_Dim1>
1194  Tensor_Dim0, Tensor_Dim1>
1195  getFTensor2N(FieldApproximationBase base, const int gg, const int bb);
1196 
1197  /** \brief Get base functions for Hdiv space
1198 
1199  Example:
1200  \code
1201  FTensor::Index<'i',3> i;
1202  FTensor::Index<'j',3> j;
1203  int nb_dofs = data.getFieldData().size();
1204  for(int gg = 0;gg!=nb_gauss_pts;++gg) {
1205  int ll = 0;
1206  auto t_n_hdiv = data.getFTensor2N<3,3>(gg,0);
1207  for(;ll!=nb_dofs;ll++) {
1208  double dot_product = t_n_hdiv(i)*t_n_hdiv(i);
1209  ++t_n_hdiv;
1210  }
1211  for(;ll!=data.getVectorN().size2()/3;ll++) {
1212  ++t_n_hdiv;
1213  }
1214  }
1215  \endcode
1216 
1217  */
1218  template <int Tensor_Dim0, int Tensor_Dim1>
1219  auto getFTensor2N(const int gg, const int bb) {
1220  return getFTensor2N<Tensor_Dim0, Tensor_Dim1>(bAse, gg, bb);
1221  }
1222 
1223  /**@}*/
1224 
1225  /** \name Auxiliary functions */
1226 
1227  /**@{*/
1228 
1229  friend std::ostream &operator<<(std::ostream &os,
1231 
1232  /**
1233  * Reset data associated with particular field name
1234  * @return error code
1235  */
1238  sPace = NOSPACE;
1239  bAse = NOBASE;
1240  iNdices.resize(0, false);
1241  localIndices.resize(0, false);
1242  dOfs.resize(0, false);
1243  fieldData.resize(0, false);
1245  }
1246 
1247  /**@}*/
1248 
1249  protected:
1250  int sEnse; ///< Entity sense (orientation)
1251  ApproximationOrder oRder; ///< Entity order
1252  FieldSpace sPace; ///< Entity space
1253  FieldApproximationBase bAse; ///< Field approximation base
1254  VectorInt iNdices; ///< Global indices on entity
1255  VectorInt localIndices; ///< Local indices on entity
1256  VectorDofs dOfs; ///< DoFs on entity
1257  VectorDouble fieldData; ///< Field data on entity
1258  boost::shared_ptr<MatrixDouble> N[LASTBASE]; ///< Base functions
1259  boost::shared_ptr<MatrixDouble>
1260  diffN[LASTBASE]; ///< Derivatives of base functions
1261  };
1262 
1263  std::bitset<LASTSPACE> sPace; ///< spaces on element
1264  std::bitset<LASTBASE> bAse; ///< bases on element
1265  ublas::matrix<int> facesNodes; ///< nodes on finite element faces
1266  std::bitset<LASTSPACE>
1267  spacesOnEntities[MBMAXTYPE]; ///< spaces on entity types
1268  std::bitset<LASTBASE> basesOnEntities[MBMAXTYPE]; ///< bases on entity types
1269  std::bitset<LASTBASE> basesOnSpaces[LASTSPACE]; ///< base on spaces
1270 
1271  boost::ptr_vector<EntData> dataOnEntities[MBMAXTYPE]; ///< data on nodes, base
1272  ///< function, dofs
1273  ///< values, etc.
1274 
1275  /**
1276  * Reset data associated with particular field name
1277  * @return error code
1278  */
1281  for (EntityType t = MBVERTEX; t != MBMAXTYPE; t++) {
1282  for (auto &e : dataOnEntities[t]) {
1283  ierr = e.resetFieldDependentData();
1284  CHKERRG(ierr);
1285  }
1286  }
1288  }
1289 
1290  DataForcesAndSourcesCore(const EntityType type);
1292 
1293  virtual MoFEMErrorCode setElementType(const EntityType type);
1294 
1295  friend std::ostream &operator<<(std::ostream &os,
1296  const DataForcesAndSourcesCore &e);
1297 protected:
1299 };
1300 
1301 /** \brief this class derive data form other data structure
1302  * \ingroup mofem_forces_and_sources_user_data_operators
1303  *
1304  *
1305  * It behaves like normal data structure it is used to share base functions with
1306  * other data structures. Dofs values, approx. order and
1307  * indices are not shared.
1308  *
1309  * Shape functions, senses are shared with other data structure.
1310  *
1311  */
1313 
1314  /** \brief Derived ata on single entity (This is passed as argument to
1315  * DataOperator::doWork) \ingroup mofem_forces_and_sources_user_data_operators
1316  * \nosubgrouping
1317  *
1318  * DerivedEntData share part information with EntData except infomation about
1319  * base functions.
1320  *
1321  */
1323 
1324  const boost::shared_ptr<DataForcesAndSourcesCore::EntData> entDataPtr;
1325  DerivedEntData(const boost::shared_ptr<DataForcesAndSourcesCore::EntData>
1326  &ent_data_ptr)
1327  : entDataPtr(ent_data_ptr) {}
1328 
1329  int getSense() const { return entDataPtr->getSense(); }
1330 
1331  boost::shared_ptr<MatrixDouble> &
1333  return entDataPtr->getNSharedPtr(base);
1334  }
1335  boost::shared_ptr<MatrixDouble> &
1337  return entDataPtr->getDiffNSharedPtr(base);
1338  }
1339  const boost::shared_ptr<MatrixDouble> &
1341  return entDataPtr->getNSharedPtr(base);
1342  }
1343  const boost::shared_ptr<MatrixDouble> &
1345  return entDataPtr->getDiffNSharedPtr(base);
1346  }
1347 
1348  };
1349 
1351  const boost::shared_ptr<DataForcesAndSourcesCore> &data_ptr);
1352  MoFEMErrorCode setElementType(const EntityType type);
1353 
1354 private:
1355  const boost::shared_ptr<DataForcesAndSourcesCore> dataPtr;
1356 };
1357 
1358 /** \name Specializations for H1/L2 */
1359 
1360 /**@{*/
1361 
1362 template <>
1364 DataForcesAndSourcesCore::EntData::getFTensor1FieldData<3>();
1365 
1366 template <>
1368 DataForcesAndSourcesCore::EntData::getFTensor1FieldData<2>();
1369 
1370 template <>
1372 DataForcesAndSourcesCore::EntData::getFTensor2FieldData<3, 3>();
1373 
1374 template <>
1376 DataForcesAndSourcesCore::EntData::getFTensor2SymmetricFieldData<3>();
1377 
1378 template <>
1380 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<3>(
1381  const FieldApproximationBase base);
1382 template <>
1384 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<3>();
1385 
1386 template <>
1388 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<2>(
1389  const FieldApproximationBase base);
1390 template <>
1392 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<2>();
1393 template <>
1395 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<3>(
1396  const FieldApproximationBase base, const int bb);
1397 template <>
1399 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<3>(const int bb);
1400 template <>
1402 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<2>(
1403  const FieldApproximationBase base, const int bb);
1404 template <>
1406 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<2>(const int bb);
1407 
1408 /**@}*/
1409 
1410 /** \name Specializations for HDiv/HCurl in 2d */
1411 
1412 /**@{*/
1413 
1414 template <>
1416 DataForcesAndSourcesCore::EntData::getFTensor2DiffN<3, 2>(
1417  FieldApproximationBase base);
1418 template <>
1420 DataForcesAndSourcesCore::EntData::getFTensor2DiffN<3, 2>(
1421  FieldApproximationBase base, const int gg, const int bb);
1422 
1423 /**@}*/
1424 
1425 /// \deprecated Use DataForcesAndSourcesCore
1427 
1428 /// \deprecated use DerivedDataForcesAndSourcesCore
1431 
1432 } // namespace MoFEM
1433 
1434 #endif //__DATASTRUCTURES_HPP
1435 
1436 /**
1437  * \defgroup mofem_forces_and_sources_user_data_operators User data operator
1438  * data structures \ingroup
1439  *
1440  * \brief Users data structures and operator
1441  *
1442  * Data structures passed by argument to MoFEM::DataOperator::doWork and generic
1443  * user operators operating on those structures.
1444  *
1445  */
const MatrixAdaptor getDiffN(const FieldApproximationBase base, const int gg)
get derivative of base functions at Gauss pts
FTensor::Tensor2_symmetric< FTensor::PackPtr< double *,(Tensor_Dim *(Tensor_Dim+1))/2 >, Tensor_Dim > getFTensor2SymmetricFieldData()
Return symmetric FTensor rank 2, i.e. matrix from filed data coeffinects.
ublas::matrix< int > facesNodes
nodes on finite element faces
boost::shared_ptr< MatrixDouble > & getNSharedPtr(const FieldApproximationBase base)
virtual boost::shared_ptr< MatrixDouble > & getDiffNSharedPtr(const FieldApproximationBase base)
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Common.hpp:60
boost::shared_ptr< MatrixDouble > & getDiffNSharedPtr(const FieldApproximationBase base)
data structure for finite element entityIt keeps that about indices of degrees of freedom...
const VectorAdaptor getN(const FieldApproximationBase base, const int gg, const int nb_base_functions)
get base functions at Gauss pts
const MatrixAdaptor getDiffN(const int gg)
get derivative of base functions at Gauss pts
DEPRECATED FTensor::Tensor1< FTensor::PackPtr< double *, 1 >, Tensor_Dim > getTensor1FormData(MatrixDouble &data)
const VectorIntAdaptor getIndicesUpToOrder(int order)
get global indices of dofs on entity up to given order
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:483
ApproximationOrder oRder
Entity order.
const VectorInt & getLocalIndices() const
get local indices of dofs on entity
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const int gg, const int bb)
Get base function as Tensor0 (Loop by integration points)
virtual boost::shared_ptr< MatrixDouble > & getNSharedPtr(const FieldApproximationBase base)
const MatrixAdaptor getDiffN(const FieldApproximationBase base, const int gg, const int nb_base_functions)
get derivatives of base functions at Gauss pts
VectorShallowArrayAdaptor< int > VectorIntAdaptor
Definition: Common.hpp:235
const VectorDofs & getFieldDofs() const
get dofs data stature FEDofEntity
VectorInt iNdices
Global indices on entity.
FTensor::Tensor2< FTensor::PackPtr< double *, 1 >, Tensor_Dim0, Tensor_Dim1 > getFTensor2FromMat(MatrixDouble &data)
Get tensor rank 2 (matrix) form data matrix (specialization)
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:526
FTensor::Tensor2< FTensor::PackPtr< double *, Tensor_Dim0 *Tensor_Dim1 >, Tensor_Dim0, Tensor_Dim1 > getFTensor2DiffN(const int gg, const int bb)
Get derivatives of base functions for Hdiv space at integration pts.
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1FieldData()
Return FTensor of rank 1, i.e. vector from filed data coeffinects.
FTensor::Tensor2< FTensor::PackPtr< double *, Tensor_Dim0 *Tensor_Dim1 >, Tensor_Dim0, Tensor_Dim1 > getFTensor2DiffN()
Get derivatives of base functions for Hdiv space.
const MatrixAdaptor getVectorDiffN(FieldApproximationBase base, const int gg)
get DiffHdiv of base functions at Gauss pts
#define THROW_MESSAGE(a)
Throw MoFEM exception.
Definition: definitions.h:602
virtual const boost::shared_ptr< MatrixDouble > & getDiffNSharedPtr(const FieldApproximationBase base) const
DEPRECATED typedef DerivedDataForcesAndSourcesCore DerivedDataForcesAndSurcesCore
MatrixDouble & getDiffN()
Get derivatives of base functions.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:490
FTensor::Tensor1< FTensor::PackPtr< double *, 1 >, Tensor_Dim > getFTensor1FromMat(MatrixDouble &data)
Get tensor rank 1 (vector) form data matrix (specialization)
virtual const MatrixDouble & getN(const FieldApproximationBase base) const
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb. of columns is equal to number of base functions on this entity.
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
DEPRECATED FTensor::Tensor2< FTensor::PackPtr< double *, 1 >, Tensor_Dim0, Tensor_Dim1 > getTensor2FormData(MatrixDouble &data)
const MatrixAdaptor getVectorDiffN(const FieldApproximationBase base, const int dof, const int gg)
get DiffHdiv of base functions at Gauss pts
const MatrixAdaptor getVectorN(const FieldApproximationBase base, const int gg)
get Hdiv of base functions at Gauss pts
std::vector< double, std::allocator< double > > DoubleAllocator
Definition: Common.hpp:207
int getSense() const
get entity sense, need to calculate base functions with conforming approximation fields ...
FTensor::Tensor2< FTensor::PackPtr< double *, Tensor_Dim0 *Tensor_Dim1 >, Tensor_Dim0, Tensor_Dim1 > getFTensor2FieldData()
Return FTensor rank 2, i.e. matrix from filed data coeffinects.
FieldSpace in [ 0, LASTSPACE )
Definition: definitions.h:172
VectorShallowArrayAdaptor< double > VectorAdaptor
Definition: Common.hpp:234
const boost::shared_ptr< DataForcesAndSourcesCore::EntData > entDataPtr
virtual const MatrixDouble & getDiffN(const FieldApproximationBase base) const
get derivatives of base functions
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0FromVec< double, DoubleAllocator >(ublas::vector< double, DoubleAllocator > &data)
std::ostream & operator<<(std::ostream &os, const DataForcesAndSourcesCore::EntData &e)
MatrixShallowArrayAdaptor< double > MatrixAdaptor
Definition: Common.hpp:266
const VectorAdaptor getFieldDataUpToOrder(int order)
get dofs values up to given order
MatrixDouble & getDiffN(const FieldApproximationBase base)
Get derivatives of base functions.
virtual const boost::shared_ptr< MatrixDouble > & getNSharedPtr(const FieldApproximationBase base) const
VectorDouble fieldData
Field data on entity.
ublas::vector< int, IntAllocator > VectorInt
Definition: Common.hpp:210
FieldApproximationBase
approximation base
Definition: definitions.h:140
const MatrixAdaptor getVectorDiffN(const int gg)
get DiffHdiv of base functions at Gauss pts
const boost::shared_ptr< DataForcesAndSourcesCore > dataPtr
const boost::shared_ptr< MatrixDouble > & getDiffNSharedPtr(const FieldApproximationBase base) const
auto getFTensor1N(const int gg, const int bb)
Get Hdiv base functions at integration point.
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Common.hpp:212
const VectorDouble & getFieldData() const
get dofs values
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base, const int gg, const int bb)
Get base function as Tensor0 (Loop by integration points)
const MatrixDouble & getDiffN() const
const MatrixAdaptor getVectorDiffN(const int dof, const int gg)
get DiffHdiv of base functions at Gauss pts
DEPRECATED FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getTensor0FormData(ublas::vector< T, A > &data)
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Common.hpp:211
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vectorExample how to use it.
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N()
Get base function as Tensor0.
FTensor::Tensor0< double * > getFTensor0N(const int bb)
Get base function as Tensor0 (Loop by integration points)
const VectorAdaptor getN(const int gg)
get base functions at Gauss pts
FieldSpace & getSpace()
Get field space.
DataForcesAndSourcesCore::EntData EntData
ApproximationOrder getOrder() const
get approximation order
FieldSpace
approximation spaces
Definition: definitions.h:165
DEPRECATED typedef DataForcesAndSourcesCore DataForcesAndSurcesCore
FTensor::Tensor2_symmetric< FTensor::PackPtr< double *, 1 >, Tensor_Dim > getFTensor2SymmetricFromMat(MatrixDouble &data)
Get symmetric tensor rank 2 form matrix.
Derived ata on single entity (This is passed as argument to DataOperator::doWork) ...
FieldApproximationBase bAse
Field approximation base.
const VectorIntAdaptor getLocalIndicesUpToOrder(int order)
get local indices of dofs on entity up to given order
ublas::vector< boost::shared_ptr< const FEDofEntity >, DofsAllocator > VectorDofs
FieldApproximationBase & getBase()
Get approximation base.
#define DEPRECATED
Definition: definitions.h:29
DerivedEntData(const boost::shared_ptr< DataForcesAndSourcesCore::EntData > &ent_data_ptr)
ublas::unbounded_array< boost::shared_ptr< const FEDofEntity >, std::allocator< boost::shared_ptr< const FEDofEntity > > > DofsAllocator
const MatrixAdaptor getVectorN(const int gg)
get Hdiv of base functions at Gauss pts
std::bitset< LASTSPACE > sPace
spaces on element
const MatrixAdaptor getDiffN(const int gg, const int nb_base_functions)
get derivatives of base functions at Gauss pts
const boost::shared_ptr< MatrixDouble > & getNSharedPtr(const FieldApproximationBase base) const
const VectorInt & getIndices() const
Get global indices of dofs on entity.
Data on single entity (This is passed as argument to DataOperator::doWork)
auto getVectorAdaptor
Get Vector adaptor.
Definition: Common.hpp:256
auto getFTensor2N()
Get base functions for Hdiv space.
FTensor::Tensor0< double * > getFTensor0N(const FieldApproximationBase base, const int bb)
Get base function as Tensor0 (Loop by integration points)
const int N
Definition: speed_test.cpp:3
std::bitset< LASTBASE > bAse
bases on element
const VectorAdaptor getN(const int gg, const int nb_base_functions)
get base functions at Gauss pts
VectorInt localIndices
Local indices on entity.
int sEnse
Entity sense (orientation)
auto getFTensor2N(const int gg, const int bb)
Get base functions for Hdiv space.
auto getFTensor1N()
Get base functions for Hdiv space.
const VectorAdaptor getN(const FieldApproximationBase base, const int gg)
get base functions at Gauss pts
this class derive data form other data structure
MatrixDouble & getN(const FieldApproximationBase base)
Get base functions.
MatrixDouble & getN()
Get base functions.
virtual int getSense() const
get entity sense, need to calculate base functions with conforming approximation fields ...
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base)
Get base function as Tensor0.