v0.8.23
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
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
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
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 *
162  * Tensor_Dim1>, Tensor_Dim0, Tensor_Dim1>
163  */
164  template <int Tensor_Dim0, int Tensor_Dim1>
166  Tensor_Dim0, Tensor_Dim1>
168  std::stringstream s;
169  s << "Not implemented for this dimension dim0 = " << Tensor_Dim0;
170  s << " and dim1 " << Tensor_Dim1;
171  THROW_MESSAGE(s.str());
172  }
173 
174  /**
175  * @brief Return symmetric FTensor rank 2, i.e. matrix from filed data
176  * coeffinects
177  *
178  * \code
179  * auto t_mat = data.getFTensor2SymmetricFieldData<3>();
180  * \endcode
181  *
182  * @tparam Tensor_Dim dimension of the tensor
183  * @return FTensor::Tensor2_symmetric<FTensor::PackPtr<double *, (Tensor_Dim
184  * * (Tensor_Dim + 1)) / 2>, Tensor_Dim>
185  */
186  template <int Tensor_Dim>
188  FTensor::PackPtr<double *, (Tensor_Dim * (Tensor_Dim + 1)) / 2>,
189  Tensor_Dim>
191  std::stringstream s;
192  s << "Not implemented for this dimension dim = " << Tensor_Dim;
193  THROW_MESSAGE(s.str());
194  }
195 
196  /**
197  * @brief Resturn scalar files as a FTensor of rank 0
198  *
199  * @return FTensor::Tensor0<FTensor::PackPtr<double *,1> >
200  */
201  FTensor::Tensor0<FTensor::PackPtr<double *, 1>> getFTensor0FieldData();
202 
203  inline VectorDofs &getFieldDofs() { return dOfs; }
204 
205  /**@}*/
206 
207  /** \name Base and space */
208 
209  /**@{*/
210 
211  /**
212  * \brief Get approximation base
213  * @return Approximation base
214  */
215  inline FieldApproximationBase &getBase() { return bAse; }
216 
217  /**
218  * \brief Get field space
219  * @return Field space
220  */
221  inline FieldSpace &getSpace() { return sPace; }
222 
223  /**
224  * Get shared pointer to base base functions
225  */
226  virtual boost::shared_ptr<MatrixDouble> &
228  return N[base];
229  }
230 
231  /**
232  * Get shared pointer to base base functions
233  */
234  virtual const boost::shared_ptr<MatrixDouble> &
236  return N[base];
237  }
238 
239  /**
240  * Get shared pointer to derivatives of base base functions
241  */
242  virtual boost::shared_ptr<MatrixDouble> &
244  return diffN[base];
245  }
246 
247  /**
248  * Get shared pointer to derivatives of base base functions
249  */
250  virtual const boost::shared_ptr<MatrixDouble> &
252  return diffN[base];
253  }
254 
255  /**@}*/
256 
257  /** \name Get base functions for H1/L2 */
258 
259  /**@{*/
260 
261  /** \brief get base functions
262  * this return matrix (nb. of rows is equal to nb. of Gauss pts, nb. of
263  * columns is equal to number of base functions on this entity.
264  *
265  * Note that for vectorial base, like Hdiv or Hcurl, in columns are
266  * vectorial base functions. For tonsorial would be tonsorial base
267  * functions. Interpretation depends on type of base, scalar, vectorial or
268  * tonsorial and dimension fo problem.
269  *
270  */
271  virtual const MatrixDouble &getN(const FieldApproximationBase base) const {
272  return *(getNSharedPtr(base));
273  }
274 
275  inline const MatrixDouble &getN() const { return getN(bAse); }
276 
277  /** \brief get derivatives of base functions
278  *
279  * Matrix at rows has nb. of Gauss pts, at columns it has derivative of
280  * base functions. Columns are structured as follows, [ dN1/dx, dN1/dy,
281  * dN1/dz, dN2/dx, dN2/dy, dN2/dz, ... ]
282  *
283  * Scalar base functions:
284  * Note that base functions are calculated in file H1.c
285  * Above description not apply for derivatives of nodal functions, since
286  * derivative of nodal functions in case of simplexes, EDGES, TRIANGLES and
287  * TETS are constant. So that matrix rows represents nb. of base
288  * functions, columns are derivatives. Nb. of columns depend on element
289  * dimension, for EDGES is one, for TRIS is 2 and TETS is 3.
290  *
291  * Note that for node element this function make no sense.
292  *
293  * Tonsorial base functions:
294  * Note: In rows ale integration pts, columns are formatted that that
295  * components of vectors and then derivatives, for example row for given
296  * integration points is formatted in array
297  * \f[
298  * t_{0,0}, t_{1,0}, t_{1,0}, t_{0,1}, t_{1,1}, t_{1,1}, t_{0,2}, t_{1,2},
299  * t_{1,2} \f] where comma express derivative, i.e. \f$t_{2,1} =
300  * \frac{\partial t_2}{\partial \xi_1}\f$
301  *
302  */
303  virtual const MatrixDouble &
304  getDiffN(const FieldApproximationBase base) const {
305  return *(getDiffNSharedPtr(base));
306  }
307 
308  inline const MatrixDouble &getDiffN() const { return getDiffN(bAse); }
309 
310  /**
311  * \brief Get base functions
312  * @param base Approximation base
313  * @return Error code
314  */
316  return *(getNSharedPtr(base));
317  }
318 
319  /**
320  * \brief Get base functions
321  *
322  * It assumed that approximation base for given field is known and stored in
323  * this data structure
324  *
325  * @return Error code
326  */
327  inline MatrixDouble &getN() { return getN(bAse); }
328 
329  /**
330  * \brief Get derivatives of base functions
331  * @param base Approximation base
332  * @return Error code
333  */
335  return *(getDiffNSharedPtr(base));
336  }
337 
338  /**
339  * \brief Get derivatives of base functions
340  *
341  * It assumed that approximation base for given field is known and stored in
342  * this data structure
343  *
344  * @return Error code
345  */
346  inline MatrixDouble &getDiffN() { return getDiffN(bAse); }
347 
348  /// \brief get base functions at Gauss pts
349  inline const VectorAdaptor getN(const FieldApproximationBase base,
350  const int gg) {
351  int size = getN(base).size2();
352  double *data = &getN(base)(gg, 0);
353  return VectorAdaptor(size,
354  ublas::shallow_array_adaptor<double>(size, data));
355  }
356 
357  /// \brief get base functions at Gauss pts
358  inline const VectorAdaptor getN(const int gg) { return getN(bAse, gg); }
359 
360  /** \brief get derivative of base functions at Gauss pts
361 
362  * returned matrix on rows has base functions, in column its derivatives.
363  *
364  * \param base Approximation base
365  * \param gg Nb. of Gauss pts.
366  *
367  */
369  const int gg) {
370  // FIXME: That is bug, it will not work if number of integration pts is
371  // equal to number of nodes on entity. User who not implementing low
372  // level DataOperator will not experience this.
373  if (getN(base).size1() == getDiffN(base).size1()) {
374  int size = getN(base).size2();
375  int dim = getDiffN(base).size2() / size;
376  double *data = &getDiffN(base)(gg, 0);
377  return MatrixAdaptor(
378  getN(base).size2(), dim,
379  ublas::shallow_array_adaptor<double>(getDiffN(base).size2(), data));
380  } else {
381  // in some cases, f.e. for derivatives of nodal base functions at only
382  // one gauss point is needed
383  return MatrixAdaptor(
384  getN(base).size1(), getN(base).size2(),
385  ublas::shallow_array_adaptor<double>(getDiffN(base).data().size(),
386  &getDiffN(base).data()[0]));
387  }
388  }
389 
390  /** \brief get derivative of base functions at Gauss pts
391 
392  * returned matrix on rows has base functions, in column its derivatives.
393  *
394  * \param gg nb. of Gauss pts.
395  *
396  */
397  inline const MatrixAdaptor getDiffN(const int gg) {
398  return getDiffN(bAse, gg);
399  }
400 
401  /** \brief get base functions at Gauss pts
402 
403  * Note that multi field element, two different field can have different
404  * approximation orders. Since we use hierarchical approximation basis,
405  * base functions are calculated once for element, using maximal
406  * approximation order on given entity.
407  *
408  * \param base Approximation base
409  * \param gg number of Gauss point
410  * \param nb_base_functions number of of base functions returned
411 
412  */
413  inline const VectorAdaptor getN(const FieldApproximationBase base,
414  const int gg, const int nb_base_functions) {
415  (void)getN()(gg, nb_base_functions -
416  1); // throw error if nb_base_functions is to big
417  double *data = &getN(base)(gg, 0);
418  return VectorAdaptor(
419  nb_base_functions,
420  ublas::shallow_array_adaptor<double>(nb_base_functions, data));
421  }
422 
423  /** \brief get base functions at Gauss pts
424 
425  * Note that multi field element, two different field can have different
426  * approximation orders. Since we use hierarchical approximation basis,
427  * base functions are calculated once for element, using maximal
428  * approximation order on given entity.
429  *
430  * \param gg number of Gauss point
431  * \param nb_base_functions number of of base functions returned
432 
433  */
434  inline const VectorAdaptor getN(const int gg, const int nb_base_functions) {
435  return getN(bAse, gg, nb_base_functions);
436  }
437 
438  /** \brief get derivatives of base functions at Gauss pts
439  *
440  * Note that multi field element, two different field can have different
441  * approximation orders. Since we use hierarchical approximation basis,
442  * base functions are calculated once for element, using maximal
443  * approximation order on given entity.
444  *
445  * \param base Approximation base
446  * \param gg nb. of Gauss point
447  * \param nb_base_functions number of of base functions
448  *
449  */
451  const int gg,
452  const int nb_base_functions) {
453  // FIXME: That is bug, it will not work if number of integration pts is
454  // equal to number of nodes on entity. User who not implementing low
455  // level DataOperator will not experience this.
456  if (getN(base).size1() == getDiffN(base).size1()) {
457  (void)getN(base)(gg,
458  nb_base_functions -
459  1); // throw error if nb_base_functions is to big
460  int dim = getDiffN(base).size2() / getN(base).size2();
461  double *data = &getDiffN(base)(gg, 0);
462  return MatrixAdaptor(nb_base_functions, dim,
463  ublas::shallow_array_adaptor<double>(
464  dim * nb_base_functions, data));
465  } else {
466  // in some cases, f.e. for derivatives of nodal base functions only one
467  // gauss point is needed
468  return MatrixAdaptor(
469  getN(base).size1(), getN(base).size2(),
470  ublas::shallow_array_adaptor<double>(getDiffN(base).data().size(),
471  &getDiffN(base).data()[0]));
472  }
473  }
474 
475  /** \brief get derivatives of base functions at Gauss pts
476  *
477  * Note that multi field element, two different field can have different
478  * approximation orders. Since we use hierarchical approximation basis,
479  * base functions are calculated once for element, using maximal
480  * approximation order on given entity.
481  *
482  * \param gg nb. of Gauss point
483  * \param nb_base_functions number of of base functions
484  *
485  */
486  inline const MatrixAdaptor getDiffN(const int gg,
487  const int nb_base_functions) {
488  return getDiffN(bAse, gg, nb_base_functions);
489  }
490 
491  /**@}*/
492 
493  /** \name Get base functions for vectorial approximation basese, i.e.
494  * Hdiv/Hcurl */
495 
496  /**@{*/
497 
498  /** \brief get Hdiv of base functions at Gauss pts
499  *
500  * \param base Approximation base
501  * \param gg nb. of Gauss point
502  *
503  */
504  template <int DIM>
506  const int gg) {
507  if (PetscUnlikely(getN(base).size2() % DIM)) {
508  THROW_MESSAGE("Wrong dimension");
509  }
510 
511  const int nb_base_functions = getN(base).size2() / DIM;
512  double *data = &getN(base)(gg, 0);
513  return MatrixAdaptor(
514  nb_base_functions, DIM,
515  ublas::shallow_array_adaptor<double>(DIM * nb_base_functions, data));
516  }
517 
518  /** \brief get Hdiv of base functions at Gauss pts
519  *
520  * \param gg nb. of Gauss point
521  * \param number of of base functions
522  *
523  */
524  template <int DIM> inline const MatrixAdaptor getVectorN(const int gg) {
525  return getVectorN<DIM>(bAse, gg);
526  }
527 
528  /** \brief get DiffHdiv of base functions at Gauss pts
529  *
530  * \param base Approximation base
531  * \param gg nb. of Gauss point
532  * \param number of of base functions
533  *
534  */
535  template <int DIM0, int DIM1>
537  const int gg) {
538  if (PetscUnlikely(getDiffN(base).size2() % (DIM0 * DIM1))) {
539  THROW_MESSAGE("Wrong dimension");
540  }
541 
542  const int nb_base_functions = getN(base).size2() / (DIM0 * DIM1);
543  double *data = &getN(base)(gg, 0);
544  return MatrixAdaptor(nb_base_functions, DIM0 * DIM1,
545  ublas::shallow_array_adaptor<double>(
546  DIM0 * DIM1 * nb_base_functions, data));
547  }
548 
549  /** \brief get DiffHdiv of base functions at Gauss pts
550  *
551  * \param gg nb. of Gauss point
552  * \param number of of base functions
553  *
554  */
555  template <int DIM0, int DIM1>
556  inline const MatrixAdaptor getVectorDiffN(const int gg) {
557  return getVectorDiffN<DIM0, DIM1>(bAse, gg);
558  }
559 
560  /** \brief get DiffHdiv of base functions at Gauss pts
561  *
562  * \param base Approximation base
563  * \param gg nb. of Gauss point
564  * \param number of of base functions
565  *
566  */
567  template <int DIM0, int DIM1>
569  const int dof, const int gg) {
570  double *data = &getDiffN(base)(gg, DIM0 * DIM1 * dof);
571  return MatrixAdaptor(
572  DIM0, DIM1, ublas::shallow_array_adaptor<double>(DIM0 * DIM1, data));
573  }
574 
575  /** \brief get DiffHdiv of base functions at Gauss pts
576  *
577  * \param gg nb. of Gauss point
578  * \param number of of base functions
579  *
580  */
581  template <int DIM0, int DIM1>
582  inline const MatrixAdaptor getVectorDiffN(const int dof, const int gg) {
583  return getVectorDiffN<DIM0, DIM1>(bAse, dof, gg);
584  }
585 
586  /**@}*/
587 
588  /** \name Get base functions with FTensor */
589 
590  /**@{*/
591 
592  /**
593  * \brief Get base function as Tensor0
594  *
595  * \param base
596  * \return Tensor0
597  *
598  */
601  double *ptr = &*getN(base).data().begin();
603  };
604 
605  /**
606  * \brief Get base function as Tensor0
607  *
608  * Return base functions for field base
609  *
610  * \return Tensor0
611  *
612  */
614  return getFTensor0N(bAse);
615  };
616 
617  /**
618  * \brief Get base function as Tensor0 (Loop by integration points)
619  *
620  * \param base
621  * \param bb base function
622  * \return Tensor0
623 
624  Note that:
625  \code
626  t0 = data.getFTensor0N(base,bb);
627  ++t0
628  \endcode
629  Increment in above code will move pointer to base function in next
630  integration point.
631 
632  *
633  */
635  getFTensor0N(const FieldApproximationBase base, const int bb) {
636  double *ptr = &getN(base)(0, bb);
637  return FTensor::Tensor0<double *>(ptr, getN(base).size2());
638  };
639 
640  /**
641  * \brief Get base function as Tensor0 (Loop by integration points)
642  *
643  * Return base functions for field base
644  *
645  * \param bb base function
646  * \return Tensor0
647  *
648  *
649  */
651  return getFTensor0N(bAse, bb);
652  };
653 
654  /**
655  * \brief Get base function as Tensor0 (Loop by integration points)
656  *
657  * \param base
658  * \param gg integration points
659  * \param bb base function
660  * \return Tensor0
661 
662  Note that:
663  \code
664  t0 = data.getFTensor0N(base,bb);
665  ++t0
666  \endcode
667  Increment in above code will move pointer to base function in next
668  integration point.
669 
670  *
671  */
673  getFTensor0N(const FieldApproximationBase base, const int gg,
674  const int bb) {
675  double *ptr = &getN(base)(gg, bb);
677  };
678 
679  /**
680  * \brief Get base function as Tensor0 (Loop by integration points)
681  *
682  * Return base functions for field base
683  *
684  * \param bb base function
685  * \return Tensor0
686  *
687  */
689  getFTensor0N(const int gg, const int bb) {
690  return getFTensor0N(bAse, gg, bb);
691  };
692 
693  /**
694  * \brief Get derivatives of base functions
695  *
696  * For volume element like tetrahedral or prism,
697  * \code
698  * Tensor1<double*,3> diff_base = data.getFTensor1DiffN<3>();
699  * \endcode
700  *
701  * For face element like triangle or quad
702  * \code
703  * Tensor1<double*,2> diff_base = data.getFTensor1DiffN<2>();
704  * \endcode
705  *
706  * \param base functions
707  * \return Tensor rank 1 (vector)
708  *
709  */
710  template <int Tensor_Dim>
712  getFTensor1DiffN(const FieldApproximationBase base);
713 
714  /**
715  * \brief Get derivatives of base functions
716  *
717  * For volume element like tetrahedral or prism,
718  * \code
719  * Tensor1<double*,3> diff_base = data.getFTensor1DiffN<3>();
720  * \endcode
721  *
722  * For face element like triangle or quad
723  * \code
724  * Tensor1<double*,2> diff_base = data.getFTensor1DiffN<2>();
725  * \endcode
726  *
727  * \return Tensor rank 1 (vector)
728  *
729  */
730  template <int Tensor_Dim>
731  FTensor::Tensor1<double *, Tensor_Dim> getFTensor1DiffN();
732 
733  /**
734  * \brief Get derivatives of base functions (Loop by integration points)
735  *
736  * For volume element like tetrahedral or prism,
737  * \code
738  * Tensor1<double*,3> diff_base = data.getFTensor1DiffN<3>(base,bb);
739  * \endcode
740  * where bb is base function. Operator ++diff_base will move tensor pointer
741  * to next integration point.
742  *
743  * For face element like triangle or quad
744  * \code
745  * Tensor1<double*,2> diff_base = data.getFTensor1DiffN<2>(base,bb);
746  * \endcode
747  *
748  * \param base functions
749  * \return Tensor rank 1 (vector)
750  *
751  */
752  template <int Tensor_Dim>
754  getFTensor1DiffN(const FieldApproximationBase base, const int bb);
755 
756  /**
757  * \brief Get derivatives of base functions (Loop by integration points)
758  *
759  * For volume element like tetrahedral or prism,
760  * \code
761  * Tensor1<double*,3> diff_base = data.getFTensor1DiffN<3>(bb);
762  * \endcode
763  * where bb is base function. Operator ++diff_base will move tensor pointer
764  * to next integration point.
765  *
766  * For face element like triangle or quad
767  * \code
768  * Tensor1<double*,2> diff_base = data.getFTensor1DiffN<2>(bb);
769  * \endcode
770  *
771  * \return Tensor rank 1 (vector)
772  *
773  */
774  template <int Tensor_Dim>
775  FTensor::Tensor1<double *, Tensor_Dim> getFTensor1DiffN(const int bb);
776 
777  /**
778  * \brief Get derivatives of base functions (Loop by integration points)
779  *
780  * For volume element like tetrahedral or prism,
781  * \code
782  * Tensor1<double*,3> diff_base = data.getFTensor1DiffN<3>(base,gg,bb);
783  * \endcode
784  * where bb is base function and gg is integration pt. Operator ++diff_base
785  * will move tensor pointer to next integration point.
786  *
787  * For face element like triangle or quad
788  * \code
789  * Tensor1<double*,2> diff_base = data.getFTensor1DiffN<2>(base,gg,bb);
790  * \endcode
791  *
792  * \return Tensor rank 1 (vector)
793  *
794  */
795  template <int Tensor_Dim>
797  getFTensor1DiffN(const FieldApproximationBase base, const int gg,
798  const int bb);
799 
800  /**
801  * \brief Get derivatives of base functions (Loop by integration points)
802  *
803  * For volume element like tetrahedral or prism,
804  * \code
805  * Tensor1<double*,3> diff_base = data.getFTensor1DiffN<3>(gg,bb);
806  * \endcode
807  * where bb is base function and gg is integration pt. Operator ++diff_base
808  * will move tensor pointer to next base function.
809  *
810  * For face element like triangle or quad
811  * \code
812  * Tensor1<double*,2> diff_base = data.getFTensor1DiffN<2>(gg,bb);
813  * \endcode
814  *
815  * \return Tensor rank 1 (vector)
816  *
817  */
818  template <int Tensor_Dim>
819  FTensor::Tensor1<double *, Tensor_Dim> getFTensor1DiffN(const int gg,
820  const int bb);
821 
822  /** \brief Get base functions for Hdiv/Hcurl spaces
823 
824  \note You probably like to use getFTensor1N(), in typical use base is
825  set automatically based on base set to field.
826 
827  * @param base Approximation base
828 
829  Example:
830  \code
831  FTensor::Index<'i',3> i;
832  int nb_dofs = data.getFieldData().size();
833  auto t_n_hdiv = data.getFTensor1N<3>();
834  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
835  int ll = 0;
836  for(;ll!=nb_dofs;ll++) {
837  double dot_product = t_n_hdiv(i)*t_n_hdiv(i);
838  ++t_n_hdiv;
839  }
840  for(;ll!=data.getVectorN().size2()/3;ll++) {
841  ++t_n_hdiv;
842  }
843  }
844  \endcode
845 
846  */
847  template <int Tensor_Dim>
849  getFTensor1N(FieldApproximationBase base);
850 
851  /** \brief Get base functions for Hdiv space
852 
853  Example:
854  \code
855  FTensor::Index<'i',3> i;
856  int nb_dofs = data.getFieldData().size();
857  auto t_n_hdiv = data.getFTensor1N<3>();
858  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
859  int ll = 0;
860  for(;ll!=nb_dofs;ll++) {
861  double dot_product = t_n_hdiv(i)*t_n_hdiv(i);
862  ++t_n_hdiv;
863  }
864  for(;ll!=data.getVectorN().size2()/3;ll++) {
865  ++t_n_hdiv;
866  }
867  }
868  \endcode
869 
870  */
871  template <int Tensor_Dim> auto getFTensor1N() {
872  return getFTensor1N<Tensor_Dim>(bAse);
873  }
874 
875  /** \brief Get derivatives of base functions for Hdiv space
876  */
877  template <int Tensor_Dim0, int Tensor_Dim1>
879  Tensor_Dim0, Tensor_Dim1>
880  getFTensor2DiffN(FieldApproximationBase base);
881 
882  /** \brief Get derivatives of base functions for Hdiv space at integration
883  * pts
884  */
885  template <int Tensor_Dim0, int Tensor_Dim1>
887  Tensor_Dim0, Tensor_Dim1>
888  getFTensor2DiffN(FieldApproximationBase base, const int gg, const int bb);
889 
890  /** \brief Get derivatives of base functions for Hdiv space
891  */
892  template <int Tensor_Dim0, int Tensor_Dim1>
893  inline FTensor::Tensor2<
895  Tensor_Dim1>
897  return getFTensor2DiffN<Tensor_Dim0, Tensor_Dim1>(bAse);
898  }
899 
900  /** \brief Get derivatives of base functions for Hdiv space at integration
901  * pts
902  */
903  template <int Tensor_Dim0, int Tensor_Dim1>
904  inline FTensor::Tensor2<
906  Tensor_Dim1>
907  getFTensor2DiffN(const int gg, const int bb) {
908  return getFTensor2DiffN<Tensor_Dim0, Tensor_Dim1>(bAse, gg, bb);
909  }
910 
911  /**
912  * \brief Get Hdiv base functions at integration point
913 
914  \code
915  FTensor::Index<'i',3> i;
916  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
917  auto t_base = data.getFTensor1N(base,gg,bb);
918  for(int bb = 0;bb!=nb_base_functions;bb++) {
919  auto dot = t_base(i)*t_base(i);
920  }
921  }
922  \endcode
923 
924  */
925  template <int Tensor_Dim>
927  getFTensor1N(FieldApproximationBase base, const int gg, const int bb);
928 
929  /**
930  * \brief Get Hdiv base functions at integration point
931 
932  \code
933  FTensor::Index<'i',3> i;
934  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
935  auto t_base = data.getFTensor1N(gg,0);
936  for(int bb = 0;bb!=nb_base_functions;bb++) {
937  double dot = t_base(i)*t_base(i);
938  }
939  }
940  \endcode
941 
942  */
943  template <int Tensor_Dim>
944  inline auto getFTensor1N(const int gg, const int bb) {
945  return getFTensor1N<Tensor_Dim>(bAse, gg, bb);
946  }
947 
948  /** \brief Get base functions for Hdiv/Hcurl spaces
949 
950  \note You probably like to use getFTensor1N(), in typical use base is
951  set automatically based on base set to field.
952 
953  * @param base Approximation base
954 
955  Example:
956  \code
957  FTensor::Index<'i',3> i;
958  FTensor::Index<'i',3> j;
959  int nb_dofs = data.getFieldData().size();
960  auto t_n_hdiv = data.getFTensor2N<3,3>();
961  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
962  int ll = 0;
963  for(;ll!=nb_dofs;ll++) {
964  double dot_product = t_n_hdiv(i,j)*t_n_hdiv(i,j);
965  ++t_n_hdiv;
966  }
967  for(;ll!=data.getVectorN().size2()/3;ll++) {
968  ++t_n_hdiv;
969  }
970  }
971  \endcode
972 
973  */
974  template <int Tensor_Dim0, int Tensor_Dim1>
976  Tensor_Dim0, Tensor_Dim1>
977  getFTensor2N(FieldApproximationBase base);
978 
979  /** \brief Get base functions for Hdiv space
980 
981  Example:
982  \code
983  FTensor::Index<'i',3> i;
984  FTensor::Index<'j',3> j;
985 
986  int nb_dofs = data.getFieldData().size();
987  auto t_n_hdiv = data.getFTensor2N<3,3>();
988  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
989  int ll = 0;
990  for(;ll!=nb_dofs;ll++) {
991  double dot_product = t_n_hdiv(i,j)*t_n_hdiv(i,j);
992  ++t_n_hdiv;
993  }
994  for(;ll!=data.getVectorN().size2()/3;ll++) {
995  ++t_n_hdiv;
996  }
997  }
998  \endcode
999 
1000  */
1001  template <int Tensor_Dim0, int Tensor_Dim1> auto getFTensor2N() {
1002  return getFTensor2N<Tensor_Dim0, Tensor_Dim1>(bAse);
1003  }
1004 
1005  /** \brief Get base functions for tensor Hdiv/Hcurl spaces
1006 
1007  \note You probably like to use getFTensor2N(), in typical use base is
1008  set automatically based on base set to field.
1009 
1010  @param base Approximation base
1011 
1012  Example:
1013  \code
1014  FTensor::Index<'i',3> i;
1015  FTensor::Index<'j',3> i;
1016  int nb_dofs = data.getFieldData().size();
1017  for(int gg = 0;gg!=nb_gauss_pts;gg++) {
1018  auto t_n_hdiv = data.getFTensor2N<3>(base,gg,bb);
1019  int ll = 0;
1020  for(;ll!=nb_dofs;ll++) {
1021  double dot_product = t_n_hdiv(i,j)*t_n_hdiv(i,j);
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_Dim0, int Tensor_Dim1>
1033  Tensor_Dim0, Tensor_Dim1>
1034  getFTensor2N(FieldApproximationBase base, const int gg, const int bb);
1035 
1036  /** \brief Get base functions for Hdiv space
1037 
1038  Example:
1039  \code
1040  FTensor::Index<'i',3> i;
1041  FTensor::Index<'j',3> j;
1042  int nb_dofs = data.getFieldData().size();
1043  for(int gg = 0;gg!=nb_gauss_pts;++gg) {
1044  int ll = 0;
1045  auto t_n_hdiv = data.getFTensor2N<3,3>(gg,0);
1046  for(;ll!=nb_dofs;ll++) {
1047  double dot_product = t_n_hdiv(i)*t_n_hdiv(i);
1048  ++t_n_hdiv;
1049  }
1050  for(;ll!=data.getVectorN().size2()/3;ll++) {
1051  ++t_n_hdiv;
1052  }
1053  }
1054  \endcode
1055 
1056  */
1057  template <int Tensor_Dim0, int Tensor_Dim1>
1058  auto getFTensor2N(const int gg, const int bb) {
1059  return getFTensor2N<Tensor_Dim0, Tensor_Dim1>(bAse, gg, bb);
1060  }
1061 
1062  /**@}*/
1063 
1064  /** \name Auxiliary functions */
1065 
1066  /**@{*/
1067 
1068  friend std::ostream &operator<<(std::ostream &os,
1070 
1071  /**
1072  * Reset data associated with particular field name
1073  * @return error code
1074  */
1077  sPace = NOSPACE;
1078  bAse = NOBASE;
1079  iNdices.resize(0, false);
1080  localIndices.resize(0, false);
1081  dOfs.resize(0, false);
1082  fieldData.resize(0, false);
1084  }
1085 
1086  /**@}*/
1087 
1088  protected:
1089  int sEnse; ///< Entity sense (orientation)
1090  ApproximationOrder oRder; ///< Entity order
1091  FieldSpace sPace; ///< Entity space
1092  FieldApproximationBase bAse; ///< Field approximation base
1093  VectorInt iNdices; ///< Global indices on entity
1094  VectorInt localIndices; ///< Local indices on entity
1095  VectorDofs dOfs; ///< DoFs on entity
1096  VectorDouble fieldData; ///< Field data on entity
1097  boost::shared_ptr<MatrixDouble> N[LASTBASE]; ///< Base functions
1098  boost::shared_ptr<MatrixDouble>
1099  diffN[LASTBASE]; ///< Derivatives of base functions
1100  };
1101 
1102  std::bitset<LASTSPACE> sPace; ///< spaces on element
1103  std::bitset<LASTBASE> bAse; ///< bases on element
1104  ublas::matrix<int> facesNodes; ///< nodes on finite element faces
1105  std::bitset<LASTSPACE>
1106  spacesOnEntities[MBMAXTYPE]; ///< spaces on entity types
1107  std::bitset<LASTBASE> basesOnEntities[MBMAXTYPE]; ///< bases on entity types
1108  std::bitset<LASTBASE> basesOnSpaces[LASTSPACE]; ///< base on spaces
1109 
1110  boost::ptr_vector<EntData> dataOnEntities[MBMAXTYPE]; ///< data on nodes, base
1111  ///< function, dofs
1112  ///< values, etc.
1113 
1114  /**
1115  * Reset data associated with particular field name
1116  * @return error code
1117  */
1120  for (EntityType t = MBVERTEX; t != MBMAXTYPE; t++) {
1121  for (auto &e : dataOnEntities[t]) {
1122  ierr = e.resetFieldDependentData();
1123  CHKERRG(ierr);
1124  }
1125  }
1127  }
1128 
1129  DataForcesAndSourcesCore(const EntityType type);
1131 
1132  virtual MoFEMErrorCode setElementType(const EntityType type);
1133 
1134  friend std::ostream &operator<<(std::ostream &os,
1135  const DataForcesAndSourcesCore &e);
1136 
1137 protected:
1139 };
1140 
1141 /** \brief this class derive data form other data structure
1142  * \ingroup mofem_forces_and_sources_user_data_operators
1143  *
1144  *
1145  * It behaves like normal data structure it is used to share base functions with
1146  * other data structures. Dofs values, approx. order and
1147  * indices are not shared.
1148  *
1149  * Shape functions, senses are shared with other data structure.
1150  *
1151  */
1153 
1154  /** \brief Derived ata on single entity (This is passed as argument to
1155  * DataOperator::doWork) \ingroup mofem_forces_and_sources_user_data_operators
1156  * \nosubgrouping
1157  *
1158  * DerivedEntData share part information with EntData except infomation about
1159  * base functions.
1160  *
1161  */
1163 
1164  const boost::shared_ptr<DataForcesAndSourcesCore::EntData> entDataPtr;
1165  DerivedEntData(const boost::shared_ptr<DataForcesAndSourcesCore::EntData>
1166  &ent_data_ptr)
1167  : entDataPtr(ent_data_ptr) {}
1168 
1169  int getSense() const { return entDataPtr->getSense(); }
1170 
1171  boost::shared_ptr<MatrixDouble> &
1173  return entDataPtr->getNSharedPtr(base);
1174  }
1175  boost::shared_ptr<MatrixDouble> &
1177  return entDataPtr->getDiffNSharedPtr(base);
1178  }
1179  const boost::shared_ptr<MatrixDouble> &
1181  return entDataPtr->getNSharedPtr(base);
1182  }
1183  const boost::shared_ptr<MatrixDouble> &
1185  return entDataPtr->getDiffNSharedPtr(base);
1186  }
1187  };
1188 
1190  const boost::shared_ptr<DataForcesAndSourcesCore> &data_ptr);
1191  MoFEMErrorCode setElementType(const EntityType type);
1192 
1193 private:
1194  const boost::shared_ptr<DataForcesAndSourcesCore> dataPtr;
1195 };
1196 
1197 /** \name Specializations for H1/L2 */
1198 
1199 /**@{*/
1200 
1201 template <>
1203 DataForcesAndSourcesCore::EntData::getFTensor1FieldData<3>();
1204 
1205 template <>
1207 DataForcesAndSourcesCore::EntData::getFTensor1FieldData<2>();
1208 
1209 template <>
1211 DataForcesAndSourcesCore::EntData::getFTensor2FieldData<3, 3>();
1212 
1213 template <>
1215 DataForcesAndSourcesCore::EntData::getFTensor2SymmetricFieldData<3>();
1216 
1217 template <>
1219 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<3>(
1220  const FieldApproximationBase base);
1221 template <>
1223 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<3>();
1224 
1225 template <>
1227 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<2>(
1228  const FieldApproximationBase base);
1229 template <>
1231 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<2>();
1232 template <>
1234 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<3>(
1235  const FieldApproximationBase base, const int bb);
1236 template <>
1238 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<3>(const int bb);
1239 template <>
1241 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<2>(
1242  const FieldApproximationBase base, const int bb);
1243 template <>
1245 DataForcesAndSourcesCore::EntData::getFTensor1DiffN<2>(const int bb);
1246 
1247 /**@}*/
1248 
1249 /** \name Specializations for HDiv/HCurl in 2d */
1250 
1251 /**@{*/
1252 
1253 template <>
1255 DataForcesAndSourcesCore::EntData::getFTensor2DiffN<3, 2>(
1256  FieldApproximationBase base);
1257 template <>
1259 DataForcesAndSourcesCore::EntData::getFTensor2DiffN<3, 2>(
1260  FieldApproximationBase base, const int gg, const int bb);
1261 
1262 /**@}*/
1263 
1264 /// \deprecated Use DataForcesAndSourcesCore
1266 
1267 /**
1268  * @brief Assemble PETSc vector
1269  *
1270  * Function extract indices from entity data and assemble vector
1271  *
1272  * <a
1273  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Vec/VecSetValues.html>See
1274  * PETSc documentation</a>
1275  *
1276  * @param V
1277  * @param data
1278  * @param ptr
1279  * @param iora
1280  * @return MoFEMErrorCode
1281  */
1282 inline MoFEMErrorCode
1284  const double *ptr, InsertMode iora) {
1285  return VecSetValues(V, data.getIndices().size(), &*data.getIndices().begin(),
1286  ptr, iora);
1287 }
1288 
1289 /**
1290  * @brief Assemble PETSc matrix
1291  *
1292  * Function extract indices from entity data and assemble vector
1293  *
1294  * <a
1295  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatSetValues.html>See
1296  * PETSc documentation</a>
1297  *
1298  * @param M
1299  * @param row_data
1300  * @param col_data
1301  * @param ptr
1302  * @param iora
1303  * @return MoFEMErrorCode
1304  */
1305 inline MoFEMErrorCode
1307  const DataForcesAndSourcesCore::EntData &col_data,
1308  const double *ptr, InsertMode iora) {
1309  return MatSetValues(
1310  M, row_data.getIndices().size(), &*row_data.getIndices().begin(),
1311  col_data.getIndices().size(), &*col_data.getIndices().begin(), ptr, iora);
1312 }
1313 
1314 } // namespace MoFEM
1315 
1316 #endif //__DATASTRUCTURES_HPP
1317 
1318 /**
1319  * \defgroup mofem_forces_and_sources_user_data_operators User data operator
1320  * data structures \ingroup
1321  *
1322  * \brief Users data structures and operator
1323  *
1324  * Data structures passed by argument to MoFEM::DataOperator::doWork and generic
1325  * user operators operating on those structures.
1326  *
1327  */
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:500
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:543
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:619
virtual const boost::shared_ptr< MatrixDouble > & getDiffNSharedPtr(const FieldApproximationBase base) const
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:507
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....
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:174
const boost::shared_ptr< DataForcesAndSourcesCore::EntData > entDataPtr
virtual const MatrixDouble & getDiffN(const FieldApproximationBase base) const
get derivatives of base functions
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:142
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:121
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
DataForcesAndSourcesCore::EntData EntData
MoFEMErrorCode MatSetValues(Mat M, const DataForcesAndSourcesCore::EntData &row_data, const DataForcesAndSourcesCore::EntData &col_data, const double *ptr, InsertMode iora)
Assemble PETSc matrix.
auto getVectorAdaptor(T1 ptr, const size_t n)
Get Vector adaptor.
Definition: Templates.hpp:43
VectorShallowArrayAdaptor< double > VectorAdaptor
Definition: Types.hpp:104
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.
ApproximationOrder getOrder() const
get approximation order
FieldSpace
approximation spaces
Definition: definitions.h:167
MoFEMErrorCode VecSetValues(Vec V, const DataForcesAndSourcesCore::EntData &data, const double *ptr, InsertMode iora)
Assemble PETSc vector.
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:105
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
const int order
approximation order
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base)
Get base function as Tensor0.