Attitude dynamics

configs/TLI_Oct2019.json

@@ -7,7 +7,7 @@
         "vel_x": -534.084,
         "vel_y": -3792.878,
         "vel_z": -867.495,
-        "ang_vel_x": 4.5,
+        "h_x": 4.5,
         "time": 1539102600
     },
     "models" : ["pos", "gyro"]

configs/apollo12SIVB.json

@@ -7,7 +7,7 @@
         "vel_x": 5106.47906333 ,
         "vel_y": -4728.66201988,
         "vel_z": 4572.38411289 ,
-        "ang_vel_x": 4.5,
+        "h_x": 4.5,
         "time": -4058816
     },
     "models" : ["pos", "gyro"]

configs/attitude_dynamics_test.json

@@ -1,8 +1,12 @@
 {
+    "parameters": {
+        "max_iter": 1e3
+    },
     "initial_condition": {
         "x": 100000000.0,
-        "ang_vel_x": 4.5,
-        "ang_vel_y": 1.0,
+        "h_x": 4.5,
+        "h_y": 0.2,
+        "h_z": 0.0,
         "time": 0.0,
         "quat_v1": 1.0,
         "quat_v2": 0.0,

configs/freefall.json

@@ -4,7 +4,7 @@
         "x": 100000000.0,
         "y": 1000.0,
         "z": 1000.0,
-        "ang_vel_x": 4.5,
+        "h_x": 4.5,
         "time": 0.0,
         "quat_v1": 1.0,
         "quat_v2": 0.0,

configs/lro.json

@@ -7,7 +7,7 @@
         "vel_x": 0.2069511992715368,
         "vel_y": -0.2105451990331415,
         "vel_z": -1.157255066656339,
-        "ang_vel_x": 4.5,
+        "h_x": 4.5,
         "time": 1651885860
     },
     "models" : ["pos", "gyro"]

configs/roadster.json

@@ -7,7 +7,7 @@
         "vel_x": -6.837449791362452,
         "vel_y": -3.855569832386001,
         "vel_z": -2.658902141013662,
-        "ang_vel_x": 4.5,
+        "h_x": 4.5,
         "time": 1517990400
     },
     "models" : ["pos", "gyro"]

configs/schema.json

