Skip to content

Commit 71abd79

Browse files
nikpaunikpau
committed
version update
1 parent c06e5b4 commit 71abd79

File tree

4 files changed

+160
-162
lines changed

4 files changed

+160
-162
lines changed

‎src/mmgdynamics/calibrated_vessels.py‎

Lines changed: 1 addition & 133 deletions
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,3 @@
1-
import math
2-
31
# Dicts of tested and measured vessels.
42
# Data from Yoshimura, Yasuo, and Yumiko Masumoto. "Hydrodynamic database
53
# and manoeuvring prediction method with medium high-speed merchant ships and fishing vessels."
@@ -227,134 +225,4 @@
227225
"I_zG": 2e12, # Moment of inertia of ship around center of gravity (m*(0.25*Lpp)**2) (Point mass Inertia)
228226
"J_z_dash": 0.011, # Added moment of inertia coefficient
229227
"a_H": 0.312 # Rudder force increase factor
230-
}
231-
232-
# In here all the important hydrodynamic derivatives are calculated via empirical
233-
# formulas from Suras and Sakir Bal (2019)
234-
def get_coef_dict(v: dict, rho: float) -> dict:
235-
"""Calculate relevant hydrodynamic derivatives based on a minimal
236-
dict and a given water density
237-
238-
Sources:
239-
Suras and Sakir Bal (2019)
240-
241-
Args:
242-
v (dict): Minimal dict of vessel parameters
243-
rho (float): water density
244-
245-
Raises:
246-
KeyError: If your dict is incomplete, calculations will be aborted
247-
248-
Returns:
249-
dict: Dict with all relevant information to in used as input in the mmg_step() function
250-
"""
251-
252-
# Length, Breadth, draft, Block coef, mass, Rudder area
253-
mininfo = ["B","Lpp","d","C_b","m","A_R"]
254-
255-
if not all (k in v for k in mininfo):
256-
raise KeyError(f"Mandatory keys not found. Need at least ['B','Lpp','d','C_b','A_R']")
257-
258-
# Unpack initial dict
259-
L, B, d, Cb, m = v["Lpp"], v["B"], v["d"], v["C_b"], v["m"]
260-
261-
# X-Coordinate of the center of gravity (m)
262-
if "x_G" not in v:
263-
v["x_G"] = float(0)
264-
265-
# X-Coordinate of the propeller (assumed as -0.5*Lpp)
266-
v["x_P"] = -0.5*L
267-
268-
# Add current water density to dict
269-
v["rho"]= rho
270-
271-
# Masses and Moment of Intertia
272-
nondim_M = 0.5 * v["rho"] * L**2 * d
273-
nondim_N = 0.5 * v["rho"] * L**4 * d
274-
v["m"] = m * v["rho"] # Displacement * water density
275-
v["I_zG"] = v["m"] * (0.25*L)**2
276-
v["m_x_dash"] = 0.05*v["m"] / nondim_M
277-
v["m_y_dash"] = (0.882-0.54*Cb*(1-1.6*d/B)-0.156*(1-0.673*Cb)*L/B+\
278-
0.826*((d*L)/(B**2))*(1-0.678*d/B)-0.638*Cb*((d*L)/(B**2))*(1-0.669*d/B))*v["m"] / nondim_M
279-
v["J_z_dash"] = ((0.01*(33-76.85*Cb*(1-0.784*Cb)+3.43*L/B*(1-0.63*Cb)))**2)*v["m"] / nondim_N
280-
281-
# Hydrodynamic derivatives
282-
# Surge
283-
v["X_vv_dash"] = 1.15*(Cb/(L/B)) - 0.18 # Yoshimura and Masumoto (2012)
284-
v["X_vvvv_dash"]= -6.68*(Cb/(L/B)) + 1.1 # Yoshimura and Masumoto (2012)
285-
v["X_rr_dash"] = (-0.0027+0.0076*Cb*d/B)*L/d # Lee et al. (1998)
286-
v["X_vr_dash"] = v["m_y_dash"] - 1.91*(Cb/(L/B)) + 0.08 # Yoshimura and Masumoto (2012)
287-
288-
# Sway
289-
v["Y_v_dash"] = -(0.5*math.pi*(2*d/L)+1.4*(Cb/(L/B))) # Yoshimura and Masumoto (2012)
290-
v["Y_vvv_dash"] = -0.185*L/B + 0.48 # Yoshimura and Masumoto (2012)
291-
v["Y_r_dash"] = v["m_x_dash"] + 0.5*Cb*B/L # Yoshimura and Masumoto (2012)
292-
v["Y_rrr_dash"] = -0.051
293-
v["Y_vrr_dash"] = -(0.26*(1-Cb)*L/B + 0.11)
294-
v["Y_vvr_dash"] = -0.75 # TODO: Change this
295-
296-
# Yaw
297-
v["N_v_dash"] = -2*d/L # Yoshimura and Masumoto (2012)
298-
v["N_vvv_dash"] = -(-0.69*Cb+0.66)# Yoshimura and Masumoto (2012)
299-
v["N_r_dash"] = -0.54*(2*d/L) + (2*d/L)**2 # Yoshimura and Masumoto (2012)
300-
v["N_rrr_dash"] = ((0.25*Cb)/(L/B))-0.056 # Yoshimura and Masumoto (2012)
301-
v["N_vrr_dash"] = -0.075*(1-Cb)*L/B-0.098 # Yoshimura and Masumoto (2012)
302-
v["N_vvr_dash"] = ((1.55*Cb)/(L/B))-0.76 # Yoshimura and Masumoto (2012)
303-
304-
# Wake coefficient
305-
if "w_P0" not in v:
306-
v["w_P0"] = 0.5*Cb - 0.05
307-
308-
# Thrust deduction factor
309-
if "t_P" not in v:
310-
v["t_P"] = -0.27
311-
312-
# Rudder force incease factor
313-
v["a_H"] = 0.627*Cb-0.153 # Quadvlieg (2013)
314-
315-
# Longituinal acting point of longitudinal force
316-
v["x_H_dash"] = -0.45 # Aoki et. al. (2006)
317-
318-
# Steering resistance deduction factor
319-
v["t_R"] = 0.39 # Yoshimura and Masumoto (2012)
320-
321-
# Espilon [Lee and Shin (1998)]
322-
v["epsilon"] = -2.3281+8.697*Cb-3.78*((2*d)/L)+1.19*Cb**2+292*((2*d)/L)**2-81.51*Cb*((2*d)/L) # Lee and Shin (1998)
323-
#v["epsilon"] = -156.5*(Cb*B/L)**2 + 41.6*(Cb*B/L) - 1.76 # Kijima (1990)
324-
#v["epsilon"] = 2.26*1.82*(1-v["w_P0"]) # Yoshimura and Masumoto (2012)
325-
326-
# Kappa
327-
v["kappa"] = 0.55-0.8*Cb*B/L
328-
329-
# Correction of flow straightening factor to yaw-rate
330-
v["l_R"] = -0.9 # Yoshimura and Masumoto (2012)
331-
332-
# Flow straightening coefficient
333-
v["gamma_R"] = 2.06*Cb*B/L+0.14
334-
335-
# Additional assumptions
336-
v["J_int"] = 0.4
337-
v["k_0"] = 0.4
338-
v["J_slo"] = -0.5
339-
340-
# Rudder aspect ratio
341-
if "f_alpha" not in v:
342-
while True:
343-
asp = input("No rudder aspect ratio found. Please enter manually: ")
344-
try:
345-
asp = float(asp)
346-
break
347-
except ValueError as e:
348-
print("Wrong input type. Only float or int allowed. Err:",e)
349-
350-
asp = float(asp)
351-
falpha =(6.13*asp)/(asp+2.25)
352-
v["f_alpha"] = falpha
353-
354-
print(f"Aspect ration saved. Rudder lift coef calculated to be {falpha}")
355-
356-
# Frictional resistance coefficent [ARBITRARY, no calculations so far]
357-
v["R_0_dash"] = 0.025
358-
359-
return v
360-
228+
}

