We Helped With This Mechanical Engineering Assignment: Have A Similar One?

ME2602: Computing, Analytical Methods, Control and Instrumentation

Assignment for MATLAB Programming

References: All information required for completing this assignment can be found in Lecture Notes (Blackboard). It is advisable to complete all tasks in all 3 MATLAB workshops before you start the assignment tasks. For additional reference, any of the following books will be adequate.

•       A. Gilat, MATLAB: An Introduction with Applications, John Wiley & Sons; 4th Ed., 2011

•       S. Chapman, J., MATLAB Programming for Engineers, Nelson Engineering; 4th Ed., 2007 •            H. Moore, MATLAB for Engineers, Pearson Education; 3rd Edition, 2011

MATLAB Tips/Tricks

•       To get info and examples on how to use a function/command, use doc,  e.g.   >> doc interp1

•       To find information for a specific topic on internet, just Google “MATLAB topic”, e.g. Google

“MATLAB interpretation” can find the information about data interpretation with  MATLAB.

•       To return to previously typed commands in Command Window, you can use “arrow up” ↑ key.

•       Upper or Low Case: MATLAB is case‐sensitive: all functions must be in low cases, but variables can be in low cases, upper cases or mixed. NO space in function and variable names.

•       Location to Save MATLAB program or data:  You can change work directory (folder) to any other directory (folder) such as “H:mywork”. See top left immediately above command window.

•       To copy a figure into a word document use “copy figure” under “Edit”. You can change figure copy options and figure properties under “Edit”.

Submission Deadline: Week 23 Wednesday 01/03/2017, 4PM, Please Confirm with TPO

What to Submit:               The report (maximum 12 pages) together with all MATLAB programs on a CD.

Requirement for Programs

•       All tasks in each problem should be solved by one program (script file) together with function files if required.

•       Use ; to suppress all intermittent results from MATLAB command window. However, all final results should be presented on the command window very clearly using formatted texts (fprintf( ), Lecture Notes Page 37 Slides 145 & 146).

•       All plots should have titles and labels for all axes.

•       To make a MATLAB program readable for other program users and yourself, it is important to include some comments and very brief explanation throughout the program. In MATLAB % is used to add comments. On any line of program, anything after % will be ignored when the program is executed. The first contiguous comment on the program becomes the script's help file. It is a good idea to identify main and important variables used in the beginning of program. Please add sufficient comments when you write the program, this will help you to avoid wasting time later when you need to debug the program. However, it is not necessary to comment every line of program. Only comment on important lines.

•       The program should only require the minimum intervention by users. If the user is required to enter any information (e.g. enter some numbers), some clear and self explanatory prompts should be presented on the command window.

Requirements for Report (Check “Sample Assignment Report” about report format required.)

•       It is desirable to include some very brief explanations for the program development, for example, to derive mathematical equations which are used in the program, to use a flow chart to explain the program flow etc.

•       All final results and plots for each task should be presented in the report. When appropriate, some discussions on results should be included.

•       All programs (MATLAB codes) should be included in the report with some annotations which explain the function of a line of code (if it is important) or a block of codes. Do not explain every code line, only explain according to tasks or sub‐tasks.

•       Maximum page number for the report: 12.

Requirements for submitting programs on CD

•       All programs and generated data files should be provided on a CD which is submitted with the report. Before you burn files to the CD, it is advisable to create a folder on the computer and put all files in this folder and test run them within MATLAB.

•       For the report marking, all programs (script files) will be tested to check their correctness.

Marking Criteria:  All problems will be marked according to the following marking criteria.  [80‐100]: In addition to [70‐79], the additional marks will be awarded to quality of MATLAB program and explanations, for example, efficient program codes, proper in‐program explanations/comments.

[70‐79]: The MATLAB programs are correct and correct results are produced (The programs supplied on CD will be tested). Program codes are presented in the report with proper annotations to explain the functions of the program. Results are presented clearly on the report and the results are discussed (if required).

[60‐69]: The MATLAB programs related to majority of tasks are correct and correct results are produced for majority tasks. (The programs supplied on CD will be tested). Program codes are presented in the report with some explanations and comments. Results are presented clearly on the report and the results are discussed (if required).

[50‐59]: The MATLAB programs related to over half of tasks are correct and correct results are produced for over half of tasks. (The programs supplied on CD will be tested). Program codes are presented in the report clearly. Produced results are presented on the report and some discussions are presented (if required).

[40‐49]: The MATLAB programs related to 1/3 of tasks are correct and correct results are produced for over 1/3 of tasks. (The programs supplied on CD will be tested). Program codes are presented in the report. Some results are presented.

