www.mooseframework.org
GeometricCut3DUserObject.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
3 /* */
4 /* All contents are licensed under LGPL V2.1 */
5 /* See LICENSE for full restrictions */
6 /****************************************************************/
7 
9 
10 // MOOSE includes
11 #include "MooseError.h"
12 
13 #include "libmesh/string_to_enum.h"
14 
15 // XFEM includes
16 #include "XFEMFuncs.h"
17 
18 template <>
19 InputParameters
21 {
22  // Get input parameters from parent class
23  InputParameters params = validParams<GeometricCutUserObject>();
24 
25  // Class description
26  params.addClassDescription("Base class for 3D XFEM Geometric Cut UserObjects");
27  // Return the parameters
28  return params;
29 }
30 
31 GeometricCut3DUserObject::GeometricCut3DUserObject(const InputParameters & parameters)
32  : GeometricCutUserObject(parameters), _center(), _normal()
33 {
34  _cut_time_ranges.push_back(std::make_pair(0.0, 0.0));
35 }
36 
37 bool GeometricCut3DUserObject::active(Real /*time*/) const { return true; }
38 
39 bool
41  std::vector<CutEdge> & /*cut_edges*/,
42  std::vector<CutNode> & /*cut_nodes*/,
43  Real /*time*/) const
44 {
45  mooseError("Invalid method: must use vector of element faces for 3D mesh cutting");
46  return false;
47 }
48 
49 bool
51  std::vector<CutFace> & cut_faces,
52  Real /*time*/) const
53 // TODO: Time evolving cuts not yet supported in 3D (hence the lack of use of the time variable)
54 {
55  bool cut_elem = false;
56 
57  for (unsigned int i = 0; i < elem->n_sides(); ++i)
58  {
59  // This returns the lowest-order type of side.
60  std::unique_ptr<Elem> curr_side = elem->side(i);
61  if (curr_side->dim() != 2)
62  mooseError("In cutElementByGeometry dimension of side must be 2, but it is ",
63  curr_side->dim());
64  unsigned int n_edges = curr_side->n_sides();
65 
66  std::vector<unsigned int> cut_edges;
67  std::vector<Real> cut_pos;
68 
69  for (unsigned int j = 0; j < n_edges; j++)
70  {
71  // This returns the lowest-order type of side.
72  std::unique_ptr<Elem> curr_edge = curr_side->side(j);
73  if (curr_edge->type() != EDGE2)
74  mooseError("In cutElementByGeometry face edge must be EDGE2, but type is: ",
75  libMesh::Utility::enum_to_string(curr_edge->type()),
76  " base element type is: ",
77  libMesh::Utility::enum_to_string(elem->type()));
78  Node * node1 = curr_edge->get_node(0);
79  Node * node2 = curr_edge->get_node(1);
80 
81  Point intersection;
82  if (intersectWithEdge(*node1, *node2, intersection))
83  {
84  cut_edges.push_back(j);
85  cut_pos.push_back(getRelativePosition(*node1, *node2, intersection));
86  }
87  }
88 
89  if (cut_edges.size() == 2)
90  {
91  cut_elem = true;
92  CutFace mycut;
93  mycut.face_id = i;
94  mycut.face_edge.push_back(cut_edges[0]);
95  mycut.face_edge.push_back(cut_edges[1]);
96  mycut.position.push_back(cut_pos[0]);
97  mycut.position.push_back(cut_pos[1]);
98  cut_faces.push_back(mycut);
99  }
100  }
101 
102  return cut_elem;
103 }
104 
105 bool
106 GeometricCut3DUserObject::cutFragmentByGeometry(std::vector<std::vector<Point>> & /*frag_edges*/,
107  std::vector<CutEdge> & /*cut_edges*/,
108  Real /*time*/) const
109 {
110  mooseError("Invalid method: must use vector of element faces for 3D mesh cutting");
111  return false;
112 }
113 
114 bool
115 GeometricCut3DUserObject::cutFragmentByGeometry(std::vector<std::vector<Point>> & /*frag_faces*/,
116  std::vector<CutFace> & /*cut_faces*/,
117  Real /*time*/) const
118 {
119  // TODO: Need this for branching in 3D
120  mooseError("cutFragmentByGeometry not yet implemented for 3D mesh cutting");
121  return false;
122 }
123 
124 bool
125 GeometricCut3DUserObject::intersectWithEdge(const Point & p1, const Point & p2, Point & pint) const
126 {
127  bool has_intersection = false;
128  double plane_point[3] = {_center(0), _center(1), _center(2)};
129  double plane_normal[3] = {_normal(0), _normal(1), _normal(2)};
130  double edge_point1[3] = {p1(0), p1(1), p1(2)};
131  double edge_point2[3] = {p2(0), p2(1), p2(2)};
132  double cut_point[3] = {0.0, 0.0, 0.0};
133 
135  plane_point, plane_normal, edge_point1, edge_point2, cut_point) == 1)
136  {
137  Point temp_p(cut_point[0], cut_point[1], cut_point[2]);
138  if (isInsideCutPlane(temp_p) && isInsideEdge(p1, p2, temp_p))
139  {
140  pint = temp_p;
141  has_intersection = true;
142  }
143  }
144  return has_intersection;
145 }
146 
147 bool
148 GeometricCut3DUserObject::isInsideEdge(const Point & p1, const Point & p2, const Point & p) const
149 {
150  Real dotp1 = (p1 - p) * (p2 - p1);
151  Real dotp2 = (p2 - p) * (p2 - p1);
152  return (dotp1 * dotp2 <= 0.0);
153 }
154 
155 Real
157  const Point & p2,
158  const Point & p) const
159 {
160  // get the relative position of p from p1
161  Real full_len = (p2 - p1).norm();
162  Real len_p1_p = (p - p1).norm();
163  return len_p1_p / full_len;
164 }
GeometricCut3DUserObject(const InputParameters &parameters)
InputParameters validParams< GeometricCutUserObject >()
std::vector< std::pair< Real, Real > > _cut_time_ranges
int plane_normal_line_exp_int_3d(double pp[3], double normal[3], double p1[3], double p2[3], double pint[3])
Definition: XFEMFuncs.C:408
InputParameters validParams< GeometricCut3DUserObject >()
virtual bool intersectWithEdge(const Point &p1, const Point &p2, Point &pint) const
Real getRelativePosition(const Point &p1, const Point &p2, const Point &p) const
virtual bool isInsideCutPlane(Point p) const =0
std::vector< Real > position
virtual bool active(Real time) const override
virtual bool cutFragmentByGeometry(std::vector< std::vector< Point >> &frag_edges, std::vector< CutEdge > &cut_edges, Real time) const override
unsigned int face_id
std::vector< unsigned int > face_edge
virtual bool cutElementByGeometry(const Elem *elem, std::vector< CutEdge > &cut_edges, std::vector< CutNode > &cut_nodes, Real time) const override
bool isInsideEdge(const Point &p1, const Point &p2, const Point &p) const