fixed-wing-sim

授权协议 GPL-3.0 License
开发语言
所属分类 应用工具、 科研计算工具
软件类型 开源软件
地区 不详
投 递 者 公羊英达
操作系统 跨平台
开源组织
适用人群 未知
 软件概览

Simulation of a Fixed-Wing Unmanned Aerial Glider

Table of Content

Overview
Installation
Airframe
Reference Frames
How to Run the Code
Results
References

Overview

An example of a non-linear flight simulation for a unmanned aerial glider with a wingspan of 1.5m. The simulation is implemented with Matlab Simulink and uses FlightGear for visualization purposes.

In addition to existing Simulink examples from the Mathworks documentation, this implementation shows how to:

  1. Compute the required aerodynamic coefficient tables using Tornado an implementation of the vortex lattice method (VLM). For more information on the Tornado implementation, see also [1].
  2. Find the trimmed gliding state and deduce longitudinal and lateral linear time invariant systems (LTI) for the trimmed state according to text book definitions such as the one described in [2].
Simulation Real Flight
Visualization of the Simulink simulation with FlightGear Test flight with the real airframe
Lateral LTI Longitudinal LTI
Characteristics of the corresponding lateral LTI system Characteristics of the corresponding longitudinal LTI system

Installation and Configuration

Besides Matlab and Simulink, you need to install FlightGear and Tornado.

FlightGear Installation

After installing FlightGear, it is necessary to copy the aircraft visualization data from your Git working copy to the FlightGear data directory. Assuming you installed FlightGear 3.4.0 on Windows, just copy the content of the working copy folder FlightGear\Aircraft\ExperimentalCarrier to C:\Program Files\FlightGear 3.4.0\data\Aircraft\ExperimentalCarrier. For other versions or operating systems, proceed accordingly.

Edit the files runFlightGear.bat and runFlightGear.m in ExperimentalCarrierSimulink/utilities and adjust the FlightGear installation path to point to the correct location.

Tornado Installation

Download Tornado here and unzip it anywhere convenient, for example to C:\tornado\T135_export.

Install the Airframe and Airfoil Definitions

To install the airframe definition, copy ExperimentalCarrier.mat to T135_export/aircraft. To install the airfoil, copy also JR001.dat to T135_export/aircraft/airfoil.

Get Rid of Interp1 Warnings in Tornado

Edit T135_export\fLattice_setup2.m and replace in calls to interp1 (4 locations) cubic with pchip. This will fix the Matlab warning

Warning: INTERP1(...,'CUBIC') will change in a future release. Use INTERP1(...,'PCHIP') instead. 

Configure the Tornado Installation Directory

Assuming you have installed Tornado under C:\tornado\T135_export. Edit mainComputeCoefficients.m and mainComputeLTIs.m in ExperimentalCarrierSimulink/code and adjust the tornado_root_directory variable definition to point to the Tornado root directory, e.g.:

tornado_root_directory = 'C:\tornado\T135_export';

You can now run the simulation. Jump here to see how, or continue reading to learn about the airframe and the relevant reference frames first.

Airframe

The simulated airframe has a twin-boom fuselage and a wing with upward cranked tips. The total wing span is 1.5m and the take-off weight is 1.56kg (actual glider equiped with on-board computer and temporarily installed electric motor for testing / take-off). Center of gravity has been found to be at 92mm from the leading edge of the main wing. Via GPS measurements a gliding velocity of about 45km/h was confirmed (at roughly zero elevator deflection). The glider uses two actuators: elevator and rudder. The rudder is asymmetrically attached to the left of the two vertical stabilizers.

Below is the airframe as defined for the vortex lattice method computation with Tornado:

Wing partition layout VLM discretization
Airfoil JR001 Example pressure distribution computed by Tornado

The airfoil JR001 features a planar pressure side which simplifies the build procedure for the wing. The profile was designed to work well with low Reynold's numbers and to provide friendly stall characteristics. It wasn't designed with gliding performance in mind.

Further drawings related to the airframe can be found here and here. The Tornado definition of the airframe is here.

Reference Frames

When specifying forces, moments, or angles a body-fixed reference frame is used. The usual convention is shown in the figure below on the left. This is the convention as introduced in [2] and also as used in Matlab. The Tornado implementation [1] uses a slightly different reference frame, see below on the right.

Standard Body-fixed Reference Frame Tornado Body-fixed Reference Frame
Standard notation for forces and moments, and linear and rotational velocities in a body-fixed reference frame. The origin is located at the center of gravity (figure reproduced from [2]) Reference frame as used in the Tornado VLM implementation. The origin is located at the leading edge of the wing and the x-axis extends aft (figure reproduced from [1]).

Stability Axes

Another reference frame which is sometimes used is a set of axes for which the x-axis is parallel to the velocity vector for an equilibrium state (e.g. trimmed gliding). Such axes are called stability axes. Choosing the principal axes in this way simplifies some equations when computing longitudinal and lateral linear systems from given aerodynamic coefficients (see also [2] pp. 45). But since the solution implemented here, finds the corresponding LTI systems by linearizing a non-linear model around an equilibrium state, this is not really an advantage.

