MATLAB MESH GENERATION

Mesh Generation for a Curved Structure

clc
clear
M=8;
N=8;
L=8;
LX=1;
LY=1;
LZ=1;
dx=LX/M;
dy=LY/N;
dz=LZ/L;
for i=1:M;
for j=1:N;
for k=1:L;
x(i,j,k)=i*dx;
y0(i,j,k)=j*dy;
y(i,j,k)=(x(i,j,k))^2+y0(i,j,k);
z(i,j,k)=k*dz;
end
end
end
for i=1:M;
for j=1:N;
for k=1:L;
plot3(x(i,j,k),y(i,j,k),z(i,j,k),'*r')
hold on
grid on
xlabel('x')
ylabel('y')
zlabel('z')
end
end
end
for k=1:L;
for i=1:M;
for j=1:N;
if (i<M)
xx=[ x(i+1,j,k) x(i,j,k) ];
yy=[ y(i+1,j,k) y(i,j,k) ];
zz=[ z(i,j,k) z(i,j,k) ];
f=line(xx,yy,zz);
elseif (i==M)
xx=[ x(i,j,k) x(i,j,k) ];
yy=[ y(i,j,k) y(i,j,k) ];
zz=[ z(i,j,k) z(i,j,k)];
f=line(xx,yy,zz);
end
end
end
end
for k=1:L;
for i=1:M;
for j=1:N;
if (j<N)
xxx=[ x(i,j+1,k) x(i,j,k) ];
yyy=[ y(i,j+1,k) y(i,j,k) ];
zzz=[ z(i,j,k) z(i,j,k) ];
f=line(xxx,yyy,zzz);
elseif (j==N)
xxx=[ x(i,j,k) x(i,j,k) ];
yyy=[ y(i,j,k) y(i,j,k) ];
zzz=[ z(i,j,k) z(i,j,k) ];
f=line(xxx,yyy,zzz);
end
end
end
end
for k=1:L;
for i=1:M;
for j=1:N;
if (k<L)
xxxx=[ x(i,j,k) x(i,j,k) ];
yyyy=[ y(i,j,k) y(i,j,k) ];
zzzz=[ z(i,j,k+1) z(i,j,k) ];
f=line(xxxx,yyyy,zzzz);
elseif (k==L)
xxxx=[ x(i,j,k) x(i,j,k) ];
yyyy=[ y(i,j,k) y(i,j,k) ];
zzzz=[ z(i,j,k) z(i,j,k) ];
f=line(xxxx,yyyy,zzzz);
end
end
end
end

Mesh Generation in Cartesian Coordinates

clc
clear
M=7;
N=7;
L=7;
LX=1;
LY=1;
LZ=1;
dx=LX/M;
dy=LY/N;
dz=LZ/L;
for k=1:L;
if (k==1 |k==L)
c=0;
else
c=1;
end
for i=1:M;
for j=1:N;
if (j==1 | i==1)
u(i,j,k)=0*c;
elseif (j==N | i==M)
u(i,j,k)=0;
else u(i,j,k)=c*randn(1,1);
end
end
end
end
for i=1:M;
for j=1:N;
for k=1:L;
x(i,j,k)=i*dx;
y(i,j,k)=j*dy;
z(i,j,k)=k*dz;
end
end
end
figure
for i=1:M;
for j=1:N;
for k=1:L;
plot3(x(i,j,k),y(i,j,k),z(i,j,k),'*r')
end
hold on
end
end
title('Domain Points')
xlabel('x')
ylabel('y')
zlabel('z')
tic
for i=1:M;
for j=1:N;
for k=1:L;
xc=x(i,j,k);
yc=y(i,j,k);
zc=z(i,j,k);
hold on
plot3(xc,yc,zc,'*b')
axis equal
%--------------------------------------------------------
%Bottom Plane
xc1=xc+dx/2;
yc1=yc+dy/2;
zc1=zc-dz/2; xc2=xc-dx/2; yc2=yc+dy/2; zc2=zc-dz/2; xc3=xc-dx/2; yc3=yc-dy/2; zc3=zc-dz/2;
xc4=xc+dx/2; yc4=yc-dy/2; zc4=zc-dz/2; v1=[xc1,xc2]; v2=[yc1,yc2]; v3=[zc1,zc2]; f=line(v1,v2,v3); v1=[xc2,xc3];
v2=[yc2,yc3]; v3=[zc2,zc3]; f=line(v1,v2,v3); v1=[xc3,xc4]; v2=[yc3,yc4]; v3=[zc3,zc4]; f=line(v1,v2,v3); v1=[xc4,xc1];
v2=[yc4,yc1]; v3=[zc4,zc1]; f=line(v1,v2,v3); %--------------------------------------------------------
%Top Plane xc5=xc+dx/2;
yc5=yc+dy/2;
zc5=zc+dz/2; xc6=xc-dx/2; yc6=yc+dy/2; zc6=zc+dz/2; xc7=xc-dx/2; yc7=yc-dy/2;
zc7=zc+dz/2;
xc8=xc+dx/2;
yc8=yc-dy/2;
zc8=zc+dz/2;
v1=[xc5,xc6];
v2=[yc5,yc6];
v3=[zc5,zc6];
f=line(v1,v2,v3);
v1=[xc6,xc7];
v2=[yc6,yc7];
v3=[zc6,zc7];
f=line(v1,v2,v3);
v1=[xc7,xc8];
v2=[yc7,yc8];
v3=[zc7,zc8];
f=line(v1,v2,v3);
v1=[xc8,xc5];
v2=[yc8,yc5];
v3=[zc8,zc5];
f=line(v1,v2,v3);
%--------------------------------------------------------
% Side Planes
v1=[xc1,xc5];
v2=[yc1,yc5];
v3=[zc1,zc5];
f=line(v1,v2,v3);
v1=[xc2,xc6];
v2=[yc2,yc6];
v3=[zc2,zc6];
f=line(v1,v2,v3);
v1=[xc3,xc7];
v2=[yc3,yc7];
v3=[zc3,zc7];
f=line(v1,v2,v3);
v1=[xc4,xc8];
v2=[yc4,yc8];
v3=[zc4,zc8];
f=line(v1,v2,v3);
%--------------------------------------------------------
end
end
end
toc

