跳转至

EarthMoon

user_EarthMoon1.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
%make a box model, which will be put on a slope
clear;
fs.randSeed(1);%build random model
B=obj_Box;%build a box object
B.name='EarthMoon';
B.GPUstatus='auto';
B.ballR=7e5;
B.isClump=0;
B.distriRate=0.2;
B.sampleW=16e6;
B.sampleL=16e6;
B.sampleH=18e6;
B.setType('topPlaten');

B.buildInitialModel();%B.show();
d=B.d;
d.mo.aR(1:d.mNum)=d.mo.aR(1:d.mNum)*0.98;%reduce the element size a bit
%--------------end initial model------------
B.gravitySediment();
B.compactSample(2);
mfs.reduceGravity(d,5);%reduce the gravity of element
%------------return and save result--------------
d.status.dispEnergy();%display the energy of the model
%d.showFilter('Group',{'sample'});
d.show('mGZ');

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();%because data is clear, it will be re-calculated
user_EarthMoon2.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
%set the material of the model
clear
load('TempModel/EarthMoon1.mat');
spaceSize=B.sampleW*6;%side length of square area
%----------------make objects----------------
fastGroupModel=0;%@@default value is 0, fast when the value is 1
earthR=(6.371e6-B.ballR)*1.125;
moonR=3.476e6-B.ballR;
packBallR=B.ballR;%record the ballR of the small box block
packBoxObj=B.d.group2Obj('sample');%make struct object from a group
packBoxObj=mfs.moveObj2Origin(packBoxObj);
earthObj=mfs.cutSphereObj(packBoxObj,earthR);
sphereObj2=mfs.cutSphereObj(packBoxObj,moonR);
boxObj=mfs.cutBoxObj(packBoxObj,B.sampleW*0.3,B.sampleW*0.1,B.sampleW*0.1);

rate=1;
sphereObj2=mfs.move(sphereObj2,B.sampleW*rate,0,0);
boxObj=mfs.move(boxObj,0,-B.sampleW*rate*2,0);
%----------------end make objects----------------

fs.randSeed(1);%build random model
B=obj_Box;%build a box object
B.name='EarthMoon';
B.GPUstatus='auto';
B.ballR=packBallR;
B.isClump=0;
B.isSample=0;
B.distriRate=0.0;
B.sampleW=B.ballR*10;
B.sampleL=B.ballR*10;
B.sampleH=B.ballR*10;
B.boundaryStatus=[0,0,0,0,1,0];
B.setType('botPlaten');

B.SET.fastGroupModel=fastGroupModel;%@@default value is 0, fast when the value is 1
B.SET.spaceSize=spaceSize;
B.buildInitialModel();%B.show();
d=B.d;



%please change the material in other simulations
matTxt=[1.7e12  0.15    20e6    200e6   0.6 5500];%load material file
Mats{1,1}=material('Soil1',matTxt,B.ballR);
Mats{1,1}.Id=1;
d.Mats=Mats;%you may remove this line in other simulations

[earthId,sphere2Id,boxId]=d.addElement(1,{earthObj,sphereObj2,boxObj});
d.addGroup('earth',earthId);
d.addGroup('sphere2',sphere2Id);
d.addGroup('box',boxId);
d.delElement('botPlaten');
d.GROUP.groupProtect=[];
d.delElement('botB');
d.mo.aX(end)=0;d.mo.aY(end)=0;d.mo.aY(end)=0;
%d.showB=2;d.show('aR');return
frame.minX=-spaceSize/2;
frame.minY=-spaceSize/2;
frame.minZ=-spaceSize/2;
frame.maxX=spaceSize/2;
frame.maxY=spaceSize/2;
frame.maxZ=spaceSize/2;
d.mo.frame=frame;
d.mo.mAZ(:)=0;%no earth gravitation

%----------test error between two methods, can be removed
tic
planetfs.resetmGXYZ(d);
planetfs.setModelGravitation(d);
time1=toc;
totalGX1=sum(d.mo.mGX(d.GROUP.earth));
%d.show('aR');return;

