MATLAB: How do i find geometry properties of triangles

centroidconnectivitydelaunaygeometrymoment of inertiatriangulationvertices

I want to write a code to find the polar moment of inertia of an airfoil. Being a complex shape, a friend suggested a numerical estimate. I wrote a code to create the airfoil geometry and split it up into triangles using delaunay. I now need to find the moment of inertia about x and y for each triangle and use the parallel axis theorem to relate those to the centroid.

I have the delaunay triangulation complete but how can i find the Ix and Iy of each triangle? I don't even know where to start.

here's the code i have:

 chord = 0.25; %chord of airfoil
nacanum = 12; % NACA XX## number for the airfoil
t = nacanum/100; %max thickness of airfoil
x=0:0.005:chord; %points along x axis from 0 to chord length
distratio = x/chord; %ratio of points along x to the chord
% function generating top half of airfoil
yt = 5*t*(0.2969*sqrt(distratio)-0.126*distratio-0.3516*distratio.^2+0.2843*distratio.^3-0.1015*distratio.^4);
ytneg=yt*-1; %generatingg bottom half of airfoil
yfull= [yt ytneg]; % full points for y axis of airfoil
y=yfull';
x2=flipud(x);
xfull=[x x2];%full points for x axis of airfoil
x=xfull'; 
TR=delaunayTriangulation(x,y); %triangulation of the airfoil
triplot(TR,x,y) %plotting the triangulation

Best Answer

You can access x,y coordinates of each triangle like below:

nodes = TR.ConnectivityList ;   % nodal connectivity matrix 
coor = TR.Points ;              % coordinates 
ntri = size(nodes,1) ;          % number of triangles 
% loop for each triangle 
for i = 1:ntri
    xi = coor(nodes(i,:),1) ;   % x coordinates of i'th rriangle 
    yi = coor(nodes(i,:),2) ;   % y coordinates of i'th triangle 
    % do what you want
end

But let me tell you, your mesh of air foil is very bad. I suggest you to go through the following meshing package: https://in.mathworks.com/matlabcentral/fileexchange/25555-mesh2d-automatic-mesh-generation

Related Solutions

MATLAB: Prevent nodes from seeing each other across a circle.

Your problem is you are trying to use a delaunay triangulation tool to create the result. A general delaunay triangulation describes a CONVEX domain. A region with a hole inside is not convex by definition, so it will fail.

(It WOULD have helped if you provided the data, so I could easily give you some examples. Instead, I will need to cook some up. And, of course, I won't spend the time to make a very good example.) I'll even try constraining the edges of the circle to be edges in the tringulation.

nc = 31;
t = linspace(0,2*pi,nc)';
ps = polyshape(cos(t)+1,sin(t)+1); % used to initially exclude generated points in the circle
t(end) = [];
xyc = [cos(t),sin(t)] + 1; % A unit radius circle centered at [1,1]
tric = (1:nc-1)';
tric = [tric,tric+1];
tric(end,2) = 1;
xyr = rand(500,2)*6 - 3;
k = isinterior(ps,xyr(:,1),xyr(:,2));
xyr(k,:) = [];
xy = [xyc;xyr];
tri = delaunayTriangulation(xy,tric)
triplot(tri)
axis equal

So, did it work? Of course not. A Delaunay trinagulation describes a convex set. So even though I constrained the edges of the circle to be edges of the triangulation, the interior of the circle is still triangulated too.

Could I have done it differently? Well, yes, that is not difficult. For example, start with the delaunay triangulation, as above, wherein I constrained it to have edges on the circle. Then remove any triangles from the triangulation that have the property that ALL 3 vertices of the triangle lie on the circular boundary. This will be sufficient to create a triangulation with the desired property - an interior "circular" hole.

Thus, starting with the constrained Delaunay triangulation, do this:

triCL = tri.ConnectivityList;
triCL(all(triCL <= 30,2),:) = [];
triplot(triCL,tri.Points(:,1),tri.Points(:,2))
axis equal

While that was not so difficult, we can use a simpler scheme yet - an alpha shape, even though it lacks an explicit constraint option.

stri = alphaShape(xy(:,1),xy(:,2),.75,'HoleThreshhold',1);
plot(stri)
axis equal

This did work of couse, since an alpha shape is able to describe a non-convex domain. Note the difference of the two triangulations around the perimeter of the domain. A delaunay triangulation will sometimes tend to have many very long thin triangles around the perimeter. You can see them in the first two figures I show. This does depends of the way the points were generated, of course. However, an alpha shape will not have that failing.

MATLAB: How to make delaunayTriangulation object from points and connectivity in 3D

You should use the function triangulation for this. Read about this.

https://in.mathworks.com/help/matlab/ref/triangulation.html

Related Question