[30‐39]: The programs fail to run, however some good efforts have been made to complete all tasks. (The programs supplied on CD will be tested). Some program codes are presented in the report.

[0‐29]: A suitable marks will be given based on the programs produced and the effort made. 

The final mark will be calculated according to the following weighting factors.

All MATLAB commands should be entered into script files, NO command should be entered in command window >> unless it is for testing purpose.

ME2602 MATLAB Assignment Problems

P1:      Use MATLAB to solve linear simultaneous equations (Estimate time

required, 20 minutes.)

Consider the statically determinate truss shown in Fig.P1. The applied force has a magnitude of P=950(1+0.2S) N at an angle of 30° from the horizontal. The inner angles, θ₁ and θ₂ are 45° and 65° respectively.    The      following         linear   simultaneous equations can be used to solve for 6 unknown forces F₁, F2, F3, Fhx, Fhy, Fry.

Note that:  The parameter S in the applied force P is a parameter related to each student’s student id and this parameter must be generated before you start to solve this problem using the following MATLAB codes.

RandStream.setGlobalStream(RandStream('mt19937ar','seed',1406833));

S=rand(1);

                                                                                                replace it with your student id

P1.3     Present the final results on command window using formatted text.

Ideas for solving this problem. P1.1

P1.3:   To present formatted texts using command  fprintf, See Lecture Notes page 37, Slides 145, 146.

P2: Use MATLAB to solve constraint optimization problems in engineering design (Estimate time required 30 minutes)

Consider the two-bar truss shown in Figure P2. The objective is to minimize the volume V of the two bars AC and BC. The cross-sectional areas of the bars AC and BC are denoted as S₁ and S₂ respectively, and the vertical position of joint is denoted as h. The constraints of the system are that the tensile stresses on the two bars must be less than or equal to σ=10⁵kPa and that the vertical position remains between 1m and 3m. Finally, S₁ and S₂ are nonnegative. (fmincon is the command for this type of question.)

The optimization objective is to minimize V: 

                        subject to constraints:                                                                          

Ideas for solving this problem.

•        The above objective function is a function of several variables. To solve this type of problem, the objective function must be first converted into a function of a single vector variable. The same should be done for constraint functions. Please see examples in Workshop 2 W2P3 and Lecture Notes Pages 26, 27, Slides 104, 105. 

•        For the examples of solving constraint optimization, please see examples in Lecture Notes Pages 29‐31, slides 113‐122. Basically, you need to define TWO functions, one function for the cost (objective), another for describing constraints. Both functions can be defined as “in program” anonymous functions or by separate function files. 

•        You need to do some preparation work to convert all constraints to the standard

format so that they can be described by the constraint function.. 

P3: Curve Fitting/Model Identification: find a suitable formula/model to fit experiment data. (Estimate time required 45 minutes)

In a wind tunnel study of dense gas dispersion in a neutral boundary layer over a rough surface, the

⎛⁻z2 ⎞

Gaussian function c = exp_⎜ 2s2 ⎟_⎠ is used to model relationship amongst the average plume _⎝

concentration ratio (c), height (z) and vertical plume speed (s). The task is to find the parameter s for the above model using collected data. This is a curve fitting problem (find parameter s so that the above mathematical model fit the collected data).

Convert a nonlinear model into a linear model so that the linear least square method can be used: 

•       Nonlinear equations such as the one we have in here must be first converted into linear equations. 

                                        c = exp⎜⎛−2sz22 ⎠⎟⎞ ⇒ ln( )c =− 21s2 z2

⎝

                                                 By defining     y= ln( ),c       x z and a= ²       =−  

                The nonlinear equation has been converted to the linear equation              y = ax

The least square method is first used to estimate the parameter a. Once a is estimated, the corresponding estimate for s can be determined.

Linear Least Curve Fitting Method: 

•       For a linear equation y=ax and experiment data set (x(1), y(1)), (x(2), y(2)), …, (x(N), y(N)), the least

square estimation for a is given by 

N

 ^aˆ ^=∑k∑₌^{N1 x k y k}( ( ))^{( ) ( ) =} x y(1) (1)( (1))x +² +x(2) (2)( (2))x y ² + ++ +"" ( (x N y Nx N( ) ())² ) _→ sum x y sum x( .* ) / ( .^ 2) _{x k} 2

k=1

Program steps: 

(1)    Calculate x and y according to the definition above.

•       In theory, the values for c should all be positive. However, some experiment values collected become negative because “measurement noise”. This may cause problem for the calculation of ln(c). This problem can be overcome by using ln(|c|) instead of ln(c).     log(abs(c)) in MATLAB.

