Added the GZ curves computation tool to the console interface

This commit is contained in:
Jose Luis Cercos Pita 2016-01-25 15:15:04 +01:00
parent f0770a7456
commit 499c685668
5 changed files with 248 additions and 178 deletions

View File

@ -35,4 +35,5 @@ from shipHydrostatics.Tools import floatingArea, BMT, mainFrameCoeff
from shipCreateWeight.Tools import createWeight
from shipCreateTank.Tools import createTank
from shipCapacityCurve.Tools import tankCapacityCurve
from shipCreateLoadCondition.Tools import createLoadCondition
from shipCreateLoadCondition.Tools import createLoadCondition
from shipGZ.Tools import gz

View File

@ -154,7 +154,8 @@ class Tank:
return ret_value
def getCoG(self, fp, vol, roll=0.0, trim=0.0):
def getCoG(self, fp, vol, roll=Units.parseQuantity("0 deg"),
trim=Units.parseQuantity("0 deg")):
"""Return the fluid volume center of gravity, provided the volume of
fluid inside the tank.
@ -162,15 +163,14 @@ class Tank:
Keyword arguments:
fp -- Part::FeaturePython object affected.
vol -- Volume of fluid (in m^3).
roll -- Ship roll angle (in degrees).
trim -- Ship trim angle (in degrees).
vol -- Volume of fluid.
roll -- Ship roll angle.
trim -- Ship trim angle.
If the fluid volume is bigger than the total tank one, it will be
conveniently clamped.
"""
# Change the units of the volume, and clamp the value
vol = vol * Units.Metre.Value**3
if vol <= 0.0:
return Vector()
if vol >= fp.Shape.Volume:
@ -187,19 +187,19 @@ class Tank:
return cog
# Get a first estimation of the level
level = vol / fp.Shape.Volume
level = vol.Value / fp.Shape.Volume
# Transform the tank shape
current_placement = fp.Placement
m = current_placement.toMatrix()
m.rotateX(radians(roll))
m.rotateY(-radians(trim))
fp.Placement = m
m.rotateX(roll.getValueAs("rad"))
m.rotateY(-trim.getValueAs("rad"))
fp.Placement = Placement(m)
# Iterate to find the fluid shape
for i in range(COMMON_BOOLEAN_ITERATIONS):
shape = self.getVolume(fp, level, return_shape=True)
error = (vol - shape.Volume) / fp.Shape.Volume
error = (vol.Value - shape.Volume) / fp.Shape.Volume
if abs(error) < 0.01:
break
level += error
@ -222,8 +222,8 @@ class Tank:
fp.Placement = current_placement
p = Part.Point(cog)
m = Matrix()
m.rotateY(radians(trim))
m.rotateX(-radians(roll))
m.rotateY(trim.getValueAs("rad"))
m.rotateX(-roll.getValueAs("rad"))
p.rotate(Placement(m))
return Vector(p.X, p.Y, p.Z)

View File

@ -48,18 +48,22 @@ class TaskPanel:
form.n_points = self.widget(QtGui.QSpinBox, "NumPoints")
form.var_trim = self.widget(QtGui.QCheckBox, "VariableTrim")
rolls = []
roll = Units.Quantity(Locale.fromString(
form.angle.text())).getValueAs('deg').Value
roll = Units.Quantity(Locale.fromString(form.angle.text()))
n_points = form.n_points.value()
var_trim = form.var_trim.isChecked()
rolls = []
for i in range(n_points):
rolls.append(roll * i / float(n_points - 1))
gzs, drafts, trims = Tools.solve(self.ship,
self.weights,
self.tanks,
rolls,
form.var_trim.isChecked())
points = Tools.gz(self.lc, rolls, var_trim)
gzs = []
drafts = []
trims = []
for p in points:
gzs.append(p[0].getValueAs('m').Value)
drafts.append(p[1].getValueAs('m').Value)
trims.append(p[2].getValueAs('deg').Value)
PlotAux.Plot(rolls, gzs, drafts, trims)
@ -175,99 +179,6 @@ class TaskPanel:
continue
except ValueError:
continue
# Extract the weights and the tanks
weights = []
index = 6
while True:
try:
ws = doc.getObjectsByLabel(obj.get('A{}'.format(index)))
except ValueError:
break
index += 1
if len(ws) != 1:
if len(ws) == 0:
msg = QtGui.QApplication.translate(
"ship_console",
"Wrong Weight label! (no instances labeled as"
"'{}' found)",
None,
QtGui.QApplication.UnicodeUTF8)
App.Console.PrintError(msg + '\n'.format(
obj.get('A{}'.format(index - 1))))
else:
msg = QtGui.QApplication.translate(
"ship_console",
"Ambiguous Weight label! ({} instances labeled as"
"'{}' found)",
None,
QtGui.QApplication.UnicodeUTF8)
App.Console.PrintError(msg + '\n'.format(
len(ws),
obj.get('A{}'.format(index - 1))))
continue
w = ws[0]
try:
if w is None or not w.PropertiesList.index("IsWeight"):
msg = QtGui.QApplication.translate(
"ship_console",
"Invalid Weight! (the object labeled as"
"'{}' is not a weight)",
None,
QtGui.QApplication.UnicodeUTF8)
App.Console.PrintError(msg + '\n'.format(
len(ws),
obj.get('A{}'.format(index - 1))))
continue
except ValueError:
continue
weights.append(w)
tanks = []
index = 6
while True:
try:
ts = doc.getObjectsByLabel(obj.get('C{}'.format(index)))
dens = float(obj.get('D{}'.format(index)))
level = float(obj.get('E{}'.format(index)))
except ValueError:
break
index += 1
if len(ts) != 1:
if len(ts) == 0:
msg = QtGui.QApplication.translate(
"ship_console",
"Wrong Tank label! (no instances labeled as"
"'{}' found)",
None,
QtGui.QApplication.UnicodeUTF8)
App.Console.PrintError(msg + '\n'.format(
obj.get('C{}'.format(index - 1))))
else:
msg = QtGui.QApplication.translate(
"ship_console",
"Ambiguous Tank label! ({} instances labeled as"
"'{}' found)",
None,
QtGui.QApplication.UnicodeUTF8)
App.Console.PrintError(msg + '\n'.format(
len(ts),
obj.get('C{}'.format(index - 1))))
continue
t = ts[0]
try:
if t is None or not t.PropertiesList.index("IsTank"):
msg = QtGui.QApplication.translate(
"ship_console",
"Invalid Tank! (the object labeled as"
"'{}' is not a tank)",
None,
QtGui.QApplication.UnicodeUTF8)
App.Console.PrintError(msg + '\n'.format(
len(ws),
obj.get('C{}'.format(index - 1))))
continue
except ValueError:
continue
tanks.append((t, dens, level))
# Let's see if several loading conditions have been selected (and
# prompt a warning)
if self.lc:
@ -281,8 +192,6 @@ class TaskPanel:
break
self.lc = obj
self.ship = ship
self.weights = weights
self.tanks = tanks
if not self.lc:
msg = QtGui.QApplication.translate(
"ship_console",

View File

@ -22,72 +22,79 @@
#***************************************************************************
import math
import FreeCAD as App
import FreeCADGui as Gui
from FreeCAD import Vector, Matrix, Placement
import Part
import Units
import FreeCAD as App
import FreeCADGui as Gui
import Instance as ShipInstance
import WeightInstance
import TankInstance
from shipHydrostatics import Tools as Hydrostatics
G = 9.81
G = Units.parseQuantity("9.81 m/s^2")
MAX_EQUILIBRIUM_ITERS = 10
DENS = 1.025 # [tons/m3], salt water
DENS = Units.parseQuantity("1025 kg/m^3")
TRIM_RELAX_FACTOR = 10.0
def solve(ship, weights, tanks, rolls, var_trim=True):
""" Compute the ship GZ curve.
@param ship Ship instance.
@param weights Considered weights.
@param tanks Considered tanks.
@param rolls List of considered roll angles.
@param var_trim True if the trim angle should be recomputed at each roll
angle, False otherwise.
@return GZ values, drafts and trim angles, for each roll angle (in 3
separated lists)
"""Compute the ship GZ stability curve
Position arguments:
ship -- Ship object
weights -- List of weights to consider
tanks -- List of tanks to consider (each one should be a tuple with the
tank instance, the density of the fluid inside, and the filling level ratio)
rolls -- List of roll angles
Keyword arguments:
var_trim -- True if the equilibrium trim should be computed for each roll
angle, False if null trim angle can be used instead.
Returned value:
List of GZ curve points. Each point contains the GZ stability length, the
equilibrium draft, and the equilibrium trim angle (0 deg if var_trim is
False)
"""
# Get the unloaded weight (ignoring the tanks for the moment).
W = 0.0
COG = Vector()
W = Units.parseQuantity("0 kg")
mom_x = Units.parseQuantity("0 kg*m")
mom_y = Units.parseQuantity("0 kg*m")
mom_z = Units.parseQuantity("0 kg*m")
for w in weights:
W += w.Proxy.getMass(w).getValueAs('kg').Value
W += w.Proxy.getMass(w)
m = w.Proxy.getMoment(w)
COG.x += m[0].getValueAs('kg*m').Value
COG.y += m[1].getValueAs('kg*m').Value
COG.z += m[2].getValueAs('kg*m').Value
COG = COG.multiply(1.0 / W)
mom_x += m[0]
mom_y += m[1]
mom_z += m[2]
COG = Vector(mom_x / W, mom_y / W, mom_z / W)
W = W * G
# Get the tanks weight
TW = 0.0
TW = Units.parseQuantity("0 kg")
VOLS = []
for t in tanks:
# t[0] = tank object
# t[1] = load density
# t[2] = filling level
vol = t[0].Proxy.getVolume(t[0], t[2]).getValueAs('m^3').Value
vol = t[0].Proxy.getVolume(t[0], t[2])
VOLS.append(vol)
TW += vol * t[1]
TW = TW * G
gzs = []
drafts = []
trims = []
points = []
for i,roll in enumerate(rolls):
App.Console.PrintMessage("{0} / {1}\n".format(i + 1, len(rolls)))
gz, draft, trim = solve_point(W, COG, TW, VOLS,
ship, tanks, roll, var_trim)
if gz is None:
point = solve_point(W, COG, TW, VOLS,
ship, tanks, roll, var_trim)
if point is None:
return []
gzs.append(gz)
drafts.append(draft)
trims.append(trim)
points.append(point)
return points
return gzs, drafts, trims
def solve_point(W, COG, TW, VOLS, ship, tanks, roll, var_trim=True):
""" Compute the ship GZ value.
@ -101,56 +108,209 @@ def solve_point(W, COG, TW, VOLS, ship, tanks, roll, var_trim=True):
angle, False otherwise.
@return GZ value, equilibrium draft, and equilibrium trim angle (0 if
variable trim has not been requested)
"""
gz = 0.0
"""
# Look for the equilibrium draft (and eventually the trim angle too)
max_draft = ship.Shape.BoundBox.ZMax / Units.Metre.Value
draft = ship.Draft.getValueAs('m').Value
max_disp = ship.Shape.Volume / Units.Metre.Value**3 * DENS * 1000.0 * G
max_draft = Units.Quantity(ship.Shape.BoundBox.ZMax, Units.Length)
draft = ship.Draft
max_disp = Units.Quantity(ship.Shape.Volume, Units.Volume) * DENS * G
if max_disp < W + TW:
msg = QtGui.QApplication.translate(
"ship_console",
"Too much weight! The ship will never displace water enough",
None,
QtGui.QApplication.UnicodeUTF8)
App.Console.PrintError(msg + ' ({} tons vs. {} tons)\n'.format(
max_disp / 1000.0 / G, (W + TW) / 1000.0 / G))
App.Console.PrintError(msg + ' ({} vs. {})\n'.format(
(max_disp / G).UserString, ((W + TW) / G).UserString))
return None
trim = 0.0
trim = Units.parseQuantity("0 deg")
for i in range(MAX_EQUILIBRIUM_ITERS):
# Get the displacement, and the bouyance application point
disp, B, _ = Hydrostatics.displacement(ship,
draft * Units.Metre,
roll * Units.Degree,
trim * Units.Degree)
disp = disp.getValueAs("kg").Value * G
B.multiply(1.0 / Units.Metre.Value)
draft,
roll,
trim)
disp *= G
# Add the tanks effect on the center of gravity
cog = Vector(COG.x * W, COG.y * W, COG.z * W)
mom_x = Units.Quantity(COG.x, Units.Length) * W
mom_y = Units.Quantity(COG.y, Units.Length) * W
mom_z = Units.Quantity(COG.z, Units.Length) * W
for i,t in enumerate(tanks):
tank_weight = VOLS[i] * t[1] * G
cog += t[0].Proxy.getCoG(t[0], VOLS[i], roll, trim).multiply(
tank_weight / Units.Metre.Value)
cog = cog.multiply(1.0 / (W + TW))
tank_cog = t[0].Proxy.getCoG(t[0], VOLS[i], roll, trim)
mom_x += Units.Quantity(tank_cog.x, Units.Length) * tank_weight
mom_y += Units.Quantity(tank_cog.y, Units.Length) * tank_weight
mom_z += Units.Quantity(tank_cog.z, Units.Length) * tank_weight
cog_x = mom_x / (W + TW)
cog_y = mom_y / (W + TW)
cog_z = mom_z / (W + TW)
# Compute the errors
draft_error = -(disp - W - TW) / max_disp
R = cog - B
draft_error = -((disp - W - TW) / max_disp).Value
R_x = cog_x - Units.Quantity(B.x, Units.Length)
R_y = cog_y - Units.Quantity(B.y, Units.Length)
R_z = cog_z - Units.Quantity(B.z, Units.Length)
if not var_trim:
trim_error = 0.0
else:
trim_error = -TRIM_RELAX_FACTOR * R.x / ship.Length.getValueAs('m').Value
trim_error = -TRIM_RELAX_FACTOR * R_x / ship.Length
# Check if we can tolerate the errors
if abs(draft_error) < 0.01 and abs(trim_error) < 0.05:
if abs(draft_error) < 0.01 and abs(trim_error) < 0.1:
break
# Get the new draft and trim
draft += draft_error * max_draft
trim += trim_error
trim += trim_error * Units.Degree
# GZ should be provided in the Free surface oriented frame of reference
c = math.cos(math.radians(roll))
s = math.sin(math.radians(roll))
return c * R.y - s * R.z, draft, trim
c = math.cos(roll.getValueAs('rad'))
s = math.sin(roll.getValueAs('rad'))
return c * R_y - s * R_z, draft, trim
def gz(lc, rolls, var_trim=True):
"""Compute the ship GZ stability curve
Position arguments:
lc -- Load condition spreadsheet
rolls -- List of roll angles to compute
Keyword arguments:
var_trim -- True if the equilibrium trim should be computed for each roll
angle, False if null trim angle can be used instead.
Returned value:
List of GZ curve points. Each point contains the GZ stability length, the
equilibrium draft, and the equilibrium trim angle (0 deg if var_trim is
False)
"""
# B1 cell must be a ship
# B2 cell must be the loading condition itself
doc = lc.Document
try:
if lc not in doc.getObjectsByLabel(lc.get('B2')):
return[]
ships = doc.getObjectsByLabel(lc.get('B1'))
if len(ships) != 1:
if len(ships) == 0:
msg = QtGui.QApplication.translate(
"ship_console",
"Wrong Ship label! (no instances labeled as"
"'{}' found)",
None,
QtGui.QApplication.UnicodeUTF8)
App.Console.PrintError(msg + '\n'.format(
lc.get('B1')))
else:
msg = QtGui.QApplication.translate(
"ship_console",
"Ambiguous Ship label! ({} instances labeled as"
"'{}' found)",
None,
QtGui.QApplication.UnicodeUTF8)
App.Console.PrintError(msg + '\n'.format(
len(ships),
lc.get('B1')))
return[]
ship = ships[0]
if ship is None or not ship.PropertiesList.index("IsShip"):
return[]
except ValueError:
return[]
# Extract the weights and the tanks
weights = []
index = 6
while True:
try:
ws = doc.getObjectsByLabel(lc.get('A{}'.format(index)))
except ValueError:
break
index += 1
if len(ws) != 1:
if len(ws) == 0:
msg = QtGui.QApplication.translate(
"ship_console",
"Wrong Weight label! (no instances labeled as"
"'{}' found)",
None,
QtGui.QApplication.UnicodeUTF8)
App.Console.PrintError(msg + '\n'.format(
lc.get('A{}'.format(index - 1))))
else:
msg = QtGui.QApplication.translate(
"ship_console",
"Ambiguous Weight label! ({} instances labeled as"
"'{}' found)",
None,
QtGui.QApplication.UnicodeUTF8)
App.Console.PrintError(msg + '\n'.format(
len(ws),
lc.get('A{}'.format(index - 1))))
continue
w = ws[0]
try:
if w is None or not w.PropertiesList.index("IsWeight"):
msg = QtGui.QApplication.translate(
"ship_console",
"Invalid Weight! (the object labeled as"
"'{}' is not a weight)",
None,
QtGui.QApplication.UnicodeUTF8)
App.Console.PrintError(msg + '\n'.format(
len(ws),
lc.get('A{}'.format(index - 1))))
continue
except ValueError:
continue
weights.append(w)
tanks = []
index = 6
while True:
try:
ts = doc.getObjectsByLabel(lc.get('C{}'.format(index)))
dens = float(lc.get('D{}'.format(index)))
level = float(lc.get('E{}'.format(index)))
dens = Units.parseQuantity("{} kg/m^3".format(dens))
except ValueError:
break
index += 1
if len(ts) != 1:
if len(ts) == 0:
msg = QtGui.QApplication.translate(
"ship_console",
"Wrong Tank label! (no instances labeled as"
"'{}' found)",
None,
QtGui.QApplication.UnicodeUTF8)
App.Console.PrintError(msg + '\n'.format(
lc.get('C{}'.format(index - 1))))
else:
msg = QtGui.QApplication.translate(
"ship_console",
"Ambiguous Tank label! ({} instances labeled as"
"'{}' found)",
None,
QtGui.QApplication.UnicodeUTF8)
App.Console.PrintError(msg + '\n'.format(
len(ts),
lc.get('C{}'.format(index - 1))))
continue
t = ts[0]
try:
if t is None or not t.PropertiesList.index("IsTank"):
msg = QtGui.QApplication.translate(
"ship_console",
"Invalid Tank! (the object labeled as"
"'{}' is not a tank)",
None,
QtGui.QApplication.UnicodeUTF8)
App.Console.PrintError(msg + '\n'.format(
len(ws),
lc.get('C{}'.format(index - 1))))
continue
except ValueError:
continue
tanks.append((t, dens, level))
return solve(ship, weights, tanks, rolls, var_trim)

View File

@ -261,12 +261,12 @@ def displacement(ship, draft=None,
B = Part.Point(Vector(cog.x, cog.y, cog.z))
m = Matrix()
m.move(Vector(0.0, 0.0, draft))
m.move(Vector(-draft * math.sin(math.radians(trim)), 0.0, 0.0))
m.rotateY(math.radians(trim))
m.move(Vector(-draft * math.sin(trim.getValueAs("rad")), 0.0, 0.0))
m.rotateY(trim.getValueAs("rad"))
m.move(Vector(0.0,
-draft * math.sin(math.radians(roll)),
-draft * math.sin(roll.getValueAs("rad")),
base_z))
m.rotateX(-math.radians(roll))
m.rotateX(-roll.getValueAs("rad"))
B.transform(m)
try: