how do i rotate an airfoil points? and plot.
Show older comments
How do I rotate an airfoil points? I have all the new chord lengths and the angles in which are needed to rotate. I need to plot all these airfoils 16 times together in respect to their new angle and chord length. I know I need a sort of for loop but evetrytime I plot, the figure is funky looking.
HELP ME PLEASE!!!
NEW Chord lengths and twist angle, R vaules
Chord length= [ 0 0.0082328 0.0164657 0.0246985 0.0329313 0.0411642 0.049397 0.0576298 0.0658626 0.07409548 0.08232831 0.09056114 0.09879397 0.1070268 0.11525963 0.12349246 ]
r =[ 0 0.03 0.06 0.09 0.12 0.15 0.18 0.21 0.24 0.27 0.3 0.33 0.36 0.39 0.42 0.45]
twist = [1.0471976 0.7561128 0.5465782 0.4134997 0.3278728 0.2699279 0.2286826 0.1980428 0.1744757 0.15582879 0.14072889 0.12826373 0.1178059 0.10891077 0.10125497 0.09459804]
BELOW
THIS IS THE CODE I STATRED WITH FOR NACA
THE PLOT SHOULD LOOK SOMETHING LIKE THIS BELOW 
c1= 1;
%c1= [0,0.008232831,0.016465661,0.024698492,0.032931323,0.041164153,0.049396984,0.057629815,0.065862645,0.074095476,0.082328307,0.090561137,0.098793968,0.107026798,0.115259629,0.12349246];
angle =[1.047197551 0.756112778 0.546578176 0.413499657 0.327872784 0.269927858 0.228682627 0.198042808 0.174475668 0.155828787 0.140728889 0.128263734 0.117805904 0.108910765 0.101254972 0.094598036];
s=num2str(5410);
%4 digits
NACA=s;
% m represents first digit out of the scalar
m =str2double(s(1))/100;
%p pulls the second digit out of the scalar
p=str2double(s(2))/10;
%t pulls and represents the third and fourth digit out of the scalar
t=str2double(s(3:4))/100;
% define x
r=linspace(0, c1, 250);
%yt fuction including the coeffiencients
yt =5*t*c1*(.2969*(sqrt(r/c1))+-.1260*(r/c1)+-.3516*(r/c1).^2+.2843*(r/c1).^3+-.1015*(r/c1).^4);
%yt2 =5*t*c2*(.2969*(sqrt(x/c2))+-.1260*(x/c2)+-.3516*(x/c2).^2+.2843*(x/c2).^3+-.1015*(x/c2).^4);
for k = 1:length(r)
if r(k) <= p*c1
yc(k)=m*(r(k)/p^2)*(2*p-(r(k)/c1));
dx(k)=(2*m)/p^2*(p-(r(k)/c1));
elseif r(k) > p*c1
yc(k)=m*((c1-r(k))/(1-p)^2)*(1+(r(k)/c1)-(2*p));
dx(k)=((2*m)/(1-p)^2)*(p-(r(k)/c1));
end
%upper and lower limits of the airfoil (xu,yu) ; (xl,yl)
angle=atan(dx(k));
xu(k)=r(k)-yt(k)*-sin(angle);
yu(k)=yc(k)+yt(k)*cos(angle);
xl(k)=r(k)+yt(k)*sin(angle);
yl(k)=yc(k)-yt(k)*cos(angle);
end
%plot of airfoil
plot(xu,yu,xlabel('x'),ylabel('y'));
title( ['NACA Airfoil ', NACA]);
hold on
plot(xl,yl,'r')
plot(r,yc,'g')
grid on
axis equal
%--------------------------------------------------------------------------------------%
% Given data vectors X and Y.
% Want to rotate the data through angle "ang" about rotation center Xc, Yc
X = [xl,xu];
Y = [yl,yu];
ang = [254]; % 0.756112778 0.546578176 0.413499657 0.327872784 0.269927858 0.228682627 0.198042808 0.174475668 0.155828787 0.140728889 0.128263734 0.117805904 0.108910765 0.101254972 0.094598036];
% Specify the coordinates of the center of rotation
Xc = 0.25 ; % Rotate about the 1/4 chord point
Yc = 0 ;
% The data is roated in a three-step process
% Step 1) Shift the data to the rotation center
Xs = X - Xc; % shifted data
Ys = Y - Yc;
% Step 2) Rotate the data
Xsr = Xs.*cos(ang) + Ys.*sin(ang); % shifted and rotated data
Ysr = -Xs*sin(ang) + Ys*cos(ang); %
% Step 3) Un-shift the data (back to the original coordinate system)
Xr = Xsr + Xc; % Rotated data
Yr = Ysr + Yc;
hold on
Accepted Answer
David Wilson
on 18 Apr 2019
Edited: David Wilson
on 18 Apr 2019
Does this look better?
Some (small) changes to your code. I swapped your rotation direction to better follow your example plot. Probably not a good idea to call a variable angle.
c1= 1;
%c1= [0,0.008232831,0.016465661,0.024698492,0.032931323,0.041164153,0.049396984,0.057629815,0.065862645,0.074095476,0.082328307,0.090561137,0.098793968,0.107026798,0.115259629,0.12349246];
s=num2str(5410);
%4 digits
NACA=s;
% m represents first digit out of the scalar
m =str2double(s(1))/100;
%p pulls the second digit out of the scalar
p=str2double(s(2))/10;
%t pulls and represents the third and fourth digit out of the scalar
t=str2double(s(3:4))/100;
% define x
r=linspace(0, c1, 250);
%yt fuction including the coeffiencients
yt =5*t*c1*(.2969*(sqrt(r/c1))+-.1260*(r/c1)+-.3516*(r/c1).^2+.2843*(r/c1).^3+-.1015*(r/c1).^4);
%yt2 =5*t*c2*(.2969*(sqrt(x/c2))+-.1260*(x/c2)+-.3516*(x/c2).^2+.2843*(x/c2).^3+-.1015*(x/c2).^4);
for k = 1:length(r)
if r(k) <= p*c1
yc(k)=m*(r(k)/p^2)*(2*p-(r(k)/c1));
dx(k)=(2*m)/p^2*(p-(r(k)/c1));
elseif r(k) > p*c1
yc(k)=m*((c1-r(k))/(1-p)^2)*(1+(r(k)/c1)-(2*p));
dx(k)=((2*m)/(1-p)^2)*(p-(r(k)/c1));
end
%upper and lower limits of the airfoil (xu,yu) ; (xl,yl)
angle=atan(dx(k));
xu(k)=r(k)-yt(k)*-sin(angle);
yu(k)=yc(k)+yt(k)*cos(angle);
xl(k)=r(k)+yt(k)*sin(angle);
yl(k)=yc(k)-yt(k)*cos(angle);
end
%% plot of airfoil
clf;
plot(xu,yu,xl,yl,r,yc);
xlabel('x'),ylabel('y')
title( ['NACA Airfoil ', NACA]);
grid on; axis equal
%% --------------------------------------------------------------------------------------%
% Given data vectors X and Y.
% Want to rotate the data through angle "ang" about rotation center Xc, Yc
X = [xl,fliplr(xu)]; % watch this
Y = [yl,fliplr(yu)];
% Specify the coordinates of the center of rotation
Xc = 0.25 ; % Rotate about the 1/4 chord point
Yc = 0 ;
% The data is roated in a three-step process
% Step 1) Shift the data to the rotation center
Xs = X - Xc; % shifted data
Ys = Y - Yc;
% Step 2) Rotate the data
angles =[1.047197551 0.756112778 0.546578176 0.413499657 0.327872784 ...
0.269927858 0.228682627 0.198042808 0.174475668 0.155828787 0.140728889 ...
0.128263734 0.117805904 0.108910765 0.101254972 0.094598036];
hold on
for i=1:length(angles)
ang = -angles(i); % in radians ??
Xsr = Xs.*cos(ang) + Ys.*sin(ang); % shifted and rotated data
Ysr = -Xs*sin(ang) + Ys*cos(ang); %
% Step 3) Un-shift the data (back to the original coordinate system)
Xr = Xsr + Xc; % Rotated data
Yr = Ysr + Yc;
plot(Xr, Yr)
end % for
hold off
2 Comments
Jorge Cabrera
on 18 Apr 2019
David Wilson
on 19 Apr 2019
Edited: David Wilson
on 19 Apr 2019
Like this?
I've rotated each of your 16 profiles by the corresponding 16 angles. I assume this is what you want. By the way, you do have a chord length of zero! Seems a bit strange.
I've rotated the aerofoils about 1/4 chord length of the actual chord. Looking at your plot above, I see that you have probably rotated about the point 1/4 max chord length. It's an easy change to make if that's what you want.
I've also lightly coloured the foils to make them look nicer.
Code is now as below. I've added a structure to contain each of the profiles. This could be in a rectangular matrix, but you might in the future want to add other fields. (BTW, I typically use English spelling, simply to piss some people off.)
%% Rotate NACA aerofoils with adjusting chord lengths
s=num2str(5410); NACA=s; %4 digits
m =str2double(s(1))/100; % m represents first digit out of the scalar
p=str2double(s(2))/10; %p pulls the second digit out of the scalar
%t pulls and represents the third and fourth digit out of the scalar
t=str2double(s(3:4))/100;
% Chord length, but dropping first (zero) entry
c2= [0,0.008232831,0.016465661,0.024698492,0.032931323...
0.041164153,0.049396984,0.057629815,0.065862645...
0.074095476,0.082328307,0.090561137,0.098793968,0.107026798...
0.115259629,0.12349246];
%c2(1) = []; % drop the zero chord length !!
aerofoil = struct('chord',[],'x',[],'y',[],'r',[],'yc',[]); % place holder
for i=1:length(c2)
c1 = c2(i);
r=linspace(0, c1, 250)'; % define x
%yt function including the coefficients
yt =5*t*c1*(.2969*(sqrt(r/c1))+-.1260*(r/c1)+-.3516*(r/c1).^2+.2843*(r/c1).^3+-.1015*(r/c1).^4);
%yt2 =5*t*c2*(.2969*(sqrt(x/c2))+-.1260*(x/c2)+-.3516*(x/c2).^2+.2843*(x/c2).^3+-.1015*(x/c2).^4);
for k = 1:length(r)
if r(k) <= p*c1
yc(k)=m*(r(k)/p^2)*(2*p-(r(k)/c1));
dx(k)=(2*m)/p^2*(p-(r(k)/c1));
elseif r(k) > p*c1
yc(k)=m*((c1-r(k))/(1-p)^2)*(1+(r(k)/c1)-(2*p));
dx(k)=((2*m)/(1-p)^2)*(p-(r(k)/c1));
end
%upper and lower limits of the airfoil (xu,yu) ; (xl,yl)
angle=atan(dx(k));
xu(k)=r(k)-yt(k)*-sin(angle);
yu(k)=yc(k)+yt(k)*cos(angle);
xl(k)=r(k)+yt(k)*sin(angle);
yl(k)=yc(k)-yt(k)*cos(angle);
end % for k
aerofoil(i).chord = c1;
aerofoil(i).x = [xl,fliplr(xu)]; aerofoil(i).y = [yl,fliplr(yu)];
aerofoil(i).r = r; aerofoil(i).yc = yc;
end
%% plot all airfoils (un-rotated)
clf; hold on
for i=1:length(aerofoil)
x = aerofoil(i).x; y = aerofoil(i).y;
r = aerofoil(i).r; yc = aerofoil(i).yc;
plot(x,y, r, yc, '--')
% plot(xu,yu,xl,yl,r,yc);
end
hold off
xlabel('x'),ylabel('y')
title( ['NACA Airfoil ', NACA]);
grid on; axis equal
%% --------------------------------------------------------------------------------------%
Yc = 0 ; % always rotate about ycenter = 0;
% Step 2) Rotate the data
angles =[1.047197551 0.756112778 0.546578176 0.413499657 0.327872784 ...
0.269927858 0.228682627 0.198042808 0.174475668 0.155828787 0.140728889 ...
0.128263734 0.117805904 0.108910765 0.101254972 0.094598036];
clf;
cmap = winter(length(angles));
hold on
for i=length(aerofoil):-1:1
% Given data vectors X and Y.
X = aerofoil(i).x; Y = aerofoil(i).y;
r = aerofoil(i).r; yc = aerofoil(i).yc;
% Want to rotate the data through angle "ang" about rotation center Xc, Yc
% Specify the coordinates of the center of rotation
Xc = 0.25*r(end) ; % Rotate about the 1/4 chord point
% The data is roated in a three-step process
% Step 1) Shift the data to the rotation center
Xs = X - Xc; % shifted data
Ys = Y - Yc;
ang = -angles(i); % in radians ??
Xsr = Xs.*cos(ang) + Ys.*sin(ang); % shifted and rotated data
Ysr = -Xs*sin(ang) + Ys*cos(ang); %
% Step 3) Un-shift the data (back to the original coordinate system)
Xr = Xsr + Xc; % Rotated data
Yr = Ysr + Yc;
plot(Xr, Yr, 'color',cmap(i,:),'linewidth',1);
h=patch(Xr, Yr, cmap(i,:),'FaceAlpha',0.2);
end % for
hold off
axis equal
grid
More Answers (0)
Categories
Find more on Airfoil tools in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