(2)    Find the estimation of a using the least square formula. 

(3)    Determine s using a.  (use the variable s1 to store the value of s in order to differentiate with the following parts).

                                     s z c dz

Use plot command to plot five graphs (c versus z, c₁ versus z, c₂ versus z, c₃ versus z, c₄ versus z) on the same figure window. 

                                 →            plot(z, c, z, c1, ….)

P4: Use MATLAB to solve complex engineering problem (Function files, program branches, conditional and logic operations, solving differential equations, numerical integration, 2-d plots) (Estimate time required 30 minutes)

Fig. P4.1 represents the “quarter-car” model of a car suspension system. The car travels at a constant horizontal velocity V_H and the road surface has a chuck hole and a bump shown in Fig.P4.2. The equations for describing the dynamics of this system can be derived as

                                   displacement y(t)                               M=350kg:  mass of quarter car

K₂:  suspension spring constant B:  shock absorber damping coefficient m=21kg:  unsprung tire mass

K₁:  tyre compliance spring constant

Figure P4.1

Figure P4.2

P4.1 The first step to solve differential equations is to create a function to describe derivative on the right hand side of differential equations. Write a MATLAB function m-file to describe the dynamic equations with dotx as the function output. In addition to usual inputs of time (t) and x as inputs, the variable parameters V_H, B, K₁ and K₂ will be used as input variables as well  

The structure of this function will be

function dotx=funcname(t, x, VH, B, K1, K2)

M=350; m=21; L=VH*t; if L>=0 & L<2     u=0; elseif L>=2 & L<2.5     u=(-0.06/0.5)*(L-2)

    ????????     ???????     ????????

end y=x(1);z=x(2);vy=x(3);vz=? dotx=[vy;?;(K2/M)*(z-y)+??? ; ???]; end

To test the correctness of this function, try the following to produce a 4x1 vector

>> funcname(3, [1;2;0;3], 19, 1100,20000, 10000)

P4.2 Given parameters B=1140Ns/m, K₂=20900N/m and K₁=10800N/m, V_H=60mph (convert to m/s), use a MATLAB differential equation solver (ode15s, ode23s or ode23t) to simulate the response of the car suspension system for the time period [0, 20] seconds. The initial conditions are all as zeros. 

•       Examples of solving multiple 1^st order differential equations, see Workshop 2 W2P4 and Lecture Notes pages 33 & 34 Slides 129‐135.  

•       For solving differential equations with additional parameters, see Lecture Notes pages 34 & 35 Slides 136 & 137.  

•       The ODE solver ode45 will NOT produce correct results for this problem, therefore ode15s, ode23s or ode23t will be used. General the ODE solver can automatically determine the best step size to use to solve ODEs. However, for this particular problem, the maximum step size must be specified to get a solution. See Below  VERY IMPORTANT.  

(a)               Plot y, z, v_y, v_z versus t in 4 subplots.

o   y, z, v_y, v_z can be picked out from X, i.e.  y=X(:,1),..., vz=X(:,4).  o For using subplot command, see Lecture Notes page 19 Slides 74 & 75. 

(b)               The car’s passenger comfort can be evaluated using the car body’s displacement (y) and vertical acceleration (a = vy ). From simulation results, determine and plot the vertical acceleration (a = vy ) and displacement (y) all against the time t. Find the maximal vertical displacement and the root mean square of the acceleration

a_rms

o   The      vertical            acceleration    can       be        determined     using    the       formula,

