• /
• # Fracture

on 20 Nov 2023
• 37
• 110
• 1
• 6
• 1303
drawframe(1);
function drawframe(f)
% Done with my kids! Explosion using voronoi cells.
persistent V T T2 Vg Xo R C L L2
rng default
N=150; % Number of voronoi domains
v=@(x)vecnorm(x); % This is used a lot...
if f==1
% Distribute points in the unit sphere, biased toward the center
Xo=randn(3,N);
Xo=Xo./v(Xo).*rand(1, N);
% Bounding layer of points that will create our outer surface. Need
% lots of them...
NA=100;
ps=randn(3,NA);
ps=1.3*ps./v(ps);
% Concatenate
X=[Xo,ps]';
% Voronoi diagram
[V,R]=voronoin(X);
mnR=cellfun(@min, R)~=1; % Which cells have inf's
ginds=unique(cell2mat(R(~mnR)')); % Get bordering nodes
Iinds=setdiff(1:size(V,1), ginds); % Get interior nodes
mxr=max(v(V(Iinds,:)'));
ginds(1)=[]; % Get rid of inf
% Make non-inf outer-nodes have unit radius * some small scale factor
V(ginds,:)=1.3*V(ginds,:)./v(V(ginds,:)')';
% Glow
Vg=ones(size(V,1),1);
Vg(ginds)=0;
TS=@(k,x,y,z,C)trisurf(k,x,y,z,'FaceC',C,'EdgeC','none');
cnt = 1;
for n = 1:length(mnR)
if mnR(n) == 1
xt=V(R{n},1);
yt=V(R{n},2);
zt=V(R{n},3);
C=[1,1,1];
k = convhull(xt,yt,zt);
T{n}=TS(k,xt,yt,zt,C/2);
hold on;
material(T{n},[0,1,0,3]);
s=1.1;
T2{n}=TS(k,xt*s,yt*s,zt*s,C);
material(T2{n},[1,0,0,3]);
if cnt == 1
set(gca, 'color', 'k');
axis equal off
axis([-1,1,-1,1,-1,1]*6);
cnt = cnt + 1;
camproj p
camva(70);
campos([-55-5 -71 52]/30);
set(gcf,'color','k');
L2=light;
L{1}=light('position',[0,0,0],'style','local');
L{2}=light('position',[0.1,0,0], 'style','local');
end
end
end
elseif f<10
for n=1:N
end
L2.Color = C/f;
elseif f >= 10
% Loop over fragments and expand
for n = 1:N
T2{n}.Vertices=1.3*V(R{n},:)*4000/f.^3;
T{n}.Vertices=T{n}.Vertices+2.5*Xo(:,n)'/f;
end
if f > 20
for n=1:2
L{n}.Color=L{n}.Color*.96;
end
L2.Color = C*f/48;
end
end
end
Tim
Tim