Example Mesh Generation

clc
clear
x0=-1;
x00=1;
y0=1;
y00=-1;
M=30;
N=30;
dx=(x00-x0)/M;
dy=(y00-y0)/N;
[x,y] = meshgrid(x0:dx:x00,y0:dy:y00);
for i=1:M+1;
for j=1:N+1;
z(i,j)=0;
end
end
tri = delaunay(x,y);
trisurf(tri,x,y,z)

Example Generating a Mesh

clc
clear
x=[1:1:10];
y=[1:1:10];
M=10;
tt=12;
d=0;
for tt=1:M;
    d=d+0.01;
    %for i=1:M
        yy=d*y;
        h=plot(x,yy,'.');
        grid on
        hold on
   % end
    set(h,'EraseMode');
end

Example Generation a Surface

clc
clear
x0=0;
x00=6;
M=16;
dz=x00-x0;
dx=dz/M;
k=2;
m=1;
y0=0
for j=1:k;
    z=j*dz;
for i=1:M;
    x(i)=i*dx
    y(i)=m*(x(i)-z)+y0;
    grid on
    plot(x,y)
grid on
xlabel('x')
ylabel('y')
hold on
end
end

Plotting Matrix Nodes, Nodes Coordinates, Mesh Function

clc
clear
LX=1
LY=1
N=12
M=12
DX=LX/N
DY=LY/M
for I=1:N
for J=1:M
x(I,J)=I*DX;
y(I,J)=J*DY;
A(I,J)=sin(x(I,J)*y(I,J))/(x(I,J)*y(I,J))+0.01*randn(1,1);
end
end

for I=1:N
for J=1:M
plot (I,J,'.')
hold on
grid on
end
end
title('Relationship Between M and N')
xlabel('N Number of Nodes')
ylabel('M Number of Nodes')

figure(2)
for I=1:N
for J=1:M
plot (x(I,J),y(I,J),'.')
hold on
grid on
end
end
title('Relationship Between M and N')
xlabel('Node Coordinate in x direction')
ylabel('Node Coordinate in y direction')

figure(3)
mesh(x,y,A)

Staggered Grid

clc
clear
%Absolute Referance Mesh Input Data x0=0; x00=1; y0=0; y00=1; lx=x00-x0; ly=y00-y0; M=10; N=10; dx=lx/M; dy=ly/N; %offseting Distance ddx=+0.5; ddy=+0.5; %Relative Referance Mesh Input Data x10=0; x11=1; y10=0; y11=1; lx1=(x11-x10)+ddx;
ly1=(y11-y10)+ddy;
M1=10;
N1=10;
dx1=lx1/M1;
dy1=ly1/N1;
% Matrix Mesh Generation matlab Command
[x,y]=meshgrid(x0:dx:x00,y0:dy:y00);
[x1,y1]=meshgrid(x10:dx1:x11,y10:dy1:y11);
% Output Data Visulization
% figure(1)
% plot(x,y,'*b')
% %figure(2)
% hold on
% plot(x1,y1,'*r')
% hold off
% grid on
figure(1)
subplot(2,2,1)
plot(x,y,'*b')
xlabel('x')
ylabel('y')
grid on
title('First Grid ')
subplot(2,2,2)
plot(x1,y1,'*r')
xlabel('x')
ylabel('y')
grid on
title('Second Grid ')
subplot(2,2,[3 4])
plot(x,y,'*b',x1,y1,'*r')
xlabel('x')
ylabel('y')
grid(gca,'minor')
title('Staggered Grid ')

Generation of 3D Grids Points

Unless otherwise noted, all content on this site is @Copyright by Ahmed Al Makky 2012-2013 - http://cfd2012.com