ME563-HW8

code/constitutive laws/Cauchy_type_selection.m

@@ -19,7 +19,9 @@
     case 7
          Cauchy = stress7(kinematics,properties,dim);
     case 8
-         Cauchy = stress8(kinematics,properties,dim);        
+         Cauchy = stress8(kinematics,properties,dim);  
+    case 9
+         Cauchy = stress9(kinematics,properties,cons);
     case 17
          %-----------------------------------------------------------------
          % Select internal variables at a particular Gauss point (igauss)

code/constitutive laws/elasticity_modulus_selection.m

@@ -19,7 +19,9 @@
     case 7
          c_tensor = ctens7(kinematics,properties,dimension);
     case 8
-         c_tensor = ctens8(kinematics,properties,dimension);            
+         c_tensor = ctens8(kinematics,properties,dimension);  
+    case 9
+         c_tensor = ctens9(kinematics,properties,cons);
     case 17
          %-----------------------------------------------------------------
          % Select internal variables at a particular Gauss point (igauss)

code/constitutive laws/elasticity_tensor/ctens9.m

@@ -0,0 +1,13 @@
+%--------------------------------------------------------------------------
+% Evaluates the constitutive tensor (in Voigt notation) for material type 1.
+%--------------------------------------------------------------------------
+function c = ctens9(kinematics,properties,cons)
+mu         = properties(2);
+lambda     = properties(3);
+J          = kinematics.J;
+lambda_    = lambda/J;
+mu_        = (mu - lambda*log(J))/J;
+c          = lambda_*cons.IDENTITY_TENSORS.c1 + mu_*cons.IDENTITY_TENSORS.c2;
+end
+
+

code/constitutive laws/stress/stress9.m

@@ -0,0 +1,16 @@
+%--------------------------------------------------------------------------
+% Evaluates the Cauchy stress tensor for material type 1.
+%--------------------------------------------------------------------------
+function Cauchy = stress9(kinematics,properties,cons)
+mu              = properties(2);
+lambda          = properties(3);
+J               = kinematics.J;
+b               = kinematics.b;
+%isochoric Mooney-Rivilin
+k = lambda+2*mu/3;
+k/mu;
+C               =kinematics.F'*kinematics.F;
+C_bar           =(J^(-2/3)*C);
+S               =J^(-2/3)*((-1/3)*((mu*cons.I+2*mu*cons.I)*C_bar^-1+(mu+mu*cons.I)*cons.I-mu*C_bar));
+Cauchy          =kinematics.F*S*kinematics.F';
+end

job_folder/1elt_3D_incompressible_neohookean/1elt_3D_incompressible_neohookean.dat

@@ -12,11 +12,11 @@ hexa8
 1
 1 1 1 4 3 2 5 8 7 6
 1
-1 5 1.0 .3356 16.67
+1 9 1.0 30 16.67
 0 4 0 0.0 0.0 0.0
 5 2 0.005
 6 2 0.005
 7 2 0.005
 8 2 0.005
-1 20.0 1.0 25 1.e-10 0.0 0.0 1 1 3 2
+1 20.0 1.0 25 1.e-3 0.0 0.0 1 1 3 2