Xfem
The Xfem Plugin
xfemHMCoupled.h
Go to the documentation of this file.
1 /************************************************************************
2 **
3 ** Copyright (C) 2017 by Geomechanic group - Tecgraf/Puc-Rio
4 ** All rights reserved.
5 **
6 ** This file is part of the "GeMA" software. It's use should respect
7 ** the terms in the license agreement that can be found together
8 ** with this source code.
9 ** It is provided AS IS, with NO WARRANTY OF ANY KIND,
10 ** INCLUDING THE WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR
11 ** A PARTICULAR PURPOSE.
12 **
13 ************************************************************************/
14 
24 #ifndef _GEMA_XFEM_HMCOUPLED_H_
25 #define _GEMA_XFEM_HMCOUPLED_H_
26 
27 #include <gmpMechanicalPhysics.h>
28 #include <gmpMechanicPoint.h>
29 
30 #include "xfemMechanic.h"
31 
32 class GmGaussAccessor;
33 
37 {
38 public:
39 
40  // constructor
41  XfemHMCoupled(const char* pluginType, GmSimulationData* simulation, QString id, QString description,
42  const GmpFemPhysicsCommonMaterialFactory* matFactory, const GmLogCategory& logger);
43 
44  // destructor
45  virtual ~XfemHMCoupled();
46 
47  // set physic map
48  virtual const QVariantMap* physicsMetaDataMap();
49  // Comments on the base class
50  virtual const char* pluginName() const { return "Xfem"; }
51  // Comments on the base class
52  virtual const GmElementDof* dofMapping(const GmElement* e) const;
53  // Comments on the base class
54  virtual bool fixedNodalForcesBc(QList<int>& nodes, QList<int>&dof, QList<double>& values) const;
55  // Comments on the base class
56  virtual bool fixedNodalDofsBc(QList<int>& nodes, QList<int>&dof, QList<double>& values, bool* constantValues) const;
57 
58  // Comments on the base class
59  virtual FemResultType fillElementDataForBc(const GmElement* e, const GmBoundaryCondition* bc, int bcIndex, int bcListIndex,
60  int border, GmpFemMatrixSet& elemMatrices, GmpFemVectorSet& elemVectors);
61 
62  // Comments on the base class
63  virtual FemResultType fillElementData(const GmElement* e, GmpFemMatrixSet& elemMatrices, GmpFemVectorSet& elemVectors);
64 
65 protected:
68  {
75  T_ID,
78 
79  // ------ NO ADDING BELOW THIS LINE
81  };
82 
85  {
87 
88  // ------ NO ADDING BELOW THIS LINE
90  };
91 
94  {
111 
112  // ------ NO ADDING BELOW THIS LINE
114  };
115 
118  {
138 
139  // ------ NO ADDING BELOW THIS LINE
141  };
142 
143  virtual bool checkAndLoadPrivateData(LuaTable& table);
144  virtual bool checkAndLoadAttributeAccessors(LuaTable& nodeTable, LuaTable& gaussTable);
145  virtual bool fillElementDataMPNIT(const GmElement* e, GmpFemMatrixSet& elemMatrices, GmpFemVectorSet& elemVectors);
146  virtual bool fillElementDataNIT(const GmElement* e, GmpFemMatrixSet& elemMatrices, GmpFemVectorSet& elemVectors);
147  virtual bool fillElementDataMP(const GmElement* e, GmpFemMatrixSet& elemMatrices, GmpFemVectorSet& elemVectors);
148  void fillElementDisplacements(const GmElement* e, GmVector& ue, const int& nEnrich);
149  void fillElementPorePressure(const GmElement* e, GmVector& pe, GmVector& pf, bool ghostNodes, const QList<int>& pfList);
150  void fillElementTemperature(const GmElement* e, GmVector& pe, GmVector& pf, bool ghostNodes, const QList<int>& pfList);
151  void fillElementCapillaryPressure(const GmElement* e, GmVector& pe, GmVector& pf, bool ghostNodes, const QList<int>& pfList);
152  void fillEnhancedBuMatrix(const GmElement* e, const GmShape* shape, const GmVector& ncoord, const GmMatrix& X,
153  const GmMatrix& J, GmMatrix& Bu, const QList<int>& gaussLevelSet, const int& nEnrich, const GmVector& nodes, GmElementDof* dofMap, XfemEnrichedElementData* exData, const int& gaussIP);
154  void fillEnhancedBuJumpMatrix(const GmElement* e, const GmShape* shape, const GmVector& ncoord, const GmMatrix& J, GmMatrix& Bu, const QList<int>& TopLevelSet, const QList<int>& BottomLevelSet, const int& nEnrich, const GmVector& nodes, GmElementDof* dofMap, XfemEnrichedElementData* exData, const int& iSubCrack);
155  void fillEnhancedBpJumpMatrix(const GmElement* e, const GmShape* shape, const GmVector& ncoord, const GmMatrix& J, GmMatrix& Bp, const QList<int>& TopLevelSet, const QList<int>& BottomLevelSet, const int& nEnrich, GmElementDof* dofMap, XfemEnrichedElementData* exData, const int& iSubCrack);
156  void fillEnhancedNpMatrix(const GmElement* e, GmVector& Np, const GmVector& N, const QList<int>& gaussLevelSet, const int& nEnrich, GmElementDof* dofMap, XfemEnrichedElementData* exData, const int& gaussIP);
157  void fillEnhancedNptopMatrix(const GmElement* e, GmVector& Np, const GmVector& N, const QList<int>& gaussLevelSet, const int& nEnrich, GmElementDof* dofMap, XfemEnrichedElementData* exData, const int& iSubCrack);
158  void fillEnhancedNpbotMatrix(const GmElement* e, GmVector& Np, const GmVector& N, const QList<int>& gaussLevelSet, const int& nEnrich, GmElementDof* dofMap, XfemEnrichedElementData* exData, const int& iSubCrack);
159  void fillEnhancedNutopMatrix(const GmElement* e, GmMatrix& Nu, const GmVector& N, const QList<int>& gaussLevelSet, const int& nEnrich, const GmVector& nodes, GmElementDof* dofMap, XfemEnrichedElementData* exData, const int& iSubCrack);
160  void fillEnhancedNubotMatrix(const GmElement* e, GmMatrix& Nu, const GmVector& N, const QList<int>& gaussLevelSet, const int& nEnrich, const GmVector& nodes, GmElementDof* dofMap, XfemEnrichedElementData* exData, const int& iSubCrack);
161  void fillEnhancedNuJumpMatrix(const GmElement* e, GmMatrix& Nu, const GmVector& N, const QList<int>& BottomLevelSet, const QList<int>& TopLevelSet, const int& nEnrich, const GmVector& nodes, GmElementDof* dofMap, XfemEnrichedElementData* exData, const int& iSubCrack);
162  void fillEnhancedNpJumpMatrix(const GmElement* e, GmMatrix& Np, const GmVector& N, const QList<int>& BottomLevelSet, const QList<int>& TopLevelSet, const int& nEnrich, GmElementDof* dofMap, XfemEnrichedElementData* exData, const int& iSubCrack);
163  void tipDOFTreat(const int& nelem, GmVector& nodes);
164  void shapeValuesBarElem(double xi, GmVector& N);
165  void unitVectSegment(const GmMatrix& Xg, GmMatrix& uVec, GmMatrix& uVec2, GmMatrix& uVec3, GmMatrix& ROT);
166  virtual void fillStandardDof(const GmElement* e, QList<int>& iUe, QList<int>& iPe);
167  virtual void fillEnhancedDof(const XfemElement* xe, QList<int>& iUe, QList<int>& iPe, QList<int>& iPf, const int& np);
168  virtual void fillStandardDofMP(const GmElement* e, QList<int>& iUe, QList<int>& iPw, QList<int>& iPc);
169  virtual void fillEnhancedDofMP(const XfemElement* xe, QList<int>& iUe, QList<int>& iPw, QList<int>& iPc, QList<int>& iPwd, QList<int>& iPcd, const int& np);
170  virtual void fillStandardDofNIT(const GmElement* e, QList<int>& iUe, QList<int>& iPw, QList<int>& iPc);
171  virtual void fillEnhancedDofNIT(const XfemElement* xe, QList<int>& iUe, QList<int>& iPw, QList<int>& iPc, QList<int>& iPwd, QList<int>& iPcd, const int& np);
172  virtual void fillStorageCompMatrix(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Sc, const GmVector& Np, const double& c);
173  virtual void fillCrackStorageCompMatrix(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Sc, const GmVector& Np, const double& c);
174  virtual void fillHeatCapacityMatrix(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Sc, const GmVector& Np, const double& c);
175  virtual void fillCrackHeatCapacityMatrix(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Sc, const GmVector& Np, const double& c, double& wn);
176  virtual void fillThermalExpansionMatrix(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Sc, const GmVector& Np, const double& c);
177  virtual void fillCrackThermalExpansionMatrix(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Sc, const GmVector& Np, const double& c, double& wn);
178  virtual void fillStorageCompMatrixMP(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Cww, GmMatrix& Cwc, GmMatrix& Ccc, GmMatrix& Ccw, const GmVector& Np, const double& c, const double& sw);
179  virtual void fillCrackStorageCompMatrixMP(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Cww, GmMatrix& Cwc, GmMatrix& Ccc, GmMatrix& Ccw, const GmVector& Np, const double& c, const double& sw, double& wn);
180  virtual void fillCouplingMecHydMatrix(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Lup, const GmVector& Np, const GmMatrix& Bu, const double& c);
181  virtual void fillCouplingCrackMecHydMatrix(const GmElement* e, const GmVector* coord, const int& ip, const double& c2, const GmVector& Np, GmMatrix& Qpfu, const GmVector& v_n, const GmMatrix& NuJump);
182  virtual void fillCouplingCrackHydMecMatrix(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Epfu, const GmVector& Np, const GmMatrix& Bu, const GmVector& v_t, const double& c2, double& wn, GmMatrix& Qpfu, const GmVector& v_n, const GmMatrix& NuJump);
183  virtual void fillCouplingHydMecMatrix(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Lup, const GmVector& Np, const GmMatrix& Bu, const double& c);
184  virtual void fillCouplingThermoMecMatrix(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Lup, const GmVector& Np, const GmMatrix& Bu, const double& c, const double& Tg);
185  virtual void fillCouplingCrackThermoMecMatrix(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Lup, const GmVector& Np, const GmMatrix& Bu, const double& c, const double& Tg, const double& wn);
186  void fillEnhancedBpMatrix(const GmElement* e, const GmShape* shape, const GmVector& ncoord,
187  const GmMatrix& X, const GmMatrix& J, const GmMatrix& B, GmMatrix& Bp, const QList<int>& gaussLevelSet, const int& nEnrich, GmElementDof* dofMap, XfemEnrichedElementData* exData, const int& gaussIP);
188  void fillFlowMatrix(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Kc, const GmMatrix& Bp, const double& c);
189  void fillConvectiveMatrix(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Kc, const GmMatrix& Bp, const GmMatrix& Np, const double& c, const GmVector& pe);
190  void fillConvectiveMatrix2(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Kc, const GmMatrix& Bp, const double& tg, const double& c);
191  void fillAdvectiveMatrix(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Kc, const GmMatrix& Bp, const double& c);
192  void fillCrackConvectiveMatrix(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Kc, const GmMatrix& Bp, const GmMatrix& Np, const double& c, double& wn, double& kfd, const GmVector& pfe);
193  void fillCrackConvectiveMatrix2(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Kc, const GmMatrix& Bp, const double& tg, const double& c, double& wn, double& kfd);
194  void fillCrackAdvectiveMatrix(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Kc, const GmMatrix& Bp, const double& c, double& wn);
195  void fillFlowMatrixMP(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Kc, const GmMatrix& Bp, const double& c, const double& krw, const double& krn);
196  void fillFlowMatrixMP2(const GmElement* e, const GmVector* coord, const int& ip, GmMatrix& Kc, const GmMatrix& Bp, const double& c, const double& krn);
197  virtual bool assemblyElemDampingMatrix(const GmElement* e, GmMatrix& elemC, const GmMatrix& Lpu, const GmMatrix& Sc, const QList<int>& iUe, const QList<int>& iPe, const int& nPf, const int& numDof, const int& np);
198  virtual bool assemblyElemDampingMatrixNIT(const GmElement* e, GmMatrix& elemC, const GmMatrix& Lpu, const GmMatrix& Sc, const GmMatrix& Rt, const GmMatrix& Ct, const GmMatrix& Qtu, const QList<int>& iUe, const QList<int>& iPe, const QList<int>& iTe, const int& nPf, const int& numDof, const int& np);
199  virtual bool assemblyElemDampingMatrixMP(const GmElement* e, GmMatrix& elemC, const GmMatrix& Lpu, const GmMatrix& Cww, const GmMatrix& Cwc, const GmMatrix& Ccc, const GmMatrix& Ccw, const QList<int>& iUe, const QList<int>& iPw, const QList<int>& iPc, const int& nPf, const int& numDof, const int& np, const double& sw);
200  virtual bool assemblyCrackDampingMatrix(const GmElement* e, GmMatrix& elemC, const GmMatrix& Qf, const QList<int>& iUe, const QList<int>& iPe, const QList<int>& iPf, const QList<int>& pfList, const int& nPf, const int& numDof, const int& np);
201  virtual bool assemblyCrackDampingMatrixNIT(const GmElement* e, GmMatrix& elemC, const GmMatrix& Qf, const QList<int>& iUe, const QList<int>& iPe, const QList<int>& iPf, const QList<int>& iTf, const QList<int>& pfList, const int& nPf, const int& numDof, const int& np, const GmMatrix& Cwdwd, const GmMatrix& Rtdtd, const GmMatrix& Cwdtd);
202  virtual bool assemblyCrackDampingMatrixMP(const GmElement* e, GmMatrix& elemC, const GmMatrix& Qf, const QList<int>& iUe, const QList<int>& iPe, const QList<int>& iPf, const QList<int>& pfList, const int& nPf, const int& numDof, const int& np, const QList<int>& iPc, const QList<int>& iPcd, const double& sw,
203  const GmMatrix& Cwdwd, const GmMatrix& Cwdcd, const GmMatrix& Ccdcd, const GmMatrix& Ccdwd);
204  virtual bool assemblyElemStiffnessMatrix(const GmElement* e, GmMatrix& elemK, const GmMatrix& Kep, const GmMatrix& Lup, const GmMatrix& Kc, const QList<int>& iUe, const QList<int>& iPe, const int& nPf, const int& numDof, const int& np);
205  virtual bool assemblyElemStiffnessMatrixNIT(const GmElement* e, GmMatrix& elemK, const GmMatrix& Kep, const GmMatrix& Lup, const GmMatrix& Kc, const GmMatrix& Qut, const GmMatrix& Qtw, const GmMatrix& Mt, const GmMatrix& Ot, const QList<int>& iUe, const QList<int>& iPe, const QList<int>& iTe, const int& nPf, const int& numDof, const int& np);
206  virtual bool assemblyElemStiffnessMatrixMP(const GmElement* e, GmMatrix& elemK, const GmMatrix& Kep, const GmMatrix& Lup, const GmMatrix& Hww, const GmMatrix& Hcc, const QList<int>& iUe, const QList<int>& iPw, const QList<int>& iPc, const int& nPf, const int& numDof, const int& np, const double& sw);
207  virtual bool assemblyCrackStiffnessMatrix(const GmElement* e, GmMatrix& elemK, const GmMatrix& Qupf, const GmMatrix& L1, const GmMatrix& L2, const GmMatrix& L3, const GmMatrix& Hf, const GmMatrix& Kfr,
208  const QList<int>& iUe, const QList<int>& iPe, const QList<int>& iPf, const QList<int>& pfList, const int& nPf, const int& numDof, const int& np);
209  virtual bool assemblyCrackStiffnessMatrixNIT(const GmElement* e, GmMatrix& elemK, const GmMatrix& Qupf, const GmMatrix& L1, const GmMatrix& L2, const GmMatrix& L3, const GmMatrix& Hf, const GmMatrix& Kfr,
210  const QList<int>& iUe, const QList<int>& iPe, const QList<int>& iTe, const QList<int>& iPf, const QList<int>& iTd, const QList<int>& pfList, const int& nPf, const int& numDof, const int& np, const GmMatrix& Mtdtd, const GmMatrix&Otdtd, const GmMatrix&Qtdwd, const GmMatrix& Ltt, const GmMatrix& Lttd, const GmMatrix& Ltdtd, const GmMatrix& Qutd);
211  virtual bool assemblyCrackStiffnessMatrixMP(const GmElement* e, GmMatrix& elemK, const GmMatrix& Qupf, const GmMatrix& L1, const GmMatrix& L2, const GmMatrix& L3, const GmMatrix& Hf, const GmMatrix& Kfr, const GmMatrix& Hcdcd,
212  const QList<int>& iUe, const QList<int>& iPe, const QList<int>& iPf, const QList<int>& pfList, const int& nPf, const int& numDof, const int& np, const QList<int>& iPc, const QList<int>& iPcd, const double& swd);
213  virtual bool getPorePressureGaussPoint(const GmElement* e, const GmVector* coord, const int& ip, GmVector& Pg, const GmVector& Np, const GmVector& pe);
214  virtual bool getTemperatureGaussPoint(const GmElement* e, const GmVector* coord, const int& ip, double& Pg, const GmVector& Np, const GmVector& pe);
215  virtual bool getSaturationAndRelativePermeability(const GmElement* e, const GmVector* coord, const int& ip, GmVector& Pg, double& Np, double& pe, double& pa);
216  virtual void assemblyInternalForceVector(const GmElement* e, GmMatrix& elemFi, const GmVector& Fint_St, const GmVector& Fint_Pp, const QList<int>& iUe, const QList<int>& iPe, const int& nPf, const int& numDof, const int& np);
217  virtual void assemblyInternalForceVectorNIT(const GmElement* e, GmMatrix& elemFi, const GmVector& Fint_St, const GmVector& Fint_Pp, const GmVector& Fint_T, const QList<int>& iUe, const QList<int>& iPe, const QList<int>& iTe, const int& nPf, const int& numDof, const int& np);
218  virtual void assemblyInternalForceVectorMP(const GmElement* e, GmMatrix& elemFi, const GmVector& Fint_St, const GmVector& Fint_Pw, const GmVector& Fint_Pc, const QList<int>& iUe, const QList<int>& iPw, const QList<int>& iPc, const int& nPf, const int& numDof, const int& np);
219  virtual void assemblyCrackInternalForceVector(const GmElement* e, GmMatrix& elemFi, const GmVector& Fint_St, const GmVector& Fint_Pp, const GmVector& Fint_Hyd,
220  const QList<int>& iUe, const QList<int>& iPe, const QList<int>& iPf, const QList<int>& pfList, const int& nPf, const int& numDof, const int& np);
221  virtual void assemblyCrackInternalForceVectorNIT(const GmElement* e, GmMatrix& elemFi, const GmVector& Fint_St, const GmVector& Fint_Pp, const GmVector& Fint_Te, const GmVector& Fint_pf, const GmVector& Fint_Tf,
222  const QList<int>& iUe, const QList<int>& iPe, const QList<int>& iTe, const QList<int>& iPf, const QList<int>& iTf, const QList<int>& pfList, const int& nPf, const int& numDof, const int& np);
223  virtual void assemblyCrackInternalForceVectorMP(const GmElement* e, GmMatrix& elemFi, const GmVector& Fint_St, const GmVector& Fint_Pw, const GmVector& Fint_Pc, const QList<int>& iUe, const QList<int>& iPe, const QList<int>& iPf, const QList<int>& pfList, const int& nPf,
224  const int& numDof, const int& np, const QList<int>& iPc, const QList<int>& iPcd, const GmVector& Fint_Pwd, const GmVector& Fint_Pcd);
225  //Return a list with the pair of ghost nodes of each fracture segment contained in this element
226  virtual bool getCrackNodePair(XfemEnrichedElementData* exData, QList<QList<int>>& pfNodes);
227  // remove repeted values
228  virtual bool insertionSortWORepeatedint(QList<int>& a, const double& tol);
229  // Sort vector values
230  virtual bool insertionSortint(QList<int>& tVec);
231  //Return the fracture rotation matrix
232  void fillFractureRotationMatrix(const GmElement* e, const GmShape* shp, GmMatrix Je, GmMatrix &R);
233  // Returns the elastic constitutive matrix for the fracture segment
234  virtual bool fillNaturalCrackStiffnessMatrix(const GmElement* e, GmMatrix& Dep, GmVector coord, const int& ip, const GmVector& S_tn, const GmVector& Etn_old, const GmVector& v_w_tn, const GmMatrix& ROT, GmMatrix& D_xy, bool open);
235  // Returns the Cohesive constitutive matrix for the fracture segment
236  virtual bool fillCohesiveCrackStiffnessMatrix(const GmElement* e, GmMatrix& Dep, GmVector coord, const int& ip, const GmVector& nrd, const GmVector& trd, const GmVector& Etn_old, const GmVector& v_w_tn, const GmVector& sold,
237  const GmMatrix& ROT, GmMatrix& D_xy, const double& rockStrength, bool open);
238  // Returns fracture segment stresses
239  virtual bool naturalFractureReturnMapping(const GmElement* e, XfemEnrichedElementData* exData, const GmVector& sold, const GmMatrix& Dep, const GmVector& enew, const GmVector& eold, GmVector& Stn_new, const GmMatrix& ROT, GmVector& S_w_xy, bool open, const int& iFracSubElement, const GmVector& ip, const int& k);
240  // Returns fracture cohesive segment stresses
241  virtual bool cohesiveFractureReturnMapping(const GmElement* e, const GmVector& sold, const GmMatrix& Dep, const GmVector& enew, const GmVector& eold, GmVector& Stn_new, const GmMatrix& ROT, GmVector& S_w_xy, const double& rockStrength, bool open, const GmVector& ip, const int& k);
242  // Split tensor to interface or Join interface to tensor
243  virtual bool splitJoinTensorToInterface(GmVector& Tensor, GmVector& Cohesive, bool type);
244 
245 };
246 
247 #endif
virtual const GmElementDof * dofMapping(const GmElement *e) const
Definition: xfemHMCoupled.cpp:168
The number of state var ids above.
Definition: xfemHMCoupled.h:80
Id for fixed node fracture flow boundary condition of the wetting fluid.
Definition: xfemHMCoupled.h:102
Nodal pore flow for fixed node pore flow boundary condition of the wetting fluid.
Definition: xfemHMCoupled.h:120
virtual void fillCrackThermalExpansionMatrix(const GmElement *e, const GmVector *coord, const int &ip, GmMatrix &Sc, const GmVector &Np, const double &c, double &wn)
Calculates Thermal Expansion for element 'e' at the given integration point (coord,...
Definition: xfemHMCoupled.cpp:7343
Id for fixed node flow boundary condition of the non-wetting fluid.
Definition: xfemHMCoupled.h:99
Id for convective boundary condition.
Definition: xfemHMCoupled.h:109
Heat flux direction for fixed heat flux boundary condition.
Definition: xfemHMCoupled.h:133
virtual void fillStorageCompMatrix(const GmElement *e, const GmVector *coord, const int &ip, GmMatrix &Sc, const GmVector &Np, const double &c)
Calculates the compressibility matrix of "storage / compressibility" for element 'e' at the given int...
Definition: xfemHMCoupled.cpp:7093
Nodal temperature for fixed enriched temperature boundary condition.
Definition: xfemHMCoupled.h:130
virtual bool fillElementDataNIT(const GmElement *e, GmpFemMatrixSet &elemMatrices, GmpFemVectorSet &elemVectors)
Fills the element matrices and vectors for a generic elements.
Definition: xfemHMCoupled.cpp:2565
Id for fixed node pore pressure boundary condition of the wetting fluid.
Definition: xfemHMCoupled.h:95
virtual void fillEnhancedDof(const XfemElement *xe, QList< int > &iUe, QList< int > &iPe, QList< int > &iPf, const int &np)
Fills the position index with the degreee of freedom, both the mechanical and hydraulic physics....
Definition: xfemHMCoupled.cpp:9748
virtual void fillCouplingCrackHydMecMatrix(const GmElement *e, const GmVector *coord, const int &ip, GmMatrix &Epfu, const GmVector &Np, const GmMatrix &Bu, const GmVector &v_t, const double &c2, double &wn, GmMatrix &Qpfu, const GmVector &v_n, const GmMatrix &NuJump)
Calculates the stiffness matrix of the mechanical- hydraulic coupling.
Definition: xfemHMCoupled.cpp:7515
void fillElementPorePressure(const GmElement *e, GmVector &pe, GmVector &pf, bool ghostNodes, const QList< int > &pfList)
Given an element, fills the vector pe with nodal pore pressure. The vector should have size equal to ...
Definition: xfemHMCoupled.cpp:5287
Nodal fracture pressure for fixed node fracture pressure boundary condition of the non-wetting fluid.
Definition: xfemHMCoupled.h:127
Basic class for the standard plane strain physics plugin object.
Definition: xfemMechanic.h:68
Id for retrieving the accessor to the pore pressure state var of the wetting fluid.
Definition: xfemHMCoupled.h:69
Id for retrieving the accessor to the enriched pressure state var of the wetting fluid.
Definition: xfemHMCoupled.h:70
virtual void fillEnhancedDofMP(const XfemElement *xe, QList< int > &iUe, QList< int > &iPw, QList< int > &iPc, QList< int > &iPwd, QList< int > &iPcd, const int &np)
Fills the position index with the degreee of freedom, both the mechanical and hydraulic physics....
Definition: xfemHMCoupled.cpp:9817
virtual void fillStandardDofNIT(const GmElement *e, QList< int > &iUe, QList< int > &iPw, QList< int > &iPc)
Fills the position index with the degreee of freedo, both the mechanical and hydraulic physics....
Definition: xfemHMCoupled.cpp:9701
void fillEnhancedNpJumpMatrix(const GmElement *e, GmMatrix &Np, const GmVector &N, const QList< int > &BottomLevelSet, const QList< int > &TopLevelSet, const int &nEnrich, GmElementDof *dofMap, XfemEnrichedElementData *exData, const int &iSubCrack)
Calculates the enhanced strain displacement matrix (B) and the Jacobian determinant in a specified (x...
Definition: xfemHMCoupled.cpp:6683
virtual bool checkAndLoadPrivateData(LuaTable &table)
Overloads default checkAndLoadPrivateData() to be able to setup materials.
Definition: xfemHMCoupled.cpp:417
void unitVectSegment(const GmMatrix &Xg, GmMatrix &uVec, GmMatrix &uVec2, GmMatrix &uVec3, GmMatrix &ROT)
Definition: xfemHMCoupled.cpp:9999
Nodal enriched pressure for fixed node enriched pressure boundary condition of the non-wetting fluid.
Definition: xfemHMCoupled.h:124
Declaration of the xfemMechanic classes.
virtual void fillCouplingCrackThermoMecMatrix(const GmElement *e, const GmVector *coord, const int &ip, GmMatrix &Lup, const GmVector &Np, const GmMatrix &Bu, const double &c, const double &Tg, const double &wn)
Calculates the stiffness matrix of the Thermo - mechanical coupling.
Definition: xfemHMCoupled.cpp:7714
Nodal temperature for fixed temperature boundary condition.
Definition: xfemHMCoupled.h:129
virtual void fillHeatCapacityMatrix(const GmElement *e, const GmVector *coord, const int &ip, GmMatrix &Sc, const GmVector &Np, const double &c)
Calculates Heat Capacity for element 'e' at the given integration point (coord, ip) storing the resul...
Definition: xfemHMCoupled.cpp:7213
Nodal node fracture flow for fixed node fracture flow boundary condition of the wetting fluid.
Definition: xfemHMCoupled.h:126
The number of boundary conditions value ids above.
Definition: xfemMechanic.h:149
The number of boundary conditions value ids above.
Definition: xfemHMCoupled.h:140
Basic interface for an XFEM element.
Definition: xfemElement.h:43
Id for fixed node fracture pressure boundary condition of the non-wetting fluid.
Definition: xfemHMCoupled.h:103
Heat transfer coefficient for convective boundary condition.
Definition: xfemHMCoupled.h:134
Id for retrieving the accessor to the Enriched Temperature state var.
Definition: xfemHMCoupled.h:76
xfemHMStateVarIds
IDs for Mechanical xfem state vars.
Definition: xfemHMCoupled.h:67
void fillEnhancedBuJumpMatrix(const GmElement *e, const GmShape *shape, const GmVector &ncoord, const GmMatrix &J, GmMatrix &Bu, const QList< int > &TopLevelSet, const QList< int > &BottomLevelSet, const int &nEnrich, const GmVector &nodes, GmElementDof *dofMap, XfemEnrichedElementData *exData, const int &iSubCrack)
Calculates the enhanced strain displacement matrix (B) and the Jacobian determinant in a specified (x...
Definition: xfemHMCoupled.cpp:5869
Id for fixed enriched temperature boundary condition.
Definition: xfemHMCoupled.h:106
Id for fixed node pore pressure boundary condition of the non-wetting fluid.
Definition: xfemHMCoupled.h:98
void fillEnhancedNptopMatrix(const GmElement *e, GmVector &Np, const GmVector &N, const QList< int > &gaussLevelSet, const int &nEnrich, GmElementDof *dofMap, XfemEnrichedElementData *exData, const int &iSubCrack)
Calculates the enhanced strain displacement matrix (B) and the Jacobian determinant in a specified (x...
Definition: xfemHMCoupled.cpp:6294
void fillEnhancedBpJumpMatrix(const GmElement *e, const GmShape *shape, const GmVector &ncoord, const GmMatrix &J, GmMatrix &Bp, const QList< int > &TopLevelSet, const QList< int > &BottomLevelSet, const int &nEnrich, GmElementDof *dofMap, XfemEnrichedElementData *exData, const int &iSubCrack)
Calculates the enhanced strain displacement matrix (B) and the Jacobian determinant in a specified (x...
Definition: xfemHMCoupled.cpp:6040
Point generation rate for point generation boundary condition.
Definition: xfemHMCoupled.h:136
Id for fixed node enriched pressure boundary condition of the wetting fluid.
Definition: xfemHMCoupled.h:97
Id for fixed node flow boundary condition of the wetting fluid.
Definition: xfemHMCoupled.h:96
Id for retrieving the accessor to the Temperature state var.
Definition: xfemHMCoupled.h:75
The number of state var ids above.
Definition: xfemMechanic.h:116
Point position for point generation boundary condition.
Definition: xfemHMCoupled.h:137
Nodal enriched pressure for fixed node enriched pressure boundary condition of the wetting fluid.
Definition: xfemHMCoupled.h:121
virtual void fillEnhancedDofNIT(const XfemElement *xe, QList< int > &iUe, QList< int > &iPw, QList< int > &iPc, QList< int > &iPwd, QList< int > &iPcd, const int &np)
Fills the position index with the degreee of freedom, both the mechanical and hydraulic physics....
Definition: xfemHMCoupled.cpp:9899
virtual void fillCouplingHydMecMatrix(const GmElement *e, const GmVector *coord, const int &ip, GmMatrix &Lup, const GmVector &Np, const GmMatrix &Bu, const double &c)
Calculates the stiffness matrix of the hydraulic - mechanical coupling.
Definition: xfemHMCoupled.cpp:7630
virtual bool fixedNodalDofsBc(QList< int > &nodes, QList< int > &dof, QList< double > &values, bool *constantValues) const
See comments on base class. Fills vectors with prescribed node pore pressures.
Definition: xfemHMCoupled.cpp:454
virtual void fillCrackStorageCompMatrix(const GmElement *e, const GmVector *coord, const int &ip, GmMatrix &Sc, const GmVector &Np, const double &c)
Calculates the compressibility matrix of "storage / compressibility" for element 'e' at the given int...
Definition: xfemHMCoupled.cpp:7160
Basic class for the HM coupled Xfem physics plugin object.
Definition: xfemHMCoupled.h:36
xfemHMBoundaryConditionIds
IDs for accepted boundary condition types.
Definition: xfemHMCoupled.h:93
Nodal node fracture flow for fixed node fracture flow boundary condition of the non-wetting fluid.
Definition: xfemHMCoupled.h:128
void fillElementCapillaryPressure(const GmElement *e, GmVector &pe, GmVector &pf, bool ghostNodes, const QList< int > &pfList)
Given an element, fills the vector pe with nodal capillary pressure. The vector should have size equa...
Definition: xfemHMCoupled.cpp:5506
virtual void fillCouplingThermoMecMatrix(const GmElement *e, const GmVector *coord, const int &ip, GmMatrix &Lup, const GmVector &Np, const GmMatrix &Bu, const double &c, const double &Tg)
Calculates the stiffness matrix of the Thermo - mechanical coupling.
Definition: xfemHMCoupled.cpp:7667
xfemGaussAttributeIds
IDs for physics Gauss attributes.
Definition: xfemHMCoupled.h:84
virtual void fillStorageCompMatrixMP(const GmElement *e, const GmVector *coord, const int &ip, GmMatrix &Cww, GmMatrix &Cwc, GmMatrix &Ccc, GmMatrix &Ccw, const GmVector &Np, const double &c, const double &sw)
Calculates the compressibility matrix of "storage / compressibility" for element 'e' at the given int...
Definition: xfemHMCoupled.cpp:7395
virtual void fillCouplingCrackMecHydMatrix(const GmElement *e, const GmVector *coord, const int &ip, const double &c2, const GmVector &Np, GmMatrix &Qpfu, const GmVector &v_n, const GmMatrix &NuJump)
Calculates the stiffness matrix of the mechanical- hydraulic coupling.
Definition: xfemHMCoupled.cpp:7556
void fillEnhancedNuJumpMatrix(const GmElement *e, GmMatrix &Nu, const GmVector &N, const QList< int > &BottomLevelSet, const QList< int > &TopLevelSet, const int &nEnrich, const GmVector &nodes, GmElementDof *dofMap, XfemEnrichedElementData *exData, const int &iSubCrack)
Calculates the enhanced strain displacement matrix (B) and the Jacobian determinant in a specified (x...
Definition: xfemHMCoupled.cpp:6534
Id for retrieving the accessor to the enriched pressure state var of the non-wetting fluid.
Definition: xfemHMCoupled.h:72
virtual void fillCrackHeatCapacityMatrix(const GmElement *e, const GmVector *coord, const int &ip, GmMatrix &Sc, const GmVector &Np, const double &c, double &wn)
Calculates Heat Capacity for element 'e' at the given integration point (coord, ip) storing the resul...
Definition: xfemHMCoupled.cpp:7254
Id for fixed heat flux boundary condition.
Definition: xfemHMCoupled.h:108
Id for retrieving the accessor to the Fracture Temperature state var.
Definition: xfemHMCoupled.h:77
virtual FemResultType fillElementData(const GmElement *e, GmpFemMatrixSet &elemMatrices, GmpFemVectorSet &elemVectors)
Fills the element matrices and vectors for a generic elements.
Definition: xfemHMCoupled.cpp:849
virtual bool splitJoinTensorToInterface(GmVector &Tensor, GmVector &Cohesive, bool type)
Split tensor to interface or Join interface to tensor type = true ==> split tensor to interface compo...
Definition: xfemHMCoupled.cpp:10135
void fillEnhancedBpMatrix(const GmElement *e, const GmShape *shape, const GmVector &ncoord, const GmMatrix &X, const GmMatrix &J, const GmMatrix &B, GmMatrix &Bp, const QList< int > &gaussLevelSet, const int &nEnrich, GmElementDof *dofMap, XfemEnrichedElementData *exData, const int &gaussIP)
Calculates the enhanced strain displacement matrix (B) and the Jacobian determinant in a specified (x...
Definition: xfemHMCoupled.cpp:7774
void shapeValuesBarElem(double xi, GmVector &N)
Returns the shape function evaluated at Xi.
Definition: xfemHMCoupled.cpp:9987
void fillEnhancedNubotMatrix(const GmElement *e, GmMatrix &Nu, const GmVector &N, const QList< int > &gaussLevelSet, const int &nEnrich, const GmVector &nodes, GmElementDof *dofMap, XfemEnrichedElementData *exData, const int &iSubCrack)
Calculates the enhanced strain displacement matrix (B) and the Jacobian determinant in a specified (x...
Definition: xfemHMCoupled.cpp:6953
void fillEnhancedNutopMatrix(const GmElement *e, GmMatrix &Nu, const GmVector &N, const QList< int > &gaussLevelSet, const int &nEnrich, const GmVector &nodes, GmElementDof *dofMap, XfemEnrichedElementData *exData, const int &iSubCrack)
Calculates the enhanced strain displacement matrix (B) and the Jacobian determinant in a specified (x...
Definition: xfemHMCoupled.cpp:6807
xfemHMBoundaryConditionValueIds
IDs for property values from accepted boundary condition types.
Definition: xfemHMCoupled.h:117
virtual void fillCrackStorageCompMatrixMP(const GmElement *e, const GmVector *coord, const int &ip, GmMatrix &Cww, GmMatrix &Cwc, GmMatrix &Ccc, GmMatrix &Ccw, const GmVector &Np, const double &c, const double &sw, double &wn)
Calculates the compressibility matrix of "storage / compressibility" for element 'e' at the given int...
Definition: xfemHMCoupled.cpp:7459
virtual ~XfemHMCoupled()
Destructor.
Definition: xfemHMCoupled.cpp:85
void fillEnhancedNpbotMatrix(const GmElement *e, GmVector &Np, const GmVector &N, const QList< int > &gaussLevelSet, const int &nEnrich, GmElementDof *dofMap, XfemEnrichedElementData *exData, const int &iSubCrack)
Calculates the enhanced strain displacement matrix (B) and the Jacobian determinant in a specified (x...
Definition: xfemHMCoupled.cpp:6414
void fillEnhancedNpMatrix(const GmElement *e, GmVector &Np, const GmVector &N, const QList< int > &gaussLevelSet, const int &nEnrich, GmElementDof *dofMap, XfemEnrichedElementData *exData, const int &gaussIP)
Calculates the enhanced strain displacement matrix (B) and the Jacobian determinant in a specified (x...
Definition: xfemHMCoupled.cpp:6176
Heat flux for fixed heat flux boundary condition.
Definition: xfemHMCoupled.h:132
Nodal temperature for fixed fracture temperature boundary condition.
Definition: xfemHMCoupled.h:131
Ambient temperature for convective boundary condition.
Definition: xfemHMCoupled.h:135
Id for retrieving the accessor to the pore pressure state var of the non-wetting fluid.
Definition: xfemHMCoupled.h:71
void fillEnhancedBuMatrix(const GmElement *e, const GmShape *shape, const GmVector &ncoord, const GmMatrix &X, const GmMatrix &J, GmMatrix &Bu, const QList< int > &gaussLevelSet, const int &nEnrich, const GmVector &nodes, GmElementDof *dofMap, XfemEnrichedElementData *exData, const int &gaussIP)
Calculates the enhanced strain displacement matrix (B) and the Jacobian determinant in a specified (x...
Definition: xfemHMCoupled.cpp:5697
The number of Gauss attribute ids above.
Definition: xfemHMCoupled.h:89
virtual const QVariantMap * physicsMetaDataMap()
Returns a reference for the HydroMechanical coupled Xfem physics attribute map, built when the functi...
Definition: xfemHMCoupled.cpp:91
Nodal pore flow for fixed node pore flow boundary condition of the non-wetting fluid.
Definition: xfemHMCoupled.h:123
Nodal pore pressure for fixed node pore pressure boundary condition of the non-wetting fluid.
Definition: xfemHMCoupled.h:122
The number of boundary conditions ids above.
Definition: xfemMechanic.h:138
Id for fixed temperature boundary condition.
Definition: xfemHMCoupled.h:105
virtual void fillCouplingMecHydMatrix(const GmElement *e, const GmVector *coord, const int &ip, GmMatrix &Lup, const GmVector &Np, const GmMatrix &Bu, const double &c)
Calculates the stiffness matrix of the mechanical- hydraulic coupling.
Definition: xfemHMCoupled.cpp:7592
Id for retrieving the accessor to the fracture pressure state var of the non-wetting fluid.
Definition: xfemHMCoupled.h:74
virtual bool checkAndLoadAttributeAccessors(LuaTable &nodeTable, LuaTable &gaussTable)
Reimplementation of the common function to init the _hmaterialPointAccessor structure.
Definition: xfemHMCoupled.cpp:442
Id for fixed node fracture flow boundary condition of the non-wetting fluid.
Definition: xfemHMCoupled.h:104
The number of boundary conditions ids above.
Definition: xfemHMCoupled.h:113
Id for point generation boundary condition.
Definition: xfemHMCoupled.h:110
arma::vec GmVector
virtual bool fillElementDataMPNIT(const GmElement *e, GmpFemMatrixSet &elemMatrices, GmpFemVectorSet &elemVectors)
Fills the element matrices and vectors for a generic elements.
Definition: xfemHMCoupled.cpp:1726
Id for fixed fracture temperature boundary condition.
Definition: xfemHMCoupled.h:107
A class used to store per element data needed only in enriched elements.
Definition: xfemEnrichedElementData.h:37
virtual bool fixedNodalForcesBc(QList< int > &nodes, QList< int > &dof, QList< double > &values) const
See comments on base class. Fills vectors with prescribed nodal forces - flow.
Definition: xfemHMCoupled.cpp:514
virtual void fillThermalExpansionMatrix(const GmElement *e, const GmVector *coord, const int &ip, GmMatrix &Sc, const GmVector &Np, const double &c)
Calculates Thermal Expansion for element 'e' at the given integration point (coord,...
Definition: xfemHMCoupled.cpp:7295
Id for fixed node fracture pressure boundary condition of the wetting fluid.
Definition: xfemHMCoupled.h:101
Nodal pore pressure for fixed node pore pressure boundary condition of the wetting fluid.
Definition: xfemHMCoupled.h:119
void fillElementTemperature(const GmElement *e, GmVector &pe, GmVector &pf, bool ghostNodes, const QList< int > &pfList)
Given an element, fills the vector Te with nodal temperatures. The vector should have size equal to t...
Definition: xfemHMCoupled.cpp:5395
virtual void fillStandardDofMP(const GmElement *e, QList< int > &iUe, QList< int > &iPw, QList< int > &iPc)
Fills the position index with the degreee of freedo, both the mechanical and hydraulic physics....
Definition: xfemHMCoupled.cpp:9654
arma::mat GmMatrix
virtual void fillStandardDof(const GmElement *e, QList< int > &iUe, QList< int > &iPe)
Fills the position index with the degreee of freedo, both the mechanical and hydraulic physics....
Definition: xfemHMCoupled.cpp:9611
Id for retrieving the accessor to the fracture pressure state var of the wetting fluid.
Definition: xfemHMCoupled.h:73
XfemHMCoupled(const char *pluginType, GmSimulationData *simulation, QString id, QString description, const GmpFemPhysicsCommonMaterialFactory *matFactory, const GmLogCategory &logger)
Constructor. Will be called by the plugin loading code.
Definition: xfemHMCoupled.cpp:59
Nodal fracture pressure for fixed node fracture pressure boundary condition of the wetting fluid.
Definition: xfemHMCoupled.h:125
void fillElementDisplacements(const GmElement *e, GmVector &ue, const int &nEnrich)
See comments on base class. Fills the vector ue with nodal displacements.
Definition: xfemHMCoupled.cpp:5179
Id for fixed node enriched pressure boundary condition of the non-wetting fluid.
Definition: xfemHMCoupled.h:100
virtual bool fillElementDataMP(const GmElement *e, GmpFemMatrixSet &elemMatrices, GmpFemVectorSet &elemVectors)
Fills the element matrices and vectors for a generic elements.
Definition: xfemHMCoupled.cpp:3570
Base Id for Gauss attribute(s) used to store the wetting fluid saturation.
Definition: xfemHMCoupled.h:86
void tipDOFTreat(const int &nelem, GmVector &nodes)
treat tip degrees of freedom
Definition: xfemHMCoupled.cpp:5572