The stability axes reference frame is not used in this implementation.

How to Run the Code

The codebase and the Simulink models can be used to:

  1. Compute aerodynamic properties and coefficients using the Tornado VLM implementation.
  2. Run the non-linear flight simulation using previously computed coefficient matrices.
  3. Extract the linear time-invariant systems for the trimmed gliding state.

Compute the Aerodynamic Coefficients

Aerodynamic coefficient matrices are computed by calling the script:

mainComputeCoefficients.m

Note, the Simulink models contain pre-computed coefficient matrices, and the result of above script is added here. So, you don't need to re-compute the coefficients if you just want to to run the simulation or compute the LTIs. The resulting matrices are visualized here.

Input to compute the coefficients is:

  1. Tornado airframe definition (geometry and airfoil),
  2. center of gravity position
  3. air speed [m/s]
  4. air density [kg/m^3]
  5. range of alpha and beta values for which coefficients shall be computed (alpha denotes the angle of attack, and beta the sideslip angle).

Output is for each coefficient a 2-dimensional matrix which contains the coefficient's value for all specified [alpha, beta] configurations. Computed coefficient matrices can be found here. The coefficient naming convention is summarized below.

Aerodynamic Coefficients

Datum Coefficients
  • CX: force coefficient in body-fixed longitudinal direction.
  • CY: force coefficient in body-fixed lateral direction.
  • CZ: force coefficient perpendicular to CX, CY (body-fixed 'lift').
  • Cl: roll moment coefficient.
  • Cm: pitch moment coefficient.
  • Cn: yaw moment coefficient.
Damping Coefficients

For each force and moment coefficient, a damping coefficient is computed. These coefficients specify how each coefficient changes when the aircraft rolls (P), pitches (Q), or yaws (R). The resulting matrices are CX_P, CX_Q, CX_R, CY_P, CY_Q, CY_R, CZ_P, CZ_Q, CZ_R, Cl_P, Cl_Q, Cl_R, Cm_P, Cm_Q, Cm_R, Cn_P, Cn_Q, Cn_R.

Control Surface Deflection Coefficients

For each control surface, coefficient derivatives *_d are computed. They denote how much the respective coefficient changes when the respective surface is deflected. The resulting matrices CX_d, CY_d, CZ_d, Cl_d, Cm_d, Cn_d are 3-dimensional where the index of the third dimension denotes the control surface.

Other Aerodynamic Properties

Tornado can also be used to estimate the neutral point of the entire airframe. The neutral point is the point on the vehicle's x-axis where the aerodynamic moment Cm remains constant independently of the angle of attack. The distance between the center of gravity and the neutral point is called the stability margin. In a classic fixed-wing configuration, the center of gravity has to be placed before the neutral point (in flight direction). If the the stability margin approaches zero, the airplane becomes unstable. For the considered airframe, Tornado calculated the neutral point to lie at 49% MAC (mean aerodynamic chord), which is for the given wing geometry simply 0.49*0.25[m].

Run the Non-Linear Flight Simulation

Open the ExperimentalCarrierSimulink.prj located here. This will start FlightGear and open a couple of Simulink models. Select and run the Plant model. Switch to the FlightGear window and use v to change the view perspective.

Control the glider

To control the glider the elevator and rudder surface deflection can be adjusted with two sliders as shown below.

Compute Lateral and Longitudinal Linear Systems for the Trimmed Gliding State

To find the trimmed gliding state and calculate the lateral and longitudinal linear systems, run the script:

mainComputeLTIs.m

Results are shown here.

The resulting linear systems correspond to the ones derived in [2] (see Sect. 5.2 and Sect. 5.3), except that they are calculated for the coordinate system of the standard body-fixed reference frame and not in the reference frame of the equilibrium axes.

Results

Aerodynamic Coefficients

Datum coefficients (forces and moments).

Elevator control surface derivatives. Describe how datum coefficients change when the elevator is deflected.

Rudder control surface derivatives. Describe how datum coefficients change when the rudder is deflected.

Asymmetric values appear due to the asymmetry of the vertical stabilizers (only the stabilizer on the port side carries the rudder control surface).

For remaining coefficients, see here.

Lateral and Longitudinal Linear Systems

Output of mainComputeLTIs.m:

alpha 5.2286 deg, beta 1 deg
 
Searching for trimmed operating point of lateral dynamics...
Number of iterations: 3
Operating point specifications were successfully met.
 
Trimmed lateral state: 
============================
1) v     : velocity in y-body           (m/s): -0.0069615
2) p     : rotation velocity x-body   (deg/s): 0.0014094
3) phi   : euler roll angle             (deg): -1.5419
4) r     : rotation velocity z-body    deg/s): -0.022565
 
Other parameters in trimmed lateral state: 
================================================
Beta     : sideslip                        (deg): -0.034558
Absolute velocity (longitudinal & lateral) (m/s): 11.5901
 
Short Period Mode Properties (longitudinal): 
                                                 Pole 1      Pole 2
