跳转至

BoxCompaction

user_BoxCompaction1.m
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
clear;
fs.randSeed(2);
sampleFileName='';
%sampleFileName='TempModel/BoxTakeGroup2.mat';
SET.sampleW=0.152;SET.sampleH=0.120;
%----------set of the sand sample
if isempty(sampleFileName)
    hRate=1.02;% box is bigger a bit, 1.0
    grainDensity=2650;
    sampleW=0.152;sampleL=0;totalM=250e-3;
    %minmum diameter, maximum diameter,
    grainSizeDistribution=[0.001,0.002,0.2;0.002,0.004,0.4;0.004,0.008,0.2;0.008,0.016,0.1]*2;
    grainR=mfs.getGradationDiameter(grainSizeDistribution,totalM/grainDensity)/2;
    SET=mfs.getBoxSample(grainR,sampleW,sampleL,hRate);
    SET.isSample=1;
else
    load(sampleFileName);
    SET.ballR=sqrt(max(grainObj.R)*min(grainObj.R));
    SET.isSample=0;
end

SET.grainDensity=grainDensity;
SET.sampleFileName=sampleFileName;
SET.boxDiameter=0.152;SET.boxHeight=0.120;
SET.addHeight1=0.05;SET.addHeight2=0.01;
%----------end set of the sand sample
%----------set of the equipments
SET.hammerMass=4.5;
SET.hammerDistance=0.457;
SET.hammerV=sqrt(SET.hammerDistance*2/9.8);
SET.hammerDiameter=0.051;
SET.hammerArea=pi*SET.hammerDiameter^2/4;
SET.hammer2DAreaRate=SET.ballR*2/SET.hammerDiameter;
SET.hammerMass2D=SET.hammerMass*SET.hammer2DAreaRate;
SET.hammerGZ2D=SET.hammerMass2D*-9.8;
SET.hammerX=0.0125;
%----------end set of the equipments

%----------start building the model
fs.randSeed(1);%random model seed, 1,2,3...
B=obj_Box;%declare a box object
B.name='BoxCompaction';
%--------------initial model------------
B.GPUstatus='auto';%program will test the CPU and GPU speed, and choose the quicker one
B.ballR=SET.ballR;
B.distriRate=0.2;%define distribution of ball radius, 
B.sampleW=SET.boxDiameter;%inner diameter
B.sampleL=0;%when L is zero, it is a 2-dimensional model
B.sampleH=SET.sampleH;%
B.isSample=SET.isSample;
B.SET=SET;
B.platenStatus=[0;0;0;0;0;1];%topPlaten
B.buildInitialModel();%B.show();
B.setUIoutput();
d=B.d;%d.breakGroup('sample');d.breakGroup('lefPlaten');

if isempty(sampleFileName)
    sampleNum=length(d.GROUP.sample);
    d.mo.setGPU('off');
    d.delElement(SET.grainNum+1:sampleNum);
    sId=d.GROUP.sample;
    d.aR(sId)=grainR;
    d.mo.aR(sId)=grainR;
    d.groupMat2Model({'sample'},1);%apply the material
    d.mo.setNearbyBall();
    d.breakGroup();
    B.convert2D(B.ballR);%change ball properties to 2D
else
    d.GROUP.sample=d.addElement(1,grainObj);
    d.setClump();
    d.SET.ballArea=grainObj.ballArea;
end
%d.mo.setShear('off');

d.mo.aMUp(:)=0.5;%increase friction coefficient to increase porosity
d.mo.setGPU('auto');

%---------- gravity sedimentation
B.gravitySediment();%you may use B.gravitySediment(10); to increase sedimentation time (10)
d.show('mV');
%return
%mfs.reduceGravity(d,2);
%------------return and save result--------------
porosity=mfs.get2DPorosity(d,0,B.sampleW,0,B.sampleH/2);
d.status.dispEnergy();%display the energy of the model
d.show('aR');

d.mo.setGPU('off');
d.clearData(1);%clear dependent data
d.recordCalHour('Step1Finish');
save(['TempModel/' B.name '1.mat'],'B','d');
save(['TempModel/' B.name '1R' num2str(B.ballR) '-distri' num2str(B.distriRate)  'aNum' num2str(d.aNum) '.mat']);
d.calculateData();
user_BoxCompaction2.m
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
%set the material of the model
clear
load('TempModel/BoxCompaction1.mat');
B.setUIoutput();%set the output
d=B.d;
d.calculateData();%calculate data
d.mo.setGPU('off');%close the GPU calculation
d.getModel();%get xyz from d.mo
%---------------delele elements on the top
SET=B.SET;

