10 #if !defined(KRATOS_TRIGEN_PFEM_MODELER_H_INCLUDED )
11 #define KRATOS_TRIGEN_PFEM_MODELER_H_INCLUDED
36 void triangulate(
char *,
struct triangulateio *,
struct triangulateio *,
struct triangulateio *);
110 Element const& rReferenceElement,
111 Condition const& rReferenceBoundaryCondition,
113 double my_alpha = 1.4,
double h_factor=0.5)
117 if (ThisModelPart.
NodesBegin()->SolutionStepsDataHas(IS_FREE_SURFACE)==
false )
118 KRATOS_THROW_ERROR(std::logic_error,
"Add ----IS_FREE_SURFACE---- variable!!!!!! ERROR",
"");
119 if (ThisModelPart.
NodesBegin()->SolutionStepsDataHas(IS_STRUCTURE)==
false )
120 KRATOS_THROW_ERROR(std::logic_error,
"Add ----IS_STRUCTURE---- variable!!!!!! ERROR",
"");
121 if (ThisModelPart.
NodesBegin()->SolutionStepsDataHas(IS_BOUNDARY)==
false )
122 KRATOS_THROW_ERROR(std::logic_error,
"Add ----IS_BOUNDARY---- variable!!!!!! ERROR",
"");
123 if (ThisModelPart.
NodesBegin()->SolutionStepsDataHas(IS_FLUID)==
false )
127 const auto inital_time = std::chrono::steady_clock::now();
146 unsigned int bucket_size = 20;
149 unsigned int max_results = 100;
154 Node work_point(0,0.0,0.0,0.0);
160 list_of_nodes.reserve(ThisModelPart.
Nodes().size());
161 for(ModelPart::NodesContainerType::iterator i_node = ThisModelPart.
NodesBegin() ; i_node != ThisModelPart.
NodesEnd() ; i_node++)
163 (list_of_nodes).push_back(*(i_node.base()));
166 KdtreeType nodes_tree1(list_of_nodes.begin(),list_of_nodes.end(), bucket_size);
168 RemoveCloseNodes(ThisModelPart, nodes_tree1, node_erase, h_factor);
176 struct triangulateio in_mid;
177 struct triangulateio out_mid;
178 struct triangulateio vorout_mid;
180 initialize_triangulateio(in_mid);
181 initialize_triangulateio(out_mid);
182 initialize_triangulateio(vorout_mid);
186 in_mid.numberofpoints = ThisModelPart.
Nodes().size();
187 in_mid.pointlist = (
REAL *) malloc(in_mid.numberofpoints * 2 *
sizeof(
REAL));
190 ModelPart::NodesContainerType::iterator nodes_begin = ThisModelPart.
NodesBegin();
191 for(
unsigned int i = 0;
i<ThisModelPart.
Nodes().size();
i++)
197 (nodes_begin +
i)->SetId(
i+1);
200 in_mid.pointlist[base] = (nodes_begin +
i)->
X();
201 in_mid.pointlist[base+1] = (nodes_begin +
i)->
Y();
203 auto& node_dofs = (nodes_begin +
i)->GetDofs();
204 for(
auto iii = node_dofs.begin(); iii != node_dofs.end(); iii++)
206 (**iii).SetEquationId(
i+1);
217 char options1[] =
"Pne";
218 triangulate(options1, &in_mid, &out_mid, &vorout_mid);
222 unsigned int el_number=out_mid.numberoftriangles;
225 std::vector<int> preserved_list1(el_number);
226 preserved_list1.resize(el_number);
228 int number_of_preserved_elems= PassAlphaShape(ThisModelPart, preserved_list1, el_number, my_alpha, out_mid);
232 clean_triangulateio(in_mid);
233 clean_triangulateio(vorout_mid);
237 struct triangulateio in2;
238 struct triangulateio out2;
239 struct triangulateio vorout2;
241 initialize_triangulateio(in2);
242 initialize_triangulateio(out2);
243 initialize_triangulateio(vorout2);
246 in2.numberofpoints = ThisModelPart.
Nodes().size();
247 in2.pointlist = (
REAL *) malloc(in2.numberofpoints * 2 *
sizeof(
REAL));
250 for(
unsigned int i = 0;
i<ThisModelPart.
Nodes().size();
i++)
253 in2.pointlist[base] = (nodes_begin +
i)->
X();
254 in2.pointlist[base+1] = (nodes_begin +
i)->
Y();
256 in2.numberoftriangles=number_of_preserved_elems;
258 in2.trianglelist = (
int *) malloc(in2.numberoftriangles * 3 *
sizeof(
int));
259 in2.trianglearealist = (
REAL *) malloc(in2.numberoftriangles *
sizeof(
REAL));
266 for(ModelPart::NodesContainerType::iterator i_node = ThisModelPart.
NodesBegin() ; i_node != ThisModelPart.
NodesEnd() ; i_node++)
268 if (i_node->FastGetSolutionStepValue(IS_FREE_SURFACE)!=0)
271 double& val=i_node->FastGetSolutionStepValue(NODAL_H);
279 for (
unsigned int el=0;
el<el_number;
el++)
281 if(
static_cast<bool>(preserved_list1[
el]) ==
true )
287 in2.trianglelist[new_base] = out_mid.trianglelist[old_base];
288 in2.trianglelist[new_base+1] = out_mid.trianglelist[old_base+1];
289 in2.trianglelist[new_base+2] = out_mid.trianglelist[old_base+2];
292 double prescribed_h = (nodes_begin + out_mid.trianglelist[old_base]-1)->FastGetSolutionStepValue(NODAL_H);
293 prescribed_h += (nodes_begin + out_mid.trianglelist[old_base+1]-1)->FastGetSolutionStepValue(NODAL_H);
294 prescribed_h += (nodes_begin + out_mid.trianglelist[old_base+2]-1)->FastGetSolutionStepValue(NODAL_H);
295 prescribed_h *= 0.3333;
297 in2.trianglearealist[
counter] = 0.5*(1.5*prescribed_h*1.5*prescribed_h);
303 for(ModelPart::NodesContainerType::iterator i_node = ThisModelPart.
NodesBegin() ; i_node != ThisModelPart.
NodesEnd() ; i_node++)
305 if (i_node->FastGetSolutionStepValue(IS_FREE_SURFACE)!=0)
307 double& nodal_h=i_node->FastGetSolutionStepValue(NODAL_H);
313 clean_triangulateio(out_mid);
322 char mesh_regen_opts[] =
"YJq1.4arn";
323 triangulate(mesh_regen_opts, &in2, &out2, &vorout2);
328 char mesh_regen_opts[] =
"YJrn";
329 triangulate(mesh_regen_opts, &in2, &out2, &vorout2);
337 typedef std::vector<PointType::Pointer>
PointVector;
357 int n_points_before_refinement = in2.numberofpoints;
359 if (out2.numberofpoints > n_points_before_refinement )
361 for(
int i = n_points_before_refinement;
i<out2.numberofpoints;
i++)
365 double&
x= out2.pointlist[base];
366 double&
y= out2.pointlist[base+1];
370 pnode->SetBufferSize(ThisModelPart.
NodesBegin()->GetBufferSize() );
372 list_of_new_nodes.push_back( pnode );
376 for(Node::DofsContainerType::iterator iii = reference_dofs.begin(); iii != reference_dofs.end(); iii++)
381 (p_new_dof)->FreeDof();
389 std::cout <<
"During refinement we added " << out2.numberofpoints-n_points_before_refinement<<
" nodes " <<std::endl;
404 unsigned int MaximumNumberOfResults = list_of_new_nodes.size();
410 if(out2.numberofpoints-n_points_before_refinement > 0)
413 kd_tree nodes_tree2(list_of_new_nodes.begin(),list_of_new_nodes.end(),bucket_size);
416 for(
int el = 0;
el< in2.numberoftriangles;
el++)
420 point_base = (in2.trianglelist[base] - 1)*2;
421 x1[0] = in2.pointlist[point_base];
422 x1[1] = in2.pointlist[point_base+1];
424 point_base = (in2.trianglelist[base+1] - 1)*2;
425 x2[0] = in2.pointlist[point_base];
426 x2[1] = in2.pointlist[point_base+1];
428 point_base = (in2.trianglelist[base+2] - 1)*2;
429 x3[0] = in2.pointlist[point_base];
430 x3[1] = in2.pointlist[point_base+1];
435 CalculateCenterAndSearchRadius(
x1[0],
x1[1],
441 int number_of_points_in_radius;
444 work_point.
Z() = 0.0;
446 number_of_points_in_radius = nodes_tree2.SearchInRadius(work_point,
radius*1.01, Results.begin(),
447 ResultsDistances.begin(), MaximumNumberOfResults);
450 *( (nodes_begin + in2.trianglelist[base]-1).base() ),
451 *( (nodes_begin + in2.trianglelist[base+1]-1).base() ),
452 *( (nodes_begin + in2.trianglelist[base+2]-1).base() )
456 for(
auto it_found = Results.begin(); it_found != Results.begin() + number_of_points_in_radius; it_found++)
458 bool is_inside =
false;
462 (*it_found)->X(),(*it_found)->Y(),
N);
465 if(is_inside ==
true)
477 for(
auto & list_of_new_node : list_of_new_nodes)
479 const array_1d<double,3>& disp = list_of_new_node->FastGetSolutionStepValue(DISPLACEMENT);
480 list_of_new_node->X0() = list_of_new_node->X() - disp[0];
481 list_of_new_node->Y0() = list_of_new_node->Y() - disp[1];
482 list_of_new_node->Z0() = 0.0;
492 (ThisModelPart.
Elements()).reserve(out2.numberoftriangles);
494 for(
int iii = 0; iii< out2.numberoftriangles; iii++)
499 *( (nodes_begin + out2.trianglelist[base]-1).base() ),
500 *( (nodes_begin + out2.trianglelist[base+1]-1).base() ),
501 *( (nodes_begin + out2.trianglelist[base+2]-1).base() )
507 if( *(ModelNodes).find( out2.trianglelist[base]).base() == *(ThisModelPart.
Nodes().end()).base() )
509 if( *(ModelNodes).find( out2.trianglelist[base+1]).base() == *(ThisModelPart.
Nodes().end()).base() )
511 if( *(ModelNodes).find( out2.trianglelist[base+2]).base() == *(ThisModelPart.
Nodes().end()).base() )
516 Element::Pointer p_element = rReferenceElement.
Create(
id, geom, properties);
517 (ThisModelPart.
Elements()).push_back(p_element);
523 ModelPart::ElementsContainerType::const_iterator el_begin = ThisModelPart.
ElementsBegin();
524 for(ModelPart::ElementsContainerType::const_iterator iii = ThisModelPart.
ElementsBegin();
528 int base = ( iii->Id() - 1 )*3;
530 (iii->GetValue(NEIGHBOUR_ELEMENTS)).resize(3);
532 for(
int i = 0;
i<3;
i++)
534 int index = out2.neighborlist[base+
i];
538 neighb(
i) = Element::WeakPointer();
542 IdentifyBoundary(ThisModelPart, rReferenceBoundaryCondition, properties, out2);
546 clean_triangulateio(in2);
548 clean_triangulateio(out2);
550 clean_triangulateio(vorout2);
551 KRATOS_WATCH(
"Finished remeshing with Trigen_PFEM_Refine")
572 virtual std::string
Info()
const
636 boost::numeric::ublas::bounded_matrix<double,2,2> mJ;
637 boost::numeric::ublas::bounded_matrix<double,2,2> mJinv;
651 unsigned int max_results = 100;
656 Node work_point(0,0.0,0.0,0.0);
658 unsigned int n_points_in_radius;
662 for(ModelPart::NodesContainerType::const_iterator in = ThisModelPart.
NodesBegin(); in != ThisModelPart.
NodesEnd(); in++)
664 radius=h_factor*in->FastGetSolutionStepValue(NODAL_H);
666 work_point[0]=in->X();
667 work_point[1]=in->Y();
668 work_point[2]=in->Z();
670 n_points_in_radius = nodes_tree1.SearchInRadius(work_point,
radius,
res.begin(),res_distances.begin(), max_results);
671 if (n_points_in_radius>1)
673 if (in->FastGetSolutionStepValue(IS_BOUNDARY)==0.0 && in->FastGetSolutionStepValue(IS_STRUCTURE)==0.0)
676 double erased_nodes = 0;
677 for(
auto i=
res.begin();
i!=
res.begin() + n_points_in_radius ;
i++)
678 erased_nodes += in->Is(TO_ERASE);
680 if( erased_nodes < 1)
681 in->Set(TO_ERASE,
true);
684 else if ( (in)->FastGetSolutionStepValue(IS_STRUCTURE)!=1.0)
690 for(
auto i=
res.begin();
i!=
res.begin() + n_points_in_radius ;
i++)
692 if ( (*i)->FastGetSolutionStepValue(IS_BOUNDARY,1)==1.0 && res_distances[
k] < 0.2*
radius && res_distances[
k] > 0.0 )
700 in->Set(TO_ERASE,
true);
714 int PassAlphaShape(
ModelPart& ThisModelPart,
std::vector<
int>& preserved_list1,
unsigned int & el_number,
double& my_alpha, struct triangulateio& out_mid)
722 array_1d<double,3>
x1,
x2,x3,
xc;
725 int number_of_preserved_elems=0;
728 for(
unsigned int el = 0;
el< el_number;
el++)
733 point_base = (out_mid.trianglelist[base] - 1)*2;
734 x1[0] = out_mid.pointlist[point_base];
735 x1[1] = out_mid.pointlist[point_base+1];
737 point_base = (out_mid.trianglelist[base+1] - 1)*2;
738 x2[0] = out_mid.pointlist[point_base];
739 x2[1] = out_mid.pointlist[point_base+1];
741 point_base = (out_mid.trianglelist[base+2] - 1)*2;
742 x3[0] = out_mid.pointlist[point_base];
743 x3[1] = out_mid.pointlist[point_base+1];
746 Geometry<Node >
temp;
748 temp.push_back( *((ThisModelPart.Nodes()).find( out_mid.trianglelist[base]).base() ) );
749 temp.push_back( *((ThisModelPart.Nodes()).find( out_mid.trianglelist[base+1]).base() ) );
750 temp.push_back( *((ThisModelPart.Nodes()).find( out_mid.trianglelist[base+2]).base() ) );
752 int number_of_structure_nodes =
int(
temp[0].FastGetSolutionStepValue(IS_STRUCTURE) );
753 number_of_structure_nodes +=
int(
temp[1].FastGetSolutionStepValue(IS_STRUCTURE) );
754 number_of_structure_nodes +=
int(
temp[2].FastGetSolutionStepValue(IS_STRUCTURE) );
757 int nfs =
int(
temp[0].FastGetSolutionStepValue(IS_FREE_SURFACE) );
758 nfs +=
int(
temp[1].FastGetSolutionStepValue(IS_FREE_SURFACE) );
759 nfs +=
int(
temp[2].FastGetSolutionStepValue(IS_FREE_SURFACE) );
762 int nfluid =
int(
temp[0].FastGetSolutionStepValue(IS_FLUID) );
763 nfluid +=
int(
temp[1].FastGetSolutionStepValue(IS_FLUID) );
764 nfluid +=
int(
temp[2].FastGetSolutionStepValue(IS_FLUID) );
776 if(number_of_structure_nodes!=3)
779 if (nfs != 0 || nfluid != 3)
782 if( AlphaShape(my_alpha,
temp) && number_of_structure_nodes!=3)
790 preserved_list1[
el] =
true;
791 number_of_preserved_elems += 1;
797 double bigger_alpha = my_alpha*10.0;
798 if( AlphaShape(bigger_alpha,
temp) && number_of_structure_nodes!=3)
806 preserved_list1[
el] =
true;
807 number_of_preserved_elems += 1;
813 preserved_list1[
el] =
false;
816 return number_of_preserved_elems;
823 void IdentifyBoundary(
ModelPart& ThisModelPart, Condition
const& rReferenceBoundaryCondition, Properties::Pointer& properties,
struct triangulateio& out2)
827 for(ModelPart::NodesContainerType::const_iterator in = ThisModelPart.NodesBegin(); in!=ThisModelPart.NodesEnd(); in++)
829 in->FastGetSolutionStepValue(IS_BOUNDARY) = 0;
832 ModelPart::ElementsContainerType::iterator elements_end = ThisModelPart.Elements().end();
834 ThisModelPart.Elements().Unique();
837 for(ModelPart::ElementsContainerType::iterator iii = ThisModelPart.ElementsBegin(); iii != ThisModelPart.ElementsEnd(); iii++)
839 int base = ( iii->Id() - 1 )*3;
841 ModelPart::ElementsContainerType::iterator el_neighb;
852 el_neighb = (ThisModelPart.Elements()).find( out2.neighborlist[base] );
853 if( el_neighb == elements_end )
857 iii->GetGeometry()[1].FastGetSolutionStepValue(IS_BOUNDARY) = 1;
858 iii->GetGeometry()[2].FastGetSolutionStepValue(IS_BOUNDARY) = 1;
863 temp1.push_back(iii->GetGeometry()(1));
864 temp1.push_back(iii->GetGeometry()(2));
866 Geometry< Node >::Pointer cond = Geometry< Node >::Pointer(
new Geometry< Node >(temp1) );
867 int id = (iii->Id()-1)*3;
868 Condition::Pointer p_cond = rReferenceBoundaryCondition.Create(
id, temp1, properties);
869 ThisModelPart.Conditions().push_back(p_cond);
875 el_neighb = (ThisModelPart.Elements()).find( out2.neighborlist[base+1] );
877 if( el_neighb == elements_end )
880 iii->GetGeometry()[2].FastGetSolutionStepValue(IS_BOUNDARY) = 1;
881 iii->GetGeometry()[0].FastGetSolutionStepValue(IS_BOUNDARY) = 1;
886 temp1.push_back(iii->GetGeometry()(2));
887 temp1.push_back(iii->GetGeometry()(0));
889 Geometry< Node >::Pointer cond = Geometry< Node >::Pointer(
new Geometry< Node >(temp1) );
890 int id = (iii->Id()-1)*3+1;
891 Condition::Pointer p_cond = rReferenceBoundaryCondition.Create(
id, temp1, properties);
892 ThisModelPart.Conditions().push_back(p_cond);
899 el_neighb = (ThisModelPart.Elements()).find( out2.neighborlist[base+2] );
900 if( el_neighb == elements_end )
903 iii->GetGeometry()[0].FastGetSolutionStepValue(IS_BOUNDARY) = 1;
904 iii->GetGeometry()[1].FastGetSolutionStepValue(IS_BOUNDARY) = 1;
909 temp1.push_back(iii->GetGeometry()(0));
910 temp1.push_back(iii->GetGeometry()(1));
912 Geometry< Node >::Pointer cond = Geometry< Node >::Pointer(
new Geometry< Node >(temp1) );
913 int id = (iii->Id()-1)*3+2;
914 Condition::Pointer p_cond = rReferenceBoundaryCondition.Create(
id, temp1, properties);
915 ThisModelPart.Conditions().push_back(p_cond);
925 bool AlphaShape(
double alpha_param, Geometry<Node >& pgeom)
931 double x0 = pgeom[0].X();
932 double x1 = pgeom[1].X();
933 double x2 = pgeom[2].X();
935 double y0 = pgeom[0].Y();
936 double y1 = pgeom[1].Y();
937 double y2 = pgeom[2].Y();
945 double detJ = mJ(0,0)*mJ(1,1)-mJ(0,1)*mJ(1,0);
947 mJinv(0,0) = mJ(1,1);
948 mJinv(0,1) = -mJ(0,1);
949 mJinv(1,0) = -mJ(1,0);
950 mJinv(1,1) = mJ(0,0);
952 bounded_matrix<double,2,2> check;
956 h = pgeom[0].FastGetSolutionStepValue(NODAL_H);
957 h += pgeom[1].FastGetSolutionStepValue(NODAL_H);
958 h += pgeom[2].FastGetSolutionStepValue(NODAL_H);
966 pgeom[0].GetSolutionStepValue(IS_BOUNDARY) = 1;
967 pgeom[1].GetSolutionStepValue(IS_BOUNDARY) = 1;
968 pgeom[2].GetSolutionStepValue(IS_BOUNDARY) = 1;
975 double x0_2 = x0*x0 + y0*y0;
976 double x1_2 =
x1*
x1 + y1*y1;
977 double x2_2 =
x2*
x2 + y2*y2;
983 mRhs[0] = (x1_2 - x0_2);
984 mRhs[1] = (x2_2 - x0_2);
989 double radius = sqrt(pow(mC[0]-x0,2)+pow(mC[1]-y0,2));
1006 inline void CalculateCenterAndSearchRadius(
const double x0,
const double y0,
1007 const double x1,
const double y1,
1008 const double x2,
const double y2,
1009 double&
xc,
double&
yc,
double&
R
1012 xc = 0.3333333333333333333*(x0+
x1+
x2);
1013 yc = 0.3333333333333333333*(y0+y1+y2);
1015 double R1 = (
xc-x0)*(
xc-x0) + (
yc-y0)*(
yc-y0);
1026 inline double CalculateVol(
const double x0,
const double y0,
1027 const double x1,
const double y1,
1028 const double x2,
const double y2
1031 return 0.5*( (
x1-x0)*(y2-y0)- (y1-y0)*(
x2-x0) );
1035 const double x1,
const double y1,
1036 const double x2,
const double y2,
1037 const double xc,
const double yc,
1038 array_1d<double,3>&
N
1043 if(area < 0.000000000001)
1058 double upper_limit = 1.0+
tol;
1059 double lower_limit = -
tol;
1061 if(
N[0] >= lower_limit &&
N[1] >= lower_limit &&
N[2] >= lower_limit &&
N[0] <= upper_limit &&
N[1] <= upper_limit &&
N[2] <= upper_limit)
1068 void Interpolate( Triangle2D3<Node >& geom,
const array_1d<double,3>&
N,
1069 unsigned int step_data_size,
1070 Node::Pointer pnode)
1072 unsigned int buffer_size = pnode->GetBufferSize();
1080 double* step_data = (pnode)->SolutionStepData().Data(
step);
1083 double* node0_data = geom[0].SolutionStepData().Data(
step);
1084 double* node1_data = geom[1].SolutionStepData().Data(
step);
1085 double* node2_data = geom[2].SolutionStepData().Data(
step);
1088 for(
unsigned int j= 0;
j<step_data_size;
j++)
1091 step_data[
j] =
N[0]*node0_data[
j] +
N[1]*node1_data[
j] +
N[2]*node2_data[
j];
1096 if (
N[0]==0.0 &&
N[1]==0.0 &&
N[2]==0.0)
1097 KRATOS_THROW_ERROR(std::logic_error,
"SOMETHING's wrong with the added nodes!!!!!! ERROR",
"");
1108 pnode->FastGetSolutionStepValue(IS_BOUNDARY)=0.0;
1109 pnode->FastGetSolutionStepValue(IS_STRUCTURE)=0.0;
1110 pnode->FastGetSolutionStepValue(IS_INTERFACE)=0.0;
1111 pnode->Set(TO_ERASE,
false);
1112 pnode->FastGetSolutionStepValue(IS_FREE_SURFACE)=0.0;
1113 pnode->FastGetSolutionStepValue(IS_FLUID)=1.0;
1116 void initialize_triangulateio( triangulateio& tr )
1118 tr.pointlist = (
REAL*)
nullptr;
1119 tr.pointattributelist = (
REAL*)
nullptr;
1120 tr.pointmarkerlist = (
int*)
nullptr;
1121 tr.numberofpoints = 0;
1122 tr.numberofpointattributes = 0;
1123 tr.trianglelist = (
int*)
nullptr;
1124 tr.triangleattributelist = (
REAL*)
nullptr;
1125 tr.trianglearealist = (
REAL*)
nullptr;
1126 tr.neighborlist = (
int*)
nullptr;
1127 tr.numberoftriangles = 0;
1128 tr.numberofcorners = 3;
1129 tr.numberoftriangleattributes = 0;
1130 tr.segmentlist = (
int*)
nullptr;
1131 tr.segmentmarkerlist = (
int*)
nullptr;
1132 tr.numberofsegments = 0;
1133 tr.holelist = (
REAL*)
nullptr;
1134 tr.numberofholes = 0;
1135 tr.regionlist = (
REAL*)
nullptr;
1136 tr.numberofregions = 0;
1137 tr.edgelist = (
int*)
nullptr;
1138 tr.edgemarkerlist = (
int*)
nullptr;
1139 tr.normlist = (
REAL*)
nullptr;
1140 tr.numberofedges = 0;
1143 void clean_triangulateio( triangulateio& tr )
1145 if(tr.pointlist !=
nullptr) free(tr.pointlist );
1146 if(tr.pointattributelist !=
nullptr) free(tr.pointattributelist );
1147 if(tr.pointmarkerlist !=
nullptr) free(tr.pointmarkerlist );
1148 if(tr.trianglelist !=
nullptr) free(tr.trianglelist );
1149 if(tr.triangleattributelist !=
nullptr) free(tr.triangleattributelist );
1150 if(tr.trianglearealist !=
nullptr) free(tr.trianglearealist );
1151 if(tr.neighborlist !=
nullptr) free(tr.neighborlist );
1152 if(tr.segmentlist !=
nullptr) free(tr.segmentlist );
1153 if(tr.segmentmarkerlist !=
nullptr) free(tr.segmentmarkerlist );
1154 if(tr.holelist !=
nullptr) free(tr.holelist );
1155 if(tr.regionlist !=
nullptr) free(tr.regionlist );
1156 if(tr.edgelist !=
nullptr) free(tr.edgelist );
1157 if(tr.edgemarkerlist !=
nullptr) free(tr.edgemarkerlist );
1158 if(tr.normlist !=
nullptr) free(tr.normlist );
1199 inline std::istream&
operator >> (std::istream& rIStream,
1207 rOStream << std::endl;
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
Short class definition.
Definition: bucket.h:57
Base class for all Conditions.
Definition: condition.h:59
Geometry< NodeType >::PointsArrayType NodesArrayType
definition of nodes container type, redefined from GeometryType
Definition: condition.h:86
Dof represents a degree of freedom (DoF).
Definition: dof.h:86
Base class for all Elements.
Definition: element.h:60
virtual Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, PropertiesType::Pointer pProperties) const
It creates a new element pointer.
Definition: element.h:202
This process removes the entities from a model part with the flag TO_ERASE.
Definition: entity_erase_process.h:70
void Execute() override
Execute method is used to execute the Process algorithms.
This class is a wrapper for a pointer to a data that is located in a different rank.
Definition: global_pointer.h:44
Definition: amatrix_interface.h:497
Definition: amatrix_interface.h:530
PropertiesType::Pointer pGetProperties(IndexType PropertiesId)
Definition: mesh.h:394
This class aims to manage meshes for multi-physics simulations.
Definition: model_part.h:77
ElementIterator ElementsBegin(IndexType ThisIndex=0)
Definition: model_part.h:1169
NodeType::Pointer CreateNewNode(int Id, double x, double y, double z, VariablesList::Pointer pNewVariablesList, IndexType ThisIndex=0)
Definition: model_part.cpp:270
SizeType GetNodalSolutionStepDataSize()
Definition: model_part.h:571
NodeIterator NodesBegin(IndexType ThisIndex=0)
Definition: model_part.h:487
ConditionsContainerType & Conditions(IndexType ThisIndex=0)
Definition: model_part.h:1381
ElementsContainerType & Elements(IndexType ThisIndex=0)
Definition: model_part.h:1189
MeshType::NodesContainerType NodesContainerType
Nodes container. Which is a vector set of nodes with their Id's as key.
Definition: model_part.h:128
ElementIterator ElementsEnd(IndexType ThisIndex=0)
Definition: model_part.h:1179
NodeIterator NodesEnd(IndexType ThisIndex=0)
Definition: model_part.h:497
MeshType & GetMesh(IndexType ThisIndex=0)
Definition: model_part.h:1791
NodesContainerType & Nodes(IndexType ThisIndex=0)
Definition: model_part.h:507
This class defines the node.
Definition: node.h:65
std::vector< std::unique_ptr< Dof< double > >> DofsContainerType
The DoF container type definition.
Definition: node.h:92
double Y() const
Definition: point.h:187
double Z() const
Definition: point.h:193
double X() const
Definition: point.h:181
void reserve(size_type dim)
Definition: pointer_vector.h:319
static double ElapsedSeconds(const std::chrono::steady_clock::time_point StartTime)
This method returns the resulting time.
Definition: timer.h:179
A generic tree data structure for spatial partitioning.
Definition: tree.h:190
Short class definition.
Definition: trigen_pfem_refine.h:64
Node::Pointer PointPointerType
Definition: trigen_pfem_refine.h:72
virtual void PrintInfo(std::ostream &rOStream) const
Print information about this object.
Definition: trigen_pfem_refine.h:578
Node PointType
Definition: trigen_pfem_refine.h:71
virtual void PrintData(std::ostream &rOStream) const
Print object's data.
Definition: trigen_pfem_refine.h:581
Tree< KDTreePartition< BucketType > > KdtreeType
Definition: trigen_pfem_refine.h:78
PointVector::iterator PointIterator
Definition: trigen_pfem_refine.h:74
TriGenPFEMModeler()
Default constructor.
Definition: trigen_pfem_refine.h:85
Bucket< 3, PointType, PointVector, PointPointerType, PointIterator, DistanceIterator > BucketType
Definition: trigen_pfem_refine.h:77
std::vector< double > DistanceVector
Definition: trigen_pfem_refine.h:75
virtual std::string Info() const
Turn back information as a string.
Definition: trigen_pfem_refine.h:572
KRATOS_CLASS_POINTER_DEFINITION(TriGenPFEMModeler)
Pointer definition of TriGenModeler.
virtual ~TriGenPFEMModeler()=default
Destructor.
void ReGenerateMesh(ModelPart &ThisModelPart, Element const &rReferenceElement, Condition const &rReferenceBoundaryCondition, EntitiesEraseProcess< Node > &node_erase, bool rem_nodes=true, bool add_nodes=true, double my_alpha=1.4, double h_factor=0.5)
Definition: trigen_pfem_refine.h:108
std::vector< PointType::Pointer > PointVector
Definition: trigen_pfem_refine.h:73
std::vector< double >::iterator DistanceIterator
Definition: trigen_pfem_refine.h:76
A three node 2D triangle geometry with linear shape functions.
Definition: triangle_2d_3.h:74
#define KRATOS_THROW_ERROR(ExceptionType, ErrorMessage, MoreInfo)
Definition: define.h:77
#define KRATOS_CATCH(MoreInfo)
Definition: define.h:110
#define KRATOS_WATCH(variable)
Definition: define.h:806
#define KRATOS_TRY
Definition: define.h:109
#define REAL
Definition: delaunator_utilities.cpp:16
Kratos::ModelPart ModelPart
Definition: kratos_wrapper.h:31
z
Definition: GenerateWind.py:163
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
bool CalculatePosition(Geometry< Node > &geom, const double xc, const double yc, const double zc, array_1d< double, 3 > &N)
Definition: transfer_utility.h:343
double CalculateVol(const double x0, const double y0, const double x1, const double y1, const double x2, const double y2)
Definition: transfer_utility.h:424
AMatrix::MatrixProductExpression< TExpression1Type, TExpression2Type > prod(AMatrix::MatrixExpression< TExpression1Type, TCategory1 > const &First, AMatrix::MatrixExpression< TExpression2Type, TCategory2 > const &Second)
Definition: amatrix_interface.h:568
std::istream & operator>>(std::istream &rIStream, LinearMasterSlaveConstraint &rThis)
input stream function
std::vector< PointTypePointer > PointVector
Definition: internal_variables_interpolation_process.h:53
T & noalias(T &TheMatrix)
Definition: amatrix_interface.h:484
void triangulate(char *, struct triangulateio *, struct triangulateio *, struct triangulateio *)
REACTION_CHECK_STIFFNESS_FACTOR int
Definition: contact_structural_mechanics_application_variables.h:75
std::ostream & operator<<(std::ostream &rOStream, const LinearMasterSlaveConstraint &rThis)
output stream function
Definition: linear_master_slave_constraint.h:432
std::vector< double > DistanceVector
Definition: internal_variables_interpolation_process.h:55
int step
Definition: face_heat.py:88
y
Other simbols definition.
Definition: generate_axisymmetric_navier_stokes_element.py:54
res
Definition: generate_convection_diffusion_explicit_element.py:211
h
Definition: generate_droplet_dynamics.py:91
x2
Definition: generate_frictional_mortar_condition.py:122
x1
Definition: generate_frictional_mortar_condition.py:121
int tol
Definition: hinsberg_optimization.py:138
R
Definition: isotropic_damage_automatic_differentiation.py:172
float radius
Definition: mesh_to_mdpa_converter.py:18
def Interpolate(variable, entity, sf_values, historical_value)
Definition: point_output_process.py:231
int k
Definition: quadrature.py:595
int j
Definition: quadrature.py:648
el
Definition: read_stl.py:25
float temp
Definition: rotating_cone.py:85
float xc
Definition: rotating_cone.py:77
float yc
Definition: rotating_cone.py:78
int counter
Definition: script_THERMAL_CORRECT.py:218
add_nodes
Definition: script.py:96
N
Definition: sensitivityMatrix.py:29
x
Definition: sensitivityMatrix.py:49
namespace
Definition: array_1d.h:793
integer i
Definition: TensorModule.f:17
e
Definition: run_cpp_mpi_tests.py:31