In[1]:= (* 飛行速度[m/s] *)
In[2]:= u=7.5;
In[3]:= (* 飛行高度[m] *)
In[4]:= h0=4.5;
In[5]:= (* 空気密度[kg/m^3] *)
In[6]:= rho=0.119*9.8;
In[7]:= (* 全機有害抵抗係数 *)
In[8]:= cd0=0.02134;
In[9]:= (* 飛行機効率 *)
In[10]:= ef=0.8;
In[11]:= (* 機体質量[kg] *)
In[12]:= m=100;
In[13]:= (* 2次元翼の揚力傾斜[rad^-1] *)
In[14]:= a0=6.088; (* DAE11 *)
In[15]:= (* 重力加速度[m/s^2] *)
In[16]:= g=9.8;
In[17]:= (* 重心位置(mac%) *)
In[18]:= h=0.33;
In[19]:= (* 空力中心位置(mac%) *)
In[20]:= hnw=0.25;
In[21]:= (* 慣性モーメント[kg m^2] *)
In[22]:= ixx=1729;
In[23]:= iyy=117;
In[24]:= izz=1762;
In[25]:= ixz=15;
In[26]:= (********************* 主翼 *********************)
In[27]:= (* スパン [m] *)
In[28]:= bw=29;
In[29]:= (* 平面形 *)
In[30]:=cw[y_]:=Module[{t},t=Abs[y];If[t<5.7,1.102,If[t<9.9,1.102+(t-5.7)*(0.888-1.102)/4.2, 0.888+(t-9.9)*(0.651-0.888)/4.6]]]
In[31]:= (* 面積[m^2] *)
In[32]:= sw=NIntegrate[cw[y],{y,-bw/2,bw/2}]
Out[32]= 28.0001
In[33]:= (* 平均空力翼弦[m] *)
In[34]:= macw=NIntegrate[cw[y]^2,{y,0,bw/2}]*2/sw
Out[34]= 0.988678
In[35]:= (* 平均空力翼弦のスパン方向の位置[m] *)
In[36]:=macwy=Module[{t},For[t=0,t<=bw/2,t+=0.001,If[cw[t]-macw<0,Break[]]];t]
Out[36]= 7.925
In[37]:= (* アスベクト比 *)
In[38]:= arw=bw^2/sw
Out[38]= 30.0356
In[39]:= (* 地面効果を考慮した主翼アスベクト比 *)
In[40]:= areh=(1+33*(h0/bw)^1.5)/(33*(h0/bw)^1.5)*arw
Out[40]= 44.9257
In[41]:= (* 主翼揚力傾斜[rad^-1] *)
In[42]:= aw=a0/(1+a0/(Pi*areh))
Out[42]= 5.83625
In[43]:= (* 上半角[rad] *)
In[44]:= dihedral[y_]:=If[y<7.8,0.01923,0.33734];
In[45]:= (* 抵抗係数 *)
In[46]:= cdw=0.0095;
In[47]:= (******************* 水平尾翼 *******************)
In[48]:= (* スパン[m] *)
In[49]:= bh=2.9;
In[50]:= (* 平面形 *)
In[51]:= ch[y_]:=0.744+Abs[y]*(0.558-0.744)/1.45
In[52]:= (* 面積[m^2] *)
In[53]:= sh=NIntegrate[ch[y],{y,-bh/2,bh/2}]
Out[53]= 1.8879
In[54]:= (* 重心から測った位置[m] *)
In[55]:= lh=4.4;
In[56]:= (* 水平尾翼容積 *)
In[57]:= vhr=lh*sh/(macw*sw)
Out[57]= 0.300066
In[58]:= (* 平均空力翼弦[m] *)
In[59]:= mach=NIntegrate[ch[y]^2,{y,0,bh/2}]*2/sh
Out[59]= 0.655429
In[60]:= (* アスベクト比 *)
In[61]:= arh=bh^2/sh
Out[61]= 4.45469
In[62]:= (* 水平尾翼揚力傾斜[rad^-1] *)
In[63]:= ah=a0/(1+a0/(Pi*arh))
Out[63]= 4.24245
In[64]:= (******************* 垂直尾翼 *******************)
In[65]:= (* スパン[m] *)
In[66]:= bv=2.3;
In[67]:= (* 平面形 *)
In[68]:=cv[y_]:=If[y<0,0.8+y*(0.8-0.65)/0.65,0.8+y*(0.49-0.8)/1.65]
In[69]:= (* 面積[m^2] *)
In[70]:= sv=NIntegrate[cv[y],{y,-0.65,1.65}]
Out[70]= 1.5355
In[71]:= (* 重心から測った位置[m] *)
In[72]:= lv=5.3;
In[73]:= (* 垂直尾翼容積 *)
In[74]:= vvr=lv*sv/(bw*sw)
Out[74]= 0.0100223
In[75]:= (* 平均空力翼弦[m] *)
In[76]:= macv=NIntegrate[cv[y]^2,{y,-0.65,1.65}]/sv
Out[76]= 0.678951
In[77]:= (* アスベクト比 *)
In[78]:= arv=bv^2/sv
Out[78]= 3.44513
In[79]:= (* 揚力傾斜[rad^-1] *)
In[80]:= av=a0/(1+a0/(Pi*arv))
Out[80]= 3.89633
In[81]:= (* 力の作用点の位置[m] *)
In[82]:= zvr=0.5;
In[83]:= zvd=0.5;
In[84]:= (**************************************************)
In[85]:= (* 全機揚力係数 *)
In[86]:= cl=2*m*g/(rho*u^2*sw)
Out[86]= 1.06709
In[87]:= (* 全機抵抗係数 *)
In[88]:= cd=cd0+cl^2/(Pi*ef*areh)
Out[88]= 0.0314248
In[89]:= (* 必要推力[N] *)
In[90]:= nt=rho*u^2*sw*cd/2
Out[90]= 28.8601
In[91]:= (* ダウンウォッシュ特性 *)
In[92]:= deda=2*aw/(Pi*areh)
Out[92]= 0.0827026
In[93]:= (* CLα *)
In[94]:= cla=aw*(1+ah/aw*sh/sw*(1-deda))
Out[94]= 6.09864
In[95]:= (* CDα *)
In[96]:= cda=2*cl*cla/(Pi*ef*areh)
Out[96]= 0.115273
In[97]:= (* 胴体容積[m^3] *)
In[98]:= vf=1;
In[99]:= (* 胴体容積比 *)
In[100]:= vfr=vf/(macw*sw)
Out[100]= 0.0361231
In[101]:= (* サイドウォッシュ特性 *)
In[102]:= dsdb=0;
In[103]:= dsdr=0;
In[104]:= (******************* 縦の安定微係数 *****************)
In[105]:= (* 無次元 *)
In[106]:= cxu=-2*nt/(rho*sw*u^2)-2*cd
Out[106]= -0.0942743
In[107]:= czu=0;
In[108]:= cmu=0;
In[109]:= cxa=cl*(1-2*cla/(Pi*ef*areh))
Out[109]= 0.951815
In[110]:= cza=-cla
Out[110]= -6.09864
In[111]:= cma=aw*((h-hnw)-vhr*ah/aw*(1-deda)+vfr*2/aw)
Out[111]= -0.628587
In[112]:= cma'=-2*vhr*lh/macw*ah*deda
Out[112]= -0.937088
In[113]:= cxq=0;
In[114]:= czq=-2*vhr*ah
Out[114]= -2.54603
In[115]:= cmq=-2*vhr*lh/macw*ah
Out[115]= -11.3308
In[116]:= czde=-sh/sw*ah*1
Out[116]= -0.286046
In[117]:= cmde=-vhr*ah*1
Out[117]= -1.27302
In[118]:= (* 有次元 *)
In[119]:= xu=rho*u*sw/(2*m)*cxu
Out[119]= -0.11544
In[120]:= zu=rho*u*sw/(2*m)*(czu-2*cl)
Out[120]= -2.61333
In[121]:= mu=0;
In[122]:= xa=rho*u^2*sw/(2*m)*cxa
Out[122]= 8.74135
In[123]:= za=rho*u^2*sw/(2*m)*cza
Out[123]= -56.0091
In[124]:= ma=rho*u^2*sw*macw/(2*iyy)*cma
Out[124]= -4.87821
In[125]:= ma'=rho*u*sw*macw^2/(4*iyy)*cma'
Out[125]= -0.479335
In[126]:= xq=0;
In[127]:= zq=rho*u*sw*macw/(4*m)*czq
Out[127]= -1.54118
In[128]:= mq=rho*u*sw*macw^2/(4*iyy)*cmq
Out[128]= -5.79589
In[129]:= zde=rho*u^2*sw/(2*m)*czde
Out[129]= -2.62701
In[130]:= mde=rho*u^2*sw*macw/(2*iyy)*cmde
Out[130]= -9.87935
In[131]:= (**************** 横の安定微係数 ********************)
In[132]:= (* 無次元 *)
In[133]:= cyb=-sv/sw*av*(1-dsdb)
Out[133]= -0.213671
In[134]:=clb=-2*a0/(sw*bw)*NIntegrate[dihedral[y]*cw[y]*y,{y,0,bw/2}]
Out[134]= -0.314195
In[135]:= cnb=vvr*av*(1-dsdb)-2*vfr*macw/bw
Out[135]= 0.0365872
In[136]:= cyp=0;
In[137]:= clp=-4*a0/(sw*bw^2)*NIntegrate[cw[y]*y^2,{y,0,bw/2}]
Out[137]= -0.875301
In[138]:=cnp=-4/(sw*bw^2)*(cl-cda)*NIntegrate[cw[y]*y^2,{y,0,bw/2}]
Out[138]= -0.136847
In[139]:= cyr=sv/sw*av*(2*lv/bw+dsdr)
Out[139]= 0.0781004
In[140]:=clr=8/(sw*bw^2)*cl*NIntegrate[cw[y]*y^2,{y,0,bw/2}]+zvr*sv/(bw*sw)*av*(2*lv/bw+dsdr)
Out[140]= 0.308187
In[141]:=cnr=-8/(sw*bw^2)*cdw*NIntegrate[cw[y]*y^2,{y,0,bw/2}]-vvr*av*(2*lv/bw+dsdr)
Out[141]= -0.0170052
In[142]:= clda=0;
In[143]:= cnda=0;
In[144]:= cydr=sv/sw*av
Out[144]= 0.213671
In[145]:= cyda=0;
In[146]:= cldr=zvd*sv/(bw*sw)*av
Out[146]= 0.00368398
In[147]:= cndr=-vvr*av
Out[147]= -0.0390502
In[148]:= (* 有次元 *)
In[149]:= yb=rho*u^2*sw/(2*m)*cyb
Out[149]= -1.96233
In[150]:= lb=rho*u^2*sw*bw/(2*ixx)*clb
Out[150]= -4.83981
In[151]:= nb=rho*u^2*sw*bw/(2*izz)*cnb
Out[151]= 0.553027
In[152]:= yp=0;
In[153]:= lp=rho*u*sw*bw^2/(4*ixx)*clp
Out[153]= -26.0671
In[154]:= np=rho*u*sw*bw^2/(4*izz)*cnp
Out[154]= -3.99908
In[155]:= yr=rho*u*sw*bw/(4*m)*cyr
Out[155]= 1.38671
In[156]:= lr=rho*u*sw*bw^2/(4*ixx)*clr
Out[156]= 9.17805
In[157]:= nr=rho*u*sw*bw^2/(4*izz)*cnr
Out[157]= -0.496944
In[158]:= xdt=0;
In[159]:= yda=0;
In[160]:= lda=0;
In[161]:= nda=0;
In[162]:= ydr=rho*u^2*sw/(2*m)*cydr
Out[162]= 1.96233
In[163]:= ldr=rho*u^2*sw*bw/(2*ixx)*cldr
Out[163]= 0.0567474
In[164]:= ndr=rho*u^2*sw*bw/(2*izz)*cndr
Out[164]= -0.590257
In[165]:= (* primed derivative *)
In[166]:= lbp=(lb+(ixz/ixx)*nb)/(1-ixz^2/(ixx*izz))
Out[166]= -4.83537
In[167]:= lpp=(lp+(ixz/ixx)*np)/(1-ixz^2/(ixx*izz))
Out[167]= -26.1037
In[168]:= lrp=(lr+(ixz/ixx)*nr)/(1-ixz^2/(ixx*izz))
Out[168]= 9.17441
In[169]:= ldap=(lda+(ixz/ixx)*nda)/(1-ixz^2/(ixx*izz))
Out[169]= 0
In[170]:= ldrp=(ldr+(ixz/ixx)*ndr)/(1-ixz^2/(ixx*izz))
Out[170]= 0.0516305
In[171]:= nbp=(nb+(ixz/izz)*lb)/(1-ixz^2/(ixx*izz))
Out[171]= 0.511864
In[172]:= npp=(np+(ixz/izz)*lp)/(1-ixz^2/(ixx*izz))
Out[172]= -4.2213
In[173]:= nrp=(nr+(ixz/izz)*lr)/(1-ixz^2/(ixx*izz))
Out[173]= -0.418842
In[174]:= ndap=(nda+(ixz/izz)*lda)/(1-ixz^2/(ixx*izz))
Out[174]= 0
In[175]:= ndrp=(ndr+(ixz/izz)*ldr)/(1-ixz^2/(ixx*izz))
Out[175]= -0.589818
前の章 目次に戻る 次の章