SPSU, UDAIPUR Lab Manual

Department of Mechanical Engineering Numerical Computation

School of Engineering
Udaipur

Laboratory Manual

ON

Numerical Computation with

MATLAB

DEPARTMENT OF MECHANICAL
ENGINEERING
1
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

LIST OF EXPERIMENTS
INDEX

Page No.
S.No. Objective
From to
Evaluation of Arithmetic, Exponential, Logarithmic,
1. Trigonometric, and Complex expressions. 3 4

2. Creating and working with arrays. 5 10

3. Creating, Plotting and Working with anonymous functions 11 12

4. Symbolic Computation 13 14

Programming in MATLAB using SCRIPT FILE and FUNCTION

5. 15 16
FILE

6. Curve Fitting and Interpolation 17 17

Application of MATLAB
7. 18 19

Application oriented Program -1 ( Mechanical vibrations ) 20 21

Application oriented Program -2 (Numerical Methods ) 22 22

Application oriented Program -3 ( Heat and mass transfer) 23 23

Application oriented Program -4 ( Engineering Mechanics) 24 24

2
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

EXERCISE NO.01

1. OBJECTIVE:

Evaluation of Arithmetic, Exponential, Logarithmic, Trigonometric, and

Complex expressions.

2. QUESTIONS

Compute the following:

(1.) Arithmetic operations
1. 25÷(25-1)
2. (1- 1/25)-1

3. [3{ (√5-1)÷(√5 +1)2}]- 1

4. Area = πr2 with ‘r’= π1/3-1
5. 7 + 8/2*5
6. 5^3/2
7. 27^1/3 + 32^0.2
8. 27^(1/3) + 32^0.2
9. 2^2^3
10. 2^(2^3)
11. 2^4/4
12. 4/2*2
13. 1+ 2\4

(2.) Exponential and Logarithm functions

1.
2. ln

3. √

4. 3x=17

3
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

(3.) Trigonometrical operations

1. sin

2. cos π

3. tan

4. sin2 + cos2

5. cosh2 x - sinh2 x , where ‘x’= 32π

(4.) Complex numbers
1. 1+3i/ 1-3i

2.
3. eπ/2*i
4. eπ/2i
(5.): Calculate the following.

(a) [37log(76)/(73 +546)] + √910

(b) 43[ ( √250 + 23)2/ ]

(c) cos sin + [tan 8 / √7]

4
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

EXERCISE NO. 02

1. OBJECTIVE:

Creating and working with arrays.

2. QUESTIONS

Do the following.

(a) Create a column vector that has elements 32, 4, 81, e 2.5, 63, cos , ln(4).

(b) Crate a column vector whose first element is 15, whose elements decrease in
steps of -5 and whose last element is -25.
(c) Crate a column vector whose first element is 26, whose elements decrease in
steps of -3.6 and whose last element is -10.
(d) Create a row vector with 100 equally spaced elements, in which the first element
is 0 and the last element is .

(e) Create a vector A first that has 16 elements in which the first four is 4, increment
is 3 and the last is 49. Then using the colon notation create a new vector. Call it
.A second that has 8 elements. The first four are the first four elements of A first
and the last four elements of A first.
(f) Create a row vector that has the following elements: 8, 10/4 , 12 x 1.4, 51,
tan85° , √26, and 0.15 .

(g) Create a row vector that has the following elements: √15 × 10 , , ,
.
, 129 ( /20)

(h) Create the following matrix A where

6 43 2 11 87
A= 12 6 34 0 5
34 18 7 41 9

2. Evaluate the following by hand and check using MATLAB:

5
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

(a) A=[ 1 2; 3 4]; B= ones(2); Compute A.*B, B./A, A^2, A.^2

(b) X=[ 1 2 3]; y=[2 3 4]; z =[1 2; 3 4] Compute x.^y , 2.^x, 2.^z
(c) Using element by element operations evaluate the following functions:
(i) y= x2/(x3+ 1) where x= 1,2,3,4,5.
(ii) y= sin xcosx/ (x3 + 3x + 2) where x runs from 1 to 7 in steps of 0.2.

Note: y(2) gives the second element of the vector y and not the value of y when x=2.

3. Enter following three matrices

2 6 1 2 −5 5
= = =
3 9 3 4 5 3

4. Check the following Linear algebra rules for above matrices

