Interacting with the environment during propagation

Each Body object and its constituent members is updated to the current state and time automatically during the numerical propagation. We stress that only those objects that are relevant for a given propagation are updated every time step. See A single function evaluation for a more detailed discussion on what happens during a function evaluation of the state derivative.

Even though the environment is updated automatically, in various cases a user has control over how it gets updated. This is the case when using any of the custom models in Tudat. Typical examples of such models are aerodynamic guidance, or thrust guidance. When defining such custom models, you will in many cases need to access the properties of bodies in the environment to perform your calculations. These properties are stored as member variables of Body objects. This page lists how to access such properties, and how to use/interpret them (if relevant). Some examples are:

Retrieving the translational state of a spacecraft w.r.t. a central body to calculate its thrust direction
Retrieving the current atmospheric density during entry to calculate the required vehicle orientation
Retrieving properties that cannot be saved as a dependent variable to calculate a custom dependent variable (note - this option does not influence the propagation results, but can be used to obtain more flexible output, see custom_dependent_variable())
…

Below, we list how to retrieve the relevant information of the current properties of a body in Tudat. Note that these can only be used during a propagation, and therefore are only relevant when setting up custom models.

In what follows below, we will assume that you have created a SystemOfBodies variable named bodies, from which you want to access various properties during the propagation.

Translational state

Retrieved directly from the Body object with the state function in Cartesian elements. Note that this state is always w.r.t. the global frame origin and orientation (see translational state reference frames). To retrieve the state of one body w.r.t. the other, if neither body is designated as the global frame origin, you can use:
current_relate_state = bodies.get('Vehicle').state - bodies.get('Earth').state
to retrieve the current state of a body named ‘Vehicle’ w.r.t. the Earth. Note that if the Earth is the global frame origin, the above will still work fine. However, the Earth state will be a zero-vector, and its subtraction from the vehicle state may be omitted for efficiency purposes

NOTE: If you are only interested in the position or velocity components, you can use the position or velocity functions.

Rotational state

The rotational state of a body during the propagation is defined by the orientation of its body-fixed frame w.r.t. the global frame orientation, which is always an inertial orientation (see Frames in the Environment). This orientation is defined internally by a quaternion (see Definition of rotational state), but during a simulation a user will interact with the rotation matrix if they wish to use any current rotation in their custom models. The rotation matrix \(\mathbf{R}^{B/I}\) from the inertial \(I\) to body-fixed frame \(B\) is retrieved from a Body object using the inertial_to_body_fixed_frame function. The inverse rotation matrix \(\mathbf{R}^{I/B}\) (body-fixed to inertial) is retrieved using the body_fixed_to_inertial_frame function.

The time-derivative of the orientation is provided in two formulations (with equivalent information content): the angular velocity vector of the body-fixed frame, and the time derivative of the rotation matrix. The angular velocity vector (of \(B\) w.r.t. inertial space), in inertial and body-fixed coordinates, is obtained from the inertial_angular_velocity and body_fixed_angular_velocity functions respectively. Note that the latter is the formulation that is used to represent the time-variation of the rotation when propagating rotational dynamics (see Rotational Dynamics). Alternatively, the time-derivative of the rotation matrix from inertial to body-fixed frame \(\dot{\mathbf{R}}^{B/I}\) is given by inertial_to_body_fixed_frame_derivative, while the derivative of the inverse rotation \(\dot{\mathbf{R}}^{I/B}\) is taken from body_fixed_to_inertial_frame_derivative.

Body inertial mass

Retrieved directly from a Body object with the mass function. Note that this mass is not (at least, not by definition) the mass used for calculation of gravitional interactions (the gravitational mass \(m_{g}\), as you would find it in Newton’s law of gravity (\(a=\frac{Gm_{g}}{r^{2}}\)), but the mass used to convert forces to accelerations and vice versa (the inertial mass \(m_{i}\), as you would find it in Newton’s law of motion \(F=m_{i}a\)). To the best of our knowledge the two masses are equal for all bodies, but various formulations of general relativity predict a difference between the two. Moreover, we have found it useful to not automatically define a gravity field for any body which happens to have a mass assigned to it. For instance, a spacecraft will have an (inertial) mass which is needed for computing most non-gravitational accelerations. But, it does not require its own gravity field to compute gravitational accelerations.

Spherical harmonic gravity field coefficients

These coefficients may be time variable (see gravity_field_variation). The current cosine and sine coefficients can be retrieved from a Body object through its gravity field. A piece of example code on retrieving these coefficients is given below for the case of Earth:
earth_gravity_field = bodies.get( "Earth" ).gravity_field_model
cosine_coefficients = earth_gravity_field.cosine_coefficients
sine_coefficients = earth_gravity_field.cosine_coefficients
Note the above will only work if the earth_gravity_field is of the type SphericalHarmonicGravityFieldModel(), which typically means that the body has default spherical harmonic gravity field settings (see Default environment models) or that spherical harmonic gravity field settings were defined using the spherical_harmonic() function). For safety, the above could be put inside the try block of a try/except construction, wherethe except block will be entered in case the gravity field type of the Earth is not spherical harmonic.