‎src/mmgdynamics/crossflow.py‎

Lines changed: 6 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -1,3 +1,5 @@
1+
from typing import Dict
2+
13
# Fig 1 in https://doi.org/10.1016/S0141-1187(98)00002-9
24
# Points have been digitalized using an online graph digitizer
35

@@ -36,8 +38,8 @@
3638
0.603975535168196,0.570336391437309]
3739
}
3840

39-
def interpolateY(x: float, vals: dict) -> float:
40-
"""Interpolate y values from a discrete function table
41+
def interpolateY(x: float, vals: Dict[str,list]) -> float:
42+
"""Linear interpolation of y values from a discrete function table
4143
4244
Args:
4345
x (float): input value for which an interpolation shall be returned
@@ -61,7 +63,8 @@ def interpolateY(x: float, vals: dict) -> float:
6163
break
6264

6365
if "func_ind" not in locals():
64-
raise RuntimeError(f"Value out of interpolation range. Choose a value between {min(func_x)} and {max(func_x)}")
66+
raise RuntimeError(f"Value out of interpolation range. "
67+
f"Choose a value between {min(func_x)} and {max(func_x)}")
6568

6669
# Define the two points to interpolate in between
6770
xa, xb = func_x[func_ind+1], func_x[func_ind]

‎src/mmgdynamics/dynamics.py‎