i) Is matrix addition commutative? A+B=B+A
ii) Is matrix addition associative?  (A+B)+C = A+(B+C)
iii) Is scalar multiplication with matrix distributive?  α*(A+B) = α*A+ α*B
iv) Is scalar multiplication with matrix distributive?  A*(B+C)=A*B+B*C
v) Is matrix multiplication commutative?  A*B=B*A
vi) For scalars a*b = a*c implies that b = c if a≠0. Is it true for matrices? Check
by computing A*B = A*C for the matrices given above.

5. Create the following matrices with the help of matrix generation functions zeros, eye
and ones.
5 0 0
0 0 0 3 3
= = 0 5 0 =
0 0 0 3 3
0 0 5

6. The following matrix G is created by putting matrices A, B and C given above on its
diagonal. In how many ways can you create this matrix using sub matrices A, B and C.
(that is you are not allowed to enter the non-zeros numbers explicitly)?

6
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

2 6 0 0 0 0
⎡3 9 0 0 0 0⎤
⎢0 0 1 2 0 0 ⎥⎥
=⎢
⎢0 0 3 4 0 0⎥
⎢0 0 0 0 −5 5⎥
⎣0 0 0 0 5 3⎦

7. Do the following on matrix G created above

i) Delete the last row and last column of the matrix.
ii) Extract the first 4 × 4 sub matrix from G.
iii) Replace G(5,5) with 4.
iv) What do you get when you type G(13) and hit return? Can you explain
how MATLAB got that answer?
v) What happens when you type G(12,1)=1and hit return?

8. Create a 20 × 20 matrix with the command = (20). Now replace the 10 ×

10 sub matrix between rows 6:15 and columns 6:15 with zeros. See the
structure of the matrix (in terms) of non-zero entries with the command ( ).
Set the 5 × 5 sub matrices in the top in the top right corner and bottom left
corner to zeros and see the structure again.

9. Create an upper triangular matrix with the following command

= (1: 6) + (7: 11,1) + (12: 15,2)
Now use the upper off diagonal terms of A to make Aa symmetric matrix with the
following command.
= + ( + 1)’

10. Create 10 × 10 random matrix with the command = (10). Now do the
following operations.
i) Multiply all the elements by 100 and then round off all elements of the
matrix to integers with the command = ( ).
ii) Replace all elements of A<10 with zeros.
iii) Replace all elements of A>10 with infinity (inf).

7
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

iv) Extract all 30≤ a ij≤50 in a vector b.

11.

12.

13.

8
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

14.

15.

16.

17. Create the following matrix by using vector notation for creating vectors with
constant spacing and/or the linspace command.
130 110 90 70 50 30 10
= 1 2.8333 4 6.6667 8.333 10.1667 12
12 22 32 42 52 62 72
2 5 8 11 14 17
⎡3 6 9 12 15 18 ⎤
⎢ ⎥
18. Given = ⎢4 7 10 13 16 19 ⎥
⎢5 8 11 14 17 20 ⎥
⎣6 9 12 15 18 21 ⎦
5 10 15 20 25 30
= 30 35 40 45 50 55
55 60 65 70 75 80
= [99 98 97 96 95 94 93 92 91]
Create the three arrays in the command window, and then, by writing one
command, replace the last four columns of the first and third rows of A with the first
four columns of the first two rows of B, the last four columns of the fourth row of A

9
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

with the elements 5 through 8 of C and the last four columns of the fifth row of A
with columns 2 through 5 of the third row of B.
19. Compute the output for the following commands by hand and check using
MATLAB:
(a) 4<=2
(b) X=4>2
(c) Y=(7>10) + (2<11) + (5*4= = 20)
(d) A= [ 7 9 8 4]; B= [4 11 2 3]; x=(a>=b)
(e) A= [1 2;3 4]; B= 2*ones(2); A==B
(f) A= [1 2;3 4]; B= 2*ones(2); A>2

20. : Compute the output for the following using MATLAB:

(a) 4&2 , 4|2
(b) 4&0 , 4|0
(c) 6, 0
(d) Y=[(7&10) +(2|0) + (~0)]*2
(e) A= [7 9 8 4]; b= [4 11 2 3]; x= (a&b)
(f) x= [-1 1 1]; y= [1 2 -3];
i. x>0 & y>0
ii. x>0 | y>0
iii. xor (x>0, y>0)
iv. any(x>0)
v. all(x>0)
vi. any(x>0) && any(y>0)
Note by using the command if x>0 && sin(1/x) we can avoid a possible division
by zero.

