Fixed-Wing UAV Flight Control Simulator

Aerosonde UAV — 6-DOF Flight Dynamics Modeling & Longitudinal Autopilot

A fixed-wing UAV flight dynamics and control simulator developed in MATLAB from first principles: nonlinear equations of motion -> trim -> linearization -> mode analysis -> autopilot design -> nonlinear closed-loop verification. Aircraft: Aerosonde UAV, modeled in full 6-DOF (longitudinal + lateral-directional). No Control System Toolbox or Optimization Toolbox — every numerical method (trim solver, linearization, gain selection, mode classification, simulation) runs on core MATLAB alone. I designed the trim, linearization, gain- selection, and mode-classification methodology, and independently verified every result against analytic predictions and the nonlinear plant.

---

Project Preview

Longitudinal Eigenvalue Pole Map	Lateral-Directional Eigenvalue Pole Map	Nonlinear Attitude-Hold Step Response
Nonlinear Altitude-Hold Response	Lateral Mode Responses — Dutch Roll / Roll / Spiral	Inner-Loop Root Locus (Pitch-Rate Damper)

---

V1 scope

V1 = a full 6-DOF (longitudinal + lateral-directional) flight dynamics model, plus a verified longitudinal autopilot.

	Modeled (trim, linearization, open-loop modes)	Controlled (closed-loop autopilot)
Longitudinal (pitch/altitude)	Yes — Phase 0	Yes — Phase 1 (rate + attitude + altitude hold)
Lateral-directional (roll/yaw)	Yes — Phase 2	No — explicitly out of scope for V1

Concretely: this simulator can trim, linearize, and correctly classify the dynamic modes of the full airframe (short-period, phugoid, Dutch roll, roll subsidence, spiral), and it flies a verified three-loop pitch/altitude-hold autopilot against the nonlinear longitudinal plant. It does not yet fly itself laterally — there is no aileron/ rudder control law, no heading hold, no turn coordination. That is the headline item for V2 (see Roadmap below). This distinction is stated explicitly here, in the report, and in RELEASES.md so scope is never overstated.

Project structure

matlab/       All MATLAB source (.m files): dynamics, trim, linearization,
              mode analysis, autopilot design/verification, and figure-
              generation scripts (longitudinal + lateral).
figures/      Generated plots (.png), saved automatically by the analysis
              and verification scripts.
docs/         Technical documentation of the model itself: state
              variables, nonlinear equations of motion, linearization/
              mode-analysis results (longitudinal in State_Variables.md /
              Model_Reference.md / State_Space_Model.md; lateral in
              Lateral_Dynamics.md), and the autopilot design record
              (Autopilot_Design.md).
references/   External source documentation, incl. Aerosonde_Parameters.md
              tracing every aircraft parameter back to its published source.
reports/      Standalone written engineering report (V1_Engineering_Report.md),
              evaluable without running any code.
simulink/     Simulink models (empty for now; planned for a later phase).
setup_project.m  Run once per MATLAB session to add matlab/ to the path.

Getting started

From the repo root in MATLAB:

setup_project                % adds matlab/ to the path (self-locating,
                              % works from any current folder)
analyze_modes                % Phase 0: trim -> linearize -> longitudinal open-loop modes
analyze_lateral_modes        % Phase 2: linearize -> lateral open-loop modes
                              %          (Dutch roll, roll subsidence, spiral)
verify_performance            % Phase 1: nonlinear attitude-hold verification
verify_altitude_performance   % Phase 1: nonlinear altitude-hold verification
generate_lateral_figures      % regenerates figures/lateral_*.png
generate_longitudinal_pole_map % regenerates figures/longitudinal_pole_map.png

Aircraft model

Aerosonde UAV, full 6-DOF (longitudinal + lateral-directional) dynamics. All aerodynamic and physical parameters are sourced from Beard & McLain, Small Unmanned Aircraft: Theory and Practice (Princeton University Press, 2012) — see references/Aerosonde_Parameters.md for the full parameter table, citation, and an explicit provenance note on how the lateral-directional coefficients were sourced (and a memory-recall error that was caught and corrected before use — kept on record deliberately, not swept under the rug). I designed the dynamics models, trim solver, linearization, mode classification, and all controller design specifically for this repository, and verified each stage's results before building the next stage on top of it.

Results summary

Open-loop dynamic modes

