[benchmarks] - SCMRUTH For Cell3

dG3D/benchmarks/torchSCMRUTH-Cell10/model.py

@@ -3,21 +3,12 @@ import torch
 from gmshpy import *
 from dG3Dpy import*
 # from dG3DpyDebug import*
-#script to launch PBC problem with a python script
-#DEFINE MICRO PROBLEM
-
-cell=["SC_MRU_T","SC_MRU_I","SMRU"] #available cells
-InpType =  cell[0] #choose the cell type
-factorStep= 1 #number of steps x 20 during reloading between points B and C
-
 
+#DEFINE MICRO PROBLEM
 
 # material law
-# matnum1 = 1 # unique number of law
-# rho = 7850e-9 # Bulk mass
 lnum1 = 11
 rho = 1e-9
-# Path = '../Model/Bounds_GS'
 Path = 'SCU_Bounds_GS'
 checkpoint = torch.load(Path)
 [XGmin, XGmax] = checkpoint['Input_Bounds']
@@ -25,7 +16,7 @@ checkpoint = torch.load(Path)
 [dTmin, dTmax] = checkpoint['Time_Bounds']
 print(XGmin, XGmax,YSmin, YSmax, dTmin, dTmax)
 
-if(InpType == "SC_MRU_T"):
+
 Radiusmean = (XGmax[0] + XGmin[0]) / 2.0
 Radiusstd = (XGmax[0] - XGmin[0]) / 2.0
 
@@ -47,46 +38,10 @@ if(InpType == "SC_MRU_T"):
 EYZmean = (XGmax[6] + XGmin[6]) / 2.0
 EYZstd = (XGmax[6] - XGmin[6]) / 2.0
   
-  
-  # model= '../Model/DInpFullModel.pt'
 model= 'Model.pt'
 kinematicinput=1
 extranorm=True
 doubleinput=True
-elif(InpType == "SC_MRU_I"):
-  model= "../Model/IncrementModel.pt"
-  kinematicinput=4
-  EXXmean = 0.0
-  EXXstd = 1.0
-  EXYmean = 0.0
-  EXYstd = 1.0
-  EYYmean = 0.0
-  EYYstd = 1.0
-  EYZmean = 0.0
-  EYZstd = 1.0
-  EZZmean = 0.0
-  EZZstd = 1.0
-  EZXmean = 0.0
-  EZXstd = 1.0
-  extranorm=False
-  doubleinput=False
-else:
-  EXXmean = (XGmax[0]+XGmin[0])/2.0
-  EXXstd = (XGmax[0]-XGmin[0])/2.0
-  EXYmean = (XGmax[1]+XGmin[1])/2.0
-  EXYstd = (XGmax[1]-XGmin[1])/2.0
-  EYYmean = (XGmax[2]+XGmin[2])/2.0
-  EYYstd = (XGmax[2]-XGmin[2])/2.0
-  EYZmean = 0.0
-  EYZstd = 1.0
-  EZZmean = 0.0
-  EZZstd = 1.0
-  EZXmean = 0.0
-  EZXstd = 1.0
-  model= '../Model/FullModel.pt'#'../Model/model.pt'#
-  kinematicinput=1
-  extranorm=True
-  doubleinput=False
 
 SXXmean = (YSmax[0]+YSmin[0])/2.0
 SXXstd = (YSmax[0]-YSmin[0])/2.0
@@ -109,9 +64,9 @@ SYZstd = (YSmax[5]-YSmin[5]) / 2.0
 
 
 numberOfInput = 8
-numInternalVars = 150#75
+numInternalVars = 37#75
 initial_h = 0.0
