v0.14.0
Loading...
Searching...
No Matches
ForcesAndSourcesCore.hpp
Go to the documentation of this file.
1/** \file ForcesAndSourcesCore.hpp
2 \brief Implementation of elements on entities.
3
4 Those element are inherited by user to implement specific implementation of
5 particular problem.
6
7*/
8
9
10
11#ifndef __FORCES_AND_SOURCES_CORE__HPP__
12#define __FORCES_AND_SOURCES_CORE__HPP__
13
14using namespace boost::numeric;
15
16namespace MoFEM {
17
18/** \brief structure to get information form mofem into EntitiesFieldData
19 * \ingroup mofem_forces_and_sources
20 *
21 */
23
25
27 typedef boost::function<int(int order_row, int order_col, int order_data)>
29
30 typedef boost::function<MoFEMErrorCode(ForcesAndSourcesCore *fe_raw_ptr,
31 int order_row, int order_col,
32 int order_data)>
34
35 /**
36 * \brief Hook to get rule
37 *
38 * \todo check preferred format how works with gcc and clang,
39 * see
40 * <http://www.boost.org/doc/libs/1_64_0/doc/html/function/tutorial.html#idp247873024>
41 */
43
44 /**
45 * @brief Set function to calculate integration rule
46 *
47 */
49
50 /** \brief Data operator to do calculations at integration points.
51 * \ingroup mofem_forces_and_sources
52
53 Is inherited and implemented by user to do calculations. It can be used in
54 many different ways but typically is used to integrate matrices (f.e.
55 stiffness matrix) and the right hand vector (f.e. force vector).
56
57 Note: It is assumed that operator is executed for symmetric problem. That
58 means that is executed for not repeating entities on finite element. For
59 example on triangle we have nodes, 3 edges and one face. Because of symmetry
60 calculations are for: nodes-nodes, nodes-edge0, nodes-edge_1, nodes-edge_2,
61 nodes-face,
62 edges_1-edges_1, edges_1-edge_1, edge_1-edge_2,
63 edge_1-edge_1, edge_1-edge_2,
64 edge_2-edge_2,
65 edge_1-face, edge_1-face, edges_3-face,
66 face - face
67
68 In case of non symmetric problem in addition calculations of lower off
69 diagonal terms. F.e. edge_1-edge_0, esges_3-edge_0, edge_3-edge_1,
70
71 In that case class variable UserDataOperator::sYmm = false;
72
73 NoteL: By default sYmm is set for symmetric problems
74
75 */
76 struct UserDataOperator;
77
78 /** \brief Use to push back operator for row operator
79
80 It can be used to calculate nodal forces or other quantities on the mesh.
81
82 */
83 boost::ptr_deque<UserDataOperator> &getOpPtrVector() { return opPtrVector; }
84
85 /**
86 * @brief Get the Entity Polynomial Base object
87 *
88 * @return boost::shared_ptr<BaseFunction>&&
89 */
91
92 /**
93 * @brief Get the User Polynomial Base object
94 *
95 * @return boost::shared_ptr<BaseFunction>&
96 */
98
99 /**
100 * @brief Matrix of integration points
101 *
102 * Columns is equal to number of integration points, numver of rows
103 * depends on dimension of finite element entity, for example for
104 * tetrahedron rows are x,y,z,weight. Last row is integration weight.
105 *
106 * FIXME: that should be moved to private class data and acessed only by
107 * member function
108 */
110
111 virtual MoFEMErrorCode preProcess();
112 virtual MoFEMErrorCode operator()();
113 virtual MoFEMErrorCode postProcess();
114
115public:
116 /** \brief Get max order of approximation for data fields
117
118 getMeasure getMaxDataOrder () return maximal order on entities, for
119 all data on the element. So for example if finite element is triangle, and
120 triangle base function have order 4 and on edges base function have order
121 2, this function return 4.
122
123 If finite element has for example 2 or more approximated fields, for
124 example Pressure (order 3) and displacement field (order 5), this function
125 returns 5.
126
127 */
128 int getMaxDataOrder() const;
129
130 /// \brief Get max order of approximation for field in rows
131 int getMaxRowOrder() const;
132
133 /// \brief Get max order of approximation for field in columns
134 int getMaxColOrder() const;
135
136 /**
137 * @brief Get the entity data
138 *
139 * @param space
140 * @param type
141 * @param side
142 * @return EntitiesFieldData::EntData&
143 */
144 auto &getEntData(const FieldSpace space, const EntityType type,
145 const int side) {
146 return dataOnElement[space]->dataOnEntities[type][side];
147 }
148
149 /**
150 * @brief Get data on entities and space
151 *
152 * Entities data are stored by space, by entity type, and entity side.
153 *
154 * @return std::array<boost::shared_ptr<EntitiesFieldData>, LASTSPACE>
155 */
157
158 /**
159 * @brief Get derived data on entities and space
160 *
161 * Entities data are stored by space, by entity type, and entity side. Derived
162 * data is used to store data on columns, so it shares infromatin about shape
163 * functions wih rows.
164 *
165 * @return std::array<boost::shared_ptr<EntitiesFieldData>, LASTSPACE>
166 */
168
169protected:
170 /**
171 * \brief get sense (orientation) of entity
172 * @param type type of entity
173 * @param data entity data
174 * @return error code
175 */
177 const EntityType type,
178 boost::ptr_vector<EntitiesFieldData::EntData> &data) const;
179
180 /**
181 * @brief Get the entity data order
182 *
183 * @param type
184 * @param space
185 * @param data
186 * @return MoFEMErrorCode
187 */
189 const EntityType type, const FieldSpace space,
190 boost::ptr_vector<EntitiesFieldData::EntData> &data) const;
191
192 /**
193 * @brief Get the entity sense (orientation)
194 *
195 * @tparam type
196 * @param data
197 * @return MoFEMErrorCode
198 */
199 template <EntityType type>
201 return getEntitySense(type, data.dataOnEntities[type]);
202 }
203
204 /**
205 * @brief Get the entity data order for given space
206 *
207 * @tparam type
208 * @param data
209 * @param space
210 * @return MoFEMErrorCode
211 */
212 template <EntityType type>
214 const FieldSpace space) const {
215 return getEntityDataOrder(type, space, data.dataOnEntities[type]);
216 }
217
218 /** \name Indices */
219
220 /**@{*/
221
222 /// \brief get node indices
223 template <typename EXTRACTOR>
225 getNodesIndices(const int bit_number,
226 FieldEntity_vector_view &ents_field, VectorInt &nodes_indices,
227 VectorInt &local_nodes_indices, EXTRACTOR &&extractor) const;
228
229 /// \brief get row node indices from FENumeredDofEntity_multiIndex
231 const int bit_number) const;
232
233 /// \brief get col node indices from FENumeredDofEntity_multiIndex
235 const int bit_number) const;
236
237 template <typename EXTRACTOR>
238 MoFEMErrorCode getEntityIndices(EntitiesFieldData &data, const int bit_number,
239 FieldEntity_vector_view &ents_field,
240 const EntityType type_lo,
241 const EntityType type_hi,
242 EXTRACTOR &&extractor) const;
243
245 getEntityRowIndices(EntitiesFieldData &data, const int bit_number,
246 const EntityType type_lo = MBVERTEX,
247 const EntityType type_hi = MBPOLYHEDRON) const;
248
250 getEntityColIndices(EntitiesFieldData &data, const int bit_number,
251 const EntityType type_lo = MBVERTEX,
252 const EntityType type_hi = MBPOLYHEDRON) const;
253
254 /// \brief get NoField indices
256 getNoFieldIndices(const int bit_number,
257 boost::shared_ptr<FENumeredDofEntity_multiIndex> dofs,
258 VectorInt &nodes_indices) const;
259
260 /// \brief get col NoField indices
262 const int bit_number) const;
263
264 /// \brief get col NoField indices
266 const int bit_number) const;
267
268 /**@}*/
269
270 /** \name Data */
271
272 /**@{*/
273
274 /** Get bit ref level in entities, and set it to data
275 */
277
278 /**
279 * \brief Get field data on nodes
280 */
281 MoFEMErrorCode getNoFieldFieldData(const int bit_number,
282 VectorDouble &ent_field_data,
283 VectorDofs &ent_field_dofs,
284 VectorFieldEntities &ent_field) const;
285
287 const int bit_number) const;
288
289 /**
290 * \brief Get data on nodes
291 * @param data Data structure
292 * @param field_name Field name
293 * @return Error code
294 */
296 const int bit_number) const;
297
299 getEntityFieldData(EntitiesFieldData &data, const int bit_number,
300 const EntityType type_lo = MBVERTEX,
301 const EntityType type_hi = MBPOLYHEDRON) const;
302
303 /**@}*/
304
305 /// \brief Get nodes on faces
307
308 /// \brief Get field approximation space and base on entities
311
312 /** \name Data form NumeredDofEntity_multiIndex */
313
314 /**@{*/
315
316 /// \brief get indices of nodal indices which are declared for problem but
317 /// not this particular element
319 const NumeredDofEntity_multiIndex &dofs,
320 VectorInt &nodes_indices) const;
321
322 /// \brief get indices by type (generic function) which are declared for
323 /// problem but not this particular element
325 const NumeredDofEntity_multiIndex &dofs,
326 EntityType type, int side_number,
327 VectorInt &indices) const;
328
330 VectorInt &nodes_indices) const;
332 EntityType type, int side_number,
333 VectorInt &indices) const;
335 VectorInt &nodes_indices) const;
337 EntityType type, int side_number,
338 VectorInt &indices) const;
339
340 /**@}*/
341
342 /**
343 * \brief another variant of getRule
344 * @param order_row order of base function on row
345 * @param order_col order of base function on columns
346 * @param order_data order of base function approximating data
347 * @return integration rule
348 *
349 * This function is overloaded by the user. The integration rule
350 * is set such that specific operator implemented by the user is
351 * integrated accurately. For example if user implement bilinear operator
352 * \f[
353 * b(u,v) =
354 * \int_\mathcal{T}
355 * \frac{\partial u_i}{\partial x_j}\frac{\partial v_i}{\partial x_j}
356 * \textrm{d}\mathcal{T}
357 * \f]
358 * then if \f$u\f$ and \f$v\f$ are polynomial of given \em order, then
359 * exact integral would be \code int getRule(int order) { return
360 * 2*(order-1); }; \endcode
361 *
362 * The integration points and weights are set appropriately for given
363 * entity type and integration rule from \ref quad.c
364 *
365 * Method \ref ForcesAndSourcesCore::getRule takes at argument takes
366 * maximal polynomial order set on the element on all fields defined on
367 * the element. If a user likes to have more control, another variant of
368 * this function can be called which distinguishing between field orders
369 * on rows, columns and data, the i.e. first argument of a bilinear form,
370 * the second argument of bilinear form and field coefficients on the
371 * element.
372 *
373 * \note If user set rule to -1 or any other negative integer, then method
374 * \ref ForcesAndSourcesCore::setGaussPts is called. In that method user
375 * can implement own (specific) integration method.
376 *
377 * \bug this function should be const
378 */
379 virtual int getRule(int order_row, int order_col, int order_data);
380
381 /** \brief set user specific integration rule
382
383 This function allows for user defined integration rule. The key is to
384 called matrix gaussPts, which is used by other MoFEM procedures. Matrix
385 has number of rows equal to problem dimension plus one, where last index
386 is used to store weight values. %Number of columns is equal to number of
387 integration points.
388
389 \note This function is called if method \ref
390 ForcesAndSourcesCore::getRule is returning integer -1 or any other
391 negative integer.
392
393 User sets
394 \code
395 MatrixDouble gaussPts;
396 \endcode
397 where
398 \code
399 gaussPts.resize(dim+1,nb_gauss_pts);
400 \endcode
401 number rows represents local coordinates of integration points
402 in reference element, where last index in row is for integration weight.
403
404 */
405 virtual MoFEMErrorCode setGaussPts(int order_row, int order_col,
406 int order_data);
407
408 /**
409 * \brief Calculate base functions
410 * @return Error code
411 */
414
415 /**
416 * \brief Calculate base functions
417 * @return Error code
418 */
420
421 /**
422 * @brief Calculate Bernstein-Bezier base
423 *
424 * @return MoFEMErrorCode
425 */
427
428 /**
429 * @brief Create a entity data on element object
430 *
431 * @return MoFEMErrorCode
432 */
434
435 /**
436 * @brief Iterate user data operators
437 *
438 * @return MoFEMErrorCode
439 */
441
442 /**@{*/
443
444 /** \name Deprecated (do not use) */
445
446 /** \deprecated Use getRule(int row_order, int col_order, int data order)
447 */
448 virtual int getRule(int order);
449
450 /** \deprecated setGaussPts(int row_order, int col_order, int data order);
451 */
452 virtual MoFEMErrorCode setGaussPts(int order);
453
454 /**@}*/
455
456 /**
457 * @brief Entity data on element entity rows fields
458 *
459 */
460 const std::array<boost::shared_ptr<EntitiesFieldData>, LASTSPACE>
462
463 /**
464 * @brief Entity data on element entity columns fields
465 *
466 */
467 const std::array<boost::shared_ptr<EntitiesFieldData>, LASTSPACE>
469
475
476 /**
477 * @brief Vector of finite element users data operators
478 *
479 */
480 boost::ptr_deque<UserDataOperator> opPtrVector;
481
482 friend class UserDataOperator;
483
484protected:
485 /**
486 * @brief Last evaluated type of element entity
487 *
488 */
490
491private:
492 /**
493 * @brief Pointer to entity polynomial base
494 *
495 */
496 boost::shared_ptr<BaseFunction> elementPolynomialBasePtr;
497
498 /**
499 * @brief Pointer to user polynomail base
500 */
501 boost::shared_ptr<BaseFunction> userPolynomialBasePtr;
502
503 /**
504 * @brief Element to integrate on the sides
505 *
506 */
508
509 /**
510 * @brief Set the pointer to face element on the side
511 *
512 * \note Function is is used by face element, while it iterates over
513 * elements on the side
514 *
515 * @param side_fe_ptr
516 * @return MoFEMErrorCode
517 */
519
520 /**u
521 * @brief Element to integrate parent or child
522 *
523 */
525
526 /**
527 * @brief Set the pointer to face element refined
528 *
529 * @param refine_fe_ptr
530 * @return MoFEMErrorCode
531 */
533
539
540protected:
541 MatrixDouble coordsAtGaussPts; ///< coordinated at gauss points
542 double elementMeasure; ///< Depending on dimension of elements, stores length,
543 ///< area or volume of element.
544};
545
547
548 /**
549 * \brief Controls loop over entities on element
550 *
551 * - OPROW is used if row vector is assembled
552 * - OPCOL is usually used if column vector is assembled
553 * - OPROWCOL is usually used for assemble matrices.
554 * - OPSPACE no field is defined for such operator. Is usually used to modify
555 * base
556 *
557 *
558 * For typical problem like Bubnov-Galerkin OPROW and OPCOL are the same. In
559 * more general case for example for non-square matrices columns and rows
560 * could have different numeration and/or different set of fields.
561 *
562 */
563 enum OpType {
564 OPROW = 1 << 0, ///< operator doWork function is executed on FE rows
565 OPCOL = 1 << 1, ///< operator doWork function is executed on FE columns
566 OPROWCOL = 1 << 2, ///< operator doWork is executed on FE rows &columns
567 OPSPACE = 1 << 3, ///< operator do Work is execute on space data
568 OPLAST = 1 << 3 ///< @deprecated would be removed
569 };
570
571 static const char *const OpTypeNames[];
572
573 char opType;
574 std::string rowFieldName;
575 std::string colFieldName;
577
578 /**
579 * This Constructor is used typically when some modification base shape
580 * functions on some approx. space is applied. Operator is run for all data
581 * on space.
582 *
583 * User has no access to field data from this operator.
584 */
585 UserDataOperator(const FieldSpace space, const char type = OPSPACE,
586 const bool symm = true);
587
588 UserDataOperator(const std::string field_name, const char type,
589 const bool symm = true);
590
591 UserDataOperator(const std::string row_field_name,
592 const std::string col_field_name, const char type,
593 const bool symm = true);
594
595 /** \brief Return raw pointer to NumeredEntFiniteElement
596 */
597 inline boost::shared_ptr<const NumeredEntFiniteElement>
599
600 /**
601 * \brief Return finite element entity handle
602 * @return Finite element entity handle
603 */
604 inline EntityHandle getFEEntityHandle() const;
605
606 /**
607 * @brief Get dimension of finite element
608 *
609 * @return int
610 */
611 inline int getFEDim() const;
612
613 /**
614 * @brief Get dimension of finite element
615 *
616 * @return int
617 */
618 inline EntityType getFEType() const;
619
620 /**
621 * @brief Get the side number pointer
622 *
623 * \note For vertex is expectation. Side basses in argument of function doWork
624 * is zero. For other entity types side can be used as argument of this
625 * function.
626 *
627 * @param side_number
628 * @param type
629 * @return boost::weak_ptr<SideNumber>
630 */
631 inline boost::weak_ptr<SideNumber> getSideNumberPtr(const int side_number,
632 const EntityType type);
633
634 /**
635 * @brief Get the side entity
636 *
637 * \note For vertex is expectation. Side basses in argument of function
638 * doWork is zero. For other entity types side can be used as argument of
639 * this function.
640 *
641 * \code
642 * MoFEMErrorCode doWork(int side, EntityType type,
643 * EntitiesFieldData::EntData &data) {
644 * MoFEMFunctionBegin;
645 *
646 * if (type == MBVERTEX) {
647 * for (int n = 0; n != number_of_nodes; ++n)
648 * EntityHandle ent = getSideEntity(n, type);
649 *
650 * // Do somthing
651 *
652 * } else {
653 * EntityHandle ent = getSideEntity(side, type);
654 *
655 * // Do somthing
656 *
657 * }
658 *
659 * MoFEMFunctionReturn(0);
660 * }
661 * \endcode
662 *
663 * @param side_number
664 * @param type
665 */
666 inline EntityHandle getSideEntity(const int side_number,
667 const EntityType type);
668
669 /**
670 * @brief Get the number of nodes on finite element
671 *
672 * @return int
673 */
674 inline int getNumberOfNodesOnElement() const;
675
676 /** \brief Get row indices
677
678 Field could be or not declared for this element but is declared for
679 problem
680
681 \param field_name
682 \param type entity type
683 \param side side number, any number if type is MBVERTEX
684 \return indices
685
686 NOTE: Using those indices to assemble matrix will result in error if new
687 non-zero values need to be created.
688
689 */
690 MoFEMErrorCode getProblemRowIndices(const std::string filed_name,
691 const EntityType type, const int side,
692 VectorInt &indices) const;
693
694 /** \brief Get col indices
695
696 Field could be or not declared for this element but is declared for
697 problem
698
699 \param field_name
700 \param type entity type
701 \param side side number, any number if type is MBVERTEX
702 \return indices
703
704 NOTE: Using those indices to assemble matrix will result in error if new
705 non-zero values need to be created.
706
707 */
708 MoFEMErrorCode getProblemColIndices(const std::string filed_name,
709 const EntityType type, const int side,
710 VectorInt &indices) const;
711
712 /** \brief Return raw pointer to Finite Element Method object
713 */
714 inline const FEMethod *getFEMethod() const;
715
716 /**
717 * \brief Get operator types
718 * @return Return operator type
719 */
720 inline int getOpType() const;
721
722 /**
723 * \brief Set operator type
724 * @param Operator type
725 */
726 inline void setOpType(const OpType type);
727
728 /**
729 * \brief Add operator type
730 */
731 inline void addOpType(const OpType type);
732
733 /**
734 * \brief get number of finite element in the loop
735 * @return number of finite element
736 */
737 inline int getNinTheLoop() const;
738
739 /**
740 * \brief get size of elements in the loop
741 * @return loop size
742 */
743 inline int getLoopSize() const;
744
745 /** \brief Get name of the element
746 */
747 inline std::string getFEName() const;
748
749 /** \name Accessing KSP */
750
751 /**@{*/
752
753 inline const PetscData::Switches &getDataCtx() const;
754
755 inline const KspMethod::KSPContext getKSPCtx() const;
756
757 inline const SnesMethod::SNESContext getSNESCtx() const;
758
759 inline const TSMethod::TSContext getTSCtx() const;
760
761 /**@}*/
762
763 /**@{*/
764
765 inline Vec getKSPf() const;
766
767 inline Mat getKSPA() const;
768
769 inline Mat getKSPB() const;
770
771 /**@}*/
772
773 /** \name Accessing SNES */
774
775 /**@{*/
776
777 inline Vec getSNESf() const;
778
779 inline Vec getSNESx() const;
780
781 inline Mat getSNESA() const;
782
783 inline Mat getSNESB() const;
784
785 /**@}*/
786
787 /** \name Accessing TS */
788
789 /**@{*/
790
791 inline Vec getTSu() const;
792
793 inline Vec getTSu_t() const;
794
795 inline Vec getTSu_tt() const;
796
797 inline Vec getTSf() const;
798
799 inline Mat getTSA() const;
800
801 inline Mat getTSB() const;
802
803 inline int getTSstep() const;
804
805 inline double getTStime() const;
806
807 inline double getTStimeStep() const;
808
809 inline double getTSa() const;
810
811 inline double getTSaa() const;
812
813 /**@}*/
814
815 /**@{*/
816
817 /** \name Base funtions and integration points */
818
819 /** \brief matrix of integration (Gauss) points for Volume Element
820 *
821 * For triangle: columns 0,1 are x,y coordinates respectively and column
822 * 2 is a weight value for example getGaussPts()(1,13) returns y
823 * coordinate of 13th Gauss point on particular volume element
824 *
825 * For tetrahedron: columns 0,1,2 are x,y,z coordinates respectively and
826 * column 3 is a weight value for example getGaussPts()(1,13) returns y
827 * coordinate of 13th Gauss point on particular volume element
828 *
829 */
830 inline MatrixDouble &getGaussPts();
831
832 /**
833 * @brief Get integration weights
834 *
835 * \code
836 * auto t_w = getFTensor0IntegrationWeight();
837 * for(int gg = 0; gg!=getGaussPts.size2(); ++gg) {
838 * // integrate something
839 * ++t_w;
840 * }
841 * \endcode
842 *
843 * @return FTensor::Tensor0<FTensor::PackPtr<double *, 1>>
844 */
845 inline auto getFTensor0IntegrationWeight();
846
847 /**@}*/
848
849 /** \name Coordinates and access to internal data */
850
851 /**@{*/
852
853 /** \brief Gauss points and weight, matrix (nb. of points x 3)
854
855 Column 0-2 integration points coordinate x, y and z, respectively. At rows
856 are integration points.
857
858 */
860
861 /**
862 * \brief Get coordinates at integration points assuming linear geometry
863 *
864 * \code
865 * auto t_coords = getFTensor1CoordsAtGaussPts();
866 * for(int gg = 0;gg!=nb_int_ptrs;gg++) {
867 * // do something
868 * ++t_coords;
869 * }
870 * \endcode
871 *
872 */
873 inline auto getFTensor1CoordsAtGaussPts();
874
875 /**@}*/
876
877 /**@{*/
878
879 /** \name Measures (area, volume, length, etc.) */
880
881 /**
882 * \brief get measure of element
883 * @return volume
884 */
885 inline double getMeasure() const;
886
887 /**
888 * \brief get measure of element
889 * @return volume
890 */
891 inline double &getMeasure();
892
893 /**}*/
894
895 /**@{*/
896
897 /** \name Loops */
898
899 using AdjCache =
900 std::map<EntityHandle,
901 std::vector<boost::weak_ptr<NumeredEntFiniteElement>>>;
902
903 /**
904 * @brief User call this function to loop over elements on the side of
905 * face. This function calls finite element with is operator to do
906 * calculations.
907 *
908 * @param fe_name name of the side element
909 * @param side_fe pointer to the side element instance
910 * @param dim dimension the of side element
911 * @param ent_for_side entity handle for which adjacent volume or face will
912 * be accessed
913 * @param verb
914 * @param sev
915 * @return MoFEMErrorCode
916 */
917 MoFEMErrorCode loopSide(const string &fe_name, ForcesAndSourcesCore *side_fe,
918 const size_t dim, const EntityHandle ent_for_side = 0,
919 const int verb = QUIET,
920 const LogManager::SeverityLevel sev = Sev::noisy,
921 AdjCache *adj_cache = nullptr
922
923 );
924
925 /**
926 * @brief User call this function to loop over parent elements. This function
927 * calls finite element with is operator to do calculations.
928 *
929 * @param fe_name
930 * @param parent_fe
931 * @param verb
932 * @param sev
933 * @return MoFEMErrorCode
934 */
935 MoFEMErrorCode loopThis(const string &fe_name,
936 ForcesAndSourcesCore *parent_fe,
937 const int verb = QUIET,
938 const LogManager::SeverityLevel sev = Sev::noisy);
939
940 /**
941 * @brief User call this function to loop over parent elements. This function
942 * calls finite element with is operator to do calculations.
943 *
944 * @param fe_name
945 * @param parent_fe
946 * @param verb
947 * @param sev
948 * @return MoFEMErrorCode
949 */
950 MoFEMErrorCode loopParent(const string &fe_name,
951 ForcesAndSourcesCore *parent_fe,
952 const int verb = QUIET,
953 const LogManager::SeverityLevel sev = Sev::noisy);
954
955 /**
956 * @brief User call this function to loop over parent elements. This function
957 * calls finite element with is operator to do calculations.
958 *
959 * @param fe_name
960 * @param child_fe
961 * @param verb
962 * @param sev
963 * @return MoFEMErrorCode
964 */
965 MoFEMErrorCode loopChildren(const string &fe_name,
966 ForcesAndSourcesCore *child_fe,
967 const int verb = QUIET,
968 const LogManager::SeverityLevel sev = Sev::noisy);
969
970 /**@}*/
971
972 inline ForcesAndSourcesCore *getPtrFE() const;
973
974 inline ForcesAndSourcesCore *getSidePtrFE() const;
975
976 inline ForcesAndSourcesCore *getRefinePtrFE() const;
977
978
979protected:
981
983
984private:
989};
990
991/// \deprecated Used ForcesAndSourcesCore instead
993
994boost::shared_ptr<const NumeredEntFiniteElement>
995ForcesAndSourcesCore::UserDataOperator::getNumeredEntFiniteElementPtr() const {
997};
998
999EntityHandle ForcesAndSourcesCore::UserDataOperator::getFEEntityHandle() const {
1000 return getNumeredEntFiniteElementPtr()->getEnt();
1001}
1002
1003int ForcesAndSourcesCore::UserDataOperator::getFEDim() const {
1005};
1006
1007EntityType ForcesAndSourcesCore::UserDataOperator::getFEType() const {
1009};
1010
1011boost::weak_ptr<SideNumber>
1012ForcesAndSourcesCore::UserDataOperator::getSideNumberPtr(
1013 const int side_number, const EntityType type) {
1014 auto &side_table_by_side_and_type =
1015 ptrFE->numeredEntFiniteElementPtr->getSideNumberTable().get<1>();
1016 auto side_it =
1017 side_table_by_side_and_type.find(boost::make_tuple(type, side_number));
1018 if (side_it != side_table_by_side_and_type.end())
1019 return *side_it;
1020 else
1021 return boost::weak_ptr<SideNumber>();
1022}
1023
1025ForcesAndSourcesCore::UserDataOperator::getSideEntity(const int side_number,
1026 const EntityType type) {
1027 if (auto side_ptr = getSideNumberPtr(side_number, type).lock())
1028 return side_ptr->ent;
1029 else
1030 return 0;
1031}
1032
1033int ForcesAndSourcesCore::UserDataOperator::getNumberOfNodesOnElement() const {
1034 return ptrFE->getNumberOfNodes();
1035}
1036
1037const FEMethod *ForcesAndSourcesCore::UserDataOperator::getFEMethod() const {
1038 return ptrFE;
1039}
1040
1041int ForcesAndSourcesCore::UserDataOperator::getOpType() const { return opType; }
1042
1043void ForcesAndSourcesCore::UserDataOperator::setOpType(const OpType type) {
1044 opType = type;
1045}
1046
1047void ForcesAndSourcesCore::UserDataOperator::addOpType(const OpType type) {
1048 opType |= type;
1049}
1050
1051int ForcesAndSourcesCore::UserDataOperator::getNinTheLoop() const {
1052 return getFEMethod()->getNinTheLoop();
1053}
1054
1055int ForcesAndSourcesCore::UserDataOperator::getLoopSize() const {
1056 return getFEMethod()->getLoopSize();
1057}
1058
1059std::string ForcesAndSourcesCore::UserDataOperator::getFEName() const {
1060 return getFEMethod()->getFEName();
1061}
1062
1063const PetscData::Switches &
1064ForcesAndSourcesCore::UserDataOperator::getDataCtx() const {
1065 return getFEMethod()->data_ctx;
1066}
1067
1069ForcesAndSourcesCore::UserDataOperator::getKSPCtx() const {
1070 return getFEMethod()->ksp_ctx;
1071}
1072
1074ForcesAndSourcesCore::UserDataOperator::getSNESCtx() const {
1075 return getFEMethod()->snes_ctx;
1076}
1077
1079ForcesAndSourcesCore::UserDataOperator::getTSCtx() const {
1080 return getFEMethod()->ts_ctx;
1081}
1082
1083Vec ForcesAndSourcesCore::UserDataOperator::getKSPf() const {
1084#ifndef NDEBUG
1085 if (getFEMethod()->ksp_f == PETSC_NULL)
1086 THROW_MESSAGE("KSP not set F vector");
1087#endif
1088 return getFEMethod()->ksp_f;
1089}
1090
1091Mat ForcesAndSourcesCore::UserDataOperator::getKSPA() const {
1092#ifndef NDEBUG
1093 if (getFEMethod()->ksp_A == PETSC_NULL)
1094 THROW_MESSAGE("KSP not set A vector");
1095#endif
1096 return getFEMethod()->ksp_A;
1097}
1098
1099Mat ForcesAndSourcesCore::UserDataOperator::getKSPB() const {
1100#ifndef NDEBUG
1101 if (getFEMethod()->ksp_B == PETSC_NULL)
1102 THROW_MESSAGE("KSP not set B vector");
1103#endif
1104 return getFEMethod()->ksp_B;
1105}
1106
1107Vec ForcesAndSourcesCore::UserDataOperator::getSNESf() const {
1108#ifndef NDEBUG
1109 if (getFEMethod()->snes_f == PETSC_NULL)
1110 THROW_MESSAGE("SNES not set F vector");
1111#endif
1112 return getFEMethod()->snes_f;
1113}
1114
1115Vec ForcesAndSourcesCore::UserDataOperator::getSNESx() const {
1116#ifndef NDEBUG
1117 if (getFEMethod()->snes_x == PETSC_NULL)
1118 THROW_MESSAGE("SNESnot set X vector");
1119#endif
1120 return getFEMethod()->snes_x;
1121}
1122
1123Mat ForcesAndSourcesCore::UserDataOperator::getSNESA() const {
1124#ifndef NDEBUG
1125 if (getFEMethod()->snes_A == PETSC_NULL)
1126 THROW_MESSAGE("SNES not set A vector");
1127#endif
1128 return getFEMethod()->snes_A;
1129}
1130
1131Mat ForcesAndSourcesCore::UserDataOperator::getSNESB() const {
1132#ifndef NDEBUG
1133 if (getFEMethod()->snes_B == PETSC_NULL)
1134 THROW_MESSAGE("SNES not set A matrix");
1135#endif
1136 return getFEMethod()->snes_B;
1137}
1138
1139Vec ForcesAndSourcesCore::UserDataOperator::getTSu() const {
1140#ifndef NDEBUG
1141 if (getFEMethod()->ts_u == PETSC_NULL)
1142 THROW_MESSAGE("TS not set U vector");
1143#endif
1144 return getFEMethod()->ts_u;
1145}
1146
1147Vec ForcesAndSourcesCore::UserDataOperator::getTSu_t() const {
1148#ifndef NDEBUG
1149 if (getFEMethod()->ts_u_t == PETSC_NULL)
1150 THROW_MESSAGE("TS not set U_t vector");
1151#endif
1152 return getFEMethod()->ts_u_t;
1153}
1154
1155Vec ForcesAndSourcesCore::UserDataOperator::getTSu_tt() const {
1156#ifndef NDEBUG
1157 if (getFEMethod()->ts_u_tt == PETSC_NULL)
1158 THROW_MESSAGE("TS not set U_tt vector");
1159#endif
1160 return getFEMethod()->ts_u_tt;
1161}
1162
1163Vec ForcesAndSourcesCore::UserDataOperator::getTSf() const {
1164#ifndef NDEBUG
1165 if (getFEMethod()->ts_F == PETSC_NULL)
1166 THROW_MESSAGE("TS not set F vector");
1167#endif
1168 return getFEMethod()->ts_F;
1169}
1170
1171Mat ForcesAndSourcesCore::UserDataOperator::getTSA() const {
1172#ifndef NDEBUG
1173 if (getFEMethod()->ts_A == PETSC_NULL)
1174 THROW_MESSAGE("TS not set A matrix");
1175#endif
1176 return getFEMethod()->ts_A;
1177}
1178
1179Mat ForcesAndSourcesCore::UserDataOperator::getTSB() const {
1180#ifndef NDEBUG
1181 if (getFEMethod()->ts_B == PETSC_NULL)
1182 THROW_MESSAGE("TS not set B matrix");
1183#endif
1184 return getFEMethod()->ts_B;
1185}
1186
1187int ForcesAndSourcesCore::UserDataOperator::getTSstep() const {
1188#ifndef NDEBUG
1189 if ((getFEMethod()->data_ctx & PetscData::PetscData::CtxSetTime).none())
1190 THROW_MESSAGE("TS not set time");
1191#endif
1192 return getFEMethod()->ts_step;
1193}
1194
1195double ForcesAndSourcesCore::UserDataOperator::getTStime() const {
1196#ifndef NDEBUG
1197 if ((getFEMethod()->data_ctx & PetscData::PetscData::CtxSetTime).none())
1198 THROW_MESSAGE("TS not set time");
1199#endif
1200 return getFEMethod()->ts_t;
1201}
1202
1203double ForcesAndSourcesCore::UserDataOperator::getTStimeStep() const {
1204#ifndef NDEBUG
1205 if ((getFEMethod()->data_ctx & PetscData::PetscData::CtxSetTime).none())
1206 THROW_MESSAGE("TS not set time");
1207#endif
1208 return getFEMethod()->ts_dt;
1209}
1210
1211double ForcesAndSourcesCore::UserDataOperator::getTSa() const {
1212#ifndef NDEBUG
1213 if ((getFEMethod()->data_ctx & (PetscData::CtxSetA | PetscData::CtxSetB))
1214 .none() ||
1215 (getFEMethod()->data_ctx & (PetscData::CtxSetX_T)).none())
1216 THROW_MESSAGE("TS not set B matrix");
1217#endif
1218 return getFEMethod()->ts_a;
1219}
1220
1221double ForcesAndSourcesCore::UserDataOperator::getTSaa() const {
1222#ifndef NDEBUG
1223 if ((getFEMethod()->data_ctx & (PetscData::CtxSetA | PetscData::CtxSetB))
1224 .none() ||
1225 (getFEMethod()->data_ctx & (PetscData::CtxSetX_TT)).none())
1226 THROW_MESSAGE("TS not set B matrix");
1227#endif
1228 return getFEMethod()->ts_aa;
1229}
1230
1231MatrixDouble &ForcesAndSourcesCore::UserDataOperator::getGaussPts() {
1232 return static_cast<ForcesAndSourcesCore *>(ptrFE)->gaussPts;
1233}
1234
1235auto ForcesAndSourcesCore::UserDataOperator::getFTensor0IntegrationWeight() {
1237 &(getGaussPts()(getGaussPts().size1() - 1, 0)));
1238}
1239
1240// MoFEMErrorCode ForcesAndSourcesCore::UserDataOperator::getPorblemRowIndices(
1241// const std::string filed_name, const EntityType type, const int side,
1242// VectorInt &indices) const {
1243// return getProblemRowIndices(filed_name, type, side, indices);
1244// }
1245
1246ForcesAndSourcesCore *ForcesAndSourcesCore::UserDataOperator::getPtrFE() const {
1247 return ptrFE;
1248}
1249
1251ForcesAndSourcesCore::UserDataOperator::getSidePtrFE() const {
1252 return ptrFE->sidePtrFE;
1253}
1254
1256ForcesAndSourcesCore::UserDataOperator::getRefinePtrFE() const {
1257 return ptrFE->refinePtrFE;
1258}
1259
1260MatrixDouble &ForcesAndSourcesCore::UserDataOperator::getCoordsAtGaussPts() {
1261 return static_cast<ForcesAndSourcesCore *>(ptrFE)->coordsAtGaussPts;
1262}
1263
1264auto ForcesAndSourcesCore::UserDataOperator::getFTensor1CoordsAtGaussPts() {
1266 &getCoordsAtGaussPts()(0, 0), &getCoordsAtGaussPts()(0, 1),
1267 &getCoordsAtGaussPts()(0, 2));
1268}
1269
1270double ForcesAndSourcesCore::UserDataOperator::getMeasure() const {
1271 return static_cast<ForcesAndSourcesCore *>(ptrFE)->elementMeasure;
1272}
1273
1274double &ForcesAndSourcesCore::UserDataOperator::getMeasure() {
1275 return static_cast<ForcesAndSourcesCore *>(ptrFE)->elementMeasure;
1276}
1277
1278/**
1279 * @brief Element used to execute operators on side of the element
1280 *
1281 * @tparam E template for side element type
1282 *
1283 */
1284template <typename E>
1286
1288
1289 /**
1290 * @brief Construct a new Op Loop Side object
1291 *
1292 * @param m_field
1293 * @param fe_name name of side (domain element)
1294 * @param side_dim dimension
1295 */
1296 OpLoopSide(MoFEM::Interface &m_field, const std::string fe_name,
1297 const int side_dim,
1298 const LogManager::SeverityLevel sev = Sev::noisy,
1299 boost::shared_ptr<AdjCache> adj_cache = nullptr)
1300 : UserDataOperator(NOSPACE, OPSPACE), sideFEPtr(new E(m_field)),
1301 fieldName(fe_name), sideDim(side_dim), sevLevel(sev),
1302 adjCache(adj_cache) {}
1303
1307 CHKERR loopSide(fieldName, sideFEPtr.get(), sideDim, 0, VERBOSE, sevLevel,
1308 adjCache.get());
1310 };
1311
1312 boost::ptr_deque<UserDataOperator> &getOpPtrVector() {
1313 return sideFEPtr->getOpPtrVector();
1314 }
1315
1316 boost::shared_ptr<E> &getSideFEPtr() { return sideFEPtr; }
1317
1318protected:
1319 const std::string fieldName;
1320 const int sideDim;
1321 boost::shared_ptr<E> sideFEPtr;
1323 boost::shared_ptr<AdjCache> adjCache;
1324};
1325
1326} // namespace MoFEM
1327
1328#endif //__FORCES_AND_SOURCES_CORE__HPP__
1329
1330/**
1331 * \defgroup mofem_forces_and_sources Forces and sources
1332 * \ingroup mofem
1333 *
1334 * \brief Manages complexities related to assembly of vector and matrices at
1335 * single finite element level.
1336 *
1337 **/
@ QUIET
Definition: definitions.h:208
@ VERBOSE
Definition: definitions.h:209
FieldApproximationBase
approximation base
Definition: definitions.h:58
FieldSpace
approximation spaces
Definition: definitions.h:82
@ LASTSPACE
FieldSpace in [ 0, LASTSPACE )
Definition: definitions.h:89
@ NOSPACE
Definition: definitions.h:83
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define DEPRECATED
Definition: definitions.h:17
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
Definition: definitions.h:561
const int dim
multi_index_container< boost::shared_ptr< NumeredDofEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, const_mem_fun< NumeredDofEntity::interface_type_DofEntity, UId, &NumeredDofEntity::getLocalUniqueId > >, ordered_non_unique< tag< Part_mi_tag >, member< NumeredDofEntity, unsigned int, &NumeredDofEntity::pArt > >, ordered_non_unique< tag< Idx_mi_tag >, member< NumeredDofEntity, DofIdx, &NumeredDofEntity::dofIdx > >, ordered_non_unique< tag< PetscGlobalIdx_mi_tag >, member< NumeredDofEntity, DofIdx, &NumeredDofEntity::petscGloablDofIdx > >, ordered_non_unique< tag< PetscLocalIdx_mi_tag >, member< NumeredDofEntity, DofIdx, &NumeredDofEntity::petscLocalDofIdx > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< NumeredDofEntity::interface_type_DofEntity, EntityHandle, &NumeredDofEntity::getEnt > > > > NumeredDofEntity_multiIndex
MultiIndex container keeps NumeredDofEntity.
SeverityLevel
Severity levels.
Definition: LogManager.hpp:33
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
UBlasVector< int > VectorInt
Definition: Types.hpp:67
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
auto type_from_handle(const EntityHandle h)
get type from entity handle
Definition: Templates.hpp:1761
ublas::vector< FieldEntity *, FieldEntAllocator > VectorFieldEntities
DEPRECATED typedef ForcesAndSourcesCore ForcesAndSurcesCore
auto dimension_from_handle(const EntityHandle h)
get entity dimension form handle
Definition: Templates.hpp:1777
std::vector< boost::weak_ptr< FieldEntity > > FieldEntity_vector_view
ublas::vector< FEDofEntity *, DofsAllocator > VectorDofs
constexpr auto field_name
base operator to do operations at Gauss Pt. level
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
data structure for finite element entity
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
structure for User Loop Methods on finite elements
EntityHandle getFEEntityHandle() const
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
Base face element used to integrate on skeleton.
void setOpType(const OpType type)
Set operator type.
EntityHandle getFEEntityHandle() const
Return finite element entity handle.
MoFEMErrorCode getProblemColIndices(const std::string filed_name, const EntityType type, const int side, VectorInt &indices) const
Get col indices.
int getLoopSize() const
get size of elements in the loop
EntityType getFEType() const
Get dimension of finite element.
const KspMethod::KSPContext getKSPCtx() const
int getNinTheLoop() const
get number of finite element in the loop
auto getFTensor1CoordsAtGaussPts()
Get coordinates at integration points assuming linear geometry.
MoFEMErrorCode loopThis(const string &fe_name, ForcesAndSourcesCore *parent_fe, const int verb=QUIET, const LogManager::SeverityLevel sev=Sev::noisy)
User call this function to loop over parent elements. This function calls finite element with is oper...
int getNumberOfNodesOnElement() const
Get the number of nodes on finite element.
boost::weak_ptr< SideNumber > getSideNumberPtr(const int side_number, const EntityType type)
Get the side number pointer.
const SnesMethod::SNESContext getSNESCtx() const
MoFEMErrorCode loopParent(const string &fe_name, ForcesAndSourcesCore *parent_fe, const int verb=QUIET, const LogManager::SeverityLevel sev=Sev::noisy)
User call this function to loop over parent elements. This function calls finite element with is oper...
const PetscData::Switches & getDataCtx() const
EntityHandle getSideEntity(const int side_number, const EntityType type)
Get the side entity.
MatrixDouble & getCoordsAtGaussPts()
Gauss points and weight, matrix (nb. of points x 3)
auto getFTensor0IntegrationWeight()
Get integration weights.
double getMeasure() const
get measure of element
std::string getFEName() const
Get name of the element.
MoFEMErrorCode getProblemRowIndices(const std::string filed_name, const EntityType type, const int side, VectorInt &indices) const
Get row indices.
boost::shared_ptr< const NumeredEntFiniteElement > getNumeredEntFiniteElementPtr() const
Return raw pointer to NumeredEntFiniteElement.
std::map< EntityHandle, std::vector< boost::weak_ptr< NumeredEntFiniteElement > > > AdjCache
void addOpType(const OpType type)
Add operator type.
OpType
Controls loop over entities on element.
@ OPCOL
operator doWork function is executed on FE columns
@ OPROW
operator doWork function is executed on FE rows
@ OPROWCOL
operator doWork is executed on FE rows &columns
@ OPSPACE
operator do Work is execute on space data
MoFEMErrorCode loopChildren(const string &fe_name, ForcesAndSourcesCore *child_fe, const int verb=QUIET, const LogManager::SeverityLevel sev=Sev::noisy)
User call this function to loop over parent elements. This function calls finite element with is oper...
virtual MoFEMErrorCode setPtrFE(ForcesAndSourcesCore *ptr)
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
MoFEMErrorCode loopSide(const string &fe_name, ForcesAndSourcesCore *side_fe, const size_t dim, const EntityHandle ent_for_side=0, const int verb=QUIET, const LogManager::SeverityLevel sev=Sev::noisy, AdjCache *adj_cache=nullptr)
User call this function to loop over elements on the side of face. This function calls finite element...
int getFEDim() const
Get dimension of finite element.
MatrixDouble & getGaussPts()
matrix of integration (Gauss) points for Volume Element
structure to get information form mofem into EntitiesFieldData
MoFEMErrorCode loopOverOperators()
Iterate user data operators.
MoFEMErrorCode getProblemNodesRowIndices(const std::string &field_name, VectorInt &nodes_indices) const
auto & getDataOnElementBySpaceArray()
Get data on entities and space.
MoFEMErrorCode getEntityIndices(EntitiesFieldData &data, const int bit_number, FieldEntity_vector_view &ents_field, const EntityType type_lo, const EntityType type_hi, EXTRACTOR &&extractor) const
MoFEMErrorCode getEntityRowIndices(EntitiesFieldData &data, const int bit_number, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
MoFEMErrorCode getRowNodesIndices(EntitiesFieldData &data, const int bit_number) const
get row node indices from FENumeredDofEntity_multiIndex
int getMaxRowOrder() const
Get max order of approximation for field in rows.
auto & getUserPolynomialBase()
Get the User Polynomial Base object.
MoFEMErrorCode getEntitySense(EntitiesFieldData &data) const
Get the entity sense (orientation)
MoFEMErrorCode getSpacesAndBaseOnEntities(EntitiesFieldData &data) const
Get field approximation space and base on entities.
MoFEMErrorCode getProblemTypeColIndices(const std::string &field_name, EntityType type, int side_number, VectorInt &indices) const
MoFEMErrorCode getNoFieldIndices(const int bit_number, boost::shared_ptr< FENumeredDofEntity_multiIndex > dofs, VectorInt &nodes_indices) const
get NoField indices
virtual MoFEMErrorCode operator()()
function is run for every finite element
MoFEMErrorCode getFaceNodes(EntitiesFieldData &data) const
Get nodes on faces.
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACE > derivedDataOnElement
Entity data on element entity columns fields.
MoFEMErrorCode getNoFieldRowIndices(EntitiesFieldData &data, const int bit_number) const
get col NoField indices
virtual MoFEMErrorCode setGaussPts(int order_row, int order_col, int order_data)
set user specific integration rule
MoFEMErrorCode getNodesIndices(const int bit_number, FieldEntity_vector_view &ents_field, VectorInt &nodes_indices, VectorInt &local_nodes_indices, EXTRACTOR &&extractor) const
get node indices
MoFEMErrorCode getEntityFieldData(EntitiesFieldData &data, const int bit_number, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
MoFEMErrorCode setRefineFEPtr(const ForcesAndSourcesCore *refine_fe_ptr)
Set the pointer to face element refined.
auto & getElementPolynomialBase()
Get the Entity Polynomial Base object.
ForcesAndSourcesCore * sidePtrFE
Element to integrate on the sides.
MatrixDouble coordsAtGaussPts
coordinated at gauss points
boost::function< MoFEMErrorCode(ForcesAndSourcesCore *fe_raw_ptr, int order_row, int order_col, int order_data)> GaussHookFun
virtual int getRule(int order_row, int order_col, int order_data)
another variant of getRule
boost::ptr_deque< UserDataOperator > opPtrVector
Vector of finite element users data operators.
MoFEMErrorCode getEntitySense(const EntityType type, boost::ptr_vector< EntitiesFieldData::EntData > &data) const
get sense (orientation) of entity
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACE > dataOnElement
Entity data on element entity rows fields.
auto & getDerivedDataOnElementBySpaceArray()
Get derived data on entities and space.
MoFEMErrorCode getProblemNodesColIndices(const std::string &field_name, VectorInt &nodes_indices) const
MoFEMErrorCode calHierarchicalBaseFunctionsOnElement()
Calculate base functions.
MatrixDouble gaussPts
Matrix of integration points.
MoFEMErrorCode calBernsteinBezierBaseFunctionsOnElement()
Calculate Bernstein-Bezier base.
auto & getEntData(const FieldSpace space, const EntityType type, const int side)
Get the entity data.
MoFEMErrorCode getEntityDataOrder(EntitiesFieldData &data, const FieldSpace space) const
Get the entity data order for given space.
virtual MoFEMErrorCode preProcess()
function is run at the beginning of loop
MoFEMErrorCode getColNodesIndices(EntitiesFieldData &data, const int bit_number) const
get col node indices from FENumeredDofEntity_multiIndex
MoFEMErrorCode setSideFEPtr(const ForcesAndSourcesCore *side_fe_ptr)
Set the pointer to face element on the side.
int getMaxColOrder() const
Get max order of approximation for field in columns.
EntityType lastEvaluatedElementEntityType
Last evaluated type of element entity.
boost::function< int(int order_row, int order_col, int order_data)> RuleHookFun
MoFEMErrorCode getProblemNodesIndices(const std::string &field_name, const NumeredDofEntity_multiIndex &dofs, VectorInt &nodes_indices) const
get indices of nodal indices which are declared for problem but not this particular element
boost::ptr_deque< UserDataOperator > & getOpPtrVector()
Use to push back operator for row operator.
GaussHookFun setRuleHook
Set function to calculate integration rule.
MoFEMErrorCode getNodesFieldData(EntitiesFieldData &data, const int bit_number) const
Get data on nodes.
MoFEMErrorCode getEntityColIndices(EntitiesFieldData &data, const int bit_number, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
boost::shared_ptr< BaseFunction > elementPolynomialBasePtr
Pointer to entity polynomial base.
MoFEMErrorCode getEntityDataOrder(const EntityType type, const FieldSpace space, boost::ptr_vector< EntitiesFieldData::EntData > &data) const
Get the entity data order.
ForcesAndSourcesCore * refinePtrFE
Element to integrate parent or child.
MoFEMErrorCode getNoFieldColIndices(EntitiesFieldData &data, const int bit_number) const
get col NoField indices
int getMaxDataOrder() const
Get max order of approximation for data fields.
MoFEMErrorCode getNoFieldFieldData(const int bit_number, VectorDouble &ent_field_data, VectorDofs &ent_field_dofs, VectorFieldEntities &ent_field) const
Get field data on nodes.
boost::shared_ptr< BaseFunction > userPolynomialBasePtr
Pointer to user polynomail base.
MoFEMErrorCode getProblemTypeRowIndices(const std::string &field_name, EntityType type, int side_number, VectorInt &indices) const
MoFEMErrorCode getProblemTypeIndices(const std::string &field_name, const NumeredDofEntity_multiIndex &dofs, EntityType type, int side_number, VectorInt &indices) const
get indices by type (generic function) which are declared for problem but not this particular element
RuleHookFun getRuleHook
Hook to get rule.
MoFEMErrorCode createDataOnElement(EntityType type)
Create a entity data on element object.
virtual MoFEMErrorCode postProcess()
function is run at the end of loop
KSPContext
pass information about context of KSP/DM for with finite element is computed
Definition: LoopMethods.hpp:73
Element used to execute operators on side of the element.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
boost::shared_ptr< E > & getSideFEPtr()
boost::shared_ptr< AdjCache > adjCache
boost::shared_ptr< E > sideFEPtr
boost::ptr_deque< UserDataOperator > & getOpPtrVector()
const std::string fieldName
const LogManager::SeverityLevel sevLevel
OpLoopSide(MoFEM::Interface &m_field, const std::string fe_name, const int side_dim, const LogManager::SeverityLevel sev=Sev::noisy, boost::shared_ptr< AdjCache > adj_cache=nullptr)
Construct a new Op Loop Side object.
static constexpr Switches CtxSetA
Definition: LoopMethods.hpp:37
static constexpr Switches CtxSetX_TT
Definition: LoopMethods.hpp:41
std::bitset< 8 > Switches
Definition: LoopMethods.hpp:33
static constexpr Switches CtxSetX_T
Definition: LoopMethods.hpp:40
static constexpr Switches CtxSetB
Definition: LoopMethods.hpp:38
Vec & snes_f
residual
Vec & snes_x
state vector
Mat & snes_B
preconditioner of jacobian matrix
Mat & snes_A
jacobian matrix
Vec & ts_F
residual vector
Vec & ts_u_tt
second time derivative of state vector
Vec & ts_u_t
time derivative of state vector
Vec & ts_u
state vector
Mat & ts_B
Preconditioner for ts_A.
Base volume element used to integrate on contact surface (could be extended to other volume elements)...
Base volume element used to integrate on skeleton.