Showing posts with label Mesh Generation. Show all posts
Showing posts with label Mesh Generation. Show all posts

Friday, December 2, 2011

GID

There is another mesh generator GID which is quite a powerful and complete mesh tool. This software is not free, but you can download it for free and use it with some limitations such as number of elements or nodes. Here I post the download link, give it a try, you won't be disappointed.
http://www.gidhome.com/download
There are versions for Windows, Mac and Linux for 32bit and 64bit.

Friday, November 18, 2011

Transfer your mesh from GMSH to Matlab

I believe examples are on of the best ways to learn, so I start from an example.

1- build your geometry in GMSH (you can easily do that by studying the first tutorial of the software)


2- Mesh the geometry by going to the Mesh tab and clicking on the 2D. This will build a triangular non structured mesh

Sunday, November 13, 2011

GMSH

Here is the link to a very useful free mesh generator and post processor software (GMSH) which is distributed under the terms of the GNU General Public License (GPL). http://geuz.org/gmsh/ There are versions available for Win, Mac and Linux.
You can find tutorials in the tutorial folder of the downloaded file.

Below is an example of unstructured mesh for a circular surface with a squared hole inside.



Friday, November 11, 2011

4 node quadrilateral element for non symmetrical rectangles

Below you may find the Matlab function for 4 node quadrilateral element which is capable of meshing non symmetric rectangles. I will try to do the same for triangles.


function [Edof,Ex,Ey]=quadmesh2(xcorner,ycorner,nx,ny)


% function [Edof,Ex,Ey]=quadmesh2(xcorner,ycorner,nx,ny)
%----------- Description ---------
% Generates quadrilateral 4 node elements
% ------------ INPUT -------------
% xcorner=[lower left, lower right, upper right, upper left] 
% ycorner=[lower left, lower right, upper right, upper left]
% nx= No. elements in x direction
% ny= No. elements in y direction
%       4         3
%       o----------o
%      |          |
%     |          | 
%    |          |
%    o-------- o
%    1         2
% ----------- OUTPUT -------------
% Ex= element x coordinates matrix [nrelem X 4]
% Ey= element y coordinates matrix [nrelem X 4]
% Edof= element DOF [nrelem X 9]
% Written by 
% Shah, 11 November 2011


ld=abs(xcorner(2)-xcorner(1));
lu=abs(xcorner(3)-xcorner(4));
hr=abs(ycorner(3)-ycorner(2));
hl=abs(ycorner(4)-ycorner(1));
for i=1:ny
        for j=1:nx
            dl=(i-1)*(nx+1)+j;
            ul=dl+(nx+1); 
            dr=dl+1;
            ur=ul+1;
            l1=(i-1)*(lu-ld)/ny+ld; l2=i*(lu-ld)/ny+ld; 
            h1=(j-1)*(hr-hl)/nx+hl; h2=j*(hr-hl)/nx+hl;
            x_cd=(i-1)*(xcorner(4)-xcorner(1))/ny+xcorner(1);
            x_cu=i*(xcorner(4)-xcorner(1))/ny+xcorner(1);
            y_cl=(j-1)*(ycorner(2)-ycorner(1))/nx+ycorner(1);
            y_cr=j*(ycorner(2)-ycorner(1))/nx+ycorner(1);
            xy_dl=[(j-1)*(l1/nx)+x_cd (i-1)*(h1/ny)+y_cl];
            xy_dr=[j*(l1/nx)+x_cd     (i-1)*(h2/ny)+y_cr];
            xy_ul=[(j-1)*(l2/nx)+x_cu  i*(h1/ny)+y_cl];
            xy_ur=[j*(l2/nx)+x_cu      i*(h2/ny)+y_cr];
            el=(i-1)*nx+j;
            Edof(el,:)=[el,2*dl-1,2*dl,2*dr-1,2*dr,2*ur-1,2*ur,2*ul-1,2*ul];
            Ex(el,:)=[xy_dl(1) xy_dr(1) xy_ur(1) xy_ul(1)];
            Ey(el,:)=[xy_dl(2) xy_dr(2) xy_ur(2) xy_ul(2)];
        end
end




figure(1);
hold on;
for i=1:length(Edof(:,1));
    plot(Ex(i,[1:end 1]),Ey(i,[1:end 1]));
end
axis equal

3 node symmetric triangular elements

Below you may find a Matlab function to generate 3 node symmetric triangular elements as I promised. 
Note that nx and ny are the number of rectangular elements in each direction.



function [Edof,Ex,Ey]=symtrimesh(xcorner,ycorner,nx,ny)

% function [Edof,Ex,Ey]=symtrimesh(xcorner,ycorner,nx,ny)
%----------- Description ---------
% Generates triangular 3 node symmetric elements
% ------------ INPUT -------------
% xcorner=[lower left, lower right, upper right, upper left] 
% ycorner=[lower left, lower right, upper right, upper left]
% nx= No. of rectangulat elements in x direction
% ny= No. of rectangulat elements in y direction
%    7         8         9        
%    o---------o---------o
%    | +       |       + |
%    |   +     |     +   |
%    |     +   |   +     |
%    |       + | +       |
%   4o---------5---------o6
%    |       + | +       |
%    |     +   |   +     |  
%    |   +     |     +   |
%    | +       |       + |
%    o---------o---------o
%    1         2         3
% ----------- OUTPUT -------------
% Ex= element x coordinates matrix [nrelem X 3]
% Ey= element y coordinates matrix [nrelem X 3]
% Edof= element DOF [nrelem X 7]
% Written by 
% Shah, 11 November 2011

l=xcorner(2)-xcorner(1);
h=ycorner(3)-ycorner(2);
for i=1:ny
        for j=1:nx
            dl=(i-1)*(2*nx+1)+2*j-1;
            dm=dl+1;
            dr=dl+2;
            ul=dl+(2*nx+1); 
            um=ul+1;
            ur=ul+2;
            xy_dl=[(j-1)*(l/nx)+xcorner(1) (i-1)*(h/ny)+ycorner(1)];
            xy_dr=[j*(l/nx)+xcorner(1) (i-1)*(h/ny)+ycorner(1)];
            xy_dm=(xy_dl+xy_dr)./2;
            xy_ul=[(j-1)*(l/nx)+xcorner(1) i*(h/ny)+ycorner(1)];
            xy_ur=[j*(l/nx)+xcorner(1) i*(h/ny)+ycorner(1)];
            xy_um=(xy_ul+xy_ur)./2;
            el=4*nx*(i-1)+4*(j-1)+1;
            if floor(i/2)-i/2~=0
                Edof(el,:)=[el,2*dl-1,2*dl,2*um-1,2*um,2*ul-1,2*ul];
                Edof(el+1,:)=[el+1,2*dl-1,2*dl,2*dm-1,2*dm,2*um-1,2*um];
                Edof(el+2,:)=[el+2,2*dm-1,2*dm,2*dr-1,2*dr,2*um-1,2*um];
                Edof(el+3,:)=[el+3,2*dr-1,2*dr,2*ur-1,2*ur,2*um-1,2*um];
                Ex(el,:)=[xy_dl(1) xy_um(1) xy_ul(1)];
                Ex(el+1,:)=[xy_dl(1) xy_dm(1) xy_um(1)];
                Ex(el+2,:)=[xy_dm(1) xy_dr(1) xy_um(1)];
                Ex(el+3,:)=[xy_dr(1) xy_ur(1) xy_um(1)];
                Ey(el,:)=[xy_dl(2) xy_um(2) xy_ul(2)];
                Ey(el+1,:)=[xy_dl(2) xy_dm(2) xy_um(2)];
                Ey(el+2,:)=[xy_dm(2) xy_dr(2) xy_um(2)];
                Ey(el+3,:)=[xy_dr(2) xy_ur(2) xy_um(2)];
            else
                Edof(el,:)=[el,2*dl-1,2*dl,2*dm-1,2*dm,2*ul-1,2*ul];
                Edof(el+1,:)=[el+1,2*dm-1,2*dm,2*um-1,2*um,2*ul-1,2*ul];
                Edof(el+2,:)=[el+2,2*dm-1,2*dm,2*ur-1,2*ur,2*um-1,2*um];
                Edof(el+3,:)=[el+3,2*dm-1,2*dm,2*dr-1,2*dr,2*ur-1,2*ur];
                Ex(el,:)=[xy_dl(1) xy_dm(1) xy_ul(1)];
                Ex(el+1,:)=[xy_dm(1) xy_um(1) xy_ul(1)];
                Ex(el+2,:)=[xy_dm(1) xy_ur(1) xy_um(1)];
                Ex(el+3,:)=[xy_dm(1) xy_dr(1) xy_ur(1)];
                Ey(el,:)=[xy_dl(2) xy_dm(2) xy_ul(2)];
                Ey(el+1,:)=[xy_dm(2) xy_um(2) xy_ul(2)];
                Ey(el+2,:)=[xy_dm(2) xy_ur(2) xy_um(2)];
                Ey(el+3,:)=[xy_dm(2) xy_dr(2) xy_ur(2)];       

             end            
        end
end


figure(1);
hold on;
for i=1:length(Edof(:,1));
    plot(Ex(i,[1:end 1]),Ey(i,[1:end 1]));
end
axis equal

3 node triangular element

Below you may find a Matlab function to generate 3 node triangular elements. The output will not be symmetric, so I will try to write one for symmetric triangular elements soon.
Note that nx and ny are the number of rectangular elements in each direction.

function [Edof,Ex,Ey]=trimesh(xcorner,ycorner,nx,ny)

% function [Edof,Ex,Ey]=trimesh(xcorner,ycorner,nx,ny)
%----------- Description ---------
% Generates triangular 3 node elements
% ------------ INPUT -------------
% xcorner=[lower left, lower right, upper right, upper left] 
% ycorner=[lower left, lower right, upper right, upper left]
% nx= No. elements in x direction
% ny= No. elements in y direction
%    3         4
%    o---------o
%    | +       |
%    |   +     |
%    |     +   | 
%    |       + |
%    o---------o
%    1         2
% ----------- OUTPUT -------------
% Ex= element x coordinates matrix [nrelem X 3]
% Ey= element y coordinates matrix [nrelem X 3]
% Edof= element DOF [nrelem X 7]
% Written by 
% Shah, 11 November 2011

l=xcorner(2)-xcorner(1);
h=ycorner(3)-ycorner(2);
for i=1:ny
        for j=1:nx
            dl=(i-1)*(nx+1)+j;
            ul=dl+(nx+1); 
            dr=dl+1;
            ur=ul+1;
            xy_dl=[(j-1)*(l/nx)+xcorner(1) (i-1)*(h/ny)+ycorner(1)];
            xy_dr=[j*(l/nx)+xcorner(1) (i-1)*(h/ny)+ycorner(1)];
            xy_ul=[(j-1)*(l/nx)+xcorner(1) i*(h/ny)+ycorner(1)];
            xy_ur=[j*(l/nx)+xcorner(1) i*(h/ny)+ycorner(1)];
            el=2*(i-1)*nx+2*j-1;
            Edof(el,:)=[el,2*dl-1,2*dl,2*dr-1,2*dr,2*ul-1,2*ul];
            Edof(el+1,:)=[el+1,2*dr-1,2*dr,2*ur-1,2*ur,2*ul-1,2*ul];
            Ex(el,:)=[xy_dl(1) xy_dr(1) xy_ul(1)];
            Ex(el+1,:)=[xy_dr(1) xy_ur(1) xy_ul(1)];
            Ey(el,:)=[xy_dl(2) xy_dr(2) xy_ul(2)];
            Ey(el+1,:)=[xy_dr(2) xy_ur(2) xy_ul(2)];
        end
end

figure(1);
hold on;
for i=1:length(Edof(:,1));
    plot(Ex(i,[1:end 1]),Ey(i,[1:end 1]));
end
axis equal

Thursday, November 10, 2011

4 node quadrilateral element

Below you may find a Matlab function written by me to generate 4 node quadrilateral elements. It is very easy and straight forward.


function [Edof,Ex,Ey]=quadmesh(xcorner,ycorner,nx,ny)


% function [Edof,Ex,Ey]=quadmesh(xcorner,ycorner,nx,ny)
%----------- Description ---------
% Generates quadrilateral 4 node elements
% ------------ INPUT -------------
% xcorner=[lower left, lower right, upper right, upper left] 
% ycorner=[lower left, lower right, upper right, upper left]
% nx= No. elements in x direction
% ny= No. elements in y direction
%    4         3
%    o---------o
%    |         |
%    |         | 
%    |         |
%    o---------o
%    1         2
% ----------- OUTPUT -------------
% Ex= element x coordinates matrix [nrelem X 4]
% Ey= element y coordinates matrix [nrelem X 4]
% Edof= element DOF [nrelem X 9]
% Written by 
% Shah, 10 November 2011


l=xcorner(2)-xcorner(1);
h=ycorner(3)-ycorner(2);
for i=1:ny
        for j=1:nx
            dl=(i-1)*(nx+1)+j;
            ul=dl+(nx+1); 
            dr=dl+1;
            ur=ul+1;
            xy_dl=[(j-1)*(l/nx)+xcorner(1) (i-1)*(h/ny)+ycorner(1)];
            xy_dr=[j*(l/nx)+xcorner(1) (i-1)*(h/ny)+ycorner(1)];
            xy_ul=[(j-1)*(l/nx)+xcorner(1) i*(h/ny)+ycorner(1)];
            xy_ur=[j*(l/nx)+xcorner(1) i*(h/ny)+ycorner(1)];
            el=(i-1)*nx+j;
            Edof(el,:)=[el,2*dl-1,2*dl,2*dr-1,2*dr,2*ur-1,2*ur,2*ul-1,2*ul];
            Ex(el,:)=[xy_dl(1) xy_dr(1) xy_ur(1) xy_ul(1)];
            Ey(el,:)=[xy_dl(2) xy_dr(2) xy_ur(2) xy_ul(2)];
        end
end


figure(1);
hold on;
for i=1:length(Edof(:,1));
    plot(Ex(i,[1:end 1]),Ey(i,[1:end 1]));
end
axis equal