📄 mainprogram.m
字号:
% REAL MT(0:10),Jar(0:10),Jmin,Jmax,Kgold,MTUC(0:10)
% DIMENSION ALPHA(0:25,0:25), RE(0:25,0:25), Z(0:1,0:25,0:25)
% DIMENSION GK(0:10,0:10), S(0:10), GP(0:10),CbyR(0:9)
% DIMENSION X(0:9),R(0:9),BETA(0:9),C(0:9),T(0:9),DUM(0:9)
% DIMENSION DT(0:11,0:11),DQ(0:11,0:11),F(0:10,0:10)
% DIMENSION CT(0:10),CP(0:10),ETA(0:10)
% CHARACTER*20 AFds(10)
% CHARACTER*20 PROPNAME,AIRFOILNAME,AFNAME,DUMMYNAME
%
% CHARACTER*3 S28
% CHARACTER*2 Goldstr
% CHARACTER*1 A1
% INTEGER AF,RA,Ufixed,ROWS,COLS,F1,F2,F3,Fzero
% COMMON /logunits/imon,ikey
% common /propel/af,pitch
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear
global MT Jar Jmin Jmax Kgold MTUC
global ALPHA RE Z
global GK S GP CbyR
global X R BETA C T DUM
global DT DQ F
global CT CP ETA
global AFds
global PROPNAME AIRFOILNAME AFNAME DUMMYNAME
global S28
global Goldstr
global A1
global AF RA LDA Ufixed ROWS COLS F1 F2 F3 F4 F5 Fzero
global af pitch
global Azero
global DR RD MF Czero C1 C2 P1 P2 P3%单位变换系数
global H A MPH RPM N %输入高度，音速，速度转速
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
X=[0.30,0.45,0.60,0.70,0.75,0.80,0.85,0.90,0.95,1.0];
%
%************************* 变量初始化 **********************************
%----1101
for I=1:26
for J=1:26
ALPHA(I,J)=0.0;%攻角矩阵
RE(I,J)=0.0;%雷诺数矩阵
for JJ=1:2
Z(JJ,I,J)=0.0;%翼型气动系数矩阵
end
end
end
%----1101
%----1102
for I=1:12
for J=1:12
DT(I,J)=0.0;
DQ(I,J)=0.0;
end
end
%----1102
%----1103
for I=1:11
MT(I)=0.0;
Jar(I)=0.0;
S(I)=0.0;
GP(I)=0.0;
CT(I)=0.0;
CP(I)=0.0;
ETA(I)=0.0;
for J=1:11
F(I,J)=0.0;%拟合系数矩阵
end
end
%----1103
%----1104
for I=1:10
R(I)=0.0;
BETA(I)=0.0;
C(I)=0.0;
T(I)=0.0;
DUM(I)=0.0;
CbyR(I)=0.0;
end
%----1104
%
%***************************** 常数赋值 *********************************
imon=6;
ikey=5;
Fzero=0; %（厚度……标志）
F1=0;
F2=0;
F3=0;
F4=0;
F5=0;
pi;%PI=4.*ATAN(1.0)
DR=pi/180.0;%degree-->radian
RD=180.0/pi;%radian-->degree
Czero=(pi^3.0)/8.0;
C1=360.0;
C2=5.E-4;
MF=1.46667;%mile per hour-->feet per second
P1 = 10.0;
P2 = 100.0;
P3 = 1000.0;
%
%************************* 输入可用翼型信息**************************
AIRFOILS;%已是调用了文件AIRFOILS，OPEN(UNIT=1,FILE='AIRFOILS.d',STATUS = 'OLD')
NAIRFOILS = 5;
AFds(1,:);
AFds(2,:);
AFds(3,:);
AFds(4,:)';
AFds(5,:)';
%
%************************** 从文件中输入数据 **************************************
%---1__continue
propeller;%***********OPEN(23, FILE = 'propeler.dat', STATUS = 'OLD')
fprintf('****************** Openning the file :propeller.m********************\n');
PROPNAME
N;
NB;
%
%---2循环
for I=1:N
R(I)=propdata(I,1); % DO 2 I = 0, N-1
CbyR(I)=propdata(I,2);
T(I)=propdata(I,3); %READ(23,*) R(I), CbyR(I), T(I)
BETA(I) = atan(pitch/(2*Radius)/(3.14158*R(I)))*180./3.14159;
C(I)=CbyR(I)*Radius;
end
fprintf('The number of blades is: %4i \n',NB);
fprintf('*********************矩阵R*****************************\n');
R
fprintf('*********************请敲回车确认*****************************\n');
pause;
fprintf('*********************矩阵CbyR*****************************\n');
CbyR
fprintf('*********************请敲回车确认*****************************\n');
pause;
fprintf('*********************矩阵Beta*****************************\n');
BETA
fprintf('*********************请敲回车确认*****************************\n');
pause;
%write(*,*) 'R = ',R(i),' CbyR = ',CbyR(i),' beta = ',beta(i)
% read(*,*)
F6=0;
AF=0;
while F6<=0
fprintf('******************Which airfoil section do you want to use****************\n');
fprintf('\n');
fprintf(' (1) flatplat\n');
fprintf(' (2) ClarkY\n');
fprintf(' (3) RAF6\n');
fprintf(' (4) NACA44Lo\n');
fprintf(' (5) Symmetri\n');
AF=input('******************请输入翼型标号****************\n');
if ((AF==1)|(AF==2)|(AF==3)|(AF==4)|(AF==5)) %如果正确输入翼型，F6=1
F6=1;
fprintf('你选择的翼型标号为%4i\n',AF);
end
if F6==0
fprintf('*********************** 错误输入，请重新输入！！！********************\n');
end
end
AIRFOILNAME=AFds(AF,:)
fprintf('*********************请敲回车确认*****************************\n');
pause;
%*******************************选择分析模型,CL/CD ADJUSTMENTS *********************************************
%---17
FLAGRA=0;
while FLAGRA<=0
fprintf('******************Select desired refinement of analysis: ****************\n');
fprintf('\n');
fprintf(' (1) Analysis by simple blade element theory\n');
fprintf(' (2) Analysis including induced velocity\n');
fprintf(' (3) Analysis including induced velocity and tip losses\n');
fprintf(' (1-3)?\n');
RA=input('******************请输入分析模型标号****************\n');
if ((RA==1)|(RA==2)|(RA==3)) %如果正确分析模型， FLAGRA=1
FLAGRA=1;
fprintf('你选择的分析模型标号为%4i\n',RA);
end
if ((RA~=1)&(RA~=2)&(RA~=3))
fprintf('********************** Incorrect Entry. Please Try Again.！！！*******************\n');
FLAGRA=0;
end
end
%---17
%---16
FLAGLDA=0;
while FLAGLDA<=0
fprintf('******************Select one of these Lift/Drag coefficient adjustments: ****************\n');
fprintf('\n');
fprintf(' (1) No Cl/Cd adjustments\n');
fprintf(' (2) Mach number adjustments\n');
fprintf(' (3) Reynolds number adjustments\n');
fprintf(' (4) Mach and Reynolds number adjustments\n');
fprintf(' (1-4)?\n');
LDA=input('******************请输入Lift/Drag coefficient调整模型标号****************\n');
if ((LDA==1)|(LDA==2)|(LDA==3)|(LDA==4)) %如果正确分析模型， FLAGRA=1
FLAGLDA=1;
fprintf('你选择的调整模型标号为%4i\n',LDA);
end
if ((LDA~=1)&(LDA~=2)&(LDA~=3)&(LDA~=4))
fprintf('********************** Incorrect Entry. Please Try Again.！！！*******************\n');
FLAGLDA=0;
end
end
%---16
%*******************************INPUT ALTITUDE*******************************
%---15
FLAGH=0;
while FLAGH<=0
H=input('******************Select altitude in thousands of feet: (0 - 75): ****************\n');
if ((H<(75.0))&(H>(0.0)))
FLAGH=1;
else
fprintf( '********************Number must be between 0.0 and 75.0 ！！！*************\n');
FLAGH=0;
end
altitude=H;
end
%---15
%调用气体函数
[RHO, NU, A]=ATMOSPROPS(H);
%
% fprintf( 'RPM =%6f\n ', RPM);
% pause;
[LDA,MPH,RPM,S28,Ufixed]=GETMPHRPM(C1,LDA, A, MF, Radius, H);
%
%**************************************************************************
GETJNJX;
%**************************************************************************
%*****************************曲线拟合初始数据*****************************
%
%******** Determine if the input data file includes the tip **********
if (R(N)>0.9999)%该点是桨尖
X(10)=0.99;
N=N-1;%共有10个点,带桨尖点, 抛弃桨尖
else%该点不是桨尖,则添加桨尖
C(N+1)=0.0;%桨尖弦长记为0
T(N+1)=T(N);%桨尖厚度记为T(N-1)
R(N+1)=1.0;%桨尖无量纲半径记为1
BETA(N+1)=(BETA(N)-BETA(N-1))/(R(N)-R(N-1))*(R(N+1)-R(N))+BETA(N);
end %至此，共有10个点,带桨尖点
%
% ******************** Curve Fit Chord ***********************
%---40头
for I=1:N+1 % Do 40 I = 0,N
F(1,I) = R(I); % F(0,I) = R(I)
F(2,I) = C(I); % F(1,I) = C(I)
end % 40 Continue
%---40尾，至此，第一行为无量纲半径，第二行为弦长，下面送给函数 spline 求得各拟合系数
%**************************************************************************
SPLINE;
%**************************************************************************
%---41头
for I=1:10 %Do 41 I = 0,9
FX = X(I);
[FY]=FNEVAL(FX, P3, F);
C(I) = FY;
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -