diff --git a/rocket_twin/systems/tank/tank.py b/rocket_twin/systems/tank/tank.py index bd621e3..c4a32b8 100644 --- a/rocket_twin/systems/tank/tank.py +++ b/rocket_twin/systems/tank/tank.py @@ -1,4 +1,9 @@ +import numpy as np from cosapp.base import System +from OCC.Core.BRepGProp import brepgprop +from OCC.Core.gp import gp_Vec, gp_XYZ +from OCC.Core.GProp import GProp_GProps +from pyoccad.create import CreateCylinder class Tank(System): @@ -24,9 +29,17 @@ class Tank(System): def setup(self): # Geometry + self.add_inward("radius", 1.0, desc="Tank base radius", unit="m") + self.add_inward("height", 1.0, desc="Tank height", unit="m") self.add_inward("weight_s", 1.0, desc="Structure weight", unit="kg") self.add_inward("weight_max", 5.0, desc="Maximum fuel capacity", unit="kg") + # pyoccad model + shape = CreateCylinder.from_base_and_dir( + np.zeros(3), gp_Vec(gp_XYZ(0, 0, self.height)), self.radius + ) + self.add_inward("shape", shape, desc="pyoccad tank model") + # Inputs self.add_inward("w_in", 0.0, desc="Fuel income rate", unit="kg/s") @@ -38,6 +51,7 @@ def setup(self): self.add_outward("weight", 1.0, desc="Weight", unit="kg") self.add_outward("cg", 1.0, desc="Center of gravity", unit="m") self.add_outward("w_out", 0.0, desc="Fuel output rate", unit="kg/s") + self.add_outward("I", np.zeros((3, 3)), desc="Tank inertia tensor", unit="kg*m**2") # Transient self.add_transient("weight_p", der="w_in - w_out", desc="Propellant weight") @@ -45,3 +59,9 @@ def setup(self): def compute(self): self.w_out = self.w_out_max * self.w_command self.weight = self.weight_s + self.weight_p + + v_prop = GProp_GProps() + brepgprop.VolumeProperties(self.shape, v_prop) + inertia = v_prop.MatrixOfInertia() + for i, j in zip(range(3), range(3)): + self.I[i, j] = inertia.Value(i + 1, j + 1) diff --git a/rocket_twin/tests/test_tank.py b/rocket_twin/tests/test_tank.py index eeccbc3..b0535c1 100644 --- a/rocket_twin/tests/test_tank.py +++ b/rocket_twin/tests/test_tank.py @@ -19,6 +19,9 @@ def test_fuel(self): sys.run_drivers() np.testing.assert_allclose(sys.weight, 16.0, atol=10 ** (-10)) + np.testing.assert_allclose( + sys.I.diagonal(), np.array([np.pi / 3, np.pi / 3, np.pi / 2]), atol=10 ** (-10) + ) def test_flight(self): sys = Tank("sys")