-macromat1 = torchANNBasedDG3DMaterialLaw(lnum1, rho, numberOfInput, numInternalVars, model, EXXmean, EXXstd, EXYmean, EXYstd, EYYmean, EYYstd, EYZmean, EYZstd, EZZmean, EZZstd, EZXmean, EZXstd, SXXmean, SXXstd, SXYmean, SXYstd, SYYmean, SYYstd, SYZmean, SYZstd, SZZmean,  SZZstd, SZXmean, SZXstd,True,1e-3)
+macromat1 = torchANNBasedDG3DMaterialLaw(lnum1, rho, numberOfInput, numInternalVars, model, EXXmean, EXXstd, EXYmean, EXYstd, EYYmean, EYYstd, EYZmean, EYZstd, EZZmean, EZZstd, EZXmean, EZXstd, SXXmean, SXXstd, SXYmean, SXYstd, SYYmean, SYYstd, SYZmean, SYZstd, SZZmean,  SZZstd, SZXmean, SZXstd,True,1e-5)
 
 
 macromat1.setInitialHValue(initial_h)
@@ -125,7 +80,7 @@ macromat1.setNumExtraInput(2)
 macromat1.setTimeArg(1)
 
 # Set initial vcalues for extra inputs in order Radius, cellsize, time
-macromat1.setInitialExtraInput(0.6)
+macromat1.setInitialExtraInput(0.5)
 macromat1.setInitialExtraInput(0.0)
 
 # Set normalization values for extra inputs in order Radius, cellsize, time
@@ -141,8 +96,8 @@ sol = 2  # Gmm=0 (default) Taucs=1 PETsc=2
 soltype =1 # StaticLinear=0 (default) StaticNonLinear=1
 nstep = 100  # number of step (used only if soltype=1)
 ftime =1.   # Final time (used only if soltype=1)
-tol=1.e-4  # relative tolerance for NR scheme (used only if soltype=1)
-nstepArch=1 # Number of step between 2 archiving (used only if soltype=1)
+tol=1.e-5  # relative tolerance for NR scheme (used only if soltype=1)
+nstepArch=5 # Number of step between 2 archiving (used only if soltype=1)
 
 
 # creation of  macro part Domain
@@ -150,7 +105,6 @@ nstepArch=1 # Number of step between 2 archiving (used only if soltype=1)
 dim =3
 beta1 = 100;
 fullDG = False;
-
 averageStrainBased = True
 
 # non DG domain
@@ -167,42 +121,18 @@ mysolver.addDomain(macrodomain1)
 mysolver.addMaterialLaw(macromat1)
 mysolver.Scheme(soltype)
 mysolver.Solver(sol)
-mysolver.snlData(nstep,ftime,tol,tol/10.)
-mysolver.snlManageTimeStep(50,5,2.,10)
+mysolver.snlData(nstep,ftime,tol)
+# mysolver.snlManageTimeStep(50,5,2.,10)
 mysolver.stepBetweenArchiving(nstepArch)
 
 
 # boundary condition
-# mysolver.displacementBC("Face",11,2,0.0)
-
-# mysolver.displacementBC("Edge",2,1,0.0)
-# mysolver.displacementBC("Edge",5,0,0.0)
 mysolver.displacementBC("Face",551,1,0.0)
 mysolver.displacementBC("Face",554,0,0.0)
 mysolver.displacementBC("Face",555,2,0.0)
 mysolver.displacementBC("Face",553,1,-5.0)
 
-# fct = PiecewiseLinearFunction()
-# fct.put(0.,0.)
-# fct.put(0.1,0.025)
-# fct.put(0.2,0.02)
-# fct.put(0.4,0.035)
-# fct.put(0.6,0.03)
-# fct.put(0.8,0.05)
-# fct.put(1.,0.0)
-# mysolver.displacementBC("Edge",4,1,fct)
-# unit = 1
-# mysolver.setNumStepTimeInterval(0.1,20*unit*1)
-# mysolver.setNumStepTimeInterval(0.2,3*unit)
-# mysolver.setNumStepTimeInterval(0.4,20*unit*factorStep)
-# mysolver.setNumStepTimeInterval(0.6,3*unit)
-# mysolver.setNumStepTimeInterval(0.8,50*unit)
-# mysolver.setNumStepTimeInterval(1.,20*unit)
-
-
 # archivage
-
-
 mysolver.internalPointBuildView("Green-Lagrange_xx",IPField.GL_XX, 1, nstepArch);
 mysolver.internalPointBuildView("Green-Lagrange_yy",IPField.GL_YY, 1, nstepArch);
 mysolver.internalPointBuildView("Green-Lagrange_zz",IPField.GL_ZZ, 1, nstepArch);