@@ -42,13 +42,13 @@
 		"initial_condition": {
 			"type": "object",
 			"properties": {
-				"ang_vel_x": {
+				"h_x": {
 					"type": "number"
 				},
-				"ang_vel_y": {
+				"h_y": {
 					"type": "number"
 				},
-				"ang_vel_z": {
+				"h_z": {
 					"type": "number"
 				},
 				"quat_v1": {

configs/test_angles.json

@@ -1,9 +1,9 @@
 {
     "initial_condition": {
 
-        "ang_vel_x" : 5.0, 
-        "ang_vel_y" : 5.0, 
-        "ang_vel_z" : 5.0, 
+        "h_x" : 5.0, 
+        "h_y" : 5.0, 
+        "h_z" : 5.0, 
 
         "quat_v1" : 10.0, 
         "quat_v2" : 10.0, 

configs/test_zeroes.json

@@ -10,9 +10,9 @@
     },
     "initial_condition": {
 
-        "ang_vel_x" : 0.0, 
-        "ang_vel_y" : 0.0, 
-        "ang_vel_z" : 0.0, 
+        "h_x" : 0.0, 
+        "h_y" : 0.0, 
+        "h_z" : 0.0, 
 
         "quat_v1" : 0.0, 
         "quat_v2" : 0.0, 

configs/tli.json

@@ -7,7 +7,7 @@
         "vel_x": -534.084,
         "vel_y": -3792.878,
         "vel_z": -867.495,
-        "ang_vel_x": 4.5,
+        "h_x": 4.5,
         "time": 1539102600
     },
     "models" : ["pos", "gyro"]

src/core/config.py

@@ -7,7 +7,7 @@
 from utils.constants import ModelEnum
 
 from core.parameters import Parameters
-from core.state.statetime import StateTime
+from core.state import statetime
 
 
 class MutationException(Exception):
@@ -34,7 +34,8 @@ def __init__(
     ):
 
         self.param = Parameters(param_dict=parameters)
-        self.init_cond = StateTime.from_dict(initial_condition)
+        self.init_cond = statetime.StateTime.from_dict(initial_condition)
+        statetime.apply_coupled_initial_conditions(self.init_cond)
 
         # convert string list to model list
         model_objs = []
@@ -57,9 +58,12 @@ def __delattr__(self, __name) -> None:
     def make_config(cls, path_str: str):
         """make_config creates a config object from json file in the proposed location.
 
-        It is dependent on state.py, and parameters.py. Changes in these files will affect this method. It will be used by main.py. Changes in this method might affect the functionality of Main.
+        It is dependent on state.py, and parameters.py. Changes in these files will affect this method. It will be used
+        by main.py. Changes in this method might affect the functionality of Main.
 
-        To construct a json file, consult schema.json and example.json in the 'configs' folder. All properties are optional, default values will be inserted if a field is not specified. However, if a property is specified its type and format has to be correct.
+        To construct a json file, consult schema.json and example.json in the 'configs' folder. All properties are
+        optional, default values will be inserted if a field is not specified. However, if a property is specified
+        its type and format has to be correct.
 
         Raises: `JsonError` if json file is not well defined.
 
@@ -74,7 +78,8 @@ def make_config(cls, path_str: str):
             except ValidationError:
                 raise JsonError("Schema validation failed.")
 
-            # Checking if gyro_bias ans gyro_noise are in the correct format if thery are specified. (other type verification is done by schema.json)
+            # Checking if gyro_bias ans gyro_noise are in the correct format if thery are specified. 
+            # (other type verification is done by schema.json)
             try:  # validate gyro_bias is a list of length 3
                 gyro_bias = data.get("parameters").get("gyro_bias")
                 if len(gyro_bias) != 3:

src/core/models/derived_models.py

@@ -1,11 +1,10 @@
 from abc import abstractmethod
 from core.state.state import State
 from utils.astropy_util import get_body_position
-from typing import Dict
 import numpy as np
-from utils.constants import BodyEnum
+from utils.constants import BodyEnum, State_Type
 from core.models.model_base import Model
-from typing import Union, Tuple, List
+from typing import Union, Tuple, List, Dict
 
 
 class DerivedStateModel(Model):
@@ -101,4 +100,101 @@ def evaluate(self, t: float, state: State) -> Dict[str, np.ndarray]:
         }
 
 
-DERIVED_MODEL_LIST: List[DerivedStateModel] = [DerivedPosition(), DerivedAttitude()]
+class InertiaModel(DerivedStateModel):
+    def evaluate(self, _: float, state: State) -> Dict[str, State_Type]:
+        """the _b postfix means the value is in the spacecraft's body frame"""
+        fill_frac = state.fill_frac
+
+        dcm = np.array([[0, 1, 0], [0, 0, -1], [-1, 0, 0]], dtype=np.int32)
+        dcmT = np.transpose(dcm)
+        # not sure why we do the above. Did i do that?
+        # probably to get the (incorrect) solidworks frame into the spacecraft
+        # body frame
+
+        # Inertia tensor when full. Structure is:
+        # [[Ixx, Ixy, Ixz],
+        #  [Iyx, Iyy, Iyz],
+        #  [Izx, Izy, Izz]].
+        # Units are (kg * m^2).
+        idf = (
+            np.array(
+                [
+                    [933513642.20, 260948256.18, 430810000.30],
+                    [260948256.18, 1070855457.07, 387172545.62],
+                    [430810000.30, 387172545.62, 629606813.62],
+                ],
+                dtype=np.float64,
+            )
+            * 1e-9
+        )
+        # inertia tensor in correct body frame
+        idf_b = np.matmul(dcm, idf)
+
+        # Inertia tensor at 125 mL. Structure is:
+        # [[Ixx, Ixy, Ixz],
+        #  [Iyx, Iyy, Iyz],
+        #  [Izx, Izy, Izz]].
+        # Units are (kg * m^2).
+        idi = (
+            np.array(
+                [
+                    [855858994.14, 229481961.55, 377087149.13],
+                    [229481961.55, 963124288.81, 353943859.15],
+                    [377087149.13, 353943859.15, 559805590.96],
+                ],
+                dtype=np.float64,
+            )
+            * 1e-9
+        )
+        # inertia tensor in correct body frame
+        idi_b = np.matmul(dcm, idi)
+
+        # Determine inertia tensor for Oxygen via linear interpolation as a function of fill
+        # fraction.
+
+        # TODO: determine whether linear fit is accurate enough
+        ineria_matrix_b = (idf_b - idi_b) * fill_frac + idi_b
+
+        # now that we know our moment of inertia (I matrix), we can calculate our angular
+        # velocities:
+        # page 14 of the 4060 Attitude Dynamics handout:
+        # h = I dot omega
+        # --> omega = inverse(I) matmul h
+
+        angular_momentum = np.array([[state.h_x, state.h_y, state.h_z]])
+
+        angular_velocity = np.matmul(np.linalg.inv(ineria_matrix_b), angular_momentum.T)
+        return {
+            "I": ineria_matrix_b,
+            "ang_vel": angular_velocity,
+        }
+
+
+class KaneModel(DerivedStateModel):
+    """Calculates the Kane damping coefficient from 2016 simulation data by K. Doyle."""
+
+    def evaluate(self, _: float, state: State) -> Dict[str, State_Type]:
+        # Coefficients below are from Kyle's work.
+        # TODO: Update them when we conduct a new Ansys analysis.
+        k = 0.00085
+        factor = 1.2
+
+        N = 50
+        kane = fill = tau1 = tau2 = np.zeros(N)
+        for i in range(N):
+            fill[i] = (i - 1) / N
+            tau1[i] = k * fill[i]
+            tau2[i] = factor * k * (1 - fill[i])
+            kane[i] = -np.sqrt(tau1[i] ** 2 + tau2[i] ** 2)
+
+        kane = kane - np.max(kane) + k
+        kane = kane - np.min(kane)
+        kane = kane * k / np.max(kane)
+
+        q = int(np.round(state.fill_frac / 0.02 + 1))
+        c = kane[q]
+
+        return {"kane_c": c}
+
+
+DERIVED_MODEL_LIST: List[DerivedStateModel] = [DerivedPosition(), DerivedAttitude(), InertiaModel(), KaneModel()]