Flight conditions

The FlightConditions class, and its derived class AtmosphericFlightConditions stores data relating to altitude, flight angles, local atmospheric properties, etc. The FlightConditions class is atypical, in the sense that a user does not provide settings for the flight conditions when creating a Body object. The reason is that the FlightConditions does not contain any ‘new’ information. Instead, it is responsible for using the existing properties of the environment and the propagation to calculate various properties related to the current state.

The reason is that FlightConditions` are related to a central body, and the object is created automatically whenever the code identifies that it is required for any of its calculations (state derivative; dependent variables, etc.). A user may also create the class themselves by using the add_flight_conditions() function. The choice between the two classes (FlightConditions and AtmosphericFlightConditions, with the latter derived from the former) is made based on the central body: if this has an atmosphere model defined, AtmosphericFlightConditions are created, if it does not, then FlightConditions are created.

Below are some examples of information that can be retrieved from the flight conditions (base class):
current_altitude = bodies.get( "Earth" ).flight_conditions.altitude
current_longitude = bodies.get( "Earth" ).flight_conditions.longiude
current_latitude = bodies.get( "Earth" ).flight_conditions.latitude
as well as its derived class that also incorporates atmospheric properties
current_airspeed = bodies.get( "Earth" ).flight_conditions.airspeed
current_freestream_density = bodies.get( "Earth" ).flight_conditions.density
current_mach_number = bodies.get( "Earth" ).flight_conditions.mach_number
The FlightConditions class also contains an object of type AerodynamicAngleCalculator, which handles the calculation of angles (latitude, longitude, flight path angle, heading angle, angle of attack, sidelip angle, bank angle) and transformations between reference frames (inertial, central-body-fixed, vertical, trajectory, aerodynamic and body-fixed frames; see this reference for details) typically used in flight dynamics. The angles and frames are listed in the tudatpy enums AerodynamicsReferenceFrameAngles and AerodynamicsReferenceFrames, respectively. Each of the angles, and the rotation between each of the frames, can be retrieved as follows (for two representative examples):
angle_calculator = bodies.get( "Earth" ).flight_conditions.aerodynamic_angle_calculator
bank_angle = angle_calculator.get_angle( environment.bank_angle )
rotation_matrix_vertical_to_body_fixed = angle_calculator.get_rotation_matrix_between_frames( environment.vertical_frame, environment.body_frame )

Aerodynamic coefficients

Aerodynamic coefficients in Tudat can be a function of a number of independent variables, such as angle of attack, Mach number, etc (see AerodynamicCoefficientsIndependentVariables for comprehensive list of options). During the propagation, the AtmosphericFlightConditions object (see above) automatically calculates the values of the independent variables, and passes the list of independent variables to an AerodynamicCoefficientInterface of the Body object (if it possesses any) to update the aerodynamic coefficients to the current state/time. The current values can be extracted from the aero_coefficient_independent_variables attribute. The current force and moment coefficients can be extracted from the coefficient interface using the current_force_coefficients and current_moment_coefficients attributes, respectively.

It may happen that a custom model influences the values of the independent variables, for instance when specifying a custom function for the angle of attack using the aerodynamic_angle_based() rotation model. If the algorithm itself depends on these angles, it may be necessary to update the aerodynamic coefficients in the guidance algorithm. One example is shown in the entry example page
# Extract Mach number from fliht conditions
mach_number = vehicle_flight_conditions.mach_number
# Compute angle attach attack according to user-defined guidance law
angle_of_attack = np.deg2rad(30 / (1 + np.exp(-2*(mach_number-9))) + 10)
# Update the variables on which the aerodynamic coefficients are based (AoA and Mach)
current_aerodynamics_independent_variables = [self.angle_of_attack, mach_number]
# Update the aerodynamic coefficients
aerodynamic_coefficient_interface.update_coefficients(
            current_aerodynamics_independent_variables, current_time)
# Extract the current force coefficients (in order: C_D, C_S, C_L)
current_force_coefficients = aerodynamic_coefficient_interface.current_force_coefficients
# Compute bank angle using guidance law requiring current_force_coefficients as input
bank_angle = ... #=f(current_force_coefficients)
In the above example, the aerodynamic coefficients are a function of angle of attack and Mach number (in that order). For an arbitrary coefficient interface, the independent variable types may be extracted using the independent_variable_names attribute.

Note that the current_force_coefficients may represent the set \(\pm[C_{D}, C_{S}, C_{L}]\) (in the aerodynamic frame) or \(\pm[C_{X}, C_{Y}, C_{Z}]\) (in the body-fixed frame). This information can be determined using the are_coefficients_in_aerodynamic_frame (for aerodynamic or body frame) and are_coefficients_in_negative_direction (for plus or minus sign).