tic
planetfs.resetmGXYZ(d);
planetfs.setGroupOuterGravitation(d);
time2=toc;
totalGX2=sum(d.mo.mGX(d.GROUP.earth));
disp(['Model:' num2str(totalGX1) '; Time:' num2str(time1)]);
disp(['Group:' num2str(totalGX2) '; Time:' num2str(time2)]);
%----------end test error between two methods, can be removed
d.show('aR');return;


%------------balance the force on each planet
d.mo.setShear('off');
d.mo.dT=d.mo.dT*4;
planetfs.resetmGXYZ(d);
planetfs.setGroupInnerGravitation(d,{'earth','sphere2','box'});
d.balance('Standard',0.2);
if B.SET.fastGroupModel==1
    d.setClump('earth');
    d.setClump('sphere2');
    d.setClump('box');
    d.mo.aBF=1e100+d.mo.aBF*1e6;
    planetfs.resetmGXYZ(d);
    d.balance('Standard',0.5);
end
d.mo.dT=d.mo.dT/4;

[escapeV,earthM,earthR]=planetfs.getEscapeSpeed(d,'earth');

d.mo.setGPU('off');
d.clearData(1);%clear dependent data
d.recordCalHour('Step1Finish');
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();%because data is clear, it will be re-calculated
user_EarthMoon3.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
clear;
load('TempModel/EarthMoon2.mat');
d.calculateData();
d.mo.setGPU('off');
B.setUIoutput();
d=B.d;
d.getModel();%reset the initial status of the model
d.resetStatus();%initialize model status, which records running information
d.mo.isHeat=1;%calculate heat in the model

visRate=0.0000;
d.mo.mVis=d.mo.mVis*visRate;
d.mo.mVY(d.GROUP.earth)=0;
d.mo.mVY(d.GROUP.sphere2)=5e3;
d.mo.mVX(d.GROUP.box)=-3e3;
planetfs.setZeroMomentum(d);
%----------set the drawing of result during iterations
showType='mV';
d.showB=4;%only frame is shown

groupCenter=planetfs.getGroupCenter(d);%record the center of each group
earthCenter=groupCenter.earth;
sphere2Center=groupCenter.sphere2;
boxCenter=groupCenter.box;
%----------end set the drawing of result during iterations
gpuStatus=d.mo.setGPU('auto');
totalCircle=100;

calGInterval=2;%gravitation calculation interval
recordInterval=2;%recordStatus interval
d.mo.dT=d.mo.dT*4;
d.tic(totalCircle);
fName=['data/step/' B.name num2str(B.ballR) '-' num2str(B.distriRate) 'loopNum'];
save([fName '0.mat']);%return;
for i=1:totalCircle
    for j=1:100
        %---------computing
        d.balance();
        if mod(j,calGInterval)==0
            planetfs.resetmGXYZ(d);
            if B.SET.fastGroupModel==1
                planetfs.setGroupOuterGravitation(d,{'earth','sphere2','box'});
            else
                planetfs.setModelGravitation(d);
            end
        end

        %---------limit space
        sphereR=B.SET.spaceSize*sqrt(2)/2;
        planetfs.limitElementInSphere(d,sphereR);

        %---------record data
        if mod(j,recordInterval)==0
            groupCenter=planetfs.getGroupCenter(d);
            earthCenter=[earthCenter;groupCenter.earth];
            sphere2Center=[sphere2Center;groupCenter.sphere2];
            boxCenter=[boxCenter;groupCenter.box];
            d.recordStatus();
        end
    end
    %---------show the result
    d.figureNumber=d.show(showType);
    view(0,90);
    plot3(earthCenter(:,1),earthCenter(:,2),earthCenter(:,3),'-b');
    plot3(sphere2Center(:,1),sphere2Center(:,2),sphere2Center(:,3),'-r');
    plot3(boxCenter(:,1),boxCenter(:,2),boxCenter(:,3),'-k');

    save([fName num2str(i) '.mat']);
    pause(0.05);
end

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();