new testY

dG3D/benchmarks/CMakeLists.txt

@@ -42,6 +42,7 @@ add_subdirectory(TransverseAnisotropicCZM)
 add_subdirectory(TransverseAnisotropicThermoMec)
 add_subdirectory(twoHole)
 add_subdirectory(crystalPlasticity)
+add_subdirectory(crystalPlasticityTaylor)
 add_subdirectory(incrementalSecantTest)
 add_subdirectory(incrementalSecantTest2)
 add_subdirectory(CrackCompositeFailure)

dG3D/benchmarks/crystalPlasticityTaylor/CMakeLists.txt

+# test file
+
+set(PYFILE cube.py)
+
+set(FILES2DELETE 
+  *.msh
+  *.csv
+)
+
+add_cm3python_test(${PYFILE} "${FILES2DELETE}")

dG3D/benchmarks/crystalPlasticityTaylor/cube.geo

+// Test case a SCB with a vertical load at its free extremity
+// Size
+
+//definition of unit
+mm = 1e-03;
+
+// volum fraction
+
+f=0.45;
+
+// characteristic size for fibe
+R=0.05*mm;             // fiber radius
+ly = 2.0*R*Sqrt(Pi/(2.0*Sqrt(3.0)*f));  //distance between fiber
+lx = ly*Sqrt(3.0)/2.0;
+
+
+x=4*ly;
+y=4*lx;
+z=0.5*mm;
+
+// Characteristic length
+Lc1=R/2.5;
+
+// definition of points
+Point(1) = { 0.0 , 0.0 , 0.0 , Lc1};
+Point(2) = {  x  , 0.0 , 0.0 , Lc1};
+Point(3) = {  x  ,  y  , 0.0 , Lc1};
+Point(4) = { 0.0 ,  y  , 0.0 , Lc1};
+Point(5) = { 0.0 , 0.0 , z , Lc1};
+Point(6) = {  x  , 0.0 , z , Lc1};
+Point(7) = {  x  ,  y  , z , Lc1};
+Point(8) = { 0.0 ,  y  , z , Lc1};
+
+// Line between points
+Line(1) = {1,2};
+Line(2) = {2,3};
+Line(3) = {3,4};
+Line(4) = {4,1};
+Line(5) = {5,6};
+Line(6) = {6,7};
+Line(7) = {7,8};
+Line(8) = {8,5};
+Line(9) = {1,5};
+Line(10)= {2,6};
+Line(11)= {3,7};
+Line(12)= {4,8};
+
+// Surface definition
+Line Loop(1) = {1,2,3,4};
+Line Loop(2) = {5,6,7,8};
+Line Loop(3) = {1,10,-5,-9};
+Line Loop(4) = {2,11,-6,-10};
+Line Loop(5) = {3,12,-7,-11};
+Line Loop(6) = {4,9,-8,-12};
+
+Plane Surface(1) = {1};
+Plane Surface(2) = {2};
+Plane Surface(3) = {3};
+Plane Surface(4) = {4};
+Plane Surface(5) = {5};
+Plane Surface(6) = {6};
+
+//VOlume
+
+Surface Loop(7) = {1,2,3,4,5,6};
+Volume(1) = {7};
+
+// Physical objects to applied BC and material
+Physical Surface(1234) = {1};
+Physical Surface(5678) = {2};
+Physical Surface(1265) = {3};
+Physical Surface(2376) = {4};
+Physical Surface(3487) = {5};
+Physical Surface(4158) = {6};
+Physical Line(12) = {1};
+Physical Line(23) = {2};
+Physical Line(34) = {3};
+Physical Line(41) = {4};
+Physical Line(56) = {5};
+Physical Line(67) = {6};
+Physical Line(78) = {7};
+Physical Line(85) = {8};
+Physical Line(15) = {9};
+Physical Line(26) = {10};
+Physical Line(37) = {11};
+Physical Line(48) = {12};
+
+Physical Point(1) ={1};
+Physical Point(2) ={2};
+Physical Point(3) ={3};
+Physical Point(4) ={4};
+Physical Point(5) ={5};
+Physical Point(6) ={6};
+Physical Point(7) ={7};
+Physical Point(8) ={8};
+
+Physical Volume(10) ={1};
+
+// define transfinite mesh
+Transfinite Line {1,2,3,4,5,6,7,8} = 2;
+Transfinite Line {9,10,11,12} = 2;
+Transfinite Surface {1,2,3,4,5,6} ;
+//Recombine Surface {1,2,3,4,5,6} ;
+Transfinite Volume {1};
+

dG3D/benchmarks/crystalPlasticityTaylor/cube.py

