DNN Gel Constitutive Law

SetDirectory[NotebookDirectory[]]

Define exact constitutive law

Import 2d flux definitions

NotebookEvaluate[Directory[] <> "/../../FluxDefinitions2d.nb"];

Import averaged Hessian

<<Hav.mx;

Define reduced Jacobian that gives the momentum flux

M=Table[H[[k]].U[[k]],{k,d}];
Map[TableForm[#]&,M]

Replacement of (Subscript[∇, 𝔽])ℙ with averaged Hessian components

Prules = Flatten[Table[(Subscript[P, i,j])->Hav[[i,j]],{i,(d^2)},{j,(d^2)}]];

MP =Table[M[[k,i,j]] /. Prules,{k,1,d,1},{i,1,d,1},{j,1,d,1}];

Manipulate[MP[[k,i,j]] ,{k,1,d,1},{i,1,d,1},{j,1,d,1}]

Import microscopic data

Ψdata1=Flatten[Transpose[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-1/StrainEnergyData.mat"],{1,2,3,5,4}],3];
BadData1=Table[Position[Ψdata1[[i]],0.],{i,1,Length[Ψdata1]}];
Fdata1=Flatten[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-1/DeformationGradientData.mat"],3];
Pdata1=Flatten[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-1/FirstPiolaKirchhoffData.mat"],3];
Fdata1=Flatten[Table[Delete[Fdata1[[i]],BadData1[[i]]],{i,1,Length[Fdata1]}],1];
Pdata1=Flatten[Table[Delete[Pdata1[[i]],BadData1[[i]]],{i,1,Length[Pdata1]}],1];
Ψdata1=Flatten[Table[Delete[Ψdata1[[i]],BadData1[[i]]],{i,1,Length[Ψdata1]}],1];

Ψdata2=Flatten[Transpose[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-2/StrainEnergyData.mat"],{1,2,3,5,4}],3];
BadData2=Table[Position[Ψdata2[[i]],0.],{i,1,Length[Ψdata2]}];
Fdata2=Flatten[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-2/DeformationGradientData.mat"],3];
Pdata2=Flatten[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-2/FirstPiolaKirchhoffData.mat"],3];
Fdata2=Flatten[Table[Delete[Fdata2[[i]],BadData2[[i]]],{i,1,Length[Fdata2]}],1];
Pdata2=Flatten[Table[Delete[Pdata2[[i]],BadData2[[i]]],{i,1,Length[Pdata2]}],1];
Ψdata2=Flatten[Table[Delete[Ψdata2[[i]],BadData2[[i]]],{i,1,Length[Ψdata2]}],1];

Ψdata4=Flatten[Transpose[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-4/StrainEnergyData.mat"],{1,2,3,5,4}],3];
BadData4=Table[Position[Ψdata4[[i]],0.],{i,1,Length[Ψdata4]}];
Fdata4=Flatten[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-4/DeformationGradientData.mat"],3];
Pdata4=Flatten[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-4/FirstPiolaKirchhoffData.mat"],3];
Fdata4=Flatten[Table[Delete[Fdata4[[i]],BadData4[[i]]],{i,1,Length[Fdata4]}],1];
Pdata4=Flatten[Table[Delete[Pdata4[[i]],BadData4[[i]]],{i,1,Length[Pdata4]}],1];
Ψdata4=Flatten[Table[Delete[Ψdata4[[i]],BadData4[[i]]],{i,1,Length[Ψdata4]}],1];

Ψdata5=Flatten[Transpose[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-5/StrainEnergyData.mat"],{1,2,3,5,4}],3];
BadData5=Table[Position[Ψdata5[[i]],0.],{i,1,Length[Ψdata5]}];
Fdata5=Flatten[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-5/DeformationGradientData.mat"],3];
Pdata5=Flatten[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-5/FirstPiolaKirchhoffData.mat"],3];
Fdata5=Flatten[Table[Delete[Fdata5[[i]],BadData5[[i]]],{i,1,Length[Fdata5]}],1];
Pdata5=Flatten[Table[Delete[Pdata5[[i]],BadData5[[i]]],{i,1,Length[Pdata5]}],1];
Ψdata5=Flatten[Table[Delete[Ψdata5[[i]],BadData5[[i]]],{i,1,Length[Ψdata5]}],1];

Ψdata6=Flatten[Transpose[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-6/StrainEnergyData.mat"],{1,2,3,5,4}],3];
BadData6=Table[Position[Ψdata6[[i]],0.],{i,1,Length[Ψdata6]}];
Fdata6=Flatten[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-6/DeformationGradientData.mat"],3];
Pdata6=Flatten[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-6/FirstPiolaKirchhoffData.mat"],3];
Fdata6=Flatten[Table[Delete[Fdata6[[i]],BadData6[[i]]],{i,1,Length[Fdata6]}],1];
Pdata6=Flatten[Table[Delete[Pdata6[[i]],BadData6[[i]]],{i,1,Length[Pdata6]}],1];
Ψdata6=Flatten[Table[Delete[Ψdata6[[i]],BadData6[[i]]],{i,1,Length[Ψdata6]}],1];

Ψdata7=Flatten[Transpose[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-7/StrainEnergyData.mat"],{1,2,3,5,4}],3];
BadData7=Table[Position[Ψdata7[[i]],0.],{i,1,Length[Ψdata7]}];
Fdata7=Flatten[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-7/DeformationGradientData.mat"],3];
Pdata7=Flatten[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-7/FirstPiolaKirchhoffData.mat"],3];
Fdata7=Flatten[Table[Delete[Fdata7[[i]],BadData7[[i]]],{i,1,Length[Fdata7]}],1];
Pdata7=Flatten[Table[Delete[Pdata7[[i]],BadData7[[i]]],{i,1,Length[Pdata7]}],1];
Ψdata7=Flatten[Table[Delete[Ψdata7[[i]],BadData7[[i]]],{i,1,Length[Ψdata7]}],1];

Ψdata8=Flatten[Transpose[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-8/StrainEnergyData.mat"],{1,2,3,5,4}],3];
BadData8=Table[Position[Ψdata8[[i]],0.],{i,1,Length[Ψdata8]}];
Fdata8=Flatten[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-8/DeformationGradientData.mat"],3];
Pdata8=Flatten[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-8/FirstPiolaKirchhoffData.mat"],3];
Fdata8=Flatten[Table[Delete[Fdata8[[i]],BadData8[[i]]],{i,1,Length[Fdata8]}],1];
Pdata8=Flatten[Table[Delete[Pdata8[[i]],BadData8[[i]]],{i,1,Length[Pdata8]}],1];
Ψdata8=Flatten[Table[Delete[Ψdata8[[i]],BadData8[[i]]],{i,1,Length[Ψdata8]}],1];

Ψdata9=Flatten[Transpose[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-9/StrainEnergyData.mat"],{1,2,3,5,4}],3];
BadData9=Table[Position[Ψdata9[[i]],0.],{i,1,Length[Ψdata9]}];
Fdata9=Flatten[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-9/DeformationGradientData.mat"],3];
Pdata9=Flatten[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-9/FirstPiolaKirchhoffData.mat"],3];
Fdata9=Flatten[Table[Delete[Fdata9[[i]],BadData9[[i]]],{i,1,Length[Fdata9]}],1];
Pdata9=Flatten[Table[Delete[Pdata9[[i]],BadData9[[i]]],{i,1,Length[Pdata9]}],1];
Ψdata9=Flatten[Table[Delete[Ψdata9[[i]],BadData9[[i]]],{i,1,Length[Ψdata9]}],1];

Ψdata10=Flatten[Transpose[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-10/StrainEnergyData.mat"],{1,2,3,5,4}],3];
BadData10=Table[Position[Ψdata10[[i]],0.],{i,1,Length[Ψdata10]}];
Fdata10=Flatten[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-10/DeformationGradientData.mat"],3];
Pdata10=Flatten[Import["/home/brian/Research/BrianAdam/ActinMesh/data/LAMMPS Hyperelastic Analysis/HyperelasticData_03_08_2021/network-10/FirstPiolaKirchhoffData.mat"],3];
Fdata10=Flatten[Table[Delete[Fdata10[[i]],BadData10[[i]]],{i,1,Length[Fdata10]}],1];
Pdata10=Flatten[Table[Delete[Pdata10[[i]],BadData10[[i]]],{i,1,Length[Pdata10]}],1];
Ψdata10=Flatten[Table[Delete[Ψdata10[[i]],BadData10[[i]]],{i,1,Length[Ψdata10]}],1];

Fdata0=Join[Fdata1,Fdata2,Fdata4,Fdata5,Fdata6,Fdata7,Fdata8,Fdata9,Fdata10];
Pdata0=Join[Pdata1,Pdata2,Pdata4,Pdata5,Pdata6,Pdata7,Pdata8,Pdata9,Pdata10];
Ψdata=Join[Ψdata1,Ψdata2,Ψdata4,Ψdata5,Ψdata6,Ψdata7,Ψdata8,Ψdata9,Ψdata10];

Fdata=Flatten[#]&/@Fdata0;
Pdata=Flatten[#]&/@Pdata0;
nData=Length[Fdata];

FPdata=Table[Fdata[[i]]->Pdata[[i]],{i,nData}];
FΨdata=Table[Fdata[[i]]->Ψdata[[i]],{i,nData}];

Bootstrapping

SampleLabels=Table[i,{i,1,nData}];
Remaining=SampleLabels;
nTestData=Ceiling[0.1*nData];

Remaining=SampleLabels;
TestSamples=Table[{},{i,1,9}];
Clear[i,j]
For[j=1,j≤9,j++,
For[i=1,i≤nTestData,i++,
samplesLeft=Length[Remaining];
testIdx=RandomInteger[{1,samplesLeft}];
testSample=Remaining[[testIdx]];
AppendTo[TestSamples[[j]],{testSample}];
Remaining=Delete[Remaining,{testIdx}];
];
];
Clear[i,j]

Define training data

Define ranges of displacement gradient ∂u/∂a, and increment

minDuDa=-0.2; maxDuDa=0.2; δDuDa=0.1;

Construct 𝔽 input data set

Fdata = Flatten[Table[{FL1,FL2,FL3,FL4,FR1,FR2,FR3,FR4},
{FL1,1+minDuDa,1+maxDuDa, δDuDa},
{FL2,minDuDa,maxDuDa, δDuDa},
{FL3,minDuDa,maxDuDa, δDuDa},
{FL4,1+minDuDa,1+maxDuDa, δDuDa},
{FR1,1+minDuDa,1+maxDuDa, δDuDa},
{FR2,minDuDa,maxDuDa, δDuDa},
{FR3,minDuDa,maxDuDa, δDuDa},
{FR4,1+minDuDa,1+maxDuDa, δDuDa}],7];

Define material parameters

matRules = {AL->0.,DL->0.,muL->2.,AR->0.,DR->0.,muR->2.};

Construct M output data

kDir=1;
Mk=Flatten[MP[[kDir]] /. matRules];

Mdata=Map[ (Mk /. 
{(Subscript[F, L,1])->#[[1]], (Subscript[F, L,2])->#[[2]], (Subscript[F, L,3])->#[[3]], (Subscript[F, L,4])->#[[4]], 
(Subscript[F, R,1])->#[[5]], (Subscript[F, R,2])->#[[6]], (Subscript[F, R,3])->#[[7]], (Subscript[F, R,4])->#[[8]]})&, Fdata];

Dimensions[Mdata]

Construct training data set

nData=Length[Mdata]

FMdata = Table[Fdata[[i]]->Mdata[[i]],{i,nData}];

Networks with single-step dimension reduction

Linear, least-squares approximation

nLayersF=1;  nLayersM=1;
nlmF = Table[{nIn},{n,nLayersF}];
nlmM =Table[{nOut},{n,nLayersM}];
nlmL=NetChain[Flatten[{nlmF,nlmM}]]

nlmLTrained = NetTrain[nlmL,FMdata]

dk=1000;
errors=Table[Norm[nlmLTrained[Fdata[[k]]] - Mdata[[k]]],{k,1,nData,dk}];
ListPlot[errors,PlotLabel->dk Mean[errors]/nData]

Single layer DNN

FdataTrain=Delete[Fdata,TestSamples[[1]]];
PdataTrain=Delete[Pdata,TestSamples[[1]]];
FPdataTrain=Delete[FPdata,TestSamples[[1]]];

nLayersF=1;  nLayersP=1;
nlmF = Table[{9,LogisticSigmoid},{n,nLayersF}];
nlmM =Table[{9,LogisticSigmoid},{n,nLayersP}];
nlmSL=NetChain[Flatten[{nlmF,nlmM}]]

nlmSLTrained = NetTrain[nlmSL,FPdataTrain]

dk=100;
errors=Table[Norm[nlmSLTrained[FdataTrain[[k]]] - PdataTrain[[k]]],{k,1,Length[FdataTrain],dk}];
ListPlot[errors,PlotLabel->dk Mean[errors]/Length[FdataTrain]]

FdataTest=Extract[Fdata,TestSamples[[1]]];
PdataTest=Extract[Pdata,TestSamples[[1]]];
FPdataTest=Extract[FPdata,TestSamples[[1]]];

dk=10;
errors=Table[Norm[nlmSLTrained[FdataTest[[k]]] - PdataTest[[k]]],{k,1,Length[FdataTest],dk}];
ListPlot[errors,PlotLabel->dk Mean[errors]/Length[FdataTest]]

Two layer DNN

FdataTrain=Delete[Fdata,TestSamples[[2]]];
PdataTrain=Delete[Pdata,TestSamples[[2]]];
FPdataTrain=Delete[FPdata,TestSamples[[2]]];

nLayersF=2;  nLayersP=2;
nlmF = Table[{9,LogisticSigmoid},{n,nLayersF}];
nlmM =Table[{9,LogisticSigmoid},{n,nLayersP}];
nlmDL=NetChain[Flatten[{nlmF,nlmM}]]

nlmDLTrained = NetTrain[nlmDL,FPdataTrain]

dk=100;
errors=Table[Norm[nlmDLTrained[FdataTrain[[k]]] - PdataTrain[[k]]],{k,1,Length[PdataTrain],dk}];
ListPlot[errors,PlotLabel->dk Mean[errors]/Length[PdataTrain]]

FdataTest=Extract[Fdata,TestSamples[[2]]];
PdataTest=Extract[Pdata,TestSamples[[2]]];
FPdataTest=Extract[FPdata,TestSamples[[2]]];

dk=10;
errors=Table[Norm[nlmSLTrained[FdataTest[[k]]] - PdataTest[[k]]],{k,1,Length[FdataTest],dk}];
ListPlot[errors,PlotLabel->dk Mean[errors]/Length[FdataTest]]

Three layer DNN

FdataTrain=Delete[Fdata,TestSamples[[3]]];
ΨdataTrain=Delete[Ψdata,TestSamples[[3]]];
FΨdataTrain=Delete[FΨdata,TestSamples[[3]]];

nLayersF=3;  nLayersΨ=3;
nlmF = Table[{9,LogisticSigmoid},{n,nLayersF}];
nlmM =Table[{1,LogisticSigmoid},{n,nLayersΨ}];
nlmTL=NetChain[Flatten[{nlmF,nlmM}]]

nlmTLTrained = NetTrain[nlmTL,FΨdataTrain]

dk=100;
errors=Table[Norm[nlmTLTrained[FdataTrain[[k]]] - ΨdataTrain[[k]]]/Abs[ΨdataTrain[[k]]],{k,1,Length[FdataTrain],dk}];
ListPlot[errors,PlotLabel->dk Mean[errors]/Length[FdataTrain]]

Four layer DNN

nLayersF=4;  nLayersM=4;
nlmF = Table[{nIn,LogisticSigmoid},{n,nLayersF}];
nlmM =Table[{nOut,LogisticSigmoid},{n,nLayersM}];
nlmQL=NetChain[Flatten[{nlmF,nlmM}]]

nlmQLTrained = NetTrain[nlmQL,FMdata]

dk=500;
errors=Table[Norm[nlmQLTrained[Fdata[[k]]] - Mdata[[k]]],{k,1,nData,dk}];
ListPlot[errors,PlotLabel->dk Mean[errors]/nData]

Networks with multiple-step dimension reduction

Three layers per dimension reduction DNN

nLayersF=3;  nLayersT=3; nLayersM=3;
nlmF = Table[{nIn,LogisticSigmoid},{n,nLayersF}];
nlmT = Table[Flatten[Table[{n,LogisticSigmoid},nLayersT]],{n,nIn,nOut,-1}];
nlmM =Table[{nOut,LogisticSigmoid},{n,nLayersM}];
nlmT3L=NetChain[Flatten[{nlmF,nlmT,nlmM}]]

nlmT3LTrained = NetTrain[nlmT3L,FMdata]

dk=1000;
errors=Table[Norm[nlmT3LTrained[Fdata[[k]]] - Mdata[[k]]],{k,1,nData,dk}];
ListPlot[errors,PlotLabel->dk Mean[errors]/nData]

Two layer DNN

nLayersF=2;  nLayersM=2;
nlmF = Table[{nIn,LogisticSigmoid},{n,nLayersF}];
nlmM =Table[{nOut,LogisticSigmoid},{n,nLayersM}];
nlmDL=NetChain[Flatten[{nlmF,nlmM}]]

nlmDLTrained = NetTrain[nlmDL,FMdata]

dk=500;
errors=Table[Norm[nlmDLTrained[Fdata[[k]]] - Mdata[[k]]],{k,1,nData,dk}];
ListPlot[errors,PlotLabel->dk Mean[errors]/nData]

Three layer DNN

nLayersF=3;  nLayersM=3;
nlmF = Table[{nIn,LogisticSigmoid},{n,nLayersF}];
nlmM =Table[{nOut,LogisticSigmoid},{n,nLayersM}];
nlmTL=NetChain[Flatten[{nlmF,nlmM}]]

nlmTLTrained = NetTrain[nlmTL,FMdata]

dk=500;
errors=Table[Norm[nlmTLTrained[Fdata[[k]]] - Mdata[[k]]],{k,1,nData,dk}];
ListPlot[errors,PlotLabel->dk Mean[errors]/nData]

Four layer DNN

nLayersF=4;  nLayersM=4;
nlmF = Table[{nIn,LogisticSigmoid},{n,nLayersF}];
nlmM =Table[{nOut,LogisticSigmoid},{n,nLayersM}];
nlmQL=NetChain[Flatten[{nlmF,nlmM}]]

nlmQLTrained = NetTrain[nlmQL,FMdata]

dk=500;
errors=Table[Norm[nlmQLTrained[Fdata[[k]]] - Mdata[[k]]],{k,1,nData,dk}];
ListPlot[errors,PlotLabel->dk Mean[errors]/nData]