NeoPZ
tpzquadratictrig.cpp
Go to the documentation of this file.
1 
5 #include "tpzquadratictrig.h"
6 
7 #include "pzfmatrix.h"
8 #include "pzvec.h"
9 #include "pzgeotriangle.h"
10 #include "pzgmesh.h"
11 #include "pzgeoel.h"
12 #include "tpztriangle.h"
13 #include "pznoderep.h"
14 #include "pzshapetriang.h"
15 #include "tpzgeoelmapped.h"
16 
17 #include <iostream>
18 #include <iostream>
19 #include <cmath>
20 
21 #include "pzgeoelrefless.h.h"
22 #include "tpzgeoelrefpattern.h.h"
23 #include "pznoderep.h.h"
24 
25 #include "tpzgeomid.h"
26 
27 #ifdef _AUTODIFF
28 #include "fad.h"
29 #endif
30 
31 using namespace std;
32 using namespace pzshape;
33 using namespace pzgeom;
34 using namespace pztopology;
35 
42 template<class T>
43 void TPZQuadraticTrig::TShape(const TPZVec<T> &param,TPZFMatrix<T> &phi,TPZFMatrix<T> &dphi)
44 {
45  T qsi = param[0], eta = param[1];
46 
47  phi(0,0) = (qsi+eta-1.)*(2.*qsi+2.*eta-1.);
48  phi(1,0) = qsi*(2.*qsi-1.);
49  phi(2,0) = eta*(2.*eta-1.);
50  phi(3,0) = -4.*qsi*(qsi+eta-1.);
51  phi(4,0) = 4.*qsi*eta;
52  phi(5,0) = -4.*eta*(qsi+eta-1.);
53  dphi(0,0) = 1.-4.*(1.-qsi-eta);
54  dphi(0,1) = -1.+4.*qsi;
55  dphi(0,2) = 0.;
56  dphi(0,3) = 4.*(1.-qsi-eta)-4.*qsi;
57  dphi(0,4) = 4.*eta;
58  dphi(0,5) = -4.*eta;
59  dphi(1,0) = 1.-4.*(1.-qsi-eta);
60  dphi(1,1) = 0.;
61  dphi(1,2) = -1.+4.*eta;
62  dphi(1,3) = -4.*qsi;
63  dphi(1,4) = 4.*qsi;
64  dphi(1,5) = 4.*(1.-qsi-eta)-4.*eta;
65 }
66 
67 template<class T>
68 void TPZQuadraticTrig::X(const TPZFMatrix<REAL> &coord, TPZVec<T>& par, TPZVec< T >& result)
69 {
70  TPZManVector<T,3> parMap(2);
71  REAL spacephi[6],spacedphi[12];
72  TPZFNMatrix<6,T> phi(NNodes,1);
73  TPZFNMatrix<12,T> dphi(2,NNodes);
74  TShape(par,phi,dphi);
75  for(int i = 0; i < 3; i++)
76  {
77  result[i] = 0.0;
78  for(int j = 0; j < 6; j++) result[i] += phi(j,0)*coord.GetVal(i,j);
79  }
80 }
81 
82 
83 template<class T>
84 void TPZQuadraticTrig::GradX(const TPZFMatrix<REAL> &nodes,TPZVec<T> &loc, TPZFMatrix<T> &gradx){
85 
86  gradx.Resize(3,2);
87  gradx.Zero();
88  int nrow = nodes.Rows();
89  int ncol = nodes.Cols();
90 #ifdef PZDEBUG
91  if(nrow != 3 || ncol != 6){
92  std::cout << "Objects of incompatible lengths, gradient cannot be computed." << std::endl;
93  std::cout << "nodes matrix must be 3x6." << std::endl;
94  DebugStop();
95  }
96 
97 #endif
98 
99  TPZFNMatrix<3,T> phi(NNodes,1);
100  TPZFNMatrix<6,T> dphi(2,NNodes);
101  TShape(loc,phi,dphi);
102  for(int i = 0; i < NNodes; i++)
103  {
104  for(int j = 0; j < 3; j++)
105  {
106  gradx(j,0) += nodes.GetVal(j,i)*dphi(0,i);
107  gradx(j,1) += nodes.GetVal(j,i)*dphi(1,i);
108  }
109  }
110 
111 }
112 
113 // TPZGeoEl *TPZQuadraticTrig::CreateBCGeoEl(TPZGeoEl *orig,int side,int bc)
114 // {
115 // if(side==6)
116 // {
117 // TPZManVector<int64_t> nodes(3); int i;
118 // for (i=0;i<3;i++) nodes[i] = orig->SideNodeIndex(side,i);
119 // int64_t index;
120 // TPZGeoEl *gel = orig->Mesh()->CreateGeoBlendElement(ETriangle,nodes,bc,index);
121 // int iside;
122 // for (iside = 0; iside <6; iside++)
123 // {
124 // TPZGeoElSide(gel,iside).SetConnectivity(TPZGeoElSide(orig,TPZShapeTriang::ContainedSideLocId(side,iside)));
125 // }
126 // TPZGeoElSide(gel,6).SetConnectivity(TPZGeoElSide(orig,side));
127 // return gel;
128 // }
129 
130 // else if(side>-1 && side<3)
131 // {
132 // TPZManVector<int64_t> nodeindexes(1);
133 // nodeindexes[0] = orig->SideNodeIndex(side,0);
134 // int64_t index;
135 // TPZGeoEl *gel = orig->Mesh()->CreateGeoBlendElement(EPoint,nodeindexes,bc,index);
136 // TPZGeoElSide(gel,0).SetConnectivity(TPZGeoElSide(orig,side));
137 // return gel;
138 // }
139 
140 // else if(side > 2 && side < 6)
141 // {
142 // TPZManVector<int64_t> nodes(2);
143 // nodes[0] = orig->SideNodeIndex(side,0);
144 // nodes[1] = orig->SideNodeIndex(side,1);
145 // int64_t index;
146 // TPZGeoEl *gel = orig->Mesh()->CreateGeoBlendElement(EOned,nodes,bc,index);
147 // TPZGeoElSide(gel,0).SetConnectivity(TPZGeoElSide(orig,TPZShapeTriang::ContainedSideLocId(side,0)));
148 // TPZGeoElSide(gel,1).SetConnectivity(TPZGeoElSide(orig,TPZShapeTriang::ContainedSideLocId(side,1)));
149 // TPZGeoElSide(gel,2).SetConnectivity(TPZGeoElSide(orig,side));
150 // gel->BuildBlendConnectivity();
151 // return gel;
152 // }
153 
154 // else PZError << "TPZGeoTriangle::CreateBCGeoEl has no bc.\n";
155 // return 0;
156 // }
157 
158 
159 // /**
160 // * Creates a geometric element according to the type of the father element
161 // */
162 // TPZGeoEl *TPZQuadraticTrig::CreateGeoElement(TPZGeoMesh &mesh, MElementType type,
163 // TPZVec<int64_t>& nodeindexes,
164 // int matid,
165 // int64_t& index)
166 // {
167 // return CreateGeoElementMapped(mesh,type,nodeindexes,matid,index);
168 // }
169 
171 /* @param gmesh mesh in which the element should be inserted
172  @param matid material id of the element
173  @param lowercorner (in/out) on input lower corner o the cube where the element should be created, on exit position of the next cube
174  @param size (in) size of space where the element should be created
175  */
176 #include "tpzchangeel.h"
177 
178 void TPZQuadraticTrig::InsertExampleElement(TPZGeoMesh &gmesh, int matid, TPZVec<REAL> &lowercorner, TPZVec<REAL> &size)
179 {
180  TPZManVector<REAL,3> co(3),shift(3),scale(3);
181  TPZManVector<int64_t,3> nodeindexes(3);
182  for (int i=0; i<3; i++) {
183  scale[i] = size[i]/3.;
184  shift[i] = 1./2.+lowercorner[i];
185  }
186 
187  for (int i=0; i<NCornerNodes; i++) {
188  ParametricDomainNodeCoord(i, co);
189  for (int j=0; j<co.size(); j++) {
190  co[j] = shift[j]+scale[j]*co[j]+(rand()*0.2/RAND_MAX)-0.1;
191  }
192  nodeindexes[i] = gmesh.NodeVec().AllocateNewElement();
193  gmesh.NodeVec()[nodeindexes[i]].Initialize(co, gmesh);
194  }
195  int64_t index;
196  gmesh.CreateGeoElement(ETriangle, nodeindexes, matid, index);
197  TPZGeoEl *gel = gmesh.Element(index);
198  int nsides = gel->NSides();
199  for (int is=0; is<nsides; is++) {
200  gel->SetSideDefined(is);
201  }
202  gel = TPZChangeEl::ChangeToQuadratic(&gmesh, index);
203  for (int node = gel->NCornerNodes(); node < gel->NNodes(); node++) {
204  TPZManVector<REAL,3> co(3);
205  gel->NodePtr(node)->GetCoordinates(co);
206  for (int i=0; i<3; i++) {
207  co[i] += (0.2*rand())/RAND_MAX - 0.1;
208  }
209  gel->NodePtr(node)->SetCoord(co);
210  }
211 }
212 
213 
214 //void TPZQuadraticTrig::ParametricDomainNodeCoord(int node, TPZVec<REAL> &nodeCoord)
215 //{
216 // if(node > this->NNodes)
217 // {
218 // DebugStop();
219 // }
220 // nodeCoord.Resize(Dimension, 0.);
221 // switch (node) {
222 // case (0):
223 // {
224 // nodeCoord[0] = 0.;
225 // nodeCoord[1] = 0.;
226 // break;
227 // }
228 // case (1):
229 // {
230 // nodeCoord[0] = 1.;
231 // nodeCoord[1] = 0.;
232 // break;
233 // }
234 // case (2):
235 // {
236 // nodeCoord[0] = 0.;
237 // nodeCoord[1] = 1.;
238 // break;
239 // }
240 // case (3):
241 // {
242 // nodeCoord[0] = 0.5;
243 // nodeCoord[1] = 0.0;
244 // break;
245 // }
246 // case (4):
247 // {
248 // nodeCoord[0] = 0.5;
249 // nodeCoord[1] = 0.5;
250 // break;
251 // }
252 // case (5):
253 // {
254 // nodeCoord[0] = 0.0;
255 // nodeCoord[1] = 0.5;
256 // break;
257 // }
258 // default:
259 // {
260 // DebugStop();
261 // break;
262 // }
263 // }
264 //}
265 
267 
268 int TPZQuadraticTrig::ClassId() const{
269  return Hash("TPZQuadraticTrig") ^ TPZNodeRep<6,pztopology::TPZTriangle>::ClassId() << 1;
270 }
271 
272 template class TPZRestoreClass< TPZGeoElRefPattern<TPZQuadraticTrig>>;
273 
274 /*@orlandini : I REALLY dont know why is this here, so I have commented the following lines.
275 If it breaks something, I am sorry.*/
276 //template class pzgeom::TPZNodeRep<6,TPZQuadraticTrig>;
277 
278 namespace pzgeom {
279  template void TPZQuadraticTrig::X(const TPZFMatrix<REAL>&, TPZVec<REAL>&, TPZVec<REAL>&);
280  template void TPZQuadraticTrig::GradX(const TPZFMatrix<REAL> &nodes,TPZVec<REAL> &loc, TPZFMatrix<REAL> &gradx);
281 
282 #ifdef _AUTODIFF
283  template void TPZQuadraticTrig::X(const TPZFMatrix<REAL>&, TPZVec<Fad<REAL> >&, TPZVec<Fad<REAL> >&);
284  template void TPZQuadraticTrig::GradX(const TPZFMatrix<REAL> &nodes,TPZVec<Fad<REAL> > &loc, TPZFMatrix<Fad<REAL> > &gradx);
285 #endif
286 
287 }
TPZAdmChunkVector::AllocateNewElement
int AllocateNewElement()
Makes more room for new elements.
Definition: pzadmchunk.h:184
TPZGeoNode::SetCoord
void SetCoord(const TPZVec< REAL > &x)
Sets all coordinates into the current node. It gets the dim values from x.
Definition: pzgnode.cpp:60
TPZGeoEl::NodePtr
TPZGeoNode * NodePtr(int i) const
Returns a pointer to the ith node of the element.
Definition: pzgeoel.cpp:2566
tpzchangeel.h
Contains the TPZChangeEl class. It is a special map.
TPZGeoMesh::CreateGeoElement
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
ETriangle
Definition: pzeltype.h:56
pzvec.h
Templated vector implementation.
pzshape
groups all classes dedicated to the computation of shape functions
Definition: pzshapeextend.h:16
TPZGeoEl::NCornerNodes
virtual int NCornerNodes() const =0
Returns the number of corner nodes of the element.
Fad
Definition: fad.h:54
TPZVec
This class implements a simple vector storage scheme for a templated class T. Utility.
Definition: pzgeopoint.h:19
TPZGeoEl::NSides
virtual int NSides() const =0
Returns the number of connectivities of the element.
tpztriangle.h
Contains the TPZTriangle class which defines the topology of a triangle.
tpzgeoelmapped.h
Contains declaration of TPZGeoElMapped class which implements a geometric element using its ancestral...
TPZFMatrix::Zero
int Zero() override
Makes Zero all the elements.
Definition: pzfmatrix.h:651
pzgmesh.h
Contains declaration of TPZMesh class which defines a geometrical mesh and contains a corresponding l...
pztopology
Groups all classes defining the structure of the master element.
Definition: PrismExtend.cpp:15
TPZChangeEl::ChangeToQuadratic
static TPZGeoEl * ChangeToQuadratic(TPZGeoMesh *Mesh, int64_t ElemIndex)
Turns an linear geoelement to quadratic.
TPZGeoEl
Defines the behaviour of all geometric elements. GeometryTPZGeoEl is the common denominator for all g...
Definition: pzgeoel.h:43
TPZGeoMesh::Element
TPZGeoEl * Element(int64_t iel)
Definition: pzgmesh.h:139
pzfmatrix.h
Contains TPZMatrixclass which implements full matrix (using column major representation).
DebugStop
#define DebugStop()
Returns a message to user put a breakpoint in.
Definition: pzerror.h:20
TPZGeoMesh::NodeVec
TPZAdmChunkVector< TPZGeoNode > & NodeVec()
Definition: pzgmesh.h:140
TPZMatrix::Rows
int64_t Rows() const
Returns number of rows.
Definition: pzmatrix.h:803
pznoderep.h
Contains the TPZNodeRep class which implements ... Clase intermediaria que guarda.
co
REAL co[8][3]
Coordinates of the eight nodes.
Definition: tpzgensubstruct.cpp:56
pzgeom::TPZNodeRep
Implements ... Geometry Topology.
Definition: pznoderep.h:40
tpzquadratictrig.h
Contains the TPZQuadraticTrig class which defines a triangular geometric element with quadratic map...
Hash
int32_t Hash(std::string str)
Definition: TPZHash.cpp:10
TPZGeoEl::NNodes
virtual int NNodes() const =0
Returns the number of nodes of the element.
TPZGeoEl::SetSideDefined
virtual void SetSideDefined(int side)=0
Flags the side as defined, this means no neighbouring element was found.
TPZGeoMesh
This class implements a geometric mesh for the pz environment. Geometry.
Definition: pzgmesh.h:48
pznoderep.h.h
Contains the implementation of the TPZNodeRep methods.
TPZMatrix::Cols
int64_t Cols() const
Returns number of cols.
Definition: pzmatrix.h:809
TPZFMatrix::Resize
int Resize(const int64_t newRows, const int64_t wCols) override
Redimension a matrix, but maintain your elements.
Definition: pzfmatrix.cpp:1016
tpzgeoelrefpattern.h.h
Contains the implementation of the TPZGeoElRefPattern methods.
TPZGeoNode::GetCoordinates
void GetCoordinates(TPZVec< REAL > &co)
Fill the coordinates of the node.
Definition: pzgnode.cpp:83
pzgeom
Groups all classes which model the geometry.
Definition: pzgeopoint.cpp:18
TPZFMatrix::GetVal
const TVar & GetVal(const int64_t row, const int64_t col) const override
Get values without bounds checking This method is faster than "Get" if DEBUG is defined.
Definition: pzfmatrix.h:566
pzshapetriang.h
Contains TPZShapeTriang class which implements the shape functions of a triangular element...
TPZRestoreClass
Implements an interface to register a class id and a restore function. Persistence.
Definition: TPZSavable.h:150
TPZFNMatrix
Non abstract class which implements full matrices with preallocated storage with (N+1) entries...
Definition: pzfmatrix.h:716
Generated on Sun Aug 16 2020 11:17:10 for NeoPZ by   doxygen 1.8.13