v0.10.0
SurfacePressure.hpp
Go to the documentation of this file.
1 /* \file SurfacePressure.hpp
2  \brief Implementation of pressure and forces on triangles surface
3 
4 */
5 
6 /* This file is part of MoFEM.
7  * MoFEM is free software: you can redistribute it and/or modify it under
8  * the terms of the GNU Lesser General Public License as published by the
9  * Free Software Foundation, either version 3 of the License, or (at your
10  * option) any later version.
11  *
12  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
13  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
14  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
15  * License for more details.
16  *
17  * You should have received a copy of the GNU Lesser General Public
18  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
19 
20 #ifndef __SURFACE_PERSSURE_HPP__
21 #define __SURFACE_PERSSURE_HPP__
22 
23 /** \brief Finite element and operators to apply force/pressures applied to
24  * surfaces \ingroup mofem_static_boundary_conditions
25  */
27 
29 
30  /**
31  * \brief Analytical force method
32  */
34 
36 
37  /**
38  * User implemented analytical force
39  * @param coords coordinates of integration point
40  * @param normal normal at integration point
41  * @param force returned force
42  * @return error code
43  */
45  const VectorDouble3 &coords,
46  const VectorDouble3 &normal,
47  VectorDouble3 &force) {
49  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
50  "You need to implement this");
52  }
53  };
54 
56 
57  LinearVaringPresssure(const VectorDouble3 &p, const double c)
59 
60  MoFEMErrorCode getForce(const EntityHandle ent, const VectorDouble3 &coords,
61  const VectorDouble3 &normal, VectorDouble3 &force);
62 
63  private:
65  const double pressureShift;
66  };
67 
68  /**
69  * Definition of face element used for integration
70  */
72  int addToRule;
74  int getRule(int order) { return 2 * order + addToRule; };
75  };
76 
77  // FE for the right-hand side (spatial configuration)
79  MyTriangleFE &getLoopFe() { return fe; }
80 
81  // FE for the left-hand side (spatial configuration, ALE)
84 
85  // FE for the right-hand side (material configuration, ALE)
88 
89  // FE for the left-hand side (material configuration, ALE)
92 
94  : mField(m_field), fe(m_field), feLhs(m_field), feMatRhs(m_field),
95  feMatLhs(m_field) {}
96 
97  struct bCForce {
98  ForceCubitBcData data;
99  Range tRis;
100  };
101  std::map<int, bCForce> mapForce;
102 
103  struct bCPressure {
104  PressureCubitBcData data;
105  Range tRis;
106  };
107  std::map<int, bCPressure> mapPressure;
108 
109  boost::ptr_vector<MethodForForceScaling> methodsOp;
110  boost::ptr_vector<MethodForAnalyticalForce> analyticalForceOp;
111 
112  using UserDataOperator =
117 
118  /// Operator for force element
120 
121  Vec F;
123  boost::ptr_vector<MethodForForceScaling> &methodsOp;
124 
126 
127  OpNeumannForce(const std::string field_name, Vec _F, bCForce &data,
128  boost::ptr_vector<MethodForForceScaling> &methods_op,
129  bool ho_geometry = false);
130 
131  VectorDouble Nf; //< Local force vector
132 
133  MoFEMErrorCode doWork(int side, EntityType type,
135  };
136 
137  /// Operator for force element
139 
141  const std::string field_name, Vec f, const Range tris,
142  boost::ptr_vector<MethodForForceScaling> &methods_op,
143  boost::shared_ptr<MethodForAnalyticalForce> &analytical_force_op,
144  const bool ho_geometry = false);
145 
146  VectorDouble nF; //< Local force vector
147 
148  MoFEMErrorCode doWork(int side, EntityType type,
150 
151  Vec F;
152 
153  private:
154  const Range tRis;
155  boost::ptr_vector<MethodForForceScaling> &methodsOp;
156  boost::shared_ptr<MethodForAnalyticalForce> analyticalForceOp;
157  const bool hoGeometry;
158  };
159 
160  /**
161  * @brief RHS-operator for pressure element (spatial configuration)
162  *
163  */
165 
166  Vec F;
168  boost::ptr_vector<MethodForForceScaling> &methodsOp;
170 
171  OpNeumannPressure(const std::string field_name, Vec _F, bCPressure &data,
172  boost::ptr_vector<MethodForForceScaling> &methods_op,
173  bool ho_geometry = false);
174 
176 
177  /**
178  * @brief Integrate pressure
179  *
180  * \f[
181  * \begin{split}
182  * \mathbf{t} &= p \mathbf{n} \\
183  * \mathbf{f}^i &= \int_\mathcal{T} {\pmb\phi}^i \mathbf{t}
184  * \textrm{d}\mathcal{T} \end{split}
185  * \f]
186  * where \f$p\f$ is pressure, \f$\mathbf{n}\f$ is normal,
187  * \f$\mathbf{t}\f$ is traction, and
188  * \f$\mathbf{f}^i\f$ is local vector of external forces for ith base
189  * function \f${\pmb\phi}^i\f$.
190  *
191  *
192  * @param side
193  * @param type
194  * @param data
195  * @return MoFEMErrorCode
196  */
197  MoFEMErrorCode doWork(int side, EntityType type,
199  };
200 
202  : public boost::enable_shared_from_this<DataAtIntegrationPts> {
203 
206 
212 
213  inline boost::shared_ptr<MatrixDouble> getSmallhMatPtr() {
214  return boost::shared_ptr<MatrixDouble>(shared_from_this(), &hMat);
215  }
216 
217  inline boost::shared_ptr<MatrixDouble> getHMatPtr() {
218  return boost::shared_ptr<MatrixDouble>(shared_from_this(), &HMat);
219  }
220 
222 
223  // Pointer to arc length method DOF, used to scale pressure in LHS
224  boost::weak_ptr<DofEntity> arcLengthDof;
225 
227 
230  };
231 
232  /**
233  * @brief Operator for computing tangent vectors
234  *
235  */
236  struct OpGetTangent : public UserDataOperator {
237 
238  boost::shared_ptr<DataAtIntegrationPts> dataAtIntegrationPts;
239  OpGetTangent(const string field_name,
240  boost::shared_ptr<DataAtIntegrationPts> dataAtIntegrationPts)
241  : UserDataOperator(field_name, UserDataOperator::OPCOL),
243 
244  int ngp;
245  MoFEMErrorCode doWork(int side, EntityType type,
247  };
248 
249  /**
250  * @brief LHS-operator for pressure element (spatial configuration)
251  *
252  * Computes linearisation of the spatial component with respect to
253  * material coordinates.
254  *
255  */
257 
258  boost::shared_ptr<DataAtIntegrationPts> dataAtIntegrationPts;
259  Mat Aij;
262 
264 
265  /**
266  * @brief Compute left-hand side
267  *
268  * Computes linearisation of the spatial component with respect to
269  * material coordinates.
270  *
271  * Virtual work of the surface pressure corresponding to a test function
272  * of the spatial configuration \f$(\delta\mathbf{x})\f$:
273  * \f[
274  * \delta W^\text{spatial}_p(\mathbf{X}, \delta\mathbf{x}) =
275  * \int\limits_\mathcal{T} p\,\mathbf{N}(\mathbf{X}) \cdot
276  * \delta\mathbf{x}\, \textrm{d}\mathcal{T} =
277  * \int\limits_{\mathcal{T}_{\xi}}
278  * p\left(\frac{\partial\mathbf{X}}{\partial\xi}\times\frac{\partial
279  * \mathbf{X}} {\partial\eta}\right) \cdot \delta\mathbf{x}\,
280  * \textrm{d}\xi\textrm{d}\eta, \f] where \f$p\f$ is pressure,
281  * \f$\mathbf{N}\f$ is a normal to the face in the material configuration
282  * and \f$\xi, \eta\f$ are coordinates in the parent space
283  * \f$(\mathcal{T}_\xi)\f$.
284  *
285  * Linearisation with respect to a variation of material coordinates
286  * \f$(\Delta\mathbf{X})\f$:
287  *
288  * \f[
289  * \textrm{D} \delta W^\text{spatial}_p(\mathbf{X}, \delta\mathbf{x})
290  * [\Delta\mathbf{X}] = \int\limits_{\mathcal{T}_{\xi}} p\left[
291  * \frac{\partial\mathbf{X}}{\partial\xi} \cdot \left(\frac{\partial
292  * \Delta \mathbf{X}}{\partial\eta}\times\delta\mathbf{x}\right)
293  * -\frac{\partial\mathbf{X}}
294  * {\partial\eta} \cdot \left(\frac{\partial\Delta
295  * \mathbf{X}}{\partial\xi}\times \delta\mathbf{x}\right)\right]
296  * \textrm{d}\xi\textrm{d}\eta \f]
297  *
298  */
299  MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
300  EntityType col_type,
303 
305  const string field_name_1, const string field_name_2,
306  boost::shared_ptr<DataAtIntegrationPts> data_at_pts, Mat aij,
307  bCPressure &data, bool ho_geometry = false)
308  : UserDataOperator(field_name_1, field_name_2,
310  dataAtIntegrationPts(data_at_pts), Aij(aij), dAta(data),
311  hoGeometry(ho_geometry) {
312  sYmm = false; // This will make sure to loop over all entities
313  };
314  };
315 
316  /**
317  * @brief Operator for computing deformation gradients in side volumes
318  *
319  */
321 
322  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
324 
325  MoFEMErrorCode doWork(int side, EntityType type,
327 
328  OpCalculateDeformation(const string field_name,
329  boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
330  bool ho_geometry = false)
332  dataAtPts(data_at_pts), hoGeometry(ho_geometry) {
333  std::fill(&doEntities[MBEDGE], &doEntities[MBMAXTYPE], false);
334  sYmm = false;
335  };
336  };
337 
338  /**
339  * @brief RHS-operator for the pressure element (material configuration)
340  *
341  * Integrates pressure in the material configuration.
342  *
343  */
345 
346  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
347  boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> sideFe;
348  std::string sideFeName;
349  Vec F;
352 
354 
356 
357  int nbRows; ///< number of dofs on rows
358  int nbIntegrationPts; ///< number of integration points
359 
360  /**
361  * @brief Integrate pressure in the material configuration.
362  *
363  * Virtual work of the surface pressure corresponding to a test function
364  * of the material configuration \f$(\delta\mathbf{X})\f$:
365  *
366  * \f[
367  * \delta W^\text{material}_p(\mathbf{x}, \mathbf{X}, \delta\mathbf{X}) =
368  * -\int\limits_\mathcal{T} p\left\{\mathbf{F}^{\intercal}\cdot
369  * \mathbf{N}(\mathbf{X}) \right\} \cdot \delta\mathbf{X}\,
370  * \textrm{d}\mathcal{T} =
371  * -\int\limits_{\mathcal{T}_{\xi}} p\left\{\mathbf{F}^{\intercal}\cdot
372  * \left(\frac{\partial\mathbf{X}}{\partial\xi}\times\frac{\partial
373  * \mathbf{X}} {\partial\eta}\right) \right\} \cdot \delta\mathbf{X}\,
374  * \textrm{d}\xi\textrm{d}\eta
375  * \f]
376  *
377  * where \f$p\f$ is pressure, \f$\mathbf{N}\f$ is a normal to the face
378  * in the material configuration, \f$\xi, \eta\f$ are coordinates in the
379  * parent space
380  * \f$(\mathcal{T}_\xi)\f$ and \f$\mathbf{F}\f$ is the deformation gradient:
381  *
382  * \f[
383  * \mathbf{F} = \mathbf{h}(\mathbf{x})\,\mathbf{H}(\mathbf{X})^{-1} =
384  * \frac{\partial\mathbf{x}}{\partial\boldsymbol{\chi}}
385  * \frac{\partial\boldsymbol{\chi}}{\partial\mathbf{X}}
386  * \f]
387  *
388  * where \f$\mathbf{h}\f$ and \f$\mathbf{H}\f$ are the gradients of the
389  * spatial and material maps, respectively, and \f$\boldsymbol{\chi}\f$ are
390  * the reference coordinates.
391  *
392  */
393  MoFEMErrorCode doWork(int side, EntityType type,
395  MoFEMErrorCode iNtegrate(EntData &row_data);
396  MoFEMErrorCode aSsemble(EntData &row_data);
397 
399  const string material_field,
400  boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
401  boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> side_fe,
402  std::string &side_fe_name, Vec f, bCPressure &data,
403  bool ho_geometry = false)
404  : UserDataOperator(material_field, UserDataOperator::OPROW),
405  dataAtPts(data_at_pts), sideFe(side_fe), sideFeName(side_fe_name),
406  F(f), dAta(data), hoGeometry(ho_geometry){};
407  };
408 
409  /**
410  * @brief Base class for LHS-operators for pressure element (material
411  * configuration)
412  *
413  * Linearisation of the material component with respect to
414  * spatial and material coordinates consists of three parts, computed
415  * by operators working on the face and on the side volume:
416  *
417  * \f[
418  * \textrm{D} \delta W^\text{material}_p(\mathbf{x}, \mathbf{X},
419  * \delta\mathbf{x})
420  * [\Delta\mathbf{x}, \Delta\mathbf{X}] = \textrm{D} \delta
421  * W^\text{(face)}_p(\mathbf{x}, \mathbf{X}, \delta\mathbf{x})
422  * [\Delta\mathbf{X}] + \textrm{D} \delta
423  * W^\text{(side volume)}_p(\mathbf{x}, \mathbf{X}, \delta\mathbf{x})
424  * [\Delta\mathbf{x}] + \textrm{D} \delta W^\text{(side volume)}_p
425  * (\mathbf{x}, \mathbf{X}, \delta\mathbf{x}) [\Delta\mathbf{X}]
426  * \f]
427  *
428  */
430 
431  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
432  boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> sideFe;
433  std::string sideFeName;
434  Mat Aij;
437 
441 
445 
448 
450 
451  virtual MoFEMErrorCode doWork(int row_side, int col_side,
452  EntityType row_type, EntityType col_type,
457  }
458 
459  virtual MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data) {
462  }
463 
464  MoFEMErrorCode aSsemble(EntData &row_data, EntData &col_data);
465 
467  const string field_name_1, const string field_name_2,
468  boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
469  boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> side_fe,
470  std::string &side_fe_name, Mat aij, bCPressure &data,
471  bool ho_geometry = false)
472  : UserDataOperator(field_name_1, field_name_2,
474  dataAtPts(data_at_pts), sideFe(side_fe), sideFeName(side_fe_name),
475  Aij(aij), dAta(data), hoGeometry(ho_geometry) {
476  sYmm = false; // This will make sure to loop over all entities
477  }
478  };
479 
480  /**
481  * @brief LHS-operator for the pressure element (material configuration)
482  *
483  * Computes linearisation of the material component with respect to
484  * material coordinates (also triggers a loop over operators
485  * from the side volume).
486  *
487  */
490 
491  /**
492  * Integrates a contribution to the left-hand side and triggers a loop
493  * over side volume operators.
494  *
495  */
496  MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
497  EntityType col_type,
500 
501  /**
502  * @brief Compute part of the left-hand side
503  *
504  * Computes the linearisation of the material component
505  * with respect to a variation of material coordinates
506  * \f$(\Delta\mathbf{X})\f$:
507  *
508  * \f[
509  * \textrm{D} \delta W^\text{(face)}_p(\mathbf{x}, \mathbf{X},
510  * \delta\mathbf{x})
511  * [\Delta\mathbf{X}] = -\int\limits_{\mathcal{T}_{\xi}} p \,
512  * \mathbf{F}^{\intercal}\cdot \left[ \frac{\partial\mathbf{X}}
513  * {\partial\xi} \cdot \left(\frac{\partial\Delta
514  * \mathbf{X}}{\partial\eta}\times\delta\mathbf{x}\right)
515  * -\frac{\partial\mathbf{X}}
516  * {\partial\eta} \cdot \left(\frac{\partial\Delta
517  * \mathbf{X}}{\partial\xi}\times \delta\mathbf{x}\right)\right]
518  * \textrm{d}\xi\textrm{d}\eta
519  * \f]
520  *
521  */
522  MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data);
523 
525  const string field_name_1, const string field_name_2,
526  boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
527  boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> side_fe,
528  std::string &side_fe_name, Mat aij, bCPressure &data,
529  bool ho_geometry = false)
530  : OpNeumannPressureMaterialLhs(field_name_1, field_name_2, data_at_pts,
531  side_fe, side_fe_name, aij, data,
532  ho_geometry) {
533  sYmm = false; // This will make sure to loop over all entities
534  };
535  };
536 
537  /**
538  * @brief LHS-operator for the pressure element (material configuration)
539  *
540  * Triggers loop over operators from the side volume
541  *
542  */
545 
546  /*
547  * Triggers loop over operators from the side volume
548  *
549  */
550  MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
551  EntityType col_type,
554 
556  const string field_name_1, const string field_name_2,
557  boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
558  boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> side_fe,
559  std::string &side_fe_name, Mat aij, bCPressure &data,
560  bool ho_geometry = false)
561  : OpNeumannPressureMaterialLhs(field_name_1, field_name_2, data_at_pts,
562  side_fe, side_fe_name, aij, data,
563  ho_geometry) {
564  sYmm = false; // This will make sure to loop over all entities
565  };
566  };
567 
568  /**
569  * @brief Base class for LHS-operators (material) on side volumes
570  *
571  */
573  : public VolOnSideUserDataOperator {
574 
576 
577  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
578  Mat Aij;
581 
584 
588 
591 
593 
594  MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
595  EntityType col_type,
598  virtual MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data) {
601  }
602  MoFEMErrorCode aSsemble(EntData &row_data, EntData &col_data);
603 
605  const string field_name_1, const string field_name_2,
606  boost::shared_ptr<DataAtIntegrationPts> data_at_pts, Mat aij,
607  bCPressure &data, bool ho_geometry = false)
608  : VolOnSideUserDataOperator(field_name_1, field_name_2,
610  dataAtPts(data_at_pts), Aij(aij), dAta(data),
611  hoGeometry(ho_geometry) {
612  sYmm = false; // This will make sure to loop over all entities
613  }
614  };
615 
616  /**
617  * @brief LHS-operator (material configuration) on the side volume
618  *
619  * Computes the linearisation of the material component
620  * with respect to a variation of spatial coordinates on the side volume.
621  */
624 
625  /**
626  * @brief Integrates over a face contribution from a side volume
627  *
628  * Computes linearisation of the material component
629  * with respect to a variation of spatial coordinates:
630  *
631  * \f[
632  * \textrm{D} \delta W^\text{(side volume)}_p(\mathbf{x}, \mathbf{X},
633  * \delta\mathbf{x})
634  * [\Delta\mathbf{x}] = -\int\limits_{\mathcal{T}_{\xi}} p
635  * \left\{\left[
636  * \frac{\partial\Delta\mathbf{x}}{\partial\boldsymbol{\chi}}\,\mathbf{H}^{-1}
637  * \right]^{\intercal}\cdot\left(\frac{\partial\mathbf{X}}{\partial\xi}
638  * \times\frac{\partial\mathbf{X}}{\partial\eta}\right)\right\}
639  * \cdot \delta\mathbf{X}\, \textrm{d}\xi\textrm{d}\eta
640  * \f]
641  *
642  */
643  MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data);
644 
646  const string field_name_1, const string field_name_2,
647  boost::shared_ptr<DataAtIntegrationPts> data_at_pts, Mat aij,
648  bCPressure &data, bool ho_geometry = false)
650  field_name_1, field_name_2, data_at_pts, aij, data, ho_geometry) {
651  sYmm = false; // This will make sure to loop over all entities
652  };
653  };
654 
655  /**
656  * @brief LHS-operator (material configuration) on the side volume
657  *
658  * Computes the linearisation of the material component
659  * with respect to a variation of material coordinates on the side volume.
660  *
661  */
664 
665  /**
666  * @brief Integrates over a face contribution from a side volume
667  *
668  * Computes linearisation of the material component
669  * with respect to a variation of material coordinates:
670  *
671  * \f[
672  * \textrm{D} \delta W^\text{(side volume)}_p(\mathbf{x}, \mathbf{X},
673  * \delta\mathbf{x})
674  * [\Delta\mathbf{X}] = \int\limits_{\mathcal{T}_{\xi}} p
675  * \left\{\left[
676  * \mathbf{h}\,\mathbf{H}^{-1}\,\frac{\partial\Delta\mathbf{X}}
677  * {\partial\boldsymbol{\chi}}\,\mathbf{H}^{-1}
678  * \right]^{\intercal}\cdot\left(\frac{\partial\mathbf{X}}{\partial\xi}
679  * \times\frac{\partial\mathbf{X}}{\partial\eta}\right)\right\}
680  * \cdot \delta\mathbf{X}\, \textrm{d}\xi\textrm{d}\eta
681  * \f]
682  */
683  MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data);
684 
686  const string field_name_1, const string field_name_2,
687  boost::shared_ptr<DataAtIntegrationPts> data_at_pts, Mat aij,
688  bCPressure &data, bool ho_geometry = false)
690  field_name_1, field_name_2, data_at_pts, aij, data, ho_geometry) {
691  sYmm = false; // This will make sure to loop over all entities
692  };
693  };
694 
695  /// Operator for flux element
696  struct OpNeumannFlux : public UserDataOperator {
697 
698  Vec F;
700  boost::ptr_vector<MethodForForceScaling> &methodsOp;
702 
703  OpNeumannFlux(const std::string field_name, Vec _F, bCPressure &data,
704  boost::ptr_vector<MethodForForceScaling> &methods_op,
705  bool ho_geometry);
706 
708 
709  MoFEMErrorCode doWork(int side, EntityType type,
711  };
712 
713  /**
714  * \brief Add operator to calculate forces on element
715  * @param field_name Field name (f.e. TEMPERATURE)
716  * @param F Right hand side vector
717  * @param ms_id Set id (SideSet or BlockSet if block_set = true)
718  * @param ho_geometry Use higher order shape functions to define curved
719  * geometry
720  * @param block_set If tru get data from block set
721  * @return ErrorCode
722  */
723  MoFEMErrorCode addForce(const std::string field_name, Vec F, int ms_id,
724  bool ho_geometry = false, bool block_set = false);
725 
726  /**
727  * \brief Add operator to calculate pressure on element
728  * @param field_name Field name (f.e. TEMPERATURE)
729  * @param F Right hand side vector
730  * @param ms_id Set id (SideSet or BlockSet if block_set = true)
731  * @param ho_geometry Use higher order shape functions to define curved
732  * geometry
733  * @param block_set If true get data from block set
734  * @return ErrorCode
735  */
736  MoFEMErrorCode addPressure(const std::string field_name, Vec F, int ms_id,
737  bool ho_geometry = false, bool block_set = false);
738 
739  /**
740  * \brief Add operator to calculate pressure on element (in ALE)
741  * @param field_name_1 Field name for spatial positions
742  * @param field_name_2 Field name for material positions
743  * @param data_at_pts Common data at integration points
744  * @param side_fe_name Name of the element in the side volume
745  * @param F Right hand side vector
746  * @param aij Tangent matrix
747  * @param ms_id Set id (SideSet or BlockSet if block_set = true)
748  * @param ho_geometry Use higher order shape functions to define curved
749  * geometry
750  * @param block_set If true get data from block set
751  * @return ErrorCode
752  */
754  addPressureAle(const std::string field_name_1, const std::string field_name_2,
755  boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
756  std::string side_fe_name, Vec F, Mat aij, int ms_id,
757  bool ho_geometry = false, bool block_set = false);
758 
759  /**
760  * \brief Add operator to calculate pressure on element
761  * @param field_name Field name (f.e. TEMPERATURE)
762  * @param F Right hand side vector
763  * @param ms_id Set id (SideSet or BlockSet if block_set = true)
764  * @param ho_geometry Use higher order shape functions to define curved
765  * geometry
766  * @param block_set If tru get data from block set
767  * @return ErrorCode
768  */
769  MoFEMErrorCode addLinearPressure(const std::string field_name, Vec F,
770  int ms_id, bool ho_geometry = false);
771 
772  /// Add flux element operator (integration on face)
773  MoFEMErrorCode addFlux(const std::string field_name, Vec F, int ms_id,
774  bool ho_geometry = false);
775 
776  /// \deprecated fixed spelling mistake
778 
780 
781  DEPRECATED typedef bCPressure
782  bCPreassure; ///< \deprecated Do not use spelling mistake
783 
784  /// \deprecated function is deprecated because spelling mistake, use
785  /// addPressure instead
786  DEPRECATED MoFEMErrorCode addPreassure(const std::string field_name, Vec F,
787  int ms_id, bool ho_geometry = false,
788  bool block_set = false);
789 };
790 
791 /** \brief Set of high-level function declaring elements and setting operators
792  * to apply forces/fluxes \ingroup mofem_static_boundary_conditions
793  */
795 
796  /**
797  * \brief Declare finite element
798  *
799  * Search cubit sidesets and blocksets with pressure bc and declare surface
800  elemen
801 
802  * Block set has to have name “PRESSURE”. Can have name “PRESSURE_01” or any
803  * other name with prefix. The first attribute of block set is pressure
804  * value.
805 
806  *
807  * @param m_field Interface insurance
808  * @param field_name Field name (f.e. DISPLACEMENT)
809  * @param mesh_nodals_positions Name of field on which ho-geometry is defined
810  * @param intersect_ptr Pointer to range to interect meshset entities
811  * @return Error code
812  */
814  MoFEM::Interface &m_field, const std::string field_name,
815  const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS",
816  Range *intersect_ptr = NULL);
817 
818  /**
819  * \brief Set operators to finite elements calculating right hand side vector
820 
821  * @param m_field Interface
822  * @param neumann_forces Map of pointers to force/pressure elements
823  * @param F Right hand side vector
824  * @param field_name Field name (f.e. DISPLACEMENT)
825  * @param mesh_nodals_positions Name of field on which ho-geometry is defined
826  * @return Error code
827  *
828  */
830  MoFEM::Interface &m_field,
831  boost::ptr_map<std::string, NeumannForcesSurface> &neumann_forces, Vec F,
832  const std::string field_name,
833  const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS");
834 
836  MoFEM::Interface &m_field, const std::string field_name,
837  const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS");
838 
840  MoFEM::Interface &m_field,
841  boost::ptr_map<std::string, NeumannForcesSurface> &neumann_forces, Vec F,
842  const std::string field_name,
843  const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS");
844 };
845 
846 /**
847  * @deprecated Do not use that name it has spelling mistake
848  */
850 
851 /**
852  * @deprecated Do not use that name it has spelling mistake
853  */
855 
856 #endif //__SURFACE_PERSSURE_HPP__
857 
858 /**
859  * \defgroup mofem_static_boundary_conditions Pressure and force boundary
860  * conditions
861  *
862  * \ingroup user_modules
863  **/
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
Integrates over a face contribution from a side volume.
DEPRECATED typedef bCPressure bCPreassure
OpNeumannForceAnalytical(const std::string field_name, Vec f, const Range tris, boost::ptr_vector< MethodForForceScaling > &methods_op, boost::shared_ptr< MethodForAnalyticalForce > &analytical_force_op, const bool ho_geometry=false)
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
Deprecated interface functions.
NeumannForcesSurface(MoFEM::Interface &m_field)
OpNeumannPressureMaterialLhs_dX_dx(const string field_name_1, const string field_name_2, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > side_fe, std::string &side_fe_name, Mat aij, bCPressure &data, bool ho_geometry=false)
virtual MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
static MoFEMErrorCode addNeumannBCElements(MoFEM::Interface &m_field, const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS", Range *intersect_ptr=NULL)
Declare finite element.
boost::shared_ptr< MatrixDouble > getHMatPtr()
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
boost::ptr_vector< MethodForForceScaling > & methodsOp
OpNeumannPressureMaterialLhs_dX_dX(const string field_name_1, const string field_name_2, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > side_fe, std::string &side_fe_name, Mat aij, bCPressure &data, bool ho_geometry=false)
VolumeElementForcesAndSourcesCoreOnSideBase::UserDataOperator UserDataOperator
OpNeumannPressureMaterialLhs(const string field_name_1, const string field_name_2, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > side_fe, std::string &side_fe_name, Mat aij, bCPressure &data, bool ho_geometry=false)
MoFEMErrorCode aSsemble(EntData &row_data, EntData &col_data)
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
MoFEMErrorCode addLinearPressure(const std::string field_name, Vec F, int ms_id, bool ho_geometry=false)
Add operator to calculate pressure on element.
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &row_data)
Integrate pressure in the material configuration.
OpNeumannPressure(const std::string field_name, Vec _F, bCPressure &data, boost::ptr_vector< MethodForForceScaling > &methods_op, bool ho_geometry=false)
static Index< 'p', 3 > p
OpNeumannPressureMaterialVolOnSideLhs(const string field_name_1, const string field_name_2, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, Mat aij, bCPressure &data, bool ho_geometry=false)
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:76
Operator for computing deformation gradients in side volumes.
FTensor::Tensor1< double, 2 > normal(FTensor::Tensor1< T1, 3 > &t_coords, FTensor::Tensor1< T2, 2 > &t_disp)
Definition: ContactOps.hpp:137
OpCalculateDeformation(const string field_name, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, bool ho_geometry=false)
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
std::array< bool, MBMAXTYPE > doEntities
If true operator is executed for entity.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
Set of high-level function declaring elements and setting operators to apply forces/fluxes.
boost::ptr_vector< MethodForForceScaling > & methodsOp
std::map< int, bCPressure > mapPressure
OpNeumannPressureMaterialRhs_dX(const string material_field, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > side_fe, std::string &side_fe_name, Vec f, bCPressure &data, bool ho_geometry=false)
MyTriangleFE & getLoopFeLhs()
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
Integrates over a face contribution from a side volume.
boost::shared_ptr< MatrixDouble > getSmallhMatPtr()
static MoFEMErrorCode addNeumannFluxBCElements(MoFEM::Interface &m_field, const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
boost::ptr_vector< MethodForAnalyticalForce > analyticalForceOp
MyTriangleFE & getLoopFeMatRhs()
LHS-operator for the pressure element (material configuration)
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
boost::shared_ptr< DataAtIntegrationPts > dataAtIntegrationPts
virtual MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
MoFEMErrorCode addForce(const std::string field_name, Vec F, int ms_id, bool ho_geometry=false, bool block_set=false)
Add operator to calculate forces on element.
DEPRECATED typedef MethodForAnalyticalForce MethodForAnaliticalForce
OpNeumannFlux(const std::string field_name, Vec _F, bCPressure &data, boost::ptr_vector< MethodForForceScaling > &methods_op, bool ho_geometry)
virtual MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
MyTriangleFE & getLoopFeMatLhs()
bool sYmm
If true assume that matrix is symmetric structure.
RHS-operator for pressure element (spatial configuration)
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
boost::ptr_vector< MethodForForceScaling > & methodsOp
virtual MoFEMErrorCode getForce(const EntityHandle ent, const VectorDouble3 &coords, const VectorDouble3 &normal, VectorDouble3 &force)
ForcesAndSourcesCore::UserDataOperator UserDataOperator
Base class for LHS-operators (material) on side volumes.
DataForcesAndSourcesCore::EntData EntData
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
boost::weak_ptr< DofEntity > arcLengthDof
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
DEPRECATED typedef MetaNeumannForces MetaNeummanForces
MoFEMErrorCode addPressure(const std::string field_name, Vec F, int ms_id, bool ho_geometry=false, bool block_set=false)
Add operator to calculate pressure on element.
std::map< int, bCForce > mapForce
DEPRECATED typedef NeumannForcesSurface NeummanForcesSurface
MyTriangleFE & getLoopFe()
OpNeumannForce(const std::string field_name, Vec _F, bCForce &data, boost::ptr_vector< MethodForForceScaling > &methods_op, bool ho_geometry=false)
boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > sideFe
Operator for computing tangent vectors.
boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > sideFe
LHS-operator (material configuration) on the side volume.
MoFEMErrorCode addPressureAle(const std::string field_name_1, const std::string field_name_2, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, std::string side_fe_name, Vec F, Mat aij, int ms_id, bool ho_geometry=false, bool block_set=false)
Add operator to calculate pressure on element (in ALE)
OpNeumannPressureMaterialVolOnSideLhs_dX_dx(const string field_name_1, const string field_name_2, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, Mat aij, bCPressure &data, bool ho_geometry=false)
MoFEM::Interface & mField
boost::ptr_vector< MethodForForceScaling > methodsOp
static MoFEMErrorCode setMassFluxOperators(MoFEM::Interface &m_field, boost::ptr_map< std::string, NeumannForcesSurface > &neumann_forces, Vec F, const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
MoFEMErrorCode getForce(const EntityHandle ent, const VectorDouble3 &coords, const VectorDouble3 &normal, VectorDouble3 &force)
boost::shared_ptr< DataAtIntegrationPts > dataAtIntegrationPts
constexpr int order
ublas::vector< int, IntAllocator > VectorInt
Definition: Types.hpp:73
VectorBoundedArray< double, 3 > VectorDouble3
Definition: Types.hpp:92
#define DEPRECATED
Definition: definitions.h:29
MoFEMErrorCode aSsemble(EntData &row_data, EntData &col_data)
MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
Compute part of the left-hand side.
MoFEMErrorCode addFlux(const std::string field_name, Vec F, int ms_id, bool ho_geometry=false)
Add flux element operator (integration on face)
DataForcesAndSourcesCore::EntData EntData
LHS-operator for the pressure element (material configuration)
OpGetTangent(const string field_name, boost::shared_ptr< DataAtIntegrationPts > dataAtIntegrationPts)
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
static MoFEMErrorCode setMomentumFluxOperators(MoFEM::Interface &m_field, boost::ptr_map< std::string, NeumannForcesSurface > &neumann_forces, Vec F, const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
Set operators to finite elements calculating right hand side vector.
LHS-operator for pressure element (spatial configuration)
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
RHS-operator for the pressure element (material configuration)
Data on single entity (This is passed as argument to DataOperator::doWork)
LinearVaringPresssure(const VectorDouble3 &p, const double c)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
OpNeumannPressureMaterialVolOnSideLhs_dX_dX(const string field_name_1, const string field_name_2, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, Mat aij, bCPressure &data, bool ho_geometry=false)
boost::shared_ptr< MethodForAnalyticalForce > analyticalForceOp
DEPRECATED typedef OpNeumannPressure OpNeumannPreassure
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
OpNeumannPressureLhs_dx_dX(const string field_name_1, const string field_name_2, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, Mat aij, bCPressure &data, bool ho_geometry=false)
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:74
DEPRECATED MoFEMErrorCode addPreassure(const std::string field_name, Vec F, int ms_id, bool ho_geometry=false, bool block_set=false)
LHS-operator (material configuration) on the side volume.
FaceElementForcesAndSourcesCoreSwitch< 0 > FaceElementForcesAndSourcesCore
Face finite element default.
Base class for LHS-operators for pressure element (material configuration)
#define _F(n)
Definition: quad.c:25
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
Integrate pressure.
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
Compute left-hand side.
boost::ptr_vector< MethodForForceScaling > & methodsOp
MyTriangleFE(MoFEM::Interface &m_field)
Finite element and operators to apply force/pressures applied to surfaces.