xfemTools.h

Go to the documentation of this file.

/************************************************************************
**
** Copyright (C) 2017 by Geomechanic group - Tecgraf/Puc-Rio
** All rights reserved.
**
** This file is part of the "GeMA" software. It's use should respect
** the terms in the license agreement that can be found together
** with this source code.
** It is provided AS IS, with NO WARRANTY OF ANY KIND,
** INCLUDING THE WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR
** A PARTICULAR PURPOSE.
**
************************************************************************/
 
#ifndef _GEMA_XFEM_TOOLS_H_
#define _GEMA_XFEM_TOOLS_H_
 
#include <xfemEnrichedElementData.h>
#include <gmpFemPhysics.h>
#include "xfemOptions.h"
#include "xfemMesh.h"
#include <assert.h>
 
#ifndef _USE_MATH_DEFINES
#define _USE_MATH_DEFINES  // To enable M_PI in Visual Studio
#endif
#include <math.h>
 
class XfemTools 
{
public:
        //XfemTools(XfemMesh* xmesh, XfemOptions& xfemOptions, GmSimulationData* simulation, const QList<GmpFemPhysics*>& physics, GmLogCategory& logger);
        XfemTools(XfemMesh* xmesh, XfemOptions& xfemOptions, GmSimulationData* simulation, QVector<GmpFemPhysics*>& physics, const GmLogCategory& logger);
 
  virtual ~XfemTools();
 
