v0.8.20
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 /** \brief data structure for finite element entity
37  * \ingroup mofem_forces_and_sources_user_data_operators
38  *
39  * It keeps that about indices of degrees of freedom, dofs data, base functions
40  * functions, entity side number, type of entities, approximation order, etc.
41  *
42  */
44 
45  /** \brief Data on single entity (This is passed as argument to
46  * DataOperator::doWork) \ingroup mofem_forces_and_sources_user_data_operators
47  * \nosubgrouping
48  *
49  * \todo Hdiv and Hcurl functions should be accessed through common interface.
50  */
51  struct EntData {
52 
53  /** \name Constructor and destructor */
54 
55  /**@{*/
56 
57  EntData();
58  virtual ~EntData() {}
59 
60  /**@}*/
61 
62  /** \name Sense, order and indices */
63 
64  /**@{*/
65 
66  /// \brief get entity sense, need to calculate base functions with
67  /// conforming approximation fields
68  virtual int getSense() const { return sEnse; }
69 
70  /// \brief get approximation order
71  inline ApproximationOrder getOrder() const { return oRder; }
72 
73  /// \brief Get global indices of dofs on entity
74  inline const VectorInt &getIndices() const { return iNdices; }
75 
76  /// \brief get global indices of dofs on entity up to given order
77  inline const VectorIntAdaptor getIndicesUpToOrder(int order) {
78  unsigned int size = 0;
79  if (iNdices.size()) {
80  size = dOfs[0]->getOrderNbDofs(order) * dOfs[0]->getNbOfCoeffs();
81  size = size < iNdices.size() ? size : iNdices.size();
82  }
83  int *data = &*iNdices.data().begin();
84  return VectorIntAdaptor(size,
85  ublas::shallow_array_adaptor<int>(size, data));
86  }
87 
88  /// \brief get local indices of dofs on entity
89  inline const VectorInt &getLocalIndices() const { return localIndices; }
90 
91  /// \brief get local indices of dofs on entity up to given order
92  inline const VectorIntAdaptor getLocalIndicesUpToOrder(int order) {
93  unsigned int size = 0;
94  if (localIndices.size()) {
95  size = dOfs[0]->getOrderNbDofs(order) * dOfs[0]->getNbOfCoeffs();
96  size = size < localIndices.size() ? size : localIndices.size();
97  }
98  int *data = &*localIndices.data().begin();
99  return VectorIntAdaptor(size,
100  ublas::shallow_array_adaptor<int>(size, data));
101  }
102 
103  inline int &getSense() { return sEnse; }
104  inline ApproximationOrder &getDataOrder() { return oRder; }
105  inline VectorInt &getIndices() { return iNdices; }
106  inline VectorInt &getLocalIndices() { return localIndices; }
107 
108  /**@}*/
109 
110  /** \name Data on entity */
111 
112  /**@{*/
113 
114  /// \brief get dofs values
115  inline const VectorDouble &getFieldData() const { return fieldData; }
116 
117  /// \brief get dofs values up to given order
118  inline const VectorAdaptor getFieldDataUpToOrder(int order) {
119  unsigned int size = 0;
120  if (fieldData.size()) {
121  size = dOfs[0]->getOrderNbDofs(order) * dOfs[0]->getNbOfCoeffs();
122  size = size < fieldData.size() ? size : fieldData.size();
123  }
124  double *data = &*fieldData.data().begin();
125  return getVectorAdaptor(data, size);
126  }
127 
128  /// \brief get dofs data stature FEDofEntity
129  inline const VectorDofs &getFieldDofs() const { return dOfs; }
130 
131  inline VectorDouble &getFieldData() { return fieldData; }
132 
133  /**
134  * @brief Return FTensor of rank 1, i.e. vector from filed data coeffinects
135  *
136  * \code
137  * auto t_vec = data.getFTensor1FieldData<3>();
138  * \endcode
139  *
140  * @tparam Tensor_Dim size of vector
141  * @return FTensor::Tensor1<FTensor::PackPtr<double *, Tensor_Dim>,
142  * Tensor_Dim>
143  */
144  template <int Tensor_Dim>
147  std::stringstream s;
148  s << "Not implemented for this dimension dim = " << Tensor_Dim;
149  THROW_MESSAGE(s.str());
150  }
151 
152  /**
153  * @brief Return FTensor rank 2, i.e. matrix from filed data coeffinects
154  *
155  * \code
156  * auto t_mat = data.getFTensor2FieldData<3,3>();
157  * \endcode
158  *
159  * @tparam Tensor_Dim0
160  * @tparam Tensor_Dim1
161  * @return FTensor::Tensor2<FTensor::PackPtr<double *, Tensor_Dim0 * Tensor_Dim1>, Tensor_Dim0, Tensor_Dim1>
162  */
163  template <int Tensor_Dim0, int Tensor_Dim1>
165  Tensor_Dim0, Tensor_Dim1>
167  std::stringstream s;
168  s << "Not implemented for this dimension dim0 = " << Tensor_Dim0;
169  s << " and dim1 " << Tensor_Dim1;
170  THROW_MESSAGE(s.str());
171  }
172 
173  /**
174  * @brief Return symmetric FTensor rank 2, i.e. matrix from filed data coeffinects
175  *
176  * \code
177  * auto t_mat = data.getFTensor2SymmetricFieldData<3>();
178  * \endcode
179  *
180  * @tparam Tensor_Dim dimension of the tensor
181  * @return FTensor::Tensor2_symmetric<FTensor::PackPtr<double *, (Tensor_Dim * (Tensor_Dim + 1)) / 2>, Tensor_Dim>
182  */
183  template <int Tensor_Dim>
185  FTensor::PackPtr<double *, (Tensor_Dim * (Tensor_Dim + 1)) / 2>,
186  Tensor_Dim>
188  std::stringstream s;
189  s << "Not implemented for this dimension dim = " << Tensor_Dim;
190  THROW_MESSAGE(s.str());
191  }
192 
193  /**
194  * @brief Resturn scalar files as a FTensor of rank 0
195  *
196  * @return FTensor::Tensor0<FTensor::PackPtr<double *,1> >
197  */
198  FTensor::Tensor0<FTensor::PackPtr<double *,1> > getFTensor0FieldData();
199 
200  inline VectorDofs &getFieldDofs() { return dOfs; }
201 
202  /**@}*/
203 
204  /** \name Base and space */
205 
206  /**@{*/
207 
208  /**
209  * \brief Get approximation base
210  * @return Approximation base
211  */
212  inline FieldApproximationBase &getBase() { return bAse; }
213 
214  /**
215  * \brief Get field space
216  * @return Field space
217  */
218  inline FieldSpace &getSpace() { return sPace; }
219 
220  /**
221  * Get shared pointer to base base functions
222  */
223  virtual boost::shared_ptr<MatrixDouble> &
225  return N[base];
226  }
227 
228  /**
229  * Get shared pointer to base base functions
230  */
231  virtual const boost::shared_ptr<MatrixDouble> &
233  return N[base];
234  }
235 
236  /**
237  * Get shared pointer to derivatives of base base functions
238  */
239  virtual boost::shared_ptr<MatrixDouble> &
241  return diffN[base];
242  }
243 
244  /**
245  * Get shared pointer to derivatives of base base functions
246  */
247  virtual const boost::shared_ptr<MatrixDouble> &
249  return diffN[base];
250  }
251 
252  /**@}*/
253 
254  /** \name Get base functions for H1/L2 */
255 
256  /**@{*/
257 
258  /** \brief get base functions
259  * this return matrix (nb. of rows is equal to nb. of Gauss pts, nb. of
260  * columns is equal to number of base functions on this entity.
261  *
262  * Note that for vectorial base, like Hdiv or Hcurl, in columns are
263  * vectorial base functions. For tonsorial would be tonsorial base
264  * functions. Interpretation depends on type of base, scalar, vectorial or
265  * tonsorial and dimension fo problem.
266  *
267  */
268  virtual const MatrixDouble &getN(const FieldApproximationBase base) const {
269  return *(getNSharedPtr(base));
270  }
271 
272  inline const MatrixDouble &getN() const { return getN(bAse); }
273 
274  /** \brief get derivatives of base functions
275  *
276  * Matrix at rows has nb. of Gauss pts, at columns it has derivative of
277  * base functions. Columns are structured as follows, [ dN1/dx, dN1/dy,
278  * dN1/dz, dN2/dx, dN2/dy, dN2/dz, ... ]
279  *
280  * Scalar base functions:
281  * Note that base functions are calculated in file H1.c
282  * Above description not apply for derivatives of nodal functions, since
283  * derivative of nodal functions in case of simplexes, EDGES, TRIANGLES and
284  * TETS are constant. So that matrix rows represents nb. of base
285  * functions, columns are derivatives. Nb. of columns depend on element
286  * dimension, for EDGES is one, for TRIS is 2 and TETS is 3.
287  *
288  * Note that for node element this function make no sense.
289  *
290  * Tonsorial base functions:
291  * Note: In rows ale integration pts, columns are formatted that that
292  * components of vectors and then derivatives, for example row for given
293  * integration points is formatted in array
294  * \f[
295  * t_{0,0}, t_{1,0}, t_{1,0}, t_{0,1}, t_{1,1}, t_{1,1}, t_{0,2}, t_{1,2},
296  * t_{1,2} \f] where comma express derivative, i.e. \f$t_{2,1} =
297  * \frac{\partial t_2}{\partial \xi_1}\f$
298  *
299  */
300  virtual const MatrixDouble &
301  getDiffN(const FieldApproximationBase base) const {
302  return *(getDiffNSharedPtr(base));
303  }
304 
305  inline const MatrixDouble &getDiffN() const { return getDiffN(bAse); }
306 
307  /**
308  * \brief Get base functions
309  * @param base Approximation base
310  * @return Error code
311  */
313  return *(getNSharedPtr(base));
314  }
315 
316  /**
317  * \brief Get base functions
318  *
319  * It assumed that approximation base for given field is known and stored in
320  * this data structure
321  *
322  * @return Error code
323  */
324  inline MatrixDouble &getN() { return getN(bAse); }
325 
326  /**
327  * \brief Get derivatives of base functions
328  * @param base Approximation base
329  * @return Error code
330  */
332  return *(getDiffNSharedPtr(base));
333  }
334 
335  /**
336  * \brief Get derivatives of base functions
337  *
338  * It assumed that approximation base for given field is known and stored in
339  * this data structure
340  *
341  * @return Error code
342  */
343  inline MatrixDouble &getDiffN() { return getDiffN(bAse); }
344 
345  /// \brief get base functions at Gauss pts
346  inline const VectorAdaptor getN(const FieldApproximationBase base,
347  const int gg) {
348  int size = getN(base).size2();
349  double *data = &getN(base)(gg, 0);
350  return VectorAdaptor(size,
351  ublas::shallow_array_adaptor<double>(size, data));
352  }
353 
354  /// \brief get base functions at Gauss pts
355  inline const VectorAdaptor getN(const int gg) { return getN(bAse, gg); }
356 
357  /** \brief get derivative of base functions at Gauss pts
358 
359  * returned matrix on rows has base functions, in column its derivatives.
360  *
361  * \param base Approximation base
362  * \param gg Nb. of Gauss pts.
363  *
364  */
366  const int gg) {
367  // FIXME: That is bug, it will not work if number of integration pts is
368  // equal to number of nodes on entity. User who not implementing low level
369  // DataOperator will not experience this.
370  if (getN(base).size1() == getDiffN(base).size1()) {
371  int size = getN(base).size2();
372  int dim = getDiffN(base).size2() / size;
373  double *data = &getDiffN(base)(gg, 0);
374  return MatrixAdaptor(
375  getN(base).size2(), dim,
376  ublas::shallow_array_adaptor<double>(getDiffN(base).size2(), data));
377  } else {
378  // in some cases, f.e. for derivatives of nodal base functions at only
379  // one gauss point is needed
380  return MatrixAdaptor(
381  getN(base).size1(), getN(base).size2(),
382  ublas::shallow_array_adaptor<double>(getDiffN(base).data().size(),
383  &getDiffN(base).data()[0]));
384  }
385  }
386 
387  /** \brief get derivative of base functions at Gauss pts
388 
389  * returned matrix on rows has base functions, in column its derivatives.
390  *
391  * \param gg nb. of Gauss pts.
392  *
393  */
394  inline const MatrixAdaptor getDiffN(const int gg) {
395  return getDiffN(bAse, gg);
396  }
397 
398  /** \brief get base functions at Gauss pts
399 
400  * Note that multi field element, two different field can have different
401  * approximation orders. Since we use hierarchical approximation basis,
402  * base functions are calculated once for element, using maximal
403  * approximation order on given entity.
404  *
405  * \param base Approximation base
406  * \param gg number of Gauss point
407  * \param nb_base_functions number of of base functions returned
408 
409  */
410  inline const VectorAdaptor getN(const FieldApproximationBase base,
411  const int gg, const int nb_base_functions) {
412  (void)getN()(gg, nb_base_functions -
413  1); // throw error if nb_base_functions is to big
414  double *data = &getN(base)(gg, 0);
415  return VectorAdaptor(
416  nb_base_functions,
417  ublas::shallow_array_adaptor<double>(nb_base_functions, data));
418  }
419 
420  /** \brief get base functions at Gauss pts
421 
422  * Note that multi field element, two different field can have different
423  * approximation orders. Since we use hierarchical approximation basis,
424  * base functions are calculated once for element, using maximal
425  * approximation order on given entity.
426  *
427  * \param gg number of Gauss point
428  * \param nb_base_functions number of of base functions returned
429 
430  */
431  inline const VectorAdaptor getN(const int gg, const int nb_base_functions) {
432  return getN(bAse, gg, nb_base_functions);
433  }
434 
435  /** \brief get derivatives of base functions at Gauss pts
436  *
437  * Note that multi field element, two different field can have different
438  * approximation orders. Since we use hierarchical approximation basis,
439  * base functions are calculated once for element, using maximal
440  * approximation order on given entity.
441  *
442  * \param base Approximation base
443  * \param gg nb. of Gauss point
444  * \param nb_base_functions number of of base functions
445  *
446  */
448  const int gg,
449  const int nb_base_functions) {
450  // FIXME: That is bug, it will not work if number of integration pts is
451  // equal to number of nodes on entity. User who not implementing low level
452  // DataOperator will not experience this.
453  if (getN(base).size1() == getDiffN(base).size1()) {
454  (void)getN(base)(gg,
455  nb_base_functions -
456  1); // throw error if nb_base_functions is to big
457  int dim = getDiffN(base).size2() / getN(base).size2();
458  double *data = &getDiffN(base)(gg, 0);
459  return MatrixAdaptor(nb_base_functions, dim,
460  ublas::shallow_array_adaptor<double>(
461  dim * nb_base_functions, data));
462  } else {
463  // in some cases, f.e. for derivatives of nodal base functions only one
464  // gauss point is needed
465  return MatrixAdaptor(
466  getN(base).size1(), getN(base).size2(),
467  ublas::shallow_array_adaptor<double>(getDiffN(base).data().size(),
468  &getDiffN(base).data()[0]));
469  }
470  }
471 
472  /** \brief get derivatives of base functions at Gauss pts
473  *
474  * Note that multi field element, two different field can have different
475  * approximation orders. Since we use hierarchical approximation basis,
476  * base functions are calculated once for element, using maximal
477  * approximation order on given entity.
478  *
479  * \param gg nb. of Gauss point
480  * \param nb_base_functions number of of base functions
481  *
482  */
483  inline const MatrixAdaptor getDiffN(const int gg,
484  const int nb_base_functions) {
485  return getDiffN(bAse, gg, nb_base_functions);
486  }
487 
488  /**@}*/
489 
490  /** \name Get base functions for vectorial approximation basese, i.e.
491  * Hdiv/Hcurl */
492 
493  /**@{*/
494 
495  /** \brief get Hdiv of base functions at Gauss pts
496  *
497  * \param base Approximation base
498  * \param gg nb. of Gauss point
499  *
500  */
501  template <int DIM>
503  const int gg) {
504  if (PetscUnlikely(getN(base).size2() % DIM)) {
505  THROW_MESSAGE("Wrong dimension");
506  }
507 
508  const int nb_base_functions = getN(base).size2() / DIM;
509  double *data = &getN(base)(gg, 0);
510  return MatrixAdaptor(
511  nb_base_functions, DIM,
512  ublas::shallow_array_adaptor<double>(DIM * nb_base_functions, data));
513  }
514 
515  /** \brief get Hdiv of base functions at Gauss pts
516  *
517  * \param gg nb. of Gauss point
518  * \param number of of base functions
519  *
520  */
521  template<int DIM>
522  inline const MatrixAdaptor getVectorN(const int gg) {
523  return getVectorN<DIM>(bAse, gg);
524  }
525 
526  /** \brief get DiffHdiv of base functions at Gauss pts
527  *
528  * \param base Approximation base
529  * \param gg nb. of Gauss point
530  * \param number of of base functions
531  *
532  */
533  template <int DIM0, int DIM1>
535  const int gg) {
536  if (PetscUnlikely(getDiffN(base).size2() % (DIM0*DIM1))) {
537  THROW_MESSAGE("Wrong dimension");
538  }
539 
540  const int nb_base_functions = getN(base).size2() / (DIM0 * DIM1);
541  double *data = &getN(base)(gg, 0);
542  return MatrixAdaptor(nb_base_functions, DIM0 * DIM1,
543  ublas::shallow_array_adaptor<double>(
544  DIM0 * DIM1 * nb_base_functions, data));
545  }
546 
547  /** \brief get DiffHdiv of base functions at Gauss pts
548  *
549  * \param gg nb. of Gauss point
550  * \param number of of base functions
551  *
552  */
553  template <int DIM0, int DIM1>
554  inline const MatrixAdaptor getVectorDiffN(const int gg) {
555  return getVectorDiffN<DIM0,DIM1>(bAse, gg);
556  }
557 
558  /** \brief get DiffHdiv of base functions at Gauss pts
559  *
560  * \param base Approximation base
561  * \param gg nb. of Gauss point
562  * \param number of of base functions
563  *
564  */
565  template <int DIM0, int DIM1>
567  const int dof, const int gg) {
568  double *data = &getDiffN(base)(gg, DIM0 * DIM1 * dof);
569  return MatrixAdaptor(
570  DIM0, DIM1, ublas::shallow_array_adaptor<double>(DIM0 * DIM1, data));
571  }
572 
573  /** \brief get DiffHdiv of base functions at Gauss pts
574  *
575  * \param gg nb. of Gauss point
576  * \param number of of base functions
577  *
578  */
579  template <int DIM0, int DIM1>
580  inline const MatrixAdaptor getVectorDiffN(const int dof, const int gg) {
581  return getVectorDiffN<DIM0,DIM1>(bAse, dof, gg);
582  }
583 
584  /**@}*/
585 
586  /** \name Get base functions with FTensor */
587 
588  /**@{*/
589 
590  /**
591  * \brief Get base function as Tensor0
592  *
593  * \param base
594  * \return Tensor0
595  *
596  */
599  double *ptr = &*getN(base).data().begin();
601  };
602 
603  /**
604  * \brief Get base function as Tensor0
605  *
606  * Return base functions for field base
607  *
608  * \return Tensor0
609  *
610  */
612  return getFTensor0N(bAse);
613  };
614 
615  /**
616  * \brief Get base function as Tensor0 (Loop by integration points)
617  *
618  * \param base
619  * \param bb base function
620  * \return Tensor0
621 
622  Note that:
623  \code
624  t0 = data.getFTensor0N(base,bb);
625  ++t0
626  \endcode
627  Increment in above code will move pointer to base function in next
628  integration point.
629 
630  *
631  */
633  getFTensor0N(const FieldApproximationBase base, const int bb) {
634  double *ptr = &getN(base)(0, bb);
635  return FTensor::Tensor0<double *>(ptr, getN(base).size2());
636  };
637 
638  /**
639  * \brief Get base function as Tensor0 (Loop by integration points)
640  *
641  * Return base functions for field base
642  *
643  * \param bb base function
644  * \return Tensor0
645  *
646  *
647  */
649  return getFTensor0N(bAse, bb);
650  };
651 
652  /**
653  * \brief Get base function as Tensor0 (Loop by integration points)
654  *
655  * \param base
656  * \param gg integration points
657  * \param bb base function
658  * \return Tensor0
659 
660  Note that:
661  \code
662  t0 = data.getFTensor0N(base,bb);
663  ++t0
664  \endcode
665  Increment in above code will move pointer to base function in next
666  integration point.
667 
668  *
669  */
671  getFTensor0N(const FieldApproximationBase base, const int gg,
672  const int bb) {
673  double *ptr = &getN(base)(gg, bb);
675  };
676 
677  /**
678  * \brief Get base function as Tensor0 (Loop by integration points)
679  *
680  * Return base functions for field base
681  *
682  * \param bb base function
683  * \return Tensor0
684  *
685  */
687  getFTensor0N(const int gg, const int bb) {
688  return getFTensor0N(bAse, gg, bb);
689  };
690 
691  /**
692  * \brief Get derivatives of base functions
693  *
694  * For volume element like tetrahedral or prism,
695  * \code
696  * Tensor1<double*,3> diff_base = data.getFTensor1DiffN<3>();
697  * \endcode
698  *
699  * For face element like triangle or quad
700  * \code
701  * Tensor1<double*,2> diff_base = data.getFTensor1DiffN<2>();
702  * \endcode
703  *
704  * \param base functions
705  * \return Tensor rank 1 (vector)
706  *
707  */
708  template <int Tensor_Dim>
710  getFTensor1DiffN(const FieldApproximationBase base);
711 
712  /**
713  * \brief Get derivatives of base functions
714  *
715  * For volume element like tetrahedral or prism,
716  * \code
717  * Tensor1<double*,3> diff_base = data.getFTensor1DiffN<3>();
718  * \endcode
719  *
720  * For face element like triangle or quad
721  * \code
722  * Tensor1<double*,2> diff_base = data.getFTensor1DiffN<2>();
723  * \endcode
724  *
725  * \return Tensor rank 1 (vector)
726  *
727  */
728  template <int Tensor_Dim>
729  FTensor::Tensor1<double *, Tensor_Dim> getFTensor1DiffN();
730 
731  /**
732  * \brief Get derivatives of base functions (Loop by integration points)
733  *
734  * For volume element like tetrahedral or prism,
735  * \code
736  * Tensor1<double*,3> diff_base = data.getFTensor1DiffN<3>(base,bb);
737  * \endcode
738  * where bb is base function. Operator ++diff_base will move tensor pointer
739  * to next integration point.
740  *
741  * For face element like triangle or quad
742  * \code
743  * Tensor1<double*,2> diff_base = data.getFTensor1DiffN<2>(base,bb);
744  * \endcode
745  *
746  * \param base functions
747  * \return Tensor rank 1 (vector)
748  *
749  */
750  template <int Tensor_Dim>
752  getFTensor1DiffN(const FieldApproximationBase base, const int bb);
753 
754  /**
755  * \brief Get derivatives of base functions (Loop by integration points)
756  *
757  * For volume element like tetrahedral or prism,
758  * \code
759  * Tensor1<double*,3> diff_base = data.getFTensor1DiffN<3>(bb);
760  * \endcode
761  * where bb is base function. Operator ++diff_base will move tensor pointer
762  * to next integration point.
763  *
764  * For face element like triangle or quad
765  * \code
766  * Tensor1<double*,2> diff_base = data.getFTensor1DiffN<2>(bb);
767  * \endcode
768  *
769  * \return Tensor rank 1 (vector)
770  *
771  */
772  template <int Tensor_Dim>
773  FTensor::Tensor1<double *, Tensor_Dim> getFTensor1DiffN(const int bb);
774 
775  /**
776  * \brief Get derivatives of base functions (Loop by integration points)
777  *
778  * For volume element like tetrahedral or prism,
779  * \code
780  * Tensor1<double*,3> diff_base = data.getFTensor1DiffN<3>(base,gg,bb);
781  * \endcode
782  * where bb is base function and gg is integration pt. Operator ++diff_base
783  * will move tensor pointer to next integration point.
784  *
785  * For face element like triangle or quad
786  * \code
787  * Tensor1<double*,2> diff_base = data.getFTensor1DiffN<2>(base,gg,bb);
788  * \endcode
789  *
790  * \return Tensor rank 1 (vector)
791  *
792  */
793  template <int Tensor_Dim>
795  getFTensor1DiffN(const FieldApproximationBase base, const int gg,
796  const int bb);
797 
798  /**
799  * \brief Get derivatives of base functions (Loop by integration points)
800  *
801  * For volume element like tetrahedral or prism,
802  * \code
803  * Tensor1<double*,3> diff_base = data.getFTensor1DiffN<3>(gg,bb);
804  * \endcode
805  * where bb is base function and gg is integration pt. Operator ++diff_base
806  * will move tensor pointer to next base function.
807  *
808  * For face element like triangle or quad
809  * \code
810  * Tensor1<double*,2> diff_base = data.getFTensor1DiffN<2>(gg,bb);
811  * \endcode
812  *
813  * \return Tensor rank 1 (vector)
814  *
815  */
816  template <int Tensor_Dim>
817  FTensor::Tensor1<double *, Tensor_Dim> getFTensor1DiffN(const int gg,
818  const int bb);
819 
820  /** \brief Get base functions for Hdiv/Hcurl spaces
821 
822  \note You probably like to use getFTensor1N(), in typical use base is
823  set automatically based on base set to field.
824 
825  * @param base Approximation base
826 
827  Example:
828  \code
829  FTensor::Index<'i',3> i;
830  int nb_dofs = data.getFieldData().size();
831  auto t_n_hdiv = data.getFTensor1N<3>();
832  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
833  int ll = 0;
834  for(;ll!=nb_dofs;ll++) {
835  double dot_product = t_n_hdiv(i)*t_n_hdiv(i);
836  ++t_n_hdiv;
837  }
838  for(;ll!=data.getVectorN().size2()/3;ll++) {
839  ++t_n_hdiv;
840  }
841  }
842  \endcode
843 
844  */
845  template <int Tensor_Dim>
847  getFTensor1N(FieldApproximationBase base);
848 
849  /** \brief Get base functions for Hdiv space
850 
851  Example:
852  \code
853  FTensor::Index<'i',3> i;
854  int nb_dofs = data.getFieldData().size();
855  auto t_n_hdiv = data.getFTensor1N<3>();
856  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
857  int ll = 0;
858  for(;ll!=nb_dofs;ll++) {
859  double dot_product = t_n_hdiv(i)*t_n_hdiv(i);
860  ++t_n_hdiv;
861  }
862  for(;ll!=data.getVectorN().size2()/3;ll++) {
863  ++t_n_hdiv;
864  }
865  }
866  \endcode
867 
868  */
869  template <int Tensor_Dim> auto getFTensor1N() {
870  return getFTensor1N<Tensor_Dim>(bAse);
871  }
872 
873  /** \brief Get derivatives of base functions for Hdiv space
874  */
875  template <int Tensor_Dim0, int Tensor_Dim1>
877  Tensor_Dim0, Tensor_Dim1>
878  getFTensor2DiffN(FieldApproximationBase base);
879 
880  /** \brief Get derivatives of base functions for Hdiv space at integration
881  * pts
882  */
883  template <int Tensor_Dim0, int Tensor_Dim1>
885  Tensor_Dim0, Tensor_Dim1>
886  getFTensor2DiffN(FieldApproximationBase base, const int gg,
887  const int bb);
888 
889  /** \brief Get derivatives of base functions for Hdiv space
890  */
891  template <int Tensor_Dim0, int Tensor_Dim1>
892  inline FTensor::Tensor2<
894  Tensor_Dim1>
896  return getFTensor2DiffN<Tensor_Dim0, Tensor_Dim1>(bAse);
897  }
898 
899  /** \brief Get derivatives of base functions for Hdiv space at integration
900  * pts
901  */
902  template <int Tensor_Dim0, int Tensor_Dim1>
903  inline FTensor::Tensor2<
905  Tensor_Dim1>
906  getFTensor2DiffN(const int gg, const int bb) {
907  return getFTensor2DiffN<Tensor_Dim0, Tensor_Dim1>(bAse, gg, bb);
908  }
909 
910  /**
911  * \brief Get Hdiv base functions at integration point
912 
913  \code
914  FTensor::Index<'i',3> i;
915  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
916  auto t_base = data.getFTensor1N(base,gg,bb);
917  for(int bb = 0;bb!=nb_base_functions;bb++) {
918  auto dot = t_base(i)*t_base(i);
919  }
920  }
921  \endcode
922 
923  */
924  template <int Tensor_Dim>
926  getFTensor1N(FieldApproximationBase base, const int gg, const int bb);
927 
928  /**
929  * \brief Get Hdiv base functions at integration point
930 
931  \code
932  FTensor::Index<'i',3> i;
933  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
934  auto t_base = data.getFTensor1N(gg,0);
935  for(int bb = 0;bb!=nb_base_functions;bb++) {
936  double dot = t_base(i)*t_base(i);
937  }
938  }
939  \endcode
940 
941  */
942  template <int Tensor_Dim>
943  inline auto getFTensor1N(const int gg, const int bb) {
944  return getFTensor1N<Tensor_Dim>(bAse, gg, bb);
945  }
946 
947  /** \brief Get base functions for Hdiv/Hcurl spaces
948 
949  \note You probably like to use getFTensor1N(), in typical use base is
950  set automatically based on base set to field.
951 
952  * @param base Approximation base
953 
954  Example:
955  \code
956  FTensor::Index<'i',3> i;
957  FTensor::Index<'i',3> j;
958  int nb_dofs = data.getFieldData().size();
959  auto t_n_hdiv = data.getFTensor2N<3,3>();
960  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
961  int ll = 0;
962  for(;ll!=nb_dofs;ll++) {
963  double dot_product = t_n_hdiv(i,j)*t_n_hdiv(i,j);
964  ++t_n_hdiv;
965  }
966  for(;ll!=data.getVectorN().size2()/3;ll++) {
967  ++t_n_hdiv;
968  }
969  }
970  \endcode
971 
972  */
973  template <int Tensor_Dim0, int Tensor_Dim1>
975  Tensor_Dim0, Tensor_Dim1>
976  getFTensor2N(FieldApproximationBase base);
977 
978  /** \brief Get base functions for Hdiv space
979 
980  Example:
981  \code
982  FTensor::Index<'i',3> i;
983  FTensor::Index<'j',3> j;
984 
985  int nb_dofs = data.getFieldData().size();
986  auto t_n_hdiv = data.getFTensor2N<3,3>();
987  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
988  int ll = 0;
989  for(;ll!=nb_dofs;ll++) {
990  double dot_product = t_n_hdiv(i,j)*t_n_hdiv(i,j);
991  ++t_n_hdiv;
992  }
993  for(;ll!=data.getVectorN().size2()/3;ll++) {
994  ++t_n_hdiv;
995  }
996  }
997  \endcode
998 
999  */
1000  template <int Tensor_Dim0, int Tensor_Dim1> auto getFTensor2N() {
1001  return getFTensor2N<Tensor_Dim0, Tensor_Dim1>(bAse);
1002  }
1003 
1004  /** \brief Get base functions for tensor Hdiv/Hcurl spaces
1005 
1006  \note You probably like to use getFTensor2N(), in typical use base is
1007  set automatically based on base set to field.
1008 
1009  @param base Approximation base
1010 
1011  Example:
1012  \code
1013  FTensor::Index<'i',3> i;
1014  FTensor::Index<'j',3> i;
1015  int nb_dofs = data.getFieldData().size();
1016  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
1017  auto t_n_hdiv = data.getFTensor2N<3>(base,gg,bb);
1018  int ll = 0;
1019  for(;ll!=nb_dofs;ll++) {
1020  double dot_product = t_n_hdiv(i,j)*t_n_hdiv(i,j);
1021  ++t_n_hdiv;
1022  }
1023  for(;ll!=data.getVectorN().size2()/3;ll++) {
1024  ++t_n_hdiv;
1025  }
1026  }
1027  \endcode
1028 
1029  */
1030  template <int Tensor_Dim0, int Tensor_Dim1>
1032  Tensor_Dim0, Tensor_Dim1>
1033  getFTensor2N(FieldApproximationBase base, const int gg, const int bb);
1034 
1035  /** \brief Get base functions for Hdiv space
1036 
1037  Example:
1038  \code
1039  FTensor::Index<'i',3> i;
1040  FTensor::Index<'j',3> j;
1041  int nb_dofs = data.getFieldData().size();
1042  for(int gg = 0;gg!=nb_gauss_pts;++gg) {
1043  int ll = 0;
1044  auto t_n_hdiv = data.getFTensor2N<3,3>(gg,0);
1045  for(;ll!=nb_dofs;ll++) {
1046  double dot_product = t_n_hdiv(i)*t_n_hdiv(i);
1047  ++t_n_hdiv;
1048  }
1049  for(;ll!=data.getVectorN().size2()/3;ll++) {
1050  ++t_n_hdiv;
1051  }
1052  }
1053  \endcode
1054 
1055  */
1056  template <int Tensor_Dim0, int Tensor_Dim1>
1057  auto getFTensor2N(const int gg, const int bb) {
1058  return getFTensor2N<Tensor_Dim0, Tensor_Dim1>(bAse, gg, bb);
1059  }
1060 
1061  /**@}*/
1062 
1063  /** \name Auxiliary functions */
1064 
1065  /**@{*/
1066 
1067  friend std::ostream &operator<<(std::ostream &os,
1069 
1070  /**
1071  * Reset data associated with particular field name
1072  * @return error code
1073  */
1076  sPace = NOSPACE;
1077  bAse = NOBASE;
1078  iNdices.resize(0, false);
1079  localIndices.resize(0, false);
1080  dOfs.resize(0, false);
1081  fieldData.resize(0, false);
1083  }
1084 
1085  /**@}*/
1086 
1087  protected:
1088  int sEnse; ///< Entity sense (orientation)
1089  ApproximationOrder oRder; ///< Entity order
1090  FieldSpace sPace; ///< Entity space
1091  FieldApproximationBase bAse; ///< Field approximation base
1092  VectorInt iNdices; ///< Global indices on entity
1093  VectorInt localIndices; ///< Local indices on entity
1094  VectorDofs dOfs; ///< DoFs on entity
1095  VectorDouble fieldData; ///< Field data on entity
1096  boost::shared_ptr<MatrixDouble> N[LASTBASE]; ///< Base functions
1097  boost::shared_ptr<MatrixDouble>
1098  diffN[LASTBASE]; ///< Derivatives of base functions
1099  };
1100 
1101  std::bitset<LASTSPACE> sPace; ///< spaces on element
1102  std::bitset<LASTBASE> bAse; ///< bases on element
1103  ublas::matrix<int> facesNodes; ///< nodes on finite element faces
1104  std::bitset<LASTSPACE>
1105  spacesOnEntities[MBMAXTYPE]; ///< spaces on entity types
1106  std::bitset<LASTBASE> basesOnEntities[MBMAXTYPE]; ///< bases on entity types
1107  std::bitset<LASTBASE> basesOnSpaces[LASTSPACE]; ///< base on spaces
1108 
1109  boost::ptr_vector<EntData> dataOnEntities[MBMAXTYPE]; ///< data on nodes, base
1110  ///< function, dofs
1111  ///< values, etc.
1112 
1113  /**
1114  * Reset data associated with particular field name
1115  * @return error code
1116  */
1119  for (EntityType t = MBVERTEX; t != MBMAXTYPE; t++) {
1120  for (auto &e : dataOnEntities[t]) {
1121  ierr = e.resetFieldDependentData();
1122  CHKERRG(ierr);
1123  }
1124  }
1126  }
1127 
1128  DataForcesAndSourcesCore(const EntityType type);
1130 
1131  virtual MoFEMErrorCode setElementType(const EntityType type);
1132 
1133  friend std::ostream &operator<<(std::ostream &os,
1134  const DataForcesAndSourcesCore &e);
1135 protected:
1137 };
1138 
1139 /** \brief this class derive data form other data structure
1140  * \ingroup mofem_forces_and_sources_user_data_operators
1141  *
1142  *
1143  * It behaves like normal data structure it is used to share base functions with
1144  * other data structures. Dofs values, approx. order and
1145  * indices are not shared.
1146  *
1147  * Shape functions, senses are shared with other data structure.
1148  *
1149  */
1151 
1152  /** \brief Derived ata on single entity (This is passed as argument to
1153  * DataOperator::doWork) \ingroup mofem_forces_and_sources_user_data_operators
1154  * \nosubgrouping
1155  *
1156  * DerivedEntData share part information with EntData except infomation about
1157  * base functions.
1158  *
1159  */
1161 
1162  const boost::shared_ptr<DataForcesAndSourcesCore::EntData> entDataPtr;
1163  DerivedEntData(const boost::shared_ptr<DataForcesAndSourcesCore::EntData>
1164  &ent_data_ptr)
1165  : entDataPtr(ent_data_ptr) {}
1166 
1167  int getSense() const { return entDataPtr->getSense(); }
1168 
1169  boost::shared_ptr<MatrixDouble> &
1171  return entDataPtr->getNSharedPtr(base);
1172  }
1173  boost::shared_ptr<MatrixDouble> &
1175  return entDataPtr->getDiffNSharedPtr(base);
1176  }
1177  const boost::shared_ptr<MatrixDouble> &
1179  return entDataPtr->getNSharedPtr(base);
1180  }
1181  const boost::shared_ptr<MatrixDouble> &
1183  return entDataPtr->getDiffNSharedPtr(base);
1184  }
1185 
1186  };
1187 
1189  const boost::shared_ptr<DataForcesAndSourcesCore> &data_ptr);
1190  MoFEMErrorCode setElementType(const EntityType type);
1191 
1192 private:
1193  const boost::shared_ptr<DataForcesAndSourcesCore> dataPtr;
1194 };
1195 
1196 /** \name Specializations for H1/L2 */
1197 
1198 /**@{*/
1199 
1200 template <>
1202 DataForcesAndSourcesCore::EntData::getFTensor1FieldData<3>();
1203 
1204 template <>
1206 DataForcesAndSourcesCore::EntData::getFTensor1FieldData<2>();
1207 
1208 template <>
1210 DataForcesAndSourcesCore::EntData::getFTensor2FieldData<3, 3>();
1211 
1212 template <>
1214 DataForcesAndSourcesCore::EntData::getFTensor2SymmetricFieldData<3>();
1215 
1216 template <>
1218 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<3>(
1219  const FieldApproximationBase base);
1220 template <>
1222 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<3>();
1223 
1224 template <>
1226 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<2>(
1227  const FieldApproximationBase base);
1228 template <>
1230 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<2>();
1231 template <>
1233 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<3>(
1234  const FieldApproximationBase base, const int bb);
1235 template <>
1237 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<3>(const int bb);
1238 template <>
1240 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<2>(
1241  const FieldApproximationBase base, const int bb);
1242 template <>
1244 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<2>(const int bb);
1245 
1246 /**@}*/
1247 
1248 /** \name Specializations for HDiv/HCurl in 2d */
1249 
1250 /**@{*/
1251 
1252 template <>
1254 DataForcesAndSourcesCore::EntData::getFTensor2DiffN<3, 2>(
1255  FieldApproximationBase base);
1256 template <>
1258 DataForcesAndSourcesCore::EntData::getFTensor2DiffN<3, 2>(
1259  FieldApproximationBase base, const int gg, const int bb);
1260 
1261 /**@}*/
1262 
1263 /// \deprecated Use DataForcesAndSourcesCore
1265 
1266 /// \deprecated use DerivedDataForcesAndSourcesCore
1269 
1270 } // namespace MoFEM
1271 
1272 #endif //__DATASTRUCTURES_HPP
1273 
1274 /**
1275  * \defgroup mofem_forces_and_sources_user_data_operators User data operator
1276  * data structures \ingroup
1277  *
1278  * \brief Users data structures and operator
1279  *
1280  * Data structures passed by argument to MoFEM::DataOperator::doWork and generic
1281  * user operators operating on those structures.
1282  *
1283  */
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)
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
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:499
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)
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:77
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
const VectorDofs & getFieldDofs() const
get dofs data stature FEDofEntity
VectorInt iNdices
Global indices on entity.
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:542
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:618
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:506
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
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
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:173
const boost::shared_ptr< DataForcesAndSourcesCore::EntData > entDataPtr
virtual const MatrixDouble & getDiffN(const FieldApproximationBase base) const
get derivatives of base functions
auto getVectorAdaptor
Get Vector adaptor.
Definition: Types.hpp:121
std::ostream & operator<<(std::ostream &os, const DataForcesAndSourcesCore::EntData &e)
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.
FieldApproximationBase
approximation base
Definition: definitions.h:141
const MatrixAdaptor getVectorDiffN(const int gg)
get DiffHdiv of base functions at Gauss pts
const boost::shared_ptr< DataForcesAndSourcesCore > dataPtr
MatrixShallowArrayAdaptor< double > MatrixAdaptor
Matrix adaptor.
Definition: Types.hpp:141
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.
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)
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
const MatrixDouble & getDiffN() const
const MatrixAdaptor getVectorDiffN(const int dof, const int gg)
get DiffHdiv of base functions at Gauss pts
VectorShallowArrayAdaptor< double > VectorAdaptor
Definition: Types.hpp:99
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:166
DEPRECATED typedef DataForcesAndSourcesCore DataForcesAndSurcesCore
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< int, IntAllocator > VectorInt
Definition: Types.hpp:75
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)
VectorShallowArrayAdaptor< int > VectorIntAdaptor
Definition: Types.hpp:100
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)
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:76
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.