K² (z − +y⁾ B (v_z −v_y). a_y = =^vy _{M M}

o   The command trapz(t, ay) can be used to calculate the integration.

P5: Use MATLAB to solve complex engineering problems (functions, differential equation solutions, data interpretation, loop, conditional and logic operations, combining arrays, 2-d plots and subplots) (Estimate time required 90 minutes)

A ball is thrown into the air with the initial velocity v₀ and the angle α to horizontal direction (see Figure P5a). The air is assumed to have friction which is proportional to the velocity. The horizontal and vertical velocities and positions can be determined by solving the following differential equations (derived from the Newton’s law).

⎧⎪dvdtx =− mB vx

⎪

                           ⎪⎪dvdty =− −g     mB vy      

                         ⎨

⎪dx ₌ vx ⎪dt ^⎪dy

                                  ⎪ = vy

⎩dt

=

                        Note that S is the student related parameter generated in the same way as P1.

P5.1 Use a differential equation solver in MATLAB (i.e. ode45, ode23 or ode15s) to solve the above equations to find the solution for v_x, v_y, x, y and corresponding time history t for the time period of [0, 6] seconds. Use subplot command to plot four graphs (v_x versus t, v_y versus t, x versus t, y versus t) on the same figure window. (Don’t present all date of v_x, v_y, x and y on command window screen and the Report. Only graphs need to be presented on the report.)

•       Examples of solving multiple 1^st order differential equations, see Workshop 2 W2P4 and Lecture Notes pages 33 & 34 Slides 129‐135. 

•       The first step is to define a MATLAB function (anonymous function or function file) to describe the ball throw dynamic equations with dotz as the function output, t (time) and z as the input. Note that the only task for this function is to describe four differential equations in vector format, nothing else should be put into this function!! 

⎡dv_x ⎤

                                                ⎢ dt ⎥     dotz(1)                    ⎡v_x⎤z(1)

                                               ⎢      ⎥

⎢dv_y ⎥     ⎢vy⎥⎥z(2) dotz =⎢⎢ dt ⎥⎥    dotz(2) z =⎢⎢ x ⎥z(3) ⎢ dx ⎥       dotz(3) ^⎢_⎣ y ^⎥_⎦z(4)

⎢ dt ⎥

                                                ⎢ dy ⎥     dotz(4)

                                               ⎢      ⎥

                                                ⎣ dt ⎦                                                  

•       Before you move to next step, it is a good idea to test this function, 

>>funcname(2, [1  2  3  4]) should produce a 4x1 column vector. 

•       How to use ODE solver?  

[th, zh]=ode54(@funcname, [0 6], initial condition vector z0). 

Here, th is the history of time t for the period of [0  6], v_x, v_y, x, y can be picked out from zh, i.e.  vxh=z(:,1)..., yh=z(:,4). 

•       For using subplot command, see Lecture Notes page 19 Slides 74 & 75.

P5.2 Find the time (t_f) when the ball strikes the ground, i.e. y=0. The way to find the time t_f is to simulate (solve) the system differential equations for a longer period time (done in P5.1) and then use “interp1” to find t_f when y(t_f)=0.

•       This can be done by the find command (Lecture Notes page 47) and/or the interpretation (Lecture Notes page 49 Slides 189 & 190). From the graph on the right, it can be seen that no data point at exactly y(t)=0, so that some approximations must be made. We can find the first data point when y(t)<0 and the last data point when y(t)>0 and then approximately determine t_f.

Method 1: Very rough approximation 

Pindex=find(y<0); tf=(t(Pindex(1))+t(Pindex(1)-1))/2; 0

Method 2: Approximate with more data points

                Pindex=find(y<0);                     

yselect=y(Pindex(1)-3:Pindex(1)+3); 0 _{Time (seconds)}t _f tselect=t(Pindex(1)-3:Pindex(1)+3); tf=interp1(yselect, tselect,0);

Figure P5b

                Note:    while-end loop is used to complete this task, but for-end loop can also be used

V(1)=v0;     %Vector V to store all bounce velocities v1, v2, .....

Alpha(1)=??; %Vector Alpha to store bounce angles alpha1, alpha2, ..... bn=1;        % The number of bounce

fprintf(' ')

fprintf('Original projection velocity and angle are %g and %g ',v0,alpha0)

fprintf('The original projection velocity is %g ',v0)

while (V(bn) >= 0.5*v0) & (bn <100)  %Restrict maximum loop number to avoid infinite loop     bn=bn+1;

    V(bn)=???*V(bn-1);

    Alpha(bn)=????;

    fprintf('The velocity and angle of bounce %g is %g and %g ',bn-1,V(bn),Alpha(bn)) end

fprintf('Bounce number before bounce velocity less than 0.5 of original velocity is:   %g ',bn-1)

T=[t1; t2+tf1; t3+tf1+tf2;.............];

X=[x1; x2+x1(end); x3+x2(end)+x1(end) .................];

Y=[y1; y2; y3;..............];

VX=[ vx1; vx2; vx3......];

VY=[ vy1; vy2; vy3.......];

NOTE: The commands given above combine COLUMN vectors together. To combine row vectors together, replace ; with ,  e.g. Y=[y1, y2, y3,..............]

P6 MATLAB program branches (Estimate time required 20 minutes, example see Workshop 3 W3P4 )

               Write a MATLAB program to call all programs written for all problems for this assignment. 

When this program is run, a selection menu similar to the following one should be presented.

Select Problem to solve:

  1: Problem P1

  2: Problem P2

……………………………………

  5: Problem 5

Once the selection is made, the relevant program should be executed to solve the problem selected. (Use MATLAB command “switch-case-end”). The program should be able to return to the selection menu before all problems have been solved.