@@ -215,16 +145,7 @@ mysolver.internalPointBuildView("sig_zz",IPField.SIG_ZZ, 1, nstepArch);
 mysolver.internalPointBuildView("sig_xy",IPField.SIG_XY, 1, nstepArch);
 mysolver.internalPointBuildView("sig_yz",IPField.SIG_YZ, 1, nstepArch);
 mysolver.internalPointBuildView("sig_xz",IPField.SIG_XZ, 1, nstepArch);
-#mysolver.internalPointBuildView("sig_VM",IPField.SVM, 1, 1);
-#mysolver.internalPointBuildView("Green-Lagrange equivalent strain",IPField.GL_EQUIVALENT_STRAIN, 1, 1);
-#mysolver.internalPointBuildView("Equivalent plastic strain",IPField.PLASTICSTRAIN, 1, 1);
-
-#mysolver.internalPointBuildView("Deformation energy",IPField.DEFO_ENERGY, 1, 1);
-#mysolver.internalPointBuildView("Plastic energy",IPField.PLASTIC_ENERGY, 1, 1);
-
-# mysolver.archivingForceOnPhysicalGroup("Edge",2,1)
-# mysolver.archivingForceOnPhysicalGroup("Edge",2,0)
-# mysolver.archivingNodeDisplacement(6,1)
+mysolver.internalPointBuildView("sig_VM",IPField.SVM, 1, 1);
 
 mysolver.archivingForceOnPhysicalGroup("Face",551, 1)
 mysolver.archivingNodeDisplacementOnPhysical(2, 553, 1, 1)

dG3D/benchmarks/torchSCMRUTH-Cell3/model.py

