Added stability hydrostatics

This commit is contained in:
Jose Luis Cercós pita 2012-03-09 13:38:55 +01:00
parent 58f55b4a2b
commit 223e49af2f
2 changed files with 387 additions and 4 deletions

View File

@ -59,6 +59,7 @@ class Plot(object):
if self.execute():
return
ImageGui.open(self.path + 'volume.png')
ImageGui.open(self.path + 'stability.png')
def createDirectory(self):
""" Create needed folder to write data and scripts.
@ -97,13 +98,16 @@ class Plot(object):
Output.write(" # 3: Wetted surface [m2]\n")
Output.write(" # 4: 1cm triming ship moment [ton m]\n")
Output.write(" # 5: Bouyance center x coordinate\n")
Output.write(" # 6: Floating area\n")
Output.write(" # 7: KBt\n")
Output.write(" # 8: BMt\n")
Output.write(" #\n")
Output.write(" #################################################################\n")
# Print data
for i in range(0,len(drafts)):
draft = drafts[i]
point = Tools.Point(ship,draft,trim)
string = "%f %f %f %f %f\n" % (point.disp, point.draft, point.wet, point.mom, point.xcb)
string = "%f %f %f %f %f %f %f %f\n" % (point.disp, point.draft, point.wet, point.mom, point.xcb, point.farea, point.KBt, point.BMt)
Output.write(string)
# Close file
Output.close()
@ -160,13 +164,28 @@ class Plot(object):
Output.write("set style 1 line linetype 1 linewidth 1 colour rgb (0):(0):(0)\n")
Output.write("set style 2 line linetype 1 linewidth 1 colour rgb (1):(0):(0)\n")
Output.write("set style 3 line linetype 1 linewidth 1 colour rgb (0):(0):(1)\n")
Output.write("set style 4 line linetype 1 linewidth 1 colour rgb (1):(0):(1)\n")
Output.write("set style 4 line linetype 1 linewidth 1 colour rgb (0.1):(0.5):(0.1)\n")
# Write plot call
Output.write("# Plot\n")
Output.write("plot '%s' using 2:1 title '$\\bigtriangleup$' axes x1y1 with lines style 1, \\\n" % (self.dataFile))
Output.write("plot '%s' using 2:1 title 'Draft' axes x1y1 with lines style 1, \\\n" % (self.dataFile))
Output.write(" '' using 3:1 title 'Wetted area' axes x2y1 with lines style 2, \\\n")
Output.write(" '' using 4:1 title '1cm trim moment' axes x3y1 with lines style 3, \\\n")
Output.write(" '' using 5:1 title 'XCB' axes x4y1 with lines style 4\n")
# Prepare second plot
Output.write("set output '%s'\n" % (self.path + 'stability.eps'))
Output.write("# X axis\n")
Output.write("set x2label '\\textit{Floating area} / $\\mathrm{m}^2$'\n")
Output.write("set x2tic\n")
Output.write("set x3label '$KB_{T}$ / $\\mathrm{m}$'\n")
Output.write("set x3tic\n")
Output.write("set x4label '$BM_{T}$ / $\\mathrm{m}$'\n")
Output.write("set x4tic\n")
# Write plot call
Output.write("# Plot\n")
Output.write("plot '%s' using 2:1 title 'Draft' axes x1y1 with lines style 1, \\\n" % (self.dataFile))
Output.write(" '' using 6:1 title 'Floating area' axes x2y1 with lines style 2, \\\n")
Output.write(" '' using 7:1 title '$KB_{T}$' axes x3y1 with lines style 3, \\\n")
Output.write(" '' using 8:1 title '$BM_{T}$' axes x4y1 with lines style 4\n")
# Close file
self.layoutFile = filename
Output.close()
@ -176,6 +195,7 @@ class Plot(object):
""" Calls pyxplot in order to plot an save an image.
@return True if error happens.
"""
# Plot
filename = self.path + 'volume'
comm = "pyxplot %s" % (self.layoutFile)
if os.system(comm):
@ -184,6 +204,16 @@ class Plot(object):
msg = Translator.translate("Plot will not generated\n")
FreeCAD.Console.PrintError(msg)
return True
# Convert volume
comm = "gs -r300 -dEPSCrop -dTextAlphaBits=4 -sDEVICE=png16m -sOutputFile=%s.png -dBATCH -dNOPAUSE %s.eps" % (filename,filename)
if os.system(comm):
msg = Translator.translate("Can't execute ghostscript. Maybe is not installed?\n")
FreeCAD.Console.PrintError(msg)
msg = Translator.translate("Generated image will not converted to png\n")
FreeCAD.Console.PrintError(msg)
return True
# Convert stability
filename = self.path + 'stability'
comm = "gs -r300 -dEPSCrop -dTextAlphaBits=4 -sDEVICE=png16m -sOutputFile=%s.png -dBATCH -dNOPAUSE %s.eps" % (filename,filename)
if os.system(comm):
msg = Translator.translate("Can't execute ghostscript. Maybe is not installed?\n")

View File

@ -355,6 +355,350 @@ def Moment(ship, draft, trim, disp, xcb):
mom1 = -data[1]*data[2]
return mom1 - mom0
def FloatingArea(ship, draft, trim):
""" Calculate ship floating area.
@param ship Selected ship instance
@param draft Draft.
@param trim Trim in degrees.
@return Ship floating area.
"""
angle = math.radians(trim)
sections = Instance.sections(ship)
xCoord = ship.xSection[:]
lines = []
area = 0.0
if not sections:
return 0.0
for i in range(0, len(sections)):
# Get the section
section = sections[i]
if len(section) < 2: # Empty section
lines.append(0.0)
continue
# Get the position of the section
x = xCoord[i]
# Get the maximum Z value
Z = draft - x*math.tan(angle)
# Format section
section = convertSection(section, x, Z)
if not section:
lines.append(0.0)
continue
# Get floating line length
line = 0.0
flag = True # Even lines compute for floating areas, odd no
j = len(section)-1
k = len(section[j])-1
while k>0:
k = k-1
if flag:
y0 = abs(section[j][k-1].y)
y1 = abs(section[j][k].y)
line = line + (y1 - y0)
flag = not flag
if flag: # Central body computation lefts
y = abs(section[j][0].y)
line = line + y
lines.append(2.0*line) # 2x because only half ship is represented
# Add area if proceed
if i > 0:
dx = xCoord[i] - xCoord[i-1]
x = 0.5*(xCoord[i] + xCoord[i-1])
line = 0.5*(lines[i] + lines[i-1])
area = area + line*dx
return area
def KBT(ship, draft, trim, roll=0.0):
""" Calculate ship Keel to Bouyance center transversal distance.
@param ship Selected ship instance
@param draft Draft.
@param trim Trim in degrees.
@param roll Roll angle in degrees.
@return [KBTx, KBTy]: \n
KBTy : TRansversal KB y coordinate \n
KBTz : TRansversal KB z coordinate
"""
angle = math.radians(trim)
rAngle = math.radians(roll)
sections = Instance.sections(ship)
xCoord = ship.xSection[:]
areas = []
vol = 0.0
vMy = 0.0
vMz = 0.0
My = []
Mz = []
kb = [0.0,0.0]
if not sections:
return [[],0.0,0.0]
for i in range(0, len(sections)):
# ------------------------------------
# Board
# ------------------------------------
# Get the section
section = sections[i]
if len(section) < 2: # Empty section
areas.append(0.0)
continue
# Get the position of the section
x = xCoord[i]
# Get the maximum Z value
Z = draft - x*math.tan(angle)
# Format section
aux = convertSection(section, x, Z)
if not aux:
areas.append(0.0)
My.append(0.0)
Mz.append(0.0)
continue
# Correct by roll angle
pts = aux[len(aux)-1]
beam = pts[0].y
for i in range(1,len(pts)):
beam = max(beam, pts[i].y)
Z = Z - beam*math.tan(rAngle)
# Format section
section = convertSection(section, x, Z)
if not section:
areas.append(0.0)
My.append(0.0)
Mz.append(0.0)
continue
# Integrate area
area = 0.0
momy = 0.0
momz = 0.0
for j in range(0, len(section)-1):
for k in range(0, min(len(section[j])-1, len(section[j+1])-1)):
# y11,z11 ------- y01,z01
# | |
# | |
# | |
# y10,z10 ------- y00,z00
y00 = abs(section[j][k].y)
z00 = section[j][k].z
y10 = abs(section[j][k+1].y)
z10 = section[j][k+1].z
y01 = abs(section[j+1][k].y)
z01 = section[j+1][k].z
y11 = abs(section[j+1][k+1].y)
z11 = section[j+1][k+1].z
dy = 0.5*((y00 - y10) + (y01 - y11))
dz = 0.5*((z01 - z00) + (z11 - z10))
y = 0.25*(y00 + y10 + y01 + y11)
z = 0.25*(z01 + z00 + z11 + z10)
area = area + dy*dz
momy = momy + y*dy*dz
momz = momz + z*dy*dz
if(len(section[j]) < len(section[j+1])):
# y01,z01 ------- y11,z11
# | __/
# | __/
# | /
# y00,z00
k = len(section[j])-1
y00 = abs(section[j][k].y)
z00 = section[j][k].z
y01 = abs(section[j+1][k].y)
z01 = section[j+1][k].z
y11 = abs(section[j+1][k+1].y)
z11 = section[j+1][k+1].z
dy = y01 - y11
dz = z01 - z00
y = 0.33*(y00 + y01 + y11)
z = 0.33*(z01 + z00 + z11)
area = area + 0.5*dy*dz
momy = momy + y*0.5*dy*dz
momz = momz + z*0.5*dy*dz
elif(len(section[j]) > len(section[j+1])):
# y01,z01
# | \__
# | \__
# | \
# y00,z00 ------- y10,z10
k = len(section[j+1])-1
y00 = abs(section[j][k].y)
z00 = section[j][k].z
y10 = abs(section[j][k+1].y)
z10 = section[j][k+1].z
y01 = abs(section[j+1][k].y)
z01 = section[j+1][k].z
dy = y00 - y10
dz = z01 - z00
y = 0.33*(y00 + y01 + y10)
z = 0.33*(z01 + z00 + z10)
area = area + 0.5*dy*dz
momy = momy + y*0.5*dy*dz
momz = momz + z*0.5*dy*dz
elif(len(section[j]) == 1):
# y1,z1 -------
# |
# |
# |
# y0,z0 -------
k = 0
y0 = abs(section[j][k].y)
z0 = section[j][k].z
y1 = abs(section[j+1][k].y)
z1 = section[j+1][k].z
dy = 0.5 * (y0 + y1)
dz = z1 - z0
y = 0.25*(y1 + y0)
z = 0.25*(z1 + z0)
area = area + dy*dz
momy = momy + y*dy*dz
momz = momz + z*dy*dz
# ------------------------------------
# StarBoard
# ------------------------------------
# Get the section
section = sections[i]
# Get the maximum Z value
Z = draft - x*math.tan(angle)
# Format section
aux = convertSection(section, x, Z)
if not aux:
areas.append(0.0)
My.append(0.0)
Mz.append(0.0)
continue
# Correct by roll angle
pts = aux[len(aux)-1]
beam = pts[0].y
for i in range(1,len(pts)):
beam = max(beam, pts[i].y)
Z = Z + beam*math.tan(rAngle)
# Format section
section = convertSection(section, x, Z)
if not section:
areas.append(0.0)
My.append(0.0)
Mz.append(0.0)
continue
# Integrate area
for j in range(0, len(section)-1):
for k in range(0, min(len(section[j])-1, len(section[j+1])-1)):
# y11,z11 ------- y01,z01
# | |
# | |
# | |
# y10,z10 ------- y00,z00
y00 = abs(section[j][k].y)
z00 = section[j][k].z
y10 = abs(section[j][k+1].y)
z10 = section[j][k+1].z
y01 = abs(section[j+1][k].y)
z01 = section[j+1][k].z
y11 = abs(section[j+1][k+1].y)
z11 = section[j+1][k+1].z
dy = 0.5*((y00 - y10) + (y01 - y11))
dz = 0.5*((z01 - z00) + (z11 - z10))
y = 0.25*(y00 + y10 + y01 + y11)
z = 0.25*(z01 + z00 + z11 + z10)
area = area + dy*dz
momy = momy - y*dy*dz
momz = momz + z*dy*dz
if(len(section[j]) < len(section[j+1])):
# y01,z01 ------- y11,z11
# | __/
# | __/
# | /
# y00,z00
k = len(section[j])-1
y00 = abs(section[j][k].y)
z00 = section[j][k].z
y01 = abs(section[j+1][k].y)
z01 = section[j+1][k].z
y11 = abs(section[j+1][k+1].y)
z11 = section[j+1][k+1].z
dy = y01 - y11
dz = z01 - z00
y = 0.33*(y00 + y01 + y11)
z = 0.33*(z01 + z00 + z11)
area = area + 0.5*dy*dz
momy = momy - y*0.5*dy*dz
momz = momz + z*0.5*dy*dz
elif(len(section[j]) > len(section[j+1])):
# y01,z01
# | \__
# | \__
# | \
# y00,z00 ------- y10,z10
k = len(section[j+1])-1
y00 = abs(section[j][k].y)
z00 = section[j][k].z
y10 = abs(section[j][k+1].y)
z10 = section[j][k+1].z
y01 = abs(section[j+1][k].y)
z01 = section[j+1][k].z
dy = y00 - y10
dz = z01 - z00
y = 0.33*(y00 + y01 + y10)
z = 0.33*(z01 + z00 + z10)
area = area + 0.5*dy*dz
momy = momy - y*0.5*dy*dz
momz = momz + z*0.5*dy*dz
elif(len(section[j]) == 1):
# y1,z1 -------
# |
# |
# |
# y0,z0 -------
k = 0
y0 = abs(section[j][k].y)
z0 = section[j][k].z
y1 = abs(section[j+1][k].y)
z1 = section[j+1][k].z
dy = 0.5 * (y0 + y1)
dz = z1 - z0
y = 0.25*(y1 + y0)
z = 0.25*(z1 + z0)
area = area + dy*dz
momy = momy - y*dy*dz
momz = momz + z*dy*dz
areas.append(area)
My.append(momy)
Mz.append(momz)
# Add volume & moment if proceed
if i > 0:
dx = xCoord[i] - xCoord[i-1]
x = 0.5*(xCoord[i] + xCoord[i-1])
area = 0.5*(areas[i] + areas[i-1])
momy = 0.5*(My[i] + My[i-1])
momz = 0.5*(Mz[i] + Mz[i-1])
vol = vol + area*dx
vMy = vMy + momy*dx
vMz = vMz + momz*dx
# Compute KBT
kb[0] = vMy / vol
kb[1] = vMz / vol
return kb
def BMT(ship, draft, trim):
""" Calculate ship Bouyance center transversal distance.
@param ship Selected ship instance
@param draft Draft.
@param trim Trim in degrees.
@return BM Bouyance to metacenter height.
"""
nRoll = 2
maxRoll = 7.0
kb0 = KBT(ship,draft,trim)
BM = 0.0
for i in range(0,nRoll):
roll = (maxRoll/nRoll)*(i+1)
kb = KBT(ship,draft,trim,roll)
# * M
# / \
# / \ BM ==|> BM = (BB/2) / tan(alpha/2)
# / \
# *-------*
# BB
BB = [kb[0] - kb0[0], kb[1] - kb0[1]]
BB = math.sqrt(BB[0]*BB[0] + BB[1]*BB[1])
BM = BM + 0.5*BB/math.tan(math.radians(0.5*roll)) / nRoll # nRoll is the weight function
return BM
class Point:
""" Hydrostatics point, that conatins: \n
draft Ship draft [m]. \n
@ -362,7 +706,10 @@ class Point:
disp Ship displacement [ton]. \n
xcb Bouyance center X coordinate [m].
wet Wetted ship area [m2].
mom triming 1cm ship moment [ton m].
mom Triming 1cm ship moment [ton m].
farea Floating area [m2].
KBt Transversal KB height [m].
BMt Transversal BM height [m].
@note Moment is positive when produce positive trim.
"""
def __init__(self, ship, draft, trim):
@ -376,10 +723,16 @@ class Point:
areasData = Displacement(ship,draft,trim)
wettedArea = WettedArea(ship,draft,trim)
moment = Moment(ship,draft,trim,areasData[1],areasData[2])
farea = FloatingArea(ship,draft,trim)
kb = KBT(ship,draft,trim)
bm = BMT(ship,draft,trim)
# Store final data
self.draft = draft
self.trim = trim
self.disp = areasData[1]
self.xcb = areasData[2]
self.wet = wettedArea
self.farea = farea
self.mom = moment
self.KBt = kb[1]
self.BMt = bm