%--------------assign new material
matTxt=load('Mats\Soil1.txt');%load material file
Mats{1,1}=material('Soil1',matTxt,B.ballR);
Mats{1,1}.Id=1;
d.Mats=Mats;%assign new material
d.groupMat2Model({'sample'},1);%apply the new material
if isempty(SET.sampleFileName)
d.mo.mM=B.SET.grainDensity*(4/3*pi*d.mo.aR(1:d.mNum).^3);
B.convert2D(B.ballR);%change ball properties to 2D
end

%------------set of compaction @@@@@@@@@@@@@@
platenId=[d.GROUP.lefPlaten;d.GROUP.rigPlaten;d.GROUP.botPlaten;d.GROUP.topPlaten];
d.mo.aMUp(platenId)=0;

topX=d.mo.aX(d.GROUP.topPlaten);
topFilter=topX>SET.hammerX&topX<SET.hammerX+SET.hammerDiameter;
d.addGroup('hammer1',d.GROUP.topPlaten(topFilter));

topFilter=topX<SET.boxDiameter-SET.hammerX&topX>SET.boxDiameter-SET.hammerX-SET.hammerDiameter;
d.addGroup('hammer2',d.GROUP.topPlaten(topFilter));

hammer1Rate=ones(size(d.GROUP.hammer1))/length(d.GROUP.hammer1);
hammer2Rate=ones(size(d.GROUP.hammer2))/length(d.GROUP.hammer2);
SET.hammer1M0=d.mo.mM(d.GROUP.hammer1);
SET.hammer2M0=d.mo.mM(d.GROUP.hammer2);
SET.hammer1GZ0=d.mo.mGZ(d.GROUP.hammer1);
SET.hammer2GZ0=d.mo.mGZ(d.GROUP.hammer2);

SET.hammer1M=SET.hammerMass2D*hammer1Rate;
SET.hammer2M=SET.hammerMass2D*hammer2Rate;
SET.hammer1GZ=SET.hammerGZ2D*hammer1Rate;
SET.hammer2GZ=SET.hammerGZ2D*hammer2Rate;
B.SET=SET;
%------------end set of compaction @@@@@@@@@@@@@@

d.mo.setGPU('auto');
d.mo.aMUp(d.GROUP.topPlaten)=0;
%d.balanceBondedModel(0.5);
%d.breakGroup();%break all connections and glue the sample
d.balance('Standard',4);
porosity=mfs.get2DPorosity(d,0,B.sampleW,0,B.sampleH/2);
d.show('StressZZ');

d.clearData(1);%clear dependent data
d.recordCalHour('Step2Finish');
save(['TempModel/' B.name '2.mat'],'B','d');
save(['TempModel/' B.name '2R' num2str(B.ballR) '-distri' num2str(B.distriRate)  'aNum' num2str(d.aNum) '.mat']);
d.calculateData();
user_BoxCompaction3.m
  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
clear;
load('TempModel/BoxCompaction2.mat');
d.calculateData();
d.mo.setGPU('off');
B.setUIoutput();
d=B.d;
d.getModel();%d.setModel();%reset the initial status of the model
d.resetStatus();%initialize model status, which records running information

showType='mV';
figureNumber=d.show(showType);
d.figureNumber=figureNumber;

SET=B.SET;
%reduce stiffness of topPlaten
d.mo.aKN(d.GROUP.topPlaten)=d.mo.aKN(d.GROUP.topPlaten)/10;
d.mo.aKS(d.GROUP.topPlaten)=d.mo.aKS(d.GROUP.topPlaten)/10;

topPlatenVis0=d.mo.mVis(d.GROUP.topPlaten)/10;
sampleVis0=d.mo.mVis(d.GROUP.sample)/10;
topPlatenVis1=d.mo.mVis(d.GROUP.topPlaten)/100;
sampleVis1=d.mo.mVis(d.GROUP.sample)/100;

d.setStandarddT();
%d.mo.dT=d.mo.dT*4;

fName=['data/step/' B.name num2str(B.ballR) '-' num2str(B.distriRate) 'loopNum'];
save([fName '0.mat']);%return;
gpuStatus=d.mo.setGPU('auto');
compactionNum=2;
compactionCircle=2;
resetCircle=1;
standardBalanceNum=0.5;
boxHeight=mean(d.mo.aZ(d.GROUP.topPlaten))-6*B.ballR;
porosity=mfs.get2DPorosity(d,0,B.SET.boxDiameter,0,boxHeight);
d.tic(compactionCircle*compactionNum*2+resetCircle*compactionNum*2);

