suzanne.soy/FreeCAD
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Finished the GZ curve implementation

Browse Source
This commit is contained in:
Jose Luis Cercos Pita 2016-01-20 15:13:00 +01:00
parent 6315084f58
commit 9770cca561
4 changed files with 178 additions and 113 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

103
src/Mod/Ship/TankInstance.py
View File

@ -26,7 +26,7 @@ from math import *
from PySide import QtGui, QtCore
import FreeCAD as App
import FreeCADGui as Gui
from FreeCAD import Base, Vector, Placement, Rotation
from FreeCAD import Base, Vector, Matrix, Placement, Rotation
import Part
import Units
from shipUtils import Paths, Math
@ -79,14 +79,23 @@ class Tank:
"""
pass
def getVolume(self, fp, level):
def getVolume(self, fp, level, return_shape=False):
"""Return the fluid volume inside the tank, provided the filling level.
Keyword arguments:
fp -- Part::FeaturePython object affected.
level -- Percentage of filling level (interval [0, 1]).
return_shape -- False if the tool should return the fluid volume value,
True if the tool should return the volume shape.
"""
max(min(level, 1.0), 0.0)
if level <= 0.0:
if return_shape:
return Part.Vertex()
return Units.Quantity(0.0, Units.Volume)
if level >= 1.0:
if return_shape:
return fp.Shape.copy()
return Units.Quantity(fp.Shape.Volume, Units.Volume)
# Build up the cutting box
bbox = fp.Shape.BoundBox
@ -97,7 +106,7 @@ class Tank:
box = App.ActiveDocument.addObject("Part::Box","Box")
length_format = USys.getLengthFormat()
box.Placement = Placement(Vector(bbox.XMin - dx,
bbox.XMin - dy,
bbox.YMin - dy,
bbox.ZMin - dz),
Rotation(App.Vector(0,0,1),0))
box.Length = length_format.format(3.0 * dx)
@ -133,7 +142,91 @@ class Tank:
random_bounds))
App.ActiveDocument.recompute()
return Units.Quantity(common.Shape.Volume, Units.Volume)
if return_shape:
ret_value = common.Shape.copy()
else:
ret_value = Units.Quantity(common.Shape.Volume, Units.Volume)
App.ActiveDocument.removeObject(common.Name)
App.ActiveDocument.removeObject(tank.Name)
App.ActiveDocument.removeObject(box.Name)
App.ActiveDocument.recompute()
return ret_value
def getCoG(self, fp, vol, roll=0.0, trim=0.0):
"""Return the fluid volume center of gravity, provided the volume of
fluid inside the tank.
The returned center of gravity is refered to the untransformed ship.
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).
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:
vol = 0.0
for solid in fp.Shape.Solids:
vol += solid.Volume
sCoG = solid.CenterOfMass
cog.x = cog.x + sCoG.x * solid.Volume
cog.y = cog.y + sCoG.y * solid.Volume
cog.z = cog.z + sCoG.z * solid.Volume
cog.x = cog.x / vol
cog.y = cog.y / vol
cog.z = cog.z / vol
return cog
# Get a first estimation of the level
level = vol / 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
# 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
if abs(error) < 0.01:
break
level += error
# Get the center of gravity
vol = 0.0
cog = Vector()
if len(shape.Solids) > 0:
for solid in shape.Solids:
vol += solid.Volume
sCoG = solid.CenterOfMass
cog.x = cog.x + sCoG.x * solid.Volume
cog.y = cog.y + sCoG.y * solid.Volume
cog.z = cog.z + sCoG.z * solid.Volume
cog.x = cog.x / vol
cog.y = cog.y / vol
cog.z = cog.z / vol
# Untransform the object to retrieve the original position
fp.Placement = current_placement
p = Part.Point(cog)
m = Matrix()
m.rotateY(radians(trim))
m.rotateX(-radians(roll))
p.rotate(Placement(m))
return Vector(p.X, p.Y, p.Z)
class ViewProviderTank:

139
src/Mod/Ship/shipGZ/PlotAux.py
View File