==============================================================================
Damping ratio                                 : 1  1
Undampted, natural frequency             (1/s): 6.8464      4.8203
Period                                     (s): Inf  Inf
Num. cycles to damp to half the amplitude     : 0  0
 
Phugoid Mode Properties (longitudinal): 
                                                 Pole 1      Pole 2
==============================================================================
Damping ratio                                 : 0.10757     0.10757
Undampted, natural frequency             (1/s): 0.43627     0.43627
Period                                     (s): 14.4863      14.4863
Num. cycles to damp to half the amplitude     : 1.0196      1.0196
 
Roll Mode Properties (lateral): 
                                                 Pole 1
==============================================================================
Damping ratio                                 : 1
Undampted, natural frequency             (1/s): 2.1251
Period                                     (s): Inf
Num. cycles to damp to half the amplitude     : 0
 
Spiral Mode Properties (lateral): 
                                                 Pole 1
==============================================================================
Damping ratio                                 : 1
Undampted, natural frequency             (1/s): 0.1132
Period                                     (s): Inf
Num. cycles to damp to half the amplitude     : 0
 
Dutch Roll Mode Properties (lateral): 
                                                 Pole 1      Pole 2
==============================================================================
Damping ratio                                 : 0.018588    0.018588
Undampted, natural frequency             (1/s): 2.3622      2.3622
Period                                     (s): 2.6603      2.6603
Num. cycles to damp to half the amplitude     : 5.9339      5.9339
 
Number of unobservable states (longitudinal)  : 0
Number of unobservable states (lateral)       : 0
 
Number of uncontrollable states (longitudinal): 0
Number of uncontrollable states (lateral)     : 0

Longitudinal system:

ltiLongitudinal =
 
  a = 

         0.004876         1.101        -1.056        -9.759
           -1.009        -7.448         2.742       -0.6098
           0.0347       -0.3791        -4.318             0
                0             0             1             0
 
  b = 
                
                0
            10.21
           -9.064
                0

Lateral system:

ltiLateral =
 
  a = 

          -0.7838        0.6121         9.755        -11.15
          -0.7469        -1.227             0        0.1967
                0             1    -3.794e-05       0.06246
           0.3551      -0.02095             0       -0.3154
 
  b = 

           -2.526
                0
                0
           -2.519

Short period and phugoid modes (poles) of the longitudinal system. Dirac and step response to elevator deflection. A negative pitch moment curve for increasing alpha values (angle of attack) confirms the glider's positive stability margin.

Spiral, rolling, and dutch role modes (poles) of the lateral system. Dirac and step response to rudder deflection.

References

[1] Melin, Tomas. Tornado, a vortex lattice MATLAB implementation for Linear Aerodynamic Wing applications, Masters thesis, Royal Institute of Technology (KTH), Sweden, December 2000.
[2] Caughey, David A. Introduction to Aircraft Stability and Control, Course Notes for Mechanical & Aerospace Engineering, Cornell University, New York, USA, 2011.

 相关资料
  • Wing意为翅膀,该框架将提供一套完整的开发方法论的支持。主要为实现TDD,在项目实践中不断提炼,才产生了该框架。

  • WingIDE是个相当优秀的 IDE;其编辑器包括大量语言的语法标签高亮显示,虽然它只是个面向 Python 的工具。源代码浏览器对浏览项目或模块非常实用(表现在可导航源代码和文档行摘要中)。虽然没有监视器,但调试器设计得很好。编辑器有优秀的命令自动完成和函数跳转列表,但是没有代码合并。面向项目风格的 IDE 对于大型产品非常有用(在这方面,除了 Komodo 以外,它是大多免费 IDE 中较好的

  • 前端生态系统有着各种各样的开源框架,这些框架都能让开发者的工作更简单一些。一般情况下,每一个人都有自己不同的需求,有的人看重设计美学,而其他一些人看重的是速度和可用度。 Wing是一个新的框架,也是极简主义爱好者的好朋友。它提供了默认格栅布局、自定义元素以及各种组件,而且它的体积只有5KB。 Wing最好的地方,在于它针对默认元素提供了自动样式(auto-styling)功能。当你将Wing CS

  • fixed 方法 把字符串显示为打字机字体。 语法: stringObject.fixed(); 示例: var Str = "graybobo", s = Str.fixed(); console.log( s ); 结果: >>> <tt>graybobo</tt>

  • 描述 (Description) 此方法使字符串以固定间距字体显示,就像它在标记中一样。 语法 (Syntax) 其语法如下 - string.fixed( ) 返回值 (Return Value) 返回带有《tt》标记的字符串。 例子 (Example) 请尝试以下示例。 <html> <head> <title>JavaScript String fixed() Metho

  • Wing FTP Server 是一个专业的跨平台FTP服务器端,它拥有不错的速度、可靠性和一个友好的配置界面。它 除了能提供FTP的基本服务功能以外,还能提供管理员终端、任务计划、基于Web的管理端,基于Web的客户端和Lua脚本扩展等,它还支持虚拟文件夹、 上传下载比率分配、磁盘容量分配,ODBC/Mysql存储账户等特性,支持 Windows, Linux, Mac OS和Solaris。