10
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

EXERCISE NO.03
1. OBJECTIVE:
Creating, Plotting and Working with anonymous functions
2. QUESTIONS:
1.) Plot the following functions for the range given using plot, fplot and ezplot in
MATLAB
i) sin 0 2
ii) cos  0 2
iii) tan −
iv) ( , )= −  = −5 5
v) ( )= −3 + log( − 1) + 100 = −10 10
Note: Use proper labels and legends wherever necessary.

2.) Create the function ( )= − sin + 1 . Do the following computation on the

function using MATLAB:
i) Find (0), (1), ( 2).
ii) Vectorize f and evaluate ( )where = 0 1 2 .
iii) Create = (−1,1), evaluate ( ), and plot . ( ).
iv) Combine the following three commands into a single command to
produce the plot that you will get at the end of third command.
= (−1,1); = ( ); ( , );
v) Use fplot to graph ( ) over form – .

3.) Create three anonymous functions corresponding to the following three

expressions:
( )= +8 + 17 − 4 − 20
( )= −4 +4
ℎ( ) = −4 −5
i) Evaluate ( ) − ( )ℎ( ) at = 3.
ii) Evaluate ( ) − ( )ℎ( ) at = [1 2 3 4 5]

11
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

( )
iii) Evaluate ( ) ( )
for over ∈ [−5,5].
( )
iv) Evaluate ( )
− ℎ( ) for any .

4.) Use the same function ( ) defined previously i.e.

( )= +8 + 17 − 4 − 20 and perform following tasks:
i) Plot the function using fplot over the appropriate interval of and locate all
the zeros of the function ( where ( ) = 0).
ii) Repeat question (i) using function .
iii) Use the function quad to integrate ( ) from = 0 to = 1 and verify resuts
by direct integration.
5.) Plot the function =3 − 26 + 10, and its first and second derivative for
−2 ≤ ≤ 4 , all in same plot.
6.) Plot the function = 2( . )
for 0.1 ≤ ≤ 60 through four different plots with
i) Both the axes linear
ii) Y-axis in log scale and X-axis linear
iii) X-axis in log scale and Y-axis linear
iv) Both the axes in log scale

7.) Plot the parametric curve ( ) = , = sin for 0 < </ 2using ezplot.
8.) Plot cardioid ( ) = 1 + cos for 0 < < 2 using ezpolar.
9.) Plot the following parametric functions using plot, fplot and ezplot for the domain
i) = × (3 )0 2
ii) = × (3 )0 2

iii) = × × − −2 2

10.) Plot cardioid ( ) = 1 + 2 × (2 ) for 0 < < 2 using ezpolar.

11.) Plot using ezplot3(x,y,z), where = × (3 ), = × (3 ) and = for
0 2
12.) Plot using ezsurf = −5/(1 + + ) for −3 ≤ ≤3 −3≤ ≤3

12
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

EXERCISE NO.4

1. OBJECTIVE:
Symbolic Computation
2. QUESTIONS:

Expand the following

(a) (x+y)3 (b) sin(x+y) (c) tan(x+y)

2 2
(d) (x -4x)(x -4x+1)- 20

2. Show that (x2-4x-4) is a factor of (x2-4x)(x2-4x+1)- 20. Find the factors of above
expression. Also find the roots of the above factors.

3. Solve the linear equations:

x+3y-z=2
x-y+z=3
3x-5y=4

4. Solve the following non linear equations:

3x3+ x2-1=0 and x4-10x2+2=0.

5. Find the following integral:

I.(a) ∫
(b) ∫
(c) ∫
II. r=[x y]T x and y: real numbers
(a)
(b)∫ .
(c) ∫ ∫ .

6. Differentiate: cos(x+y), sin(x+y).

7. Solve the following:

(a) + x2 =0; x(0)= x0

(b) + w2y=0; y(0)=y0 and |y=0 = V0.

(c) 4 × +3 = cos ; x(0)= 3

13
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

(d) 2 −6 + 4 = 1; x(0)=1 and |x=0 = 4. Solve without initial condition and

with initial condition. Also plot the solution in the range -10 to 10.

