NeoPZ
pzgmesh.cpp
Go to the documentation of this file.
1 
6 #include "pzgmesh.h"
7 
8 #ifdef HAVE_CONFIG_H
9 #include <pz_config.h>
10 #endif
11 
12 #include "pzvec.h"
13 #include "pzvec_extras.h"
14 #include "pzadmchunk.h"
15 #include "pzcmesh.h"
16 #include "pzcompel.h"
17 #include "pzgnode.h"
18 #include "TPZMaterial.h"
19 #include "pzerror.h"
20 #include "pzgeoel.h"
21 #include "pzmatrix.h"
22 
23 #include "TPZRefPattern.h"
24 #include "TPZRefPatternDataBase.h"
25 #include "tpzgeoelrefpattern.h"
26 #include "tpzgeoblend.h"
27 
28 #ifdef BORLAND
29 #include <io.h>
30 #include <fcntl.h>
31 #endif
32 
33 #include <sstream>
34 #include <string>
35 
36 #include "TPZStream.h"
37 
38 #include "pzlog.h"
39 
40 #ifdef LOG4CXX
41 static LoggerPtr logger(Logger::getLogger("pz.mesh.tpzgeomesh"));
42 #endif
43 
44 
45 using namespace std;
46 
48 fName(), fElementVec(0), fNodeVec(0)
49 {
50  fReference = 0;
51  fNodeMaxId = -1;
52  fElementMaxId = -1;
53  fDim = -1;
54 }
55 
57 TPZSavable(cp)
58 {
59  this->operator =(cp);
60 }
61 
63 {
64  this->CleanUp();
65 
66  this->fName = cp.fName;
67  int64_t i, n = cp.fNodeVec.NElements();
68  this->fNodeVec.Resize( n );
69  for(i = 0; i < n; i++)
70  {
71  this->fNodeVec[i] = cp.fNodeVec[i];
72  }//for
73  n = cp.fElementVec.NElements();
74  this->fElementVec.Resize( n );
75  for(i = 0; i < n; i++)
76  {
77  if (cp.fElementVec[i])
78  {
79  this->fElementVec[i] = cp.fElementVec[i]->Clone(*this);
80  }
81  else
82  {
83  this->fElementVec[i] = NULL;
84  }
85  }
86 
87  this->fNodeMaxId = cp.fNodeMaxId;
88  this->fElementMaxId = cp.fElementMaxId;
90 
91  this->fReference = NULL;
92 
93  this->fDim = cp.fDim;
94  return *this;
95 }
96 
98 {
99 #ifdef PZDEBUG2
100  std::cout << "Deleting TPZGeoMesh " << (void *) this << std::endl;
101 #endif
102  CleanUp();
103 }
104 
107  int64_t i, nel = fElementVec.NElements();
108  for (i = 0; i < nel; i++) {
109  TPZGeoEl *el = fElementVec[i];
110  if (el) {
111  if(el->HasSubElement()) el->ResetSubElements();
112  if(el->FatherIndex() != -1) el->SetFatherIndex(-1);
113  delete el;
114  fElementVec[i] = 0;
115  }
116  }
117  fElementVec.Resize(0);
119  fNodeVec.Resize(0);
121  this->fNodeMaxId = -1;
122  this->fElementMaxId = -1;
123 }
124 
125 void TPZGeoMesh::SetName (const std::string &nm)
126 {
127  fName = nm;
128 }
129 
130 void TPZGeoMesh::Print (std::ostream & out) const
131 {
132  out << "\n\t\t GEOMETRIC TPZGeoMesh INFORMATIONS:\n\n";
133  out << "TITLE-> " << fName << "\n\n";
134  out << "number of nodes = " << fNodeVec.NElements() << "\n";
135  out << "number of elements = " << fElementVec.NElements()-fElementVec.NFreeElements() << "\n";
136  out << "number of free elements = " << fElementVec.NFreeElements() << "\n";
137 
138  out << "\n\tGeometric Node Information:\n\n";
139  int64_t i;
140  int64_t nnodes = fNodeVec.NElements();
141  for(i=0; i<nnodes; i++)
142  {
143  out << "Index: " << i << " ";
144  fNodeVec[i].Print(out);
145  }
146  out << "\n\tGeometric Element Information:\n\n";
147  int64_t nelem = fElementVec.NElements();
148  for(i=0; i<nelem; i++)
149  {
150  if(fElementVec[i]) fElementVec[i]->Print(out);
151  out << "\n";
152  }
153 
154  out << "\nInterface materials : \n\n";
155  InterfaceMaterialsMap::const_iterator w, e = this->fInterfaceMaterials.end();
156 
157  const int n = this->fInterfaceMaterials.size();
158  int l, r, m;
159  out << "number = " << n << "\n";
160  out << "left material / right material -> interface material\n";
161  for(w = this->fInterfaceMaterials.begin(); w != e; w++)
162  {
163  l = w->first.first;
164  r = w->first.second;
165  m = w->second;
166  out << l << " / " << r << " -> " << m << "\n";
167  }
168 }
169 void TPZGeoMesh::PrintTopologicalInfo(std::ostream & out) const
170 {
171  out << "TITLE-> " << fName << "\n";
172  out << "Number of nodes = " << fNodeVec.NElements() << "\n";
173  out << "Number of elements = " << fElementVec.NElements()-fElementVec.NFreeElements() << "\n";
174 
175  out << "\n\tGeometric Node Information:\n";
176  int64_t i;
177  int64_t nnodes = fNodeVec.NElements();
178  for(i=0; i<nnodes; i++)
179  {
180  fNodeVec[i].Print(out);
181  }
182  out << "\n\tGeometric Element Information:\n";
183  int64_t nelem = fElementVec.NElements();
184  for(i=0; i<nelem; i++)
185  {
186  if(fElementVec[i]) fElementVec[i]->PrintTopologicalInfo(out);
187  out << "\n";
188  }
189 }
190 
192 {
193  int64_t i,nnodes=nos.NElements();
194  for(i = 0; i < nnodes; i++) nodep[i] = &fNodeVec[nos[i]];
195 }
196 
198 {
199  TPZGeoEl *elp;
200  int64_t i,nelements=fElementVec.NElements();
201  for(i=0;i<nelements;i++)
202  {
203  elp = fElementVec[i];
204  if(elp) elp->ResetReference();
205  }
206  fReference = 0;
207 }
208 
210 {
211  ResetReference();
212  fReference = cmesh;
213  TPZGeoEl *gel;
214  TPZCompEl *cel;
215  int64_t i,nelem = cmesh->ElementVec().NElements();
216  for(i=0;i<nelem;i++)
217  {
218  cel = cmesh->ElementVec()[i];
219  if(cel)
220  {
221  gel = cel->Reference();
222  if(!gel)
223  {
224  PZError << "RestoreReference incomplete. Exist computational element with geometrical\n";
225  PZError << "element not belongs to the current geometrical mesh.\n";
226  return;
227  }
228  gel->SetReference(cel);
229  }
230  }
231 }
232 
233 // GetBoundaryElements returns all elements beweeen NodFrom and NodTo counterclock wise
234 // this method uses the connectivity of the elements
235 // BuildConnectivity should be called to initialize the connectivity information
236 // this method will only work for grid with 2-D topology
237 // the current version will only work for a grid with only one level
239 {
240  //Find a first element whose first node on the side is NodFrom
241  //TPZGeoEl *def = 0;
242  //TPZAVLMap<int,TPZGeoEl *> elmap(def);
243  map<int64_t,TPZGeoEl *> elmap;
244  int64_t i,nelements=NElements();
245  for(i=0;i<nelements;i++)
246  {
247  TPZGeoEl *el = fElementVec[i];
248  if(el) elmap[el->Id()]=fElementVec[i];
249  }
250 
251  int64_t currentnode = NodFrom;
252  TPZGeoEl *candidate = 0;
253  int candidateside = 0;
254  while(currentnode != NodTo)
255  {
256  // put all elements connected to currentnode in elmap, eliminate the elements
257  // from elmap which do not contain the node
258  BuildElementsAroundNode(currentnode,elmap);
259 #ifdef LOG4CXX
260  if (logger->isDebugEnabled())
261  {
262  std::stringstream sout;
263  std::map<int64_t, TPZGeoEl *>::iterator it;
264  sout << "Elements around node " << currentnode << " : ";
265  for(it=elmap.begin(); it!=elmap.end(); it++)
266  {
267  sout << it->second->Index() << "|";
268  int in;
269  for(in=0; in<it->second->NNodes(); in++)
270  {
271  sout << it->second->NodeIndex(in) << ":";
272  }
273  sout << " -- ";
274  }
275  LOGPZ_DEBUG(logger,sout.str());
276  }
277 #endif
278 
279  // find, within elmap the element which has currentnode as its first boundary side
280  // node
281  FindElement(elmap, currentnode, candidate, candidateside);
282  // if the element found is already contained in the list, we have a circular list
283  // if no element was found, the topology may not be two dimensional
284  if(!candidate)
285  {
286  LOGPZ_WARN(logger,"GetBoundaryElements no adjacent element found");
287  break;
288  }
289  int64_t index = 0;
290  int64_t nelvec = ElementVec.NElements();
291  while(index<nelvec && ElementVec[index] != candidate) index++;
292  if(index <nelvec && Sides[index]==candidateside) break;
293  ElementVec.Push(candidate);
294  Sides.Push(candidateside);
295  elmap.erase(elmap.begin(), elmap.end());//CleanUp();
296  elmap[candidate->Id()] = candidate;
297  // initialize the list in which to look for connected elements
298  currentnode = candidate->SideNodeIndex(candidateside,1);
299  }
300 }
301 
302 // Find all elements in elmap or neighbour of elements in elmap which contain a node
303 //void TPZGeoMesh::BuildElementsAroundNode(int currentnode,TPZAVLMap<int,TPZGeoEl*> &elmap){
304 void TPZGeoMesh::BuildElementsAroundNode(int64_t currentnode,map<int64_t,TPZGeoEl*> &elmap)
305 {
306  // first eliminate all elements which do not contain currentnode
307  //TPZPix iel = elmap.First();
308  map<int64_t, TPZGeoEl *>::iterator ielprev,iel=elmap.begin();
309  TPZGeoEl *el;
310  int64_t i;
311  while(iel!=elmap.end())
312  {
313  el = iel->second;
314  ielprev=iel;
315  iel++;//elmap.Next(iel);
316  int numnode = el->NNodes();
317  for(i=0; i< numnode; i++)
318  {
319  if(el->NodeIndex(i) == currentnode) break;
320  }
321  if(i == numnode)
322  {
323  elmap.erase(ielprev);
324  }
325  }
326  iel = elmap.begin();
327  while(iel!=elmap.end())
328  {
329  el = iel->second;//elmap.Contents(iel);
330  iel++;//elmap.Next(iel);
331  int nside = el->NSides();
332  for(int is=0; is<nside; is++)
333  {
334  TPZGeoElSide neigh = el->Neighbour(is);
335  if(!neigh.Exists()) continue;
336  int numnode = neigh.Element()->NNodes();
337  for(i=0; i< numnode; i++)
338  {
339  if(neigh.Element()->NodeIndex(i) == currentnode)
340  {
341  if(elmap.find(neigh.Element()->Id())==elmap.end())
342  {
343  // this should be implemented as a stack, so that we dont have to
344  // go through the list again each time
345  elmap[neigh.Element()->Id()] = neigh.Element();
346  iel = elmap.begin();
347  }
348  break; // get out of the loop over the nodes
349  }
350  }
351  }
352  }
353 }
354 
355 // find, within elmap the element which has currentnode as its first boundary side
356 // node
357 //void TPZGeoMesh::FindElement(TPZAVLMap<int,TPZGeoEl *> &elmap,int currentnode,TPZGeoEl* &candidate,int &candidateside) {
358 void TPZGeoMesh::FindElement(std::map<int64_t,TPZGeoEl *> &elmap,int64_t currentnode,TPZGeoEl* &candidate,int &candidateside)
359 {
360  candidate = 0;
361  //TPZPix iel = elmap.First();
362  map<int64_t , TPZGeoEl *>::iterator ielprev, iel = elmap.begin();
363  while(iel!=elmap.end()) {
364  TPZGeoEl *el = iel->second;//elmap.Contents(iel);
365  ielprev=iel;
366  iel++;//elmap.Next(iel);
367  int ns = el->NSides();
368  int is = el->NCornerNodes();
369  for(; is < ns; is++) {
370  TPZGeoElSide thisside(el,is);
371  if (thisside.Dimension() > 1) continue;
372  TPZGeoElSide neigh = el->Neighbour(is);
373  TPZGeoElSide father = el->Father2(is);
374 
375 #ifdef PZDEBUG2
376  std::stringstream sout;
377  sout << __PRETTY_FUNCTION__ << " for elidx " << el->Index() << std::endl;
378  sout << "\tthiside " << el->Index() << "/" << is << std::endl;
379  if (neigh.Element()) sout << "\tneighsd " << neigh.Element()->Index() << "/" << neigh.Side() << std::endl;
380  if (father.Element())sout << "\tfathers " << father.Element()->Index() << "/" << father.Side() << std::endl;
381 #endif
382 
383  if(neigh == thisside && !father.Exists() && el->SideNodeIndex(is,0) == currentnode) {
384  candidate = el;
385  candidateside = is;
386 
387 #ifdef PZDEBUG2
388  sout << "\t\t\tNew Candidate found el/side = " << el << "/" << is << std::endl;
389  LOGPZ_DEBUG (logger,sout.str().c_str());
390 #endif
391 
392  return;
393  }
394 
395 #ifdef PZDEBUG2
396  else
397  {
398  sout << "candidate doesn't match...";
399  LOGPZ_DEBUG(logger,sout.str().c_str());
400  }
401 #endif
402  }
403  }
404 
405 #ifdef LOG4CXX
406  {
407  std::stringstream sout;
408  sout << "No neighbour found for node " << currentnode << " elmap \n";
409  for(iel=elmap.begin(); iel!=elmap.end(); iel++)
410  {
411  iel->second->Print(sout);
412  }
413  LOGPZ_WARN(logger,sout.str());
414  }
415 #endif
416 
417 }
418 
420 {
421  int64_t i=0, in, nnodes = fNodeVec.NElements();
422  while(i<nnodes && fNodeVec[i].Id() == -1) i++;
423  if(i == nnodes) return 0;
424  TPZGeoNode *gnkeep = &fNodeVec[i];
425  REAL distkeep = 0.;
426  for(in=0;in<3;in++)
427  distkeep += (co[in]-(gnkeep->Coord(in)))*(co[in]-(gnkeep->Coord(in)));
428  while(i< nnodes) {
429  TPZGeoNode *gn = &fNodeVec[i];
430  REAL dist = 0.;
431  for(in=0;in<3;in++) dist += (co[in]-gn->Coord(in))*(co[in]-gn->Coord(in));
432  if(dist < distkeep)
433  {
434  gnkeep = gn;
435  distkeep = dist;
436  }
437  i++;
438  while(i<nnodes && fNodeVec[i].Id() == -1) i++;
439  }
440  return gnkeep;
441 }
442 
444 {
445  int i=0, in, nnodes = fNodeVec.NElements();
446  while(i<nnodes && fNodeVec[i].Id() == -1) i++;
447  if(i == nnodes) return 0;
448  TPZGeoNode *gnkeep = &fNodeVec[i];
449  nodeFoundIndex = i;
450 
451  REAL distkeep = 0.;
452  for(in=0;in<3;in++)
453  distkeep += (co[in]-(gnkeep->Coord(in)))*(co[in]-(gnkeep->Coord(in)));
454  while(i< nnodes) {
455  TPZGeoNode *gn = &fNodeVec[i];
456  REAL dist = 0.;
457  for(in=0;in<3;in++) dist += (co[in]-gn->Coord(in))*(co[in]-gn->Coord(in));
458  if(dist < distkeep)
459  {
460  gnkeep = gn;
461  distkeep = dist;
462  nodeFoundIndex = i;
463  }
464  i++;
465  while(i<nnodes && fNodeVec[i].Id() == -1) i++;
466  }
467  return gnkeep;
468 }
469 
472 TPZGeoEl * TPZGeoMesh::FindCloseElement(TPZVec<REAL> &x, int64_t & InitialElIndex, int targetDim) const
473 {
474  // this method will not work if targetDim == 0 because it navigates through elements of dimension targetDim
475  TPZManVector<REAL,3> xcenter(3);
476  TPZManVector<TPZManVector<REAL, 3>, 8> cornercenter;
477  // what happens if gelnext == 0 or if InitialIndex is out of scope??
478  if (InitialElIndex >= fElementVec.NElements()) {
479  DebugStop();
480  }
481  TPZGeoEl *gelnext = ElementVec()[InitialElIndex];
482 
483  if (gelnext == 0) {
484  DebugStop();
485  }
486 
487  TPZGeoEl *gel = 0;
488  REAL geldist;
489  std::map<REAL,int> cornerdist;
490  while (gel != gelnext) {
491  gel = gelnext;
492  gelnext = 0;
493  cornerdist.clear();
494  // compute the corner coordinates
495  for (int ic=0; ic<gel->NCornerNodes(); ic++) {
496  TPZManVector<REAL,3> xcorner(3);
497  gel->NodePtr(ic)->GetCoordinates(xcorner);
498  REAL dist = 0.;
499  for (int i=0; i<3; i++) {
500  dist += (x[i]-xcorner[i])*(x[i]-xcorner[i]);
501  }
502  dist = sqrt(dist);
503  cornerdist[dist]=ic;
504  }
505  // compute the distance of the center of the element
506  {
507  geldist = 0.;
508  TPZManVector<REAL,3> xcenter(3);
509  TPZGeoElSide gelside(gel,gel->NSides()-1);
510  gelside.CenterX(xcenter);
511  geldist = dist(x,xcenter);
512  }
513 #ifdef LOG4CXX
514  if (logger->isDebugEnabled()) {
515  std::stringstream sout;
516  sout << "FindClosestElement Tried element index " << gel->Index() << std::endl;
517  sout << "Distance from the center " << geldist << std::endl;
518  LOGPZ_DEBUG(logger, sout.str())
519  }
520 #endif
521  // find the closest corner node
522  REAL closestcorner = cornerdist.begin()->first;
523  // if the center node is closer than the cornernode, return the element
524  if (geldist < closestcorner || closestcorner < 1.e-15) {
525 #ifdef LOG4CXX
526  if (logger->isDebugEnabled()) {
527  std::stringstream sout;
528  sout << "Distance from the closest corner " << closestcorner << "bailing out " << std::endl;
529  LOGPZ_DEBUG(logger, sout.str())
530  }
531 #endif
532  InitialElIndex = gel->Index();
533  return gel;
534  }
535  // look for all neighbours of the corner side
536  TPZGeoElSide gelside(gel,cornerdist.begin()->second);
537  TPZGeoElSide neighbour = gelside.Neighbour();
538  std::map<REAL,int> distneigh;
539  while (neighbour != gelside)
540  {
541  if (neighbour.Element()->Dimension() == targetDim)
542  {
543  TPZManVector<REAL,3> center(3);
544  TPZGeoElSide centerneigh(neighbour.Element(),neighbour.Element()->NSides()-1);
545  centerneigh.CenterX(center);
546  REAL distcenter = dist(center,x);
547  distneigh[distcenter] = neighbour.Element()->Index();
548  }
549  neighbour = neighbour.Neighbour();
550  }
551  // choose the element whose center is closest to the coordinate
552  REAL gelnextdist = 0.;
553  if (distneigh.size() == 0) {
554  gelnext = gel;
555  gelnextdist = geldist;
556  }
557  else {
558  gelnext = ElementVec()[distneigh.begin()->second];
559  gelnextdist = distneigh.begin()->first;
560  }
561 #ifdef LOG4CXX
562  if (logger->isDebugEnabled()) {
563  std::stringstream sout;
564  sout << "Closest element index " << gelnext->Index() << std::endl;
565  sout << "Distance from the center " << gelnextdist << std::endl;
566  LOGPZ_DEBUG(logger, sout.str())
567  }
568 #endif
569  // return if its center distance is larger than the distance of the current element
570  if (geldist <= gelnextdist) {
571  InitialElIndex = gel->Index();
572  return gel;
573  }
574  }
575  InitialElIndex = gel->Index();
576  return gel;
577 }
578 
579 TPZGeoEl * TPZGeoMesh::FindElement(TPZVec<REAL> &x, TPZVec<REAL> & qsi, int64_t & InitialElIndex, int targetDim) const {
580  TPZGeoEl *res = FindApproxElement(x, qsi, InitialElIndex, targetDim);
581  TPZManVector<REAL,3> xaprox(3);
582  res->X(qsi, xaprox);
583  REAL dist = 0.;
584  for (int i=0; i<3; i++) {
585  dist += (xaprox[i]-x[i])*(xaprox[i]-x[i]);
586  }
587  dist = sqrt(dist);
588  REAL zero;
589  ZeroTolerance(zero);
590  if (dist > zero*100)
591  {
592  std::stringstream sout;
593  sout << "Element not found for coordinate x " << x << std::endl;
594  sout << "Coordinate found " << xaprox << std::endl;
595  sout << "Distance error " << dist << std::endl;
596  sout << "Closest element index " << res->Index() << " El param " << qsi << std::endl;
597  LOGPZ_ERROR(logger, sout.str())
598 // DebugStop();
599  }
600  return res;
601 }
602 
603 TPZGeoEl * TPZGeoMesh::FindElementCaju(TPZVec<REAL> &x, TPZVec<REAL> & qsi, int64_t & InitialElIndex, int targetDim)
604 {
605  FindCloseElement(x, InitialElIndex, targetDim);
606  TPZGeoEl * gel = this->ElementVec()[InitialElIndex]->LowestFather();
607 
608  qsi.Resize(targetDim, 0.);
609 
610  REAL Tol;
611  ZeroTolerance(Tol);
612 
613  if(gel->Dimension() != targetDim)
614  {
615  bool foundDim = false;
616  for(int s = 0; s < gel->NSides(); s++)
617  {
618  TPZGeoElSide gelSide(gel,s);
619  TPZGeoElSide neighSide(gelSide.Neighbour());
620  while(neighSide != gelSide)
621  {
622  if(neighSide.Element()->Dimension() == targetDim)
623  {
624  gel = neighSide.Element();
625  foundDim = true;
626  break;
627  }
628  neighSide = neighSide.Neighbour();
629  }
630  if(foundDim)
631  {
632  break;
633  }
634  }
635  }
636  if(gel->Dimension() != targetDim)
637  {
638  return NULL;
639  }
640  else if(gel->ComputeXInverse(x, qsi, Tol) == true)
641  {
642  gel = FindSubElement(gel, x, qsi, InitialElIndex);
643  return gel;
644  }
645 
646  TPZManVector<REAL,3> projection(gel->Dimension());
647  int64_t count = 0;
648  bool mustStop = false;
649  bool projectOrthogonal = true;
650  int bissectionCalled = 0;
651  while(gel->ComputeXInverse(x, qsi, Tol) == false &&
652  mustStop == false)
653  {
654  int side = -1;
655  if(projectOrthogonal)
656  {
657  side = gel->ProjectInParametricDomain(qsi, projection);
658  }
659  else
660  {
661  side = gel->ProjectBissectionInParametricDomain(qsi, projection);
662  projectOrthogonal = true;
663  }
664  TPZGeoElSide mySide(gel,side);
665  TPZGeoElSide neighSide(mySide.Neighbour());
666  bool targetNeighFound = false;
667  while(mySide != neighSide && targetNeighFound == false)
668  {
669  if(neighSide.Element()->LowestFather() == mySide.Element()->LowestFather())
670  { // it means that neighSide->element is my own Lowest father, so
671  // is not the targetNeigh... its better keep searching!!!
672  neighSide = neighSide.Neighbour();
673  }
674  else if(neighSide.Element()->Dimension() != targetDim)
675  {
676  // it means that neighSide->element() dont have the targetDim,
677  // is not the targetNeigh... its better keep searching!!!
678  neighSide = neighSide.Neighbour();
679  }
680  else
681  {
682  // it means that this neighbour is a good candidate!!!
683  // (because is not my ancestor, not even with dimension != targetDim)
684  targetNeighFound = true;
685  }
686  }
687  TPZGeoEl * neighgel = neighSide.Element()->LowestFather();
688  if(neighgel != gel)
689  {
690  gel = neighgel;
691  if(qsi.NElements() != gel->Dimension())
692  {
693  qsi.Resize(gel->Dimension(), 0.);
694  }
695  bissectionCalled = 0;
696  }
697  else
698  {
699  projectOrthogonal = false;
700  bissectionCalled++;
701  if(bissectionCalled == 2)
702  {
703  mustStop = true;
704  }
705  }
706 
707 #ifdef LOG4CXX
708  if(logger->isDebugEnabled())
709  {
710  std::stringstream sout;
711  TPZManVector<REAL,3> par(neighSide.Dimension()),xloc(3,0.);
712  neighSide.CenterPoint(par);
713  neighSide.X(par, xloc);
714  REAL dist = 0.;
715  for (int i=0; i<3; i++) {
716  dist += (x[i]-xloc[i])*(x[i]-xloc[i]);
717  }
718  dist = sqrt(dist);
719  sout << "within element " << gel->Index() << " parameter is " << qsi << " projected parameter " << projection;
720  sout << "\nBest guess is on side " << side << " which gives neighSide " << neighSide << std::endl;
721  sout << "Distance of the center point is " << dist;
722  LOGPZ_DEBUG(logger, sout.str())
723  }
724 #endif
725 
726  count++;
727  if(count > NElements())
728  {
729  mustStop = true;
730  }
731  }
732 
733  if(mustStop)
734  {
735 #ifdef PZDEBUG
736 // DebugStop();
737 #endif
738  return NULL;//not found...
739  }
740 
741  gel = FindSubElement(gel, x, qsi, InitialElIndex);
742  return gel;
743 }
744 
745 
747 TPZGeoEl *TPZGeoMesh::FindApproxElement(TPZVec<REAL> &x, TPZVec<REAL> & qsi, int64_t & InitialElIndex, int targetDim) const {
748  FindCloseElement(x, InitialElIndex, targetDim);
749  TPZGeoEl * gel = this->ElementVec()[InitialElIndex]->LowestFather();
750 
751 #ifdef LOG4CXX
752  if (logger->isDebugEnabled()) {
753  std::stringstream sout;
754  sout << "x coordinate " << x << std::endl;
755  sout << "element index " << gel->Index() << std::endl;
756  sout << "coordinates of the corner nodes of the element\n";
757  for (int i=0; i<gel->NNodes(); i++) {
758  TPZGeoNode *ptr = gel->NodePtr(i);
759  for (int ic=0; ic<3; ic++) {
760  sout << ptr->Coord(ic) << " ";
761  }
762  sout << std::endl;
763  }
764  LOGPZ_DEBUG(logger, sout.str())
765  }
766 #endif
767  qsi.Resize(targetDim, 0.);
768  qsi.Fill(0.);
769 
770  if(gel->Dimension() != targetDim)
771  {
772  // the targetDim should be respected because it is a parameter of FindCloseElement
773  DebugStop();
774  }
775  if(gel->Dimension() != targetDim)
776  {
777  DebugStop();
778  }
779  REAL zero;
780  ZeroTolerance(zero);
781 
782  std::set<TPZGeoEl *> tested;
783  // this method will call ComputeXInverse if the element dimension != 3
784  bool memberQ = gel->ComputeXInverse(x, qsi,zero);
785  if(memberQ)
786  {
787 #ifdef LOG4CXX
788  if (logger->isDebugEnabled())
789  {
790  LOGPZ_DEBUG(logger, "Going into the FindSubElement alternative")
791  }
792 #endif
793  gel = FindSubElement(gel, x, qsi, InitialElIndex);
794  return gel;
795  }
796  tested.insert(gel);
797 #ifdef LOG4CXX
798  if (logger->isDebugEnabled())
799  {
800  std::stringstream sout;
801  sout << "Looking for x = " << x << std::endl;
802  sout << "Tried out gel index " << gel->Index() << std::endl;
803  sout << "found qsi " << qsi << std::endl;
804  TPZManVector<REAL> locx(3);
805  gel->X(qsi, locx);
806  sout << "found x co " << locx << std::endl;
807  LOGPZ_DEBUG(logger, sout.str())
808  }
809 #endif
810 
811  TPZManVector<REAL,3> projection(gel->Dimension());
812  int side = -1;
813  side = gel->ProjectInParametricDomain(qsi, projection);
814 
815 
816  TPZGeoElSide mySide(gel,side);
817  TPZManVector<REAL,3> xclose(3);
818  gel->X(projection, xclose);
819  REAL mindist = 0.;
820  for (int i=0; i<3; i++) {
821  mindist += (xclose[i]-x[i])*(xclose[i]-x[i]);
822  }
823  TPZGeoElSide bestside(mySide);
824  TPZManVector<REAL,3> bestproj(projection);
825 
826 
827  TPZStack<TPZGeoElSide> allneigh;
828  TPZGeoElSide neighSide(mySide.Neighbour());
829  while (neighSide != mySide) {
830  if (neighSide.Element()->Dimension() == targetDim && tested.find(neighSide.Element()) == tested.end()) {
831  allneigh.Push(neighSide);
832  }
833  neighSide = neighSide.Neighbour();
834  }
835 
836  while (allneigh.size())
837  {
838  TPZGeoElSide gelside = allneigh.Pop();
839  TPZGeoEl *locgel = gelside.Element();
840  if (tested.find(locgel) != tested.end()) {
841  continue;
842  }
843  qsi.Fill(0.);
844  if(locgel->ComputeXInverse(x, qsi, zero*100.) == true)
845  {
846 #ifdef LOG4CXX
847  if (logger->isDebugEnabled())LOGPZ_DEBUG(logger, "FOUND ! Going into the FindSubElement alternative")
848 #endif
849  gel = FindSubElement(locgel, x, qsi, InitialElIndex);
850  return gel;
851  }
852 
853 #ifdef LOG4CXX
854  if (logger->isDebugEnabled())
855  {
856  std::stringstream sout;
857  sout << "Looking for x = " << x << std::endl;
858  sout << "Tried out gel index " << locgel->Index() << std::endl;
859  sout << "found qsi " << qsi << std::endl;
860  TPZManVector<REAL> locx(3);
861  gel->X(qsi, locx);
862  sout << "found x co " << locx << std::endl;
863  LOGPZ_DEBUG(logger, sout.str())
864  }
865 #endif
866 
867 
868  tested.insert(locgel);
869  int side = -1;
870  side = locgel->ProjectInParametricDomain(qsi, projection);
871  TPZManVector<REAL,3> xclose(3);
872  TPZGeoElSide locgelside(locgel,side);
873  locgelside.X(projection, xclose);
874  REAL dist = 0.;
875  for (int i=0; i<3; i++) {
876  dist += (xclose[i]-x[i])*(xclose[i]-x[i]);
877  }
878  if (dist < mindist) {
879  mindist = dist;
880  bestside = locgelside;
881  bestproj = projection;
882  }
883  mySide = TPZGeoElSide(locgel,side);
884  TPZGeoElSide neighSide(mySide.Neighbour());
885  while (neighSide != mySide) {
886  if (neighSide.Element()->Dimension() == targetDim && tested.find(neighSide.Element()) == tested.end()) {
887  allneigh.Push(neighSide);
888  }
889  neighSide = neighSide.Neighbour();
890  }
891  }
892  TPZGeoEl *bestgel = bestside.Element();
893  qsi = bestproj;
894  gel = FindSubElement(bestgel, x, qsi, InitialElIndex);
895  return gel;
896 
897  return bestgel;
898 }
899 
900 TPZGeoEl * TPZGeoMesh::FindSubElement(TPZGeoEl * gel, TPZVec<REAL> &x, TPZVec<REAL> & qsi, int64_t & InitialElIndex) const {
901  REAL Tol;
902  ZeroTolerance(Tol);
903 
904 #ifdef PZDEBUG
905  TPZManVector<REAL,3> locqsi(qsi);
906  if(gel->ComputeXInverse(x,locqsi,Tol*100.) == false)
907  {
908  //The given gel does NOT contains the given x!!!
909  std::cout << "FindSubElement called for non conforming point\n";
910  }
911 #endif
912  TPZManVector<REAL,3> locx(x);
913 
914  if(gel->HasSubElement() == false)
915  {
916  InitialElIndex = gel->Index();
917  return gel;
918  }
919  else
920  {
921  TPZStack<TPZGeoEl*> subElements;
922  gel->GetAllSiblings(subElements);
923 
924  int nsons = subElements.NElements();
925  TPZGeoEl * son = NULL;
926  TPZVec< TPZVec<REAL> > qsiSonVec(nsons);
927  for(int s = 0; s < nsons; s++)
928  {
929  son = subElements[s];
930  qsiSonVec[s].Resize(son->Dimension(),0.);
931  if(son->ComputeXInverse(locx, qsiSonVec[s],Tol*100.))
932  {
933  qsi = qsiSonVec[s];
934  InitialElIndex = son->Index();
935  return FindSubElement(son, locx, qsi, InitialElIndex);
936  }
937  }
938 
939  //no son found (accuracy problems)
940  {
941  std::map<REAL,int> dist;
942  for(int s = 0; s < nsons; s++)
943  {
944  son = subElements[s];
945  TPZVec<REAL> qsiSonProj(son->Dimension());
946  son->ProjectInParametricDomain(qsiSonVec[s], qsiSonProj);
947  REAL distToProj = 0.;
948  for(int c = 0; c < son->Dimension(); c++)
949  {
950  distToProj += (qsiSonProj[c] - qsiSonVec[s][c]) * (qsiSonProj[c] - qsiSonVec[s][c]);
951  }
952  dist[sqrt(distToProj)] = s;
953  }
954  int sonPosition = dist.begin()->second;//smaller distance to parametric domain
955 
956  qsi = qsiSonVec[sonPosition];
957  son = subElements[sonPosition];
958  InitialElIndex = son->Index();
959 
960  return FindSubElement(son, locx, qsi, InitialElIndex);
961  }
962  }
963 
964  return NULL;
965 }
966 
968 {
969 
970  TPZVec<int> SideNum(NNodes(),-1);
971  TPZVec<TPZGeoEl *> NeighNode(NNodes(),0);
972  int64_t nelem = NElements();
973  int64_t iel = 0;
974  for(iel=0; iel<nelem; iel++)
975  {
976  TPZGeoEl *gel = fElementVec[iel];
977  if(!gel) continue;
978  int ncor = gel->NCornerNodes();
979  int in;
980  for(in=0; in<ncor; in++) {
981  int64_t nod = gel->NodeIndex(in);
982  if(SideNum[nod] == -1)
983  {
984  NeighNode[nod] = gel;
985  SideNum[nod] = in;
986  if(gel->SideIsUndefined(in)) gel->SetSideDefined(in);
987  }
988  else
989  {
990  TPZGeoElSide neigh(NeighNode[nod],SideNum[nod]);
991  TPZGeoElSide gelside(gel,in);
992  if(!neigh.NeighbourExists(gelside))
993  {
994  neigh.SetConnectivity(gelside);
995  }
996  }
997  }
998  }
999  for(iel=0; iel<nelem; iel++)
1000  {
1001  TPZGeoEl *gel = fElementVec[iel];
1002  if(!gel) continue;
1003  int ncor = gel->NCornerNodes();
1004  int nsides = gel->NSides();
1005  int is;
1006  for(is=ncor; is<nsides; is++)
1007  {
1008  if( gel->SideIsUndefined(is))
1009  {
1010  gel->SetSideDefined(is);
1011  TPZGeoElSide gelside(gel,is);
1012  TPZStack<TPZGeoElSide> neighbours;
1013  gelside.ComputeNeighbours(neighbours);
1014  int64_t nneigh = neighbours.NElements();
1015  int64_t in;
1016  for(in=0; in<nneigh; in++) {
1017  if(neighbours[in].Side() == -1)
1018  {
1019  std::cout << "TPZGeoMesh::BuildConnectivity : Inconsistent mesh detected analysing element/side:" ;
1020  std::cout << gelside ;
1021  std::cout << std::endl;
1022  continue;
1023  }
1024  if(neighbours[in].Element()->SideIsUndefined(neighbours[in].Side()))
1025  {
1026  neighbours[in].Element()->SetSideDefined(neighbours[in].Side());
1027  }
1028  gelside.SetConnectivity(neighbours[in]);
1029  }
1030  }
1031  }
1032  }
1033 
1034  //Verify node coordinates for curved elements
1035 #ifdef PZDEBUG
1036  const int64_t nel = this->NElements();
1037  for(int64_t el = 0; el < nel; el++)
1038  {
1039  TPZGeoEl * gel = this->ElementVec()[el];
1040  if(!gel) continue;
1041  gel->VerifyNodeCoordinates();
1042  }//for el
1043 #endif
1044 
1045  //Build the data structure of blend elements
1046  int64_t Qelem = this->NElements();
1047  for(int64_t el = 0; el < Qelem; el++)
1048  {
1049  TPZGeoEl * gel = this->ElementVec()[el];
1050  if(!gel) continue;
1051  gel->BuildBlendConnectivity();
1052  }//for el
1053 
1054 }
1055 
1057 
1058  TPZVec<int> SideNum(NNodes(),-1);
1059  TPZVec<TPZGeoEl *> NeighNode(NNodes(),0);
1060  int64_t nelem = NElements();
1061  int64_t iel = 0;
1062  while(iel<nelem && fElementVec[iel] == 0) iel++;
1063 
1064  int64_t numsearch =1;
1065  // if there are no elements, do nothing
1066  while(iel < nelem) {
1067  TPZGeoEl *el = fElementVec[iel];
1068  int numsides = el->NSides();
1069  int side;
1070  for(side = 0;side<numsides;side++)
1071  {
1072  // check whether all entries in NeighNode are equal
1073 
1074  int equalnode = 1;
1075  int numsidenodes = el->NSideNodes(side);
1076  int64_t sidenode = el->SideNodeIndex(side,0);
1077  TPZGeoEl *neigh = NeighNode[sidenode];
1078  int sidenumber = SideNum[sidenode];
1079  for(int64_t sn = 0;sn < numsidenodes; sn++)
1080  {
1081  sidenode = el->SideNodeIndex(side,sn);
1082  if (neigh != NeighNode[sidenode])
1083  {
1084  equalnode=0;
1085  break;
1086  }
1087  }
1088 
1089  if(equalnode && neigh == 0)
1090  {
1091  if(el->SideIsUndefined(side))
1092  {
1093  int elloaded = 0;
1094  for(int in=0; in<el->NNodes(); in++)
1095  {
1096  if(NeighNode[el->NodeIndex(in)] == el) elloaded = 1;
1097  }
1098  // this element is not loaded and its side is undefined
1099 
1100  // load the element side in the NeighNode vector
1101  for(int sn=0;!elloaded && sn < numsidenodes; sn++) \
1102  {
1103  sidenode = el->SideNodeIndex(side,sn);
1104  NeighNode[sidenode] = el;
1105  SideNum[sidenode] = side;
1106  }
1107  numsearch++;
1108  }
1109  } else if(equalnode && side == sidenumber && neigh == el)
1110  {
1111  // unload the element side
1112  for(int sn=0;sn < numsidenodes; sn++)
1113  {
1114  sidenode = el->SideNodeIndex(side,sn);
1115  NeighNode[sidenode] = 0;
1116  SideNum[sidenode] = -1;
1117  }
1118  // if no neighbouring element was detected during the loop
1119  // define the element side as undefined
1120  TPZGeoElSide neighbour = el->Neighbour(side);
1121  if(!neighbour.Exists()) el->SetSideDefined(side);
1122  numsearch++;
1123  }
1124  else if(equalnode && neigh != el)
1125  {
1126  // we found a neigbour
1127  TPZManVector<int64_t> SideNodes(numsidenodes);
1128  // detect which side of the neigbour is loaded witin NeighNode
1129  for(int sn=0;sn < numsidenodes; sn++)
1130  {
1131  sidenode = el->SideNodeIndex(side,sn);
1132  SideNodes[sn] = sidenode;
1133  }
1134  int neighside = neigh->WhichSide(SideNodes);
1135  TPZGeoElSide neighbour(neigh,neighside);
1136  // WhichSide will tell the side number which contains the vector
1137  // of node ids SideNodes
1138  if(neighbour.Side() != -1 && !el->NeighbourExists(side,neighbour))
1139  {
1140  TPZGeoElSide(el,side).SetConnectivity(neighbour);
1141  numsearch++;
1142  }
1143  }
1144  } // loop over the sides
1145  iel++;
1146  while(iel<nelem && fElementVec[iel] == 0) iel++;
1147  if(iel==nelem && numsearch)
1148  {
1149  numsearch = 0;
1150  iel = 0;
1151  while(iel<nelem && fElementVec[iel] == 0) iel++;
1152  }
1153  }
1154 }
1155 
1156 //Cedric : 03/03/99
1158 {
1159  int64_t nel = fElementVec.NElements();
1160  TPZGeoEl *gel = 0;
1161  for(int64_t i=0;i<nel;i++)
1162  {
1163  gel = fElementVec[i];
1164  if(gel && gel->Id()==elid) break;
1165  }
1166  return gel;
1167 }
1168 
1170 {
1171  int64_t i=0;
1172  int64_t index = gel->Index();
1173  if (ElementVec()[index] == gel) return index;
1174 
1175  int64_t numel = ElementVec().NElements();
1176  while ( i < numel )
1177  {
1178  if (ElementVec()[i] == gel) break;
1179  i++;
1180  }
1181  if(i<numel) return i;
1182  else return -1;
1183 }
1184 
1186 {
1187  int64_t i=0;
1188  int64_t numel = NodeVec().NElements();
1189  while ( i < numel )
1190  {
1191  if (&NodeVec()[i] == nod) break;
1192  i++;
1193  }
1194  return i;
1195 }
1196 
1197 #include "TPZGeoElement.h"
1198 #include "TPZGeoCube.h"
1199 #include "pzshapecube.h"
1200 #include "TPZRefCube.h"
1201 #include "pzshapelinear.h"
1202 #include "TPZGeoLinear.h"
1203 #include "TPZRefLinear.h"
1204 #include "pzrefquad.h"
1205 #include "pzshapequad.h"
1206 #include "pzgeoquad.h"
1207 #include "pzshapetriang.h"
1208 #include "pzreftriangle.h"
1209 #include "pzgeotriangle.h"
1210 #include "pzshapeprism.h"
1211 #include "pzrefprism.h"
1212 #include "pzgeoprism.h"
1213 #include "pzshapetetra.h"
1214 #include "pzreftetrahedra.h"
1215 #include "pzgeotetrahedra.h"
1216 #include "pzshapepiram.h"
1217 #include "pzrefpyram.h"
1218 #include "pzgeopyramid.h"
1219 #include "pzgeopoint.h"
1220 #include "pzrefpoint.h"
1221 #include "pzshapepoint.h"
1222 #include "Hash/TPZHash.h"
1223 #include "tpzgeoelmapped.h"
1224 
1225 using namespace pzgeom;
1226 using namespace pzrefine;
1227 using namespace pzshape;
1228 
1229 TPZGeoEl *TPZGeoMesh::CreateGeoElementMapped(MElementType type, TPZVec<int64_t>& nodeindexes, int matid, int64_t& index)
1230 {
1231  switch( type ){
1232  case 0://point
1233  {
1234  TPZGeoEl * gel =
1235  new TPZGeoElMapped<TPZGeoElRefPattern<TPZGeoPoint> > (nodeindexes, matid, *this, index);
1236  return gel;
1237  }
1238  case 1://line
1239  {
1240  TPZGeoEl *gel =
1242  (nodeindexes, matid, *this, index);
1243  return gel;
1244  }
1245  case 2://triangle
1246  {
1247  TPZGeoEl *gel =
1249  (nodeindexes, matid, *this, index);
1250  return gel;
1251  }
1252  case 3://quadrilatera
1253  {
1254  TPZGeoEl* gel =
1256  (nodeindexes, matid, *this, index);
1257  return gel;
1258  }
1259  case 4://tetraedra
1260  {
1261  TPZGeoEl*gel =
1263  (nodeindexes, matid, *this, index);
1264  return gel;
1265  }
1266  case 5://pyramid
1267  {
1268  TPZGeoEl *gel =
1270  (nodeindexes, matid, *this, index);
1271  return gel;
1272  }
1273  case 6://prism
1274  {
1275  TPZGeoEl*gel =
1277  (nodeindexes, matid, *this, index);
1278  return gel;
1279  }
1280  case 7://cube
1281  {
1282  TPZGeoEl*gel =
1284  (nodeindexes, matid, *this, index);
1285  return gel;
1286  }
1287  default:
1288  {
1289  PZError << "TPZGeoMesh::CreateGeoElement type element not exists:"
1290  << " type = " << type << std::endl;
1291  return NULL;
1292  }
1293  }
1294 }
1295 
1297  TPZVec<int64_t>& nodeindexes,
1298  int matid,
1299  int64_t& index,
1300  int reftype){
1301 
1302 #ifdef PZDEBUG
1303  {
1304  for (int i=0; i<nodeindexes.size(); i++) {
1305  if (nodeindexes[i] < 0) {
1306  DebugStop();
1307  }
1308  }
1309  }
1310 #endif
1311  if (reftype == 0)
1312  {
1313  switch( type )
1314  {
1315  case 0://point
1316  return new TPZGeoElement< TPZGeoPoint, TPZRefPoint>(nodeindexes, matid, *this, index );
1317  case 1://line
1318  return new TPZGeoElement< TPZGeoLinear, TPZRefLinear>(nodeindexes, matid, *this, index );
1319  case 2://triangle
1320  return new TPZGeoElement< TPZGeoTriangle, TPZRefTriangle >(nodeindexes, matid, *this, index );
1321  case 3://quadrilatera
1322  return new TPZGeoElement< TPZGeoQuad, TPZRefQuad >(nodeindexes, matid, *this, index );
1323  case 4://tetraedra
1324  return new TPZGeoElement< TPZGeoTetrahedra, TPZRefTetrahedra >(nodeindexes, matid, *this, index );
1325  case 5:
1326  return new TPZGeoElement< TPZGeoPyramid, TPZRefPyramid >(nodeindexes, matid, *this, index );
1327  case 6:
1328  return new TPZGeoElement< TPZGeoPrism, TPZRefPrism >(nodeindexes, matid, *this, index );
1329  case 7:
1330  return new TPZGeoElement< TPZGeoCube, TPZRefCube >(nodeindexes, matid, *this, index );
1331  default:
1332  PZError << "TPZGeoMesh::CreateGeoElement type element not exists:"
1333  << " type = " << type << std::endl;
1334  return NULL;
1335  }
1336  }
1337  else
1338  {
1339  //TPZAutoPointer<TPZRefPattern> ref = GetUniformPattern(type);
1340  switch( type )
1341  {
1342  case 0://point
1343  {
1344  TPZGeoElRefPattern<TPZGeoPoint> * gel = new TPZGeoElRefPattern<TPZGeoPoint> (nodeindexes, matid, *this, index);
1345  return gel;
1346  }
1347  case 1://line
1348  {
1349  TPZGeoElRefPattern < TPZGeoLinear > *gel = new TPZGeoElRefPattern < TPZGeoLinear >(nodeindexes, matid, *this, index);
1350  //gel->SetRefPattern (ref);
1351  return gel;
1352  }
1353  case 2://triangle
1354  {
1355  TPZGeoElRefPattern < TPZGeoTriangle > *gel = new TPZGeoElRefPattern < TPZGeoTriangle >(nodeindexes, matid, *this, index);
1356  //gel->SetRefPattern (ref);
1357  return gel;
1358  }
1359  case 3://quadrilatera
1360  {
1361  TPZGeoElRefPattern < TPZGeoQuad > * gel = new TPZGeoElRefPattern < TPZGeoQuad >(nodeindexes, matid, *this, index);
1362  //gel->SetRefPattern (ref);
1363  return gel;
1364  }
1365  case 4://tetraedra
1366  {
1367  TPZGeoElRefPattern < TPZGeoTetrahedra > *gel = new TPZGeoElRefPattern < TPZGeoTetrahedra >(nodeindexes, matid, *this, index);
1368  //gel->SetRefPattern (ref);
1369  return gel;
1370  }
1371  case 5://pyramid
1372  {
1373  TPZGeoElRefPattern < TPZGeoPyramid > *gel = new TPZGeoElRefPattern < TPZGeoPyramid >(nodeindexes, matid, *this, index);
1374  //gel->SetRefPattern (ref);
1375  return gel;
1376  }
1377  case 6://prism
1378  {
1379  TPZGeoElRefPattern < TPZGeoPrism > *gel = new TPZGeoElRefPattern < TPZGeoPrism >(nodeindexes, matid, *this, index);
1380  //gel->SetRefPattern (ref);
1381  return gel;
1382  }
1383  case 7://cube
1384  {
1385  TPZGeoElRefPattern < TPZGeoCube > *gel = new TPZGeoElRefPattern < TPZGeoCube >(nodeindexes, matid, *this, index);
1386  //gel->SetRefPattern (ref);
1387  return gel;
1388  }
1389  default:
1390  {
1391  PZError << "TPZGeoMesh::CreateGeoElement type element not exists:"
1392  << " type = " << type << std::endl;
1393  return NULL;
1394  }
1395  }
1396  }
1397  //return NULL;
1398 }
1399 
1400 TPZGeoEl *TPZGeoMesh::CreateGeoBlendElement(MElementType type, TPZVec<int64_t>& nodeindexes, int matid, int64_t& index)
1401 {
1402  switch( type ){
1403  case 0://point
1404  {
1405  TPZGeoEl *gel = new TPZGeoElRefPattern<TPZGeoBlend<pzgeom::TPZGeoPoint> > (nodeindexes,matid,*this);
1406  return gel;
1407  }
1408  case 1://line
1409  {
1410  TPZGeoEl *gel = new TPZGeoElRefPattern<TPZGeoBlend<pzgeom::TPZGeoLinear> > (nodeindexes,matid,*this);
1411  return gel;
1412  }
1413  case 2://triangle
1414  {
1415  TPZGeoEl *gel = new TPZGeoElRefPattern<TPZGeoBlend<pzgeom::TPZGeoTriangle> > (nodeindexes,matid,*this);
1416  return gel;
1417  }
1418  case 3://quadrilateral
1419  {
1420  TPZGeoEl *gel = new TPZGeoElRefPattern<TPZGeoBlend<pzgeom::TPZGeoQuad> > (nodeindexes,matid,*this);
1421  return gel;
1422  }
1423  case 4://tetraedron
1424  {
1425  TPZGeoEl *gel = new TPZGeoElRefPattern<TPZGeoBlend<pzgeom::TPZGeoTetrahedra> > (nodeindexes,matid,*this);
1426  return gel;
1427  }
1428  case 5://pyramid
1429  {
1430  TPZGeoEl *gel = new TPZGeoElRefPattern<TPZGeoBlend<pzgeom::TPZGeoPyramid> > (nodeindexes,matid,*this);
1431  return gel;
1432  }
1433  case 6://prism
1434  {
1435  TPZGeoEl *gel = new TPZGeoElRefPattern<TPZGeoBlend<pzgeom::TPZGeoPrism> > (nodeindexes,matid,*this);
1436  return gel;
1437  }
1438  case 7://cube
1439  {
1440  TPZGeoEl *gel = new TPZGeoElRefPattern<TPZGeoBlend<pzgeom::TPZGeoCube> > (nodeindexes,matid,*this);
1441  return gel;
1442  }
1443  default:
1444  {
1445  PZError << "TPZGeoMesh::CreateGeoElementRefPattern type element not exists:"
1446  << " type = " << type << std::endl;
1447  return NULL;
1448  }
1449  }
1450 }
1451 
1453  return Hash("TPZGeoMesh");
1454 }
1455 
1456 void TPZGeoMesh::DeleteElement(TPZGeoEl *gel,int64_t index)
1457 {
1458  if(!gel) DebugStop();
1459 
1460  if(index < 0 || gel != fElementVec[index])
1461  {
1462  index = ElementIndex(gel);
1463  if(index < 0)
1464  {
1465  PZError << "TPZGeoMesh::DeleteElement index error\n";
1466  return;
1467  }
1468  }
1469  if (gel->HasSubElement())
1470  {
1471  for (int i=0;i<gel->NSubElements();i++)
1472  {
1473  TPZGeoEl* son = gel->SubElement(i);
1474  DeleteElement(son,son->Index());
1475  }
1476  }
1477  gel->RemoveConnectivities();
1478  if(gel) delete gel;
1479  fElementVec[index] = NULL; //this is already called in the ~TPZGeoEl(), but twice is not a problem
1480  //fElementVec.SetFree(index); this is already called in the ~TPZGeoEl()
1481 }
1482 
1483 #ifndef BORLAND
1484 template class TPZRestoreClass<TPZGeoMesh>;
1485 #endif
1486 
1487 void TPZGeoMesh::Read(TPZStream &buf, void *context) { //ok
1488  buf.Read(&fName, 1);
1491  fNodeVec.Read(buf, context);
1492  buf.Read(&fNodeMaxId);
1493  buf.Read(&fElementMaxId);
1494  buf.Read(&fDim);
1495  int64_t ninterfacemaps;
1496  buf.Read(&ninterfacemaps);
1497  int64_t c;
1498  for (c = 0; c < ninterfacemaps; c++) {
1499  int vals[3];
1500  buf.Read(vals, 3);
1501  fInterfaceMaterials[pair<int, int>(vals[0], vals[1])] = vals[2];
1502  }
1503 }
1504 
1505 void TPZGeoMesh::Write(TPZStream &buf, int withclassid) const { //ok
1506 #ifdef LOG4CXX
1507  if (logger->isDebugEnabled()) {
1508  LOGPZ_DEBUG(logger, __PRETTY_FUNCTION__);
1509  }
1510 #endif
1511  buf.Write(&fName);
1512 // if (fReference) {
1513 // std::cout << __PRETTY_FUNCTION__ << "TPZGeoMesh::Write Trying to write a gmesh with non-null reference. Call ResetReference() and try again." << std::endl;
1514 // std::cout.flush();
1515 // DebugStop();
1516 // }
1519  fNodeVec.Write(buf, withclassid);
1520  buf.Write(&fNodeMaxId);
1521  buf.Write(&fElementMaxId);
1522  buf.Write(&fDim);
1523  int64_t ninterfacemaps = fInterfaceMaterials.size();
1524  buf.Write(&ninterfacemaps);
1525  for (auto elem : fInterfaceMaterials) {
1526  int vals[3];
1527  vals[0] = elem.first.first;
1528  vals[1] = elem.first.second;
1529  vals[2] = elem.second;
1530  buf.Write(vals, 3);
1531  }
1532 }
1533 
1534 int TPZGeoMesh::AddInterfaceMaterial(int leftmaterial, int rightmaterial, int interfacematerial)
1535 {
1536  std::pair<int, int> leftright(leftmaterial, rightmaterial);
1537  InterfaceMaterialsMap::iterator w, e;
1538  e = fInterfaceMaterials.end();
1539  w = fInterfaceMaterials.find(leftright);
1540  if (w == e)
1541  {
1542  fInterfaceMaterials[leftright] = interfacematerial;
1543  return 1;
1544  }
1545  return 0;
1546 }
1547 
1548 int TPZGeoMesh::InterfaceMaterial(int leftmaterial, int rightmaterial)
1549 {
1550  std::pair<int, int> leftright(leftmaterial, rightmaterial);
1551  InterfaceMaterialsMap::iterator w, e;
1552  e = fInterfaceMaterials.end();
1553 
1554  //trying to find an interface material associated to left and right materials
1555  w = fInterfaceMaterials.find(leftright);
1556  if (w != e) return w->second;
1557 
1558  //when left and right are equal and no exception case was inserted in interfacematerialmap return the same material of left and right
1559  if (leftmaterial == rightmaterial) return leftmaterial;
1560 
1561  //error message
1562  std::stringstream mess;
1563  mess << "\nTPZGeoMesh::InterfaceMaterial - Interface material not found left " << leftmaterial << " right " << rightmaterial ;
1564  PZError << mess.str() << std::endl;
1565  return GMESHNOMATERIAL;
1566 }
1567 
1569 {
1570  InterfaceMaterialsMap::iterator b, e;
1571  b = fInterfaceMaterials.begin();
1572  e = fInterfaceMaterials.end();
1573  fInterfaceMaterials.erase(b, e);
1574 }
1575 
1577 {
1578  TPZGeoElSide side;
1579  int64_t iel;
1580  const int64_t nelem = this->NElements();
1581  for(iel = 0; iel < nelem; iel++)
1582  {
1583  TPZGeoEl * gel = this->ElementVec()[iel];
1584  if (!gel) continue;
1585  const int nsides = gel->NSides();
1586  int is;
1587  for(is = 0; is < nsides; is++)
1588  {
1589  gel->SetNeighbour(is, side);
1590  }
1591  }
1592  this->SetName("reset");
1593 }
1594 
1596 REAL TPZGeoMesh::Area(int matid)
1597 {
1598  std::set<int> matids;
1599  matids.insert(matid);
1600  return Area(matids);
1601 }
1602 
1605 {
1606  int64_t nel = NElements();
1607  int meshdim = Dimension();
1608  std::set<int> matids;
1609  for (int64_t el=0; el<nel; el++) {
1610  TPZGeoEl *gel = Element(el);
1611  if (!gel || gel->Dimension() != meshdim) {
1612  continue;
1613  }
1614  matids.insert(gel->MaterialId());
1615  }
1616  return Area(matids);
1617 }
1618 
1620 REAL TPZGeoMesh::Area(std::set<int> &matids)
1621 {
1622  TPZVec<int> NeedsComputing(NElements(),1);
1623  int meshdim = Dimension();
1624  int64_t nel = NElements();
1625  REAL result = 0.;
1626  for (int64_t el=0; el<nel; el++) {
1627  TPZGeoEl *gel = ElementVec()[el];
1628  if (!gel || !NeedsComputing[el]) {
1629  continue;
1630  }
1631  if (gel->Dimension() != meshdim || matids.find(gel->MaterialId()) == matids.end() || gel->HasSubElement()) {
1632  NeedsComputing[el] = 0;
1633  continue;
1634  }
1635  TPZGeoElSide gelside(gel,gel->NSides()-1);
1636  TPZGeoElSide neighbour = gelside.Neighbour();
1637  while (neighbour != gelside) {
1638  int64_t neighindex = neighbour.Element()->Index();
1639  NeedsComputing[neighindex] = 0;
1640  neighbour = neighbour.Neighbour();
1641  }
1642  result += gel->SideArea(gel->NSides()-1);
1643  }
1644  return result;
1645 }
virtual int64_t SideNodeIndex(int side, int nodenum) const =0
Returns the index of the nodenum node of side.
void GetNodePtr(TPZVec< int64_t > &nos, TPZVec< TPZGeoNode *> &nodep)
Fills the nodep vector with pointers to the nodes identified by their indexes.
Definition: pzgmesh.cpp:191
Contains TPZShapeTetra class which implements the shape functions of a tetrahedral element...
int fDim
dimension of the geometric domain
Definition: pzgmesh.h:70
virtual TPZGeoElSide Father2(int side) const
Returns the father/side of the father which contains the side of the sub element. ...
Definition: pzgeoel.cpp:376
const int64_t numel
Number of elements to test.
Definition: pzsubcmesh.cpp:47
TPZGeoNode * NodePtr(int i) const
Returns a pointer to the ith node of the element.
Definition: pzgeoel.cpp:2566
TPZGeoEl * FindSubElement(TPZGeoEl *gel, TPZVec< REAL > &x, TPZVec< REAL > &qsi, int64_t &InitialElIndex) const
Returns the subelement that contains the given point x and it respective point in parametric domain q...
Definition: pzgmesh.cpp:900
Contains definitions to LOGPZ_DEBUG, LOGPZ_INFO, LOGPZ_WARN, LOGPZ_ERROR and LOGPZ_FATAL, and the implementation of the inline InitializePZLOG(string) function using log4cxx library or not. It must to be called out of "#ifdef LOG4CXX" scope.
void BuildElementsAroundNode(int64_t currentnode, std::map< int64_t, TPZGeoEl *> &elmap)
Find all elements in elmap or neighbour of elements in elmap which contain a node.
Definition: pzgmesh.cpp:304
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
TPZAdmChunkVector< TPZGeoEl * > fElementVec
List of pointers to finite elements.
Definition: pzgmesh.h:58
virtual void SetNeighbour(int side, const TPZGeoElSide &neighbour)=0
Fill in the data structure for the neighbouring information.
Contains declaration of TPZGeoNode class which defines a geometrical node.
REAL Coord(int i) const
Returns i-th coordinate of the current node.
Definition: pzgnode.h:99
void CleanUp()
Deletes all items in the TPZGeoMesh.
Definition: pzgmesh.cpp:106
Contains declaration of TPZCompEl class which defines the interface of a computational element...
Templated vector implementation.
void WritePointers(const TPZVec< TPZAutoPointer< T >> &vec)
Definition: TPZStream.h:214
Implements a generic geometric element with a uniform refinement pattern. Geometry.
Definition: TPZGeoElement.h:22
void Write(TPZStream &buf, int withclassid) const override
Writes this object to the TPZStream buffer. Include the classid if withclassid = true.
Definition: pzgmesh.cpp:1505
groups all classes dedicated to the computation of shape functions
Definition: pzshapeextend.h:16
Defines PZError.
virtual int NCornerNodes() const =0
Returns the number of corner nodes of the element.
Contains the TPZRefQuad class which implements the uniform refinement of a geometric quadrilateral el...
int64_t NElements() const
Number of elements of the mesh.
Definition: pzgmesh.h:129
iterator begin()
int64_t fNodeMaxId
Maximum id used by all nodes of this mesh.
Definition: pzgmesh.h:64
Declarates the TPZBlock<REAL>class which implements block matrices.
int64_t fElementMaxId
Maximum id used by all elements of this mesh.
Definition: pzgmesh.h:67
Utility class which represents an element with its side. The Geometric approximation classes Geometry...
Definition: pzgeoelside.h:83
InterfaceMaterialsMap fInterfaceMaterials
Datastructure which indicates the index of the interface material which needs to be created between t...
Definition: pzgmesh.h:78
int64_t FatherIndex()
Definition: pzgeoel.h:349
int64_t NElements() const
Access method to query the number of elements of the vector.
#define LOGPZ_WARN(A, B)
Define log for warnings.
Definition: pzlog.h:91
TPZGeoElSide Neighbour() const
Definition: pzgeoel.h:754
int ClassId() const override
Define the class id associated with the class.
Definition: pzgmesh.cpp:1452
virtual int NSides() const =0
Returns the number of connectivities of the element.
static TPZSavable * GetInstance(const int64_t &objId)
virtual int64_t NodeIndex(int i) const =0
Returns the index of the ith node the index is the location of the node in the nodevector of the mesh...
void ReadPointers(TPZVec< TPZAutoPointer< T >> &vec)
Definition: TPZStream.h:305
This class implements a geometric element which uses its ancestral to compute its jacobian...
virtual TPZGeoEl * CreateGeoBlendElement(MElementType type, TPZVec< int64_t > &nodeindexes, int matid, int64_t &index)
Creates a geometric element in same fashion of CreateGeoElement but here the elements are blend...
Definition: pzgmesh.cpp:1400
static TPZGraphNode gn
Definition: pzgraphmesh.cpp:76
int WhichSide(TPZVec< int64_t > &SideNodeIds)
Returns the side number which is connected to the SideNodes returns -1 if no side is found...
Definition: pzgeoel.cpp:165
void Read(TPZStream &buf, void *context) override
read objects from the stream
Definition: pzgmesh.cpp:1487
Contains TPZShapeLinear class which implements the shape functions of a linear one-dimensional elemen...
Contains the TPZRefPyramid class which implements the uniform refinement of a geometric hexahedral el...
int64_t NodeIndex(TPZGeoNode *nod)
Returns the index of the given node into the fNodeVec.
Definition: pzgmesh.cpp:1185
TPZGeoNode * FindNode(TPZVec< REAL > &co)
Returns the nearest node to the coordinate. This method is VERY INEFFICIENT.
Definition: pzgmesh.cpp:419
void RestoreReference(TPZCompMesh *cmesh)
Restore all reference in elements from computational mesh criated from current geometrical mesh previ...
Definition: pzgmesh.cpp:209
Contains declaration of TPZGeoElMapped class which implements a geometric element using its ancestral...
Contains TPZShapeCube class which implements the shape functions of a hexaedral element.
virtual TPZGeoEl * SubElement(int is) const =0
Returns a pointer to the subelement is.
TPZGeoEl * LowestFather()
Returns a pointer to the higher level father.
Definition: pzgeoel.h:337
void ClearInterfaceMaterialsMap()
Delete all interface materials in map.
Definition: pzgmesh.cpp:1568
int64_t size() const
Returns the number of elements of the vector.
Definition: pzvec.h:196
virtual void Resize(const int64_t newsize, const T &object)
Resizes the vector object reallocating the necessary storage, copying the existing objects to the new...
Definition: pzvec.h:373
Contains the TPZRefPrism class which implements the uniform refinement of a geometric prism element...
virtual void GetAllSiblings(TPZStack< TPZGeoEl *> &unrefinedSons)
[deprecated] use YoungestChildren
Definition: pzgeoel.cpp:410
void Resize(const int newsize)
Increase the size of the chunk vector.
Definition: pzadmchunk.h:280
Contains declaration of TPZMesh class which defines a geometrical mesh and contains a corresponding l...
void Push(const T object)
Pushes a copy of the object on the stack.
Definition: pzstack.h:80
std::string fName
TPZGeoMesh name for model identification.
Definition: pzgmesh.h:52
void BuildConnectivityOld()
Alternative method for computing the connectivity.
Definition: pzgmesh.cpp:1056
virtual void Write(const bool val)
Definition: TPZStream.cpp:8
virtual void Print(std::ostream &out=std::cout) const
Print the information of the grid to an ostream.
Definition: pzgmesh.cpp:130
Implements a generic geometric element which is refined according to a generic refinement pattern...
Definition: pzgmesh.h:35
Defines the behaviour of all geometric elements. GeometryTPZGeoEl is the common denominator for all g...
Definition: pzgeoel.h:43
Contains the TPZGeoBlend class which implements a blending map from curved boundaries to the interior...
TPZGeoEl * Element(int64_t iel)
Definition: pzgmesh.h:139
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 Index() const
Returns the index of the element within the element vector of the mesh.
Definition: pzgeoel.h:730
Contains the TPZGeoTetrahedra class which implements the geometry of a tetrahedral element...
#define DebugStop()
Returns a message to user put a breakpoint in.
Definition: pzerror.h:20
int64_t NNodes() const
Number of nodes of the mesh.
Definition: pzgmesh.h:126
void CenterX(TPZVec< REAL > &Xcenter) const
return the coordinates of the center of the side in real space
#define LOGPZ_DEBUG(A, B)
Define log for debug info.
Definition: pzlog.h:87
void X(TPZVec< REAL > &loc, TPZVec< REAL > &result) const
X coordinate of a point loc of the side.
int64_t NFreeElements() const
Access method to return the number of free elements.
Definition: pzadmchunk.h:87
TPZAdmChunkVector< TPZGeoNode > & NodeVec()
Definition: pzgmesh.h:140
void BuildBlendConnectivity()
Set connectivity information elements with blend geometric map.
Definition: pzgeoel.cpp:1959
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
TPZCompMesh * fReference
Computational mesh associated.
Definition: pzgmesh.h:55
string res
Definition: test.py:151
TPZGeoEl * FindElement(TPZVec< REAL > &x, TPZVec< REAL > &qsi, int64_t &InitialElIndex, int targetDim) const
Returns the element that contains the given point x and it respective point in parametric domain qsi...
Definition: pzgmesh.cpp:579
Contains declaration of TPZGeoElRefPattern class which implements a generic geometric element which i...
virtual void PrintTopologicalInfo(std::ostream &out=std::cout) const
Definition: pzgmesh.cpp:169
void SetConnectivity(const TPZGeoElSide &neighbour) const
Groups all classes which model the h refinement These classes are used as template arguments of...
Definition: pzrefpoint.cpp:15
REAL co[8][3]
Coordinates of the eight nodes.
TPZGeoMesh()
Constructors and destructor.
Definition: pzgmesh.cpp:47
Contains TPZShapePoint class which implements the shape function associated with a point...
virtual TPZGeoElSide Neighbour(int side)=0
Returns a pointer to the neighbour and the neighbourside along side of the current element...
#define LOGPZ_ERROR(A, B)
Define log for errors (cout)
Definition: pzlog.h:93
Contains declaration of TPZCompMesh class which is a repository for computational elements...
Contains the TPZRefPattern class which defines the topology of the current refinement pattern to a me...
TPZGeoMesh & operator=(const TPZGeoMesh &cp)
Operator of copy.
Definition: pzgmesh.cpp:62
void CompactDataStructure(CompactScheme type=CompactScheme::ALWAYS)
Sets the method to compact the data structure based on the.
Definition: pzadmchunk.h:222
Contains the TPZGeoCube class which implements the geometry of hexahedra element. ...
virtual int NSideNodes(int side) const =0
Returns the number of nodes for a particular side.
int32_t Hash(std::string str)
Definition: TPZHash.cpp:10
Contains the TPZGeoPoint class which implements the geometry of a point element or 0-D element...
Contains TPZMatrix<TVar>class, root matrix class.
virtual int HasSubElement() const =0
Return 1 if the element has subelements.
void ComputeNeighbours(TPZStack< TPZGeoElSide > &compneigh)
Returns the set of neighbours as computed by the intersection of neighbours along vertices...
virtual int NNodes() const =0
Returns the number of nodes of the element.
Contains TPZShapePrism class which implements the shape functions of a prism element.
TPZGeoEl * Element() const
Definition: pzgeoelside.h:162
int NeighbourExists(int side, const TPZGeoElSide &gel)
Returns 1 if gel is a neighbour of the element along side.
Definition: pzgeoel.cpp:226
Contains the TPZRefTriangle class which implements the uniform refinement of a geometric triangular e...
void GetBoundaryElements(int64_t IndexNodeFrom, int64_t IndexNodeTo, TPZStack< TPZGeoEl *> &ElementVec, TPZStack< int > &Sides)
GetBoundaryElements returns all elements beweeen NodFrom and NodTo counterclock wise this method uses...
Definition: pzgmesh.cpp:238
void BuildConnectivity()
Build the connectivity of the grid.
Definition: pzgmesh.cpp:967
virtual int SideIsUndefined(int side)=0
Returns 1 if the side has not been defined by buildconnectivity.
virtual void SetSideDefined(int side)=0
Flags the side as defined, this means no neighbouring element was found.
Contains the TPZRefTetrahedra class which implements the uniform refinement of a geometric tetrahedra...
void DeleteElement(TPZGeoEl *gel, int64_t index=-1)
Centralized method to delete elements.
Definition: pzgmesh.cpp:1456
T Pop()
Retrieve an object from the stack.
Definition: pzstack.h:87
Contains the TPZRefPatternDataBase class which defines data base of patterns.
virtual int Dimension() const =0
Returns the dimension of the element.
REAL dist(TPZVec< T1 > &vec1, TPZVec< T1 > &vec2)
Definition: pzvec_extras.h:124
int AddInterfaceMaterial(int leftmaterial, int rightmaterial, int interfacematerial)
Add an interface material associated to left and right element materials.
Definition: pzgmesh.cpp:1534
int64_t Id() const
Returns the Id of the element.
Definition: pzgeoel.h:223
int Dimension()
Get Dimension.
Definition: pzgmesh.h:291
bool ComputeXInverse(TPZVec< REAL > &XD, TPZVec< REAL > &ksi, REAL Tol)
Computes the XInverse and returns true if ksi belongs to master element domain.
Definition: pzgeoel.cpp:686
TPZGeoEl * Reference() const
Return a pointer to the corresponding geometric element if such exists, return 0 otherwise.
Definition: pzcompel.cpp:1137
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...
virtual REAL SideArea(int side)
Returns the area from the face.
Definition: pzgeoel.cpp:1064
This class implements a geometric mesh for the pz environment. Geometry.
Definition: pzgmesh.h:48
MElementType
Define the element types.
Definition: pzeltype.h:52
virtual int ProjectInParametricDomain(TPZVec< REAL > &qsi, TPZVec< REAL > &qsiInDomain)=0
Ortogonal projection from given qsi to a qsiInDomain (all in the element parametric domain) ...
This class implements a stack object. Utility.
Definition: pzcheckmesh.h:14
int Dimension() const
the dimension associated with the element/side
void ResetReference()
Reset the element referenced by the geometric element to NULL.
Definition: pzgeoel.h:431
Implements computational mesh. Computational Mesh.
Definition: pzcmesh.h:47
TPZAdmChunkVector< TPZCompEl * > & ElementVec()
Returns a reference to the element pointers vector.
Definition: pzcmesh.h:183
void RemoveConnectivities()
It removes the connectivities of the element.
Definition: pzgeoel.cpp:1098
void SetReference(TPZCompEl *elp)
Make the current element reference to the computational element.
Definition: pzgeoel.h:746
int Side() const
Definition: pzgeoelside.h:169
virtual TPZGeoEl * CreateGeoElementMapped(MElementType type, TPZVec< int64_t > &nodeindexes, int matid, int64_t &index)
Generic method for creating a geometric element. Putting this method centrally facilitates the modifi...
Definition: pzgmesh.cpp:1229
REAL Area()
Compute the area of the domain.
Definition: pzgmesh.cpp:1604
Contains the TPZRefLinear class which implements the uniform refinement of a geometric linear element...
void ResetConnectivities()
Reset all connectivities.
Definition: pzgmesh.cpp:1576
void Fill(const T &copy, const int64_t from=0, const int64_t numelem=-1)
Will fill the elements of the vector with a copy object.
Definition: pzvec.h:460
virtual ~TPZGeoMesh()
Destructor.
Definition: pzgmesh.cpp:97
int Exists() const
Definition: pzgeoelside.h:175
Contains TPZShapePiram class which implements the shape functions of a pyramid element.
Contains declaration of the abstract TPZStream class. TPZStream defines the interface for saving and ...
Contains the TPZGeoPyramid class which implements the geometry of pyramid element.
int InterfaceMaterial(int leftmaterial, int rightmaterial)
Returns the interface material associated to left and right element materials. If no interface materi...
Definition: pzgmesh.cpp:1548
bool VerifyNodeCoordinates(REAL tol=1e-1)
Verify coordinate of element nodes checking if they are coincident to the X mapping of the corner nod...
Definition: pzgeoel.cpp:1722
Defines the interface for saving and reading data. Persistency.
Definition: TPZStream.h:50
int NeighbourExists(const TPZGeoElSide &neighbour) const
Returns 1 if neighbour is a neighbour of the element along side.
int64_t NElements() const
Returns the number of elements of the vector.
Definition: pzvec.h:190
#define GMESHNOMATERIAL
Identifier indicating the no material is associated.
Definition: pzgmesh.h:24
TPZGeoEl * FindApproxElement(TPZVec< REAL > &x, TPZVec< REAL > &qsi, int64_t &InitialElIndex, int targetDim) const
find an element/parameter close to the point
Definition: pzgmesh.cpp:747
Contains the TPZGeoPrism class which implements the geometry of a prism element.
TPZGeoEl * FindCloseElement(TPZVec< REAL > &x, int64_t &InitialElIndex, int targetDim) const
Returns the element that is close to the given point x.
Definition: pzgmesh.cpp:472
Contains declaration of TPZGeoElement class which implements a generic geometric element with a unifo...
Extra utilities for TPZVec. Implementations of the saxpy, sscal, sdot, intercept, max and min functio...
clarg::argString m("-m", "input matrix file name (text format)", "matrix.txt")
void GetCoordinates(TPZVec< REAL > &co)
Fill the coordinates of the node.
Definition: pzgnode.cpp:83
TPZAdmChunkVector< TPZGeoNode > fNodeVec
List of nodes.
Definition: pzgmesh.h:61
void Read(TPZStream &buf, void *context) override
read objects from the stream
Definition: pzadmchunk.h:118
Groups all classes which model the geometry.
Definition: pzgeopoint.cpp:18
REAL ZeroTolerance()
Returns the tolerance to Zero value. Actually: .
Definition: pzreal.h:633
void ResetReference()
Resets all load references in elements and nodes.
Definition: pzgmesh.cpp:197
This class defines the interface to save and restore objects from TPZStream objects. Persistency.
Definition: TPZSavable.h:67
Defines the interface of a computational element. Computational Element.
Definition: pzcompel.h:59
virtual int NSubElements() const =0
Returns the number of subelements of the element independent of the fact whether the element has alr...
Contains TPZShapeQuad class which implements the shape functions of a quadrilateral element...
int64_t ElementIndex(TPZGeoEl *gel)
Returns the index of the given element into the fElementVec.
Definition: pzgmesh.cpp:1169
virtual void ResetSubElements()=0
Reset all subelements to NULL.
Contains the TPZRefCube class which implements the uniform refinement of a geometric hexahedral eleme...
Contains the TPZRefPoint class which implements the uniform refinement of a geometric point element...
virtual int ProjectBissectionInParametricDomain(TPZVec< REAL > &qsi, TPZVec< REAL > &qsiInDomain)=0
Projection from given qsi to a qsiInDomain (in the element boundary) using bissection method from giv...
TPZGeoEl * FindElementCaju(TPZVec< REAL > &x, TPZVec< REAL > &qsi, int64_t &InitialElIndex, int targetDim)
Returns the element that contains the given point x and it respective point in parametric domain qsi...
Definition: pzgmesh.cpp:603
static void WritePointer(const TPZSavable *obj, TPZStream *stream)
Contains TPZShapeTriang class which implements the shape functions of a triangular element...
void SetName(const std::string &nm)
Definition: pzgmesh.cpp:125
Implements an interface to register a class id and a restore function. Persistence.
Definition: TPZSavable.h:150
#define PZError
Defines the output device to error messages and the DebugStop() function.
Definition: pzerror.h:15
virtual void Read(bool &val)
Definition: TPZStream.cpp:91
TPZAdmChunkVector< TPZGeoEl * > & ElementVec()
Methods for handling pzlists.
Definition: pzgmesh.h:138
void Write(TPZStream &buf, int withclassid) const override
Writes this object to the TPZStream buffer. Include the classid if withclassid = true.
Definition: pzadmchunk.h:132