v0.9.1
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 
93  // Pointer to arc length method DOF, used to scale pressure in LHS
94  boost::shared_ptr<DofEntity> arcLengthDof;
95 
97  boost::shared_ptr<DofEntity> arc_length_dof = nullptr)
98  : mField(m_field), fe(m_field), feLhs(m_field), feMatRhs(m_field),
99  feMatLhs(m_field), arcLengthDof(arc_length_dof) {}
100 
101  struct bCForce {
102  ForceCubitBcData data;
103  Range tRis;
104  };
105  std::map<int, bCForce> mapForce;
106 
107  struct bCPressure {
108  PressureCubitBcData data;
109  Range tRis;
110  };
111  std::map<int, bCPressure> mapPressure;
112 
113  boost::ptr_vector<MethodForForceScaling> methodsOp;
114  boost::ptr_vector<MethodForAnalyticalForce> analyticalForceOp;
115 
116  using UserDataOperator =
121 
122  /// Operator for force element
124 
125  Vec F;
127  boost::ptr_vector<MethodForForceScaling> &methodsOp;
128 
130 
131  OpNeumannForce(const std::string field_name, Vec _F, bCForce &data,
132  boost::ptr_vector<MethodForForceScaling> &methods_op,
133  bool ho_geometry = false);
134 
135  VectorDouble Nf; //< Local force vector
136 
137  MoFEMErrorCode doWork(int side, EntityType type,
139  };
140 
141  /// Operator for force element
143 
145  const std::string field_name, Vec f, const Range tris,
146  boost::ptr_vector<MethodForForceScaling> &methods_op,
147  boost::shared_ptr<MethodForAnalyticalForce> &analytical_force_op,
148  const bool ho_geometry = false);
149 
150  VectorDouble nF; //< Local force vector
151 
152  MoFEMErrorCode doWork(int side, EntityType type,
154 
155  Vec F;
156 
157  private:
158  const Range tRis;
159  boost::ptr_vector<MethodForForceScaling> &methodsOp;
160  boost::shared_ptr<MethodForAnalyticalForce> analyticalForceOp;
161  const bool hoGeometry;
162  };
163 
164  /**
165  * @brief RHS-operator for pressure element (spatial configuration)
166  *
167  */
169 
170  Vec F;
172  boost::ptr_vector<MethodForForceScaling> &methodsOp;
174 
175  OpNeumannPressure(const std::string field_name, Vec _F, bCPressure &data,
176  boost::ptr_vector<MethodForForceScaling> &methods_op,
177  bool ho_geometry = false);
178 
180 
181  /**
182  * @brief Integrate pressure
183  *
184  * \f[
185  * \begin{split}
186  * \mathbf{t} &= p \mathbf{n} \\
187  * \mathbf{f}^i &= \int_\mathcal{T} {\pmb\phi}^i \mathbf{t}
188  * \textrm{d}\mathcal{T} \end{split}
189  * \f]
190  * where \f$p\f$ is pressure, \f$\mathbf{n}\f$ is normal,
191  * \f$\mathbf{t}\f$ is traction, and
192  * \f$\mathbf{f}^i\f$ is local vector of external forces for ith base
193  * function \f${\pmb\phi}^i\f$.
194  *
195  *
196  * @param side
197  * @param type
198  * @param data
199  * @return MoFEMErrorCode
200  */
201  MoFEMErrorCode doWork(int side, EntityType type,
203  };
204 
206 
207  vector<vector<VectorDouble>> tangent;
208 
209  boost::shared_ptr<MatrixDouble> hMat;
210  boost::shared_ptr<MatrixDouble> FMat;
211  boost::shared_ptr<MatrixDouble> HMat;
212  boost::shared_ptr<VectorDouble> detHVec;
213  boost::shared_ptr<MatrixDouble> invHMat;
214 
216 
218  hMat = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
219  FMat = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
220  HMat = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
221  detHVec = boost::shared_ptr<VectorDouble>(new VectorDouble());
222  invHMat = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
223 
224  faceRowData = nullptr;
225  }
226 
229  };
230 
231  /**
232  * @brief Operator for computing tangent vectors
233  *
234  */
235  struct OpGetTangent : public UserDataOperator {
236 
237  boost::shared_ptr<DataAtIntegrationPts> dataAtIntegrationPts;
238  OpGetTangent(const string field_name,
239  boost::shared_ptr<DataAtIntegrationPts> dataAtIntegrationPts)
240  : UserDataOperator(field_name, UserDataOperator::OPCOL),
242 
243  int ngp;
244  MoFEMErrorCode doWork(int side, EntityType type,
246  };
247 
248  /**
249  * @brief LHS-operator for pressure element (spatial configuration)
250  *
251  * Computes linearisation of the spatial component with respect to
252  * material coordinates.
253  *
254  */
256 
257  boost::shared_ptr<DataAtIntegrationPts> dataAtIntegrationPts;
258  Mat Aij;
260  boost::shared_ptr<NeumannForcesSurface> surfacePressure;
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,
308  boost::shared_ptr<NeumannForcesSurface> surface_pressure,
309  bool ho_geometry = false)
310  : UserDataOperator(field_name_1, field_name_2,
312  dataAtIntegrationPts(data_at_pts), Aij(aij), dAta(data),
313  surfacePressure(surface_pressure), hoGeometry(ho_geometry) {
314  sYmm = false; // This will make sure to loop over all entities
315  };
316  };
317 
318  /**
319  * @brief Operator for computing deformation gradients in side volumes
320  *
321  */
323 
324  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
326 
327  MoFEMErrorCode doWork(int side, EntityType type,
329 
330  OpCalculateDeformation(const string field_name,
331  boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
332  bool ho_geometry = false)
334  dataAtPts(data_at_pts), hoGeometry(ho_geometry) {
335  std::fill(&doEntities[MBEDGE], &doEntities[MBMAXTYPE], false);
336  sYmm = false;
337  };
338  };
339 
340  /**
341  * @brief RHS-operator for the pressure element (material configuration)
342  *
343  * Integrates pressure in the material configuration.
344  *
345  */
347 
348  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
349  boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> sideFe;
350  std::string sideFeName;
351  Vec F;
354 
356 
358 
359  int nbRows; ///< number of dofs on rows
360  int nbIntegrationPts; ///< number of integration points
361 
362  /**
363  * @brief Integrate pressure in the material configuration.
364  *
365  * Virtual work of the surface pressure corresponding to a test function
366  * of the material configuration \f$(\delta\mathbf{X})\f$:
367  *
368  * \f[
369  * \delta W^\text{material}_p(\mathbf{x}, \mathbf{X}, \delta\mathbf{X}) =
370  * -\int\limits_\mathcal{T} p\left\{\mathbf{F}^{\intercal}\cdot
371  * \mathbf{N}(\mathbf{X}) \right\} \cdot \delta\mathbf{X}\,
372  * \textrm{d}\mathcal{T} =
373  * -\int\limits_{\mathcal{T}_{\xi}} p\left\{\mathbf{F}^{\intercal}\cdot
374  * \left(\frac{\partial\mathbf{X}}{\partial\xi}\times\frac{\partial
375  * \mathbf{X}} {\partial\eta}\right) \right\} \cdot \delta\mathbf{X}\,
376  * \textrm{d}\xi\textrm{d}\eta
377  * \f]
378  *
379  * where \f$p\f$ is pressure, \f$\mathbf{N}\f$ is a normal to the face
380  * in the material configuration, \f$\xi, \eta\f$ are coordinates in the
381  * parent space
382  * \f$(\mathcal{T}_\xi)\f$ and \f$\mathbf{F}\f$ is the deformation gradient:
383  *
384  * \f[
385  * \mathbf{F} = \mathbf{h}(\mathbf{x})\,\mathbf{H}(\mathbf{X})^{-1} =
386  * \frac{\partial\mathbf{x}}{\partial\mathbf{\chi}}
387  * \frac{\partial\mathbf{\chi}}{\partial\mathbf{X}}
388  * \f]
389  *
390  * where \f$\mathbf{h}\f$ and \f$\mathbf{H}\f$ are the gradients of the
391  * spatial and material maps, respectively, and \f$\mathbf{\chi}\f$ are
392  * the reference coordinates.
393  *
394  */
395  MoFEMErrorCode doWork(int side, EntityType type,
397  MoFEMErrorCode iNtegrate(EntData &row_data);
398  MoFEMErrorCode aSsemble(EntData &row_data);
399 
401  const string material_field,
402  boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
403  boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> side_fe,
404  std::string &side_fe_name, Vec f, bCPressure &data,
405  bool ho_geometry = false)
406  : UserDataOperator(material_field, UserDataOperator::OPROW),
407  dataAtPts(data_at_pts), sideFe(side_fe), sideFeName(side_fe_name),
408  F(f), dAta(data), hoGeometry(ho_geometry){};
409  };
410 
411  /**
412  * @brief Base class for LHS-operators for pressure element (material
413  * configuration)
414  *
415  * Linearisation of the material component with respect to
416  * spatial and material coordinates consists of three parts, computed
417  * by operators working on the face and on the side volume:
418  *
419  * \f[
420  * \textrm{D} \delta W^\text{material}_p(\mathbf{x}, \mathbf{X},
421  * \delta\mathbf{x})
422  * [\Delta\mathbf{x}, \Delta\mathbf{X}] = \textrm{D} \delta
423  * W^\text{(face)}_p(\mathbf{x}, \mathbf{X}, \delta\mathbf{x})
424  * [\Delta\mathbf{X}] + \textrm{D} \delta
425  * W^\text{(side volume)}_p(\mathbf{x}, \mathbf{X}, \delta\mathbf{x})
426  * [\Delta\mathbf{x}] + \textrm{D} \delta W^\text{(side volume)}_p
427  * (\mathbf{x}, \mathbf{X}, \delta\mathbf{x}) [\Delta\mathbf{X}]
428  * \f]
429  *
430  */
432 
433  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
434  boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> sideFe;
435  std::string sideFeName;
436  Mat Aij;
438  boost::shared_ptr<NeumannForcesSurface> surfacePressure;
440 
444 
448 
451 
453 
454  virtual MoFEMErrorCode doWork(int row_side, int col_side,
455  EntityType row_type, EntityType col_type,
460  }
461 
462  virtual MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data) {
465  }
466 
467  MoFEMErrorCode aSsemble(EntData &row_data, EntData &col_data);
468 
470  const string field_name_1, const string field_name_2,
471  boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
472  boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> side_fe,
473  std::string &side_fe_name, Mat aij, bCPressure &data,
474  boost::shared_ptr<NeumannForcesSurface> surface_pressure,
475  bool ho_geometry = false)
476  : UserDataOperator(field_name_1, field_name_2,
478  dataAtPts(data_at_pts), sideFe(side_fe), sideFeName(side_fe_name),
479  Aij(aij), dAta(data), surfacePressure(surface_pressure),
480  hoGeometry(ho_geometry) {
481  sYmm = false; // This will make sure to loop over all entities
482  }
483  };
484 
485  /**
486  * @brief LHS-operator for the pressure element (material configuration)
487  *
488  * Computes linearisation of the material component with respect to
489  * material coordinates (also triggers a loop over operators
490  * from the side volume).
491  *
492  */
495 
496  /**
497  * Integrates a contribution to the left-hand side and triggers a loop
498  * over side volume operators.
499  *
500  */
501  MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
502  EntityType col_type,
505 
506  /**
507  * @brief Compute part of the left-hand side
508  *
509  * Computes the linearisation of the material component
510  * with respect to a variation of material coordinates
511  * \f$(\Delta\mathbf{X})\f$:
512  *
513  * \f[
514  * \textrm{D} \delta W^\text{(face)}_p(\mathbf{x}, \mathbf{X},
515  * \delta\mathbf{x})
516  * [\Delta\mathbf{X}] = -\int\limits_{\mathcal{T}_{\xi}} p \,
517  * \mathbf{F}^{\intercal}\cdot \left[ \frac{\partial\mathbf{X}}
518  * {\partial\xi} \cdot \left(\frac{\partial\Delta
519  * \mathbf{X}}{\partial\eta}\times\delta\mathbf{x}\right)
520  * -\frac{\partial\mathbf{X}}
521  * {\partial\eta} \cdot \left(\frac{\partial\Delta
522  * \mathbf{X}}{\partial\xi}\times \delta\mathbf{x}\right)\right]
523  * \textrm{d}\xi\textrm{d}\eta
524  * \f]
525  *
526  */
527  MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data);
528 
530  const string field_name_1, const string field_name_2,
531  boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
532  boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> side_fe,
533  std::string &side_fe_name, Mat aij, bCPressure &data,
534  boost::shared_ptr<NeumannForcesSurface> surface_pressure,
535  bool ho_geometry = false)
536  : OpNeumannPressureMaterialLhs(field_name_1, field_name_2, data_at_pts,
537  side_fe, side_fe_name, aij, data,
538  surface_pressure, ho_geometry) {
539  sYmm = false; // This will make sure to loop over all entities
540  };
541  };
542 
543  /**
544  * @brief LHS-operator for the pressure element (material configuration)
545  *
546  * Triggers loop over operators from the side volume
547  *
548  */
551 
552  /*
553  * Triggers loop over operators from the side volume
554  *
555  */
556  MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
557  EntityType col_type,
560 
562  const string field_name_1, const string field_name_2,
563  boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
564  boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> side_fe,
565  std::string &side_fe_name, Mat aij, bCPressure &data,
566  boost::shared_ptr<NeumannForcesSurface> surface_pressure,
567  bool ho_geometry = false)
568  : OpNeumannPressureMaterialLhs(field_name_1, field_name_2, data_at_pts,
569  side_fe, side_fe_name, aij, data,
570  surface_pressure, ho_geometry) {
571  sYmm = false; // This will make sure to loop over all entities
572  };
573  };
574 
575  /**
576  * @brief Base class for LHS-operators (material) on side volumes
577  *
578  */
580  : public VolOnSideUserDataOperator {
581 
583 
584  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
585  Mat Aij;
587  boost::shared_ptr<NeumannForcesSurface> surfacePressure;
589 
592 
596 
599 
601 
602  MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
603  EntityType col_type,
606  virtual MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data) {
609  }
610  MoFEMErrorCode aSsemble(EntData &row_data, EntData &col_data);
611 
613  const string field_name_1, const string field_name_2,
614  boost::shared_ptr<DataAtIntegrationPts> data_at_pts, Mat aij,
615  bCPressure &data,
616  boost::shared_ptr<NeumannForcesSurface> surface_pressure,
617  bool ho_geometry = false)
618  : VolOnSideUserDataOperator(field_name_1, field_name_2,
620  dataAtPts(data_at_pts), Aij(aij), dAta(data),
621  surfacePressure(surface_pressure), hoGeometry(ho_geometry) {
622  sYmm = false; // This will make sure to loop over all entities
623  }
624  };
625 
626  /**
627  * @brief LHS-operator (material configuration) on the side volume
628  *
629  * Computes the linearisation of the material component
630  * with respect to a variation of spatial coordinates on the side volume.
631  */
634 
635  /**
636  * @brief Integrates over a face contribution from a side volume
637  *
638  * Computes linearisation of the material component
639  * with respect to a variation of spatial coordinates:
640  *
641  * \f[
642  * \textrm{D} \delta W^\text{(side volume)}_p(\mathbf{x}, \mathbf{X},
643  * \delta\mathbf{x})
644  * [\Delta\mathbf{x}] = -\int\limits_{\mathcal{T}_{\xi}} p
645  * \left\{\left[
646  * \frac{\partial\Delta\mathbf{x}}{\partial\mathbf{\chi}}\,\mathbf{H}^{-1}
647  * \right]^{\intercal}\cdot\left(\frac{\partial\mathbf{X}}{\partial\xi}
648  * \times\frac{\partial\mathbf{X}}{\partial\eta}\right)\right\}
649  * \cdot \delta\mathbf{X}\, \textrm{d}\xi\textrm{d}\eta
650  * \f]
651  *
652  */
653  MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data);
654 
656  const string field_name_1, const string field_name_2,
657  boost::shared_ptr<DataAtIntegrationPts> data_at_pts, Mat aij,
658  bCPressure &data,
659  boost::shared_ptr<NeumannForcesSurface> surface_pressure,
660  bool ho_geometry = false)
661  : OpNeumannPressureMaterialVolOnSideLhs(field_name_1, field_name_2,
662  data_at_pts, aij, data,
663  surface_pressure, ho_geometry) {
664  sYmm = false; // This will make sure to loop over all entities
665  };
666  };
667 
668  /**
669  * @brief LHS-operator (material configuration) on the side volume
670  *
671  * Computes the linearisation of the material component
672  * with respect to a variation of material coordinates on the side volume.
673  *
674  */
677 
678  /**
679  * @brief Integrates over a face contribution from a side volume
680  *
681  * Computes linearisation of the material component
682  * with respect to a variation of material coordinates:
683  *
684  * \f[
685  * \textrm{D} \delta W^\text{(side volume)}_p(\mathbf{x}, \mathbf{X},
686  * \delta\mathbf{x})
687  * [\Delta\mathbf{X}] = \int\limits_{\mathcal{T}_{\xi}} p
688  * \left\{\left[
689  * \mathbf{h}\,\mathbf{H}^{-1}\,\frac{\partial\Delta\mathbf{X}}
690  * {\partial\mathbf{\chi}}\,\mathbf{H}^{-1}
691  * \right]^{\intercal}\cdot\left(\frac{\partial\mathbf{X}}{\partial\xi}
692  * \times\frac{\partial\mathbf{X}}{\partial\eta}\right)\right\}
693  * \cdot \delta\mathbf{X}\, \textrm{d}\xi\textrm{d}\eta
694  * \f]
695  */
696  MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data);
697 
699  const string field_name_1, const string field_name_2,
700  boost::shared_ptr<DataAtIntegrationPts> data_at_pts, Mat aij,
701  bCPressure &data,
702  boost::shared_ptr<NeumannForcesSurface> surface_pressure,
703  bool ho_geometry = false)
704  : OpNeumannPressureMaterialVolOnSideLhs(field_name_1, field_name_2,
705  data_at_pts, aij, data,
706  surface_pressure, ho_geometry) {
707  sYmm = false; // This will make sure to loop over all entities
708  };
709  };
710 
711  /// Operator for flux element
712  struct OpNeumannFlux : public UserDataOperator {
713 
714  Vec F;
716  boost::ptr_vector<MethodForForceScaling> &methodsOp;
718 
719  OpNeumannFlux(const std::string field_name, Vec _F, bCPressure &data,
720  boost::ptr_vector<MethodForForceScaling> &methods_op,
721  bool ho_geometry);
722 
724 
725  MoFEMErrorCode doWork(int side, EntityType type,
727  };
728 
729  /**
730  * \brief Add operator to calculate forces on element
731  * @param field_name Field name (f.e. TEMPERATURE)
732  * @param F Right hand side vector
733  * @param ms_id Set id (SideSet or BlockSet if block_set = true)
734  * @param ho_geometry Use higher order shape functions to define curved
735  * geometry
736  * @param block_set If tru get data from block set
737  * @return ErrorCode
738  */
739  MoFEMErrorCode addForce(const std::string field_name, Vec F, int ms_id,
740  bool ho_geometry = false, bool block_set = false);
741 
742  /**
743  * \brief Add operator to calculate pressure on element
744  * @param field_name Field name (f.e. TEMPERATURE)
745  * @param F Right hand side vector
746  * @param ms_id Set id (SideSet or BlockSet if block_set = true)
747  * @param ho_geometry Use higher order shape functions to define curved
748  * geometry
749  * @param block_set If true get data from block set
750  * @return ErrorCode
751  */
752  MoFEMErrorCode addPressure(const std::string field_name, Vec F, int ms_id,
753  bool ho_geometry = false, bool block_set = false);
754 
755  /**
756  * \brief Add operator to calculate pressure on element (in ALE)
757  * @param field_name_1 Field name for spatial positions
758  * @param field_name_2 Field name for material positions
759  * @param data_at_pts Common data at integration points
760  * @param side_fe_name Name of the element in the side volume
761  * @param F Right hand side vector
762  * @param aij Tangent matrix
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 true get data from block set
767  * @return ErrorCode
768  */
770  addPressureAle(const std::string field_name_1, const std::string field_name_2,
771  boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
772  std::string side_fe_name, Vec F, Mat aij, int ms_id,
773  boost::shared_ptr<NeumannForcesSurface> surface_pressure,
774  bool ho_geometry = false, bool block_set = false);
775 
776  /**
777  * \brief Add operator to calculate pressure on element
778  * @param field_name Field name (f.e. TEMPERATURE)
779  * @param F Right hand side vector
780  * @param ms_id Set id (SideSet or BlockSet if block_set = true)
781  * @param ho_geometry Use higher order shape functions to define curved
782  * geometry
783  * @param block_set If tru get data from block set
784  * @return ErrorCode
785  */
786  MoFEMErrorCode addLinearPressure(const std::string field_name, Vec F,
787  int ms_id, bool ho_geometry = false);
788 
789  /// Add flux element operator (integration on face)
790  MoFEMErrorCode addFlux(const std::string field_name, Vec F, int ms_id,
791  bool ho_geometry = false);
792 
793  /// \deprecated fixed spelling mistake
795 
797 
798  DEPRECATED typedef bCPressure
799  bCPreassure; ///< \deprecated Do not use spelling mistake
800 
801  /// \deprecated function is deprecated because spelling mistake, use
802  /// addPressure instead
803  DEPRECATED MoFEMErrorCode addPreassure(const std::string field_name, Vec F,
804  int ms_id, bool ho_geometry = false,
805  bool block_set = false);
806 };
807 
808 /** \brief Set of high-level function declaring elements and setting operators
809  * to apply forces/fluxes \ingroup mofem_static_boundary_conditions
810  */
812 
813  /**
814  * \brief Declare finite element
815  *
816  * Search cubit sidesets and blocksets with pressure bc and declare surface
817  elemen
818 
819  * Block set has to have name “PRESSURE”. Can have name “PRESSURE_01” or any
820  * other name with prefix. The first attribute of block set is pressure
821  * value.
822 
823  *
824  * @param m_field Interface insurance
825  * @param field_name Field name (f.e. DISPLACEMENT)
826  * @param mesh_nodals_positions Name of field on which ho-geometry is defined
827  * @param intersect_ptr Pointer to range to interect meshset entities
828  * @return Error code
829  */
831  MoFEM::Interface &m_field, const std::string field_name,
832  const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS",
833  Range *intersect_ptr = NULL);
834 
835  /**
836  * \brief Set operators to finite elements calculating right hand side vector
837 
838  * @param m_field Interface
839  * @param neumann_forces Map of pointers to force/pressure elements
840  * @param F Right hand side vector
841  * @param field_name Field name (f.e. DISPLACEMENT)
842  * @param mesh_nodals_positions Name of field on which ho-geometry is defined
843  * @return Error code
844  *
845  */
847  MoFEM::Interface &m_field,
848  boost::ptr_map<std::string, NeumannForcesSurface> &neumann_forces, Vec F,
849  const std::string field_name,
850  const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS");
851 
853  MoFEM::Interface &m_field, const std::string field_name,
854  const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS");
855 
857  MoFEM::Interface &m_field,
858  boost::ptr_map<std::string, NeumannForcesSurface> &neumann_forces, Vec F,
859  const std::string field_name,
860  const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS");
861 };
862 
863 /**
864  * @deprecated Do not use that name it has spelling mistake
865  */
867 
868 /**
869  * @deprecated Do not use that name it has spelling mistake
870  */
872 
873 #endif //__SURFACE_PERSSURE_HPP__
874 
875 /**
876  * \defgroup mofem_static_boundary_conditions Pressure and force boundary
877  * conditions
878  *
879  * \ingroup user_modules
880  **/
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
NeumannForcesSurface(MoFEM::Interface &m_field, boost::shared_ptr< DofEntity > arc_length_dof=nullptr)
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.
boost::shared_ptr< NeumannForcesSurface > surfacePressure
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< DataAtIntegrationPts > dataAtPts
boost::ptr_vector< MethodForForceScaling > & methodsOp
VolumeElementForcesAndSourcesCoreOnSideBase::UserDataOperator UserDataOperator
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:506
MoFEMErrorCode addLinearPressure(const std::string field_name, Vec F, int ms_id, bool ho_geometry=false)
Add operator to calculate pressure on element.
boost::shared_ptr< MatrixDouble > HMat
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)
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:74
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, boost::shared_ptr< NeumannForcesSurface > surface_pressure, bool ho_geometry=false)
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:128
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:482
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:513
Set of high-level function declaring elements and setting operators to apply forces/fluxes.
boost::ptr_vector< MethodForForceScaling > & methodsOp
std::map< int, bCPressure > mapPressure
vector< vector< VectorDouble > > tangent
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.
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)
boost::shared_ptr< NeumannForcesSurface > surfacePressure
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)
Base class for LHS-operators (material) on side volumes.
OpNeumannPressureMaterialVolOnSideLhs_dX_dx(const string field_name_1, const string field_name_2, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, Mat aij, bCPressure &data, boost::shared_ptr< NeumannForcesSurface > surface_pressure, bool ho_geometry=false)
boost::shared_ptr< MatrixDouble > hMat
constexpr int order
DataForcesAndSourcesCore::EntData EntData
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
boost::shared_ptr< DofEntity > arcLengthDof
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
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, boost::shared_ptr< NeumannForcesSurface > surface_pressure, bool ho_geometry=false)
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
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, boost::shared_ptr< NeumannForcesSurface > surface_pressure, bool ho_geometry=false, bool block_set=false)
Add operator to calculate pressure on element (in ALE)
boost::shared_ptr< MatrixDouble > FMat
Operator for computing tangent vectors.
boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > sideFe
boost::shared_ptr< NeumannForcesSurface > surfacePressure
LHS-operator (material configuration) on the side volume.
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
ForcesAndSourcesCore::UserDataOperator UserDataOperator
ublas::vector< int, IntAllocator > VectorInt
Definition: Types.hpp:71
VectorBoundedArray< double, 3 > VectorDouble3
Definition: Types.hpp:87
#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
boost::shared_ptr< VectorDouble > detHVec
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
boost::shared_ptr< MethodForAnalyticalForce > analyticalForceOp
DEPRECATED typedef OpNeumannPressure OpNeumannPreassure
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, boost::shared_ptr< NeumannForcesSurface > surface_pressure, bool ho_geometry=false)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:72
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.
Base class for LHS-operators for pressure element (material configuration)
#define _F(n)
Definition: quad.c:25
boost::shared_ptr< MatrixDouble > invHMat
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
OpNeumannPressureLhs_dx_dX(const string field_name_1, const string field_name_2, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, Mat aij, bCPressure &data, boost::shared_ptr< NeumannForcesSurface > surface_pressure, bool ho_geometry=false)
OpNeumannPressureMaterialVolOnSideLhs_dX_dX(const string field_name_1, const string field_name_2, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, Mat aij, bCPressure &data, boost::shared_ptr< NeumannForcesSurface > surface_pressure, bool ho_geometry=false)
MyTriangleFE(MoFEM::Interface &m_field)
Finite element and operators to apply force/pressures applied to surfaces.
OpNeumannPressureMaterialVolOnSideLhs(const string field_name_1, const string field_name_2, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, Mat aij, bCPressure &data, boost::shared_ptr< NeumannForcesSurface > surface_pressure, bool ho_geometry=false)