Trimming and Linearization of the HL-20 Airframe

This example uses:

This is Part 1 of a five-part example series on design and tuning of the flight control system for the HL-20 vehicle. This part deals with trimming and linearization of the airframe.

HL-20 Model

The HL-20 model is adapted from the model described in NASA HL-20 Lifting Body Airframe (Aerospace Blockset). This is a 6-DOF model of the vehicle during the final descent and landing phase of the flight. No thrust is used during this phase and the airframe is gliding to the landing strip.

open_system('csthl20_trim')

This version of the model includes the equations of motion (EOM), the force and moment calculation from the aerodynamic tables, the environment model, and the "Controls Selector" block which maps aileron, elevator, and rudder demands to deflections of the six control surfaces.

Batch Trimming

Trimming consists of calculating aileron, elevator, and rudder deflections that zero out the forces and moments on the airframe, or equivalently, keep the body velocities ub,vb,wb and angular rates p,q,r steady. Because thrust is not used during descent, one degree-of-freedom is lost and the trim condition must be relaxed to let ub to vary. The trim values of the deflections da,de,dr depend on the airframe orientation relative to the wind. This orientation is characterized by the angle-of-attack (AoA) alpha and the sideslip angle (AoS) beta.

With the operspec and findop functions, you can efficiently compute the trim deflections over a grid of (alpha,beta) values covering the operating range of the vehicle. Here we trim the model for 8 values of alpha ranging from -10 to 25 degrees, and 5 values of beta ranging from -10 to +10 degrees. The nominal altitude and speed are set to 10,000 feet and Mach 0.6.

d2r = pi/180;            % degrees to radians
m2ft = 3.28084;          % meter to feet
Altitude = 10000/m2ft;   % Nominal altitude
Mach = 0.6;              % Nominal Mach
alpha_vec = -10:5:25;	 % Alpha Range
beta_vec = -10:5:10;     % Beta Range
[alpha,beta] = ndgrid(alpha_vec,beta_vec);   % (Alpha,Beta) grid

Use operspec to create an array of operating point specifications.

opspec = operspec('csthl20_trim',size(alpha));

opspec(1)

ans = 


 Operating point specification for the Model csthl20_trim.
 (Time-Varying Components Evaluated at time t=0)

States: 
----------
     x         Known    SteadyState     Min         Max        dxMin       dxMax   
___________ ___________ ___________ ___________ ___________ ___________ ___________
                                                                                   
(1.) csthl20_trim/HL20 Airframe/6DOF (Euler Angles)/Calculate DCM & Euler Angles/phi theta psi
     0         false       true        -Inf         Inf        -Inf         Inf    
 -0.19945      false       true        -Inf         Inf        -Inf         Inf    
     0         false       true        -Inf         Inf        -Inf         Inf    
(2.) csthl20_trim/HL20 Airframe/6DOF (Euler Angles)/p,q,r 
     0         false       true        -Inf         Inf        -Inf         Inf    
     0         false       true        -Inf         Inf        -Inf         Inf    
     0         false       true        -Inf         Inf        -Inf         Inf    
(3.) csthl20_trim/HL20 Airframe/6DOF (Euler Angles)/ub,vb,wb
  202.67       false       true        -Inf         Inf        -Inf         Inf    
     0         false       true        -Inf         Inf        -Inf         Inf    
  23.257       false       true        -Inf         Inf        -Inf         Inf    
(4.) csthl20_trim/HL20 Airframe/6DOF (Euler Angles)/xe,ye,ze
-12071.9115    false       true        -Inf         Inf        -Inf         Inf    
     0         false       true        -Inf         Inf        -Inf         Inf    
-3047.9999     false       true        -Inf         Inf        -Inf         Inf    

Inputs: 
----------
  u   Known  Min   Max 
_____ _____ _____ _____
                       
(1.) csthl20_trim/da
  0   false -Inf   Inf 
(2.) csthl20_trim/de
  0   false -Inf   Inf 
(3.) csthl20_trim/dr
  0   false -Inf   Inf 

Outputs: 
----------
  y   Known  Min   Max 
_____ _____ _____ _____
                       
(1.) csthl20_trim/p;q;r (1-3)
  0   false -Inf   Inf 
  0   false -Inf   Inf 
  0   false -Inf   Inf 
(2.) csthl20_trim/phi;theta;psi (4-6)
  0   false -Inf   Inf 
  0   false -Inf   Inf 
  0   false -Inf   Inf 
(3.) csthl20_trim/alpha (7)
  0   false -Inf   Inf 
(4.) csthl20_trim/beta (8)
  0   false -Inf   Inf 
(5.) csthl20_trim/Mach (9)
  0   false -Inf   Inf 
(6.) csthl20_trim/Ax,Ay,Az (10-12)
  0   false -Inf   Inf 
  0   false -Inf   Inf 
  0   false -Inf   Inf

Specify the equilibrium conditions for each orientation of the airframe. To do this:

Specify the orientation by fixing the outputs alpha and beta to their desired values.
Specify the airframe speed by fixing the Mach output to 0.6.
Mark the angular rates p,q,r as steady.
Mark the velocities vb and wb as steady.

for ct=1:40
   % Specify alpha angle
   opspec(ct).Outputs(3).y = alpha(ct);
   opspec(ct).Outputs(3).Known = true;
   % Specify beta angle
   opspec(ct).Outputs(4).y = beta(ct);
   opspec(ct).Outputs(4).Known = true;
   % Specify Mach speed
   opspec(ct).Outputs(5).y = Mach;
   opspec(ct).Outputs(5).Known = true;
   % Mark p,q,r as steady
   opspec(ct).States(2).SteadyState = true(3,1);
   % Mark vb,wb as steady
   opspec(ct).States(3).SteadyState = [false;true;true];
   % (phi,theta,psi) and (Xe,Ye,Ze) are not steady
   opspec(ct).States(1).SteadyState = false(3,1);
   opspec(ct).States(4).SteadyState = false(3,1);