+#coding-Utf-8-*-
+import torch
+from gmshpy import *
+from dG3Dpy import*
+# from dG3DpyDebug import*
+
+#DEFINE MICRO PROBLEM
+
+# material law
+lnum1 = 11
+rho = 1e-9
+Path = 'SCU_Bounds_GS'
+checkpoint = torch.load(Path)
+[XGmin, XGmax] = checkpoint['Input_Bounds']
+[YSmin, YSmax] = checkpoint['Output_Bounds']
+[dTmin, dTmax] = checkpoint['Time_Bounds']
+print(XGmin, XGmax,YSmin, YSmax, dTmin, dTmax)
+
+
+Radiusmean = (XGmax[0] + XGmin[0]) / 2.0
+Radiusstd = (XGmax[0] - XGmin[0]) / 2.0
+
+EXXmean = (XGmax[1]+XGmin[1])/2.0
+EXXstd = (XGmax[1]-XGmin[1])/2.0
+
+EYYmean = (XGmax[2]+XGmin[2])/2.0
+EYYstd = (XGmax[2]-XGmin[2])/2.0
+
+EZZmean = (XGmax[3]+XGmin[3])/2.0
+EZZstd = (XGmax[3]-XGmin[3])/2.0
+
+EXYmean = (XGmax[4]+XGmin[4])/2.0
+EXYstd = (XGmax[4]-XGmin[4])/2.0
+
+EZXmean = (XGmax[5] + XGmin[5]) / 2.0
+EZXstd = (XGmax[5] - XGmin[5]) / 2.0
+
+EYZmean = (XGmax[6] + XGmin[6]) / 2.0
+EYZstd = (XGmax[6] - XGmin[6]) / 2.0
+
+model= 'Model.pt'
+kinematicinput=1
+extranorm=True
+doubleinput=True
+
+SXXmean = (YSmax[0]+YSmin[0])/2.0
+SXXstd = (YSmax[0]-YSmin[0])/2.0
+
+SYYmean = (YSmax[1]+YSmin[1])/2.0
+SYYstd = (YSmax[1]-YSmin[1])/2.0
+
+SZZmean = (YSmax[2]+YSmin[2])/2.0
+SZZstd = (YSmax[2]-YSmin[2])/2.0
+
+
+SXYmean = (YSmax[3]+YSmin[3])/2.0
+SXYstd = (YSmax[3]-YSmin[3])/2.0
+
+SZXmean = (YSmax[4]+YSmin[4]) / 2.0
+SZXstd = (YSmax[4]-YSmin[4]) / 2.0
+
+SYZmean = (YSmax[5]+YSmin[5]) / 2.0
+SYZstd = (YSmax[5]-YSmin[5]) / 2.0
+
+
+numberOfInput = 8
+numInternalVars = 37
+initial_h = 0.0
+macromat1 = torchANNBasedDG3DMaterialLaw(lnum1, rho, numberOfInput, numInternalVars, model, EXXmean, EXXstd, EXYmean, EXYstd, EYYmean, EYYstd, EYZmean, EYZstd, EZZmean, EZZstd, EZXmean, EZXstd, SXXmean, SXXstd, SXYmean, SXYstd, SYYmean, SYYstd, SYZmean, SYZstd, SZZmean,  SZZstd, SZXmean, SZXstd,True,1e-5)
+
+
+macromat1.setInitialHValue(initial_h)
+macromat1.setKinematicInput(kinematicinput)
+macromat1.setNeedExtraNorm(extranorm)
+macromat1.setDoubleInput(doubleinput)
+
+
+macromat1.setNumExtraInput(2)
+# Set the location of time in the extra array
+macromat1.setTimeArg(1)
+
+# Set initial vcalues for extra inputs in order Radius, cellsize, time
+macromat1.setInitialExtraInput(0.6)
+macromat1.setInitialExtraInput(0.0)
+
+# Set normalization values for extra inputs in order Radius, cellsize, time
+macromat1.setNormExtraInp(Radiusmean, Radiusstd)
+macromat1.setNormExtraInp(dTmin, dTmax)
+
+
+macrogeofile="rubics.geo"
+macromeshfile="rubics.msh" # name of mesh file
+
+# solver
+sol = 2  # Gmm=0 (default) Taucs=1 PETsc=2
+soltype =1 # StaticLinear=0 (default) StaticNonLinear=1
+nstep = 100  # number of step (used only if soltype=1)
+ftime =1.   # Final time (used only if soltype=1)
+tol=1.e-5  # relative tolerance for NR scheme (used only if soltype=1)
+nstepArch=5 # Number of step between 2 archiving (used only if soltype=1)
+
+
+# creation of  macro part Domain
+
+dim =3
+beta1 = 100
+fullDG = False
+averageStrainBased = True
+
+# non DG domain
+nfield1 = 11
+macrodomain1 = dG3DDomain(10,nfield1,0,lnum1,fullDG,dim)
+macrodomain1.stabilityParameters(beta1)
+
+# creation of Solver
+mysolver = nonLinearMechSolver(10)
+mysolver.loadModel(macromeshfile)
+
+mysolver.addDomain(macrodomain1)
+mysolver.addMaterialLaw(macromat1)
+mysolver.Scheme(soltype)
+mysolver.Solver(sol)
+mysolver.snlData(nstep,ftime,tol)
+# mysolver.snlManageTimeStep(50,5,2.,10)
+mysolver.stepBetweenArchiving(nstepArch)
+
+
+# boundary condition
+mysolver.displacementBC("Face",554,0,0.0)
+mysolver.displacementBC("Face",551,1,0.0)
+mysolver.displacementBC("Face",555,2,0.0)
+mysolver.displacementBC("Face",553,1,-5.0)
+
+
+# archivage
+mysolver.internalPointBuildView("Green-Lagrange_xx",IPField.GL_XX, 1, nstepArch);
+mysolver.internalPointBuildView("Green-Lagrange_yy",IPField.GL_YY, 1, nstepArch);
+mysolver.internalPointBuildView("Green-Lagrange_zz",IPField.GL_ZZ, 1, nstepArch);
+mysolver.internalPointBuildView("Green-Lagrange_xy",IPField.GL_XY, 1, nstepArch);
+mysolver.internalPointBuildView("Green-Lagrange_yz",IPField.GL_YZ, 1, nstepArch);
+mysolver.internalPointBuildView("Green-Lagrange_xz",IPField.GL_XZ, 1, nstepArch);
+mysolver.internalPointBuildView("sig_xx",IPField.SIG_XX, 1, nstepArch);
+mysolver.internalPointBuildView("sig_yy",IPField.SIG_YY, 1, nstepArch);
+mysolver.internalPointBuildView("sig_zz",IPField.SIG_ZZ, 1, nstepArch);
+mysolver.internalPointBuildView("sig_xy",IPField.SIG_XY, 1, nstepArch);
+mysolver.internalPointBuildView("sig_yz",IPField.SIG_YZ, 1, nstepArch);
+mysolver.internalPointBuildView("sig_xz",IPField.SIG_XZ, 1, nstepArch);
+mysolver.internalPointBuildView("sig_VM",IPField.SVM, 1, 1);
+
+mysolver.archivingForceOnPhysicalGroup("Face",551, 1)
+mysolver.archivingNodeDisplacementOnPhysical(2, 553, 1, 1)
+
+# solve
+mysolver.solve()
+
+# check = TestCheck()
+# check.equal(34.22952279032552,mysolver.getArchivedForceOnPhysicalGroup("Edge", 2, 1),1.e-4)
+
+
+
+

