85 lines
2.4 KiB
Python
85 lines
2.4 KiB
Python
import numpy as np
|
|
from input.user_input import *
|
|
|
|
|
|
def cd_Palumbo(fr, re, a_top):
|
|
k_cd = 19.8 # +/- 1.2
|
|
k1 = 1.50 * 10 ** 7 # +/- 0.78
|
|
k2 = 0.86 * 10 ** (-7) # +/- 0.43
|
|
a_top0 = np.pi / 4 * 1.383 ** 2 * V_design ** (2/3) # a_top0 for zero-pressure shape
|
|
|
|
value = 0.2 * k_cd/fr * (k1/re + k2 * re) * a_top / (a_top0 ** (3/2))
|
|
|
|
if value > 1.6:
|
|
value = 1.6
|
|
elif value <= 0.1:
|
|
value = 0.1
|
|
|
|
return value
|
|
|
|
|
|
def cd_PalumboLow(fr, re, a_top):
|
|
k_cd = 19.8 - 1.2
|
|
k1 = (1.50 - 0.78) * 10 ** 7
|
|
k2 = (0.86 - 0.43) * 10 ** (-7)
|
|
a_top0 = np.pi / 4 * 1.383 ** 2 * V_design ** (2/3) # a_top0 for zero-pressure shape
|
|
|
|
value = 0.2 * k_cd/fr * (k1/re + k2 * re) * a_top / (a_top0 ** (3/2))
|
|
|
|
if value > 1.6:
|
|
value = 1.6
|
|
elif value <= 0.1:
|
|
value = 0.1
|
|
|
|
return value
|
|
|
|
|
|
def cd_PalumboHigh(fr, re, a_top):
|
|
k_cd = 19.8 + 1.2
|
|
k1 = (1.50 + 0.78) * 10 ** 7
|
|
k2 = (0.86 + 0.43) * 10 ** (-7)
|
|
a_top0 = np.pi / 4 * 1.383 ** 2 * V_design ** (2/3) # a_top0 for zero-pressure shape
|
|
|
|
value = 0.2 * k_cd/fr * (k1/re + k2 * re) * a_top / (a_top0 ** (3/2))
|
|
|
|
if value > 1.6:
|
|
value = 1.6
|
|
elif value <= 0.1:
|
|
value = 0.4
|
|
|
|
return value
|
|
|
|
|
|
def drag(c_d, rho_air, A_proj, v):
|
|
return 0.5 * c_d * rho_air * A_proj * v ** 2
|
|
|
|
|
|
def cd_sphere(Re):
|
|
"This function computes the drag coefficient of a sphere as a function of the Reynolds number Re."
|
|
# Curve fitted after fig . A -56 in Evett and Liu: "Fluid Mechanics and Hydraulics"
|
|
|
|
from numpy import log10, array, polyval
|
|
|
|
if Re <= 0.0:
|
|
CD = 0.0
|
|
elif Re > 8.0e6:
|
|
CD = 0.2
|
|
elif Re > 0.0 and Re <= 0.5:
|
|
CD = 24.0 / Re
|
|
elif Re > 0.5 and Re <= 100.0:
|
|
p = array([4.22, -14.05, 34.87, 0.658])
|
|
CD = polyval(p, 1.0 / Re)
|
|
elif Re > 100.0 and Re <= 1.0e4:
|
|
p = array([-30.41, 43.72, -17.08, 2.41])
|
|
CD = polyval(p, 1.0 / log10(Re))
|
|
elif Re > 1.0e4 and Re <= 3.35e5:
|
|
p = array([-0.1584, 2.031, -8.472, 11.932])
|
|
CD = polyval(p, log10(Re))
|
|
elif Re > 3.35e5 and Re <= 5.0e5:
|
|
x1 = log10(Re / 4.5e5)
|
|
CD = 91.08 * x1 ** 4 + 0.0764
|
|
else:
|
|
p = array([-0.06338, 1.1905, -7.332, 14.93])
|
|
CD = polyval(p, log10(Re))
|
|
return CD
|