Skip to content
Snippets Groups Projects
Commit d80af810 authored by Van Dung NGUYEN's avatar Van Dung NGUYEN
Browse files

add new cas test

parent 3435b970
No related branches found
No related tags found
1 merge request!20Master
......@@ -111,3 +111,4 @@ add_subdirectory(multiscaleCohesiveTest2D_fullDG_rotateRVE)
add_subdirectory(multiscaleEnhanceStrain)
add_subdirectory(nonLinearMixedBC_2D_DirectionFollowing)
add_subdirectory(interpolationPBC_2DShifted)
add_subdirectory(DG_PRI6_resetBC)
# test file
set(PYFILE twoHole.py)
set(FILES2DELETE
*.csv
disp*
stress*
previousScheme*
)
add_cm3python_mpi_test(4 ${PYFILE} "${FILES2DELETE}")
mm=1.0e-3;
n=1;
L=10*mm;
sl1=0.4*L/n;
Lx = 2*L;
Ly = 4*L;
Point(1)={0,0,0,sl1};
Point(2)={Lx,0,0,sl1};
Point(3)={Lx,Ly,0,sl1};
Point(4)={0,Ly,0,sl1};
Line(1)={1,2};
Line(2)={2,3};
Line(3)={3,4};
Line(4)={4,1};
x0=0.5*Lx/n; y0=0.5*Ly/n; r=0.2*L/n;
x=0.3*Lx;
y=0.4*Ly;
r = 0.2*L;
sl2=0.7*sl1;
p1=newp; Point(p1)={x-r,y,0,sl2};
p2=newp; Point(p2)={x,y+r,0,sl2};
p3=newp; Point(p3)={x+r,y,0,sl2};
p4=newp; Point(p4)={x,y-r,0,sl2};
pc=newp; Point(pc)={x,y,0,sl2};
c1 = newreg; Circle(c1) = {p1,pc,p2};
c2 = newreg; Circle(c2) = {p2,pc,p3};
c3 = newreg; Circle(c3) = {p3,pc,p4};
c4 = newreg; Circle(c4) = {p4,pc,p1};
l[1]=newreg; Line Loop(l[1]) = {c1,c2,c3,c4};
x=0.7*Lx;
y=0.6*Ly;
r = 0.2*L;
sl2=0.4*sl1;
p1=newp; Point(p1)={x-r,y,0,sl2};
p2=newp; Point(p2)={x,y+r,0,sl2};
p3=newp; Point(p3)={x+r,y,0,sl2};
p4=newp; Point(p4)={x,y-r,0,sl2};
pc=newp; Point(pc)={x,y,0,sl2};
c1 = newreg; Circle(c1) = {p1,pc,p2};
c2 = newreg; Circle(c2) = {p2,pc,p3};
c3 = newreg; Circle(c3) = {p3,pc,p4};
c4 = newreg; Circle(c4) = {p4,pc,p1};
l[2]=newreg; Line Loop(l[2]) = {c1,c2,c3,c4};
l[0]=newreg;
Line Loop(l[0])={1,2,3,4};
Plane Surface(11)={l[]};
//Recombine Surface{11};
Physical Surface(16) = {11};
Physical Line(17) = {1};
Physical Line(18) = {3};
Physical Point(19) = {4};
Physical Point(20) = {1};
//+
//Recombine Surface {11};
//+
Extrude {0, 0, 1e-3} {
Surface{11}; Layers{2}; Recombine;
}
//+
Physical Volume(83) = {1};
//+
Physical Surface(84) = {37};
//+
Physical Surface(85) = {45};
This diff is collapsed.
#coding-Utf-8-*-
from gmshpy import *
from dG3Dpy import*
#script to launch beam problem with a python script
lawnum1 = 1 # unique number of law
rho = 7850
young = 28.9e9
nu = 0.3
sy0 = 150.e6
h = young/50.
# geometry
meshfile="twoHole.msh" # name of mesh file
# solver
sol = 2 # Gmm=0 (default) Taucs=1 PETsc=2
soltype = 1 # StaticLinear=0 (default) StaticNonLinear=1
nstep = 25 # 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=1 # Number of step between 2 archiving (used only if soltype=1)
harden = LinearExponentialJ2IsotropicHardening(1, sy0, h, 0., 10.)
law1 = J2LinearDG3DMaterialLaw(lawnum1,rho,young,nu,harden)
# creation of ElasticField
beta1 = 1e2
fullDG = True;
averageStrainBased = False
myfield1 = dG3DDomain(1000,83,0,lawnum1,fullDG,3)
myfield1.stabilityParameters(beta1)
myfield1.averageStrainBased(averageStrainBased)
# creation of Solver
mysolver = nonLinearMechSolver(1000)
mysolver.loadModel(meshfile)
mysolver.addDomain(myfield1)
mysolver.addMaterialLaw(law1)
mysolver.Scheme(soltype)
mysolver.Solver(sol)
mysolver.snlData(nstep,ftime,tol)
mysolver.stepBetweenArchiving(nstepArch)
mysolver.options("-ksp_type preonly -pc_type lu -pc_factor_mat_solver_package mumps")
# BC
#mysolver.displacementBC("Volume",83,2,0.)
mysolver.displacementBC("Face",84,0,0.)
mysolver.displacementBC("Face",84,1,0.)
mysolver.displacementBC("Face",84,2,0.)
mysolver.displacementBC("Face",85,0,0.)
mysolver.displacementBC("Face",85,1,5e-4)
mysolver.displacementBC("Face",85,2,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.internalPointBuildView("Damage",IPField.DAMAGE, 1, 1)
mysolver.archivingForceOnPhysicalGroup("Face", 84, 1)
mysolver.archivingNodeDisplacement(19,1,1)
t1 = mysolver.solve()
check = TestCheck()
check.equal(-2.839032e+03,mysolver.getArchivedForceOnPhysicalGroup("Face", 84, 1),1.e-4)
nstep = 5 # 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)
mysolver.snlData(nstep,ftime,tol)
fct = LinearFunctionTime(0,5e-4,ftime,0.)
mysolver.resetBoundaryConditions();
mysolver.displacementBC("Face",84,0,0.)
mysolver.displacementBC("Face",84,1,0.)
mysolver.displacementBC("Face",84,2,0.)
mysolver.displacementBC("Face",85,0,0.)
mysolver.displacementBC("Face",85,1,fct)
mysolver.displacementBC("Face",85,2,0.)
mysolver.solve()
check = TestCheck()
check.equal(2.570988e+03,mysolver.getArchivedForceOnPhysicalGroup("Face", 84, 1),1.e-4)
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment