Keyboard shortcuts

Press or to navigate between chapters

Press S or / to search in the book

Press ? to show this help

Press Esc to hide this help

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

ConstructMeaning 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
equationIntroduce 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

CheckpointEvidence to look forIf it fails
Source selectionThe generated TOML points to your edited file and intended classCorrect the active model, source roots, or overlay
Modelica compilationThe DAE command above finishesRead the parser, name-resolution, or equation diagnostic
FMI export and native buildAn FMU and its host library are producedDiagnose exporter support or the C build before trying controller changes
InitializationValid parameters, finite stationary sensors, correct gravity and actuator directionsCheck units, axes, initial conditions, and parameter bindings
MissionReadiness, arming, waypoint progress, and the required final conditionUse the console and telemetry to distinguish startup errors from flight behavior
ComparisonBaseline and changed-model logs, with the same firmware and gainsVerify 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.