how do i rotate an airfoil points? and plot.

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

Answer 1

[David Wilson]

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 

[Jorge Cabrera]

Thank you David! plot works fine. But how would i get it to plot with new chord lenght values?

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];

to plot with the same rotated angles

[David Wilson]

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