(e) Use ode45 and plot the equation below in given time span.

= ; (0) = 3, ∈ [0 2 × ]

(f) 2 +2 +6 + 4 = 1; ̇ (0) = 1, ̈ (0) = 1, (0) = 1. Plot the

given differential equation.

(g) + 10 + 35 + 50 + 24 = 5 +4 + 24, where

= cos(2 + 1) + 5

(h) Solve the above equation with initial conditions and plot the result. y(0) =
3, ̇ (0) = 2, ̈ (0) = ⃛(0) = 0

8. (a) Find the Laplace transform of f(t) = −1.25 + 3.5t + 1.25

(b) Find the Laplace transform of f(t) = 4 − 4 cos + 2 sin

( )
(c) Find the Inverse Laplace transform of F(s) =
( )

( )
(d) Find the Inverse Laplace transform of F(s) =
( )

9. Find the following integrals:

(a) ∫
( )( )( )

(b) ∫

(c) ∫ cos

14
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

EXPERIMENT NO. 05

1. OBJECTIVE:
Programming in MATLAB using SCRIPT FILE and FUNCTION FILE

2. QUESTIONS:

1. Three forces are applied to a bracket as shown. Determine the total (equivalent)
force applied to the bracket. Given F 1=400N, F2=500N and F3=700N.

2. Find the coefficient of friction μ associated with the block as shown in figure.

The results of six tests are as follows:

m(kg) 2 4 5 10 20 50

F(N) 12.5 23.5 30 61 117 294

3. Find the determinant and division of two matrix shown below.

1 2 3 2 3
= 2 2 = 1 1
1 2 1 3 2
4. Compute simple and compound interest

15
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

= × ×

= [(1 + / ) − 1]

Take amount A=10000 , percentage interest = 10%, time t=10years and n=4.

5. Write a program in MATLAB using script file and function file

i) To determine factorial of a given no.
ii) To evaluate the sum of following series
0 1 2
+ + + ⋯+
0! 1! 2! !
iii) To evaluate
a) SINE SERIES

till n terms

b) COSINE SERIES

c) EXPONENTIAL SERIES

d) FIBONACCI SERIES

1+1 +2+3 +5 +8+⋯+

iv) Solve the given system of linear equations for various values of ‘r’:
5 + 2( ) + ( ) = 2
3 + 6( ) + (2 − 1) = 5
2 + ( − 1) + (3 ) = 5

6. Using the concept of function as input to another function solve the following:
i) Determine root of the following function ( ) = − 32 + ( − 22) +
100 for described value of r (dynamically).
ii) Evaluate the function = sin( ) + and ℎ = 20 − 3 − 2 + 10
for desired values of x,y. (dynamically)

16
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

EXPERIMENT NO. 06

1. OBJECTIVE:

Curve Fitting and Interpolation

2. QUESTIONS:

1. Given the following data:

Mass (kg) 5 10 20 50 100

Deflection (mm) 15.5. 33.07 53.39 140.24 301.03

Find the spring constant assuming linear relationship between face and
deflection.

2. Fit the linear and quadratic fits to given data and also plot residuals and show
norm of residuals

x 0 0.2 0.3 0.4 0.52 0.6 0.7 0.8 0.9 1.0

y 0 2.2 3.45 4.8 6.55 7.8 9.45 11.2 13.05 15

/
3. Fit a curve ( ) = where and are unknown constants to the following
data which shows time vs pressure reading for vacuum pump.

T (sec) 0 0.5 1.0 5.0 10.0 20.0

P (Pa) 760 620 528 85 14 0.16

Also determine the constants and .

4. Consider two sample steps involved in interpolation

x 0 1.570 2.356 3.141 3.927 4.712 5.497 6.283

y 0 1.000 0.707 0.000 -0.707 -1.000 -0.707 0.000

17
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

EXPERIMENT NO.07
1. OBJECTIVE:
Application of MATLAB
2. QUESTIONS:
1. A mass is subjected to two harmonic motions given by 1 ( ) = cos and
2 ( ) = cos( + ) with =1 , = 20 / and = 1 / . Plot the
resulting motion using MATLAB and identify the beat frequency.

2. Plot the variation of natural frequency and time period with the static deflection.
Take range of static deflection.

3. Develop MATLAB Program to find free vibration response of viscously damped