+#coding-Utf-8-*-
+from gmshpy import *
+from dG3Dpy import*
+#from dG3DpyDebug import*
+
+#script to launch beam problem with a python script
+
+# material law
+lawnum   = 1 # unique number of law
+rho      = 7850.
+temperature = 273.
+# geometry
+geofile="cube.geo"
+meshfile="cube.msh" # name of mesh file
+
+# solver
+sol = 2  # Gmm=0 (default) Taucs=1 PETsc=2
+soltype = 1 # StaticLinear=0 (default) StaticNonLinear=1
+nstep = 1000   # number of step (used only if soltype=1)
+ftime =1.   # Final time (used only if soltype=1)
+tol=1.e-6   # relative tolerance for NR scheme (used only if soltype=1)
+nstepArch=50 # Number of step between 2 archiving (used only if soltype=1)
+fullDg = 0 #O = CG, 1 = DG
+space1 = 0 # function space (Lagrange=0)
+
+
+#  compute solution and BC (given directly to the solver
+# creation of law
+law1 = crystalPlasticityDG3DMaterialLaw(lawnum,rho,temperature,"grains.inp","material_test.i01")
+
+# creation of ElasticField
+nfield = 10 # number of the field (physical number of surface)
+myfield1 = dG3DDomain(1000,nfield,space1,lawnum,fullDg)
+myfield1.stabilityParameters(30.)
+#myfield1.matrixByPerturbation(1,1,1,1e-8)
+# creation of Solver
+mysolver = nonLinearMechSolver(1000)
+mysolver.createModel(geofile,meshfile,3,2)
+#mysolver.loadModel(meshfile)
+mysolver.addDomain(myfield1)
+mysolver.addMaterialLaw(law1)
+mysolver.Scheme(soltype)
+mysolver.Solver(sol)
+mysolver.snlData(nstep,ftime,tol)
+mysolver.stepBetweenArchiving(nstepArch)
+# BC
+#shearing
+#mysolver.displacementBC("Face",1234,0,0.)
+#mysolver.displacementBC("Face",1234,1,0.)
+#mysolver.displacementBC("Face",1234,2,0.)
+#mysolver.displacementBC("Face",5678,1,0.00)
+#mysolver.displacementBC("Face",5678,0,0.00001)
+#mysolver.forceBC("Face",2376,2,100000000)
+#mysolver.forceBC("Face",4158,2,-100000000)
+#mysolver.forceBC("Face",5678,0,100000000)
+#tension along z
+mysolver.displacementBC("Face",1234,2,0.)
+d1=0.000005
+cyclicFunction1=cycleFunctionTime(0.,0.,ftime/4., d1/2., ftime/2., 0., 3.*ftime/4., d1/2., ftime, d1);
+mysolver.displacementBC("Face",5678,2,cyclicFunction1)
+mysolver.displacementBC("Face",2376,0,0.)
+mysolver.displacementBC("Face",1265,1,0.)
+
+
+mysolver.internalPointBuildView("svm",IPField.SVM, 1, 1)
+mysolver.internalPointBuildView("sig_xx",IPField.SIG_XX, 1, 1)
+mysolver.internalPointBuildView("sig_yy",IPField.SIG_YY, 1, 1)
+mysolver.internalPointBuildView("sig_zz",IPField.SIG_ZZ, 1, 1)
+mysolver.internalPointBuildView("sig_xy",IPField.SIG_XY, 1, 1)
+mysolver.internalPointBuildView("sig_yz",IPField.SIG_YZ, 1, 1)
+mysolver.internalPointBuildView("sig_xz",IPField.SIG_XZ, 1, 1)
+mysolver.internalPointBuildView("epl",IPField.PLASTICSTRAIN, 1, 1)
+mysolver.archivingForceOnPhysicalGroup("Face", 1234, 2)
+mysolver.archivingForceOnPhysicalGroup("Face", 5678, 2)
+
+mysolver.solve()
+
+check = TestCheck()
+check.equal(1.873222e-04,mysolver.getArchivedForceOnPhysicalGroup("Face", 5678, 2),1.e-6)
+
+

dG3D/benchmarks/crystalPlasticityTaylor/grains.inp

