Skip to content
Snippets Groups Projects
Commit edf56129 authored by Ludovic Noels's avatar Ludovic Noels
Browse files

new test

parent 957c3119
Branches
No related tags found
No related merge requests found
Showing
with 1525 additions and 0 deletions
......@@ -256,6 +256,7 @@ add_subdirectory(torchSMRU)
add_subdirectory(torchUniaxial)
add_subdirectory(torchNonLocal)
add_subdirectory(testPBCWithReducedRVE)
add_subdirectory(stochDMN)
add_subdirectory(switchSystem)
add_subdirectory(switchSystem_quadDispLinNonlocal)
add_subdirectory(GenericPerturbationMaterialLaw)
......
# test file
set(PYFILE SimulationDMN.py)
set(FILES2DELETE
rve1_Level4_ShearRes.csv
)
add_cm3python_test(${PYFILE} "${FILES2DELETE}")
import torch.nn as nn
import torch.nn.functional as F
import pandas as pd
import random
from PIL import Image as im
import torch
from os import listdir
from os.path import isfile, join
import numpy as np
######## Strain (E_00,E_11,E_22, 2*E_01,2*E_02,2*E_12) ##############
class BasicMat(nn.Module):
def __init__(self,level,pos):
super(BasicMat,self).__init__()
self.N_mat = 2
self.level = level
self.pos = pos
self.N_Interface = 1
def FillMatrix_Ns_NT(self,angle):
Ns = torch.zeros(6*self.N_Interface, 3*self.N_Interface)
angle1 = torch.zeros(1) # for 2D
N0 = torch.cos(angle1) * torch.sin(2.0*np.pi*angle)
N1 = torch.cos(angle1) * torch.cos(2.0*np.pi*angle)
N2 = torch.sin(angle1) # for 2D
Ns[0, 0] = N0
Ns[1, 1] = N1
Ns[2, 2] = N2
Ns[3, 0] = N1
Ns[3, 1] = N0
Ns[4, 0] = N2
Ns[4, 2] = N0
Ns[5, 1] = N2
Ns[5, 2] = N1
return Ns
def forward(self, C, W1, Para_N):
W0 = 1.0 - W1
Para = Para_N[self.pos]
angle = Para_N[self.pos][0]
Ns = self.FillMatrix_Ns_NT(angle)
NT = torch.t(Ns)
K = C[0]/W0 + C[1]/W1
M_KNs = torch.einsum('ik,kj->ij', K, Ns)
M_NKN = torch.einsum('ik,kj->ij', NT, M_KNs)
invM = torch.linalg.inv(M_NKN)
M = torch.einsum('ik,kj->ij', invM, NT)
M = torch.einsum('ik,kj->ij', M, (C[1]-C[0]))
M = torch.einsum('ik,kj->ij', Ns, M)
M = torch.einsum('ik,kj->ij', (C[0]-C[1]), M)
Ch = C[0]*W0 + C[1]*W1 + M
return Ch
##################################################
##################################################
class Mat_Net_VF(nn.Module):
def __init__(self, Total_level:int, level:int, pos:int, void=False):
super(Mat_Net_VF,self).__init__()
self.N_mat = 2
self.Total_level = Total_level
self.level = level
self.N_Interface = 1
self.Dim = 2
self.void=void
if(self.level == 0):
self.pos = 0
self.locPos = 0
self.NTotal_Interface = np.power(self.N_mat, self.Total_level) -1
self.N_W = np.power(2, self.Total_level-1)-1
# create control matrix
self.N_entry = 9
else:
self.pos = pos
self.locPos = pos + 1- np.power(self.N_mat,self.level)
self.mat = []
level = self.level+1
for i in range(self.N_mat):
if(level==self.Total_level-1):
if(self.void):
self.mat.append(VoidMat(self.Total_level,level,2*self.pos+1+i))
else:
self.mat.append(BasicMat(level,2*self.pos+1+i))
else:
self.mat.append(Mat_Net_VF(self.Total_level,level,2*self.pos+1+i, self.void))
if(self.level == 0):
self.Para_N = nn.Parameter(6.0*torch.rand((self.NTotal_Interface,self.Dim-1))-3.0)
if(self.void):
self.Para_W = nn.Parameter(6.0*torch.rand((self.NTotal_Interface,self.N_mat))-3.0)
else:
self.Para_W = nn.Parameter(torch.rand((self.N_W,self.N_mat)))
else:
self.Para_N = None
self.Para_W = None
self.act = nn.Sigmoid()
def FillMatrix_Ns_NT(self,angle):
Ns = torch.zeros(6*self.N_Interface, 3*self.N_Interface)
angle1 = torch.zeros(1) # for 2D
N0 = torch.cos(angle1) * torch.sin(2.0*np.pi*angle)
N1 = torch.cos(angle1) * torch.cos(2.0*np.pi*angle)
N2 = torch.sin(angle1) # for 2D
Ns[0, 0] = N0
Ns[1, 1] = N1
Ns[2, 2] = N2
Ns[3, 0] = N1
Ns[3, 1] = N0
Ns[4, 0] = N2
Ns[4, 2] = N0
Ns[5, 1] = N2
Ns[5, 2] = N1
return Ns
def forward(self, C, Vf, Para_N=None, Para_W=None):
if(self.Para_N != None):
Para_N = self.act(self.Para_N)
Para_W = self.act(self.Para_W)
if(self.level == 0):
CHom = []
for i in range(C.size()[0]):
V1 = Vf[i]
W0tmp = Para_W[self.pos][0]
W1tmp = Para_W[self.pos][1]
angle = Para_N[self.pos][0]
Ns= self.FillMatrix_Ns_NT(angle)
NT = torch.t(Ns)
Vm = 1.0 - V1
W = (Vm*W0tmp + V1*W1tmp, Vm*(1.0-W0tmp) + V1*(1.0-W1tmp))
W1 = (V1*W1tmp/W[0], V1*(1.0-W1tmp)/W[1])
Cmat = []
for k in range(self.N_mat):
Cmat.append(self.mat[k](C[i],W1[k],Para_N, Para_W))
K = Cmat[0]/W[0] + Cmat[1]/W[1]
M_KNs = torch.einsum('ik,kj->ij', K, Ns)
M_NKN = torch.einsum('ik,kj->ij', NT, M_KNs)
invM = torch.linalg.inv(M_NKN)
M = torch.einsum('ik,kj->ij', invM, NT)
M = torch.einsum('ik,kj->ij', M, (Cmat[1]-Cmat[0]))
M = torch.einsum('ik,kj->ij', Ns, M)
M = torch.einsum('ik,kj->ij', (Cmat[0]-Cmat[1]), M)
Ch = Cmat[0]*W[0] + Cmat[1]*W[1] + M
CHom.append(Ch)
CHom = torch.stack(CHom, dim=0)
return CHom
else:
V1 = Vf
W0tmp = Para_W[self.pos][0]
W1tmp = Para_W[self.pos][1]
angle = Para_N[self.pos][0]
Ns = self.FillMatrix_Ns_NT(angle)
NT = torch.t(Ns)
Vm = 1.0 - V1
W = (Vm*W0tmp + V1*W1tmp, Vm*(1.0-W0tmp) + V1*(1.0-W1tmp))
W1 = (V1*W1tmp/W[0], V1*(1.0-W1tmp)/W[1])
Cmat = []
for k in range(self.N_mat):
if((not self.void) and self.level==self.Total_level-2):
Cmat.append(self.mat[k](C,W1[k],Para_N))
else:
Cmat.append(self.mat[k](C,W1[k],Para_N, Para_W))
K = Cmat[0]/W[0] + Cmat[1]/W[1]
M_KNs = torch.einsum('ik,kj->ij', K, Ns)
M_NKN = torch.einsum('ik,kj->ij', NT, M_KNs)
invM = torch.linalg.inv(M_NKN)
M = torch.einsum('ik,kj->ij', invM, NT)
M = torch.einsum('ik,kj->ij', M, (Cmat[1]-Cmat[0]))
M = torch.einsum('ik,kj->ij', Ns, M)
M = torch.einsum('ik,kj->ij', (Cmat[0]-Cmat[1]), M)
Ch = Cmat[0]*W[0] + Cmat[1]*W[1] + M
return Ch
def my_loss(output, target):
mask = torch.zeros((6,6))
mask[0,0] = mask[1,1]= mask[1,0]= mask[0,1]= 1.0
mask[3,3] = mask[0,3]= mask[3,0]= mask[1,3]= mask[3,1]= 1.0
norm = torch.norm((mask*target),dim=(1, 2))
loss = torch.norm((mask*(output - target)),dim=(1, 2))
loss = torch.mean(torch.div(loss,norm))
return loss
def getFiles(folder, index):
onlyfiles = [f for f in listdir(folder) if isfile(join(folder, f))]
allFiles = list()
for fname in onlyfiles:
if fname[:len(index)]==index:
allFiles.append(fname)
return allFiles
def TrainData(fn,N_add=0):
names = ['F_XX','F_XY','F_XZ','F_YX','F_YY','F_YZ','F_ZX','F_ZY','F_ZZ','P_XX','P_XY','P_XZ','P_YX','P_YY','P_YZ','P_ZX','P_ZY','P_ZZ']
dat = pd.read_csv(fn, sep=';')
strainPath = []
target = []
for i in range(len(dat['F_XX'])):
F = np.zeros(9)
P = np.zeros(9)
for j in range(9):
F[j] = dat[names[j]][i]
strainPath.append(F)
for j in range(9):
P[j] = dat[names[j+9]][i]
target.append(P)
if(N_add == 0):
return strainPath, target
else:
L_strainPath = []
L_target = []
for i in range(1,len(strainPath)):
d_strainPath = (strainPath[i] - strainPath[i-1])/(N_add+1)
d_target = (target[i] - target[i-1])/(N_add+1)
for j in range(N_add+1):
L_strainPath.append(strainPath[i-1] + j*d_strainPath)
L_target.append(target[i-1] + j*d_target)
L_strainPath.append(strainPath[-1])
L_target.append(target[-1])
return L_strainPath, L_target
File added
File added
File added
This diff is collapsed.
#coding-Utf-8-*-
import torch,csv
from gmshpy import *
from dG3Dpy import*
import pandas as pd
import numpy as np
from NNW_Tool import *
StressState = ['XZ']
Load = 'Shear'#'Tensile'#'UniStrain'#
prefix = './'
level = 4
Para = './rve1_Level'+str(level)+'.dat'
############################################
lawnum1 = 11 # unique number of law
lawnum2 = 12 # unique number of law
E = 4.668E3 #Hexcel Corporation, HexPly 8552, Epoxy matrix (180C/356F curing matrix), Product Data Sheet (2016).
nu = 0.39
K = E/3./(1.-2.*nu) # Bulk mudulus
mu =E/2./(1.+nu) # Shear mudulus
rho = 7850e-9 # Bulk mass
sy0 = 100. #MPa
h1 = 0.05*E
h2 = 0.0
harden1 = LinearExponentialJ2IsotropicHardening(1, sy0, h1, h2, 0.0)
law1 = J2LinearDG3DMaterialLaw(lawnum1,rho,E,nu,harden1)
law1.setStrainOrder(11)
#--------------------
rhoCF = 1750.e-9
youngCF = 12.99e3
youngACF = 231.e3
nuCF = 0.46
nu_minorCF = 0.3*youngCF/youngACF
GACF = 11.3e3
Ax = 0.0 # direction of anisotropy
Ay = 0.0
Az = 1.0
law2 = TransverseIsotropicDG3DMaterialLaw(lawnum2,rhoCF,youngCF,nuCF,youngACF,GACF,nu_minorCF,Ax,Ay,Az)
lawnum3 = 13
law3 = StochDMNDG3DMaterialLaw(lawnum3, rho, E, nu, Para, False)
law3.addLaw(law1)
law3.addLaw(law2)
comMat = computeMaterialLaw(law3)
########## Para_DMN ##################
names = ['F_XX','F_XY','F_XZ','F_YX','F_YY','F_YZ','F_ZX','F_ZY','F_ZZ','P_XX','P_XY','P_XZ','P_YX','P_YY','P_YZ','P_ZX','P_ZY','P_ZZ']
##############################################################################################################
def makeTest(NewPara, Stress_State, strainPath, Id):
comMat.reset_Parameter(law3,NewPara)
comMat.resetState(law3)
time=0.
TimeStep = 1.
N = int(len(strainPath))
ResArray = np.zeros([N,len(names)])
for ite in range(0,N):
# print(f'step = {ite}')
time = time + TimeStep
comMat.setTime(time,TimeStep)
comMat.nextStep()
# current plane strain
FxxCur = strainPath[ite][0]
FxyCur = strainPath[ite][1]
FxzCur = strainPath[ite][2]
FyxCur = strainPath[ite][3]
FyyCur = strainPath[ite][4]
FyzCur = strainPath[ite][5]
FzxCur = strainPath[ite][6]
FzyCur = strainPath[ite][7]
FzzCur = strainPath[ite][8]
# set new value
comMat.setDeformationGradient(0,0,FxxCur)
comMat.setDeformationGradient(0,1,FxyCur)
comMat.setDeformationGradient(0,2,FxzCur)
comMat.setDeformationGradient(1,0,FyxCur)
comMat.setDeformationGradient(1,1,FyyCur)
comMat.setDeformationGradient(1,2,FyzCur)
comMat.setDeformationGradient(2,0,FzxCur)
comMat.setDeformationGradient(2,1,FzyCur)
comMat.setDeformationGradient(2,2,FzzCur)
if(Load=='Tensile'):
comMat.computeStressState_PlaneStrainUniaxialStress(Stress_State)
# compute uniaxial tension in xy plane
else:
comMat.computeStressState()
#-----------------------------
Pstress = vectorDouble()
Fstrain = vectorDouble()
comMat.getStress(Pstress)
comMat.getDeformationGradient(Fstrain)
for k in range(9):
ResArray[ite][k] = [*Fstrain][k]
ResArray[ite][9+k] = [*Pstress][k]
df = pd.DataFrame(ResArray, columns=names)
df.to_csv(Id,sep=';', index=False)
return
######################################################################
if(Load=='Tensile'):
strain_file = './TensileTest_rve1.csv'
elif(Load=='Shear'):
strain_file = './ShearTest_rve1.csv'
else:
strain_file = './UniAxialStrain_rve1.csv'
strainPath, target = TrainData(strain_file, 0)
for r in range(1,2):
#f_Para = './DMNPara/Vf_rve'+str(r)+'_Level'+str(level)+'.dat'
#Id = prefix+'Vf_rve'+str(r)+'_Level'+str(level)+'_'+Load+'Res.csv'
f_Para = './rve'+str(r)+'_Level'+str(level)+'.dat'
Id = prefix+'rve'+str(r)+'_Level'+str(level)+'_'+Load+'Res.csv'
makeTest(f_Para, StressState[0], strainPath, Id)
## only used for rve6 with 300 training data
#f_Para = './DMNPara/rve6_Level'+str(level)+'_300Data.dat'
#Id = prefix+'rve6_Level'+str(level)+'_'+Load+'Res300.csv'
#makeTest(f_Para, StressState[0], strainPath, Id)
def checkEqual(ref, cur, tol):
if(abs(ref-cur)>abs(tol*ref)):
print("Error : reference current value ",ref," current value ",cur," relative error ", (ref-cur)/ref," tolerance ",tol)
fileres = csv.reader(open('rve1_Level4_ShearRes.csv'), delimiter=';')
res = list(fileres)
checkEqual(104.42366995900333,float(res[-1][10]),1e-5)
This diff is collapsed.
This diff is collapsed.
2 4
0.340573
0.587945
0.534264
0.712348
0.412997
0.282536
0.214489
0.424969
0.416157
0.761119
0.406089
0.760345
0.586586
0.540862
0.553134
0.360250
0.348736 0.841227
0.532382 0.701098
0.755990 0.038612
0.767671 0.000509
0.071756 0.986992
0.855184 0.010136
0.201100 0.987834
2 5
0.340573
0.435054
0.598946
0.293477
0.714262
0.610611
0.584946
0.750831
0.453267
0.558729
0.503317
0.702507
0.585985
0.394929
0.559591
0.632166
0.810744
0.288432
0.248234
0.762603
0.281320
0.315645
0.643404
0.584055
0.416766
0.507472
0.408590
0.425881
0.601115
0.489706
0.711104
0.376614
0.365614 0.775721
0.341979 0.812477
0.260299 0.928265
0.635131 0.717336
0.255349 0.920762
0.757339 0.548806
0.019977 0.949399
0.800015 0.402815
0.305050 0.687248
0.261533 0.980192
0.962335 0.009232
0.215508 0.998362
0.266344 0.786783
0.451069 0.157134
0.940725 0.067286
2 6
0.340573
0.670392
0.349202
0.205869
0.458847
0.567929
0.601584
0.368198
0.498060
0.354232
0.361268
0.516431
0.405086
0.463373
0.365639
0.289337
0.331127
0.616288
0.540821
0.487142
0.514288
0.492431
0.600788
0.370790
0.607889
0.364331
0.315987
0.591228
0.458033
0.459838
0.432231
0.803561
0.392439
0.345416
0.360358
0.267370
0.613706
0.588223
0.346349
0.606858
0.718428
0.511169
0.501894
0.663934
0.507077
0.508546
0.576033
0.401050
0.774856
0.389313
0.271908
0.747693
0.673013
0.432352
0.749366
0.659967
0.545777
0.492803
0.355272
0.368044
0.547305
0.547341
0.547728
0.555978
0.217304 0.342548
0.860130 0.277656
0.751733 0.325143
0.903281 0.635444
0.601140 0.217737
0.272466 0.928900
0.631753 0.127694
0.070277 0.996763
0.149054 0.965170
0.792052 0.395897
0.803078 0.661394
0.680806 0.013371
0.930098 0.032065
0.267104 0.972153
0.552167 0.466321
0.298368 0.995249
0.411056 0.294686
0.562516 0.700027
0.465957 0.798457
0.295143 0.717964
0.656197 0.517031
0.580394 0.576568
0.733499 0.375565
0.770791 0.600426
0.652680 0.313815
0.056106 0.953685
0.160853 0.937369
0.351752 0.961326
0.779485 0.035860
0.527706 0.527426
0.762613 0.423474
2 4
0.384944
0.584948
0.503443
0.748252
0.403994
0.218615
0.421130
0.498576
0.394711
0.791111
0.394519
0.705776
0.611666
0.502434
0.557018
0.397448
0.417159 0.856860
0.354686 0.605716
0.795007 0.040695
0.812561 0.000641
0.051145 0.987588
0.744331 0.014509
0.166711 0.971239
2 5
0.384944
0.442340
0.592711
0.315481
0.654925
0.651525
0.515867
0.724513
0.394723
0.512124
0.460559
0.700907
0.621376
0.266288
0.554892
0.576257
0.751609
0.310890
0.262603
0.735274
0.308016
0.298461
0.487173
0.625587
0.468074
0.513051
0.456222
0.516671
0.640112
0.436472
0.636921
0.451228
0.379362 0.747296
0.335752 0.838511
0.321838 0.885125
0.593690 0.711188
0.312680 0.881226
0.706591 0.660522
0.016644 0.943154
0.862911 0.333316
0.233325 0.738741
0.365329 0.977270
0.969599 0.012367
0.184044 0.993768
0.254350 0.833926
0.522969 0.245831
0.892519 0.070574
2 6
0.384944
0.622959
0.367865
0.205250
0.465293
0.665043
0.564167
0.407302
0.497519
0.254002
0.278925
0.544034
0.350614
0.473079
0.390938
0.284807
0.344515
0.613649
0.529853
0.491916
0.501545
0.489536
0.664025
0.425482
0.599941
0.422066
0.314330
0.638029
0.459067
0.418971
0.432285
0.712912
0.456606
0.399109
0.409412
0.297045
0.635883
0.613138
0.233662
0.591024
0.776835
0.491496
0.427045
0.743610
0.507230
0.507260
0.622402
0.298698
0.760550
0.451397
0.224799
0.768761
0.625199
0.452343
0.728577
0.639181
0.569914
0.576443
0.419528
0.312724
0.548105
0.548151
0.567024
0.515996
0.320343 0.266729
0.813113 0.332679
0.786955 0.269893
0.849106 0.667875
0.710728 0.177380
0.312604 0.904648
0.592230 0.099342
0.052768 0.998326
0.170623 0.944943
0.811451 0.362302
0.764127 0.669773
0.686951 0.012862
0.908512 0.049300
0.330653 0.959447
0.523422 0.440859
0.340094 0.994271
0.481473 0.393931
0.616591 0.707259
0.360578 0.759413
0.331040 0.710498
0.615746 0.622626
0.541163 0.610656
0.697229 0.374025
0.646286 0.507966
0.684406 0.303024
0.057231 0.929852
0.203867 0.910260
0.492533 0.931706
0.818467 0.030542
0.529233 0.526257
0.816741 0.341225
2 4
0.448363
0.601386
0.522366
0.705649
0.430365
0.263281
0.280493
0.458865
0.420845
0.764256
0.434218
0.662898
0.529725
0.514405
0.588303
0.335330
0.359378 0.847412
0.454168 0.676351
0.721434 0.035296
0.808670 0.000703
0.051023 0.977724
0.759887 0.015693
0.221034 0.985565
2 5
0.448363
0.441716
0.601546
0.309494
0.678817
0.630782
0.569725
0.709014
0.438990
0.563445
0.458912
0.693604
0.613675
0.329557
0.579493
0.485943
0.773215
0.281789
0.253076
0.746205
0.262404
0.313679
0.579597
0.609589
0.465401
0.479081
0.451130
0.464725
0.622266
0.458995
0.635754
0.526263
0.372182 0.758078
0.287791 0.856468
0.303629 0.900711
0.518116 0.714419
0.285078 0.899070
0.742119 0.615701
0.018485 0.925162
0.847603 0.336845
0.220292 0.752303
0.368034 0.979576
0.966267 0.010943
0.186813 0.996636
0.272853 0.798884
0.514752 0.208270
0.921091 0.066582
2 6
0.448363
0.650682
0.373386
0.187813
0.470795
0.723474
0.582636
0.395412
0.512342
0.277682
0.295795
0.520427
0.372881
0.480250
0.390171
0.256876
0.314378
0.616522
0.515842
0.489354
0.514458
0.467291
0.604457
0.419576
0.580524
0.365092
0.312485
0.640030
0.441331
0.414942
0.432281
0.758169
0.424229
0.408516
0.389888
0.289453
0.632769
0.604215
0.248385
0.596528
0.750837
0.507077
0.488018
0.722022
0.506563
0.498729
0.611439
0.327243
0.775952
0.464698
0.226144
0.742956
0.631090
0.454175
0.733479
0.607374
0.570759
0.586544
0.390893
0.353932
0.552341
0.552281
0.566596
0.547715
0.292740 0.288477
0.798518 0.338644
0.780193 0.282303
0.861793 0.648240
0.691739 0.171844
0.279414 0.918872
0.584330 0.119308
0.054498 0.997804
0.157268 0.945609
0.786840 0.399303
0.786956 0.691248
0.674279 0.013634
0.916743 0.038704
0.302997 0.978155
0.475024 0.469055
0.368793 0.993978
0.431485 0.326794
0.578640 0.724868
0.434641 0.746113
0.380913 0.693636
0.572693 0.546226
0.548040 0.613030
0.695708 0.362225
0.706616 0.438709
0.669329 0.307022
0.060054 0.944691
0.184004 0.908325
0.333349 0.954122
0.789351 0.032883
0.528066 0.528344
0.792140 0.366826
2 4
0.471991
0.633669
0.594408
0.725253
0.412541
0.269981
0.320577
0.494458
0.413248
0.753671
0.464664
0.728443
0.532299
0.601512
0.584737
0.402462
0.367571 0.834628
0.345244 0.703220
0.747307 0.053190
0.850459 0.000715
0.043072 0.990091
0.828717 0.012447
0.202405 0.987841
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment