NeoPZ
pzmatplaca2.cpp
Go to the documentation of this file.
1 
6 #include "pzmatplaca2.h"
7 #include "TPZMaterial.h"
8 
9 #include "pzbndcond.h"
10 #include <math.h>
11 #include <fstream>
12 
13 using namespace std;
14 
15 #include "pzvec.h"
16 #include "pzerror.h"
17 
18 
19 TPZMatPlaca2::TPZMatPlaca2(int num, STATE h, STATE f, STATE E1 , STATE E2 ,
20  STATE ni1 , STATE ni2 , STATE G12 , STATE G13 ,
21  STATE G23 , TPZFMatrix<STATE> &naxes, TPZVec<STATE> &xf) :
22 TPZRegisterClassId(&TPZMatPlaca2::ClassId),
23 TPZMaterial(num),
24 fIdfMax(6),fE1(E1), fE2(E2), fG12(G12), fG13(G13), fG23(G23),
25 fh(h),ff(f),fmi(1./(-1.+ni1*ni2)),fni1(ni1),fni2(ni2),
26 fnaxes(naxes),
27 fRmat(6,6,0.),fRmatT(6,6,0.),
28 fKxxR(6,6,0.),fKyyR(6,6,0.),
29 fKxyR(6,6,0.),fKyxR(6,6,0.),
30 fBx0R(6,6,0.),fB0xR(6,6,0.),
31 fBy0R(6,6,0.),fB0yR(6,6,0.),
32 fB00R(6,6,0.),
33 fKxx(6,6,0.),fKyy(6,6,0.),
34 fKxy(6,6,0.),fKyx(6,6,0.),
35 fBx0(6,6,0.),fB0x(6,6,0.),
36 fBy0(6,6,0.),fB0y(6,6,0.),
37 fB00(6,6,0.), fXF(xf)
38 {
39  STATE Small , k, mi;
40  Small = E1*STATE(1.E-5);
41  k = 5./6.; // coeficiente de cisalhamento
42  mi = 1.0/(-1.0 + ni1 * ni2);
43 
44  fRmat(0,0) = fnaxes(0,0); fRmat(0,1) = fnaxes(0,1); fRmat(0,2) = fnaxes(0,2);
45  fRmat(1,0) = fnaxes(1,0); fRmat(1,1) = fnaxes(1,1); fRmat(1,2) = fnaxes(1,2);
46  fRmat(2,0) = fnaxes(2,0); fRmat(2,1) = fnaxes(2,1); fRmat(2,2) = fnaxes(2,2);
47 
48  fRmat(3,3) = fnaxes(0,0); fRmat(3,4) = fnaxes(0,1); fRmat(3,5) = fnaxes(0,2);
49  fRmat(4,3) = fnaxes(1,0); fRmat(4,4) = fnaxes(1,1); fRmat(4,5) = fnaxes(1,2);
50  fRmat(5,3) = fnaxes(2,0); fRmat(5,4) = fnaxes(2,1); fRmat(5,5) = fnaxes(2,2);
51 
52  fRmat.Transpose(&fRmatT);
53 
54 
55 
56  fKxx(0,0) = -E1*h*mi;
57  fKxx(0,4) = -E1*f*h*mi;
58  fKxx(1,1) = G12*h/2. + h*Small/4.;
59  fKxx(1,3) = -f*G12*h/2.;
60  fKxx(2,2) = G13*h*k/2.;
61  fKxx(3,1) = -f*G12*h/2.;
62  fKxx(3,3) = f*f*G12*h/2. + G12*h*h*h/24.;
63  fKxx(4,0) = -E1*f*h*mi;
64  fKxx(4,4) = -E1*f*f*h*mi - E1*h*h*h*mi/12.;
65 
66  fKxy(0,1) = -E1*h*mi*ni2;
67  fKxy(0,3) = E1*f*h*mi*ni2;
68  fKxy(1,0) = G12*h/2. - h*Small/4.;
69  fKxy(1,4) = f*G12*h/2.;
70  fKxy(3,0) = -f*G12*h/2.;
71  fKxy(3,4) = -f*f*G12*h/2. - G12*h*h*h/24.;
72  fKxy(4,1) = -E1*f*h*mi*ni2;
73  fKxy(4,3) = E1*f*f*h*mi*ni2 + E1*h*h*h*mi*ni2/12.;
74 
75  fKxy.Transpose(&fKyx);
76 
77  fKyy(0,0) = G12*h/2. + h*Small/4.;
78  fKyy(0,4) = f*G12*h/2.;
79  fKyy(1,1) = -E2*h*mi;
80  fKyy(1,3) = E2*f*h*mi;
81  fKyy(2,2) = G23*h*k/2.;
82  fKyy(3,1) = E2*f*h*mi;
83  fKyy(3,3) = -E2*f*f*h*mi - E2*h*h*h*mi/12.;
84  fKyy(4,0) = f*G12*h/2.;
85  fKyy(4,4) = f*f*G12*h/2. + G12*h*h*h/24.;
86 
87  fB0x(4,2) = G13*h*k/2.;
88  fB0x(5,1) = -h*Small/2.;
89 
90  fB0x.Transpose(&fBx0);
91 
92  fB0y(3,2) = -G23*h*k/2.;
93  fB0y(5,0) = h*Small/2.;
94 
95  fB0y.Transpose(&fBy0);
96 
97  fB00(3,3) = G23*h*k/2.;
98  fB00(4,4) = G13*h*k/2.;
99  fB00(5,5) = h*Small;
100 
101  fKxxR = fRmatT * (fKxx * fRmat);
102  fKyxR = fRmatT * (fKyx * fRmat);
103  fKxyR = fRmatT * (fKxy * fRmat);
104  fKyyR = fRmatT * (fKyy * fRmat);
105  fB0xR = fRmatT * (fB0x * fRmat);
106  fB0yR = fRmatT * (fB0y * fRmat);
107  fBx0R = fRmatT * (fBx0 * fRmat);
108  fBy0R = fRmatT * (fBy0 * fRmat);
109  fB00R = fRmatT * (fB00 * fRmat);
110 
111 }
112 
113 void TPZMatPlaca2::SetNAxes(TPZFMatrix<STATE> &n) {
114 
115  fnaxes = n;
116  int numl = fIdfMax/3;
117  int i;
118  for(i=0;i<numl;i++) {
119  fRmat(3*i+0,3*i+0) = fnaxes(0,0);
120  fRmat(3*i+0,3*i+1) = fnaxes(0,1);
121  fRmat(3*i+0,3*i+2) = fnaxes(0,2);
122  fRmat(3*i+1,3*i+0) = fnaxes(1,0);
123  fRmat(3*i+1,3*i+1) = fnaxes(1,1);
124  fRmat(3*i+1,3*i+2) = fnaxes(1,2);
125  fRmat(3*i+2,3*i+0) = fnaxes(2,0);
126  fRmat(3*i+2,3*i+1) = fnaxes(2,1);
127  fRmat(3*i+2,3*i+2) = fnaxes(2,2);
128  }
129 
130  fRmat.Transpose(&fRmatT);
131 
132  TPZFMatrix<STATE> tmp;
133  fKxx.Multiply(fRmat,tmp);
134  fRmatT.Multiply(tmp,fKxxR);
135  fKyxR = fRmatT * (fKyx * fRmat);
136  fKxyR = fRmatT * (fKxy * fRmat);
137  fKyyR = fRmatT * (fKyy * fRmat);
138  fB0xR = fRmatT * (fB0x * fRmat);
139  fB0yR = fRmatT * (fB0y * fRmat);
140  fBx0R = fRmatT * (fBx0 * fRmat);
141  fBy0R = fRmatT * (fBy0 * fRmat);
142  fB00R = fRmatT * (fB00 * fRmat);
143 }
144 
146 ofstream placatest("placatest.dat");
147 void TPZMatPlaca2::Contribute(TPZMaterialData &data,
148  REAL weight,
149  TPZFMatrix<STATE> &ek,
150  TPZFMatrix<STATE> &ef) {
151  // this method adds the contribution of the material to the stiffness
152  // matrix and right hand side
153 
154  // check on the validity of the arguments
155  //rows x cols
156  TPZFMatrix<REAL> &dphi = data.dphix;
157  TPZFMatrix<REAL> &phi = data.phi;
158  TPZManVector<REAL,3> &x = data.x;
159  TPZFMatrix<STATE> axes(data.axes.Rows(),data.axes.Cols());
160  for (int r=0; r<axes.Rows(); r++) {
161  for (int c=0; c<axes.Cols(); c++) {
162  axes(r,c) = data.axes(r,c);
163  }
164  }
165 
166  if(phi.Cols() != 1 || dphi.Rows() != 2 || phi.Rows() != dphi.Cols()){
167  PZError << "TPZMatPlaca2.contr, inconsistent input data : phi.Cols() = "
168  << phi.Cols() << " dphi.Cols + " << dphi.Cols() <<
169  " phi.Rows = " << phi.Rows() << " dphi.Rows = " <<
170  dphi.Rows() << "\n";
171  }
172  if(fForcingFunction) {
173  fForcingFunction->Execute(x,fXF);//fXf = xfloat
174  }
175 
176 
177  STATE Dax1n1, Dax1n2, Dax2n1, Dax2n2;
178 
179  Dax1n1 = axes(0,0)* fnaxes(0,0) + axes(0,1)* fnaxes(0,1) + axes(0,2)* fnaxes(0,2);
180  Dax1n2 = axes(0,0)* fnaxes(1,0) + axes(0,1)* fnaxes(1,1) + axes(0,2)* fnaxes(1,2);
181  Dax2n1 = axes(1,0)* fnaxes(0,0) + axes(1,1)* fnaxes(0,1) + axes(1,2)* fnaxes(0,2);
182  Dax2n2 = axes(1,0)* fnaxes(1,0) + axes(1,1)* fnaxes(1,1) + axes(1,2)* fnaxes(1,2);
183 // intern33 = axes(2,0)* fnaxes(2,0) + axes(2,1)* fnaxes(2,1) + axes(2,2)* fnaxes(2,2);
184 // if(fabs(fabs(intern33)-1.) > 1.e-6) {
185 // PZError << "third axes of the material not parallel to element axis " << intern33 << endl;
186 // PZError << axes(2,0) << ' ' << axes(2,1) << ' ' << axes(2,2) << endl;
187 // PZError << fnaxes(2,0) << ' ' << fnaxes(2,1) << ' ' << fnaxes(2,2) << endl;
188 // }
189 
190  TPZFMatrix<STATE> Kn1n1(fIdfMax,fIdfMax,0.),Kn2n2(fIdfMax,fIdfMax,0.),
191  Kn1n2(fIdfMax,fIdfMax,0.),Kn2n1(fIdfMax,fIdfMax,0.),
192  Bn10(fIdfMax,fIdfMax,0.) ,B0n1(fIdfMax,fIdfMax,0.),
193  Bn20(fIdfMax,fIdfMax,0.),B0n2(fIdfMax,fIdfMax,0.),
194  B000(fIdfMax,fIdfMax,0.);
195 
196 
197  Kn1n1 = fKxxR * Dax1n1 * Dax1n1 + fKyyR * Dax1n2 * Dax1n2 +
198  fKxyR * Dax1n1 * Dax1n2 + fKyxR * Dax1n1 * Dax1n2 ;
199 
200  Kn2n2 = fKxxR * Dax2n1 * Dax2n1 + fKyyR * Dax2n2 * Dax2n2 +
201  fKxyR * Dax2n2 * Dax2n1 + fKyxR * Dax2n2 * Dax2n1;
202 
203  Kn1n2 = fKxxR * Dax1n1 * Dax2n1 + fKyyR * Dax1n2 * Dax2n2 +
204  fKxyR * Dax1n1 * Dax2n2 + fKyxR * Dax1n2 * Dax2n1;
205 
206  Kn2n1 = fKxxR * Dax2n1 * Dax1n1 + fKyyR * Dax2n2 * Dax1n2 +
207  fKxyR * Dax2n1 * Dax1n2 + fKyxR * Dax2n2 * Dax1n1;
208 
209 
210  B0n1 = fB0xR * Dax1n1 + fB0yR * Dax1n2 ;
211 
212  Bn10 = fBx0R * Dax1n1 + fBy0R * Dax1n2 ;
213 
214  B0n2 = fB0xR * Dax2n1 + fB0yR * Dax2n2;
215 
216  Bn20 = fBx0R * Dax2n1 + fBy0R * Dax2n2;
217 
218  B000 = fB00R ;
219 
220  /*
221  Kn1n1.Print("Kn1n1 = ",placatest);
222  Kn1n2.Print("Kn1n2 = ",placatest);
223  Kn2n2.Print("Kn2n2 = ",placatest);
224  Kn2n1.Print("Kn2n1 = ",placatest);
225  B0n1.Print("B0n1 = ",placatest);
226  Bn10.Print("Bn10 = ",placatest);
227  B0n2.Print("B0n2 = ",placatest);
228  Bn20.Print("B0n1 = ",placatest);
229  B000.Print("B00 = ",placatest);
230  placatest.flush();
231  */
232  int nshape = phi.Rows();
233  TPZFMatrix<STATE> KIJ(fIdfMax,fIdfMax);
234 
235 
236 
237  int idf,jdf,i,j;
238  STATE contrib[3];
239  for(i=0;i<3;i++) {
240  contrib[i]=0.;
241  for(j=0;j<3;j++) {
242  contrib[i] += fXF[j]*fnaxes(j,i);
243  }
244  }
245  for(i=0; i<nshape; i++) {
246  for(idf=0; idf<3; idf++) {
247  ef(fIdfMax*i+idf) += weight*contrib[idf]*phi(i,0);
248  }
249  for(idf=3; idf<fIdfMax; idf++) ef(fIdfMax*i+idf) += weight*fXF[idf]*phi(i,0);
250  for(j=0; j<nshape; j++) {
251  STATE dphi_0i(dphi(0,i)),dphi_1i(dphi(1,i)),dphi_0j(dphi(0,j)),dphi_1j(dphi(1,j));
252  STATE phi_i(phi(i,0)),phi_j(phi(j,0));
253  KIJ = ((STATE)weight)*(dphi_0i*dphi_0j*Kn1n1+
254  dphi_0i*dphi_1j*Kn1n2+
255  dphi_1i*dphi_0j*Kn2n1+
256  dphi_1i*dphi_1j*Kn2n2+
257  dphi_0i*phi_j *Bn10 +
258  dphi_1i*phi_j *Bn20 +
259  phi_i *dphi_0j*B0n1 +
260  phi_i *dphi_1j*B0n2 +
261  phi_i *phi_j *B000 );
262  for(idf=0; idf<fIdfMax; idf++) for(jdf=0; jdf<fIdfMax; jdf++)
263  ek(i*fIdfMax+idf,j*fIdfMax+jdf) += KIJ(idf,jdf);
264  }
265  }
266 }
267 
268 void TPZMatPlaca2::ContributeBC(TPZMaterialData &data,
269  REAL weight,
270  TPZFMatrix<STATE> &ek,
271  TPZFMatrix<STATE> &ef,
272  TPZBndCond &bc) {
273  TPZFMatrix<REAL> &phi = data.phi;
274 
275 // if(bc.Material().operator ->() != this){
276 // PZError << "TPZMatPlaca2.ContributeBC warning : this material didn't create the boundary condition!\n";
277 // }
278 
279  if(bc.Type() < 0 && bc.Type() > 2){
280  PZError << "TPZMatPlaca2.aplybc, unknown boundary condition type :" <<
281  bc.Type() << " boundary condition ignored\n";
282  }
283 
284  int numdof = NStateVariables();
285  int numnod = ek.Rows()/numdof;
286  int r = numdof;
287 
288  int idf,jdf,in,jn;
289  switch(bc.Type()){
290 
291  case 0:
292  for(in=0 ; in<numnod ; ++in){
293  for(idf = 0;idf<r;idf++) {
294  (ef)(in*r+idf,0) += gBigNumber*phi(in,0)*bc.Val2()(idf,0)*weight;
295  }
296  for(jn=0 ; jn<numnod ; ++jn) {
297  for(idf = 0;idf<r;idf++) {
298  ek(in*r+idf,jn*r+idf) += gBigNumber*phi(in,0)*phi(jn,0)*weight;
299  }
300  }
301  }
302  break;
303 
304  case 1:
305  for(in=0 ; in<numnod ; ++in){
306  for(idf = 0;idf<r;idf++) {
307  (ef)(in*r+idf,0) += phi(in,0)*bc.Val2()(idf,0)*weight;
308  }
309  }
310  break;
311 
312  case 2:
313  for(in=0 ; in<numnod ; ++in){
314  for(idf = 0;idf<r;idf++) {
315  (ef)(in*r+idf,0) += phi(in,0)*bc.Val2()(idf,0)*weight;
316  }
317  for(jn=0 ; jn<numnod ; ++jn) {
318  for(idf = 0;idf<r;idf++) {
319  for(jdf = 0;jdf<r;jdf++) {
320  ek(in*r+idf,jn*r+jdf) += bc.Val1()(idf,jdf)*phi(in,0)*phi(jn,0)*weight;
321  }
322  }
323  }
324  }//fim switch
325  }
326 }
327 
328 int TPZMatPlaca2::NFluxes() {return 1;}
329 
330 void TPZMatPlaca2::Flux(TPZVec<REAL> &/*x*/,TPZVec<STATE> &/*u*/,TPZFMatrix<STATE> &/*dudx*/,TPZFMatrix<REAL> &/*axes*/,TPZVec<STATE> &/*fl*/) {
331  PZError << "TPZMatPlaca2::Flux is called\n";
332 }
333 
334 void TPZMatPlaca2::Errors(TPZVec<REAL> &/*x*/,TPZVec<STATE> &/*u*/,TPZFMatrix<STATE> &/*dudx*/,TPZFMatrix<REAL> &/*axes*/,TPZVec<STATE> &/*flux*/,
335  TPZVec<STATE> &/*u_exact*/,TPZFMatrix<STATE> &/*du_exact*/,TPZVec<REAL> &/*values*/) {
336  PZError << "TPZMatPlaca2::Errors is called\n";
337 }
338 
339 void TPZMatPlaca2::Print(std::ostream & out) {
340  //out << "Material type TPZMatPlaca2 -- number = " << Id() << "\n";
341  //out << "Matrix xk -> "; fXk.Print("fXk",out);
342  //out << "Matrix xC -> "; fCk.Print("fCf",out);
343  //out << "Matrix xf -> "; fXf.Print("fXf",out);
344 }
345 
347 int TPZMatPlaca2::VariableIndex(const std::string &name){
348  if(!strcmp(name.c_str(),"State")) return 0;
349  if(!strcmp(name.c_str(),"Deslocx")) return 2;// Desloc. eixo x global
350  if(!strcmp(name.c_str(),"Deslocy")) return 3;// Desloc. eixo y global
351  if(!strcmp(name.c_str(),"Deslocz")) return 4;// Desloc. eixo z global
352  if(!strcmp(name.c_str(),"MnVec")) return 5;// Momentos nas direcoes dos eixos n1 e n2
353  if(!strcmp(name.c_str(),"MaVec")) return 50;// Momentos nas direcoes dos eixos a1 e a2
354  if(!strcmp(name.c_str(),"Vn1")) return 8;// forca cortante Vn1 (positiva se antihorario)
355  if(!strcmp(name.c_str(),"Vn2")) return 9;// forca cortante Vn2 (positiva se antihorario)
356  if(!strcmp(name.c_str(),"Sign1")) return 10;// tens� normal na dire�o n1
357  if(!strcmp(name.c_str(),"Sign2")) return 11;// tens� normal na dire�o n2
358  if(!strcmp(name.c_str(),"Taun1n2")) return 12;// tens� cisalhamento eixos n1 e n2
359  if(!strcmp(name.c_str(),"NaVec")) return 54;//Tensoes normais nas direcoes dos eixos a1,a2
360  if(!strcmp(name.c_str(),"Taun1n3")) return 13;// tens� cisalhamento eixos n1 e n3
361  if(!strcmp(name.c_str(),"Taun2n3")) return 14;// tens� cisalhamento eixos n2 e n3
362  if(!strcmp(name.c_str(),"Displacement")) return 15;// translacoes u,v,w
363 
364 
365 
366  return TPZMaterial::VariableIndex(name);
367 }
368 
371 int TPZMatPlaca2::NSolutionVariables(int var){
372  if(var == 2) return 1; // Desloc. na superficie de referencia na eixo x global
373  if(var == 3) return 1; // Desloc. na superficie de referencia na eixo y global
374  if(var == 4) return 1; // Desloc. na superficie de referencia na eixo z global
375  if(var == 5) return 3; // Momentos na superficie de referencia nas direcoes dos eixos n1 e n2
376  if(var== 50) return 3; // Momentos na superficie de referencia nas direcoes dos eixos a1 e a2
377  if(var == 8) return 1;
378  if(var == 9) return 1;
379  if(var == 10) return 1;
380  if(var == 11) return 1;
381  if(var == 12) return 1;
382  if(var == 54) return 3;
383  if(var == 13) return 1;
384  if(var == 14) return 1;
385  if(var == 15) return 3;
386 
387  return TPZMaterial::NSolutionVariables(var); // todos os Deslocamentos nodais e suas derivadas nas
388  // direcoes dos eixos axes (a1 e a2)
389 }
390 
392 void TPZMatPlaca2::Solution(TPZVec<STATE> &Sol,TPZFMatrix<STATE> &DSol,
393  TPZFMatrix<REAL> &axes,int var,TPZVec<STATE> &Solout){
394 
395  REAL k = 5./6.;
396 
397  if(var == 2) {
398  Solout.Resize(1);
399  Solout[0] = Sol[0];
400  return;
401  }
402  if(var == 3) {
403  Solout.Resize(1);
404  Solout[0] = Sol[1];
405  return;
406  }
407  if(var == 4) {
408  Solout.Resize(1);
409  Solout[0] = Sol[2];
410  return;
411  }
412  if(var == 15) {
413  Solout.Resize(3);
414  Solout[0] = Sol[0];
415  Solout[1] = Sol[1];
416  Solout[2] = Sol[2];
417  return;
418  }
419 
420 
421  TPZVec<REAL> Soln(6);
422 
423  TPZFMatrix<REAL> DSolnax(2,6),DSolnn(2,6);
424 
425 
426  int idf,jdf;
427  for(idf=0; idf<6; idf++) {
428  Soln[idf] = 0;
429  DSolnax(0,idf) = 0;
430  DSolnax(1,idf) = 0.;
431  for(jdf=0; jdf<6; jdf++) {
432  Soln[idf] += fRmat(idf,jdf)*Sol[jdf];
433  DSolnax(0,idf) += fRmat(idf,jdf)*DSol(0,jdf);
434  DSolnax(1,idf) += fRmat(idf,jdf)*DSol(1,jdf);
435  }
436  }
437 
438  TPZFMatrix<REAL> Rmatan(2,2);
439  Rmatan(0,0)= axes(0,0)* fnaxes(0,0) + axes(0,1)* fnaxes(0,1) + axes(0,2)* fnaxes(0,2);
440  Rmatan(1,0)= axes(0,0)* fnaxes(1,0) + axes(0,1)* fnaxes(1,1) + axes(0,2)* fnaxes(1,2);
441  Rmatan(0,1)= axes(1,0)* fnaxes(0,0) + axes(1,1)* fnaxes(0,1) + axes(1,2)* fnaxes(0,2);
442  Rmatan(1,1)= axes(1,0)* fnaxes(1,0) + axes(1,1)* fnaxes(1,1) + axes(1,2)* fnaxes(1,2);
443 
444  for(idf=0;idf<6;idf++) {
445  DSolnn(0,idf) = Rmatan(0,0)*DSolnax(0,idf)+Rmatan(0,1)*DSolnax(1,idf);
446  DSolnn(1,idf) = Rmatan(1,0)*DSolnax(0,idf)+Rmatan(1,1)*DSolnax(1,idf);
447  }
448 
449  if (var==5){
450  REAL Mn1;
451  Mn1 = -fE1*fh*fh*fh*fmi*(-fni2*DSolnn(1,3)+DSolnn(0,4))/12.;
452  Mn1 += fE1*ff*fh*fmi*(ff*fni2*DSolnn(1,3)-fni2*DSolnn(1,1)
453  -ff*DSolnn(0,4)-DSolnn(0,0));
454 
455 
456  REAL Mn2;
457  Mn2 =-(fE2*fh*fh*fh*fmi*(-DSolnn(1,3) + fni1*DSolnn(0,4)))/12.;
458  Mn2 += fE2*ff*fh*fmi*(ff*DSolnn(1,3) - DSolnn(1,1) -
459  ff*fni1*DSolnn(0,4) - fni1*DSolnn(0,0));
460 
461  REAL Mn1n2;
462  Mn1n2 = fG12*fh*fh*fh*(DSolnn(1,4) - DSolnn(0,3))/24.;
463  Mn1n2 += (ff*fG12*fh*(ff*DSolnn(1,4) + DSolnn(1,0) -
464  ff*DSolnn(0,3) + DSolnn(0,1)))/2.;
465 
466  Solout.Resize(3);
467  Solout[0]=Mn1;
468  Solout[1]=Mn2;
469  Solout[2]=Mn1n2;
470  return;
471  }
472 
473  if (var==50){
474 
475  REAL Mn1;
476  Mn1 = -fE1*fh*fh*fh*fmi*(-fni2*DSolnn(1,3)+DSolnn(0,4))/12.;
477  Mn1 += fE1*ff*fh*fmi*(ff*fni2*DSolnn(1,3)-fni2*DSolnn(1,1)
478  -ff*DSolnn(0,4)-DSolnn(0,0));
479 
480 
481  REAL Mn2;
482  Mn2 =-(fE2*fh*fh*fh*fmi*(-DSolnn(1,3) + fni1*DSolnn(0,4)))/12.;
483  Mn2 += fE2*ff*fh*fmi*(ff*DSolnn(1,3) - DSolnn(1,1) -
484  ff*fni1*DSolnn(0,4) - fni1*DSolnn(0,0));
485 
486  REAL Mn1n2;
487  Mn1n2 = fG12*fh*fh*fh*(DSolnn(1,4) - DSolnn(0,3))/24.;
488  Mn1n2 += (ff*fG12*fh*(ff*DSolnn(1,4) + DSolnn(1,0) -
489  ff*DSolnn(0,3) + DSolnn(0,1)))/2.;
490  REAL Ma1;
491  Ma1 = Mn1 * Rmatan(0,0) * Rmatan(0,0) + Mn2*Rmatan(1,0) * Rmatan(1,0) + 2.* Rmatan(0,0) * Rmatan(1,0) * Mn1n2;
492  REAL Ma2;
493  Ma2 = Mn1 * Rmatan(0,1) * Rmatan(0,1) + Mn2*Rmatan(1,1) * Rmatan(1,1) + 2.* Rmatan(1,1) * Rmatan(0,1) * Mn1n2;
494 
495  REAL Ma1a2;
496  Ma1a2 = Mn1n2 * (Rmatan(0,0)*Rmatan(0,0) + Rmatan(1,0) *Rmatan(0,1))+ (Mn1-Mn2)*Rmatan(0,0)*Rmatan(0,1);
497 
498  Solout.Resize(3);
499  Solout[0]=Ma1;
500  Solout[1]=Ma2;
501  Solout[2]=Ma1a2;
502  return;
503  }
504 
505  if (var == 8) {
506  REAL Vn1;
507  Vn1 = (fG12*fh*k*(Soln[4]+DSolnn(0,2)));
508  Solout[0] = Vn1;
509  return;
510  }
511  if (var == 9) {
512  REAL Vn2;
513  Vn2 = (fG12*fh*k*(-Soln[3]+ DSolnn(1,2)))/2.;
514  Solout[0] = Vn2;
515  return;
516  }
517 
518  REAL Sign1 = 0.;
519 
520  if (var == 10 || var == 54) {
521  REAL z;
522  z=0.;
523  Sign1 = -fE1*fmi*(fni2*(-(ff + z)*DSolnn(1,3) + DSolnn(1,1))+
524  (ff + z)*DSolnn(0,4) + DSolnn(0,0));
525  if(var == 10) {
526  Solout[0] = Sign1;
527  return;
528  }
529  }
530  REAL Sign2 = 0.;
531  if (var == 11 || var == 54) {
532  REAL z;
533  z=0.;
534  Sign2=-fE2*fmi*(-(ff + z)*DSolnn(1,3) + DSolnn(1,1) +
535  fni1*((ff + z)*DSolnn(0,4) + DSolnn(0,0)));
536  if(var == 11) {
537  Solout[0] = Sign2;
538  return;
539  }
540  }
541  REAL Taun1n2 = 0.;
542  if (var == 12 || var == 54) {
543  REAL z;
544  z=0.;
545  Taun1n2=(fG12*((ff + z)*DSolnn(1,4) + DSolnn(1,0) -
546  (ff + z)*DSolnn(0,3) + DSolnn(0,1)))/2.;
547  if(var == 12) {
548  Solout[0] = Taun1n2;
549  return;
550  }
551  }
552 
553  if(var == 54) {
554  REAL Siga1;
555  Siga1 = Sign1 * Rmatan(0,0) * Rmatan(0,0) + Sign2*Rmatan(1,0) * Rmatan(1,0) + 2.* Rmatan(0,0) * Rmatan(1,0) * Taun1n2;
556  REAL Siga2;
557  Siga2 = Sign1 * Rmatan(0,1) * Rmatan(0,1) + Sign2*Rmatan(1,1) * Rmatan(1,1) + 2.* Rmatan(1,1) * Rmatan(0,1) * Taun1n2;
558 
559  REAL Siga1a2;
560  Siga1a2 = Taun1n2 * (Rmatan(0,0)*Rmatan(0,0) + Rmatan(1,0) *Rmatan(0,1))+ (Sign1-Sign2)*Rmatan(0,0)*Rmatan(0,1);
561 
562  Solout.Resize(3);
563  Solout[0]=Siga1*fh;
564  Solout[1]=Siga2*fh;
565  Solout[2]=Siga1a2*fh;
566  return;
567  }
568  if (var == 13 ) {
569  REAL Taun1n3;
570  Taun1n3=fG13*(Soln[4] + DSolnn(0,2))/2.;
571  Solout[0] = Taun1n3;
572  return;
573  }
574  if (var == 14 ) {
575  REAL Taun2n3;
576  Taun2n3=fG23*(-Soln[3] + DSolnn(1,2))/2.;
577  Solout[0] = Taun2n3;
578  return;
579  }
580 
581  TPZMaterial::Solution(Sol,DSol,axes,var,Solout);
582 }
583 
584 int TPZMatPlaca2::ClassId() const{
585  return Hash("TPZMatPlaca2") ^ TPZMaterial::ClassId() << 1;
586 }
TPZMaterial::Solution
virtual void Solution(TPZMaterialData &data, int var, TPZVec< STATE > &Solout)
Returns the solution associated with the var index based on the finite element approximation.
Definition: TPZMaterial.cpp:200
TPZFunction::Execute
virtual void Execute(const TPZVec< REAL > &x, TPZVec< TVar > &f, TPZFMatrix< TVar > &df)
Performs function computation.
Definition: pzfunction.h:38
TPZMatPlaca2::Solution
virtual void Solution(TPZVec< STATE > &Sol, TPZFMatrix< STATE > &DSol, TPZFMatrix< REAL > &axes, int var, TPZVec< STATE > &Solout) override
Returns the solution associated with the var index based on the finite element approximation.
Definition: pzmatplaca2.cpp:392
TPZMatPlaca2::fB0x
TPZFMatrix< STATE > fB0x
Definition: pzmatplaca2.h:28
TPZMatPlaca2::fB00R
TPZFMatrix< STATE > fB00R
Definition: pzmatplaca2.h:27
TPZMatPlaca2::fRmat
TPZFMatrix< STATE > fRmat
Definition: pzmatplaca2.h:26
TPZMaterialData::x
TPZManVector< REAL, 3 > x
value of the coordinate at the integration point
Definition: pzmaterialdata.h:71
bc
clarg::argBool bc("-bc", "binary checkpoints", false)
TPZMaterial.h
TPZBndCond::Type
int Type()
Definition: pzbndcond.h:249
pzvec.h
Templated vector implementation.
TPZMatPlaca2::fKyyR
TPZFMatrix< STATE > fKyyR
Definition: pzmatplaca2.h:27
TPZMatPlaca2::fRmatT
TPZFMatrix< STATE > fRmatT
Definition: pzmatplaca2.h:26
TPZMatPlaca2
DESCRIBE PLEASE.
Definition: pzmatplaca2.h:20
TPZMatPlaca2::fnaxes
TPZFMatrix< STATE > fnaxes
Definition: pzmatplaca2.h:25
TPZMaterial::VariableIndex
virtual int VariableIndex(const std::string &name)
Returns the variable index associated with the name.
Definition: TPZMaterial.cpp:141
pzerror.h
Defines PZError.
TPZMatPlaca2::NStateVariables
virtual int NStateVariables() const override
Returns the number of state variables associated with the material.
Definition: pzmatplaca2.h:40
h
clarg::argBool h("-h", "help message", false)
pzmatplaca2.h
Contains the TPZMatPlaca2 class.
tmp
AutoPointerMutexArrayInit tmp
Definition: tpzautopointer.cpp:42
TPZMatPlaca2::VariableIndex
virtual int VariableIndex(const std::string &name) override
Returns the variable index associated with the name.
Definition: pzmatplaca2.cpp:347
TPZMatPlaca2::fBy0R
TPZFMatrix< STATE > fBy0R
Definition: pzmatplaca2.h:27
TPZMatPlaca2::fKxxR
TPZFMatrix< STATE > fKxxR
Definition: pzmatplaca2.h:27
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
TPZMaterialData::dphix
TPZFNMatrix< 660, REAL > dphix
values of the derivative of the shape functions
Definition: pzmaterialdata.h:59
TPZMaterial
This abstract class defines the behaviour which each derived class needs to implement.
Definition: TPZMaterial.h:39
TPZMatPlaca2::fKyx
TPZFMatrix< STATE > fKyx
Definition: pzmatplaca2.h:28
TPZMatPlaca2::SetNAxes
void SetNAxes(TPZFMatrix< STATE > &n)
Modify the direction of the fibres for the plate.
Definition: pzmatplaca2.cpp:113
pzbndcond.h
Contains the TPZBndCond class which implements a boundary condition for TPZMaterial objects...
TPZMatPlaca2::fBx0
TPZFMatrix< STATE > fBx0
Definition: pzmatplaca2.h:28
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
TPZMatPlaca2::Print
void Print(std::ostream &out) override
Prints out the data associated with the material.
Definition: pzmatplaca2.cpp:339
test.f
f
Definition: test.py:287
TPZMatPlaca2::ClassId
virtual int ClassId() const override
Define the class id associated with the class.
Definition: pzmatplaca2.cpp:584
TPZMatPlaca2::fKxx
TPZFMatrix< STATE > fKxx
Definition: pzmatplaca2.h:28
TPZMatPlaca2::fXF
TPZVec< STATE > fXF
Definition: pzmatplaca2.h:29
TPZBndCond
This class defines the boundary condition for TPZMaterial objects.
Definition: pzbndcond.h:29
placatest
ofstream placatest("placatest.dat")
Output file to write data of the test over shell (placa)
TPZMatPlaca2::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: pzmatplaca2.cpp:268
TPZMatrix::Rows
int64_t Rows() const
Returns number of rows.
Definition: pzmatrix.h:803
TPZMatPlaca2::NFluxes
virtual int NFluxes() override
Returns the number of components which form the flux function.
Definition: pzmatplaca2.cpp:328
TPZMaterialData::axes
TPZFNMatrix< 9, REAL > axes
axes indicating the directions of the derivatives of the shapefunctions
Definition: pzmaterialdata.h:63
TPZMatPlaca2::fB0yR
TPZFMatrix< STATE > fB0yR
Definition: pzmatplaca2.h:27
TPZBndCond::Val1
TPZFMatrix< STATE > & Val1()
Definition: pzbndcond.h:253
TPZMatPlaca2::fB00
TPZFMatrix< STATE > fB00
Definition: pzmatplaca2.h:28
TPZMatPlaca2::fBy0
TPZFMatrix< STATE > fBy0
Definition: pzmatplaca2.h:28
TPZMatPlaca2::fKxyR
TPZFMatrix< STATE > fKxyR
Definition: pzmatplaca2.h:27
TPZMatPlaca2::Flux
virtual void Flux(TPZVec< REAL > &x, TPZVec< STATE > &u, TPZFMatrix< STATE > &dudx, TPZFMatrix< REAL > &axes, TPZVec< STATE > &fl) override
Computes the value of the flux function to be used by ZZ error estimator.
Definition: pzmatplaca2.cpp:330
TPZMaterial::gBigNumber
static REAL gBigNumber
Big number to penalization method, used for Dirichlet conditions.
Definition: TPZMaterial.h:83
Hash
int32_t Hash(std::string str)
Definition: TPZHash.cpp:10
TPZMaterial::ClassId
int ClassId() const override
Unique identifier for serialization purposes.
Definition: TPZMaterial.cpp:382
TPZMatPlaca2::fB0xR
TPZFMatrix< STATE > fB0xR
Definition: pzmatplaca2.h:27
idf
static int idf[6]
Definition: TPZCompMeshTools.cpp:135
TPZMaterial::NSolutionVariables
virtual int NSolutionVariables(int var)
Returns the number of variables associated with the variable indexed by var.
Definition: TPZMaterial.cpp:176
TPZMatrix::Multiply
virtual void Multiply(const TPZFMatrix< TVar > &A, TPZFMatrix< TVar > &res, int opt=0) const
It mutiplies itself by TPZMatrix<TVar>A putting the result in res.
Definition: pzmatrix.cpp:1916
TPZMaterial::fForcingFunction
TPZAutoPointer< TPZFunction< STATE > > fForcingFunction
Pointer to forcing function, it is the right member at differential equation.
Definition: TPZMaterial.h:47
TPZMatPlaca2::TPZMatPlaca2
TPZMatPlaca2(int num, STATE h, STATE f, STATE E1, STATE E2, STATE ni1, STATE ni2, STATE G12, STATE G13, STATE G23, TPZFMatrix< STATE > &naxes, TPZVec< STATE > &xf)
Definition: pzmatplaca2.cpp:19
TPZMatPlaca2::fKxy
TPZFMatrix< STATE > fKxy
Definition: pzmatplaca2.h:28
TPZMatPlaca2::fB0y
TPZFMatrix< STATE > fB0y
Definition: pzmatplaca2.h:28
TPZMatPlaca2::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: pzmatplaca2.cpp:147
TPZMatrix::Cols
int64_t Cols() const
Returns number of cols.
Definition: pzmatrix.h:809
TPZMatPlaca2::fBx0R
TPZFMatrix< STATE > fBx0R
Definition: pzmatplaca2.h:27
TPZMatPlaca2::fKyxR
TPZFMatrix< STATE > fKyxR
Definition: pzmatplaca2.h:27
TPZMatPlaca2::NSolutionVariables
virtual int NSolutionVariables(int var) override
Returns the number of variables associated with the variable indexed by var.
Definition: pzmatplaca2.cpp:371
TPZFMatrix::Transpose
void Transpose(TPZMatrix< TVar > *const T) const override
It makes *T the transpose of current matrix.
Definition: pzfmatrix.cpp:1073
TPZMatPlaca2::fKyy
TPZFMatrix< STATE > fKyy
Definition: pzmatplaca2.h:28
PZError
#define PZError
Defines the output device to error messages and the DebugStop() function.
Definition: pzerror.h:15
TPZMatPlaca2::Errors
virtual 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: pzmatplaca2.cpp:334
Generated on Sun Aug 16 2020 11:17:03 for NeoPZ by   doxygen 1.8.13