dG3D/benchmarks/torchSCMRUTH-Cell3/rubics.geo

+SetFactory('OpenCASCADE');
+
+lx = 20;
+ly = 20;
+lz = 20;
+
+nx = 1;
+ny = 1;
+nz = 1;
+
+//nx11 = 2;
+//nx12 = 2;
+//nx13 = 2;
+
+//nz11 = 1;
+
+toll = (lx/nx)/2;
+e = 2;
+Point(1)={0, 0, 0, e};
+Point(2)={lx/nx, 0, 0, e};
+Point(3)={lx/nx, ly/ny, 0, e};
+Point(4)={0, ly/ny, 0, e};
+
+
+Line(1)={1, 2};
+Line(2)={2, 3};
+Line(3)={3, 4};
+Line(4)={4, 1};
+
+
+Line Loop(100)={1, 2, 3, 4};
+Surface(100)={100};
+
+Extrude {0, 0, lz/nz} {
+  Surface{100}; Layers{nz}; Recombine;
+}
+
+
+
+Transfinite Line{1,3}=nx+1 Using Progression 1;
+Transfinite Line{2,4}=ny+1 Using Progression 1;
+
+Transfinite Surface {100};
+Recombine Surface{100};
+Recombine Volume{1};
+
+
+Physical Surface(551) = {101};
+Physical Surface(552) = {102};
+Physical Surface(553) = {103};
+Physical Surface(554) = {104};
+Physical Surface(555) = {100};
+Physical Surface(556) = {105};
+
+Physical Volume(11)={1};

dG3D/benchmarks/torchSCMRUTH-Cell3/rubics.msh

+$MeshFormat
+4.1 0 8
+$EndMeshFormat
+$Entities
+8 12 6 1
+1 0 0 0 0 
+2 20 0 0 0 
+3 20 20 0 0 
+4 0 20 0 0 
+5 0 0 20 0 
+6 20 0 20 0 
+7 20 20 20 0 
+8 0 20 20 0 
+1 -9.99999993922529e-08 -1e-07 -1e-07 20.0000001 1e-07 1e-07 0 2 1 -2 
+2 19.9999999 -9.99999993922529e-08 -1e-07 20.0000001 20.0000001 1e-07 0 2 2 -3 
+3 -9.99999993922529e-08 19.9999999 -1e-07 20.0000001 20.0000001 1e-07 0 2 3 -4 
+4 -1e-07 -9.99999993922529e-08 -1e-07 1e-07 20.0000001 1e-07 0 2 4 -1 
+5 -1e-07 -1e-07 -9.99999993922529e-08 1e-07 1e-07 20.0000001 0 2 1 -5 
+6 19.9999999 -1e-07 -9.99999993922529e-08 20.0000001 1e-07 20.0000001 0 2 2 -6 
+7 -9.99999993922529e-08 -1e-07 19.9999999 20.0000001 1e-07 20.0000001 0 2 5 -6 
+8 19.9999999 19.9999999 -9.99999993922529e-08 20.0000001 20.0000001 20.0000001 0 2 3 -7 
+9 19.9999999 -9.99999993922529e-08 19.9999999 20.0000001 20.0000001 20.0000001 0 2 6 -7 
+10 -1e-07 19.9999999 -9.99999993922529e-08 1e-07 20.0000001 20.0000001 0 2 4 -8 
+11 -9.99999993922529e-08 19.9999999 19.9999999 20.0000001 20.0000001 20.0000001 0 2 7 -8 
+12 -1e-07 -9.99999993922529e-08 19.9999999 1e-07 20.0000001 20.0000001 0 2 8 -5 
+100 -1.000000100000001 -1.000000100000001 -1e-07 21.0000001 21.0000001 1e-07 1 555 4 1 2 3 4 
+101 -9.99999993922529e-08 -1e-07 -9.99999993922529e-08 20.0000001 1e-07 20.0000001 1 551 4 5 7 -6 -1 
+102 19.9999999 -9.99999993922529e-08 -9.99999993922529e-08 20.0000001 20.0000001 20.0000001 1 552 4 6 9 -8 -2 
+103 -9.99999993922529e-08 19.9999999 -9.99999993922529e-08 20.0000001 20.0000001 20.0000001 1 553 4 8 11 -10 -3 
+104 -1e-07 -9.99999993922529e-08 -9.99999993922529e-08 1e-07 20.0000001 20.0000001 1 554 4 10 12 -5 -4 
+105 -1.000000100000001 -1.000000100000001 19.9999999 21.0000001 21.0000001 20.0000001 1 556 4 7 9 11 12 
+1 -1.000000100000001 -1.000000100000001 -9.99999993922529e-08 21.0000001 21.0000001 20.0000001 1 11 6 -101 -102 -103 -104 -100 105 
+$EndEntities
+$Nodes
+27 27 1 27
+0 1 0 1
+1
+0 0 0
+0 2 0 1
+2
+20 0 0
+0 3 0 1
+3
+20 20 0
+0 4 0 1
+4
+0 20 0
+0 5 0 1
+5
+0 0 20
+0 6 0 1
+6
+20 0 20
+0 7 0 1
+7
+20 20 20
+0 8 0 1
+8
+0 20 20
+1 1 0 1
+9
+10 0 0
+1 2 0 1
+10
+20 10 0
+1 3 0 1
+11
+10 20 0
+1 4 0 1
+12
+0 10 0
+1 5 0 1
+13
+0 0 10
+1 6 0 1
+14
+20 0 10
+1 7 0 1
+15
+10 0 20
+1 8 0 1
+16
+20 20 10
+1 9 0 1
+17
+20 10 20
+1 10 0 1
+18
+0 20 10
+1 11 0 1
+19
+10 20 20
+1 12 0 1
+20
+0 10 20
+2 100 0 1
+21
+10 10 0
+2 101 0 1
+22
+10 0 10
+2 102 0 1
+23
+20 10 10
+2 103 0 1
+24
+10 20 10
+2 104 0 1
+25
+0 10 10
+2 105 0 1
+26
+10 10 20
+3 1 0 1
+27
+10 10 10
+$EndNodes
+$Elements
+7 7 1 7
+2 100 10 1
+1 1 2 3 4 9 10 11 12 21 
+2 101 10 1
+2 1 2 6 5 9 14 15 13 22 
+2 102 10 1
+3 2 3 7 6 10 16 17 14 23 
+2 103 10 1
+4 3 4 8 7 11 18 19 16 24 
+2 104 10 1
+5 4 1 5 8 12 13 20 18 25 
+2 105 10 1
+6 5 6 7 8 15 17 19 20 26 
+3 1 12 1
+7 1 2 3 4 5 6 7 8 9 12 13 10 14 11 16 18 15 20 17 19 21 22 25 23 24 26 27 
+$EndElements