@ -22,36 +22,37 @@
#***************************************************************************
import os
import math
from PySide import QtGui, QtCore
import FreeCAD
import FreeCADGui
from FreeCAD import Base
import Spreadsheet
from shipUtils import Paths
class Plot(object):
def __init__(self, x, y, disp, xcb, ship):
""" Constructor. performs the plot and shows it.
@param x X coordinates.
@param y Transversal computed areas.
@param disp Ship displacement.
@param xcb Bouyancy center length.
@param ship Active ship instance.
"""
self.plot(x, y, disp, xcb, ship)
self.spreadSheet(x, y, ship)
def __init__(self, roll, gz, draft, trim):
""" Plot the GZ curve
def plot(self, x, y, disp, xcb, ship):
""" Perform the areas curve plot.
@param x X coordinates.
@param y Transversal areas.
@param disp Ship displacement.
@param xcb Bouyancy center length.
@param ship Active ship instance.
@return True if error happens.
Position arguments:
roll -- List of roll angles (in degrees).
gz -- List of GZ values (in meters).
draft -- List of equilibrium drafts (in meters).
trim -- List of equilibrium trim angles (in degrees).
"""
self.plot(roll, gz)
self.spreadSheet(roll, gz, draft, trim)
def plot(self, roll, gz):
""" Plot the GZ curve.
Position arguments:
roll -- List of roll angles (in degrees).
gz -- List of GZ values (in meters).
"""
try:
import Plot
plt = Plot.figure('Areas curve')
plt = Plot.figure('GZ')
except ImportError:
msg = QtGui.QApplication.translate(
"ship_console",
@ -60,90 +61,46 @@ class Plot(object):
QtGui.QApplication.UnicodeUTF8)
FreeCAD.Console.PrintWarning(msg + '\n')
return True
# Plot areas curve
areas = Plot.plot(x, y, 'Transversal areas')
areas.line.set_linestyle('-')
areas.line.set_linewidth(2.0)
areas.line.set_color((0.0, 0.0, 0.0))
# Get perpendiculars data
Lpp = ship.Length.getValueAs('m').Value
FPx = 0.5 * Lpp
APx = -0.5 * Lpp
maxArea = max(y)
# Plot perpendiculars
FP = Plot.plot([FPx, FPx], [0.0, maxArea])
FP.line.set_linestyle('-')
FP.line.set_linewidth(1.0)
FP.line.set_color((0.0, 0.0, 0.0))
AP = Plot.plot([APx, APx], [0.0, maxArea])
AP.line.set_linestyle('-')
AP.line.set_linewidth(1.0)
AP.line.set_color((0.0, 0.0, 0.0))
# Add annotations for prependiculars
gz_plot = Plot.plot(roll, gz, 'GZ curve')
gz_plot.line.set_linestyle('-')
gz_plot.line.set_linewidth(1.0)
gz_plot.line.set_color((0.0, 0.0, 0.0))
ax = Plot.axes()
ax.annotate('AP', xy=(APx + 0.01 * Lpp, 0.01 * maxArea), size=15)
ax.annotate('AP', xy=(APx + 0.01 * Lpp, 0.95 * maxArea), size=15)
ax.annotate('FP', xy=(FPx + 0.01 * Lpp, 0.01 * maxArea), size=15)
ax.annotate('FP', xy=(FPx + 0.01 * Lpp, 0.95 * maxArea), size=15)
# Add some additional data
addInfo = ("$XCB = {0} \\; \\mathrm{{m}}$\n"
"$Area_{{max}} = {1} \\; \\mathrm{{m}}^2$\n"
"$\\bigtriangleup = {2} \\; \\mathrm{{tons}}$".format(
xcb,
maxArea,
disp))
ax.text(0.0,
0.01 * maxArea,
addInfo,
verticalalignment='bottom',
horizontalalignment='center',
fontsize=20)
# Write axes titles
Plot.xlabel(r'$x \; \mathrm{m}$')
Plot.ylabel(r'$Area \; \mathrm{m}^2$')
Plot.xlabel(r'$\phi \; [\mathrm{deg}]$')
Plot.ylabel(r'$GZ \; [\mathrm{m}]$')
ax.xaxis.label.set_fontsize(20)
ax.yaxis.label.set_fontsize(20)
# Show grid
Plot.grid(True)
# End
plt.update()
return False
def spreadSheet(self, x, y, ship):
""" Write the output data file.
@param x X coordinates.
@param y Transversal areas.
@param ship Active ship instance.
def spreadSheet(self, roll, gz, draft, trim):
""" Create a Spreadsheet with the results
Position arguments:
roll -- List of roll angles (in degrees).
gz -- List of GZ values (in meters).
draft -- List of equilibrium drafts (in meters).
trim -- List of equilibrium trim angles (in degrees).
"""
s = FreeCAD.activeDocument().addObject('Spreadsheet::Sheet',
'Areas curve')
'GZ')
# Print the header
s.set("A1", "x [m]")
s.set("B1", "area [m^2]")
s.set("C1", "FP x")
s.set("D1", "FP y")
s.set("E1", "AP x")
s.set("F1", "AP y")
s.set("A1", "roll [deg]")
s.set("B1", "GZ [m]")
s.set("C1", "draft [m]")
s.set("D1", "trim [deg]")
# Print the perpendiculars data
Lpp = ship.Length.getValueAs('m').Value
FPx = 0.5 * Lpp
APx = -0.5 * Lpp
maxArea = max(y)
s.set("C2", str(FPx))
s.set("D2", str(0.0))
s.set("C3", str(FPx))
s.set("D3", str(maxArea))
s.set("E2", str(APx))
s.set("F2", str(0.0))
s.set("E3", str(APx))
s.set("F3", str(maxArea))
# Print the data
for i in range(len(x)):
s.set("A{}".format(i + 2), str(x[i]))
s.set("B{}".format(i + 2), str(y[i]))
for i in range(len(roll)):
s.set("A{}".format(i + 2), str(roll[i]))
s.set("B{}".format(i + 2), str(gz[i]))
s.set("C{}".format(i + 2), str(draft[i]))
s.set("D{}".format(i + 2), str(trim[i]))
# Recompute
FreeCAD.activeDocument().recompute()

12
src/Mod/Ship/shipGZ/TaskPanel.py
View File

@ -55,11 +55,13 @@ class TaskPanel:
for i in range(n_points):
rolls.append(roll * i / float(n_points - 1))
gz = Tools.solve(self.ship,
self.weights,
self.tanks,
rolls,
form.var_trim.isChecked())
gzs, drafts, trims = Tools.solve(self.ship,
self.weights,
self.tanks,
rolls,
form.var_trim.isChecked())
PlotAux.Plot(rolls, gzs, drafts, trims)
return True

37
src/Mod/Ship/shipGZ/Tools.py
View File

@ -47,7 +47,8 @@ def solve(ship, weights, tanks, rolls, var_trim=True):
@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 for each roll angle
@return GZ values, drafts and trim angles, for each roll angle (in 3
separated lists)
"""
# Get the unloaded weight (ignoring the tanks for the moment).
W = 0.0
@ -63,33 +64,43 @@ def solve(ship, weights, tanks, rolls, var_trim=True):
# Get the tanks weight
TW = 0.0
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
VOLS.append(vol)
TW += vol * t[1]
TW = TW * G
gzs = []
drafts = []
trims = []
for i,roll in enumerate(rolls):
App.Console.PrintMessage("{0} / {1}\n".format(i + 1, len(rolls)))
gz = solve_point(W, COG, TW, ship, tanks, roll, var_trim)
gz, draft, trim = solve_point(W, COG, TW, VOLS,
ship, tanks, roll, var_trim)
if gz is None:
return []
gzs.append(gz)
drafts.append(draft)
trims.append(trim)
return gzs
return gzs, drafts, trims
def solve_point(W, COG, TW, ship, tanks, roll, var_trim=True):
def solve_point(W, COG, TW, VOLS, ship, tanks, roll, var_trim=True):
""" Compute the ship GZ value.
@param W Empty ship weight.
@param COG Empty ship Center of mass.
@param TW Tanks weights.
@param VOLS List of tank volumes.
@param tanks Considered tanks.
@param roll Roll angle.
@param var_trim True if the trim angle should be recomputed at each roll
angle, False otherwise.
@return GZ value
@return GZ value, equilibrium draft, and equilibrium trim angle (0 if
variable trim has not been requested)
"""
gz = 0.0
@ -107,17 +118,21 @@ def solve_point(W, COG, TW, ship, tanks, roll, var_trim=True):
max_disp / 1000.0 / G, (W + TW) / 1000.0 / G))
return None
trim = 0.0
for i in range(MAX_EQUILIBRIUM_ITERS):
# Get the displacement, and the bouyance application point
disp, B, Cb = Hydrostatics.displacement(ship, draft, roll, trim)
disp *= 1000.0 * G
# Add the tanks effect on the center of gravity
# TODO
# ---
cog = Vector(COG.x * W, COG.y * W, COG.z * 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))
# Compute the errors
draft_error = -(disp - W - TW) / max_disp
R = COG - B
R = cog - B
if not var_trim:
trim_error = 0.0
else:
@ -131,9 +146,7 @@ def solve_point(W, COG, TW, ship, tanks, roll, var_trim=True):
draft += draft_error * max_draft
trim += trim_error
App.Console.PrintMessage("draft and trim: {}, {}\n".format(draft, trim))
# 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
return c * R.y - s * R.z, draft, trim