东南大学Matlab作业3

1 Matlab Worksheet 3 Part A 1 Using function conv m m to make convolution between the following to functions x and h x 3 11 7 0 1 7 5 0 2 h 11 9 0 7 3 2 0 1 nx 2 6 nh 0 7 Plot the functions and convolution results Answer home work 3 A 1 m x 3 11 7 0 1 7 5 0 2 h 11 9 0 7 3 2 0 1 nx 2 6 nh 0 7 y ny conv m x nx h nh subplot 3 1 1 stem nx x ylabel x n axis 6 10 20 20 subplot 3 1 2 stem nh h ylabel h n axis 6 10 20 20 subplot 3 1 3 stem ny y xlabel n ylabel y n axis 6 10 50 50 conv m m function y ny conv m x nx h nh convolution function conv m m nyb nx 1 nh 1 nye nx length x nh length h ny nyb nye y conv x h end 6 4 20246810 20 0 20 x n 6 4 20246810 20 0 20 h n 6 4 20246810 50 0 50 n y n 2 2 Plot the frequency response over for the following transfer function by letting 0 where is the frequency rad sample with appropriate labels and title j ez 9 06 1 2 zz z zH Answer home work 3 A 2 m delta 0 01 Omega 0 delta pi z exp j Omega H z z 2 1 6 z 0 9 subplot 2 1 1 plot Omega abs H m LineWidth 2 xlabel 0 Omega pi ylabel H Omega axis 0 pi 0 max abs H subplot 2 1 2 plot Omega atan2 imag H real H k LineWidth 2 xlabel 0 Omega pi ylabel pi Phi H pi title frequency response axis 0 pi pi pi 00 511 522 53 0 5 10 15 0 H 00 511 522 53 2 0 2 0 H frequency response 3 Use fft to analyse following signal by plotting the original signal and its spectrum 1024 4672sin 1024 1372sin 1024 322sin nnn nx 3 Answer home work 3 A 3 m Examine signal components sampling rate 1024Hz N 1024 dt 1 N t 0 dt 1 dt define the signal x sin 2 pi 32 t sin 2 pi 137 t sin 2 pi 467 t subplot 2 1 1 plot t x axis 0 1 1 1 xlabel nT seconds ylabel x n X fft x df 1 f 0 df N 1 subplot 2 1 2 plot f abs X axis 0 N 2 0 120 xlabel k Hz ylabel X k 00 10 20 30 40 50 60 70 80 91 1 0 5 0 0 5 1 nT seconds x n 050100150200250300350400450500 0 50 100 k Hz X k 4 Use the fast Fourier transform function fft to analyse following signal Plot the original signal and the magnitude of its spectrum linearly and logarithmically Apply Hamming window to reduce the leakage 1024 7 4672sin 1024 4 1372sin 1024 5 322sin nnn nx The hamming window can be coded in Matlab as for n 1 N hamming n 0 54 0 46 cos 2 n N 1 pi N end 4 where N is the data length in the FFT Answer home work 3 A 4 m Windowing for reducing leakage sampling rate 1024Hz N 1024 dt 1 N t 0 dt 1 dt define the signal x sin 2 pi 32 5 t sin 2 pi 137 4 t sin 2 pi 467 7 t Apply rectangular window subplot 3 2 1 plot t x axis 0 1 5 5 xlabel nT seconds ylabel x n rectang win X fft x df 1 f 0 df N 1 subplot 3 2 2 plot f 20 log10 abs X axis 0 N 2 30 60 xlabel k Hz ylabel X k db Apply triangular window for n 1 N w n 1 abs 2 n N 1 N end x1 x w subplot 3 2 3 plot t x1 axis 0 1 5 5 xlabel nT seconds ylabel x n triang win X fft x1 df 1 f 0 df N 1 subplot 3 2 4 plot f 20 log10 abs X axis 0 N 2 30 60 xlabel k Hz ylabel X k db Apply Hamming window for n 1 N w n 0 54 0 46 cos 2 n N 1 pi N end x1 x w subplot 3 2 5 plot t x1 axis 0 1 5 5 xlabel nT seconds ylabel x n Hamming win X fft x1 df 1 f 0 df N 1 subplot 3 2 6 plot f 20 log10 abs X 5 axis 0 N 2 30 60 xlabel k Hz ylabel X k db 00 51 5 0 5 nT seconds x n rectang win 0200400 20 0 20 40 60 k Hz X k db 00 51 5 0 5 nT seconds x n triang win 0200400 20 0 20 40 60 k Hz X k db 00 51 5 0 5 nT seconds x n Hamming win 0200400 20 0 20 40 60 k Hz X k db 6 Part B Simulation Using SIMULINK INTRODUCTION The objective of this laboratory is to learn about various properties of signals and systems by doing simulations in SIMULINK PART 1 BASICS OF SIMULINK 1 1 UNDERSTANDING SIMPLE WAVEFORMS AND INTEGRATION Create a pulse of height 2 units from time 0 to 4 seconds by subtracting two unit steps and adding a gain Connect this pulse to an integrator with a gain of 0 5 and a zero initial condition Connect oscilloscopes to show the pulse and the output of the integrator You may wish to name your simulation block diagram to do so use the save as feature under edit Your block diagram should be similar as below u Before simulating you need to pull down the simulation header and double click on parameters Unless told otherwise always assume that you can use the ode 45 integration algorithm shown in this window and the other default parameters Typically you will only alter the start and stop times but for this first simulation you can use the default values of 0 and 10 seconds Double click on the oscilloscopes to get the windows in which the traces will appear pull down the simulate menu and click on run Plot below the integrator input and output waveforms 7 Repeat the experiment but with an initial condition of 2 on the integrator Again draw the results 1 2 FIRST ORDER SYSTEM A single time constant may be simulated using the transfer function block in which you enter the coefficients of the numerator and denominator polynomials Set up the configuration in a new SIMULINK window to realise the transfer function 1 s 2 with the input unit step and an oscilloscope connected to the output of the transfer function block Plot the block diagram in the space below Simulate the system for 5 seconds and plot the response 8 1 3 SECOND ORDER SYSTEM For the second order all pole transfer function you will recall that if a time scale of 2 2 00 22 0 ss 0t is used for plotting the step response the response shape will only be affected by changes in the damping ratio This can also be shown if we normalise the transfer function by replacing s 0 by sn to give To study the effect of varying on the step response we will therefore use the 12 1 2 nn ss transfer function 1 s2 2 s 1 Set up the following configuration for a simulation study i A unit step input ii Connect the unit step to a transfer function of 1 s2 2 s 1 with 1 0 iii Take a summing block and connect the input step to a ve input and the transfer function output to a ve input iv Connect the output of this summer to a square function This is obtained by using f u in the nonlinear blocks v Connect the output of this squarer to an integrator vi Connect two oscilloscopes to the circuit one to the transfer function output and the other to the output of the integrator vii Also connect simout blocks to the same signals as the oscilloscopes Rename the one connected to the integrator output simout1 Simulate the system with the values of listed in the table below and fill in the other figures from the simulation results To get an accurate value at the end of each run you can type simout1 in the MATLAB window You can also measure the overshoot by making use of the maximum command simply type max simout in the MATLAB window Table 1 2 01 00 90 70 50 30 10 01 overshoot 2 10761 2500 1 17781 0571 1 0000 1 1285 2 1431 3 9651 home work 3 B 1 3 slx 9 when 1 0 PART 2 SIMULATION OF AIRCRAFT PITCH ANGLE AND ALTITUTE pitch 倾斜 倾斜 The purpose IS to use SIMULINK to simulate a much simplified model of the longitudinal motion of a fighter aircraft The angle of attack is the angle between the direction a plane is pointing and the direction in which it actually moves through the air For a plane flying at approximately constant altitude this is equivalent to the pitch angle as illustrated in Fig 2 1 This angle is important because it produces a lift force perpendicular 垂直地 to the axis of the plane and hence a normal acceleration also shown in an the figure v Normal Accel an Pitch Angle Fig 2 1 Schematic of aircraft attitude 10 The pilot wants to be able to control the pitch angle and does so ultimately by rotating the front fins and tail elevators of the aircraft shown in Fig 2 2 The first task is therefore to model the effect of these movements on the pitch rate and normal acceleration q an Fig 2 2 Illustration of control surfaces 11 2 1 MODELLING THE AIRCRAFT AND ACTUATOR DYNAMICS i Normal Acceleration an The acceleration of the aircraft in a direction perpendicular to its axis the normal accel is an determined mainly by the angle of the tail elevators of the aircraft shown in Fig 2 2 Indeed 1 aerodynamic modelling shows that this relationship can be described by the differential equation 111 16342 36 4 1982 8 1 3 nnn aaa Convert this relationship into a transfer function form 转移函数 1 2 2 19 436 21634 s 3 18 82 ss H ss ii Pitch Rate q The rate at which the pitch angle changes the pitch rate q is determined mainly by the angle of 2 the front fins of the aircraft shown in Fig 2 2 Indeed aerodynamic modelling shows that this relationship can be described by the differential equation 22 27 1125 782 81 3 qqq Convert this relationship into a transfer function form 转移函数 2 2 7 2511 27 s 3 18 82 s H ss iii Actuator Dynamics Gears 齿轮齿轮 The tail elevators of the aircraft are driven directly by a hydraulic actuator which has a transfer function 3 G ss 1414 To a unit step input U s 1 s the response R s G s U s What form of 14 11 14 14 ssss step response r t would you expect from this actuator and what is its time constant 进行 L T 得 1 14 tetr t The front fins and the tail elevators are both driven from the actuator by the same drive shaft through a gear box The same gear on this drive shaft connects to the tail elevator gear wheel which has 500 teeth and the front fin gear wheel which has 100 teeth the relationship between 1 and 2 is 4 5 2 1 1 2 t t 12 2 2 SIMULATING THE AIRCRAFT AND ACTUATOR DYNAMICS i For the aircraft model you will need from the Continuous or Maths Library a Transfer Fcn block to represent the tail elevator angle normal accel relationship eqn 1 a Transfer Fcn block to represent front fin angle pitch rate relationship see eqn 2 a Transfer Fcn block to represent the hydraulic 液压的 actuator dynamics see eqn 3 a Gain block to represent the gears see eqn 4 from the Sources Library a Step Input block to act as a test input set Step time 5sec and Final value 0 01 rad from the Sinks Library a Scope block to monitor the pitch rate output set Horizontal Range 20sec ii Now connect the blocks together using the step input to drive the actuator and the actuator to drive the tail elevators The actuator is also used to drive the front fins via the gearbox The Scope should be connected to the output of the front fin angle pitch rate block Note you can take a branch from a signal line by dragging with the right mouse button from the point on the line at which you want the branch to start Plot the block diagram you have created in the box below iii From the pull down menus of your simulation window select Simulation Parameters in order to define how the simulation is to be performed Set 13 Stop Time 20 sec Initial Step Size 0 001 Max Step Size 0 01 iv Finally double click on your scope so that you can watch the progress of the simulation and then select Simulation Start from the pull down menu of your simulation window You may want to adjust the Vertical Range of your Scope and re run the simulation to get a good picture Plot the response of the aircraft pitch rate below and comment on these results Note that it can be shown that if a transfer function has a denominator polynomial 分母多项式 with a zero or negative coefficient then it must have a root with a zero or positive real part Plot of pitch rate response to a step change in control input the tail elevatorsthe front fins test input 14 Comments on response Would you fly this plane What do the results mean physically 2 3 SIMULATING THE PITCH RATE CONTROL SYSTEM The designers are quite glad they simulated the response of the aircraft before trying it for real They now decide to improve the response by measuring the output ie what is actually happening to the pitch rate and subtracting this from the input ie what they would like to happen to produce an error signal The error signal will be amplified and then used to drive the actuator see Fig 2 3 In this way the actuator will automatically act so as to reduce any differences between the input demand and the output response Fig 2 3 Feedback control of pitch rate i Before you simulate the control system described above you will probably find that you are short of space in your simulation window Besides it would be nice to package the aircraft actuator dynamics as a single block so that they are distinct from anything added later as shown in Fig 2 3 If you press the left mouse button in an empty space of the window you get a rubber band box which you can drag to surround a group of elements When you release the button all elements inside are selected Use this technique to select all the aircraft actuator dynamics but leave out the step input and scope monitor From the pull down menu under Edit select Create a Subsystem and your diagram will suddenly look as shown in Fig 2 4 Step Input Scope Subsystem Fig 2 4 New diagram ii Rename the Subsystem as Actuator Aircraft Dynamics as shown in Fig 2 3 by editing the text below the block You can see what is inside by double clicking on the block itself try it Notice that Pitch Rate Signal Generator Sum Actuator Aircraft Dynamics 0 5 Gain 15 the input and output connections are labelled in 1 out 1 etc Edit these to give them more meaningful names iii Now construct the planned control system as illustrated in Fig 2 3 Notice the pilot is now represented by a Signal Generator from the Sources library which should be set to produce a Square Wave of frequency 0 1 Hz and peak amplitude 0 5 iv The designers are unsure what value to use for the Gain Open the Pitch Rate Scope and run the simulation first for the Gain 0 5 as shown above and then for the Gain 5 Plot the two responses you obtain and comment on the relative advantages and disadvantages of each Plot of pitch rate response to a square wave in control input a with Gain 0 5 b with Gain 5 16 Comments on responses in terms of overshot and response time etc PART 3 SIMILATION OF NONLINEAR SHIP ROLL DYNAMICS The rolling motion of ships is of considerable interest to naval architects because even today around 50 of ships lost at sea sink as a result of a capsize The aim of this laboratory is to study such behaviour by constructing a nonlinear differential equation model for the system converting this model into a phase variable block diagram form suitable for simulation simulating the nonlinear response to sinusoidal excitation calculating the steady state response for a ship whose cargo has shifted dangerously to one side u t Fig 3 1 Ship schematic 17 3 1 MODELLING AND SIMULATION OF SHIP DYNAMICS i Consider the ship sketch shown in Fig 3 1 Given that The effective inertia 惯性 of the ship about its roll axis J 1 The damping moment 阻尼力矩 torque 扭矩 due to friction 摩擦 between the body of the ship and the water and due to turbulence 湍流 round the bilge keels 2 cc 12 3 The restoring moment torque due to the buoyancy 浮力 and shape of the ship as it is pushed to one side 3 kkk 12 3 3 5 The input forcing moment due to the wave forces acting on the ship u t Using Newton s second law of motion the nonlinear differential TorqueMomentJ equation which describes the rolling motion of the ship can be written as 4 5 3 3 21 3 21 kkkcctuJ ii Re arrange equation 4 to obtain an expression for the highest derivative of the output signal 5 5 3 3 21 3 21 1 kkkcctu J iii Hence complete the drawing of the SIMULINK block diagram for the ship roll system 18 3 2 SIMULATING THE SHIP DYNAMICS i For the ship model you will need from the Continuous and Maths Library a few Integrator blocks assume zero initial conditions for now a few Gain blocks not absolutely necessary from the Nonlinear Library a couple of Fcn blocks to implement the damping and stiffness terms eqns 2 and 3 from the Sources Library a Sine Wave block to simulate the sea waves u t set Ampl 0 3 and Freq 0 6 rad s from the Sinks Library a To Workspace block for the output ship roll angle of the simulation This stores the results in a MATLAB variable which we can plot later set the Variable name theta and the Maximum number of timesteps ie points stored 5000 an oscilloscope if you wish to see a waveform whilst the simulation is taking place ii Drag the appropriate blocks over into your simulation window using the mouse and double click on them to enter the appropriate parameter values The numerical values for the coefficients of eqn 4 are 6 150 20 8 0 8 0 02 0 1 32121 kkkccJ Note these values are given in consis