end

To fully characterize the trim condition, also

Set p=0 to prevent rolling.
Set the roll/pitch/yaw angles (phi,theta,psi) to (0,alpha,beta) to align the wind and earth frames.
Specify the airframe position (Xe,Ye,Ze) as (0,0,-Altitude).

for ct=1:40
   % Set (phi,theta,psi) to (0,alpha,beta)
   opspec(ct).States(1).x = [0 ; alpha(ct)*d2r ; beta(ct)*d2r];
   opspec(ct).States(1).Known = true(3,1);
   % Set p=0 (no rolling)
   opspec(ct).States(2).x(1) = 0;
   opspec(ct).States(2).Known(1) = true;
   % Set (Xe,Ye,Ze) to (0,0,-Altitude)
   opspec(ct).States(4).x = [0 ; 0 ; -Altitude];
   opspec(ct).States(4).Known = true(3,1);
end

Now use findop to compute the trim conditions for all 40 (alpha,beta) combinations in one go. This batch mode approach involves a single compilation of the model. FINDOP uses optimization to solve the nonlinear equations characterizing each equilibrium. Here we use the "SQP" algorithm for this task.

% Set options for FINDOP solver
TrimOptions = findopOptions;
TrimOptions.OptimizationOptions.Algorithm = 'sqp';
TrimOptions.DisplayReport = 'off';

% Trim model
[ops,rps] = findop('csthl20_trim',opspec,TrimOptions);

This returns 8-by-5 arrays OPS (operating conditions) and RPS (optimization reports). You can use RPS to verify that each trim condition was successfully calculated. Results for the first (alpha,beta) pair are shown below.

[alpha(1) beta(1)]

ans =

   -10   -10

ops(1)

ans = 


 Operating point for the Model csthl20_trim.
 (Time-Varying Components Evaluated at time t=0)

States: 
----------
    x     
__________
          
(1.) csthl20_trim/HL20 Airframe/6DOF (Euler Angles)/Calculate DCM & Euler Angles/phi theta psi
    0     
 -0.17453 
 -0.17453 
(2.) csthl20_trim/HL20 Airframe/6DOF (Euler Angles)/p,q,r 
    0     
 -0.15825 
 0.008004 
(3.) csthl20_trim/HL20 Airframe/6DOF (Euler Angles)/ub,vb,wb
 191.0911 
 -34.2143 
 -33.6945 
(4.) csthl20_trim/HL20 Airframe/6DOF (Euler Angles)/xe,ye,ze
    0     
    0     
-3047.9999

Inputs: 
----------
   u    
________
        
(1.) csthl20_trim/da
-23.9841
(2.) csthl20_trim/de
-6.4896 
(3.) csthl20_trim/dr
 4.0858

rps(1).TerminationString

ans =

    'Operating point specifications were successfully met.'

Batch Linearization

The gains of the flight control system are typically scheduled as a function of alpha and beta, see Part 2 (Angular Rate Control in the HL-20 Autopilot) for more details. To tune these gains, you need linearized models of the HL-20 airframe at the 40 trim conditions. Use linearize to compute these models from the trim operating conditions ops.

% Linearize airframe dynamics at each trim condition
G = linearize('csthl20_trim','csthl20_trim/HL20 Airframe',ops);

size(G)

8x5 array of state-space models.
Each model has 34 outputs, 9 inputs, and 12 states.

The linear equivalent of the "Controls Selector" block depends on the amount of elevator deflection and should be computed for qbar_inv=1 (nominal dynamic pressure at Mach=0.6). For convenience, also linearize this block at the 40 trim conditions.

CS = linearize('csthl20_trim','csthl20_trim/Controls Selector',ops);

% Zero out a/b and qbar_inv channels
CS = [CS(:,1:3) zeros(6,2)];

Linear Model Simplification

The linearized airframe models have 12 states:

xG = G.StateName

xG =

  12x1 cell array

    {'phi theta psi(1)'}
    {'phi theta psi(2)'}
    {'phi theta psi(3)'}
    {'p,q,r (1)'       }
    {'p,q,r (2)'       }
    {'p,q,r (3)'       }
    {'ub,vb,wb(1)'     }
    {'ub,vb,wb(2)'     }
    {'ub,vb,wb(3)'     }
    {'xe,ye,ze(1)'     }
    {'xe,ye,ze(2)'     }
    {'xe,ye,ze(3)'     }

Some states are not under the authority of the roll/pitch/yaw autopilot and other states contribute little to the design of this autopilot. For control purposes, the most important states are the roll angle phi, the body velocities ub,vb,wb, and the angular rates p,q,r. Accordingly, use modred to obtain a 7th-order model that only retains these states.

G7 = G;
xKeep = {...
   'phi theta psi(1)'
   'ub,vb,wb(1)'
   'ub,vb,wb(2)'
   'ub,vb,wb(3)'
   'p,q,r(1)'
   'p,q,r(2)'
   'p,q,r(3)'};
[~,xElim] = setdiff(xG,xKeep);
for ct=1:40
   G7(:,:,ct) = modred(G(:,:,ct),xElim,'truncate');
end

With these linearized models in hand, you can move to the task of tuning and scheduling the flight control system gains. See Angular Rate Control in the HL-20 Autopilot for Part 2 of this example.

Trimming and Linearization of the HL-20 Airframe

HL-20 Model

Batch Trimming

Batch Linearization

Linear Model Simplification

Related Topics