  int pointDimension() const { return _ndim; }
  //Compute Unit Vector
  virtual bool unitVectSegment(GmMatrix Xg, GmMatrix& uVec) const;
  //Compute Unit Vector of a plane
  virtual bool unitVectPlane(const GmMatrix & Xg, GmMatrix & uVec, GmMatrix & uVec2, GmMatrix & uVec3, GmMatrix & ROT) const;
  //Compute intersection of two planes
  virtual bool intersect3D_2Planes(const GmMatrix& Xg, GmMatrix& uVec, const GmMatrix& uVec2, GmMatrix& uVec3, GmMatrix & ROT, int& inter, const double& characteristicLength) const;
  //Compute intersection of a line with a plane
  virtual bool intersect3D_SegmentPlane(const GmMatrix & P0, const GmMatrix & P1, const GmMatrix & Pn2, GmMatrix & n2, GmMatrix & Line, int & inter) const;
  //Compute intersection of two planes
  virtual bool testIntersect3D_2Planes(const GmMatrix & Pn1, GmMatrix & n1, const GmMatrix & Pn2, GmMatrix & n2, GmMatrix & Line, int & inter, XfemMesh* xmesh, XfemElement* xe, XfemElement* xe2) const;
  //  returns a vector containing the i
  virtual bool intersectVolume(GmMatrix& nodes, QList<GmVector>& element, QList<GmVector>& fracData, GmMatrix& iNodes, QList<GmVector>& iData) const;
  //Create Rectangular Region at the Crack Tip
  virtual bool RectangularRegion(double Lch, GmMatrix Xg, GmMatrix uVec, GmMatrix& polygonVert) const;
  //Checks point relative position with convex polygon
  virtual bool pointInConvexPolygon(double tol, GmMatrix polygonVert, GmVector Xip, int& nIn) const;
  virtual bool polygonCentroid(const GmMatrix& polygonVert, GmVector& Xc) const;
  //Defines the integration points inside region to check propagation criterium
  virtual bool defIPinRegion(XfemMesh* xmesh, GmMatrix polygonVert, QList<GmVector>& IPInfo, GmVector crackInfo) const;
  //check and divide polygon if is cut by a segment
  virtual bool CutPolygonWithPolyline(XfemMesh* xmesh, GmVector crackInfo, GmMatrix& polygonVert, int idCrack) const;
  //Defines the integration points inside semicircle to check propagation criterion
  virtual bool defIPinCircle(XfemMesh* xmesh, GmVector crackInfo, QList<GmVector>& IPInfo) const;
  //Compute the cycle parameter
  virtual bool cycInc(int iVert, int nVert, int& nCycInc) const;
  //Compute the parametric position of line C (Seg1(xc1,yc1), Seg2(xc2,yc2)) on intersection with segment S(Node1(xs1, ys1), Node2(xs2, ys2))
  virtual bool parametricDist(GmVector S1, GmVector S2, GmVector N1, GmVector N2, double tol, double& tValue, bool& bIntersects) const;
  //Checks the intersection between Polyline A and B and cuts B
  virtual bool CutPLineWithPLine(XfemMesh* xmesh, double tol, GmVector Seg1node1, GmVector Seg1node2, GmVector Seg2node1, GmVector Seg2node2, 
          GmVector& NewCrackTip, bool& bCut, bool& bCutHeadVertexA, bool& bCutTailVertexA, bool& bCutHeadVertexB, bool& bCutTailVertexB) const;
  //Returns cartesian coordinates of parametric position (tValue) in segment(x1, y1), (x2, y2)
  virtual bool CoordsFromParam(GmVector Seg2node1, GmVector Seg2node2, double tValue, GmVector& NewCrackTip) const;
  //Principal stresses determination
  virtual bool principalStresses2D(GmVector stress, GmVector& Ps) const;
  //principal Stress Direction determination
  virtual bool principalStressDirection2D(GmVector aStress, GmVector aPStress, int nPS, GmVector& uPs) const;
  //For each region with propagation, computes new fracture paths
  virtual bool compPropagExtension(XfemMesh* xmesh, XfemFracture* fracture, int nside, GmMatrix Xnodes, GmVector Xg, GmVector aPropDirection, QList<GmVector>& aCrack, bool isHead, GmMatrix crackInfo) const;
  //For each region with propagation, computes new fracture paths
  virtual bool testCompPropagExtension3D(XfemMesh* xmesh, XfemElement* xe, QList<int> elemList, GmVector v_n, GmMatrix planeCoord, GmMatrix Xg, QList<GmVector>& aCrack, int& elemId, bool& next) const;
  //For each region with propagation, computes new fracture paths
  virtual bool compPropagExtension3D(XfemMesh* xmesh, XfemElement* xe, QList<int> elemList, GmVector v_n, GmMatrix& planeCoord, GmMatrix Xg, QList<GmVector>& aCrack, int& elemId, bool& next) const;
  //get caracteristic length
  virtual bool caracteristiclength(GmMatrix Xnodes, double& Lch) const;
  //set zero verifying tolerance
  virtual bool setZeroVeriTol(GmMatrix & Xnodes, double & Lch) const;
  //Checks the intersection between segment C {(xc1, yc1), (xc2, yc2)} and a regular xfemMesh
  virtual bool parametricSegmentInMesh(XfemMesh* xmesh, XfemFracture* fracture, double tol) const;
  //Define crack tips in tail and head values
  virtual bool checkAndInsertTips(QList<double>& tVec) const;
  // remove repeted values
  virtual bool insertionSortWORepeated(QList<double>& a, double tol) const;
  // Sort vector values
  virtual bool insertionSort(QList<double>& tVec) const;
  // remove repeted values
  virtual bool insertionSortWORepeatedint(QList<int>& a, double tol) const;
  // Sort vector values
  virtual bool insertionSortint(QList<int>& tVec) const;
  //Computes area of triangle defined by segment a,b and point c,  The signal of the area indicates the position of the point relative to the segment
  virtual bool orient2D(double xa, double ya, double xb, double yb, double xc, double yc, double& det) const ;
  // define elements crossed by a crack
  virtual bool defElemCuts(XfemMesh* xmesh) const;
  // define elemCutInfo based on the number of fractures
  virtual bool defElemSets(XfemMesh* xmesh) const;
  // set integration points for each subelement
  virtual bool setIntPoints(XfemMesh* xmesh, GmCell* e, QList<GmMatrix> subElem) const;
  // create sub elements with TetGen Library
  virtual bool createSubElem(XfemMesh* xmesh, GmCell* e, QList<GmMatrix>& subElem) const;
  // set data for integration points for each subelement
  virtual bool setEnrichedData(XfemMesh* xmesh, GmCell* e) const;
  // set ghost nodes into xfemMesh
  virtual bool setGhostNodes(XfemMesh* xmesh) const;
  // set ghost nodes into xfemMesh
  virtual bool setGhostNodes3D(XfemMesh* xmesh) const;
  // loads subelements 
  virtual bool loadSubElement(XfemMesh* xmesh, GmCell* e, GmMatrix subElemNode, int type, int rule, bool isFracture) const;
  //Set integration rule
  void setIntegrationRule(int rule) { _ruleSet = rule; };
  // Check the position of the crack points into the element
  virtual bool checkSegmentInRegularMesh(XfemMesh* xmesh, GmMatrix segCrackCoord, GmVector& mElemCutPos,double tol) const;
  //Compute the intersection between point crack (xc, yc) and segment S = node1(xs1, ys1), Node2(xs2, ys2)
  virtual bool intPointLine(GmVector N1, GmVector N2, GmVector crack, double tol, int& nInter) const;
  //Compute the intersection between point crack (xc, yc) and segment S = node1(xs1, ys1), Node2(xs2, ys2)
  virtual bool intSegment(GmVector N1, GmVector N2, GmVector crack, double tol, int& nInter) const;
  //define if a point is inside a sphere
  virtual bool inSphere(GmMatrix crack, GmVector Ip, double characteristicLength, double tol, int& nInter) const;
  // Divides polygon X by a segment S into two different polygons
  virtual bool dividePolygonWithSegment(GmCell* e, GmMatrix S, QList<GmMatrix>& subElemData, double tol) const;
  // Divides several Polygons by several segments for cases of multiple fractures in a element
  virtual bool dividePolygonWithSeveralSegments(XfemMesh* xmesh, GmCell* e, QList<GmMatrix>& subElemData, double tol) const;
  // Add enriched DOF due to cracks interaction
  virtual bool addEnrichedDof(XfemMesh* xmesh, int numDof, const QString type) const;
  // set Crack Front Data 
  virtual bool setCrackFrontData(XfemMesh * xmesh) const;
  // Define DOF sets for intersected fractures
  virtual bool defDofSets(XfemMesh* xmesh)const;
  //Define DOF sets for Mechanic physics
  virtual bool mechanicDofSets(XfemMesh* xmesh, int old, int& numDof, QList<GmVector> aCrack)const;
  //Define DOF sets for HydroMechanic physics
  virtual bool hydroMecDofSets(XfemMesh* xmesh, int old, int& numDof, QList<GmVector> aCrack)const;
  // divide and add sub elements
  virtual bool subElemDiv(XfemMesh* xmesh) const;
  // divide and add sub elements
  virtual bool subElemDiv3D(XfemMesh* xmesh) const;
  //test with voro++ library
  void test() const;
  //test with TetGen library
  void test2() const;
  //Return the fracture rotation matrix
  void fillFractureRotationMatrix(GmCell* e, const GmShape * shp, GmMatrix J, GmMatrix & R) const;
  //void setAcc(GmValueAccessor* hnAcc) { _HNodeAcc = hnAcc; };
  virtual bool pointOnSide(GmCell* e, GmMatrix xGhost, GmVector& xIp, int side) const;
  //Indicates the relative location of the secondary crack respect to the main one
  virtual bool fillInterCrackNodalLevelSet(GmVector Seg1node1, GmVector Seg1node2, GmVector Seg2node1, GmVector& interactingCraks) const;
  //Compute nodal levelset function (H)
  void fillNodalLevelSetMatrix(GmMatrix X, GmMatrix Xg, GmMatrix& H)const;
  //Compute Gauss levelset function (H)
  void fillGaussLevelSetMatrix(GmMatrix X, GmMatrix Xg, GmMatrix& H)const;
  void fillInteractingCraksNodalLevelSets(GmCell* e, XfemEnrichedElementData* exData, int nEnrich, int nCracks, QList<double>crackIndex, QList<GmMatrix>& nodalLevelSet)const;
  void fillInteractingCraksGaussLevelSets(GmCell* e, XfemEnrichedElementData* exData, int nEnrich, int nCracks, QList<double>crackIndex, QList<GmMatrix>& gaussLevelSet)const;
  //Mix list with enrichment functions
  //virtual bool mixList(XfemEnrichedElementData* exDatac, QList<GmMatrix>& nodalLevelSet, QList<double> crackIndex, double value, bool add, bool test, int type, int iSubCrack) const;
  //Add DOFs to adjacent elements to an intersection
  //virtual bool correctDOF(XfemMesh* xmesh, QList<GmVector>& aCrack, int& old) const;
  virtual bool correctDOF(XfemMesh* xmesh, GmCell* e, XfemElement* xe, XfemEnrichedElementData* exData, QList<XfemFracture*> fractureData, QList<GmVector>& aCrack, int& old, bool& add) const;
  //Delete old specific node DOF to be modified by the following subroutines
  virtual bool deleteNodeDOF(XfemMesh* xmesh, double value) const;
  //check DOFs to adjacent elements to an intersection
  virtual bool checkDOF(XfemMesh* xmesh, double value, QList<int>& enrichList) const;
  //Add adjacent elements to an intersection to the list to be updated for the assembler
  virtual bool addChangedElem(XfemMesh* xmesh, double value, QList<int>& enrichList) const;
  //fill Enrichment Functions
  virtual bool saveEnrichmentFunctions(XfemMesh* xmesh, double value) const;
  //alter Enrichment Functions
  virtual bool alterEnrichmentFunctions(XfemMesh* xmesh, double value) const;
  // Create global DOF list
  virtual bool createDofList(XfemMesh* xmesh, const QString type) const;
  // check DOFs and Enrichment Functions of adjacent elements to an intersection 
  virtual bool checkData(XfemMesh* xmesh) const;
  //Transform stresses in global XY coordinate system to local crack plane coordinate system
  virtual bool getCrackPlaneStresses(GmVector Sxy, GmMatrix ROT, GmVector& Stn) const;
  //Transform Strains in global XY coordinate system to local crack plane coordinate system
  virtual bool getCrackPlaneStrains(GmVector Exy, GmMatrix ROT, GmVector& Etn) const;
  // Reorganize from vector to tensorial representation and vice versa
  virtual bool vectorialToTensorial(GmVector& Av, GmMatrix& Am, unsigned sc, QString type = "strain", bool mode = true) const;
 
protected:
  GmValueAccessor*       _nodeAcc;           
  XfemGaussAccessor*     _stressAcc;         
  XfemGaussAccessor*     _strainAcc;         
  GmCellAccessor*        _strengthAcc;       
  int                    _ruleSet;           
  XfemOptions            _xfemOptions;       
  GmSimulationData*      _SimulationData;    
  QVector<GmpFemPhysics*>  _physics;           // physics data
 
  //GmValueAccessor*       _HNodeAcc;
 
//private:
  int     _ndim;                    
  const GmLogCategory& _logger;     
                                                       
};
 
#endif