Lines changed: 132 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -1,12 +1,12 @@
11
import math
2-
from typing import Optional
32
import copy
4-
from warnings import warn
53
import numpy as np
6-
from scipy.integrate import solve_ivp
74
import matplotlib.pyplot as plt
85
import mmgdynamics.crossflow as cf
9-
from typing import Any, List
6+
7+
from typing import Optional, Any, List, Dict
8+
from scipy.integrate import solve_ivp
9+
from warnings import warn
1010

1111
"""
1212
Set up the System of ODEs for vessel maneuvering prediction after
@@ -446,6 +446,134 @@ def C_6c(psi: float, p: dict) -> float:
446446

447447
return val1 - val2 - val3
448448

449+
# In here all the important hydrodynamic derivatives are calculated via empirical
450+
# formulas from Suras and Sakir Bal (2019)
451+
def calibrate(v: Dict[str,float], rho: float) -> Dict[str,float]:
452+
"""Calculate relevant hydrodynamic derivatives based on a minimal
453+
dict and a given water density
454+
455+
Sources:
456+
Suras and Sakir Bal (2019)
457+
458+
Args:
459+
v (dict): Minimal dict of vessel parameters
460+
rho (float): water density
461+
462+
Raises:
463+
KeyError: If your dict is incomplete, calculations will be aborted
464+
465+
Returns:
466+
dict: Dict with all relevant information to be used as input in the step() function
467+
"""
468+
469+
# Length, Breadth, draft, Block coef, mass, Rudder area
470+
mininfo = ["B","Lpp","d","C_b","m","A_R"]
471+
472+
if not all (k in v for k in mininfo):
473+
raise KeyError(f"Mandatory keys not found. Need at least ['B','Lpp','d','C_b','A_R']")
474+
475+
# Unpack initial dict
476+
L, B, d, Cb, m = v["Lpp"], v["B"], v["d"], v["C_b"], v["m"]
477+
478+
# X-Coordinate of the center of gravity (m)
479+
if "x_G" not in v:
480+
v["x_G"] = float(0)
481+
482+
# X-Coordinate of the propeller (assumed as -0.5*Lpp)
483+
v["x_P"] = -0.5*L
484+
485+
# Add current water density to dict
486+
v["rho"]= rho
487+
488+
# Masses and Moment of Intertia
489+
nondim_M = 0.5 * v["rho"] * L**2 * d
490+
nondim_N = 0.5 * v["rho"] * L**4 * d
491+
v["m"] = m * v["rho"] # Displacement * water density
492+
v["I_zG"] = v["m"] * (0.25*L)**2
493+
v["m_x_dash"] = 0.05*v["m"] / nondim_M
494+
v["m_y_dash"] = (0.882-0.54*Cb*(1-1.6*d/B)-0.156*(1-0.673*Cb)*L/B+\
495+
0.826*((d*L)/(B**2))*(1-0.678*d/B)-0.638*Cb*((d*L)/(B**2))*(1-0.669*d/B))*v["m"] / nondim_M
496+
v["J_z_dash"] = ((0.01*(33-76.85*Cb*(1-0.784*Cb)+3.43*L/B*(1-0.63*Cb)))**2)*v["m"] / nondim_N
497+
498+
# Hydrodynamic derivatives
499+
# Surge
500+
v["X_vv_dash"] = 1.15*(Cb/(L/B)) - 0.18 # Yoshimura and Masumoto (2012)
501+
v["X_vvvv_dash"]= -6.68*(Cb/(L/B)) + 1.1 # Yoshimura and Masumoto (2012)
502+
v["X_rr_dash"] = (-0.0027+0.0076*Cb*d/B)*L/d # Lee et al. (1998)
503+
v["X_vr_dash"] = v["m_y_dash"] - 1.91*(Cb/(L/B)) + 0.08 # Yoshimura and Masumoto (2012)
504+
505+
# Sway
506+
v["Y_v_dash"] = -(0.5*math.pi*(2*d/L)+1.4*(Cb/(L/B))) # Yoshimura and Masumoto (2012)
507+
v["Y_vvv_dash"] = -0.185*L/B + 0.48 # Yoshimura and Masumoto (2012)
508+
v["Y_r_dash"] = v["m_x_dash"] + 0.5*Cb*B/L # Yoshimura and Masumoto (2012)
509+
v["Y_rrr_dash"] = -0.051
510+
v["Y_vrr_dash"] = -(0.26*(1-Cb)*L/B + 0.11)
511+
v["Y_vvr_dash"] = -0.75 # TODO: Change this
512+
513+
# Yaw
514+
v["N_v_dash"] = -2*d/L # Yoshimura and Masumoto (2012)
515+
v["N_vvv_dash"] = -(-0.69*Cb+0.66)# Yoshimura and Masumoto (2012)
516+
v["N_r_dash"] = -0.54*(2*d/L) + (2*d/L)**2 # Yoshimura and Masumoto (2012)
517+
v["N_rrr_dash"] = ((0.25*Cb)/(L/B))-0.056 # Yoshimura and Masumoto (2012)
518+
v["N_vrr_dash"] = -0.075*(1-Cb)*L/B-0.098 # Yoshimura and Masumoto (2012)
519+
v["N_vvr_dash"] = ((1.55*Cb)/(L/B))-0.76 # Yoshimura and Masumoto (2012)
520+
521+
# Wake coefficient
522+
if "w_P0" not in v:
523+
v["w_P0"] = 0.5*Cb - 0.05
524+
525+
# Thrust deduction factor
526+
if "t_P" not in v:
527+
v["t_P"] = -0.27
528+
529+
# Rudder force incease factor
530+
v["a_H"] = 0.627*Cb-0.153 # Quadvlieg (2013)
531+
532+
# Longituinal acting point of longitudinal force
533+
v["x_H_dash"] = -0.45 # Aoki et. al. (2006)
534+
535+
# Steering resistance deduction factor
536+
v["t_R"] = 0.39 # Yoshimura and Masumoto (2012)
537+
538+
# Espilon [Lee and Shin (1998)]
539+
v["epsilon"] = -2.3281+8.697*Cb-3.78*((2*d)/L)+1.19*Cb**2+292*((2*d)/L)**2-81.51*Cb*((2*d)/L) # Lee and Shin (1998)
540+
#v["epsilon"] = -156.5*(Cb*B/L)**2 + 41.6*(Cb*B/L) - 1.76 # Kijima (1990)
541+
#v["epsilon"] = 2.26*1.82*(1-v["w_P0"]) # Yoshimura and Masumoto (2012)
542+
543+
# Kappa
544+
v["kappa"] = 0.55-0.8*Cb*B/L
545+
546+
# Correction of flow straightening factor to yaw-rate
547+
v["l_R"] = -0.9 # Yoshimura and Masumoto (2012)
548+
549+
# Flow straightening coefficient
550+
v["gamma_R"] = 2.06*Cb*B/L+0.14
551+
552+
# Additional assumptions
553+
v["J_int"] = 0.4
554+
v["k_0"] = 0.4
555+
v["J_slo"] = -0.5
556+
557+
# Rudder aspect ratio
558+
if "f_alpha" not in v:
559+
while True:
560+
asp = input("No rudder aspect ratio found. Please enter manually: ")
561+
try:
562+
asp = float(asp)
563+
break
564+
except ValueError as e:
565+
print("Wrong input type. Only float or int allowed. Err:",e)
566+
567+
asp = float(asp)
568+
falpha =(6.13*asp)/(asp+2.25)
569+
v["f_alpha"] = falpha
570+
571+
print(f"Aspect ration saved. Rudder lift coef calculated to be {falpha}")
572+
573+
# Frictional resistance coefficent [ARBITRARY, no calculations so far]
574+
v["R_0_dash"] = 0.025
575+
576+
return v
449577

450578
def turning_maneuver(ivs: np.ndarray, vessel: dict,
451579
time: int, dir: str = "starboard",

‎src/mmgdynamics/example.py‎

Lines changed: 21 additions & 22 deletions
Original file line numberDiff line numberDiff line change
@@ -1,8 +1,13 @@
1-
import numpy as np
21
import json
3-
from mmgdynamics.dynamics import turning_maneuver, plot_r, plot_trajecory, plot_zigzag, zigzag_maneuver
2+
import numpy as np
43
import mmgdynamics.calibrated_vessels as cvs
54

5+
from mmgdynamics.dynamics import (
6+
turning_maneuver, plot_r,
7+
plot_trajecory, plot_zigzag,
8+
zigzag_maneuver, calibrate
9+
)
10+
611
"""In here you find some basic tests to validate the
712
MMG model.
813
@@ -11,16 +16,11 @@
1116
- ZigZag test
1217
1318
Note:
14-
See the docs for each function for addidtional information
19+
See the docs for each function for additional information
20+
or use help(somefunction)
1521
1622
"""
1723

18-
# Example river dict
19-
some_river = {
20-
"wd_avg": 14.,
21-
"B_K": 200.,
22-
}
23-
2424
# Example minimal dict for a German Großmotorgüterschiff (inland transport vessel)
2525
fully_loaded_GMS = {
2626
"m": 3614.89, # Displacement
@@ -44,24 +44,23 @@
4444
5.0] # Propeller revs [s⁻¹]
4545
)
4646

47-
# Vessel parameter dictionary calculated from minimal dict
48-
#vs = cvs.get_coef_dict(fully_loaded_GMS,1000)
49-
50-
# Uncomment next line to use a pre-calibrated vessel
51-
vs = cvs.GMS1
52-
53-
# Print the vessel dict to output
54-
print(json.dumps(vs,sort_keys=True, indent=4))
47+
# Uncomment to calibrate a vessel from the minimal dict above
48+
#vs = calibrate(fully_loaded_GMS,rho = 1000)
5549

50+
# Use a pre-calibrated vessel
51+
vessel = cvs.GMS1
5652

5753
ZIGZAG: bool = False
5854

5955
iters = 600
60-
s = turning_maneuver(ivs, vs, iters, "starboard")
61-
p = turning_maneuver(ivs, vs, iters, "starboard", water_depth=150.)
62-
plot_trajecory([s, p], vs)
56+
s = turning_maneuver(ivs, vessel, iters, "starboard")
57+
p = turning_maneuver(ivs, vessel, iters, "starboard", water_depth=150.)
58+
plot_trajecory([s, p], vessel)
6359
plot_r(p)
6460

6561
if ZIGZAG:
66-
z, l = zigzag_maneuver(ivs, vs, rise_time=30, max_deg=20)
67-
plot_zigzag(z, l)
62+
z, l = zigzag_maneuver(ivs, vessel, rise_time=30, max_deg=20)
63+
plot_zigzag(z, l)
64+
65+
# Print the vessel dict to output
66+
print(json.dumps(vessel,sort_keys=True, indent=4))

0 commit comments

Comments
 (0)