Load Analysis - Wing Bending moment , Wing Shear force - Schrenk approximation

아래 파이썬 코드는 아래(Following link)의 문서를 따라 하면서 코드로 구성해본 내용이다.

Following http://www.lightaircraftassociation.co.uk/2010/Engineering/Design/schrenk%20approximation.pdf

 

1. The Schrenk approximation

(1) Construct a quarter ellipse with a length equal to the semi-span. The height of the ellipse is given by the equation:

Ellipse Height = 4𝑆/𝜋𝑏 * √(1−(2𝑦 / 𝑏)^2) (At the aircraft center line  Height = 4𝑆/𝜋𝑏 ) 
where  𝑆 = wing area,  𝑏 = wing span,  𝑦 = distance along span from aircraft center-line

(2) Plot this against span location, 𝑦, on a graph, and add to this a plot of the wing chord.

(3) Draw a line averaging the two plots.

 

"이 포스팅은 쿠팡 파트너스 활동의 일환으로, 이에 따른 일정액의 수수료를 제공받습니다."

 

This will generate a curve showing the parameter 𝑐𝐶𝑙𝑎, where 𝑐 = local wing chord ,  𝐶𝑙𝑎 = local lift coefficient for a global lift coefficient ,  𝐶𝐿 = 1.0

아래 코드는 구글 코랩을 이용하여 작성함

import pandas as pd
import matplotlib.pyplot as plt
import math

WingStation = [114, 110, 105, 100, 95, 90, 85, 80, 75, 70, 60, 50, 40, 30, 25, 20, 15,15, 0]
df1 = pd.DataFrame(WingStation, columns=['Wing Station'])

b = 19.0  #span
S = 66.5  #Wing area

a=0
b_y2_list = []
while a < len(df1):
    b_y = df1.iloc[a,0]
    b_y2 = 2 * b_y / (12 * b)
    b_y2_list.append(b_y2)
    a+=1

df2 = pd.DataFrame(b_y2_list, columns = ['b_y2'])
df3 = pd.concat([df1,df2],axis=1)

#Wing chord at each span location
c_len = [2.00, 2.11, 2.24, 2.37,2.50,2.63,2.76,2.89,3.03,3.16,3.42,3.68,3.95,4.21,4.34,4.47,4.61,4.61,5.00]

df4 = pd.DataFrame(c_len, columns = ['CHORD'])
df5 = pd.concat([df3,df4],axis=1)


#Ellipse
i = 0
temp = []

while i < len(df5):
    
  Ellipse = ((4 * S) / (math.pi * b)) * math.sqrt(1 - ((df5.iloc[i,1])**2))
    
  temp.append(Ellipse)
  df6 = pd.DataFrame(temp, columns = ['ELLIPSE'])
    
  i+=1
    
df7 = pd.concat([df5,df6],axis=1)


i = 0
temp = []
temp1 = []
while i < len(df7):
    cCla = ((df7.iloc[i,2]) + (df7.iloc[i,3])) / 2
    Cla = cCla / (df7.iloc[i,2])
    temp.append(cCla)  
    temp1.append(Cla)
    df8 = pd.DataFrame(temp, columns = ['SCRENK'])
    df9 = pd.DataFrame(temp1, columns = ['Cla'])
    i+=1
df8 = pd.concat([df7,df8],axis=1)
df9 = pd.concat([df8,df9],axis=1)


df9.plot(x='Wing Station',y=['CHORD','ELLIPSE', 'SCRENK'],style=['*','--','o'],grid=True)
plt.title("SCHRENK LIFT DISTRIBUTION")
plt.xlabel("SPANWISE STATION(INCHES FROM CENTRE-LINE)")
plt.ylabel("cCla")

df9.plot(x='Wing Station',y=['Cla'],style=['-.'],grid=True)
plt.title("LOCAL LIFT COEFFICIENT Cla FOR UNIT CL VS SPAN")
plt.xlabel("SPANWISE STATION (INCHES FROM CENTRE-LINE)")
plt.ylabel("Cla")

#DELTA Y(IN)
delY = [0,4.00,5.00,5.00,5.00,5.00,5.00,5.00,5.00,5.00,10.00,10,10,10,5,5,5,0,15]
#ELEMENT CHORD(FT)
c_elem = [0,2.05,2.17,2.3,2.43,2.57,2.7,2.83,2.96,3.09,3.29,3.55,3.82,4.08,4.28,4.41,4.54,4.61,4.8]

