From 59aeeaea4368418e027076770bb06ac5d510fbb0 Mon Sep 17 00:00:00 2001 From: Thomas Dickson Date: Wed, 17 Jun 2026 20:37:13 +0100 Subject: [PATCH 1/2] feat: scipy 5-DOF optimizer, depowering, and input validation - Replace NLopt with scipy optimizer for Python 3.13 compatibility - Add iterative depowering with heel limit enforcement (phi_max) - Fix leeway clamp and radians conversion in force balance - Fix depowering solver progression past first phi_max solution - Add input validation and error handling across VPP stack - Add comprehensive test suite (hydro, aero, resistance, API, utils) Co-Authored-By: Claude Sonnet 4.6 --- runVPP.py | 4 +- src/AeroMod.py | 93 ++++++--- src/HydroMod.py | 110 ++++++++++- src/UtilsMod.py | 63 ++++-- src/VPPMod.py | 416 ++++++++++++++++++++++++++------------- src/api.py | 74 +++++-- tests/conftest.py | 2 + tests/test_aero.py | 124 ++++++++++++ tests/test_api.py | 192 +++++++++++++++--- tests/test_hydro.py | 129 ++++++++++++ tests/test_resistance.py | 12 ++ tests/test_utils.py | 5 +- tests/test_vpp.py | 245 ++++++++++++++++++++++- 13 files changed, 1243 insertions(+), 226 deletions(-) create mode 100644 tests/conftest.py create mode 100644 tests/test_aero.py create mode 100644 tests/test_hydro.py diff --git a/runVPP.py b/runVPP.py index ca737e6..a06ed97 100755 --- a/runVPP.py +++ b/runVPP.py @@ -1,4 +1,4 @@ -#!/opt/miniconda3/bin/python +#!/usr/bin/env python3 # -*- coding: utf-8 -*- import numpy as np @@ -36,7 +36,7 @@ vpp = VPP(Yacht=YD41) vpp.set_analysis( - tws_range=np.arange(4.0, 22.0, 2.0), twa_range=np.linspace(30.0, 180.0, 31) + tws_range=np.arange(4.0, 22.0, 2.0), twa_range=np.linspace(28.0, 180.0, 39) ) vpp.run(verbose=False) diff --git a/src/AeroMod.py b/src/AeroMod.py index 88101c7..44f0a99 100644 --- a/src/AeroMod.py +++ b/src/AeroMod.py @@ -7,18 +7,46 @@ __version__ = "1.0.1" __email__ = "M.Lauber@soton.ac.uk" -import numpy as np -from scipy.interpolate import interp1d -from scipy.optimize import fsolve -from scipy.optimize import root +import functools + import matplotlib.pyplot as plt +import numpy as np + from src.UtilsMod import build_interp_func +@functools.lru_cache(maxsize=512) +def wind_triangle(tws, twa, vb): + """Analytical wind triangle: TWS/TWA/VB → (AWA, AWS) in degrees. + + Uses the law of cosines on the velocity triangle formed by the true + wind, boat speed, and apparent wind vectors. + """ + twa_rad = np.radians(twa) + aws = np.sqrt(tws**2 + vb**2 + 2 * tws * vb * np.cos(twa_rad)) + if aws < 1e-12: + return twa, 0.0 + cos_awa = np.clip((tws * np.cos(twa_rad) + vb) / aws, -1.0, 1.0) + awa = np.degrees(np.arccos(cos_awa)) + return awa, aws + + class AeroMod(object): def __init__(self, Yacht, rho=1.225, mu=0.0000181): """ - Initializes an Aero Model, given a set of sails + Aerodynamic force model. + + Computes sail drive force, side force, and heeling moment from the + yacht's sail plan using ORC aerodynamic coefficients. + + Parameters + ---------- + Yacht : Yacht + Yacht object containing sail definitions and hull geometry. + rho : float, optional + Air density (kg/m^3). Default is 1.225 (ISA sea level). + mu : float, optional + Dynamic viscosity of air (Pa.s). Default is 1.81e-5. """ # physical params self.rho = rho @@ -74,7 +102,31 @@ def _measure_sails(self): # prototype top function in hydro mod def update(self, vb, phi, tws, twa, flat, RED): """ - Update the aero model for current iter + Update aerodynamic forces for current sailing state. + + Solves the wind triangle, computes sail coefficients, and projects + forces into the boat reference frame. + + Parameters + ---------- + vb : float + Boat speed (m/s). + phi : float + Heel angle (degrees). + tws : float + True wind speed (m/s). + twa : float + True wind angle (degrees). + flat : float + Sail flattening factor (0.62 to 1.0). Reduces lift and drag. + RED : float + Reef/reduction factor. Values > 1 apply jib furling (ftj = RED - 1), + values <= 1 apply mainsail reefing (rfm = RED). + + Returns + ------- + tuple of float + (Fx, Fy, Mx) — drive force (N), side force (N), heeling moment (N.m). """ self.vb = max(0, vb) self.phi = max(0, phi) @@ -112,7 +164,7 @@ def _compute_forces(self): # side-force is horizontal component of Fh self.Fy *= np.cos(np.radians(self.phi)) - + # heeling moment self.Mx = self.Fy * self._vce() @@ -136,12 +188,15 @@ def _get_coeffs(self): self.cl = 0.0 self.cd = 0.0 kpp = 0.0 + sail_cd = {} for sail in self.sails: - - self.cl += sail.cl(self.awa) * sail.area * sail.bk - self.cd += sail.cd(self.awa) * sail.area * sail.bk - kpp += sail.cl(self.awa) ** 2 * sail.area * sail.bk * sail.kp + cl_i = sail.cl(self.awa) + cd_i = sail.cd(self.awa) + sail_cd[id(sail)] = cd_i + self.cl += cl_i * sail.area * sail.bk + self.cd += cd_i * sail.area * sail.bk + kpp += cl_i ** 2 * sail.area * sail.bk * sail.kp self.cl /= self.area self.cd /= self.area @@ -156,7 +211,7 @@ def _get_coeffs(self): for sail in self.sails: if sail.type == "jib": self.fcdj = ( - sail.bk * sail.cd(self.awa) * sail.area / (self.cd * self.area) + sail.bk * sail_cd[id(sail)] * sail.area / (self.cd * self.area) ) # final lift and drag @@ -170,17 +225,7 @@ def _update_windTriangle(self): """ find AWS and AWA for a given TWS, TWA and VB """ - _awa_ = lambda awa: self.vb * np.sin(awa / 180.0 * np.pi) - self.tws * np.sin( - (self.twa - awa) / 180.0 * np.pi - ) - self.awa = fsolve(_awa_, self.twa)[0] - self.aws = np.sqrt( - (self.tws * np.sin(self.twa / 180.0 * np.pi)) ** 2 - + (self.tws * np.cos(self.twa / 180.0 * np.pi) + self.vb) ** 2 - ) - # self.awa = np.arccos((self.tws*np.cos(np.radians(self.twa)) + self.vb) / np.sqrt((self.tws**2) + (self.vb**2) + - # 2*self.tws*self.vb * np.cos(np.radians(self.twa)))) - # self.aws = (self.tws * np.sin(np.radians(self.twa))) / np.sin(self.awa) + self.awa, self.aws = wind_triangle(self.tws, self.twa, self.vb) def _area(self): @@ -194,7 +239,7 @@ def _area(self): def _vce(self): """ - Vectical centre of effort lift/drag weigted + Vertical centre of effort, lift/drag weighted. """ sum = 0.0 for sail in self.sails: diff --git a/src/HydroMod.py b/src/HydroMod.py index 4901793..bcc69c0 100644 --- a/src/HydroMod.py +++ b/src/HydroMod.py @@ -7,14 +7,33 @@ __version__ = "1.0.1" __email__ = "M.Lauber@soton.ac.uk" -import numpy as np -from scipy.interpolate import RegularGridInterpolator import warnings + import matplotlib.pyplot as plt +import numpy as np +from scipy.interpolate import RegularGridInterpolator class HydroMod(object): def __init__(self, Yacht, rho=1025.0, mu=0.00119, g=9.81): + """ + Hydrodynamic resistance and righting moment model. + + Computes total resistance (viscous + residuary + induced), side force + from appendage lift, and righting moment using ORC methods and + interpolated resistance surfaces. + + Parameters + ---------- + Yacht : Yacht + Yacht object containing hull geometry and appendage definitions. + rho : float, optional + Seawater density (kg/m^3). Default is 1025.0. + mu : float, optional + Dynamic viscosity of seawater (Pa.s). Default is 1.19e-3. + g : float, optional + Gravitational acceleration (m/s^2). Default is 9.81. + """ # physical parameters self.rho = rho @@ -111,7 +130,7 @@ def _get_Ri(self): self.Teff = np.hstack((self.Teff, appendage.teff)) self.Ksfj = ( - 0.5 * self.rho * self.vb ** 2 * self.cla * self.leeway / 180.0 * np.pi + 0.5 * self.rho * self.vb ** 2 * self.cla * np.radians(self.leeway) ) self.Ksf = np.sum(self.Ksfj) @@ -122,24 +141,99 @@ def _get_Ri(self): def _cf(self, L): """ - Flate plate turbulent boudnary layer friction coefficient. - Take a length scale, such that it can be used for appendags as well + Flat plate turbulent boundary layer friction coefficient (ITTC 1957) + with ITTC 1978 roughness allowance. Takes a length scale so it can + be used for hull and appendages. """ Re = max( 1e4, self.vb * L / self.nu ) # prevents dividing by zero, lowest for turbulence on plate - return 0.066 * (np.log10(Re) - 2.03) ** (-2) + cf = 0.066 * (np.log10(Re) - 2.03) ** (-2) + ks = self.yacht.roughness + if ks > 0: + cf += (105.0 * (ks / self.l) ** (1.0 / 3.0) - 0.64) * 1e-3 + return cf + + def _added_resistance_waves(self, twa): + """Added resistance in waves (simplified Gerritsma scaling). + + Uses a Bretschneider-type spectrum with hull-form scaling to + estimate the mean added resistance from ocean waves. + + Parameters + ---------- + twa : float + True wind angle (degrees), used as wave encounter angle + unless ``yacht.wave_direction`` overrides it. + + Returns + ------- + float + Added resistance in Newtons. Returns 0 when Hs = 0. + """ + Hs = self.yacht.Hs + Ts = self.yacht.Ts + if Hs <= 0 or Ts <= 0: + return 0.0 - def update(self, vb, phi, leeway): + # wave encounter angle + if self.yacht.wave_direction is not None: + mu = np.radians(twa - self.yacht.wave_direction) + else: + mu = np.radians(twa) + + # heading correction — full drag head-on, near-zero following + cos2_mu = np.cos(mu) ** 2 + + # hull-form coefficient (empirical, typical displacement hull) + C_aw = 6.0 + + # wave encounter frequency + omega_0 = 2.0 * np.pi / Ts + omega_e = abs(omega_0 - omega_0 ** 2 * self.vb * np.cos(mu) / self.g) + + # resonance tuning factor (peak when encounter ≈ hull natural freq) + omega_n = np.sqrt(self.g / self.l) + r = omega_e / omega_n if omega_n > 0 else 0.0 + f_omega = r ** 2 * np.exp(1.0 - r ** 2) if r > 0 else 0.0 + + Raw = C_aw * (Hs ** 2 / self.l) * (self.bwl ** 2 / self.tc) * f_omega * cos2_mu + # convert to Newtons (rho * g scaling) + Raw *= self.rho * self.g / 1000.0 + + return max(0.0, Raw) + + def update(self, vb, phi, leeway, twa=0.0): + """ + Update hydrodynamic forces for current sailing state. + + Parameters + ---------- + vb : float + Boat speed (m/s). + phi : float + Heel angle (degrees). + leeway : float + Leeway angle (degrees). + twa : float, optional + True wind angle (degrees). Used for wave encounter angle. + + Returns + ------- + tuple of float + (Fx, Fy, Mx) — total resistance (N), side force (N), + total righting moment (N.m). + """ self.vb = max(0, vb) self.phi = max(0, phi) - self.leeway = max(0, leeway) + self.leeway = leeway self.lsm, self.lvr, self.btr = self.yacht.measureLSM() self.fn = self.vb / (np.sqrt(self.g * self.lsm)) # resistance self.Fx = self._get_Rr() + self._get_Rv() + self._get_Ri() + self.Fx += self._added_resistance_waves(twa) # keel side force, calculated when _get_Ri() is called self.Fy = self.Ksf * np.cos(self.phi / 180.0 * np.pi) diff --git a/src/UtilsMod.py b/src/UtilsMod.py index 57fe6b5..ac2cd1e 100644 --- a/src/UtilsMod.py +++ b/src/UtilsMod.py @@ -6,6 +6,7 @@ import matplotlib.pyplot as plt import numpy as np from scipy import interpolate +from scipy.interpolate import RegularGridInterpolator KNOTS_TO_MPS = 0.5144 stl = [ @@ -39,11 +40,17 @@ def json_write(data, fname): def build_interp_func(fname, i=1, kind="linear"): """ - build interpolatison function and returns it in a list + build interpolation function and returns it in a list """ a = np.genfromtxt("dat/" + fname + ".dat", delimiter=",", skip_header=1) - # linear for now, this is not good, might need to polish data outside - return interpolate.interp1d(a[0, :], a[i, :], kind=kind, fill_value="extrapolate") + # Filter out NaN values for make_interp_spline compatibility + mask = ~(np.isnan(a[0, :]) | np.isnan(a[i, :])) + x = a[0, mask] + y = a[i, mask] + k = {"linear": 1, "quadratic": 2, "cubic": 3}.get(kind, 1) + spline = interpolate.make_interp_spline(x, y, k=k) + spline.extrapolate = True + return spline def _polar(n) -> plt.Figure: @@ -94,11 +101,24 @@ def _get_vmg(dat, twa_range): return np.array([ix[sup], iy[sdn]], dtype=int), np.array([sup, sdn]) -def _get_cross(dat, n): +def _get_cross(dat, n, pad=2): + """Return [start, end) TWA indices where sail *n* is fastest. + + Parameters + ---------- + dat : ndarray + Shape ``(ntwa, nsails, nvars)``. + n : int + Sail index. + pad : int + Number of extra TWA indices to include on each side for visual + overlap. Use 0 for a tight (no-overlap) range. + """ max_ = np.where(dat[:, n, 0] >= np.max(dat[:, :, 0], axis=1))[0] if len(max_) > 0: idx = np.array( - [max(min(max_) - 2, 0), min(max(max_) + 2, len(dat[:, n, 0]))], dtype=int + [max(min(max_) - pad, 0), min(max(max_) + pad, len(dat[:, n, 0]))], + dtype=int, ) return idx else: @@ -125,8 +145,12 @@ def polar_plot(VPP_list, n, save, fname="Polars.png") -> None: for i in range(len(VPP.tws_range)): vmg, ids = _get_vmg(VPP.store[i, :, :, :], VPP.twa_range) for k in range(VPP.Nsails): - idx = _get_cross(VPP.store[i, :, :, :], k) for j in range(n): + # Speed plot (j=0) uses overlap padding so sail curves + # connect visually; heel/leeway (j>0) use tight range + # to avoid discontinuities at sail crossover points. + pad = 2 if j == 0 else 0 + idx = _get_cross(VPP.store[i, :, :, :], k, pad=pad) lab = "_nolegend_" if k == 0: lab = name + " " + f"{VPP.tws_range[i]/KNOTS_TO_MPS:.1f}" @@ -149,9 +173,23 @@ def polar_plot(VPP_list, n, save, fname="Polars.png") -> None: markersize=4, mfc="None", ) - # add legend only on first axis - ax[0].legend(title=r"TWS (knots)", loc=1, bbox_to_anchor=(1.05, 1.05)) + # TWS legend on every axis + for j in range(n): + ax[j].legend(title=r"TWS (knots)", loc=1, bbox_to_anchor=(1.05, 1.05)) + + # Sail colour legend along the bottom of the figure + from matplotlib.lines import Line2D + VPP = VPP_list[-1] + sail_handles = [ + Line2D([0], [0], color=cols[k % 7], lw=2, label=VPP.sail_name[k]) + for k in range(VPP.Nsails) + ] + fig.legend( + handles=sail_handles, title="Sail set", + loc="lower center", ncol=VPP.Nsails, frameon=True, + ) plt.tight_layout() + fig.subplots_adjust(bottom=0.12) if save: plt.savefig(fname, dpi=96) else: @@ -178,9 +216,12 @@ def sail_chart(VPP, save, fname="SailChart.png"): for j in range(ntwa): if sailset[i, j] == id: sail[i + 1, j + 1] = 1.0 - func = interpolate.interp2d(twas, twss, sail, kind="cubic") - data = func(xnew, ynew) - data = np.where(data > 1.0, 1.0, data) + func = RegularGridInterpolator( + (twss, twas), sail, method="cubic", bounds_error=False, fill_value=0.0 + ) + yy, xx = np.meshgrid(ynew, xnew, indexing="ij") + data = func((yy, xx)) + data = np.clip(data, 0.0, 1.0) ax[0].contour( np.radians(xnew), ynew, data, levels=[0.4], colors=cols[id], alpha=0.8 ) diff --git a/src/VPPMod.py b/src/VPPMod.py index 1232d5b..0abb112 100644 --- a/src/VPPMod.py +++ b/src/VPPMod.py @@ -9,10 +9,10 @@ import logging import warnings +from concurrent.futures import ProcessPoolExecutor -import nlopt import numpy as np -from scipy.optimize import root +from scipy.optimize import least_squares, minimize, root from tqdm import trange from src.AeroMod import AeroMod @@ -24,6 +24,91 @@ debug_mode = logging.getLogger().getEffectiveLevel() == logging.DEBUG +def _solve_point(yacht, tws, twa, sail_index, phi_max, lim_up, lim_dn): + """Solve one (tws, twa, sail) grid point independently. + + Creates its own AeroMod/HydroMod so there is no shared mutable state. + This is a module-level function so it can be pickled by ProcessPoolExecutor. + + Returns + ------- + tuple + (i, j, n, result_array) where result_array is [vb_kts, phi, leeway, flat, red], + or None if the point was skipped. + """ + i, j, n = sail_index + + aero = AeroMod(yacht) + hydro = HydroMod(yacht) + aero.sails[1] = yacht.sails[n + 1] + aero.up = aero.sails[1].up + + # Skip invalid sail/TWA combinations + Nsails = len(yacht.sails) - 1 + dn_limit = 135.0 if Nsails != 1 else 200.0 + if aero.up and twa >= dn_limit: + return None + if not aero.up and twa <= lim_up: + return None + + # Initial guesses + vb0 = 0.8 * tws + phi0 = 0.0 + leeway0 = 100.0 / twa if (twa > 1.0 and 100.0 / twa < 2 * tws) else 2 * tws + + def resid(x0, twa_, tws_, flat=1.0, red=2.0): + vb_, phi_, leeway_ = x0 + Fxh, Fyh, Mxh = hydro.update(vb_, phi_, leeway_, twa_) + Fxa, Fya, Mxa = aero.update(vb_, phi_, tws_, twa_, flat, red) + return [(Fxh - Fxa) ** 2, (Mxh - Mxa) ** 2, (Fyh - Fya) ** 2] + + flat = 1.0 + red = 2.0 + + sol = root(resid, [vb0, phi0, leeway0], + args=(twa, tws, flat, red), method="lm") + vb, phi, leeway = sol.x + + if phi <= phi_max: + res = np.array([vb, phi, leeway, flat, red]) + res[0] /= KNOTS_TO_MPS + return (i, j, n, res) + + # Depowering: flatten then reef + lo = [0, 0, -2] + hi = [np.inf, phi_max, 6] + margin = phi_max - 1.0 + + def _clamp(vb_, phi_, leeway_): + return [max(vb_, 0), min(max(phi_, 0), phi_max), + min(max(leeway_, -2), 6)] + + for flat in np.arange(0.98, 0.60, -0.02): + sol = least_squares( + resid, _clamp(vb, phi_max, leeway), + args=(twa, tws, flat, red), bounds=(lo, hi), + ) + vb, phi, leeway = sol.x + if phi <= margin: + res = np.array([vb, phi, leeway, flat, red]) + res[0] /= KNOTS_TO_MPS + return (i, j, n, res) + + flat = 0.62 + for red in np.arange(1.9, 0.45, -0.1): + sol = least_squares( + resid, _clamp(vb, phi_max, leeway), + args=(twa, tws, flat, red), bounds=(lo, hi), + ) + vb, phi, leeway = sol.x + if phi <= margin: + break + + res = np.array([vb, phi, leeway, flat, red]) + res[0] /= KNOTS_TO_MPS + return (i, j, n, res) + + class VPP(object): """A VPP Class that run an analysis on a given Yacht.""" @@ -49,8 +134,9 @@ def __init__(self, Yacht): warnings.filterwarnings( "ignore", "The iteration is not making good progress" ) + self.upToDate = False - def set_analysis(self, tws_range, twa_range): + def set_analysis(self, tws_range, twa_range, phi_max=35.0): """ Sets the analysis range. Parameters @@ -61,19 +147,28 @@ def set_analysis(self, tws_range, twa_range): A numpy.array with the different TWA to run the analysis at. """ - if tws_range.max() <= 35.0 and tws_range.min() >= 2.0: - logging.debug("Analysis set for TWS: ", tws_range) - self.tws_range = tws_range * KNOTS_TO_MPS - else: - logging.debug("Analysis only valid for TWS range : 2. < TWS < 35. knots.") - - if twa_range.max() <= 180.0 and twa_range.min() >= 0.0: - self.twa_range = twa_range - logging.debug("Analysis set for TWA: ", self.twa_range) - else: - logging.debug( - "Analysis only valid for TWA range : 0. < TWA < 180. degrees." + self.phi_max = phi_max + + if tws_range.size == 0: + raise ValueError("TWS range is empty. Ensure min and max TWS are not equal.") + if twa_range.size == 0: + raise ValueError("TWA range is empty. Ensure min and max TWA are not equal.") + + if tws_range.min() < 2.0 or tws_range.max() > 35.0: + raise ValueError( + f"TWS range [{tws_range.min():.1f}, {tws_range.max():.1f}] " + f"is outside valid bounds [2.0, 35.0] knots." ) + self.tws_range = tws_range * KNOTS_TO_MPS + logger.debug("Analysis set for TWS: %s", tws_range) + + if twa_range.min() < 0.0 or twa_range.max() > 180.0: + raise ValueError( + f"TWA range [{twa_range.min():.1f}, {twa_range.max():.1f}] " + f"is outside valid bounds [0.0, 180.0] degrees." + ) + self.twa_range = twa_range + logger.debug("Analysis set for TWA: %s", self.twa_range) # prepare storage array self.Nsails = len(self.yacht.sails) - 1 # main not counted @@ -95,51 +190,33 @@ def set_analysis(self, tws_range, twa_range): # flag for later self.upToDate = True - def Vb(self, x, grad): - # this should not be used - if grad.size > 0: - grad = 0.0 - return self.vb0 - - def SumForce(self, res, x, grad, twa_, tws_): - # this should not be used - if grad.size > 0: - grad[:, :] = 0.0 - - vb0 = x[0] - phi0 = x[1] - leeway = x[2] - flat = x[3] - red = x[4] - - Fxh, Fyh, Mxh = self.hydro.update(vb0, phi0, leeway) - Fxa, Fya, Mxa = self.aero.update(vb0, phi0, tws_, twa_, flat, red) - - res[0] = Fxh - Fxa - res[1] = Mxh - Mxa - res[2] = Fyh - Fya - - return None + def run(self, verbose=False, method="iterative", progress_callback=None): + """ + Run the analysis for the given analysis range. + Parameters + ---------- + verbose + A logical, if True, prints results of equilibrium at each TWA/TWS. + method + Solver method: "iterative" (3-DOF sequential with depowering loop), + "parallel" (same solver, multiprocessing across grid points), or + "5dof" (scipy SLSQP 5-DOF constrained optimizer). + progress_callback + Optional callable(tws_index, tws_kts, n_tws) invoked after each + TWS wind speed is completed. + """ - def run_NLopt(self, verbose=False): - logging.info("Optimisation start") + if method == "5dof": + return self._run_5dof(verbose) + elif method == "parallel": + return self._run_parallel() + elif method != "iterative": + raise ValueError(f"Unknown method '{method}'. Use 'iterative', 'parallel', or '5dof'.") if not self.upToDate: - raise "VPP run stop: no analysis set!" - - # gradient-free optimization because the gradient of our - # objective function cannot be evaluated - opt = nlopt.opt(nlopt.LN_COBYLA, 5) + raise RuntimeError("VPP run stop: no analysis set!") - # out three parameters are x = [v_b, hell, leeway, flat, red] - opt.set_lower_bounds([0.0, 0.0, 0.0, 0.0, 0.0]) - opt.set_upper_bounds([float("inf"), self.phi_max, 6.0, 1.0, 2.0]) - - # the function we want to maximise - opt.set_max_objective(self.Vb) - - # set solver tolerance - opt.set_xtol_rel(1e-6) + n_tws = len(self.tws_range) for i, tws in enumerate(self.tws_range): logging.debug("Sailing in TWS : %.1f" % (tws / KNOTS_TO_MPS)) @@ -164,6 +241,8 @@ def run_NLopt(self, verbose=False): if (twa > 1.0 and 100.0 / twa < 2 * tws) else 2 * tws ) + self.flat = 1.0 + self.red = 2.0 # don't do low twa with downwind sails if (self.aero.up == True) and (twa >= self.lim_dn): @@ -171,110 +250,181 @@ def run_NLopt(self, verbose=False): if (self.aero.up == False) and (twa <= self.lim_up): continue - # vector-valued constraint - constrain = lambda res, x, grad: self.SumForce( - res, x, grad, twa_=twa, tws_=tws - ) - opt.add_equality_mconstraint(constrain, np.full(5, 1e-8)) + vb, phi, leeway, flat, red = self._depower_solve(twa, tws) + self.vb0, self.phi0, self.leeway0 = vb, phi, leeway - x0 = np.array([self.vb0, self.phi0, self.leeway0, 1.0, 2.0]) - res = opt.optimize(x0) + logging.debug( + "Equilibrium residuals (Fx, Fy, Mx): ", + self.resid([vb, phi, leeway], twa, tws, flat, red), + ) # store data for later - self.store[i, j, n, :] = res[:] * np.array( - [1.0 / KNOTS_TO_MPS, 1, 1, 1, 1] - ) + res = np.array([vb, phi, leeway, flat, red]) + self.store[i, j, n, :] = res * np.array([1.0 / KNOTS_TO_MPS, 1, 1, 1, 1]) - # clean up - opt.remove_equality_constraints() + if progress_callback: + progress_callback(i, tws / KNOTS_TO_MPS, n_tws) logging.info("Optimization successful.") - def run(self, verbose=False): - """ - Run the analysis for the given analysis range. - Parameters - ---------- - verbose - A logical, if True, prints results of equilibrium at each TWA/TWS. + def _run_parallel(self): + """Run iterative solver in parallel across all grid points.""" + if not self.upToDate: + raise RuntimeError("VPP run stop: no analysis set!") + + # Build work items: (yacht, tws, twa, (i,j,n), phi_max, lim_up, lim_dn) + work = [] + for i, tws in enumerate(self.tws_range): + for n in range(self.Nsails): + for j, twa in enumerate(self.twa_range): + work.append(( + self.yacht, tws, twa, (i, j, n), + self.phi_max, self.lim_up, self.lim_dn + )) + + with ProcessPoolExecutor() as pool: + results = pool.map(_solve_point, *zip(*work)) + + for result in results: + if result is not None: + i, j, n, res = result + self.store[i, j, n, :] = res + + logging.info("Parallel optimization successful.") + + def _depower_solve(self, twa, tws): + """Solve 3-DOF equilibrium, depowering iteratively if heel exceeds phi_max. + + Depowering follows real sailing practice: + 1. Flatten sails (flat: 1.0 -> 0.62) + 2. Reef main / furl jib (RED: 2.0 -> 0.5) + + The bounded solver pins phi at phi_max, so we check for genuine + margin (1 degree) before accepting a solution — otherwise the boat + is still overpowered and needs more depowering. """ + flat = 1.0 + red = 2.0 + + sol = root(self.resid, [self.vb0, self.phi0, self.leeway0], + args=(twa, tws, flat, red), method="lm") + vb, phi, leeway = sol.x + + if phi <= self.phi_max: + return vb, phi, leeway, flat, red + + # Heel exceeds limit — use bounded solver to enforce phi <= phi_max + lo = [0, 0, -2] + hi = [np.inf, self.phi_max, 6] + margin = self.phi_max - 1.0 + + def _clamp_guess(vb_, phi_, leeway_): + return [max(vb_, 0), min(max(phi_, 0), self.phi_max), + min(max(leeway_, -2), 6)] + + # Stage 1: Flatten sails (1.0 -> 0.62) + for flat in np.arange(0.98, 0.60, -0.02): + sol = least_squares( + self.resid, _clamp_guess(vb, self.phi_max, leeway), + args=(twa, tws, flat, red), bounds=(lo, hi), + ) + vb, phi, leeway = sol.x + if phi <= margin: + return vb, phi, leeway, flat, red + + # Stage 2: Reef main / furl jib (RED: 2.0 -> 0.5) + flat = 0.62 + for red in np.arange(1.9, 0.45, -0.1): + sol = least_squares( + self.resid, _clamp_guess(vb, self.phi_max, leeway), + args=(twa, tws, flat, red), bounds=(lo, hi), + ) + vb, phi, leeway = sol.x + if phi <= margin: + return vb, phi, leeway, flat, red + # If still over, return best we got + return vb, phi, leeway, flat, red + + def _run_5dof(self, verbose=False): + """Run 5-DOF constrained optimization using scipy SLSQP. + + Simultaneously optimizes [vb, phi, leeway, flat, red] to maximize + boat speed subject to force/moment equilibrium constraints. + """ if not self.upToDate: - raise "VPP run stop: no analysis set!" + raise RuntimeError("VPP run stop: no analysis set!") for i, tws in enumerate(self.tws_range): logging.debug("Sailing in TWS : %.1f" % (tws / KNOTS_TO_MPS)) for n in range(self.Nsails): self.aero.sails[1] = self.yacht.sails[n + 1] - - logging.debug( - "Sail Config : ", - self.aero.sails[0].name + " + " + self.aero.sails[1].name, - ) - self.aero.up = self.aero.sails[1].up for j in trange(len(self.twa_range), disable=not debug_mode): twa = self.twa_range[j] - self.vb0 = 0.8 * tws - self.phi0 = 0 - self.leeway0 = ( - 100.0 / twa - if (twa > 1.0 and 100.0 / twa < 2 * tws) - else 2 * tws - ) - self.flat = 1.0 - self.red = 2.0 - # don't do low twa with downwind sails if (self.aero.up == True) and (twa >= self.lim_dn): continue if (self.aero.up == False) and (twa <= self.lim_up): continue - sol = root( - self.resid, - [self.vb0, self.phi0, self.leeway0], - args=(twa, tws), - method="lm", + vb_guess = 0.8 * tws + leeway_guess = ( + 100.0 / twa + if (twa > 1.0 and 100.0 / twa < 2 * tws) + else 2 * tws ) - self.vb0, self.phi0, self.leeway0 = res = sol.x - - if verbose and not sol.success: - logger.debug(sol.message) - - # # contraints - # con1 = {'type': 'eq', 'fun': self.Fx, 'args': (twa, tws)} - # con2 = {'type': 'eq', 'fun': self.Fy, 'args': (twa, tws)} - # con3 = {'type': 'eq', 'fun': self.Mx, 'args': (twa, tws)} - # con = (con1, con2, con3) - - # # initial guess at this twa/tws - # x0 = [self.vb0, self.phi0, self.leeway0, self.flat, self.red] - - # # minimize - # sol = minimize(self.objective, args=(twa, tws), x0=x0, method='SLSQP', - # constraints=con, bounds=self.bnds, tol=1e-2, - # options={"maxiter": 100, "disp": verbose}) - - # # get result - # self.vb0, self.phi0, self.leeway0, self.flat, self.red = res = sol.x - - logging.debug( - "Equilibrium residuals (Fx, Fy, Mx): ", - self.resid(sol.x, twa, tws), + x0 = np.array([vb_guess, 0.0, leeway_guess, 1.0, 2.0]) + + def _forces(x): + vb, phi, leeway, flat, red = x + Fxh, Fyh, Mxh = self.hydro.update(vb, phi, leeway, twa) + Fxa, Fya, Mxa = self.aero.update(vb, phi, tws, twa, flat, red) + return Fxh, Fyh, Mxh, Fxa, Fya, Mxa + + constraints = [ + {"type": "eq", "fun": lambda x: _forces(x)[0] - _forces(x)[3]}, + {"type": "eq", "fun": lambda x: _forces(x)[2] - _forces(x)[5]}, + {"type": "eq", "fun": lambda x: _forces(x)[1] - _forces(x)[4]}, + ] + + # Cache forces to avoid redundant evaluations + _cache = {} + + def _cached_forces(x): + key = tuple(x) + if key not in _cache: + _cache[key] = _forces(x) + return _cache[key] + + constraints = [ + {"type": "eq", "fun": lambda x: _cached_forces(x)[0] - _cached_forces(x)[3]}, + {"type": "eq", "fun": lambda x: _cached_forces(x)[2] - _cached_forces(x)[5]}, + {"type": "eq", "fun": lambda x: _cached_forces(x)[1] - _cached_forces(x)[4]}, + ] + + result = minimize( + lambda x: -x[0], + x0, + method="SLSQP", + bounds=self.bnds, + constraints=constraints, + options={"maxiter": 200, "ftol": 1e-8}, ) - # store data for later - self.store[i, j, n, : len(res)] = ( - res[:] * np.array([1.0 / KNOTS_TO_MPS, 1, 1, 1, 1])[: len(res)] + _cache.clear() + + res = result.x + self.store[i, j, n, :] = res * np.array( + [1.0 / KNOTS_TO_MPS, 1, 1, 1, 1] ) - logging.info("Optimization successful.") + logging.info("5-DOF optimization successful.") - def resid(self, x0, twa, tws): + def resid(self, x0, twa, tws, flat=1.0, red=2.0): """ Computes the residuals of the force/moment equilibrium at the given state. Parameters @@ -285,6 +435,10 @@ def resid(self, x0, twa, tws): A float of the TWA at which to compute the residuals. tws A float of the TWs at which to compute the residuals. + flat + Sail flattening factor (0.62 to 1.0). Default 1.0 (no flattening). + red + Reef/reduction factor. Default 2.0 (full sails, no reef/furl). Returns ------- Numpy.Array @@ -292,13 +446,11 @@ def resid(self, x0, twa, tws): """ vb0 = x0[0] - phi0 = x0[1] # min(x0[1], self.phi_max) + phi0 = x0[1] leeway = x0[2] - flat = 1.0 # x0[3] - red = 1.0 # x0[4] - Fxh, Fyh, Mxh = self.hydro.update(vb0, phi0, leeway) - Fxa, Fya, Mxa = self.aero.update(vb0, phi0, tws, twa, flat, 2.0) + Fxh, Fyh, Mxh = self.hydro.update(vb0, phi0, leeway, twa) + Fxa, Fya, Mxa = self.aero.update(vb0, phi0, tws, twa, flat, red) return [(Fxh - Fxa) ** 2, (Mxh - Mxa) ** 2, (Fyh - Fya) ** 2] diff --git a/src/api.py b/src/api.py index 14ad2f6..56a3867 100644 --- a/src/api.py +++ b/src/api.py @@ -13,7 +13,7 @@ sys.path.append(os.path.realpath(".")) from src.SailMod import Jib, Kite, Main from src.VPPMod import VPP -from src.YachtMod import Keel, Rudder, Yacht +from src.YachtMod import Keel, Rudder, ShortKeel, Yacht app = Flask(__name__) @@ -27,17 +27,35 @@ def ping(): def data_to_vpp(data: Dict[str, Any]) -> VPP: - - keel = Keel( - Cu=float(data["keel"]["Cu"]), - Cl=float(data["keel"]["Cl"]), - Span=float(data["keel"]["Span"]) - ) + + keel_data = data["keel"] + keel_type = keel_data.get("type", "fin") + # Also detect from keys if type is missing or inconsistent + if keel_type == "short" or "Length" in keel_data: + keel = ShortKeel( + Length=float(keel_data["Length"]), + Depth=float(keel_data["Depth"]), + Tc_ratio=float(keel_data.get("Tc_ratio", 0.15)), + ) + else: + keel = Keel( + Cu=float(keel_data["Cu"]), + Cl=float(keel_data["Cl"]), + Span=float(keel_data["Span"]), + ) rudder = Rudder( - Cu=float(data["rudder"]["Cu"]), - Cl=float(data["rudder"]["Cu"]), + Cu=float(data["rudder"]["Cu"]), + Cl=float(data["rudder"]["Cu"]), Span=float(data["rudder"]["Span"]) ) + # Environment parameters + roughness = float(data.get("roughness", 150e-6)) + Hs = float(data.get("Hs", 0.0)) + Ts = float(data.get("Ts", 0.0)) + wave_direction = data.get("wave_direction") + if wave_direction is not None: + wave_direction = float(wave_direction) + yacht = Yacht( Name=data["yacht"]["Name"], Lwl=float(data["yacht"]["Lwl"]), @@ -52,6 +70,10 @@ def data_to_vpp(data: Dict[str, Any]) -> VPP: Fa=float(data["yacht"]["Fa"]), Boa=float(data["yacht"]["Boa"]), Loa=float(data["yacht"]["Loa"]), + roughness=roughness, + Hs=Hs, + Ts=Ts, + wave_direction=wave_direction, App=[keel, rudder], Sails=[ Main( @@ -60,6 +82,10 @@ def data_to_vpp(data: Dict[str, Any]) -> VPP: E=float(data["main"]["E"]), Roach=float(data["main"]["Roach"]), BAD=float(data["main"]["BAD"]), + data_source=data.get("data_source", "orc"), + cl_data=data["main"].get("cl_data"), + cd_data=data["main"].get("cd_data"), + sail_type=data["main"].get("sail_type"), ), Jib( name=data["jib"]["Name"], @@ -67,32 +93,48 @@ def data_to_vpp(data: Dict[str, Any]) -> VPP: J=float(data["jib"]["J"]), LPG=float(data["jib"]["LPG"]), HBI=float(data["jib"]["HBI"]), + data_source=data.get("data_source", "orc"), + cl_data=data["jib"].get("cl_data"), + cd_data=data["jib"].get("cd_data"), + sail_type=data["jib"].get("sail_type"), ), Kite( name=data["kite"]["Name"], area=float(data["kite"]["area"]), vce=float(data["kite"]["vce"]), + data_source=data.get("data_source", "orc"), + cl_data=data["kite"].get("cl_data"), + cd_data=data["kite"].get("cd_data"), + sail_type=data["kite"].get("sail_type"), ), ], ) - + vpp = VPP(Yacht=yacht) vpp.set_analysis( tws_range=np.array(data["tws_range"]), twa_range=np.array(data["twa_range"]), ) - return vpp - + return vpp, data.get("method", "iterative") + @app.route("/api/vpp/", methods=["POST"]) def makevppresults(): data = request.get_json() + if data is None: + return jsonify({"error": "Request body must be valid JSON."}), 400 - # TODO: Support multiple implementations of different sails: require API design - # TODO: Error handling incorrect ranges + try: + vpp, method = data_to_vpp(data) + except (KeyError, TypeError, ValueError) as e: + logging.warning("Invalid VPP input: %s", e) + return jsonify({"error": f"Invalid input: {e}"}), 400 - vpp = data_to_vpp(data) - vpp.run(verbose=True) + try: + vpp.run(verbose=True, method=method) + except Exception as e: + logging.exception("VPP simulation failed") + return jsonify({"error": f"Simulation failed: {e}"}), 500 return jsonify(vpp.results()) diff --git a/tests/conftest.py b/tests/conftest.py new file mode 100644 index 0000000..da2faa6 --- /dev/null +++ b/tests/conftest.py @@ -0,0 +1,2 @@ +import matplotlib +matplotlib.use("Agg") diff --git a/tests/test_aero.py b/tests/test_aero.py new file mode 100644 index 0000000..4003109 --- /dev/null +++ b/tests/test_aero.py @@ -0,0 +1,124 @@ +"""Tests for AeroMod — wind triangle and coefficient caching.""" + +import numpy as np +import pytest + + +class TestWindTriangleCache: + """Verify that wind_triangle uses lru_cache.""" + + def test_wind_triangle_is_cached(self): + from src.AeroMod import wind_triangle + + # lru_cache-decorated functions have cache_info() + assert hasattr(wind_triangle, "cache_info"), "wind_triangle should be lru_cache-decorated" + wind_triangle.cache_clear() + wind_triangle(10.0, 90.0, 5.0) + wind_triangle(10.0, 90.0, 5.0) + info = wind_triangle.cache_info() + assert info.hits == 1 + assert info.misses == 1 + + def test_cache_returns_same_result(self): + from src.AeroMod import wind_triangle + + wind_triangle.cache_clear() + r1 = wind_triangle(10.0, 45.0, 6.0) + r2 = wind_triangle(10.0, 45.0, 6.0) + assert r1 == r2 + + +class TestSailCoefficientCalls: + """Verify that _get_coeffs doesn't call sail.cl() redundantly.""" + + def test_get_coeffs_calls_cl_once_per_sail(self): + """Each sail's cl(awa) should be called at most once per _get_coeffs().""" + from unittest.mock import MagicMock, patch + from src.AeroMod import AeroMod + + # Create a mock yacht with minimal structure + mock_sail = MagicMock() + mock_sail.cl = MagicMock(return_value=1.0) + mock_sail.cd = MagicMock(return_value=0.1) + mock_sail.area = 30.0 + mock_sail.bk = 1.0 + mock_sail.kp = 1.0 + mock_sail.type = "main" + + aero = object.__new__(AeroMod) + aero.sails = [mock_sail] + aero.area = 30.0 + aero.awa = 30.0 + aero.flat = 1.0 + aero.fcdmult = lambda flat: 1.0 + + # Provide _heff so CE computation works + aero._heff = lambda awa: 10.0 + + aero._get_coeffs() + + # cl should be called exactly once per sail (not twice as before) + assert mock_sail.cl.call_count == 1, ( + f"sail.cl() called {mock_sail.cl.call_count} times, expected 1" + ) + assert mock_sail.cd.call_count == 1, ( + f"sail.cd() called {mock_sail.cd.call_count} times, expected 1" + ) + + +class TestWindTriangle: + """Verify analytical wind triangle: given TWS, TWA, VB → AWS, AWA.""" + + @pytest.mark.parametrize( + "tws, twa, vb, expected_aws, expected_awa", + [ + # Beam reach: TWA=90°, VB=5, TWS=10 + # AWS = sqrt(10² + 5² + 2·10·5·cos(90°)) = sqrt(125) ≈ 11.18 + # AWA = arccos((10·cos(90°) + 5) / 11.18) = arccos(5/11.18) ≈ 63.43° + (10.0, 90.0, 5.0, np.sqrt(125), np.degrees(np.arccos(5.0 / np.sqrt(125)))), + # Close-hauled: TWA=45°, VB=6, TWS=12 + ( + 12.0, + 45.0, + 6.0, + np.sqrt(12**2 + 6**2 + 2 * 12 * 6 * np.cos(np.radians(45))), + np.degrees( + np.arccos( + (12 * np.cos(np.radians(45)) + 6) + / np.sqrt(12**2 + 6**2 + 2 * 12 * 6 * np.cos(np.radians(45))) + ) + ), + ), + # Dead downwind: TWA=180°, VB=4, TWS=10 + # AWS = sqrt(100 + 16 - 80) = sqrt(36) = 6 + # AWA = arccos((10·cos(180°) + 4) / 6) = arccos(-6/6) = 180° + (10.0, 180.0, 4.0, 6.0, 180.0), + # No boat speed: VB=0 → AWA=TWA, AWS=TWS + (10.0, 90.0, 0.0, 10.0, 90.0), + # Boat speed equals TWS, broad reach + (8.0, 120.0, 8.0, + np.sqrt(8**2 + 8**2 + 2 * 8 * 8 * np.cos(np.radians(120))), + np.degrees( + np.arccos( + (8 * np.cos(np.radians(120)) + 8) + / np.sqrt(8**2 + 8**2 + 2 * 8 * 8 * np.cos(np.radians(120))) + ) + )), + ], + ids=["beam_reach", "close_hauled", "dead_downwind", "no_boat_speed", "broad_reach"], + ) + def test_wind_triangle_analytical(self, tws, twa, vb, expected_aws, expected_awa): + """AWA and AWS must match law-of-cosines closed-form solution.""" + from src.AeroMod import wind_triangle + + awa, aws = wind_triangle(tws, twa, vb) + assert aws == pytest.approx(expected_aws, abs=1e-10) + assert awa == pytest.approx(expected_awa, abs=1e-10) + + def test_no_fsolve_import_used(self): + """The analytical path must not import or call fsolve.""" + from src.AeroMod import wind_triangle + import inspect + + source = inspect.getsource(wind_triangle) + assert "fsolve" not in source diff --git a/tests/test_api.py b/tests/test_api.py index 53aab60..f1e068a 100644 --- a/tests/test_api.py +++ b/tests/test_api.py @@ -4,6 +4,9 @@ from src.api import app +HEADERS = {"content-type": "application/json", "Accept-Charset": "UTF-8"} + + def test_ping_route(): client = app.test_client() response = client.get("/ping") @@ -11,29 +14,27 @@ def test_ping_route(): assert response.data.decode("utf-8") == "Pong! The server is up and running." -def make_yd41(): - yacht = dict( - { - "Name": "YD41", - "Lwl": 11.90, - "Vol": 6.05, - "Bwl": 3.18, - "Tc": 0.4, - "WSA": 28.20, - "Tmax": 2.30, - "Amax": 1.051, - "Mass": 6500, - "Ff": 1.5, - "Fa": 1.5, - "Boa": 4.2, - "Loa": 12.5, - } - ) - keel = dict({"Cu": 1.00, "Cl": 0.78, "Span": 1.90}) - rudder = dict({"Cu": 0.48, "Cl": 0.22, "Span": 1.15}) - main = dict({"Name": "MN1", "P": 16.60, "E": 5.60, "Roach": 0.1, "BAD": 1.0}) - jib = dict({"Name": "J1", "I": 16.20, "J": 5.10, "LPG": 5.40, "HBI": 1.8}) - kite = dict({"Name": "A2", "area": 150.0, "vce": 9.55}) +def make_yd41(**overrides): + yacht = { + "Name": "YD41", + "Lwl": 11.90, + "Vol": 6.05, + "Bwl": 3.18, + "Tc": 0.4, + "WSA": 28.20, + "Tmax": 2.30, + "Amax": 1.051, + "Mass": 6500, + "Ff": 1.5, + "Fa": 1.5, + "Boa": 4.2, + "Loa": 12.5, + } + keel = {"Cu": 1.00, "Cl": 0.78, "Span": 1.90} + rudder = {"Cu": 0.48, "Cl": 0.22, "Span": 1.15} + main = {"Name": "MN1", "P": 16.60, "E": 5.60, "Roach": 0.1, "BAD": 1.0} + jib = {"Name": "J1", "I": 16.20, "J": 5.10, "LPG": 5.40, "HBI": 1.8} + kite = {"Name": "A2", "area": 150.0, "vce": 9.55} tws_range = np.arange(4.0, 7.0, 2.0).tolist() twa_range = np.linspace(30.0, 180.0, 5).tolist() @@ -48,17 +49,156 @@ def make_yd41(): "tws_range": tws_range, "twa_range": twa_range, } + d.update(overrides) return d +def post_vpp(data): + client = app.test_client() + return client.post("/api/vpp/", data=json.dumps(data), headers=HEADERS) + + def test_vpp_simulation(): d = make_yd41() + response = post_vpp(d) + assert response.status_code == 200 + + +def test_empty_tws_range_returns_400(): + d = make_yd41(tws_range=[]) + response = post_vpp(d) + assert response.status_code == 400 + assert "empty" in response.json["error"].lower() + + +def test_empty_twa_range_returns_400(): + d = make_yd41(twa_range=[]) + response = post_vpp(d) + assert response.status_code == 400 + assert "empty" in response.json["error"].lower() + + +def test_tws_out_of_range_returns_400(): + d = make_yd41(tws_range=[1.0, 5.0]) + response = post_vpp(d) + assert response.status_code == 400 + assert "outside valid bounds" in response.json["error"] + - json_string = json.dumps(d) - headers = {"content-type": "application/json", "Accept-Charset": "UTF-8"} +def test_tws_above_range_returns_400(): + d = make_yd41(tws_range=[10.0, 40.0]) + response = post_vpp(d) + assert response.status_code == 400 + assert "outside valid bounds" in response.json["error"] + +def test_twa_out_of_range_returns_400(): + d = make_yd41(twa_range=[-10.0, 90.0]) + response = post_vpp(d) + assert response.status_code == 400 + assert "outside valid bounds" in response.json["error"] + + +def test_missing_field_returns_400(): + d = make_yd41() + del d["keel"] + response = post_vpp(d) + assert response.status_code == 400 + assert "error" in response.json + + +def test_invalid_json_returns_400(): client = app.test_client() - response = client.post("/api/vpp/", data=json_string, headers=headers) + response = client.post("/api/vpp/", data="not json", headers=HEADERS) + assert response.status_code == 400 + + +def test_short_keel_with_type_field(): + """Short keel payload with explicit type='short' should succeed.""" + d = make_yd41() + d["keel"] = {"type": "short", "Length": 1.2, "Depth": 0.90, "Tc_ratio": 0.15} + response = post_vpp(d) + assert response.status_code == 200 + + +def test_short_keel_without_type_field(): + """Short keel payload detected from keys alone (no type field).""" + d = make_yd41() + d["keel"] = {"Length": 1.2, "Depth": 0.90, "Tc_ratio": 0.15} + response = post_vpp(d) + assert response.status_code == 200 + + +def test_fin_keel_with_type_field(): + """Fin keel payload with explicit type='fin' should succeed.""" + d = make_yd41() + d["keel"] = {"type": "fin", "Cu": 1.00, "Cl": 0.78, "Span": 1.90} + response = post_vpp(d) + assert response.status_code == 200 + + +def test_api_5dof_method(): + """API accepts method='5dof' and returns results.""" + d = make_yd41() + d["method"] = "5dof" + response = post_vpp(d) + assert response.status_code == 200 + results = np.array(response.json["results"]) + assert np.any(results[:, :, :, 0] > 0), "5-DOF should produce non-zero speeds" + + +def test_api_data_source(): + """API accepts data_source parameter.""" + d = make_yd41() + d["data_source"] = "orc" + response = post_vpp(d) + assert response.status_code == 200 + + +def test_api_sail_type_sym_kite(): + """API accepts sail_type in kite section and produces results.""" + d = make_yd41() + d["kite"]["sail_type"] = "sym_kite" + response = post_vpp(d) + assert response.status_code == 200 + results = np.array(response.json["results"]) + assert np.any(results[:, :, :, 0] > 0) + + +def test_api_sail_type_low_performance(): + """API accepts low-performance sail_type for main and jib.""" + d = make_yd41() + d["main"]["sail_type"] = "main_low" + d["jib"]["sail_type"] = "jib_low" + response = post_vpp(d) + assert response.status_code == 200 + results = np.array(response.json["results"]) + assert np.any(results[:, :, :, 0] > 0) + + +def test_api_roughness(): + """API accepts roughness parameter and rough hull is slower.""" + d_smooth = make_yd41(roughness=0.0) + d_rough = make_yd41(roughness=500e-6) + r_smooth = post_vpp(d_smooth) + r_rough = post_vpp(d_rough) + assert r_smooth.status_code == 200 + assert r_rough.status_code == 200 + speeds_smooth = np.array(r_smooth.json["results"])[:, :, :, 0].max() + speeds_rough = np.array(r_rough.json["results"])[:, :, :, 0].max() + assert speeds_smooth > speeds_rough, "Rough hull should be slower" + + +def test_api_wave_params(): + """API accepts Hs and Ts parameters.""" + d = make_yd41(Hs=1.0, Ts=6.0) + response = post_vpp(d) + assert response.status_code == 200 + +def test_api_wave_direction(): + """API accepts wave_direction parameter.""" + d = make_yd41(Hs=1.0, Ts=6.0, wave_direction=45.0) + response = post_vpp(d) assert response.status_code == 200 diff --git a/tests/test_hydro.py b/tests/test_hydro.py new file mode 100644 index 0000000..2c57a37 --- /dev/null +++ b/tests/test_hydro.py @@ -0,0 +1,129 @@ +"""Tests for HydroMod — roughness and wave resistance.""" + +import numpy as np + +from src.HydroMod import HydroMod +from src.SailMod import Jib, Kite, Main +from src.YachtMod import Keel, Rudder, Yacht + + +def _make_yacht(**kwargs): + """Create a YD41-like yacht with optional overrides.""" + defaults = dict( + Name="TestYacht", + Lwl=11.90, + Vol=6.05, + Bwl=3.18, + Tc=0.4, + WSA=28.20, + Tmax=2.30, + Amax=1.051, + Mass=6500, + Ff=1.5, + Fa=1.5, + Boa=4.2, + Loa=12.5, + App=[Keel(Cu=1.00, Cl=0.78, Span=1.90), Rudder(Cu=0.48, Cl=0.22, Span=1.15)], + Sails=[ + Main("MN1", P=16.60, E=5.60, Roach=0.1, BAD=1.0), + Jib("J1", I=16.20, J=5.10, LPG=5.40, HBI=1.8), + Kite("A2", area=150.0, vce=9.55), + ], + ) + defaults.update(kwargs) + return Yacht(**defaults) + + +class TestRoughness: + def test_roughness_increases_resistance(self): + """Higher roughness should increase total resistance.""" + yacht_smooth = _make_yacht(roughness=0.0) + yacht_rough = _make_yacht(roughness=300e-6) + + hydro_smooth = HydroMod(yacht_smooth) + hydro_rough = HydroMod(yacht_rough) + + Fx_smooth, _, _ = hydro_smooth.update(3.0, 5.0, 2.0) + Fx_rough, _, _ = hydro_rough.update(3.0, 5.0, 2.0) + + assert Fx_rough > Fx_smooth, "Rough hull should have more resistance" + + def test_zero_roughness_gives_smooth_cf(self): + """Zero roughness should give bare ITTC 1957 Cf (no allowance).""" + yacht = _make_yacht(roughness=0.0) + hydro = HydroMod(yacht) + hydro.vb = 3.0 + cf = hydro._cf(hydro.l) + + # ITTC 1957 only + Re = 3.0 * hydro.l / hydro.nu + cf_expected = 0.066 * (np.log10(Re) - 2.03) ** (-2) + assert abs(cf - cf_expected) < 1e-10 + + def test_default_roughness_adds_allowance(self): + """Default 150 μm roughness should add positive Cf allowance.""" + yacht = _make_yacht(roughness=150e-6) + hydro = HydroMod(yacht) + hydro.vb = 3.0 + cf_rough = hydro._cf(hydro.l) + + yacht_smooth = _make_yacht(roughness=0.0) + hydro_smooth = HydroMod(yacht_smooth) + hydro_smooth.vb = 3.0 + cf_smooth = hydro_smooth._cf(hydro_smooth.l) + + assert cf_rough > cf_smooth + + +class TestWaveResistance: + def test_no_waves_zero_resistance(self): + """Hs=0 should add zero wave resistance.""" + yacht = _make_yacht(Hs=0.0, Ts=0.0) + hydro = HydroMod(yacht) + raw = hydro._added_resistance_waves(90.0) + assert raw == 0.0 + + def test_waves_add_resistance(self): + """Hs > 0 should produce positive added resistance upwind.""" + yacht = _make_yacht(Hs=1.0, Ts=6.0) + hydro = HydroMod(yacht) + hydro.vb = 3.0 + raw = hydro._added_resistance_waves(45.0) + assert raw > 0.0, "Waves should add resistance upwind" + + def test_wave_resistance_increases_with_Hs(self): + """Larger waves should produce more resistance.""" + yacht_small = _make_yacht(Hs=0.5, Ts=6.0) + yacht_large = _make_yacht(Hs=1.5, Ts=6.0) + + hydro_small = HydroMod(yacht_small) + hydro_large = HydroMod(yacht_large) + hydro_small.vb = 3.0 + hydro_large.vb = 3.0 + + raw_small = hydro_small._added_resistance_waves(45.0) + raw_large = hydro_large._added_resistance_waves(45.0) + + assert raw_large > raw_small + + def test_wave_resistance_higher_upwind(self): + """Upwind wave resistance should exceed downwind.""" + yacht = _make_yacht(Hs=1.0, Ts=6.0) + hydro = HydroMod(yacht) + hydro.vb = 3.0 + + raw_upwind = hydro._added_resistance_waves(30.0) + raw_downwind = hydro._added_resistance_waves(150.0) + + assert raw_upwind > raw_downwind + + def test_flat_water_unchanged(self): + """VPP output with Hs=0 should match no-wave behaviour.""" + yacht = _make_yacht(Hs=0.0, Ts=0.0) + hydro = HydroMod(yacht) + + Fx_flat, Fy_flat, Mx_flat = hydro.update(3.0, 5.0, 2.0, twa=45.0) + Fx_no_twa, Fy_no_twa, Mx_no_twa = hydro.update(3.0, 5.0, 2.0, twa=0.0) + + # With no waves, twa shouldn't matter for resistance + assert abs(Fx_flat - Fx_no_twa) < 1e-6 diff --git a/tests/test_resistance.py b/tests/test_resistance.py index f557ecb..f395ff4 100644 --- a/tests/test_resistance.py +++ b/tests/test_resistance.py @@ -43,3 +43,15 @@ def test_Rr_interpolation(): np_testing.assert_approx_equal(YD41._interp_Rr((0.700, 3.0, 9.0)), 357.062, 4) np_testing.assert_approx_equal(YD41._interp_Rr((0.700, 9.0, 2.5)), 38.0526, 4) np_testing.assert_approx_equal(YD41._interp_Rr((0.700, 9.0, 9.0)), 42.2353, 4) + + +def test_build_interp_func_no_deprecation(): + """Verify build_interp_func doesn't use deprecated interp1d.""" + import warnings + + from src.UtilsMod import build_interp_func + with warnings.catch_warnings(): + warnings.simplefilter("error", DeprecationWarning) + func = build_interp_func("main", i=1) + result = func(30.0) + assert isinstance(float(result), float) diff --git a/tests/test_utils.py b/tests/test_utils.py index 10d2036..80d65c6 100644 --- a/tests/test_utils.py +++ b/tests/test_utils.py @@ -1,8 +1,9 @@ from src.SailMod import Jib, Kite, Main from src.YachtMod import Keel, Rudder, Yacht + def return_YD41_particulars(): - + YD41 = Yacht( Name="YD41", Lwl=11.90, @@ -25,4 +26,4 @@ def return_YD41_particulars(): Kite("A5", area=75.0, vce=2.75), ], ) - return YD41 \ No newline at end of file + return YD41 diff --git a/tests/test_vpp.py b/tests/test_vpp.py index f1e0211..b4a5d99 100644 --- a/tests/test_vpp.py +++ b/tests/test_vpp.py @@ -1,11 +1,13 @@ +import os import numpy as np -from tests.test_utils import return_YD41_particulars +from src.SailMod import Jib, Kite, Main from src.VPPMod import VPP -from src.SailMod import Jib, Main +from tests.test_utils import return_YD41_particulars -def test_single_sail_set(): + +def test_single_sail_set(tmp_path): YD41 = return_YD41_particulars() YD41_no_kite = YD41 @@ -22,5 +24,238 @@ def test_single_sail_set(): vpp.run(verbose=False) vpp.write("results") - vpp.polar(3, False) - vpp.SailChart(False) + + polar_path = str(tmp_path / "test_polar.png") + sail_path = str(tmp_path / "test_sail.png") + vpp.polar(3, True, fname=polar_path) + vpp.SailChart(True, fname=sail_path) + assert os.path.exists(polar_path), "Polar plot was not created" + assert os.path.exists(sail_path), "Sail chart was not created" + + +def test_sail_chart_no_deprecation_warning(tmp_path): + """Verify sail_chart doesn't use deprecated interp2d.""" + import warnings + + yacht = return_YD41_particulars() + yacht.sails = [ + Main("MN1", P=16.60, E=5.60, Roach=0.1, BAD=1.0), + Jib("J1", I=16.20, J=5.10, LPG=5.40, HBI=1.8), + ] + vpp = VPP(Yacht=yacht) + vpp.set_analysis( + tws_range=np.array([6.0, 10.0]), + twa_range=np.linspace(30.0, 180.0, 16), + ) + vpp.run(verbose=False) + fname = str(tmp_path / "test_sailchart.png") + with warnings.catch_warnings(): + warnings.simplefilter("error", DeprecationWarning) + vpp.SailChart(save=True, fname=fname) + assert os.path.exists(fname), "Sail chart was not created" + + +def test_phi_max_configurable(): + yacht = return_YD41_particulars() + yacht.sails = [Main("MN1", P=16.60, E=5.60, Roach=0.1, BAD=1.0), + Jib("J1", I=16.20, J=5.10, LPG=5.40, HBI=1.8)] + vpp = VPP(Yacht=yacht) + vpp.set_analysis(tws_range=np.array([10.0]), + twa_range=np.linspace(30.0, 180.0, 3), + phi_max=25.0) + assert vpp.phi_max == 25.0 + + +def test_run_without_analysis_raises(): + """Issue #46: raise with string literal should be a proper exception.""" + import pytest + + yacht = return_YD41_particulars() + vpp = VPP(Yacht=yacht) + with pytest.raises(RuntimeError, match="no analysis set"): + vpp.run() + + +def test_set_analysis_empty_tws_raises(): + import pytest + + yacht = return_YD41_particulars() + vpp = VPP(Yacht=yacht) + with pytest.raises(ValueError, match="TWS range is empty"): + vpp.set_analysis(tws_range=np.array([]), twa_range=np.linspace(30, 180, 5)) + + +def test_set_analysis_empty_twa_raises(): + import pytest + + yacht = return_YD41_particulars() + vpp = VPP(Yacht=yacht) + with pytest.raises(ValueError, match="TWA range is empty"): + vpp.set_analysis(tws_range=np.array([6.0, 10.0]), twa_range=np.array([])) + + +def test_set_analysis_tws_below_minimum_raises(): + import pytest + + yacht = return_YD41_particulars() + vpp = VPP(Yacht=yacht) + with pytest.raises(ValueError, match="outside valid bounds"): + vpp.set_analysis(tws_range=np.array([1.0, 5.0]), twa_range=np.linspace(30, 180, 5)) + + +def test_set_analysis_tws_above_maximum_raises(): + import pytest + + yacht = return_YD41_particulars() + vpp = VPP(Yacht=yacht) + with pytest.raises(ValueError, match="outside valid bounds"): + vpp.set_analysis(tws_range=np.array([10.0, 40.0]), twa_range=np.linspace(30, 180, 5)) + + +def test_set_analysis_twa_out_of_range_raises(): + import pytest + + yacht = return_YD41_particulars() + vpp = VPP(Yacht=yacht) + with pytest.raises(ValueError, match="outside valid bounds"): + vpp.set_analysis(tws_range=np.array([6.0, 10.0]), twa_range=np.array([-5.0, 90.0])) + + +def test_5dof_solver_runs(): + """5-DOF solver produces non-zero speeds.""" + yacht = return_YD41_particulars() + yacht.sails = [Main("MN1", P=16.60, E=5.60, Roach=0.1, BAD=1.0), + Jib("J1", I=16.20, J=5.10, LPG=5.40, HBI=1.8)] + vpp = VPP(Yacht=yacht) + vpp.set_analysis(tws_range=np.array([8.0]), + twa_range=np.linspace(40.0, 160.0, 4)) + vpp.run(method="5dof") + speeds = vpp.store[0, :, 0, 0] + assert np.any(speeds > 0), "5-DOF solver should produce non-zero speeds" + + +def test_5dof_vs_iterative_comparable(): + """5-DOF and iterative solvers should produce speeds within 15%.""" + yacht = return_YD41_particulars() + yacht.sails = [Main("MN1", P=16.60, E=5.60, Roach=0.1, BAD=1.0), + Jib("J1", I=16.20, J=5.10, LPG=5.40, HBI=1.8)] + + vpp_iter = VPP(Yacht=yacht) + vpp_iter.set_analysis(tws_range=np.array([8.0]), + twa_range=np.linspace(40.0, 160.0, 4)) + vpp_iter.run(method="iterative") + + vpp_5dof = VPP(Yacht=yacht) + vpp_5dof.set_analysis(tws_range=np.array([8.0]), + twa_range=np.linspace(40.0, 160.0, 4)) + vpp_5dof.run(method="5dof") + + speeds_iter = vpp_iter.store[0, :, 0, 0] + speeds_5dof = vpp_5dof.store[0, :, 0, 0] + + # Compare only non-zero entries + mask = speeds_iter > 0.1 + if np.any(mask): + ratio = speeds_5dof[mask] / speeds_iter[mask] + assert np.all(ratio > 0.85) and np.all(ratio < 1.15), ( + f"5-DOF vs iterative speed ratio out of 15% band: {ratio}" + ) + + +def test_invalid_method_raises(): + import pytest + + yacht = return_YD41_particulars() + yacht.sails = [Main("MN1", P=16.60, E=5.60, Roach=0.1, BAD=1.0), + Jib("J1", I=16.20, J=5.10, LPG=5.40, HBI=1.8)] + vpp = VPP(Yacht=yacht) + vpp.set_analysis(tws_range=np.array([8.0]), + twa_range=np.linspace(40.0, 160.0, 3)) + with pytest.raises(ValueError, match="Unknown method"): + vpp.run(method="bogus") + + +def test_sail_type_main_low(): + """Main with sail_type='main_low' loads lower CL coefficients.""" + main_hi = Main("MN1", P=16.60, E=5.60, Roach=0.1, BAD=1.0) + main_lo = Main("MN1", P=16.60, E=5.60, Roach=0.1, BAD=1.0, sail_type="main_low") + # At peak CL angle (~28 deg), low should have less lift + assert main_lo.cl(28) < main_hi.cl(28) + + +def test_sail_type_jib_low(): + """Jib with sail_type='jib_low' loads lower CL coefficients.""" + jib_hi = Jib("J1", I=16.20, J=5.10, LPG=5.40, HBI=1.8) + jib_lo = Jib("J1", I=16.20, J=5.10, LPG=5.40, HBI=1.8, sail_type="jib_low") + assert jib_lo.cl(27) < jib_hi.cl(27) + + +def test_sail_type_sym_kite(): + """Symmetric spinnaker has higher CL than the default kite.""" + kite_default = Kite("A2", area=150.0, vce=9.55) + kite_sym = Kite("A2", area=150.0, vce=9.55, sail_type="sym_kite") + # ORC symmetric spinnaker has peak CL ~1.456 vs default ~1.08 + assert kite_sym.cl(67) > kite_default.cl(67) + + +def test_sail_type_asym_kite_variants(): + """Asymmetric spinnaker variants load and produce different coefficients.""" + kite_cl = Kite("A2", area=150.0, vce=9.55, sail_type="asym_cl_kite") + kite_pole = Kite("A2", area=150.0, vce=9.55, sail_type="asym_pole_kite") + # Both should produce non-zero CL at 67 deg + assert kite_cl.cl(67) > 1.0 + assert kite_pole.cl(67) > 1.0 + # Pole tack has slightly higher peak CL than centerline + assert kite_pole.cl(67) > kite_cl.cl(67) + + +def test_parallel_matches_iterative(): + """Parallel solver must produce identical results to sequential iterative.""" + yacht = return_YD41_particulars() + yacht.sails = [Main("MN1", P=16.60, E=5.60, Roach=0.1, BAD=1.0), + Jib("J1", I=16.20, J=5.10, LPG=5.40, HBI=1.8), + Kite("A2", area=150.0, vce=9.55)] + + tws_range = np.array([6.0, 10.0]) + twa_range = np.linspace(40.0, 160.0, 5) + + vpp_seq = VPP(Yacht=yacht) + vpp_seq.set_analysis(tws_range=tws_range, twa_range=twa_range) + vpp_seq.run(method="iterative") + + vpp_par = VPP(Yacht=yacht) + vpp_par.set_analysis(tws_range=tws_range, twa_range=twa_range) + vpp_par.run(method="parallel") + + np.testing.assert_allclose( + vpp_par.store, vpp_seq.store, atol=1e-6, + err_msg="Parallel results differ from sequential iterative" + ) + + +def test_parallel_method_accepted(): + """VPP.run(method='parallel') should not raise.""" + yacht = return_YD41_particulars() + yacht.sails = [Main("MN1", P=16.60, E=5.60, Roach=0.1, BAD=1.0), + Jib("J1", I=16.20, J=5.10, LPG=5.40, HBI=1.8)] + vpp = VPP(Yacht=yacht) + vpp.set_analysis(tws_range=np.array([8.0]), + twa_range=np.linspace(40.0, 160.0, 3)) + vpp.run(method="parallel") + assert np.any(vpp.store[0, :, 0, 0] > 0) + + +def test_sym_kite_vpp_runs(): + """VPP runs with symmetric spinnaker coefficients.""" + from src.YachtMod import Keel, Rudder, Yacht + yacht = Yacht( + Name="YD41", Lwl=11.90, Vol=6.05, Bwl=3.18, Tc=0.4, WSA=28.20, + Tmax=2.30, Amax=1.051, Mass=6500, Ff=1.5, Fa=1.5, Boa=4.2, Loa=12.5, + App=[Keel(Cu=1.00, Cl=0.78, Span=1.90), Rudder(Cu=0.48, Cl=0.22, Span=1.15)], + Sails=[Main("MN1", P=16.60, E=5.60, Roach=0.1, BAD=1.0), + Jib("J1", I=16.20, J=5.10, LPG=5.40, HBI=1.8), + Kite("A2", area=150.0, vce=9.55, sail_type="sym_kite")]) + vpp = VPP(Yacht=yacht) + vpp.set_analysis(tws_range=np.array([10.0]), twa_range=np.linspace(40.0, 160.0, 4)) + vpp.run() + assert np.any(vpp.store[0, :, :, 0] > 0) From 1c3ab042c90d6e3542779c11b6f54611383f60be Mon Sep 17 00:00:00 2001 From: Thomas Dickson Date: Wed, 17 Jun 2026 20:42:35 +0100 Subject: [PATCH 2/2] =?UTF-8?q?fix:=20revert=20api.py=20to=20master=20?= =?UTF-8?q?=E2=80=94=20ShortKeel=20belongs=20in=20sail-daring=20PR?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- src/api.py | 74 ++++-------------- tests/test_api.py | 192 +++++++--------------------------------------- 2 files changed, 42 insertions(+), 224 deletions(-) diff --git a/src/api.py b/src/api.py index 56a3867..14ad2f6 100644 --- a/src/api.py +++ b/src/api.py @@ -13,7 +13,7 @@ sys.path.append(os.path.realpath(".")) from src.SailMod import Jib, Kite, Main from src.VPPMod import VPP -from src.YachtMod import Keel, Rudder, ShortKeel, Yacht +from src.YachtMod import Keel, Rudder, Yacht app = Flask(__name__) @@ -27,35 +27,17 @@ def ping(): def data_to_vpp(data: Dict[str, Any]) -> VPP: - - keel_data = data["keel"] - keel_type = keel_data.get("type", "fin") - # Also detect from keys if type is missing or inconsistent - if keel_type == "short" or "Length" in keel_data: - keel = ShortKeel( - Length=float(keel_data["Length"]), - Depth=float(keel_data["Depth"]), - Tc_ratio=float(keel_data.get("Tc_ratio", 0.15)), - ) - else: - keel = Keel( - Cu=float(keel_data["Cu"]), - Cl=float(keel_data["Cl"]), - Span=float(keel_data["Span"]), - ) + + keel = Keel( + Cu=float(data["keel"]["Cu"]), + Cl=float(data["keel"]["Cl"]), + Span=float(data["keel"]["Span"]) + ) rudder = Rudder( - Cu=float(data["rudder"]["Cu"]), - Cl=float(data["rudder"]["Cu"]), + Cu=float(data["rudder"]["Cu"]), + Cl=float(data["rudder"]["Cu"]), Span=float(data["rudder"]["Span"]) ) - # Environment parameters - roughness = float(data.get("roughness", 150e-6)) - Hs = float(data.get("Hs", 0.0)) - Ts = float(data.get("Ts", 0.0)) - wave_direction = data.get("wave_direction") - if wave_direction is not None: - wave_direction = float(wave_direction) - yacht = Yacht( Name=data["yacht"]["Name"], Lwl=float(data["yacht"]["Lwl"]), @@ -70,10 +52,6 @@ def data_to_vpp(data: Dict[str, Any]) -> VPP: Fa=float(data["yacht"]["Fa"]), Boa=float(data["yacht"]["Boa"]), Loa=float(data["yacht"]["Loa"]), - roughness=roughness, - Hs=Hs, - Ts=Ts, - wave_direction=wave_direction, App=[keel, rudder], Sails=[ Main( @@ -82,10 +60,6 @@ def data_to_vpp(data: Dict[str, Any]) -> VPP: E=float(data["main"]["E"]), Roach=float(data["main"]["Roach"]), BAD=float(data["main"]["BAD"]), - data_source=data.get("data_source", "orc"), - cl_data=data["main"].get("cl_data"), - cd_data=data["main"].get("cd_data"), - sail_type=data["main"].get("sail_type"), ), Jib( name=data["jib"]["Name"], @@ -93,48 +67,32 @@ def data_to_vpp(data: Dict[str, Any]) -> VPP: J=float(data["jib"]["J"]), LPG=float(data["jib"]["LPG"]), HBI=float(data["jib"]["HBI"]), - data_source=data.get("data_source", "orc"), - cl_data=data["jib"].get("cl_data"), - cd_data=data["jib"].get("cd_data"), - sail_type=data["jib"].get("sail_type"), ), Kite( name=data["kite"]["Name"], area=float(data["kite"]["area"]), vce=float(data["kite"]["vce"]), - data_source=data.get("data_source", "orc"), - cl_data=data["kite"].get("cl_data"), - cd_data=data["kite"].get("cd_data"), - sail_type=data["kite"].get("sail_type"), ), ], ) - + vpp = VPP(Yacht=yacht) vpp.set_analysis( tws_range=np.array(data["tws_range"]), twa_range=np.array(data["twa_range"]), ) - return vpp, data.get("method", "iterative") - + return vpp + @app.route("/api/vpp/", methods=["POST"]) def makevppresults(): data = request.get_json() - if data is None: - return jsonify({"error": "Request body must be valid JSON."}), 400 - try: - vpp, method = data_to_vpp(data) - except (KeyError, TypeError, ValueError) as e: - logging.warning("Invalid VPP input: %s", e) - return jsonify({"error": f"Invalid input: {e}"}), 400 + # TODO: Support multiple implementations of different sails: require API design + # TODO: Error handling incorrect ranges - try: - vpp.run(verbose=True, method=method) - except Exception as e: - logging.exception("VPP simulation failed") - return jsonify({"error": f"Simulation failed: {e}"}), 500 + vpp = data_to_vpp(data) + vpp.run(verbose=True) return jsonify(vpp.results()) diff --git a/tests/test_api.py b/tests/test_api.py index f1e068a..53aab60 100644 --- a/tests/test_api.py +++ b/tests/test_api.py @@ -4,9 +4,6 @@ from src.api import app -HEADERS = {"content-type": "application/json", "Accept-Charset": "UTF-8"} - - def test_ping_route(): client = app.test_client() response = client.get("/ping") @@ -14,27 +11,29 @@ def test_ping_route(): assert response.data.decode("utf-8") == "Pong! The server is up and running." -def make_yd41(**overrides): - yacht = { - "Name": "YD41", - "Lwl": 11.90, - "Vol": 6.05, - "Bwl": 3.18, - "Tc": 0.4, - "WSA": 28.20, - "Tmax": 2.30, - "Amax": 1.051, - "Mass": 6500, - "Ff": 1.5, - "Fa": 1.5, - "Boa": 4.2, - "Loa": 12.5, - } - keel = {"Cu": 1.00, "Cl": 0.78, "Span": 1.90} - rudder = {"Cu": 0.48, "Cl": 0.22, "Span": 1.15} - main = {"Name": "MN1", "P": 16.60, "E": 5.60, "Roach": 0.1, "BAD": 1.0} - jib = {"Name": "J1", "I": 16.20, "J": 5.10, "LPG": 5.40, "HBI": 1.8} - kite = {"Name": "A2", "area": 150.0, "vce": 9.55} +def make_yd41(): + yacht = dict( + { + "Name": "YD41", + "Lwl": 11.90, + "Vol": 6.05, + "Bwl": 3.18, + "Tc": 0.4, + "WSA": 28.20, + "Tmax": 2.30, + "Amax": 1.051, + "Mass": 6500, + "Ff": 1.5, + "Fa": 1.5, + "Boa": 4.2, + "Loa": 12.5, + } + ) + keel = dict({"Cu": 1.00, "Cl": 0.78, "Span": 1.90}) + rudder = dict({"Cu": 0.48, "Cl": 0.22, "Span": 1.15}) + main = dict({"Name": "MN1", "P": 16.60, "E": 5.60, "Roach": 0.1, "BAD": 1.0}) + jib = dict({"Name": "J1", "I": 16.20, "J": 5.10, "LPG": 5.40, "HBI": 1.8}) + kite = dict({"Name": "A2", "area": 150.0, "vce": 9.55}) tws_range = np.arange(4.0, 7.0, 2.0).tolist() twa_range = np.linspace(30.0, 180.0, 5).tolist() @@ -49,156 +48,17 @@ def make_yd41(**overrides): "tws_range": tws_range, "twa_range": twa_range, } - d.update(overrides) return d -def post_vpp(data): - client = app.test_client() - return client.post("/api/vpp/", data=json.dumps(data), headers=HEADERS) - - def test_vpp_simulation(): d = make_yd41() - response = post_vpp(d) - assert response.status_code == 200 - - -def test_empty_tws_range_returns_400(): - d = make_yd41(tws_range=[]) - response = post_vpp(d) - assert response.status_code == 400 - assert "empty" in response.json["error"].lower() - - -def test_empty_twa_range_returns_400(): - d = make_yd41(twa_range=[]) - response = post_vpp(d) - assert response.status_code == 400 - assert "empty" in response.json["error"].lower() - - -def test_tws_out_of_range_returns_400(): - d = make_yd41(tws_range=[1.0, 5.0]) - response = post_vpp(d) - assert response.status_code == 400 - assert "outside valid bounds" in response.json["error"] - -def test_tws_above_range_returns_400(): - d = make_yd41(tws_range=[10.0, 40.0]) - response = post_vpp(d) - assert response.status_code == 400 - assert "outside valid bounds" in response.json["error"] + json_string = json.dumps(d) + headers = {"content-type": "application/json", "Accept-Charset": "UTF-8"} - -def test_twa_out_of_range_returns_400(): - d = make_yd41(twa_range=[-10.0, 90.0]) - response = post_vpp(d) - assert response.status_code == 400 - assert "outside valid bounds" in response.json["error"] - - -def test_missing_field_returns_400(): - d = make_yd41() - del d["keel"] - response = post_vpp(d) - assert response.status_code == 400 - assert "error" in response.json - - -def test_invalid_json_returns_400(): client = app.test_client() - response = client.post("/api/vpp/", data="not json", headers=HEADERS) - assert response.status_code == 400 - - -def test_short_keel_with_type_field(): - """Short keel payload with explicit type='short' should succeed.""" - d = make_yd41() - d["keel"] = {"type": "short", "Length": 1.2, "Depth": 0.90, "Tc_ratio": 0.15} - response = post_vpp(d) - assert response.status_code == 200 - - -def test_short_keel_without_type_field(): - """Short keel payload detected from keys alone (no type field).""" - d = make_yd41() - d["keel"] = {"Length": 1.2, "Depth": 0.90, "Tc_ratio": 0.15} - response = post_vpp(d) - assert response.status_code == 200 - - -def test_fin_keel_with_type_field(): - """Fin keel payload with explicit type='fin' should succeed.""" - d = make_yd41() - d["keel"] = {"type": "fin", "Cu": 1.00, "Cl": 0.78, "Span": 1.90} - response = post_vpp(d) - assert response.status_code == 200 - - -def test_api_5dof_method(): - """API accepts method='5dof' and returns results.""" - d = make_yd41() - d["method"] = "5dof" - response = post_vpp(d) - assert response.status_code == 200 - results = np.array(response.json["results"]) - assert np.any(results[:, :, :, 0] > 0), "5-DOF should produce non-zero speeds" - - -def test_api_data_source(): - """API accepts data_source parameter.""" - d = make_yd41() - d["data_source"] = "orc" - response = post_vpp(d) - assert response.status_code == 200 - - -def test_api_sail_type_sym_kite(): - """API accepts sail_type in kite section and produces results.""" - d = make_yd41() - d["kite"]["sail_type"] = "sym_kite" - response = post_vpp(d) - assert response.status_code == 200 - results = np.array(response.json["results"]) - assert np.any(results[:, :, :, 0] > 0) - - -def test_api_sail_type_low_performance(): - """API accepts low-performance sail_type for main and jib.""" - d = make_yd41() - d["main"]["sail_type"] = "main_low" - d["jib"]["sail_type"] = "jib_low" - response = post_vpp(d) - assert response.status_code == 200 - results = np.array(response.json["results"]) - assert np.any(results[:, :, :, 0] > 0) - - -def test_api_roughness(): - """API accepts roughness parameter and rough hull is slower.""" - d_smooth = make_yd41(roughness=0.0) - d_rough = make_yd41(roughness=500e-6) - r_smooth = post_vpp(d_smooth) - r_rough = post_vpp(d_rough) - assert r_smooth.status_code == 200 - assert r_rough.status_code == 200 - speeds_smooth = np.array(r_smooth.json["results"])[:, :, :, 0].max() - speeds_rough = np.array(r_rough.json["results"])[:, :, :, 0].max() - assert speeds_smooth > speeds_rough, "Rough hull should be slower" - - -def test_api_wave_params(): - """API accepts Hs and Ts parameters.""" - d = make_yd41(Hs=1.0, Ts=6.0) - response = post_vpp(d) - assert response.status_code == 200 - + response = client.post("/api/vpp/", data=json_string, headers=headers) -def test_api_wave_direction(): - """API accepts wave_direction parameter.""" - d = make_yd41(Hs=1.0, Ts=6.0, wave_direction=45.0) - response = post_vpp(d) assert response.status_code == 200