IStart=0;Iend=0;
for n=1:compactionNum
    %compaction of left hammer
    d.mo.mM(d.GROUP.hammer1)=SET.hammer1M;%assign load to left hammer
    d.mo.mGZ(d.GROUP.hammer1)=SET.hammer1GZ;
    d.mo.mVZ(d.GROUP.hammer1)=-SET.hammerV;
    d.mo.mVis(d.GROUP.topPlaten)=topPlatenVis1;
    d.mo.mVis(d.GROUP.sample)=sampleVis1;
    d.show('mVis');

    IStart=Iend+1;
    Iend=IStart+compactionCircle-1;
    for i=IStart:Iend
        d.mo.setGPU(gpuStatus);
        d.balance('Standard',standardBalanceNum);
        d.clearData(1);
        d.mo.setGPU('off');
        save([fName num2str(i) '.mat'],'-v7.3');
        d.calculateData();

        porosity1=mfs.get2DPorosity(d,0,B.SET.boxDiameter,0,boxHeight);
        porosity=[porosity;porosity1];
        d.toc();%show the note of time;
        d.show(showType);
    end
    %d.show('mV');return
    %reset hammer1
    d.mo.mM(d.GROUP.hammer1)=SET.hammer1M0;%reset left hammer
    d.mo.mGZ(d.GROUP.hammer1)=SET.hammer1GZ0;
    d.mo.mVis(d.GROUP.topPlaten)=topPlatenVis0;
    d.mo.mVis(d.GROUP.sample)=sampleVis0;
    IStart=Iend+1;
    Iend=IStart+resetCircle-1;
    for i=IStart:Iend
        d.mo.setGPU(gpuStatus);
        d.balance('Standard',standardBalanceNum);
        d.clearData(1);
        d.mo.setGPU('off');
        save([fName num2str(i) '.mat'],'-v7.3');
        d.calculateData();

        porosity1=mfs.get2DPorosity(d,0,B.SET.boxDiameter,0,boxHeight);
        porosity=[porosity;porosity1];
        d.toc();%show the note of time;
        d.show(showType);
    end

    %compaction of right hammer
    d.mo.mM(d.GROUP.hammer2)=SET.hammer2M;
    d.mo.mGZ(d.GROUP.hammer2)=SET.hammer2GZ;
    d.mo.mVZ(d.GROUP.hammer2)=-SET.hammerV;
    d.mo.mVis(d.GROUP.topPlaten)=topPlatenVis1;
    d.mo.mVis(d.GROUP.sample)=sampleVis1;
    IStart=Iend+1;
    Iend=IStart+compactionCircle-1;
    for i=IStart:Iend
        d.mo.setGPU(gpuStatus);
        d.balance('Standard',standardBalanceNum);
        d.clearData(1);
        d.mo.setGPU('off');
        save([fName num2str(i) '.mat'],'-v7.3');
        d.calculateData();

        porosity1=mfs.get2DPorosity(d,0,B.SET.boxDiameter,0,boxHeight);
        porosity=[porosity;porosity1];
        d.toc();%show the note of time;
        d.show(showType);
    end
    %reset hammer2
    d.mo.mM(d.GROUP.hammer2)=SET.hammer2M0;%reset right hammer
    d.mo.mGZ(d.GROUP.hammer2)=SET.hammer2GZ0;
    d.mo.mVis(d.GROUP.topPlaten)=topPlatenVis0;
    d.mo.mVis(d.GROUP.sample)=sampleVis1;
    IStart=Iend+1;
    Iend=IStart+resetCircle-1;
    for i=IStart:Iend
        d.mo.setGPU(gpuStatus);
        d.balance('Standard',standardBalanceNum);
        d.clearData(1);
        d.mo.setGPU('off');
        save([fName num2str(i) '.mat'],'-v7.3');
        d.calculateData();

        porosity1=mfs.get2DPorosity(d,0,B.SET.boxDiameter,0,boxHeight);
        porosity=[porosity;porosity1];
        d.toc();%show the note of time;
        d.show(showType);
    end
end
figure;
plot(porosity);

d.mo.setGPU('off');
d.clearData(1);
d.recordCalHour('BoxCrush3Finish');
save(['TempModel/' B.name '3.mat'],'d');
save(['TempModel/' B.name '3R' num2str(B.ballR) '-' num2str(B.distriRate)  'aNum' num2str(d.aNum) '.mat']);
d.calculateData();