NeoPZ
TPZGenSpecialGrid.cpp
Go to the documentation of this file.
1 
6 #include "pzmanvector.h"
7 #include "pzstack.h"
8 #include "pzgeoel.h"
9 
10 #include "TPZGenSpecialGrid.h"
11 
12 #include <math.h>
13 
18  // Initial mesh data
19  // In this case is a octahedron
20  // Octahedron has 6 nodes. We considerer that the nodes is on axes and radius is the distance between any node and origin.
21  const int64_t nnode = 6;
22  // Octahedron has 8 (triangular) faces
23  const int64_t nelem = 8;
24 
25  // Initial nodes and initial triangular faces of the octahedron
26  REAL initialcoord[nnode][3] = {{center[0]-radius,center[1],center[2]},{center[0],center[1]+radius,center[2]},{center[0],center[1],center[2]-radius},{center[0],center[1],center[2]+radius},{center[0],center[1]-radius,center[2]},{center[0]+radius,center[1],center[2]}};
27  int indices[nelem][nnode] = {{3,4,5},{3,5,1},{3,1,0},{3,0,4},{4,0,2},{4,2,5},{2,0,1},{5,2,1}};
28 
29  // Geometric element vector
30  TPZGeoEl *elvec[nelem], *gel;
31 
32  // Creating geometric initial mesh
33  TPZGeoMesh *gmesh = new TPZGeoMesh();
34 
35  // Initializing nodes of the polygonal initial mesh
36  int64_t node;
37  for(node=0; node<nnode; node++) {
38  int nodind = gmesh->NodeVec().AllocateNewElement();
39  TPZManVector<REAL> coord(3);
40  coord[0] = initialcoord[node][0];
41  coord[1] = initialcoord[node][1];
42  coord[2] = initialcoord[node][2];
43  gmesh->NodeVec()[nodind] = TPZGeoNode(node,coord,*gmesh);
44  }
45  // Creating triangular elements
46  int64_t el, index;
47  for(el=0; el<nelem; el++) {
49  for(node=0; node<3; node++) nodind[node]=indices[el][node];
50  elvec[el] = gmesh->CreateGeoElement(ETriangle,nodind,1,index);
51  }
52  // Building connectivity for initial mesh - OCTAHEDRON
53  gmesh->BuildConnectivity();
54 
55  //Loop making uniform refinement and changing coordinates of the nodes (projecting into the sphere) until tolerance is reached
56  TPZManVector<REAL> baryparam(3,0.), barycenter(3,0.);
57  for(int ii=0;ii<nUniformRefs;ii++) {
58  // Make a uniform refinement
59  UniformRefinement(1,gmesh,2);
60 
61  // Projecting all nodes into the sphere
62  TPZVec<REAL> coordinates(3,0.);
63  REAL norm;
64  int i;
65  for(node=0;node<gmesh->NNodes();node++) {
66  TPZGeoNode *gnode = &gmesh->NodeVec()[node];
67  gnode->GetCoordinates(coordinates);
68  for(i=0;i<3;i++)
69  coordinates[i] -= center[i];
70  norm = sqrt(coordinates[0]*coordinates[0] + coordinates[1]*coordinates[1] + coordinates[2]*coordinates[2]);
71  for(i=0;i<3;i++) coordinates[i] /= norm;
72  for(i=0;i<3;i++) coordinates[i] += center[i];
73  gnode->SetCoord(coordinates);
74  }
75  }
76 
77  // Before return, all the "father" elements are deleted, but the "son" elements are not
78  for(el=0;el<gmesh->NElements();el++) {
79  gel = gmesh->ElementVec()[el];
80  if(!gel) continue;
81  if(gel->HasSubElement()) {
82  gel->ResetSubElements();
83  gmesh->DeleteElement(gel);
84  }
85  else
86  gel->SetFatherIndex(-1);
87  }
88  // The connectivity are updated
89  gmesh->ResetConnectivities();
90  gmesh->BuildConnectivity();
91 
92  return gmesh;
93 }
94 
99  // Initial mesh data
100  // In this case is a octahedron
101  // Octahedron has 6 nodes. We considerer that the nodes is on axes and radius is the distance between any node and origin.
102  const int64_t nnode = 6;
103  // Octahedron has 8 (triangular) faces
104  const int64_t nelem = 8;
105 
106  // Initial nodes and initial triangular faces of the octahedron
107  REAL initialcoord[nnode][3] = {{center[0]-radius,center[1],center[2]},{center[0],center[1]+radius,center[2]},{center[0],center[1],center[2]-radius},{center[0],center[1],center[2]+radius},{center[0],center[1]-radius,center[2]},{center[0]+radius,center[1],center[2]}};
108  int indices[nelem][nnode] = {{3,4,5},{3,5,1},{3,1,0},{3,0,4},{4,0,2},{4,2,5},{2,0,1},{5,2,1}};
109 
110  // Geometric element vector
111  TPZGeoEl *elvec[nelem], *gel;
112 
113  // Creating geometric initial mesh
114  TPZGeoMesh *gmesh = new TPZGeoMesh();
115 
116  // Initializing nodes of the polygonal initial mesh
117  int64_t node;
118  for(node=0; node<nnode; node++) {
119  int64_t nodind = gmesh->NodeVec().AllocateNewElement();
120  TPZManVector<REAL> coord(3);
121  coord[0] = initialcoord[node][0];
122  coord[1] = initialcoord[node][1];
123  coord[2] = initialcoord[node][2];
124  gmesh->NodeVec()[nodind] = TPZGeoNode(node,coord,*gmesh);
125  }
126  // Creating triangular elements
127  int64_t el, index;
128  for(el=0; el<nelem; el++) {
130  for(node=0; node<3; node++) nodind[node]=indices[el][node];
131  elvec[el] = gmesh->CreateGeoElement(ETriangle,nodind,1,index);
132  }
133  // Building connectivity for initial mesh - OCTAHEDRON
134  gmesh->BuildConnectivity();
135 
136  //Loop making uniform refinement and changing coordinates of the nodes (projecting into the sphere) until tolerance is reached
137  TPZManVector<REAL> baryparam(3,0.), barycenter(3,0.);
138  REAL dist = 1.;
139  while(tol < dist) {
140  // Make a uniform refinement
141  UniformRefinement(1,gmesh,2);
142 
143  // Projecting all nodes into the sphere
144  TPZVec<REAL> coordinates(3,0.);
145  REAL norm;
146  int i;
147  for(node=0;node<gmesh->NNodes();node++) {
148  TPZGeoNode *gnode = &gmesh->NodeVec()[node];
149  gnode->GetCoordinates(coordinates);
150  for(i=0;i<3;i++)
151  coordinates[i] -= center[i];
152  norm = sqrt(coordinates[0]*coordinates[0] + coordinates[1]*coordinates[1] + coordinates[2]*coordinates[2]);
153  for(i=0;i<3;i++) coordinates[i] /= norm;
154  for(i=0;i<3;i++) coordinates[i] += center[i];
155  gnode->SetCoord(coordinates);
156  }
157 
158  // Find element no null to calculate distance between barycenter and projection into sphere
159  for(el=0;el<gmesh->NElements();el++) {
160  gel = gmesh->ElementVec()[el];
161  if(!gel || gel->HasSubElement()) continue;
162  gel->CenterPoint(6,baryparam);
163  gel->X(baryparam,barycenter);
164  break;
165  }
166  dist = radius - sqrt((barycenter[0]-center[0])*(barycenter[0]-center[0]) + (barycenter[1]-center[1])*(barycenter[1]-center[1]) + (barycenter[2]-center[2])*(barycenter[2]-center[2]));
167  if(dist < 0.0) {
168  delete gmesh;
169  gmesh = 0;
170  return gmesh;
171  }
172  }
173 
174  // Before return, all the "father" elements are deleted, but the "son" elements are not
175  for(el=0;el<gmesh->NElements();el++) {
176  gel = gmesh->ElementVec()[el];
177  if(!gel) continue;
178  if(gel->HasSubElement()) {
179  gel->ResetSubElements();
180  gmesh->DeleteElement(gel);
181  }
182  else
183  gel->SetFatherIndex(-1);
184  }
185  // The connectivity are updated
186  gmesh->ResetConnectivities();
187  gmesh->BuildConnectivity();
188 
189  return gmesh;
190 }
191 
195 void TPZGenSpecialGrid::UniformRefinement(int nUniformRefs,TPZGeoMesh *gmesh, const int dimelfordivision, bool allmaterial, const int matidtodivided) {
197  for(int Division=0; Division<nUniformRefs; Division++)
198  {
199  int64_t nels = gmesh->NElements();
200  for(int64_t elem = 0; elem < nels; elem++)
201  {
202  TPZGeoEl * gel = gmesh->ElementVec()[elem];
203  if(!gel || gel->HasSubElement())
204  continue;
205  if(dimelfordivision > 0 && gel->Dimension() != dimelfordivision) continue;
206  if(!allmaterial){
207  if(gel->MaterialId() == matidtodivided){
208  gel->Divide(filhos);
209  }
210  }
211  else{
212  gel->Divide(filhos);
213  }
214  }
215  }
216  gmesh->ResetConnectivities();
217  gmesh->BuildConnectivity();
218 }
int AllocateNewElement()
Makes more room for new elements.
Definition: pzadmchunk.h:184
void SetCoord(const TPZVec< REAL > &x)
Sets all coordinates into the current node. It gets the dim values from x.
Definition: pzgnode.cpp:60
int MaterialId() const
Returns the material index of the element.
Definition: pzgeoel.h:250
virtual TPZGeoEl * CreateGeoElement(MElementType type, TPZVec< int64_t > &cornerindexes, int matid, int64_t &index, int reftype=1)
Generic method for creating a geometric element. Putting this method centrally facilitates the modifi...
Definition: pzgmesh.cpp:1296
int64_t NElements() const
Number of elements of the mesh.
Definition: pzgmesh.h:129
virtual void Divide(TPZVec< TPZGeoEl *> &pv)
Divides the element and puts the resulting elements in the vector.
Definition: pzgeoel.h:742
virtual void CenterPoint(int side, TPZVec< REAL > &masscent) const =0
It returns the coordinates from the center of the side of the element in the element coordinate space...
static const double tol
Definition: pzgeoprism.cpp:23
Defines the behaviour of all geometric elements. GeometryTPZGeoEl is the common denominator for all g...
Definition: pzgeoel.h:43
virtual void SetFatherIndex(int64_t fatherindex)
Sets the father element index This method is not called SetFather in order to avoid implicit conversi...
Definition: pzgeoel.h:490
int64_t NNodes() const
Number of nodes of the mesh.
Definition: pzgmesh.h:126
Free store vector implementation.
TPZAdmChunkVector< TPZGeoNode > & NodeVec()
Definition: pzgmesh.h:140
expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ sqrt
Definition: tfadfunc.h:120
Declaration of the generator of special grids class.
static void UniformRefinement(const int nUniformRefs, TPZGeoMesh *gmesh, const int dimelfordivision, bool allmaterial=true, const int matidtodivided=0)
Make uniform refinement of the geometric mesh. This method finishs making the new connectivities for ...
A simple stack.
virtual int HasSubElement() const =0
Return 1 if the element has subelements.
void BuildConnectivity()
Build the connectivity of the grid.
Definition: pzgmesh.cpp:967
void DeleteElement(TPZGeoEl *gel, int64_t index=-1)
Centralized method to delete elements.
Definition: pzgmesh.cpp:1456
virtual int Dimension() const =0
Returns the dimension of the element.
REAL dist(TPZVec< T1 > &vec1, TPZVec< T1 > &vec2)
Definition: pzvec_extras.h:124
Implements a geometric node in the pz environment. Geometry.
Definition: pzgnode.h:31
virtual void X(TPZVec< REAL > &qsi, TPZVec< REAL > &result) const =0
Return the coordinate in real space of the point coordinate in the master element space...
This class implements a geometric mesh for the pz environment. Geometry.
Definition: pzgmesh.h:48
static TPZGeoMesh * GeneratePolygonalSphereFromOctahedron(TPZVec< REAL > &center, REAL radius, int nUniformRefs)
Static function to generate a polygonal mesh approximating a sphere from a octahedron mesh (polygonal...
void ResetConnectivities()
Reset all connectivities.
Definition: pzgmesh.cpp:1576
void GetCoordinates(TPZVec< REAL > &co)
Fill the coordinates of the node.
Definition: pzgnode.cpp:83
static int nodind[7][8]
virtual void ResetSubElements()=0
Reset all subelements to NULL.
TPZAdmChunkVector< TPZGeoEl * > & ElementVec()
Methods for handling pzlists.
Definition: pzgmesh.h:138