18 STATE lambda,STATE coef1,STATE coef2,STATE coef3) :
27 if(coef1 != -1.0)
fCoef1 = coef1;
29 if(coef2 != -1.0)
fCoef2 = coef2;
31 if(coef3 != -1.0)
fCoef3 = coef3;
52 out <<
"name of material : " <<
Name() <<
"\n";
53 out <<
"properties : \n";
77 int nshape = phi.
Rows();
79 for(ish=0; ish<nshape; ish++) {
81 ef(ish*3+i) -= (U.
val().d(i)*dphi(0,ish)+U.
val().d(3+i)*dphi(1,ish)+U.
val().d(6+i)*dphi(2,ish))*weight;
82 for(jsh=0; jsh<nshape; jsh++) {
84 ek(ish*3+i,jsh*3+j) += (
85 U.
d(i).d(j )*dphi(0,ish)*dphi(0,jsh)+U.
d(i+3).d(j )*dphi(1,ish)*dphi(0,jsh)+U.
d(i+6).d(j )*dphi(2,ish)*dphi(0,jsh)+
86 U.
d(i).d(j+3)*dphi(0,ish)*dphi(1,jsh)+U.
d(i+3).d(j+3)*dphi(1,ish)*dphi(1,jsh)+U.
d(i+6).d(j+3)*dphi(2,ish)*dphi(1,jsh)+
87 U.
d(i).d(j+6)*dphi(0,ish)*dphi(2,jsh)+U.
d(i+3).d(j+6)*dphi(1,ish)*dphi(2,jsh)+U.
d(i+6).d(j+6)*dphi(2,ish)*dphi(2,jsh)
105 STATE global[3][3][9];
106 STATE ux,uy,uz,vx,vy,vz,wx,wy,wz;
119 for(i=0; i<3; i++)
fK2[i][i] = 0.;
130 for(i=0; i<3; i++)
fK3[i][i] = 0.;
138 for(i=0; i<3; i++)
fK4[i][i] = 0.;
146 for(i=0; i<3; i++)
fK6[i][i] = 0.;
162 for(i=0; i<3; i++)
fK8[i][i] = 0.;
170 fGradDetF[0][0] = (vy+1.)*(wz+1.) - wy*vz;
174 fGradDetF[1][1] = (ux+1.)*(wz+1.) - wx*uz;
178 fGradDetF[2][2] = (ux+1.)*(vy+1.) - vx*uy;
192 fL1[0][1] = fGradDetF[0][0]*fGradDetF[0][1];
193 fL1[0][2] = fGradDetF[0][0]*fGradDetF[0][2];
194 fL1[1][0] = fGradDetF[0][1]*fGradDetF[0][0];
195 fL1[1][1] = fGradDetF[0][1]*fGradDetF[0][1];
196 fL1[1][2] = fGradDetF[0][1]*fGradDetF[0][2];
197 fL1[2][0] = fGradDetF[0][2]*fGradDetF[0][0];
198 fL1[2][1] = fGradDetF[0][2]*fGradDetF[0][1];
199 fL1[2][2] = fGradDetF[0][2]*fGradDetF[0][2];
201 fL2[0][0] = fGradDetF[0][0]*fGradDetF[1][0];
202 fL2[0][1] = fGradDetF[0][0]*fGradDetF[1][1];
203 fL2[0][2] = fGradDetF[0][0]*fGradDetF[1][2];
204 fL2[1][0] = fGradDetF[0][1]*fGradDetF[1][0];
205 fL2[1][1] = fGradDetF[0][1]*fGradDetF[1][1];
206 fL2[1][2] = fGradDetF[0][1]*fGradDetF[1][2];
207 fL2[2][0] = fGradDetF[0][2]*fGradDetF[1][0];
208 fL2[2][1] = fGradDetF[0][2]*fGradDetF[1][1];
209 fL2[2][2] = fGradDetF[0][2]*fGradDetF[1][2];
211 fL3[0][0] = fGradDetF[0][0]*fGradDetF[2][0];
212 fL3[0][1] = fGradDetF[0][0]*fGradDetF[2][1];
213 fL3[0][2] = fGradDetF[0][0]*fGradDetF[2][2];
214 fL3[1][0] = fGradDetF[0][1]*fGradDetF[2][0];
215 fL3[1][1] = fGradDetF[0][1]*fGradDetF[2][1];
216 fL3[1][2] = fGradDetF[0][1]*fGradDetF[2][2];
217 fL3[2][0] = fGradDetF[0][2]*fGradDetF[2][0];
218 fL3[2][1] = fGradDetF[0][2]*fGradDetF[2][1];
219 fL3[2][2] = fGradDetF[0][2]*fGradDetF[2][2];
221 fL4[0][0] = fGradDetF[1][0]*fGradDetF[0][0];
222 fL4[0][1] = fGradDetF[1][0]*fGradDetF[0][1];
223 fL4[0][2] = fGradDetF[1][0]*fGradDetF[0][2];
224 fL4[1][0] = fGradDetF[1][1]*fGradDetF[0][0];
225 fL4[1][1] = fGradDetF[1][1]*fGradDetF[0][1];
226 fL4[1][2] = fGradDetF[1][1]*fGradDetF[0][2];
227 fL4[2][0] = fGradDetF[1][2]*fGradDetF[0][0];
228 fL4[2][1] = fGradDetF[1][2]*fGradDetF[0][1];
229 fL4[2][2] = fGradDetF[1][2]*fGradDetF[0][2];
231 fL5[0][0] = fGradDetF[1][0]*fGradDetF[1][0];
232 fL5[0][1] = fGradDetF[1][0]*fGradDetF[1][1];
233 fL5[0][2] = fGradDetF[1][0]*fGradDetF[1][2];
234 fL5[1][0] = fGradDetF[1][1]*fGradDetF[1][0];
235 fL5[1][1] = fGradDetF[1][1]*fGradDetF[1][1];
236 fL5[1][2] = fGradDetF[1][1]*fGradDetF[1][2];
237 fL5[2][0] = fGradDetF[1][2]*fGradDetF[1][0];
238 fL5[2][1] = fGradDetF[1][2]*fGradDetF[1][1];
239 fL5[2][2] = fGradDetF[1][2]*fGradDetF[1][2];
241 fL6[0][0] = fGradDetF[1][0]*fGradDetF[2][0];
242 fL6[0][1] = fGradDetF[1][0]*fGradDetF[2][1];
243 fL6[0][2] = fGradDetF[1][0]*fGradDetF[2][2];
244 fL6[1][0] = fGradDetF[1][1]*fGradDetF[2][0];
245 fL6[1][1] = fGradDetF[1][1]*fGradDetF[2][1];
246 fL6[1][2] = fGradDetF[1][1]*fGradDetF[2][2];
247 fL6[2][0] = fGradDetF[1][2]*fGradDetF[2][0];
248 fL6[2][1] = fGradDetF[1][2]*fGradDetF[2][1];
249 fL6[2][2] = fGradDetF[1][2]*fGradDetF[2][2];
251 fL7[0][0] = fGradDetF[2][0]*fGradDetF[0][0];
252 fL7[0][1] = fGradDetF[2][0]*fGradDetF[0][1];
253 fL7[0][2] = fGradDetF[2][0]*fGradDetF[0][2];
254 fL7[1][0] = fGradDetF[2][1]*fGradDetF[0][0];
255 fL7[1][1] = fGradDetF[2][1]*fGradDetF[0][1];
256 fL7[1][2] = fGradDetF[2][1]*fGradDetF[0][2];
257 fL7[2][0] = fGradDetF[2][2]*fGradDetF[0][0];
258 fL7[2][1] = fGradDetF[2][2]*fGradDetF[0][1];
259 fL7[2][2] = fGradDetF[2][2]*fGradDetF[0][2];
261 fL8[0][0] = fGradDetF[2][0]*fGradDetF[1][0];
262 fL8[0][1] = fGradDetF[2][0]*fGradDetF[1][1];
263 fL8[0][2] = fGradDetF[2][0]*fGradDetF[1][2];
264 fL8[1][0] = fGradDetF[2][1]*fGradDetF[1][0];
265 fL8[1][1] = fGradDetF[2][1]*fGradDetF[1][1];
266 fL8[1][2] = fGradDetF[2][1]*fGradDetF[1][2];
267 fL8[2][0] = fGradDetF[2][2]*fGradDetF[1][0];
268 fL8[2][1] = fGradDetF[2][2]*fGradDetF[1][1];
269 fL8[2][2] = fGradDetF[2][2]*fGradDetF[1][2];
271 fL9[0][0] = fGradDetF[2][0]*fGradDetF[2][0];
272 fL9[0][1] = fGradDetF[2][0]*fGradDetF[2][1];
273 fL9[0][2] = fGradDetF[2][0]*fGradDetF[2][2];
274 fL9[1][0] = fGradDetF[2][1]*fGradDetF[2][0];
275 fL9[1][1] = fGradDetF[2][1]*fGradDetF[2][1];
276 fL9[1][2] = fGradDetF[2][1]*fGradDetF[2][2];
277 fL9[2][0] = fGradDetF[2][2]*fGradDetF[2][0];
278 fL9[2][1] = fGradDetF[2][2]*fGradDetF[2][1];
279 fL9[2][2] = fGradDetF[2][2]*fGradDetF[2][2];
282 STATE detF = (ux+1.)*(vy+1.)*(wz+1.) + vx*wy*uz + wx*uy*vz - wx*(vy+1.)*uz - (ux+1.)*wy*vz - vx*uy*(wz+1.);
285 cout <<
"\nDeterminante negativo!\n";
295 global[i][j][0] = fac1*
fL1[i][j];
296 global[i][j][1] = fac1*
fL2[i][j]+fac2*
fK2[i][j];
297 global[i][j][2] = fac1*
fL3[i][j]+fac2*
fK3[i][j];
298 global[i][j][3] = fac1*
fL4[i][j]+fac2*
fK4[i][j];
299 global[i][j][4] = fac1*
fL5[i][j];
300 global[i][j][5] = fac1*
fL6[i][j]+fac2*
fK6[i][j];
301 global[i][j][6] = fac1*
fL7[i][j]+fac2*
fK7[i][j];
302 global[i][j][7] = fac1*
fL8[i][j]+fac2*
fK8[i][j];
303 global[i][j][8] = fac1*
fL9[i][j];
305 global[i][i][0] += c*
fE1[i];
306 global[i][i][4] += c*
fE5[i];
307 global[i][i][8] += c*
fE9[i];
311 STATE gradJx[3],gradJy[3],gradJz[3];
312 STATE gradtrCx[3],gradtrCy[3],gradtrCz[3];
314 int nshape = phi.
Rows();
316 gradJx[i] = fGradDetF[0][i];
317 gradJy[i] = fGradDetF[1][i];
318 gradJz[i] = fGradDetF[2][i];
324 STATE *efptr = &ef(0,0);
325 REAL *dphiptr = &dphi(0,0);
326 int nrowek = ek.
Rows();
327 STATE *ekptr = &ek(0,0);
329 STATE kval[3] = {(fac2*gradJx[k]+c*gradtrCx[k]),(fac2*gradJy[k]+c*gradtrCy[k]),(fac2*gradJz[k]+c*gradtrCz[k])};
330 for(i=0; i<nshape; i++) {
331 efptr[i*3+k] += -weight*(dphiptr[3*i]*kval[0]+
332 dphiptr[3*i+1]*kval[1]+
333 dphiptr[3*i+2]*kval[2]);
338 for(i=0; i<nshape; i++) {
339 for(j=0; j<nshape; j++) {
340 ekptr[i*3+k+nrowek*(j*3+l)] += weight*(
341 dphiptr[3*i]*dphiptr[3*j]*global[k][l][0] + dphiptr[3*i]*dphiptr[1+3*j]*global[k][l][1] +
342 dphiptr[3*i]*dphiptr[2+3*j]*global[k][l][2] + dphiptr[1+3*i]*dphiptr[3*j]*global[k][l][3] +
343 dphiptr[1+3*i]*dphiptr[1+3*j]*global[k][l][4] + dphiptr[1+3*i]*dphiptr[2+3*j]*global[k][l][5] +
344 dphiptr[2+3*i]*dphiptr[3*j]*global[k][l][6] + dphiptr[2+3*i]*dphiptr[1+3*j]*global[k][l][7] +
345 dphiptr[2+3*i]*dphiptr[2+3*j]*global[k][l][8] );
356 if(!strcmp(
"Displacement6",name.c_str()))
return 1;
357 if(!strcmp(
"displacement",name.c_str()))
return 2;
358 if(!strcmp(
"Solution",name.c_str()))
return 2;
359 if(!strcmp(
"Derivate",name.c_str()))
return 3;
360 if(!strcmp(
"VonMises",name.c_str()))
return 4;
361 if(!strcmp(
"POrder",name.c_str()))
return 10;
367 if(var == 1)
return 6;
368 if(var == 2)
return 3;
369 if(var == 3)
return 9;
370 if(var == 4)
return 1;
371 if(var == 10)
return 1;
377 if(var == 1) Solout.
Resize(6,0.);
378 if(var == 2) Solout.
Resize(3,0.);
379 if(var == 1|| var == 2) {
387 for(
int i=0;i<3;i++) {
388 Solout[k++] = DSol(0,i);
389 Solout[k++] = DSol(1,i);
390 Solout[k++] = DSol(2,i);
406 STATE ux,uy,uz,vx,vy,vz,wx,wy,wz;
419 STATE J = (ux+1.)*(vy+1.)*(wz+1.) + vx*wy*uz + wx*uy*vz - wx*(vy+1.)*uz - (ux+1.)*wy*vz - vx*uy*(wz+1.);
425 STATE trsigma = sigmaF(0,0)+ sigmaF(1,1)+ sigmaF(2,2);
426 S = sigmaF - STATE(trsigma/3.0)*I;
429 for(i=0;i<3;i++)
for(j=0;j<3;j++) J2 += S(i,j)* S(i,j);
430 Solout[0] =
sqrt(3.0*J2);
450 values[1] =
pow(sol[0] - u_exact[0],(STATE)2.0);
451 values[1] +=
pow(sol[1] - u_exact[1],(STATE)2.0);
452 values[1] +=
pow(sol[2] - u_exact[2],(STATE)2.0);
456 for(
int i=0;i<3;i++) {
457 values[2] +=
pow(dsol[k++] - du_exact(0,i),(STATE)2.0);
458 values[2] +=
pow(dsol[k++] - du_exact(1,i),(STATE)2.0);
459 values[2] +=
pow(dsol[k++] - du_exact(2,i),(STATE)2.0);
462 values[0] = values[1]+values[2];
472 const STATE BIGNUMBER = 1.e12;
474 const int phr = phi.
Rows();
477 v2[0] = bc.
Val2()(0,0);
478 v2[1] = bc.
Val2()(1,0);
479 v2[2] = bc.
Val2()(2,0);
482 int numbersol = data.
sol.
size();
483 if (numbersol != 1) {
489 for(in = 0 ; in < phr; in++) {
490 ef(3*in+0,0) += BIGNUMBER * v2[0] * phi(in,0) * weight;
491 ef(3*in+1,0) += BIGNUMBER * v2[1] * phi(in,0) * weight;
492 ef(3*in+2,0) += BIGNUMBER * v2[2] * phi(in,0) * weight;
494 for (jn = 0 ; jn < phr; jn++) {
495 ek(3*in+0,3*jn+0) += BIGNUMBER * phi(in,0) * phi(jn,0) * weight;
496 ek(3*in+1,3*jn+1) += BIGNUMBER * phi(in,0) * phi(jn,0) * weight;
497 ek(3*in+2,3*jn+2) += BIGNUMBER * phi(in,0) * phi(jn,0) * weight;
503 for(in = 0 ; in < phi.
Rows(); in++) {
504 ef(3*in+0,0) += v2[0] * phi(in,0) * weight;
505 ef(3*in+1,0) += v2[1] * phi(in,0) * weight;
506 ef(3*in+2,0) += v2[2] * phi(in,0) * weight;
510 for(in = 0 ; in < phi.
Rows(); in++) {
511 ef(3*in+0,0) += v2[0] * phi(in,0) * weight;
512 ef(3*in+1,0) += v2[1] * phi(in,0) * weight;
513 ef(3*in+2,0) += v2[2] * phi(in,0) * weight;
514 for(jn=0; jn<phi.
Rows(); jn++)
516 for(idf=0; idf<3; idf++)
for(jdf=0; jdf<3; jdf++)
518 ek(3*in+idf,3*jn+jdf) += bc.
Val1()(
idf,jdf);
524 for(in = 0 ; in < phr; in++) {
525 ef(3*in+0,0) += BIGNUMBER * (0. - data.
sol[0][0]) * v2[0] * phi(in,0) * weight;
526 ef(3*in+1,0) += BIGNUMBER * (0. - data.
sol[0][1]) * v2[1] * phi(in,0) * weight;
527 ef(3*in+2,0) += BIGNUMBER * (0. - data.
sol[0][2]) * v2[2] * phi(in,0) * weight;
528 for (jn = 0 ; jn < phr; jn++) {
529 ek(3*in+0,3*jn+0) += BIGNUMBER * phi(in,0) * phi(jn,0) * weight * v2[0];
530 ek(3*in+1,3*jn+1) += BIGNUMBER * phi(in,0) * phi(jn,0) * weight * v2[1];
531 ek(3*in+2,3*jn+2) += BIGNUMBER * phi(in,0) * phi(jn,0) * weight * v2[2];
537 for(in = 0; in < 3; in ++)
538 v2[in] = - ( v1(in,0) * data.
normal[0] +
539 v1(in,1) * data.
normal[1] +
540 v1(in,2) * data.
normal[2] );
543 for(in = 0 ; in < phi.
Rows(); in++) {
544 ef(3*in+0,0) += v2[0] * phi(in,0) * weight;
545 ef(3*in+1,0) += v2[1] * phi(in,0) * weight;
546 ef(3*in+2,0) += v2[2] * phi(in,0) * weight;
552 PZError <<
"TPZElastitity3D::ContributeBC error - Wrong boundary condition type" << std::endl;
626 inline int ith(
const int i,
const int j)
631 void TPZMatHyperElastic::ContributeEnergy(
TPZVec<REAL> &x,
633 FADFADREAL &U, REAL weight)
672 void TPZMatHyperElastic::ContributeBCEnergy(
TPZVec<REAL> & x,
virtual void Solution(TPZMaterialData &data, int var, TPZVec< STATE > &Solout)
Returns the solution associated with the var index based on the finite element approximation.
virtual void Execute(const TPZVec< REAL > &x, TPZVec< TVar > &f, TPZFMatrix< TVar > &df)
Performs function computation.
TPZManVector< REAL, 3 > normal
normal to the element at the integration point
TPZManVector< REAL, 3 > x
value of the coordinate at the integration point
clarg::argBool bc("-bc", "binary checkpoints", false)
Implements a hyper elasticity material.
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.
virtual int VariableIndex(const std::string &name)
Returns the variable index associated with the name.
TPZGradSolVec dsol
vector of the derivatives of the solution at the integration point
virtual ~TPZMatHyperElastic()
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...
virtual int VariableIndex(const std::string &name) override
TPZFMatrix< STATE > & Val2(int loadcase=0)
TPZFNMatrix< 220, REAL > phi
vector of shapefunctions (format is dependent on the value of shapetype)
std::string Name() override
Returns the name of the material.
TPZFNMatrix< 660, REAL > dphix
values of the derivative of the shape functions
virtual void Print(std::ostream &out=std::cout)
Prints out the data associated with the material.
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...
This abstract class defines the behaviour which each derived class needs to implement.
Contains declaration of TPZElementMatrix struct which associates an element matrix with the coeficien...
Contains the TPZBndCond class which implements a boundary condition for TPZMaterial objects...
int64_t size() const
Returns the number of elements of the vector.
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...
Contains TPZMatrixclass which implements full matrix (using column major representation).
#define DebugStop()
Returns a message to user put a breakpoint in.
Contains the TPZMatHyperElastic class which implements a hyper elasticity material.
This class defines the boundary condition for TPZMaterial objects.
int64_t Rows() const
Returns number of rows.
virtual void Print(std::ostream &out) override
Prints out the data associated with the material.
expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ expr_ sqrt
TPZMatHyperElastic(int nummat, STATE e, STATE mu, STATE nu=-1., STATE lambda=-1., STATE coef1=-1., STATE coef2=-1., STATE coef3=-1.)
TPZFMatrix< STATE > & Val1()
int32_t Hash(std::string str)
Contains TPZMatrix<TVar>class, root matrix class.
virtual void Solution(TPZVec< STATE > &Sol, TPZFMatrix< STATE > &DSol, TPZFMatrix< REAL > &axes, int var, TPZVec< STATE > &Solout) override
int ClassId() const override
Unique identifier for serialization purposes.
TPZFlopCounter pow(const TPZFlopCounter &orig, const TPZFlopCounter &xp)
Returns the power and increments the counter of the power.
virtual int NSolutionVariables(int var)
Returns the number of variables associated with the variable indexed by var.
virtual int ClassId() const override
Define the class id associated with the class.
TPZAutoPointer< TPZFunction< STATE > > fForcingFunction
Pointer to forcing function, it is the right member at differential equation.
virtual int NSolutionVariables(int var) override
Returns the number of variables associated with the variable indexed by var.
TPZSolVec sol
vector of the solutions at the integration point
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...
#define PZError
Defines the output device to error messages and the DebugStop() function.
virtual int NStateVariables() const override
Returns the number of state variables associated with the material.