i=0
A_elem = []
dYCav2 = []
ElemUnit = []
while i < len(delY):
    A_elem.append(round(((delY[i] * 0.08333) * c_elem[i]),3))  #1in = 0.08333ft
    dYCav2.append(round((delY[i] * 0.08333) * (c_elem[i]**2) , 3))
    i+=1

df10 = pd.DataFrame(delY, columns = ['DELTA Y(IN)'])
df11 = pd.DataFrame(c_elem, columns = ['ELEMENT CHORD(FT)'])
df12 = pd.DataFrame(A_elem, columns = ['ELEMENT AREA(SQ FT)'])
df13 = pd.DataFrame(dYCav2, columns = ['DELAT Y * Elem Area'])

df14 = pd.concat([df9,df10,df11,df12,df13],axis=1)

ElemUnitCla = [0]
ElemLift = [0]
LimitShearForce = [0]

L = 5605  # Lift [unit : LB]
S = 66.5  # Wing area [unit : SQ FT]

i=0
while i < len(df14):
  ElemUnitCla.append(round((((df14.iloc[i,5]) + (df14.iloc[(i+1),5])) / 2),3))

  LimitShearForce.append(ElemLift[i] + LimitShearForce[i])  

  i+=1
  ElemLift.append ( round(((L / S) * df14.iloc[i,8] * ElemUnitCla[i]),3)) 

  if i+1 == len(df14):
    del LimitShearForce[0]
    LimitShearForce.append(ElemLift[i] - ElemLift[i])
    break

df15 = pd.DataFrame(ElemUnitCla, columns = ['ELEMENT UNIT Cla(LB)'])
df16 = pd.DataFrame(ElemLift, columns = ['ELEMENT LIFT(LB)'])
df17 = pd.DataFrame(LimitShearForce, columns = ['LIMIT SHEAR FORCE(LB)'])
df18 = pd.concat([df14,df15,df16,df17],axis=1)

df18.loc[17, 'LIMIT SHEAR FORCE(LB)'] = df18.loc[17, 'ELEMENT LIFT(LB)'] - df18.loc[18, 'ELEMENT LIFT(LB)'] #17th Limit shear force value change

i = 1
LimitBendingMoment = [0]
UltimateShearForce = []
UltimateBendingMoment = []

while i < len(df18):
  #print("i th : ", i ,"  LimitBendingMoment[i-1]",LimitBendingMoment[i-1]," | LimitShearForce[i-1]",LimitShearForce[i-1]," | LimitShearForce[i]",LimitShearForce[i])
  #print("deltyY : ",df18.iloc[i-1,0] - df18.iloc[i,0])

  deltaY = df18.iloc[i-1,0] - df18.iloc[i,0]  #delta Y from wing tip
  temp = round(((df18.iloc[i-1,12]) + (df18.iloc[i,12])),0)
  LimitBendingMoment.append (round((( LimitBendingMoment[i-1] + ( (deltaY / 2 * temp) ))),0))

  i+=1

i = 0
while i < len(df18):
  UltimateShearForce.append ( df18.iloc[i,12] * 1.5 )
  UltimateBendingMoment.append ( LimitBendingMoment[i] * 1.5 )
  i+=1


df19 = pd.DataFrame(LimitBendingMoment, columns = ['LIMIT BENDING MOMENT(LB-IN)'])
df20 = pd.DataFrame(UltimateShearForce, columns = ['Ultimate Shear Force(LB)'])
df21 = pd.DataFrame(UltimateBendingMoment, columns = ['Ultimate Bending Moment(LB-IN)'])
df22 = pd.concat([df18,df19,df20,df21],axis=1)

 

df22.plot(x='Wing Station',y='Ultimate Shear Force(LB)',style='-*',grid=True)
plt.title("Aerodynamic Ultimate Shear Force Diagram(without Inertia")
plt.xlabel("SPANWISE STATION(INCHES FROM CENTRE-LINE)")
plt.ylabel("Ultimate Shear Force(LB)")

df22.plot(x='Wing Station',y='Ultimate Bending Moment(LB-IN)',style='-.',grid=True)
plt.title("Aerodynamic Ultimate Bending Moment Diagram(without Inertia")
plt.xlabel("SPANWISE STATION(INCHES FROM CENTRE-LINE)")
plt.ylabel("Ultimate Bending Moment (LB-IN)")