NeoPZ
pzmatorthotropic.cpp
Go to the documentation of this file.
1 
6 #include "pzmatorthotropic.h"
7 #include "pzelmat.h"
8 #include "pzbndcond.h"
9 #include "pzmatrix.h"
10 #include "pzfmatrix.h"
11 #include "pzerror.h"
12 #include "pzmanvector.h"
13 #include <math.h>
14 #include <cmath>
15 #include <fstream>
16 
17 using namespace std;
18 
19 TPZMatOrthotropic::TPZMatOrthotropic(int nummat,TPZFMatrix<STATE> naxes,STATE eppx,STATE eppy,
20  STATE eppz,STATE vxy,STATE vyz,STATE vzx,
21  STATE gxy,STATE gyz,STATE gzx) :
22 TPZMaterial(nummat),
23 fKXX(3,3,0.),fKYY(3,3,0.),fKZZ(3,3,0.),
24 fKXY(3,3,0.),fKYX(3,3,0.),fKXZ(3,3,0.),
25 fKZX(3,3,0.),fKYZ(3,3,0.),
26 fKZY(3,3,0.),fXf(3,1,0.)
27 {
29  //normaliza os eixos
30  Normalize(naxes);
31  fLocAxs = naxes;
32  fGxy = gxy;
33  fGyz = gyz;
34  fGzx = gzx;
35  fEppx = eppx;
36  fEppy = eppy;
37  fEppz = eppz;
38  fVxy = vxy;
39  fVyx = fEppy*fVxy/fEppx;
40  fVyz = vyz;
41  fVzy = fEppz*fVyz/fEppy;
42  fVzx = vzx;
43  fVxz = fEppx*fVzx/fEppz;
44  fNumNom = -1.0 + fVxy*fVyx + fVxz*fVzx + fVyz*fVzy + fVxy*fVyz*fVzx + fVyx*fVzy*fVxz;
45 
46  int i,j;
47  for(i=0; i<3; i++) for(j=0; j<3; j++) fKXX(i,j) = 0.;
48  fKXX(0,0) = -(fEppx - fEppx*fVyz*fVzy)/fNumNom;
49  fKXX(1,1) = fGxy;
50  fKXX(2,2) = fGzx;
51 
52  for(i=0; i<3; i++) for(j=0; j<3; j++) fKXY(i,j) = 0.;//XY
53  fKXY(0,1) = -(fEppx*fVyx + fEppx*fVyz*fVzx)/fNumNom;
54  fKXY(1,0) = fGxy;
55 
56  for(i=0; i<3; i++) for(j=0; j<3; j++) fKXZ(i,j) = 0.;//XZ
57  fKXZ(0,2) = -(fEppx*fVzx + fEppx*fVyx*fVzy)/fNumNom;
58  fKXZ(2,0) = fGzx;
59 
60  for(i=0; i<3; i++) for(j=0; j<3; j++) fKYX(i,j) = 0.;//YX
61  fKYX(1,0) = -(fEppy*fVxy + fEppy*fVxz*fVzy)/fNumNom;
62  fKYX(0,1) = fGxy;
63 
64  for(i=0; i<3; i++) for(j=0; j<3; j++) fKYY(i,j) = 0.;
65  fKYY(0,0) = fGxy;
66  fKYY(1,1) = -(fEppy - fEppy*fVxz*fVzx)/fNumNom;
67  fKYY(2,2) = fGyz;
68 
69  for(i=0; i<3; i++) for(j=0; j<3; j++) fKYZ(i,j) = 0.;//YZ
70  fKYZ(1,2) = -(fEppy*fVzy + fEppy*fVxy*fVzx)/fNumNom;
71  fKYZ(2,1) = fGyz;
72 
73  for(i=0; i<3; i++) for(j=0; j<3; j++) fKZX(i,j) = 0.;//ZX
74  fKZX(2,0) = -fEppz*(fVxz + fVxy*fVyz)/fNumNom;
75  fKZX(0,2) = fGzx;
76 
77  for(i=0; i<3; i++) for(j=0; j<3; j++) fKZY(i,j) = 0.;//ZY
78  fKZY(2,1) = -fEppz*(fVyz + fVxz*fVyx)/fNumNom;
79  fKZY(1,2) = fGyz;
80 
81  for(i=0; i<3; i++) for(j=0; j<3; j++) fKZZ(i,j) = 0.;
82  fKZZ(0,0) = fGzx;
83  fKZZ(1,1) = fGyz;
84  fKZZ(2,2) = -(fEppz - fEppz*fVxy*fVyx)/fNumNom;
85 }
86 
87 TPZMatOrthotropic::~TPZMatOrthotropic() {
88 }
89 
90 int TPZMatOrthotropic::NStateVariables() const {
91  return 3;
92 }
93 
94 void TPZMatOrthotropic::Print(std::ostream &out) {
95  out << "name of material : " << Name() << "\n";
96  out << "properties : \n";
97  out << "Eppx = " << fEppx << "\tEppy = " << fEppy << "\tEppz = " << fEppz << endl;
98  out << "vxy = " << fVxy << "\tvyx = " << fVyx << endl;
99  out << "vzx = " << fVzx << "\tvxz = " << fVxz << endl;
100  out << "vxy = " << fVxy << "\tvyx = " << fVyx << endl;
101  out << "gxy = " << fGxy << "\tgyz = " << fGyz << "\tgzx = " << fGzx << endl;
102  out << "numnomvar = " << fNumNom << endl;
103  out << "\n>>>>>>>>>>>>>> MATRIZES KIJ <<<<<<<<<<<<<<<<<<\n\n";
104  fKXX.Print("KXX",out);
105  fKXY.Print("KXY",out);
106  fKXZ.Print("KXZ",out);
107  fKYX.Print("KYX",out);
108  fKYY.Print("KYY",out);
109  fKYZ.Print("KYZ",out);
110  fKZX.Print("KZX",out);
111  fKZY.Print("KZY",out);
112  fKZZ.Print("KZZ",out);
113 
114  TPZMaterial::Print(out);
115 }
116 
118 ofstream MatrizesK("MatrizesK.out");
119 void TPZMatOrthotropic::Contribute(TPZMaterialData &data,
120  REAL weight,
121  TPZFMatrix<STATE> &ek,
122  TPZFMatrix<STATE> &ef) {
123 
124  TPZFMatrix<REAL> &dphi = data.dphix;
125  TPZFMatrix<REAL> &phi = data.phi;
126  TPZManVector<REAL,3> &x = data.x;
127  TPZFMatrix<STATE> axes(data.axes.Rows(),data.axes.Cols());
128  for (int r=0; r<axes.Rows(); r++) {
129  for (int c=0; c<axes.Cols(); c++) {
130  axes(r,c) = data.axes(r,c);
131  }
132  }
133 
134  int phr = phi.Rows();
135  if(fForcingFunction) {
136  TPZManVector<STATE> res(3);//,&fXf(0,0),3);// phi(in, 0) = phi_in
137  fForcingFunction->Execute(x,res); // dphi(i,j) = dphi_j/dxi
138  int i;
139  for(i=0; i<3; i++) fXf(i,0) = res[i];
140  }
141 
142  TPZFMatrix<STATE> Rt(fLocAxs),R(3,3,0.),newaxes(3,3,0.);
143  TPZFMatrix<STATE> kxx(3,3,0.),kxy(3,3,0.),kxz(3,3,0.);
144  TPZFMatrix<STATE> kyx(3,3,0.),kyy(3,3,0.),kyz(3,3,0.);
145  TPZFMatrix<STATE> kzx(3,3,0.),kzy(3,3,0.),kzz(3,3,0.);
146 
147  Rt.Transpose(&R);
153  axes.Transpose(&newaxes);
154  newaxes = Rt*newaxes;
155 
156  kxx = R*(fKXX*Rt);
157  kxy = R*(fKXY*Rt);
158  kxz = R*(fKXZ*Rt);
159  kyx = R*(fKYX*Rt);
160  kyy = R*(fKYY*Rt);
161  kyz = R*(fKYZ*Rt);
162  kzx = R*(fKZX*Rt);
163  kzy = R*(fKZY*Rt);
164  kzz = R*(fKZZ*Rt);
165 
166  MatrizesK << "\n>>>>>>>>>>>>>> MATRIZES KIJ <<<<<<<<<<<<<<<<<<\n\n";
167  kxx.Print("KXX",MatrizesK);
168  kxy.Print("KXY",MatrizesK);
169  kxz.Print("KXZ",MatrizesK);
170  kyx.Print("KYX",MatrizesK);
171  kyy.Print("KYY",MatrizesK);
172  kyz.Print("KYZ",MatrizesK);
173  kzx.Print("KZX",MatrizesK);
174  kzy.Print("KZY",MatrizesK);
175  kzz.Print("KZZ",MatrizesK);
176 
177  for( int in = 0; in < phr; in++ ) {
178  ef(3*in, 0) += weight * phi(in, 0) * fXf(0,0);
179  ef(3*in+1, 0) += weight * phi(in, 0) * fXf(1,0);
180  ef(3*in+2, 0) += weight * phi(in, 0) * fXf(2,0);
181  for( int jn = 0; jn < phr; jn++ ) {
182 
184  REAL dphin0 = newaxes(0,0)*dphi(0,in)+newaxes(0,1)*dphi(1,in)+newaxes(0,2)*dphi(2,in);
185  REAL dphin1 = newaxes(1,0)*dphi(0,in)+newaxes(1,1)*dphi(1,in)+newaxes(1,2)*dphi(2,in);
186  REAL dphin2 = newaxes(2,0)*dphi(0,in)+newaxes(2,1)*dphi(1,in)+newaxes(2,2)*dphi(2,in);
188  REAL dphjn0 = newaxes(0,0)*dphi(0,jn)+newaxes(0,1)*dphi(1,jn)+newaxes(0,2)*dphi(2,jn);
189  REAL dphjn1 = newaxes(1,0)*dphi(0,jn)+newaxes(1,1)*dphi(1,jn)+newaxes(1,2)*dphi(2,jn);
190  REAL dphjn2 = newaxes(2,0)*dphi(0,jn)+newaxes(2,1)*dphi(1,jn)+newaxes(2,2)*dphi(2,jn);
192  int k,l;
193  for(k=0; k<3; k++) {
194  for(l=0; l<3; l++) {
195  ek(3*in+k,3*jn+l) += weight * ( dphin0*dphjn0*kxx(k,l) +
196  dphin0*dphjn1*kxy(k,l) +
197  dphin0*dphjn2*kxz(k,l) +
198  dphin1*dphjn0*kyx(k,l) +
199  dphin1*dphjn1*kyy(k,l) +
200  dphin1*dphjn2*kyz(k,l) +
201  dphin2*dphjn0*kzx(k,l) +
202  dphin2*dphjn1*kzy(k,l) +
203  dphin2*dphjn2*kzz(k,l)
204  );
205  }
206  }
207  }
208  }
209 }
210 
211 void TPZMatOrthotropic::ContributeBC(TPZMaterialData &data,
212  REAL weight,
213  TPZFMatrix<STATE> &ek,
214  TPZFMatrix<STATE> &ef,
215  TPZBndCond &bc) {
216  TPZFMatrix<REAL> &phi = data.phi;
217 
218  const STATE BIGNUMBER = 1.e12;
219 
220  int phr = phi.Rows();
221  short in,jn,k,l;
222  STATE v2[3];
223  v2[0] = bc.Val2()(0,0);
224  v2[1] = bc.Val2()(1,0);
225  v2[2] = bc.Val2()(2,0);
226 
227  switch (bc.Type()) {
228  case 0 :// Dirichlet condition
229  for(in = 0 ; in < phr; in++) {
230  ef(3*in,0) += BIGNUMBER * v2[0] * phi(in,0) * weight;// forced v2[0] x displacement
231  ef(3*in+1,0) += BIGNUMBER * v2[1] * phi(in,0) * weight;// forced v2[1] y displacement
232  ef(3*in+2,0) += BIGNUMBER * v2[2] * phi(in,0) * weight;// forced v2[2] y displacement
233  for (jn = 0 ; jn < phi.Rows(); jn++) {
234  for(k=0; k<3; k++) {
235  ek(3*in+k,3*jn+k) += BIGNUMBER * phi(in,0) * phi(jn,0) * weight;
236  }//o Big. deve ser s�na diagonal: OK!
237  }
238  }
239  break;
240 
241  case 1 :// Neumann condition
242  for(in = 0 ; in < phi.Rows(); in++) {// componentes da tra�o normal ao contorno
243  ef(3*in,0) += v2[0] * phi(in,0) * weight; // tra�o em x (ou press�)
244  ef(3*in+1,0) += v2[1] * phi(in,0) * weight; // tra�o em y (ou press�)
245  ef(3*in+2,0) += v2[2] * phi(in,0) * weight; // tra�o em z (ou press�)
246  }
247  break;
248 
249  case 2 :// condi�o mista
250  for(in = 0 ; in < phi.Rows(); in++) { //Sigmaij
251  ef(3*in, 0) += v2[0] * phi(in, 0) * weight;// Neumann x
252  ef(3*in+1, 0) += v2[1] * phi(in, 0) * weight;// Neumann y
253  ef(3*in+2, 0) += v2[2] * phi(in, 0) * weight;// Neumann z
254  for (jn = 0 ; jn < phi.Rows(); jn++) {
255  for(k=0; k<3; k++) {
256  for(l=0; l<3; l++) {
257  ek(3*in+k,3*jn+l) += bc.Val1()(k,l) * phi(in,0) * phi(jn,0) * weight;//integral de contorno
258  }
259  }
260  }
261  }
262  }
263 }
264 
265 
267 int TPZMatOrthotropic::VariableIndex(const std::string &name){
268  if("DisplacX" == name) return 0;
269  if(!strcmp("DisplacX",name.c_str())) return 0;
270  if(!strcmp("DisplacY",name.c_str())) return 1;
271  if(!strcmp("DisplacZ",name.c_str())) return 2;
272  if(!strcmp("DisplacN0",name.c_str())) return 3;
273  if(!strcmp("DisplacN1",name.c_str())) return 4;
274  if(!strcmp("DisplacN2",name.c_str())) return 5;
275  if(!strcmp("GlobDisplac",name.c_str())) return 6;
276  if(!strcmp("FiberDisplac",name.c_str())) return 7;
277  if(!strcmp("Tension",name.c_str())) return 8;
278  if(!strcmp("Tensor",name.c_str())) return 9;
279  if(!strcmp("SigX",name.c_str())) return 10;
280  if(!strcmp("SigY",name.c_str())) return 11;
281  if(!strcmp("SigZ",name.c_str())) return 12;
282  if(!strcmp("TauXY",name.c_str())) return 13;
283  if(!strcmp("TauXZ",name.c_str())) return 14;
284  if(!strcmp("TauYZ",name.c_str())) return 15;
285 
286 
287  cout << "TPZMatOrthotropic::VariableIndex Error\n";
288  return -1;
289 }
290 
291 int TPZMatOrthotropic::NSolutionVariables(int var){
292 
293  if(var > -1 && var < 6) return 1;//escalares
294  if(var == 6 || var == 7) return 3;//vetores
295  if(var == 8) return 6;//deve ser
296  if(var == 9) return 9;//tensor completo (sim�rico)
297  if(var > 9 && var < 16) return 1;
298  cout << "TPZMatOrthotropic::NSolutionVariables Error\n";
299  return 0;
300 }
301 
302 void TPZMatOrthotropic::Solution(TPZVec<STATE> &Sol,TPZFMatrix<STATE> &DSol,TPZFMatrix<REAL> &axespar,int var,TPZVec<STATE> &Solout){ //OBS.:acrescentado ostream
303  TPZFMatrix<STATE> axes(axespar.Rows(),axespar.Cols());
304  for (int r=0; r<axes.Rows(); r++) {
305  for (int c=0; c<axes.Cols(); c++) {
306  axes(r,c) = axespar(r,c);
307  }
308  }
309 
310  if(var == 0){
311  Solout.Resize(1);
312  Solout[0] = Sol[0];//eixo global X
313  return;
314  }
315  if(var == 1){
316  Solout.Resize(1);
317  Solout[0] = Sol[1];//eixo global Y
318  return;
319  }
320  if(var == 2){
321  Solout.Resize(1);
322  Solout[0] = Sol[2];//eixo global Z
323  return;
324  }
325  if(var == 6){
326  Solout.Resize(3);
327  Solout[0] = Sol[0];//eixos globais
328  Solout[1] = Sol[1];// { X , Y , Z }
329  Solout[2] = Sol[2];//
330  return;
331  }
332  if(var==3 || var==4 || var==5 || var==7){
333  cout << "TPZMatOrthotropic::Solution not implemented\n";
334  return;
335  TPZFMatrix<STATE> SolN(3,1,0.),sol(3,1,0.);
336  sol(0,0) = Sol[0];
337  sol(1,0) = Sol[1];
338  sol(2,0) = Sol[2];
339  SolN = fLocAxs*sol;
340  if(var == 3){
341  Solout.Resize(1);
342  Solout[0] = SolN(0,0);//eixo local n0
343  return;
344  }
345  if(var == 4){
346  Solout.Resize(1);
347  Solout[0] = SolN(1,0);//eixo local n1
348  return;
349  }
350  if(var == 5){
351  Solout.Resize(1);
352  Solout[0] = SolN(2,0);//eixo local n2
353  return;
354  }
355  if(var == 7){
356  Solout.Resize(3);
357  Solout[0] = SolN(0,0);//eixos locais
358  Solout[1] = SolN(1,0);// { n0 , n1 , n2 }
359  Solout[2] = SolN(2,0);//
360  return;
361  }
362  }
363  if(var > 7 && var < 16) {
364  Solout.Resize(6);
365  TPZFMatrix<STATE> axest(3,3,0.),floc_axt(3,3,0.),grdUGdn(3,3,0.),grdUlocdn(3,3,0.);
366  axes.Transpose(&axest);
367  floc_axt = fLocAxs*axest;
368  grdUGdn = floc_axt*DSol;
369  TPZFMatrix<STATE> grdUGdnT(3,3,0.),grdULdn(3,3,0.),grdULdnT(3,3,0.);
370  grdUGdn.Transpose(&grdUGdnT);
371  grdULdn = fLocAxs*grdUGdnT;
372  grdULdn.Transpose(&grdULdnT);
373  TPZFMatrix<STATE> FibStrain(3,3,0.);
374  FibStrain = ((STATE)0.5)*(grdULdnT + grdULdn);
375 
376  STATE epsx = FibStrain(0,0);// du/dx
377  STATE epsy = FibStrain(1,1);// dv/dy
378  STATE epsz = FibStrain(2,2);// dw/dz
379  STATE epsxy = FibStrain(0,1);
380  STATE epsyz = FibStrain(1,2);
381  STATE epszx = FibStrain(2,0);
382 
383  STATE SigX = -fEppx*(epsx*(1.-fVyz*fVzy)+epsy*(fVyx+fVyz*fVzx)+epsz*(fVzx+fVyx*fVzy))/fNumNom;
384  STATE SigY = -fEppy*(epsy*(1.-fVxz*fVzx)+epsx*(fVxy+fVxz*fVzy)+epsz*(fVzy+fVxy*fVzx))/fNumNom;
385  STATE SigZ = -fEppz*(epsz*(1.-fVxy*fVyx)+epsx*(fVxz+fVxy*fVyz)+epsy*(fVyz+fVxz*fVyx))/fNumNom;
386 
387  STATE TauXY = 2.*fGxy*epsxy;
388  STATE TauYZ = 2.*fGyz*epsyz;
389  STATE TauZX = 2.*fGzx*epszx;
390 
391  TPZFMatrix<STATE> Tensor(3,3,0.);
392  Tensor(0,0) = SigX;
393  Tensor(1,1) = SigY;
394  Tensor(2,2) = SigZ;
395  Tensor(0,1) = Tensor(1,0) = TauXY;
396  Tensor(1,2) = Tensor(2,1) = TauYZ;
397  Tensor(2,0) = Tensor(0,2) = TauZX;
398 
399  TPZFMatrix<STATE> locaxsT(3,3,0.);
400  fLocAxs.Transpose(&locaxsT);
401  TPZFMatrix<STATE> SigGlob = locaxsT*(Tensor*fLocAxs);
402 
403  if(var == 10){
404  Solout.Resize(1);
405  Solout[0] = SigX;
406  return;
407  }
408  if(var == 11){
409  Solout.Resize(1);
410  Solout[0] = SigY;
411  return;
412  }
413  if(var == 12){
414  Solout.Resize(1);
415  Solout[0] = SigZ;
416  return;
417  }
418  if(var == 13){
419  Solout.Resize(1);
420  Solout[0] = TauXY;
421  return;
422  }
423  if(var == 14){
424  Solout.Resize(1);
425  Solout[0] = TauZX;
426  return;
427  }
428  if(var == 15){
429  Solout.Resize(1);
430  Solout[0] = TauYZ;
431  return;
432  }
433 
435  Solout[0] = SigGlob(0,0);//SigX
436  Solout[1] = SigGlob(1,1);//SigY
437  Solout[2] = SigGlob(2,2);//SigZ
438  Solout[3] = SigGlob(0,1);//TauXY
439  Solout[4] = SigGlob(1,2);//TauYZ
440  Solout[5] = SigGlob(2,0);//TauZX
441 
442  if(var == 8) return;
443 
444  //aqui var �9 e o tensor �sim�rico
445  Solout.Resize(9);
447  Solout[0] = SigGlob(0,0);//00
448  Solout[1] = SigGlob(0,1);//01=10
449  Solout[2] = SigGlob(2,0);//02=20
450  Solout[3] = SigGlob(0,1);//10=01
451  Solout[4] = SigGlob(1,1);//11
452  Solout[5] = SigGlob(1,2);//12=21
453  Solout[6] = SigGlob(2,0);//20=02
454  Solout[7] = SigGlob(1,2);//21=12
455  Solout[8] = SigGlob(2,2);//22
456  return;
457  }
458 }
459 
460 void TPZMatOrthotropic::Flux(TPZVec<REAL> &/*x*/, TPZVec<STATE> &/*Sol*/, TPZFMatrix<STATE> &/*DSol*/, TPZFMatrix<REAL> &/*axes*/, TPZVec<STATE> &/*flux*/) {
461  //Flux(TPZVec<REAL> &x, TPZVec<STATE> &Sol, TPZFMatrix<STATE> &DSol, TPZFMatrix<REAL> &axes, TPZVec<STATE> &flux)
462 }
463 
464 void TPZMatOrthotropic::Errors(TPZVec<REAL> &/*x*/,TPZVec<STATE> &u,
465  TPZFMatrix<STATE> &dudx, TPZFMatrix<REAL> &axes, TPZVec<STATE> &/*flux*/,
466  TPZVec<STATE> &u_exact,TPZFMatrix<STATE> &du_exact,TPZVec<REAL> &values) {
467 
468  //TPZVec<REAL> sol(1),dsol(3);
469  TPZManVector<STATE> sol(1),dsol(3);
470  Solution(u,dudx,axes,1,sol);
471  Solution(u,dudx,axes,2,dsol);
472  if(dudx.Rows()<3) {
473  STATE dx = du_exact(0,0)*axes(0,0)+du_exact(1,0)*axes(0,1);
474  STATE dy = du_exact(0,0)*axes(1,0)+du_exact(1,0)*axes(1,1);
475  STATE parc1 = fabs(dx-dudx(0,0));
476  STATE parc2 = fabs(dy-dudx(1,0));
477  //Norma L2
478  values[1] = pow(fabs(u[0] - u_exact[0]),(STATE)2.0);
479  //seminorma
480  values[2] = pow(parc1,(STATE)2.)+pow(parc2,(STATE)2.);
481  //Norma Energia
482  values[0] = values[1]+values[2];
483  return;
484  }
485  //values[1] : eror em norma L2
486  values[1] = pow(sol[0] - u_exact[0],(STATE)2.0);
487  //values[2] : erro em semi norma H1
488  values[2] = pow(dsol[0] - du_exact(0,0),(STATE)2.0);
489  if(dudx.Rows()>1) values[2] += pow(dsol[1] - du_exact(1,0),(STATE)2.0);
490  if(dudx.Rows()>2) values[2] += pow(dsol[2] - du_exact(2,0),(STATE)2.0);
491  //values[0] : erro em norma H1 <=> norma Energia
492  values[0] = values[1]+values[2];
493 }
494 
495 void TPZMatOrthotropic::Normalize(TPZFMatrix<STATE> &naxes) {
501  int i;
503  STATE norm2 = naxes(0,0)*naxes(0,0)+naxes(0,1)*naxes(0,1)+naxes(0,2)*naxes(0,2);
504  if(norm2 == 0.){PZError << "TPZMatOrthotropic::Normalize: Eixo nulo nao e valido (eixo 1)\n"; exit(-1);}
505  for(i=0;i<3;i++) naxes(0,i) /= sqrt(norm2);
506  norm2 = naxes(1,0)*naxes(1,0)+naxes(1,1)*naxes(1,1)+naxes(1,2)*naxes(1,2);
507  if(norm2 == 0.){PZError << "TPZMatOrthotropic::Normalize: Eixo nulo nao e valido (eixo 2)\n"; exit(-1);}
508  for(i=0;i<3;i++) naxes(1,i) /= sqrt(norm2);
509  norm2 = naxes(2,0)*naxes(2,0)+naxes(2,1)*naxes(2,1)+naxes(2,2)*naxes(2,2);
510  if(norm2 != 0.) for(i=0;i<3;i++) naxes(2,i) /= sqrt(norm2);
512  STATE componente_K = naxes(0,0)*naxes(1,1) - naxes(0,1)*naxes(1,0);
513  STATE componente_J = -naxes(0,0)*naxes(1,2) + naxes(0,2)*naxes(1,0);
514  STATE componente_I = naxes(0,1)*naxes(1,2) - naxes(0,2)*naxes(1,1);
516  if(componente_I == 0. && componente_J == 0. && componente_K == 0.){
517  PZError << "TPZMatOrthotropic::Normalize: Os dois primeiros eixos nao devem ser paralelos\n";
518  PZError << "Programa abortado\n";
519  exit(-1);
520  }
523  norm2 = componente_I*componente_I+componente_J*componente_J+componente_K*componente_K;
524  norm2 = sqrt(norm2);
525  naxes(2,0) = componente_I/norm2;
526  naxes(2,1) = componente_J/norm2;
527  naxes(2,2) = componente_K/norm2;
529  naxes(1,2) = naxes(2,0)*naxes(0,1) - naxes(2,1)*naxes(0,0);//K
530  naxes(1,1) = -naxes(2,0)*naxes(0,2) + naxes(2,2)*naxes(0,0);//J
531  naxes(1,0) = naxes(2,1)*naxes(0,2) - naxes(2,2)*naxes(0,1);//I
532  naxes.Print(" * * * Eixos das fibras * * *");
533 }
TPZFunction::Execute
virtual void Execute(const TPZVec< REAL > &x, TPZVec< TVar > &f, TPZFMatrix< TVar > &df)
Performs function computation.
Definition: pzfunction.h:38
TPZMatOrthotropic::NStateVariables
virtual int NStateVariables() const override
Returns the number of state variables associated with the material.
Definition: pzmatorthotropic.cpp:90
fabs
expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ fabs
Definition: tfadfunc.h:140
TPZMatOrthotropic::fKYZ
TPZFMatrix< STATE > fKYZ
Definition: pzmatorthotropic.h:21
TPZMatOrthotropic::Print
virtual void Print(std::ostream &out) override
Prints out the data associated with the material.
Definition: pzmatorthotropic.cpp:94
TPZMaterialData::x
TPZManVector< REAL, 3 > x
value of the coordinate at the integration point
Definition: pzmaterialdata.h:71
TPZMatOrthotropic::Name
virtual std::string Name() override
Returns the name of the material.
Definition: pzmatorthotropic.h:45
bc
clarg::argBool bc("-bc", "binary checkpoints", false)
TPZMatOrthotropic::Normalize
void Normalize(TPZFMatrix< STATE > &naxes)
Verifies the consistency of the axles.
Definition: pzmatorthotropic.cpp:495
TPZMatOrthotropic::NSolutionVariables
virtual int NSolutionVariables(int var) override
Returns the number of variables associated with the variable indexed by var.
Definition: pzmatorthotropic.cpp:291
TPZBndCond::Type
int Type()
Definition: pzbndcond.h:249
TPZMatOrthotropic::TPZMatOrthotropic
TPZMatOrthotropic(int nummat, TPZFMatrix< STATE > naxes, STATE eppx, STATE eppy, STATE eppz, STATE vxy, STATE vyz, STATE vzx, STATE gxy, STATE gyz, STATE gzx)
Definition: pzmatorthotropic.cpp:19
pzerror.h
Defines PZError.
TPZMatOrthotropic::fLocAxs
TPZFMatrix< STATE > fLocAxs
Definition: pzmatorthotropic.h:22
pzmatorthotropic.h
Contains the TPZMatOrthotropic class.
TPZBndCond::Val2
TPZFMatrix< STATE > & Val2(int loadcase=0)
Definition: pzbndcond.h:255
TPZMaterialData::phi
TPZFNMatrix< 220, REAL > phi
vector of shapefunctions (format is dependent on the value of shapetype)
Definition: pzmaterialdata.h:55
TPZMatOrthotropic::fKZZ
TPZFMatrix< STATE > fKZZ
Definition: pzmatorthotropic.h:21
TPZMatOrthotropic::fKXY
TPZFMatrix< STATE > fKXY
Definition: pzmatorthotropic.h:21
TPZMaterialData::dphix
TPZFNMatrix< 660, REAL > dphix
values of the derivative of the shape functions
Definition: pzmaterialdata.h:59
TPZMaterial::Print
virtual void Print(std::ostream &out=std::cout)
Prints out the data associated with the material.
Definition: TPZMaterial.cpp:103
TPZMaterial
This abstract class defines the behaviour which each derived class needs to implement.
Definition: TPZMaterial.h:39
pzelmat.h
Contains declaration of TPZElementMatrix struct which associates an element matrix with the coeficien...
TPZMatOrthotropic::fKYY
TPZFMatrix< STATE > fKYY
Definition: pzmatorthotropic.h:21
pzbndcond.h
Contains the TPZBndCond class which implements a boundary condition for TPZMaterial objects...
TPZVec::Resize
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
TPZMatOrthotropic::~TPZMatOrthotropic
virtual ~TPZMatOrthotropic()
Definition: pzmatorthotropic.cpp:87
dx
expr_ dx(i) *cos(expr_.val())
TPZMatOrthotropic::fKZX
TPZFMatrix< STATE > fKZX
Definition: pzmatorthotropic.h:21
pzfmatrix.h
Contains TPZMatrixclass which implements full matrix (using column major representation).
TPZBndCond
This class defines the boundary condition for TPZMaterial objects.
Definition: pzbndcond.h:29
TPZMatOrthotropic::fKYX
TPZFMatrix< STATE > fKYX
Definition: pzmatorthotropic.h:21
pzmanvector.h
Free store vector implementation.
TPZMatrix::Rows
int64_t Rows() const
Returns number of rows.
Definition: pzmatrix.h:803
TPZMatOrthotropic::fKZY
TPZFMatrix< STATE > fKZY
Definition: pzmatorthotropic.h:21
sqrt
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
TPZMaterialData::axes
TPZFNMatrix< 9, REAL > axes
axes indicating the directions of the derivatives of the shapefunctions
Definition: pzmaterialdata.h:63
MatrizesK
ofstream MatrizesK("MatrizesK.out")
Output file to write stiffness matrix.
test.res
string res
Definition: test.py:151
TPZBndCond::Val1
TPZFMatrix< STATE > & Val1()
Definition: pzbndcond.h:253
pzmatrix.h
Contains TPZMatrix<TVar>class, root matrix class.
TPZMatOrthotropic::VariableIndex
virtual int VariableIndex(const std::string &name) override
Definition: pzmatorthotropic.cpp:267
pow
TPZFlopCounter pow(const TPZFlopCounter &orig, const TPZFlopCounter &xp)
Returns the power and increments the counter of the power.
Definition: pzreal.h:487
TPZMatOrthotropic::ContributeBC
virtual void ContributeBC(TPZMaterialData &data, REAL weight, TPZFMatrix< STATE > &ek, TPZFMatrix< STATE > &ef, TPZBndCond &bc) override
It computes a contribution to the stiffness matrix and load vector at one BC integration point...
Definition: pzmatorthotropic.cpp:211
TPZMaterial::fForcingFunction
TPZAutoPointer< TPZFunction< STATE > > fForcingFunction
Pointer to forcing function, it is the right member at differential equation.
Definition: TPZMaterial.h:47
TPZMatOrthotropic::Flux
virtual void Flux(TPZVec< REAL > &x, TPZVec< STATE > &Sol, TPZFMatrix< STATE > &DSol, TPZFMatrix< REAL > &axes, TPZVec< STATE > &flux) override
Computes the value of the flux function to be used by ZZ error estimator.
Definition: pzmatorthotropic.cpp:460
TPZMatrix::Cols
int64_t Cols() const
Returns number of cols.
Definition: pzmatrix.h:809
TPZMatrix::Print
virtual void Print(std::ostream &out) const
Definition: pzmatrix.h:253
rdt.values
def values
Definition: rdt.py:119
TPZMatOrthotropic::Contribute
virtual void Contribute(TPZMaterialData &data, REAL weight, TPZFMatrix< STATE > &ek, TPZFMatrix< STATE > &ef) override
It computes a contribution to the stiffness matrix and load vector at one integration point...
Definition: pzmatorthotropic.cpp:119
TPZMatOrthotropic::fXf
TPZFMatrix< STATE > fXf
Definition: pzmatorthotropic.h:25
TPZFMatrix::Transpose
void Transpose(TPZMatrix< TVar > *const T) const override
It makes *T the transpose of current matrix.
Definition: pzfmatrix.cpp:1073
TPZMatOrthotropic::fKXZ
TPZFMatrix< STATE > fKXZ
Definition: pzmatorthotropic.h:21
TPZMatOrthotropic::fKXX
TPZFMatrix< STATE > fKXX
Definition: pzmatorthotropic.h:21
PZError
#define PZError
Defines the output device to error messages and the DebugStop() function.
Definition: pzerror.h:15
TPZMatOrthotropic::Errors
void Errors(TPZVec< REAL > &x, TPZVec< STATE > &u, TPZFMatrix< STATE > &dudx, TPZFMatrix< REAL > &axes, TPZVec< STATE > &flux, TPZVec< STATE > &u_exact, TPZFMatrix< STATE > &du_exact, TPZVec< REAL > &values) override
Computes the error due to the difference between the interpolated flux and the flux computed based o...
Definition: pzmatorthotropic.cpp:464
TPZMatOrthotropic::Solution
virtual void Solution(TPZVec< STATE > &Sol, TPZFMatrix< STATE > &DSol, TPZFMatrix< REAL > &axes, int var, TPZVec< STATE > &Solout) override
Definition: pzmatorthotropic.cpp:302
Generated on Sun Aug 16 2020 11:17:03 for NeoPZ by   doxygen 1.8.13