system using following data:
M=450kg k=26519.2N/m c=1000Ns/m Xo=0.539m
̇ =1 m/s

4. Solve the first order differential equation

= + (0) = 0

5. Solve second order differential equation

̈+ sin = 0
With initial conditions (0) = 1, ̇ (0) = 0

6. Figure shows the slider crank mechanism. Calculate and plot the position,
velocity and acceleration of piston for one revolution of crank rotating at
constant speed. Crank=125 mm ; Connecting Rod=250mm

18
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

7. Plot SFD and BMD curves for the beam shown in figure:

19
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

Program -1
An analytical expression for the response of a damped single degree of freedom
system to given initial displacement and velocity is given by
x(t)= Ce  n t cos( d t   )
Where C and  represent the amplitude and phase angle of the response,
respectively having the values
2
  x   0  1   n x0   0 
C  x  n 0  ,   tan 
2

d d
0
   

 d  1   2 n

(a)Plot the response of the system using MATLAB for ω n=5rad/s, ζ = 0.05, 0.1, 0.2
subjected to the initial conditions x(0) = 0, x  0  0  60 cm/s.
(b) Plot the response of the system in above program using MATLAB for ωn = 5
rad/sec and ζ = 1.0 subjected to the initial conditions x(0) = 0, x  0  0  60 cm/s.

20
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

(a)

(b)

21
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

Program 2

The formula for the fourth order Runge-Kutta method (RK4) is given below. Consider
the problem
y  f (t , y )
y (t 0 )  

Define h to be the time step size and ti = t0 + ih. Then the following formula
y0 = 
k1 = hf(t i , yi )
h k
k 2 = h*f(t i + , yi + 1 )
2 2
h k
k 3 = h*f (t i + , yi + 2 )
2 2
k 4 = h *f(t i + h; yi + k 3 )
1
yi+1 = yi +
(k1 + 2k 2 + 2k 3 + k 4 )
6
Solve the following ordinary differential equation by using 4 th order ranga kutta
method by using Matlab

y  y  t 2  1
y(0)  0.5
The exact solution for this problem is y  t 2  2t  1  0.5et and we are interested in the
value of y for 0  t  2
Take the time step h=0.5, 0.05, 0.01

Compare the values and plot the graph to find the error difference

22
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

Program -3
Solve the unsteady 1-D heat conduction equation using the finite difference
method
– Specify the material and thermal conductivity

– Solve the heat conduction equation

– Plot the temperature distribution vs thickness with respect to time

Case 1. How much time required to cool a 50 mm thickness hot rod having
temperature 4000 C is dipped in to a water under room condition with
temperature 250C (for Aluminum and copper)

Case 2. At room temperature aluminum plate of 10 mm thickness is at heated

to 1000 C on one side and 300 C on other side. At what time the temperature
will be in equilibrium condition.

For both cases take time step 0.5 s

Material Aluminum (Al) Gold (Au) Copper (Cu)

Thermal 202.4 W/m∙K 297.73 W/m∙K 387.6 W/m∙K
Conductivity
Density 2719 kg/m3 19320 kg/m3 8978 kg/m3
Cp 871 J/kg∙K 129.81 J/kg∙K 381 J/kg∙K
Thermal 8.546e-05 m2/s 1.19e-04 m2/s 1.13e-04 m2/s
Diffusivity

23
SPSU, UDAIPUR Lab Manual
Department of Mechanical Engineering Numerical Computation

Program -4
1.Figure shows a weight W hung from the end of a horizontal pole of negligible weight.

The pole is attached to the wall by a pivot and is supported by a cable attached to the
wall at a higher point. The tension T, in the cable is given by

Wlc l p
T
d l p2  d 2
where T = tension in the cable, W = weight of the object, lc = length of the cable, lp =
length of the pole and d = distance along the pole at which the cable is attached.

Write a MATLAB program to (a) determine the distance (d) at which the cable can be
attached to the pole in order to minimize the tension in the cable, (b) plot the tension in
the cable as a function of d.
Given: W = 250 N, lp = 50 cm, lc = 40 cm.

2.Figure shows the location of the center of gravity of a 5000 N truck for the unloaded
Condition. The location of the added load WL is at a distance of x inches behind the
rear axle. Write a MATLAB program and plot WL as a function of x for x ranging from 0
to 60 mm.

24