asc/senior-design
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity
Files
main
senior-design/grimace/Grimace.m
asc 5bacd30bad reorganize file structure
2010-04-17 12:00:00 -05:00

349 lines
9.8 KiB
Matlab
Executable File
Raw Permalink Blame History

%Grimace
%all alfas in DEGREES!
%by ascorrea, with help from Jacob Huffman
%"filler" variable just makes it easier to switch between Ronald and Grimace
function Grimace(filler, pranformfile, pranformConstraint)
%LOCATIONS--------------------------------------------------------------]
rootLOCATION='\\ad.uiuc.edu\ae\correa2\Desktop\AE441\AE441 - Home Edition\Matlab code\';
airfoilLOCATION='\\ad.uiuc.edu\ae\correa2\Desktop\AE441\AE441 - Home Edition\Airfoils';
pranformLOCATION='\\ad.uiuc.edu\ae\correa2\Desktop\AE441\AE441 - Home Edition\Planforms';
pranformConstraintLOCATION='\\ad.uiuc.edu\engr\ews\homes\desktop\AE441\AE441 - Home Edition\Planform Constraints';
paramLOCATION='\\ad.uiuc.edu\engr\ews\homes\desktop\AE441\AE441 - Home Edition\Parameter Files';
%-----------------------------------------------------------------------]
n=50;
h1=13.33;
h2=10.5;
cd(pranformLOCATION)
[spans, chords, offsets, airfoil]=textread(pranformfile, '%f %f %f %s');
cd(pranformConstraintLOCATION)
[spansCONST, chordsCONST, offsetsCONST, heightCONST]=textread(pranformConstraint, '%f %f %f %f');
cd(rootLOCATION)
cd(airfoilLOCATION)
[planformLE planformTE]=plotpranform(spans, chords, offsets, n);
[constraintLE constraintTE]=plotconstraint(spansCONST, chordsCONST, offsetsCONST, n);
plotheight=plotThick(spans, chords, offsets, airfoil, h1, n);
plotheight2=plotThick(spans, chords, offsets, airfoil, h2, n);
planformLE=[flipud(planformLE); -planformLE(:,1) planformLE(:,2)];
planformTE=[flipud(planformTE); -planformTE(:,1) planformTE(:,2)];
constraintLE=[flipud(constraintLE); -constraintLE(:,1) constraintLE(:,2)];
constraintTE=[flipud(constraintTE); -constraintTE(:,1) constraintTE(:,2)];
plotheight=[flipud(plotheight); -plotheight(:,1) plotheight(:,2)];
plotheight2=[flipud(plotheight2); -plotheight2(:,1) plotheight2(:,2)];
clf
hold on
plot(planformLE (:,1), planformLE(:,2),'k','LineWidth',2)
plot(planformTE (:,1), planformTE(:,2),'k','LineWidth',2)
plot(constraintLE (:,1), constraintLE(:,2),'LineWidth',1.5)
plot(constraintTE (:,1), constraintTE(:,2),'LineWidth',1.5)
plot(plotheight(:,1),plotheight(:,2),'r','LineWidth',1)
plot(plotheight2(:,1),plotheight2(:,2),'r--','LineWidth',1)
axis equal
axis off
hold off
cd(pranformLOCATION)
filename=regexprep(pranformfile, '.plan', 'grim.tiff');
print ('-dtiff', filename)
cd(rootLOCATION)
function [plotLE plotTE]=plotpranform(spans, chords, offsets, n)
Chords = findChords(n, chords, spans);
%plot planform
for i = 1:(size(chords)-1)
dist=0;
for j = 1:n
dist = (spans(i+1)-spans(i))./n;
if i==1
xplotLE(i,j+1)=spans(1)+dist*(j);
else
xplotLE(i,1)=xplotLE(i-1,end);
xplotLE(i,j+1)=xplotLE(i-1,end)+dist*(j);
end
fitX(1)=spans(i);
fitX(2)=spans(i+1);
fitY(1)=offsets(i);
fitY(2)=offsets(i+1);
a = polyfit(fitX, fitY,1);
yplotLE(i,j) = a(1)*xplotLE(i,j)+a(2);
end
end
plotLE=zeros(1,2);
plotTE=zeros(1,2);
for i=1:(size(chords)-1)
for j=1:n
plotLE=[plotLE; xplotLE(i,j) -yplotLE(i,j)];
plotTE=[plotTE; xplotLE(i,j) -(yplotLE(i,j)+Chords(i,j))];
end
end
plotTE(1,:)=[];
plotLE(1,:)=[];
end
function [plotLE plotTE]=plotconstraint(spans, chords, offsets, n)
%plot constraint
for i = 1:(size(chords)-1)
for j = 1:n
dist = (spans(i+1)-spans(i))./n;
if i==1
xplotLE(i,j+1)=spans(1)+dist*(j);
else
xplotLE(i,1)=xplotLE(i-1,end);
xplotLE(i,j+1)=xplotLE(i-1,end)+dist*(j);
end
yplotLE(i,j) = offsets(i+1);
end
end
plotLE=zeros(1,2);
plotTE=zeros(1,2);
for i=1:(size(chords)-1)
for j=1:n
plotLE=[plotLE; xplotLE(i,j) -yplotLE(i,j)];
plotTE=[plotTE; xplotLE(i,j) (-yplotLE(i,j)-chords(i+1))];
end
end
plotTE(1,:)=[];
plotLE(1,:)=[];
end
function plotThick=plotThick(spans, chords, offsets, airfoil, h, n)
n=n/10;
Chords = findChords(n, chords, spans);
plotThickLE=zeros(1,2);
plotThickTE=zeros(1,2);
for i = 1:(size(chords)-1)
for j = 1:n
Chords(i,j);
dist = (spans(i+1)-spans(i))./n;
if i==1
xplotLE(i,j+1)=spans(1)+dist*(j);
else
xplotLE(i,1)=xplotLE(i-1,end);
xplotLE(i,j+1)=xplotLE(i-1,end)+dist*(j);
end
fitX(1)=spans(i);
fitX(2)=spans(i+1);
fitY(1)=offsets(i);
fitY(2)=offsets(i+1);
a = polyfit(fitX, fitY,1);
yplotLE(i,j) = a(1)*xplotLE(i,j)+a(2);
chordloc=linspace(0, Chords(i,j));
for k=1:length(chordloc)-1
thicknessLE(k)=findthickness(airfoil{i}, Chords(i,j), chordloc(k));
clc
fprintf('%g:%u:%u:%u/%g:%u:%u:%u',h,i,j,k,h,(length(chords)-1),n,(length(chordloc)-1))
if thicknessLE(k)>=h
plotThickLE=[plotThickLE; xplotLE(i,j) -(chordloc(k)+yplotLE(i,j))];
break
end
end
for k=1:length(chordloc)-1
thicknessTE(k)=findthickness(airfoil{i}, Chords(i,j), chordloc(length(chordloc)-k));
if thicknessTE(k)>=h
plotThickTE=[plotThickTE; xplotLE(i,j) -(chordloc(length(chordloc)-k)+yplotLE(i,j))];
break
end
end
end
end
plotLE=zeros(1,2);
plotTE=zeros(1,2);
for i=1:(size(chords)-1)
for j=1:n
plotLE=[plotLE; xplotLE(i,j) -yplotLE(i,j)];
plotTE=[plotTE; xplotLE(i,j) -(yplotLE(i,j)+Chords(i,j))];
end
end
plotTE(1,:)=[];
plotLE(1,:)=[];
plotThickLE(1,:)=[];
plotThickTE(1,:)=[];
plotThick=[plotThickLE; flipud(plotThickTE)];
end
function thickness=findthickness(airfoil, Chord, xLoc)
Airfoil = load(airfoil);
Airfoil = Airfoil*Chord;
% blah....goodnight
AirfoilX=Airfoil(:,1);
[min_diff, k]=min(abs(AirfoilX));
for j=1:k
AirfoilYupper(j)=Airfoil(j,2);
AirfoilXupper(j)=Airfoil(j,1);
end
for j=1:length(Airfoil)-k
AirfoilYlower(j)=Airfoil(length(Airfoil)-(j-1), 2);
AirfoilXlower(j)=Airfoil(length(Airfoil)-(j-1), 1);
end
AirfoilYupper=fliplr(AirfoilYupper);
AirfoilXupper=fliplr(AirfoilXupper);
AirfoilYlower=fliplr(AirfoilYlower);
AirfoilXlower=fliplr(AirfoilXlower);
AirfoilYlower=[0 AirfoilYlower];
AirfoilXlower=[0 AirfoilXlower];
[min_diff, k]=min(abs(AirfoilX));
[min_diff, kU]=min(abs(AirfoilXupper-xLoc));
[min_diff, kL]=min(abs(AirfoilXlower-xLoc));
if xLoc<AirfoilXupper(kU)
xFitU(1)=AirfoilXupper(kU-1);
xFitU(2)=AirfoilXupper(kU);
yFitU(1)=AirfoilYupper(kU-1);
yFitU(2)=AirfoilYupper(kU);
else
xFitU(1)=AirfoilXupper(kU);
xFitU(2)=AirfoilXupper(kU+1);
yFitU(1)=AirfoilYupper(kU);
yFitU(2)=AirfoilYupper(kU+1);
end
if xLoc<AirfoilXlower(kL)
xFitL(1)=AirfoilXlower(kL-1);
xFitL(2)=AirfoilXlower(kL);
yFitL(1)=AirfoilYlower(kL-1);
yFitL(2)=AirfoilYlower(kL);
else
xFitL(1)=AirfoilXlower(kL);
xFitL(2)=AirfoilXlower(kL+1);
yFitL(1)=AirfoilYlower(kL);
yFitL(2)=AirfoilYlower(kL+1);
end
aU = polyfit(xFitU, yFitU,1);
aL = polyfit(xFitL, yFitL,1);
% yplotLE(i,j) = a(1)*xplotLE(i,j)+a(2);
% thickness = AirfoilYupper(xLoc_indu)-AirfoilYlower(xLoc_indl);
thickness = (aU(1)*xLoc+aU(2))-(aL(1)*xLoc+aL(2));
end
function Chords = findChords( n, Chord, X )
for( i = 1:(length(Chord)-1) )
Xtemp = linspace(X(i),X(i+1),n);
Chords(i,:) = (Chord(i+1)-Chord(i))./(X(i+1)-X(i)).*(Xtemp-X(i))+Chord(i);
end
end
function [wakeplot]=findwake(airfoil, chords, spans, offsets, n, alt, v_mach, laminarflow, width)
n=n/10;
[q v_fps dynvisc]=findconstants(alt, v_mach);
Chords = findChords(n, chords, spans);
thicknessplot=zeros(100,1);
wakeplotX=zeros(100,1);
wakeplotY=zeros(100,1);
width=width/2;
% for i = 1:(size(chords)-1)
for i=1:1:(size(chords)-1)
for j = 1:n
% wakecoord(j,:)=linspace(0, Chords(i,j), n);
dist = (spans(i+1)-spans(i))./n;
y(j)=spans(i)+j*dist;
if i==1
xplotLE(i,j+1)=spans(1)+dist*(j);
else
xplotLE(i,1)=xplotLE(i-1,end);
xplotLE(i,j+1)=xplotLE(i-1,end)+dist*(j);
end
fitX(1)=spans(i);
fitX(2)=spans(i+1);
fitY(1)=offsets(i);
fitY(2)=offsets(i+1);
a = polyfit(fitX, fitY,1);
yplotLE(i,j) = a(1)*xplotLE(i,j)+a(2);
x=linspace(0, Chords(i,j));
for k = 1:length(x)-1
thickness(k,j)=findthickness(airfoil{i}, Chords(i,j), x(k));
end
wakeX(k,j)=x(k)+yplotLE(i,j);
wakeY(k,j)=y(j);
end
for k = 1:100
wakeX(k,j)=x(k)+yplotLE(i,j);
wakeY(k,j)=y(j);
end
thickness=[thicknessplot, wake];
wakeplotX=[wakeplotX, wakeX];
wakeplotY=[wakeplotY, wakeY];
end
end
thicknessplot(:,1)=[ ];
wakeplotX(:,1)=[ ];
wakeplotY(:,1)=[ ];
wakeplot(1,:)=[0];
flipud(thicknessplot);
surf(wakeplotY, wakeplotX, thicknessplot)
end
function [q v_fps dynvisc]=findconstants(alt, v_mach)
if alt < 36152
temp=59-.00356*alt;
press=2116*((temp+459.7)/518.6)^5.256;
elseif alt < 82345
temp=-70;
press=473.1*exp(1.73-.000048*alt);
else
print('Program Not Valid for given altitude')
end
%slpcuft: clugs ber cubic foot
density_slpcuft=press/(1718*(temp+459.7));
v_fps=v_mach*sqrt(1.4*1718*(temp+459.7));
q=finddynpress(alt, v_mach);
S = 110.4;
beta = 1.458e-6;
temp_k=[(temp-32) / (1.8)] + 273 ;
viscosity_pascalseconds = (beta*temp_k^(1.5))/(temp_k+S);
%kgpmcu:kilograms per meter cubed
%msqpse:meters squared per second
%sqftpsec:square feet per second
density_kgpmcu=density_slpcuft*515.4;
dynamic_viscosity_msqpsec=viscosity_pascalseconds/density_kgpmcu;
dynamic_viscosity_sqftpsec=dynamic_viscosity_msqpsec*10.76;
dynvisc=dynamic_viscosity_sqftpsec;
density=density_slpcuft;
end
end
Reference in New Issue View Git Blame Copy Permalink