Read the Modelica model
Know which source your run uses
The first mission uses modelica/FastDyn/Copter.mo, shown below, with the
library in third_party/common/modelica_models. The QAV-R models extend the
same wrapper. Check the active entry’s model_file, model, source_roots,
output, and [FMU].compiler in your generated TOML before an experiment.
Keep edited models in your source directory. Study scripts regenerate their
files under out/; edits to those generated copies will not persist.
A little Modelica before the full model
| Construct | Meaning in this vehicle |
|---|---|
within FastDyn; | Put the class in the FastDyn package |
model Qavr ... end Qavr; | Define the named model FastDyn.Qavr |
parameter Real mass = 0.5; | A real-valued quantity fixed during this experiment, in kg by convention |
Real force_b[3]; | A vector with three real components; Modelica indices start at 1 |
Real inertia[3,3]; | A 3-by-3 matrix; {...} constructs a vector and [...] can construct a matrix |
extends Copter(mass = bare_mass); | Reuse the parent model and change a parameter binding |
equation | Introduce relationships the compiler must satisfy, rather than a sequence of assignments |
der(omega) | The time derivative of motor speed; this introduces dynamics |
Changing tau_up adjusts the motor law already present. Replacing a motor law
or adding an external-force equation changes the physics represented by the
model. In both cases, keep the exported actuator and sensor contract intact
unless you also intend to change FastDyn’s physics backend.
1. Open the vehicle wrapper
The source below is included directly from modelica/FastDyn/Copter.mo, so it
matches the file in this checkout. It connects a reusable plant to the PWM and
sensor interface expected by the ArduPilot drivers.
Show the complete FastDyn.Copter Modelica model
within FastDyn;
model Copter
parameter Real mass = 2.5644001 "Gazebo gs_drone equivalent mass [kg]";
parameter Real inertia[3,3] = diagonal({0.02601237985, 0.02590943825, 0.045571756801})
"Inertia about the CG in body FLU [kg*m^2]";
parameter Real Ct = 8.54858e-6 "Thrust coefficient [N/(rad/s)^2]";
parameter Real Cm = 0.016 "Rotor torque/thrust ratio [m]";
parameter Real arm_length = 0.215 "Arm length [m] (430 mm motor-to-motor wheelbase)";
parameter Real Cl_p = -0.2 "Rolling moment coefficient per roll rate";
parameter Real Cm_q = -0.2 "Pitching moment coefficient per pitch rate";
parameter Real Cn_r = -0.1 "Yawing moment coefficient per yaw rate";
parameter Real motor_thrust_scale[4] = {1, 1, 1, 1}
"Per-motor thrust effectiveness multiplier for fault/robustness studies";
parameter Real tau_up = 0.0125 "Motor spin-up time constant [s]";
parameter Real tau_down = 0.025 "Motor spin-down time constant [s]";
parameter Real body_area = 0.1 "Reference area for rate damping [m^2]";
parameter Real drag_area[3] = {0.06, 0.08, 0.12} "Body drag areas [m^2]";
parameter Real linear_drag[3] = {0.12, 0.12, 0.18} "Body linear drag [N*s/m]";
parameter Real leg_x = 0.17 "Ground contact X offset [m]";
parameter Real leg_y = 0.17 "Ground contact Y offset [m]";
parameter Real leg_z = -0.10 "Ground contact Z offset in body FLU [m]";
parameter Real ground_k = 3000 "Contact stiffness [N/m]";
parameter Real ground_c = 150 "Contact normal damping [N*s/m]";
parameter Real ground_tangent_c = 25 "Contact tangential damping [N*s/m]";
FastDyn.QuadrotorWithExternalWrench plant(
external_force_b = {0, 0, 0},
external_moment_b = {0, 0, 0},
ground_z = 0.0,
vehicle_mass = mass,
J = inertia,
gravity = 9.8,
mag_world_enu = earth_mag_enu,
Ct = Ct,
Cm = Cm,
arm_length = arm_length,
Cl_p = Cl_p,
Cm_q = Cm_q,
Cn_r = Cn_r,
motor_thrust_scale = motor_thrust_scale,
tau_up = tau_up,
tau_down = tau_down,
S = body_area,
CdA = drag_area,
linear_drag = linear_drag,
leg_x = leg_x,
leg_y = leg_y,
leg_z = leg_z,
ground_k = ground_k,
ground_c = ground_c,
ground_tangent_c = ground_tangent_c);
parameter Real pwm_min = 1100.0 "Minimum motor PWM used by the Gazebo gs_drone ArduPilot control block";
parameter Real pwm_max = 1900.0 "Maximum motor PWM used by the Gazebo gs_drone ArduPilot control block";
parameter Real omega_min = 0.0 "Motor speed at minimum PWM [rad/s]";
parameter Real omega_max = 1300.0 "Aerodynamic motor speed at maximum PWM [rad/s]";
parameter Real lat0 = 40.414929 "Reference latitude [deg]";
parameter Real lon0 = -86.932387 "Reference longitude [deg]";
parameter Real ground_alt_wgs84 = 149.0 "WGS84 ellipsoid altitude of the local ground collision plane [m]";
parameter Real accel_bias[3] = {0, 0, 0} "Accelerometer bias [m/s^2]";
parameter Real gyro_bias[3] = {0, 0, 0} "Gyroscope bias [rad/s]";
parameter Real mag_bias[3] = {0, 0, 0} "Magnetometer bias [Gauss]";
parameter Real mag_motor_bias[3] = {0, 0, 0}
"Magnetometer bias at full normalized mean motor load [Gauss]";
parameter Real earth_mag_enu[3] = {0.21, 0, -0.45}
"Earth magnetic field in world ENU axes [Gauss]";
parameter Real current_idle_a = 0.0
"Empirical current proxy intercept [A]";
parameter Real current_per_motor_load_a = 0.0
"Empirical current proxy slope versus normalized mean motor load [A]";
parameter Real mag_current_slope[3] = {0, 0, 0}
"Body-FRD magnetometer bias slope versus estimated current [Gauss/A]";
parameter Real gps_bias[3] = {0, 0, 0} "GPS bias N/E/altitude [m]";
parameter Real baro_alt_bias = 0.0 "Barometer relative altitude bias [m]";
parameter Real earth_radius_m = 6378137.0 "Spherical Earth radius used for local geodetic conversion [m]";
parameter Real pi = 3.141592653589793;
input Real pwm[4](start = {1000, 1000, 1000, 1000}) "Motor PWM commands";
output Real accel[3] "Body FRD accelerometer [m/s^2]";
output Real gyro[3] "Body FRD gyroscope [rad/s]";
output Real mag[3] "Body FRD magnetometer [Gauss]";
output Real gps[3] "GPS latitude, longitude, altitude";
output Real vel_ned[3] "GPS velocity NED [m/s]";
output Real yaw_deg "Yaw [deg]";
output Real baro_altitude_m "Barometer relative altitude [m]";
output Real baro_pressure_pa "Barometer pressure [Pa]";
output Real baro_temperature_c "Barometer temperature [degC]";
output Real baro_climb_rate_mps "Barometer climb rate [m/s]";
output Real motor_cmd[4] "Motor commands after PWM scaling [rad/s]";
output Real estimated_current_a "Empirical propulsion-current proxy [A]";
protected
Real pwm_span;
Real pwm_norm[4];
Real gps_lat_lon[2];
Real geodetic_origin[3] "Reference latitude, longitude, and Earth radius";
Real yaw_rad;
Real mean_motor_load;
equation
pwm_span = pwm_max - pwm_min;
for i in 1:4 loop
pwm_norm[i] = min(1.0, max(0.0, (pwm[i] - pwm_min) / pwm_span));
motor_cmd[i] = omega_min + (omega_max - omega_min) * pwm_norm[i];
plant.omega_cmd[i] = motor_cmd[i];
end for;
accel = {plant.accel[1], -plant.accel[2], -plant.accel[3]} + accel_bias;
gyro = {plant.gyro[1], -plant.gyro[2], -plant.gyro[3]} + gyro_bias;
mean_motor_load = sum(plant.omega_m .* plant.omega_m) /
(4.0 * max(1.0, omega_max * omega_max));
estimated_current_a = max(0.0,
current_idle_a + current_per_motor_load_a * mean_motor_load);
mag = {plant.mag[1], -plant.mag[2], -plant.mag[3]} + mag_bias +
mag_motor_bias * mean_motor_load + mag_current_slope * estimated_current_a;
// Avoid collisions between the caller parameters and the function locals
// during function projection in the pinned Rumoca compiler.
geodetic_origin = {lat0, lon0, earth_radius_m};
gps_lat_lon = Geodesy.localNorthEastToLatLon(
geodetic_origin[1],
geodetic_origin[2],
plant.p[1] + gps_bias[1],
-plant.p[2] + gps_bias[2],
geodetic_origin[3]);
gps[1] = gps_lat_lon[1];
gps[2] = gps_lat_lon[2];
gps[3] = ground_alt_wgs84 + plant.p[3] + gps_bias[3];
vel_ned[1] = plant.v_w[1];
vel_ned[2] = -plant.v_w[2];
vel_ned[3] = -plant.v_w[3];
yaw_rad = atan2(2.0 * (plant.q[1] * plant.q[4] + plant.q[2] * plant.q[3]),
1.0 - 2.0 * (plant.q[3] * plant.q[3] + plant.q[4] * plant.q[4]));
yaw_deg = -yaw_rad * 180.0 / pi;
baro_altitude_m = plant.p[3] + baro_alt_bias;
baro_temperature_c = 15.0 - 0.0065 * (ground_alt_wgs84 + baro_altitude_m);
baro_pressure_pa = 101325.0 * (1.0 - 2.25577e-5 * (ground_alt_wgs84 + baro_altitude_m)) ^ 5.25588;
baro_climb_rate_mps = plant.v_w[3];
end Copter;
2. Follow the physics
The wrapper’s plant adapts the library’s
Vehicles.Templates.QuadrotorPlant.
The local QuadrotorWithExternalWrench adds external force and moment inputs
for the load exercise; those inputs default to zero.
Its actuator inputs are four motor speeds. Four motor states lag their commands,
and each motor produces thrust proportional to the square of its speed. The motor
moment map turns those thrusts into roll, pitch, and yaw moments. Body drag and
ground contacts add forces and moments before rigid-body integration.
// The motor is a dynamic state, not an instantaneous actuator.
der(omega) = tau_inv * omega_error;
thrust = Ct * omega * omega;
// Four thrusts produce a body moment through the motor geometry.
M_rotor = motor_moment_map * F_m;
Find these equations in the library file. Change tau_up and tau_down to
model a different motor response; change the equations to model another
actuator law. A new law should be checked against measured actuator data.
3. Keep arrays as arrays
The inertia tensor is a Real[3,3]; drag areas, forces, and sensor vectors are
Real[3]. The wrapper supplies J = inertia to the rigid body. Keeping this
as a matrix makes the Modelica equations match the physical notation.
The QAV-R uses its estimated inertia tensor. The payload experiment keeps the vehicle mass and inertia fixed and adds an external force and its moment instead, so these are two separate modeling experiments.
4. Inspect the compiler output
In your chosen environment, run:
rumoca compile modelica/FastDyn/Qavr.mo --model FastDyn.Qavr \
--source-root modelica --source-root third_party/common/modelica_models \
--emit dae-json --output out/qavr-dae.json
Expected output includes:
wrote Dae IR (json) to out/qavr-dae.json
This checks name resolution and lowering. It does not by itself establish that the FMI exporter supports every feature in the lowered model, or that the resulting plant flies correctly. Those are separate build and simulation checkpoints.
5. Check the frame convention
The plant’s body axes are forward, left, up. The firmware interface uses forward, right, down, so the wrapper changes the signs of Y and Z sensor components. This tutorial chooses world X=north, Y=west, Z=up and converts GPS velocity to NED. Its magnetic-field vector is configured consistently with that choice; a frame label alone cannot determine the numerical rotation.
Diagnose one stage at a time
| Checkpoint | Evidence to look for | If it fails |
|---|---|---|
| Source selection | The generated TOML points to your edited file and intended class | Correct the active model, source roots, or overlay |
| Modelica compilation | The DAE command above finishes | Read the parser, name-resolution, or equation diagnostic |
| FMI export and native build | An FMU and its host library are produced | Diagnose exporter support or the C build before trying controller changes |
| Initialization | Valid parameters, finite stationary sensors, correct gravity and actuator directions | Check units, axes, initial conditions, and parameter bindings |
| Mission | Readiness, arming, waypoint progress, and the required final condition | Use the console and telemetry to distinguish startup errors from flight behavior |
| Comparison | Baseline and changed-model logs, with the same firmware and gains | Verify that a rebuilt FMU was loaded and that no unrelated settings changed |
When moving between compiler generations, choose a fresh FMU output directory
or move your old generated output aside. Rumoca can refuse to overwrite an
archive from an older format (not a recognized previous product); that is an
artifact conflict, not a failure of your new Modelica equations.
A changed source file or a successful DAE export alone is not the final result. For your own physics experiment, retain the source, run TOML, compiler revision, FMU, and telemetry together so you or someone else can reproduce the comparison.
Next, select the QAV-R model and see which physical properties change when you move to a smaller airframe.