+*SOLID SECTION, ELSET= GRAIN1  , MATERIAL= GRAIN001
+*MATERIAL,NAME= GRAIN001
+*USER MATERIAL,CONSTANTS= 5
+      -0.879,      55.471,      57.955,         1.0,    1.000000,
+*DEPVAR
+ 22000
+*INITIAL CONDITIONS,TYPE=SOLUTION,USER
+*SOLID SECTION, ELSET= GRAIN2  , MATERIAL= GRAIN002
+*MATERIAL,NAME= GRAIN002
+*USER MATERIAL,CONSTANTS= 5
+      -0.879,      55.471,      57.955,         1.0,    1.000000,
+*DEPVAR
+ 23
+*INITIAL CONDITIONS,TYPE=SOLUTION,USER
+*SOLID SECTION, ELSET= GRAIN3  , MATERIAL= GRAIN003
+*MATERIAL,NAME= GRAIN003
+*USER MATERIAL,CONSTANTS= 5
+      -0.879,      55.471,      57.955,         1.0,    1.000000,
+*DEPVAR
+ 23
+*INITIAL CONDITIONS,TYPE=SOLUTION,USER
+*SOLID SECTION, ELSET= GRAIN4  , MATERIAL= GRAIN004
+*MATERIAL,NAME= GRAIN004
+*USER MATERIAL,CONSTANTS= 5
+      -0.879,      55.471,      57.955,         1.0,    1.000000,
+*DEPVAR
+ 23
+*INITIAL CONDITIONS,TYPE=SOLUTION,USER
+*SOLID SECTION, ELSET= GRAIN5  , MATERIAL= GRAIN005
+*MATERIAL,NAME= GRAIN005
+*USER MATERIAL,CONSTANTS= 5
+      -0.879,      55.471,      57.955,         1.0,    1.000000,
+*DEPVAR
+ 23
+*INITIAL CONDITIONS,TYPE=SOLUTION,USER
+*SOLID SECTION, ELSET= GRAIN6  , MATERIAL= GRAIN006
+*MATERIAL,NAME= GRAIN006
+*USER MATERIAL,CONSTANTS= 5
+      -0.879,      55.471,      57.955,         1.0,    1.000000,
+*DEPVAR
+ 23
+*INITIAL CONDITIONS,TYPE=SOLUTION,USER
+*SOLID SECTION, ELSET= GRAIN7  , MATERIAL= GRAIN007
+*MATERIAL,NAME= GRAIN007
+*USER MATERIAL,CONSTANTS= 5
+      -0.879,      55.471,      57.955,         1.0,    1.000000,
+*DEPVAR
+ 23
+*INITIAL CONDITIONS,TYPE=SOLUTION,USER
+*SOLID SECTION, ELSET= GRAIN8  , MATERIAL= GRAIN008
+*MATERIAL,NAME= GRAIN008
+*USER MATERIAL,CONSTANTS= 5
+      -0.879,      55.471,      57.955,         1.0,    1.000000,
+*DEPVAR
+ 23
+*INITIAL CONDITIONS,TYPE=SOLUTION,USER

dG3D/benchmarks/crystalPlasticityTaylor/material_test.i01

+     0          must be 1 if previous step exists in tmp.03
+     1					number of constitutive crystal symmetries 
+     -2					must not be 1 if cubic material. A negative value means that the Chaboche model must be activated
+    0.0     0.0    A    B
+   499.1d3  203.6d3  154.6d3		elastic coefficient C11  C12  C44
+     1					number of slip modes
+     1  					active slip modes
+     1  2  4  5  8			First basis of slip systems
+     1   Bishop-Hill slip systems for fcc
+    24   100.0  0.664e5			nb. of slip systems, exponent m, gam0
+     7					0= VOCE, 1= SWIFT
+ 2.74d-10  17.86  0.5  1.5  673  100.0 0.5311  0.129  4.5d12  1.157e8  25.92    260.7		burgers, grain-size, p, q, T (not used anymore), actVol measured in burgers^3, "beta" in F0=beta*mu*b3, "alpha" in sigma=alpha*mu*b*sqrt(rho), rho_0, k1, rho_sat/rho_0, sigmaHat
+    0  1 -1    1  1  1
+   -1  0  1    1  1  1
+    1 -1  0    1  1  1
+    0 -1 -1   -1 -1  1
+    1  0  1   -1 -1  1
+   -1  1  0   -1 -1  1
+    0  1 -1   -1  1  1
+    1  0  1   -1  1  1
+   -1 -1  0   -1  1  1
+    0 -1 -1    1 -1  1
+   -1  0  1    1 -1  1
+    1  1  0    1 -1  1
+    2 -1 -1    1  1  1
+   -1  2 -1    1  1  1
+   -1 -1  2    1  1  1
+   -2  1 -1   -1 -1  1
+    1 -2 -1   -1 -1  1
+    1  1  2   -1 -1  1
+   -2 -1 -1   -1  1  1
+    1  2 -1   -1  1  1
+    1 -1  2   -1  1  1
+    2  1 -1    1 -1  1
+   -1 -2 -1    1 -1  1
+   -1  1  2    1 -1  1
+    0
+    1.0					Must be <>1.0 if latent hardening (matrix included after this line)
+
+     0     1                 # of relaxation modes & # of grains in 1 cluster
+    -1                       type of interaction for the 2nd relaxation mode
+     0.0   1.0    0.0        relaxation tensor
+     1.0   0.0    0.0
+     0.0   0.0    0.0 
+    -1                       type of interaction for the 1st relaxation mode
+    -1.0   0.0    1.0        relaxation tensor
+     0.0   1.0    0.0
+     0.0   0.0    0.0 
+    1.d0    1.d0
+
+    -1                       type of interaction for the 2nd relaxation mode
+     0.0   0.0    0.0        relaxation tensor
+     0.0   0.0    0.0
+     0.0   0.0    1.0