HFXFEM-Multiple-Cracks/elementstressplot.m at master · neclipse/HFXFEM-Multiple-Cracks · GitHub

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
%% stress, pressure and disp in one encriched elemeent during debugging
allgauss=obj.INTELEM(1,9).EnrichGauss;
% allgauss=obj.EnrichItems(1).INTELEM(1,2).EnrichGauss;
% allgauss=obj.Mytips.NextElem.GaussPntDictM;
% elems=mesh.findelems([6.35,30],elemdict);
% elem=elemdict(elems(1));
% elem.calstress;
% allgauss=elem.GaussPntDictM;
% allgauss=obj.Elemahead.EnrichGauss;
x=[allgauss.X];
y=[allgauss.Y];
p=[allgauss.P];
u=[allgauss.U];
% phi=zeros(size(allgauss));
% for i=1:length(allgauss)
%     phi(i)=allgauss(i).Enf{2}.Phi;
% end
ux=u(1,:);
uy=u(2,:);
spn=[allgauss.Stressp];
figure
plot3(x,y,spn(1,:),'o')
figure
plot3(x,y,spn(2,:),'o')
% figure
% plot3(x,y,phi,'o')
figure
plot3(x,y,p,'*')
figure
plot3(x,y,ux,'s')
figure
plot3(x,y,uy,'s')
%% Common postprocessing
% obj.LeakoffVolume/obj.CrackVolume
% encrack.LeakoffVolume/encrack.CrackVolume
% the second parameter is mesh structure, "3" means trimesh. In fact, I can
% put any number, as long as its size is not 2, which is for quadmesh.
% Step1.Postprocess(6).plotsurf('SPN',3,0,2);
%Step1.Postprocess(end).plotsurf('SPN',3,0,2);
% Step1.Postprocess(end).plotsurf('P',3,0);
% obj.Postprocess(4).plotsurf('P',3,0);
% obj.Postprocess(4).plotsurf('SPN',3,0,2);
% obj.showme(1,'Pfrack', 'second', 'Aperture');
% obj.showme(1,'Pfrack','second','CTraction','component',2);
% obj.showme(1,'Aperture','second','CTraction','component',2);
postprocess.EnrichItems{1}.showme(1,'Aperture','second','CTraction','component',2);
obj.EnrichItems(1).showme(1,'Aperture')
obj.EnrichItems(2).showme(3,'Aperture')

%% Plot crack shape history at give time
Step1.snapshot('crackgeo',47,'timeincs',30,'basecoordinate',30);
% Step1.snapshot('crackgeo',[1,10,20,40]);
timelist=[0,5,10,15,20,25,30]; % time steps to report
timeincs=[Step1.Postprocess.Inc];
incmat=repmat(timeincs,length(timelist),1);
timediff=incmat-transpose(timelist);
[~,inclist]=min(abs(timediff),[],2);
Step1.snapshot('crackgeo',inclist,'timeincs',timelist);

%% element plot for arbitrary element

elems=obj.Preprocess.Mesher.findelems([0.08,30.04],obj.ElemDict);
elem=obj.ElemDict(elems(1));
elem.calstress(false);
allgauss=elem.GaussPntDictM;
% allgauss=obj.Elemahead.GaussPntDictM;
x=[allgauss.X];
y=[allgauss.Y];
p=[allgauss.P];
u=[allgauss.U];
ux=u(1,:);
uy=u(2,:);
% spn=[allgauss.Stressp];
% figure
% plot3(x,y,spn(2,:),'o')
figure
plot3(x,y,p,'*')
% figure
% plot3(x,y,uy,'s')
% figure
% plot3(x,y,ux,'s')
%%
elem=obj.Elemdict(9253);
elem.calstress;
allgauss=elem.GaussPntDictM;
% allgauss=obj.Elemahead.GaussPntDictM;
x=[allgauss.X];
y=[allgauss.Y];
p=[allgauss.P];
uy=[allgauss.Uy];
spn=[allgauss.Stressp];
figure
plot3(x,y,spn(2,:),'o')
figure
plot3(x,y,p,'*')
%% path plot