Mode	wn [rad/s]	zeta	Eigenvalue(s)	Stable?
Short-period (longitudinal)	15.42	0.444	—	Yes
Phugoid (longitudinal)	0.36	0.129	—	Yes
Dutch roll (lateral)	8.17	0.19	-1.5709 +/- 8.0202i	Yes
Roll subsidence (lateral)	—	— (tau = 0.030 s)	-32.8858	Yes
Spiral (lateral)	—	— (tau = 19.9 s)	-0.0502	Yes

All five open-loop modes are stable at this trim condition. Two findings are worth stating plainly rather than glossing over (full reasoning in reports/V1_Engineering_Report.md and docs/Lateral_Dynamics.md Section 4):

• A stable spiral mode is one valid outcome, not the textbook answer. Many real aircraft have a mildly unstable spiral mode (tolerated because it's slow enough to correct). This airframe's particular balance of dihedral effect (C_ell_beta) against weathercock stability (C_n_beta) lands on the stable side, verified independently against the standard sign criterion C_ell_beta*C_n_r - C_ell_r*C_n_beta = +0.0091 > 0.
• Dutch roll's damping (zeta=0.19) is low but normal for an unaugmented airframe with no yaw damper — not a modeling error. Real aircraft typically need active yaw-rate feedback (yaw damper, V2/future work) to raise this toward zeta ~= 0.35+.

Three-loop longitudinal autopilot

Loop	Gain(s)	Verified nonlinear result	Spec	Pass?
Rate (inner)	kq = 0.1490	short-period zeta: 0.444 -> 0.700	—	—
Attitude (middle)	kp = 24.70, ki = 36.3125	ts=1.915s, OS=9.69%, ess=0.0010 deg (5 deg step)	ts<2s, OS<10%, ess~0	PASS
Altitude (outer, P-only)	kh = 0.013297	ts=25.10s, OS=0%, ess=0.036 m (50 m step)	—	PASS

Full design story — sign-convention derivation, the structurally fixed pole in the 5-state augmented system, the anti-windup fix for integrator windup under saturation, the P-only altitude-loop proof of zero steady-state error, and the documented 20° stress-test failure mode — is in docs/Autopilot_Design.md.

Full figure set (root locus, pole maps, step responses, control effort, stress tests) lives in figures/.

Phase 0 — Trim, Linearization, and Open-Loop Longitudinal Modes (complete)

• longitudinal_dynamics.m — nonlinear 4-state equations of motion (x = [u; w; q; theta]), forced by elevator deflection.
• trim_analysis.m — Newton-Raphson solver for straight- and-level trim (alpha_trim, delta_e_trim, and the trim thrust that zeros udot).
• linearize_model.m — builds the state-space A (4x4) and B (4x1) matrices by central-difference numerical linearization about the trim point (no linmod, no toolbox).
• extract_modes.m / analyze_modes.m — eigenvalue-based separation of the short-period and phugoid modes, reporting natural frequency, damping ratio, and period for each.

Verified results: Lift = Weight = 107.91 N at trim; trim thrust = 18.44 N; Short-Period wn = 15.42 rad/s, zeta = 0.444; Phugoid wn = 0.36 rad/s, zeta = 0.129 — both modes stable, with the expected large frequency separation between them.

Phase 1 — Pitch Attitude Hold + Altitude Hold Autopilot (verified)

Nested (cascaded) PID architecture, matching the structure used by PX4/ArduPilot rather than full-state-feedback/LQR:

• Inner loop: pitch-rate damper on q (fast), designed by sweeping the rate-feedback gain and checking closed-loop eigenvalues directly (eig()-based "manual root locus" — no rlocus()). Chosen gain kq = 0.1490 drives short-period zeta from 0.444 to 0.700.
• Outer loop: pitch-attitude hold on theta (slower), commanding q to the inner loop, with integral action (ki) for zero steady-state error to a step attitude command. Chosen gains kp = 24.70, ki = 36.3125, selected via a linear-model gain search scored on simulated step response, not raw pole position (a 5-state augmented system has a structurally fixed real pole that makes worst-case eigenvalue position a misleading metric here — see docs/Autopilot_Design.md Section 2).
• Anti-windup: conditional-integration anti-windup in simulate_closed_loop_nonlinear.m, added after a 20 deg stress test revealed textbook integrator windup under sustained elevator saturation.
• Verification: closed-loop design validated by simulating the full nonlinear dynamics (ode45, not the linear model) against a step attitude command, with elevator saturation (±25 deg), checking settling time, overshoot, and steady-state error against spec.

Design targets: settling time < 2 s, overshoot < 10%, zero steady-state error to a step theta command.

Result (5 deg step, nonlinear plant): settling time 1.915 s, overshoot 9.69%, steady-state error 0.0010 deg — all three PASS, with ~zero linear-to-nonlinear gap. Demonstrated valid command envelope is approximately up to 10 deg; a 20 deg step exposes a sustained, phugoid-frequency limit cycle from repeated elevator saturation (an actuator/trim-authority limitation at fixed thrust, not a gain-tuning defect) — documented as a known limitation in docs/Autopilot_Design.md rather than papered over with a retune that would not have fixed it. See that document for the full inner-loop and outer-loop design story, including why the initial linear design's ~4% margin turned out not to need widening.

• Altitude loop: outer-most, slowest loop, commanding theta to the attitude loop. theta_cmd = kh*(h_cmd - h) — P-only, no integrator — provably exact zero steady-state altitude error at this project's fixed thrust (see docs/Autopilot_Design.md Section 3 for the full derivation and its caveats: fixed-thrust/no-wind command-tracking only, not a disturbance-rejection guarantee). Chosen gain kh = 0.013297; a wider gain scan confirmed genuinely faster linear responses exist but were rejected because they demand proportionally larger theta_cmd, which the (already-necessary) ±10° clamp overrides anyway. Result (50 m step, nonlinear plant, both theta_cmd and elevator saturation enforced): settling time 25.10 s (vs. 15.25 s linear prediction — slower, not unstable), 0% overshoot on altitude, 0.036 m steady-state error, elevator saturates for only 0.03 s at t=0 and never again. Clean, graceful degradation under the clamp — unlike the attitude loop's 20° failure, no limit cycle.

A systematic altitude-step envelope sweep (20/50/100/200 m, analogous to the attitude loop's 5/10/20° sequence) is deliberately deferred, not run for V1.

Phase 2 — Lateral-Directional Dynamics (modeling only, complete)

Adds the lateral-directional half of the 6-DOF model, decoupled from (but built on top of) the Phase 0/1 longitudinal trim — modeling and open-loop mode analysis only; no lateral autopilot is built here (see V1 scope above).

• lateral_dynamics.m — nonlinear 4-state equations of motion (x_lat = [v; p; r; phi]), forced by aileron and rudder deflection. Roll/yaw are coupled through the product of inertia Jxz; the coupled 2x2 Newton-Euler system is solved directly by matrix inversion rather than quoted from a textbook Gamma-coefficient table (see docs/Lateral_Dynamics.md Section 2 for the full derivation).
• linearize_lateral_model.m — builds A_lat (4x4) and B_lat (4x2) by the same central-difference method as linearize_model.m.
• analyze_lateral_modes.m — classifies eig(A_lat) into Dutch roll (complex pair), roll subsidence, and spiral (real eigenvalues, split by magnitude rather than assumed sign).
• generate_lateral_figures.m — pole map and initial-condition mode- response figures.

Verified results: Dutch roll wn = 8.17 rad/s, zeta = 0.19 (period 0.78 s); roll subsidence time constant 0.030 s (stable); spiral time constant 19.9 s (stable). All four lateral eigenvalues stable. I cross-checked several A_lat entries exactly against their analytic form (g*cos(theta_trim), -u_trim, tan(theta_trim)), confirming the linearization correctly recovers the equations of motion derived in Section 2. See docs/Lateral_Dynamics.md for the full derivation, linearization, and mode-analysis writeup, and the Results summary above for the honest sign-check discussion on spiral stability and Dutch roll damping.

Roadmap — V2 and beyond

Not started, listed here so scope stays explicit rather than implied:

• Lateral autopilot (heading hold): aileron/rudder control law, bank-angle hold, heading hold — the natural next step on top of the Phase 2 lateral model. This is the headline V2 item.
• Yaw damper (would directly address Dutch roll's low zeta=0.19 noted above).
• Wind/gust disturbance model and rejection.
• Waypoint navigation.
• A systematic altitude-step envelope sweep (20/50/100/200 m), deferred from Phase 1.
• Simulink model equivalents of the MATLAB-only design (currently an empty placeholder folder).

License

This project is licensed under the MIT License. See LICENSE for details.