Faculty of Engineering Spring 2017 Mechanical Engineering Department

Matlab Sheet 4 - Solution

Conditional Statements and Loops

1. It is desired to compute the sum of the first 10 terms of the series

14𝑘 3 − 20𝑘 2 + 5𝑘. 𝑘 = 1,2,3, …
Write and run the program to calculate the sum.
The following M-file gives the correct sum:
M-File:
total = 0;
for k = 1:10
total(k+1) =total(k)+ 14*k^3 - 20*k^2 +5*k;
end
disp('The sum for 10 terms is: ')
disp(total(k+1))
Command Window:
The sum for 10 terms is:
34925
The following M-file would not give the correct sum (why?).
M-File:
total = 0;
for k = 1:10
total = 14*k^3 - 20*k^2 +5*k;
end
disp('The sum for 10 terms is: ')
disp(total)
Command Window:
The sum for 10 terms is:
12050

Dr./ Ahmed Nagib Elmekawy 1 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

2. Create Find the results of the following operations by hand and use MATLAB to
check your results.
a. 𝑧 = 6 > 3 + 8
b. 𝑧 = 6 + 3 > 8
c. 𝑧 = 4 > (2 + 9)
d. 𝑧 = (4 < 7) + 3
e. 𝑧 = 4 < 7 + 3
f. 𝑧 = (4 < 7) − 5
g. 𝑧 = 4 < (7 ∗ 5)
h. 𝑧 = 2⁄5 >= 5
The answers are
(a) z = 0, (b) z = 1, (c) z = 0,
(d) z = 4, (e) z = 1, (f) z = 5,
(g) z = 1, (h) z = 0

3. Given: 𝑎 = −2, 𝑏 = 3, 𝑐 = 5. Evaluate the following expressions without

using MATLAB. Check the answers with MATLAB.
(a) 𝑦 = 𝑎 − 𝑏 > 𝑎 − 𝑐 < 𝑏
(b) 𝑦 = −4 < 𝑎 < 0
(c) 𝑦 = 𝑎 − 𝑐 < = 𝑏 > 𝑎 + 𝑐
M-File:
clear, clc
a=-2; b=3; c=5;
disp('Part (a)')
y=a-b>a-c<b
disp('Part (b)')
y=-4<a<0
disp('Part (c)')
y=a-c<=b>a+c
disp('Part (d)')
y=3*(c+a~=a/b-b)==(a+c)~=b
Command Window:
Part (a)
y =
1
Part (b)
y =

Dr./ Ahmed Nagib Elmekawy 2 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

0
Part (c)
y =
0
Part (d)
y =
1

4. For the arrays x and y given below, use MATLAB to find all the elements in x that
are greater than the corresponding elements in y.
𝑥 = [−3 0 0 2 6 8] 𝑦 = [−5 − 2 0 3 4 10]
M-File:
x = [-3,0,0,2,6,8];
y = [-5,-2,0,3,4,10];
n = (x>y)
z = find(x>y)
Command Window:
n=
1 1 0 0 1 0
z=
125
Thus the first, second, and fifth elements of x are greater than the
corresponding elements
in y.

5. The arrays 𝑝𝑟𝑖𝑐𝑒_𝐴 and 𝑝𝑟𝑖𝑐𝑒_𝐵 given below contain the price in dollars of two
stocks over 10 days. Use MATLAB to determine how many days the price of
stock A was above the price of stock B.
𝑝𝑟𝑖𝑐𝑒_𝐴 = [19 18 22 21 25 19 17 21 27 29]
𝑝𝑟𝑖𝑐𝑒_𝐵 = [22 17 20 19 24 18 16 25 28 27]
M-File:
price_A = [19,18,22,21,25,19,17,21,27,29];
price_B = [22,17,20,19,24,18,16,25,28,27];
length(find(price_A>price_B))
ans =
7

Dr./ Ahmed Nagib Elmekawy 3 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

6. Suppose that x = [-3, 0, 0, 2, 5, 8] and y = [-5, -2, 0, 3, 4, 10]. Find the results of
the following operations by hand and use MATLAB to check your results.
a. 𝑧 = 𝑦 < ~𝑥
b. 𝑧 = 𝑥 & 𝑦
c. 𝑧 = 𝑥 | 𝑦
M-File:
x = [-3, 0, 0, 2, 5, 8]
y = [-5, -2, 0, 3, 4, 10]
z0=~x
z1= y<~x
z2= x&y
z3=x|y
The answers are
(a) z = [1,1,1,0,0,0], (b) z = [1,0,0,1,1,1]
(c) z = [1,1,0,1,1,1]

7. Evaluate the following expressions without using MATLAB. Check the answers
with MATLAB
a) -3&3
b) ~5<4&~0> -3
c) -2&2>3|8/3
M-File:
clear, clc
disp('Part (a)')
-3&3
disp('Part (b)')
~5<4&~0>-3
disp('Part (c)')
-2&2>3|8/3
disp('Part (d)')
-3<-1<~0|5<4<3
Command Window:
Part (a)
ans =
1
Part (b)
ans =
1

Dr./ Ahmed Nagib Elmekawy 4 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

Part (c)
ans =
1
Part (d)
ans =
1
8. Create The height and speed of a projectile (such as a thrown ball) launched with
a speed of at an angle A to the horizontal are given by
ℎ(𝑡) = 𝑣0 𝑡 sin 𝐴 − 0.5g𝑡 2

𝑣(𝑡) = √𝑣02 − 2𝑣0 g𝑡 sin 𝐴 + g 2 𝑡 2

where g is the acceleration due to gravity. The projectile will strike the ground
when ℎ(𝑡) = 0, which gives the time to hit 𝑡ℎ𝑖𝑡 = 2(𝑣0 /g) sin 𝐴. Suppose that 𝐴
= 30°, 𝑣0 = 40 m/s, and g =9.8 m/s2. Use the MATLAB relational and logical
operators to find the times when
a. The height is no less than 15 m.
b. The height is no less than 15 m and the speed is simultaneously no greater
than 36 m/s.
c. The height is less than 5 m or the speed is greater than 35 m/s.
It is helpful to plot the height and speed versus time before answering the
questions, especially part (c). The plot is shown in the figure.

Dr./ Ahmed Nagib Elmekawy 5 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

The rest of the following script file can be used to compute the answers with
more precision than can be obtained by reading the plot.

M-File:
v0 = 40;g = 9.81;A = 30*pi/180;
t_hit = 2*v0*sin(A)/g;
t = 0:t_hit/100:t_hit;
h = v0*t*sin(A)-0.5*g*t.^2;
v = sqrt(v0^2-2*v0*g*sin(A)*t+g^2*t.^2);
plot(t,h,t,v),xlabel(0Time (sec)0),
gtext(0Height0),gtext(0Speed0),grid
%
% Part (a).
ua = find(h>=15);
t1a = (ua(1)-1)*(t_hit/100)
t2a = ua(length(ua)-1)*(t_hit/100)
%
% Part (b).
ub = find(h>=15&v<=36);
t1b = (ub(1)-1)*(t_hit/100)
t2b = ub(length(ub)-1)*(t_hit/100)
%
% Part (c).
uc = find(~(h<5|v>35));
t1c = (uc(1)-1)*(t_hit/100)
t2c = uc(length(uc)-1)*(t_hit/100)

When run, this file produces the results: t1a = 1.0194, t2a = 3.0581, t1b =
1.0601, t2b = 3.0173, t1c = 1.5494, t2c = 2.5280.
Thus the height is no less than 15 meters for 1.0194 ≤ t ≤3.0581 seconds.
The height is no less than 15 meters and the speed is simultaneously no
greater than 36 meters/second for 1.0601 ≤ t ≤3.0173 seconds.
The height is less than 5 meters or the speed is greater than 35
meters/second for t ≤1.5494 t ≥2.528 seconds.

Dr./ Ahmed Nagib Elmekawy 6 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

Part C
time (sec) h<5 v>35 h<5|v>35 time (sec) h<5 v>35 h<5|v>35 time (sec) h<5 v>35 h<5|v>35
0 TRUE TRUE TRUE 1.508665 FALSE TRUE TRUE 3.017329 FALSE TRUE TRUE
0.040775 TRUE TRUE TRUE 1.549439 FALSE FALSE FALSE 3.058104 FALSE TRUE TRUE
0.081549 TRUE TRUE TRUE 1.590214 FALSE FALSE FALSE 3.098879 FALSE TRUE TRUE
0.122324 TRUE TRUE TRUE 1.630989 FALSE FALSE FALSE 3.139653 FALSE TRUE TRUE
0.163099 TRUE TRUE TRUE 1.671764 FALSE FALSE FALSE 3.180428 FALSE TRUE TRUE
0.203874 TRUE TRUE TRUE 1.712538 FALSE FALSE FALSE 3.221203 FALSE TRUE TRUE
0.244648 TRUE TRUE TRUE 1.753313 FALSE FALSE FALSE 3.261978 FALSE TRUE TRUE
0.285423 FALSE TRUE TRUE 1.794088 FALSE FALSE FALSE 3.302752 FALSE TRUE TRUE
0.326198 FALSE TRUE TRUE 1.834862 FALSE FALSE FALSE 3.343527 FALSE TRUE TRUE
0.366972 FALSE TRUE TRUE 1.875637 FALSE FALSE FALSE 3.384302 FALSE TRUE TRUE
0.407747 FALSE TRUE TRUE 1.916412 FALSE FALSE FALSE 3.425076 FALSE TRUE TRUE
0.448522 FALSE TRUE TRUE 1.957187 FALSE FALSE FALSE 3.465851 FALSE TRUE TRUE
0.489297 FALSE TRUE TRUE 1.997961 FALSE FALSE FALSE 3.506626 FALSE TRUE TRUE
0.530071 FALSE TRUE TRUE 2.038736 FALSE FALSE FALSE 3.547401 FALSE TRUE TRUE
0.570846 FALSE TRUE TRUE 2.079511 FALSE FALSE FALSE 3.588175 FALSE TRUE TRUE
0.611621 FALSE TRUE TRUE 2.120285 FALSE FALSE FALSE 3.62895 FALSE TRUE TRUE
0.652396 FALSE TRUE TRUE 2.16106 FALSE FALSE FALSE 3.669725 FALSE TRUE TRUE
0.69317 FALSE TRUE TRUE 2.201835 FALSE FALSE FALSE 3.710499 FALSE TRUE TRUE
0.733945 FALSE TRUE TRUE 2.24261 FALSE FALSE FALSE 3.751274 FALSE TRUE TRUE
0.77472 FALSE TRUE TRUE 2.283384 FALSE FALSE FALSE 3.792049 FALSE TRUE TRUE
0.815494 FALSE TRUE TRUE 2.324159 FALSE FALSE FALSE 3.832824 TRUE TRUE TRUE
0.856269 FALSE TRUE TRUE 2.364934 FALSE FALSE FALSE 3.873598 TRUE TRUE TRUE
0.897044 FALSE TRUE TRUE 2.405708 FALSE FALSE FALSE 3.914373 TRUE TRUE TRUE
0.937819 FALSE TRUE TRUE 2.446483 FALSE FALSE FALSE 3.955148 TRUE TRUE TRUE
0.978593 FALSE TRUE TRUE 2.487258 FALSE FALSE FALSE 3.995923 TRUE TRUE TRUE
1.019368 FALSE TRUE TRUE 2.528033 FALSE FALSE FALSE 4.036697 TRUE TRUE TRUE
1.060143 FALSE TRUE TRUE 2.568807 FALSE TRUE TRUE 4.077472 TRUE TRUE TRUE
1.100917 FALSE TRUE TRUE 2.609582 FALSE TRUE TRUE
1.141692 FALSE TRUE TRUE 2.650357 FALSE TRUE TRUE
1.182467 FALSE TRUE TRUE 2.691131 FALSE TRUE TRUE
1.223242 FALSE TRUE TRUE 2.731906 FALSE TRUE TRUE
1.264016 FALSE TRUE TRUE 2.772681 FALSE TRUE TRUE
1.304791 FALSE TRUE TRUE 2.813456 FALSE TRUE TRUE
1.345566 FALSE TRUE TRUE 2.85423 FALSE TRUE TRUE
1.38634 FALSE TRUE TRUE 2.895005 FALSE TRUE TRUE
1.427115 FALSE TRUE TRUE 2.93578 FALSE TRUE TRUE
1.46789 FALSE TRUE TRUE 2.976555 FALSE TRUE TRUE

Dr./ Ahmed Nagib Elmekawy 7 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

9. Write a program in a M-File that finds the smallest even integer that is divisible
by 13 and by 16 whose square root is greater than 120. Use a loop in the
program. The loop should start from 1 and stop when the number is found. The
program prints the message “The required number is:” and then prints the
number.
M-File:
clear, clc
i=0;
s=0;
while s<=120
i=i+1;
if rem(i,2)==0 && rem(i,13)==0 && rem(i,16)==0
s=sqrt(i);
end
end
fprintf('The required number is: %i\n',i)

Command Window:

The required number is: 14560

Dr./ Ahmed Nagib Elmekawy 8 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

10. Write a program in a script file that determines the real roots of a quadratic
equation. Name the file “quadroots”. When the file runs, it asks the user to enter
the values of the constants a, b, and c. To calculate the roots of the equation the
program calculates the discriminant D, given by:
𝐷 = 𝑏 2 − 4𝑎𝑐
If D > 0, the program displays message “The equation has two roots,” and the
roots are displayed in the next line.
If D = 0, the program displays message “The equation has one root,” and the root
is displayed in the next line.
If D < 0, the program displays message “The equation has no real roots.”
Run the script file in the Command Window three times to obtain solutions to
the following three equations:
a. 3𝑥 2 + 6𝑥 + 3 = 0
b. −3𝑥 2 + 4𝑥 − 6 = 0
c. −3𝑥 2 + 7𝑥 + 5 = 0

M-File:
clear, clc
for k=1:3
disp('For the equation ax^2+bx+c')
a=input('Enter a: ');
b=input('Enter b: ');
c=input('Enter c: ');
D=b^2-4*a*c;
if D<0
fprintf('\nThe equation has no real roots.\n\n')
elseif D==0
root=-b/(2*a);
fprintf('\nThe equation has one root,\n')
fprintf(' %.3f\n\n',root)
else
r1=(-b+sqrt(D))/(2*a);
r2=(-b-sqrt(D))/(2*a);
fprintf('\nThe equation has two roots,\n')
fprintf(' %.3f and %.3f\n\n',r1,r2)
end
end

Dr./ Ahmed Nagib Elmekawy 9 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

Command Window:

For the equation ax^2+bx+c

Enter a: 3
Enter b: 6
Enter c: 3
The equation has one root,
-1.000
For the equation ax^2+bx+c
Enter a: -3
Enter b: 4
Enter c: -6
The equation has no real roots.
For the equation ax^2+bx+c
Enter a: -3
Enter b: 7
Enter c: 5
The equation has two roots,
-0.573 and 2.907

Dr./ Ahmed Nagib Elmekawy 10 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

11. The Pythagorean theorem states that 𝑎2 + 𝑏 2 = 𝑐 2 . Write a MATLAB program

in a script file that finds all the combinations of triples a, b, and c that are positive
integers all smaller or equal to 50 that satisfy the Pythagorean theorem. Display
the results in a three-column table in which every row corresponds to one triple.
The first three rows of the table are:
3 4 5
5 12 13
6 8 10
M-File:
clear, clc
id=1;
for k=1:50
for j=k+1:50
for i=j+1:50
if i^2==k^2+j^2
a(id)=k;
b(id)=j;
c(id)=i;
id=id+1;
end
end
end
end
table=[a' b' c']
Command Window:
table =
3 4 5
5 12 13
6 8 10
7 24 25
8 15 17
9 12 15
9 40 41
10 24 26
12 16 20
12 35 37
14 48 50
15 20 25
15 36 39
16 30 34
18 24 30
20 21 29
21 28 35

Dr./ Ahmed Nagib Elmekawy 11 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

24 32 40
27 36 45
30 40 50

12. Rewrite the following statements to use only one if statement.

𝑖𝑓 𝑥 < 𝑦
𝑖𝑓 𝑧 < 10
𝑤 =𝑥∗𝑦∗𝑧
𝑒𝑛𝑑
𝑒𝑛𝑑
M-File:
if (x<y)&(z<10)
w = x*y*z
end

Dr./ Ahmed Nagib Elmekawy 12 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

13. Figure 1 shows a mass-spring model of the type used to design packaging
systems and vehicle suspensions, for example. The springs exert a force that is
proportional to their compression, and the proportionality constant is the
spring constant k. The two side springs provide additional resistance if the
weight W is too heavy for the center spring. When the weight W is gently placed,
it moves through a distance x before coming to rest. From statics, the weight
force must balance the spring forces at this new position. Thus
𝑊 = 𝑘1 𝑥 𝑖𝑓 𝑥 < 𝑑
𝑊 = 𝑘1 𝑥 + 2𝑘2 (𝑥 − 𝑑) 𝑖𝑓 𝑥 ≥ 𝑑

Figure 1.
These relations can be used to generate the plot of 𝑥 versus 𝑊.
a. Create a function file that computes the distance x, using the input
parameters 𝑊, 𝑘1 , 𝑘2 , and d. Test your function for the following two
cases, using the values 𝑘1 = 104 N/m; 𝑘2 = 1.5 × 104 N/m; 𝑑 = 0.1 m.
𝑊 = 500 N
𝑊 = 2000 N
b. Use your program to plot 𝑥 versus 𝑊 for 0 ≤ 𝑊 ≤ 3000 N for the values
of 𝑘1 , 𝑘2 and 𝑑 given in part a.
M-File:
(a) The function file is
function package(W,k1,k2,d);
if W < k1*d
x = W/k1;

Dr./ Ahmed Nagib Elmekawy 13 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

else
x = (W+2*k2*d)/(k1+2*k2);
end
disp(‘The distance traveled is:’)
disp(x)

The session is
>>package(500,10000,15000,0.1)
‘The distance traveled is:’
0.0500
>>package(2000,10000,15000,0.1)
‘The distance traveled is:’
0.1250

(b) The function file in part (a) must be modified to provide an output
variable x and to eliminate the display statements. It is

function x = package(W,k1,k2,d);
if W < k1*d
x = W/k1;
else
x = (W+2*k2*d)/(k1+2*k2);
end

Note that the inputs must be scalars because of the if statement. Thus the
following script file would not give the correct plot. See the text for a
discussion of this pitfall.
W = 0:5:3000;
plot(W,package(W,10000,15000,0.1)
The correct script file to obtain the plot is the
following.
k1 = 10000;k2 = 15000;d = 0.1;
if k1*d <= 3000
W = [0,k1*d,3000];
x = [package(W(1),k1,k2,d),package(W(2),k1,k2,d),...
package(W(3),k1,k2,d)];plot(W,x,0- 0),...
xlabel(‘Weight (N)’),ylabel(‘Distance (m)’),...
title(‘k1 = 10,000, k2 = 15,000, d = 0.1’)
else
W = [0,3000];
x = [package(W(1),k1,k2,d),package(W(2),k1,k2,d)];
plot(W,x,0-0),xlabel(‘Weight (N)’),...
ylabel(‘Distance (m)’),title(‘k1 = 10,000, k2 =
15,000, d = 0.1’)
end

Dr./ Ahmed Nagib Elmekawy 14 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

It is easier to solve this problem using a for loop. The following script file is
the solution using such a loop.
W = linspace(0,3000,500);
for k = 1:500
x(k) = package(W(k),10000,15000,0.1);
end
plot(W,x)

14. Consider the array A.

Write a program that computes the array B by computing the natural

logarithm of all the elements of A whose value is no less than 1, and adding
20 to each element that is equal to or greater than 1. Do this in two ways:
a. By using a for loop with conditional statements.
b. By using a logical array as a mask.

a) M-File:

A = [3, 5, -4; -8, -1, 33; -17, 6, -9];

for m = 1:size(A,1)
for n = 1:size(A,2)
if A(m,n) >= 1
C(m,n) = log(A(m,n));
else
C(m,n) = A(m,n)+20;
end
end
end
A = C

b) M-File:

A = [3, 5, -4; -8, -1, 33; -17, 6, -9]

B = (A >= 1)
A(B) = log(A(B));
A( B) = A( B) + 20;
A

Dr./ Ahmed Nagib Elmekawy 15 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

15. A company has the choice of producing up to four different products with its
machinery, which consists of lathes, grinders, and milling machines. The
number of hours on each machine required to produce a product is given in the
following table, along with the number of hours available per week on each type
of machine. Assume that the company can sell everything it produces. The profit
per item for each product appears in the last line of the table.

M-File:

unit_profit = [100, 150, 90, 120];

hours = [1,2,0.5,3;0,2,4,1;3,1,5,2];
p = [0,0,0,0];
max_profit = 0;
another = [];
P = [];
for p1 = 0:15
for p2 = 0:15
for p3 = 0:8
for p4 = 0:14
p = [p1,p2,p3,p4];
if hours*p0 <= [40, 30, 45]0
profit = unit_profit*p0;
else
profit = 0;
end
if profit > max_profit
max_profit = profit;
production = [p1,p2,p3,p4];
elseif (max_profit - profit)<=1
another = [another;p];
P = [P,profit];
end
end
end
end

Dr./ Ahmed Nagib Elmekawy 16 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

end
disp(0Optimal production is:0)
disp(production)
disp(0The profit is:0)
disp(max_profit)
[m,n] = size(another);
if m>1
disp(0The number of answers giving the same profit
is:0)
disp(m-1)
disp(0The other possible answers are:0)
disp(another(2:m,:))
disp(0Other profit:0)
disp(P(2:m)0)
else
disp(0There is only one solution.0)
end
The results are that the optimal production is given by the vector [10 15 0 0].
Thus we should manufacture 10 units of product 1, 15 units of product 2, and
no units of products 3 and 4. The profit is $3250. There is only one solution.
(b) If we produce 9, 14, 0, and 0 units of each product, the profit would be
9(100) + 14(150) = $3000. If we produce 11, 16, 1, and 1 units of each
product, the profit would be 11(100) + 16(150) + 1(90) + 1(120) = $3710, but
this would require 46.5 hours on the lathe, 37 hours on the grinder, and 56
hours on the milling machine, which is more than the available hours.

Dr./ Ahmed Nagib Elmekawy 17 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

16. Create Use a while loop to determine how many terms in the series 2k, k = 1,
2, 3, . . . , are required for the sum of the terms to exceed 2000. What is the sum
for this number of terms?
M-File:
sum = 0;k = 0;
while sum <= 2000
k = k + 1;
sum = sum + 2^k;
end
k
sum
The answers are k = 10 and sum = 2046.

17. Compute Use a loop in MATLAB to determine how long it will take to accumulate
$1,000,000 in a bank account if you deposit $10,000 initially and $10,000
at the end of each year; the account pays 6 percent annual interest.
M-File:
amt = 10000;
k = 0;
while amt < 1e+6
k = k+1;
amt = amt*1.06 +10000;
end
amt
k
The result is amt = 1.0418e+006 and k = 33. Thus, after 33 years, the amount
will be
$1,041,800.

Dr./ Ahmed Nagib Elmekawy 18 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

3
18. One numerical method for calculating the cubic root of a number √𝑃 , is in
iterations. The process starts by choosing a value 𝑥1 as a first estimate of the
solution. Using this value, a second, more accurate value 𝑥2 can be calculated
with 𝑥2 = (𝑃⁄𝑥12 + 2𝑥1 )/3, which is then used for calculating a third, still more
accurate value 𝑥3 , and so on. The general equation for calculating the value of
from the value of 𝑥𝑖 is 𝑥𝑖+1 = (𝑃⁄𝑥𝑖2 + 2𝑥𝑖 )/3. Write a MATLAB program that
calculates the cubic root of a number. In the program use 𝑥1 = 𝑃 for the first
estimate of the solution. Then, by using the general equation in a loop, calculate
new, more accurate values. Stop the looping when the estimated relative error
𝑥𝑖+1 −𝑥𝑖
E defined by 𝐸 = | | is smaller than 0.00001. Use the program to calculate:
𝑥𝑖
3 3 3
a. √100 b. √53701 c. √19.35
M-File:
clear, clc
n=[100 53701 19.35];
for j=1:3
P=n(j);
x=P;
E=1;
while E>.00001
x_old=x;
x=(P/x^2+2*x)/3;
E=abs((x-x_old)/x_old);
end
fprintf('The cube root of %.0f is %.1f\n',P,x)
end
Command Window:
The cube root of 100 is 4.6
The cube root of 53701 is 37.7
The cube root of 19 is 2.7

Dr./ Ahmed Nagib Elmekawy 19 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

19. The overall grade in a course is determined from the grades of 6 quizzes, 3
midterms, and a final exam, using the following scheme: Quizzes: Quizzes are
graded on a scale from 0 to 10. The grade of the lowest quiz is dropped and the
average of the 5 quizzes with the higher grades constitutes 30% of the course
grade. Midterms and final exam: Midterms and final exams are graded on a scale
from 0 to 100. If the average of the midterm scores is higher than the score on
the final exam, the average of the midterms constitutes 50% of the course grade
and the grade of the final exam constitutes 20% of the course grade. If the final
grade is higher than the average of the midterms, the average of the midterms
constitutes 20% of the course grade and the grade of the final exam constitutes
50% of the course grade.
Write a computer program in a script file that determines the course grade for
a student. The program first asks the user to enter the six quiz grades (in a
vector), the three midterm grades (in a vector), and the grade of the final. Then
the program calculates a numerical course grade (a number between 0 and
100). Finally, the program assigns a letter grade according to the following key:
A for 𝐺𝑟𝑎𝑑𝑒 ≥ 90, B for 80 ≤ 𝐺𝑟𝑎𝑑𝑒 < 90, C for 70 ≤ 𝐺𝑟𝑎𝑑𝑒 < 80, D for 60 ≤
𝐺𝑟𝑎𝑑𝑒 < 70, and E for a grade lower than 60. Execute the program for the
following cases:
(a) Quiz grades: 6, 10, 6, 8, 7, 8. Midterm grades: 82, 95, 89. Final exam: 81.
(b) Quiz grades: 9, 5, 8, 8, 7, 6. Midterm grades: 78, 82, 75. Final exam: 81.

M-File:
clear, clc
for j=1:2
quiz=input('Please enter the quiz grades as a vector [x
x x x x x]: ');
mid=input('Please enter the midterm grades as a vector
[x x x]: ');
final=input('Please enter the final exam grade: ');
q_c=(sum(quiz)-min(quiz))/5;
if mean(mid)>final
grade=3*q_c + 0.5*mean(mid) + 0.2*final;
else
grade=3*q_c + 0.2*mean(mid) + 0.5*final;

Dr./ Ahmed Nagib Elmekawy 20 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

end
if grade>=90
letter='A';
elseif grade>=80
letter='B';
elseif grade>=70
letter='C';
elseif grade>=60
letter='D';
else
letter='E';
end
fprintf('\nThe overall course grade is %.1f for a
letter grade of
%s\n\n',grade,letter)
end
Command Window:
Please enter the quiz grades as a vector [x x x x x x]:
[6 10 6 8 7 8]
Please enter the midterm grades as a vector [x x x]:
[82 95 89]
Please enter the final exam grade: 81
The overall course grade is 83.9 for a letter grade of
B
Please enter the quiz grades as a vector [x x x x x x]:
[9 5 8 8 7 6]
Please enter the midterm grades as a vector [x x x]:
[78 82 75]
Please enter the final exam grade: 81
The overall course grade is 79.0 for a letter grade of
C

Dr./ Ahmed Nagib Elmekawy 21 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

20. Given Cam is a mechanical device that transforms rotary motion into linear
motion. The shape of the disc is designed to produce a specified displacement
profile. A displacement profile is a plot of the displacement of the follower as a
function of the angle of rotation of the cam. The motion of a certain cam is
given by the following equations:

M-File:
theta=linspace(0,2*pi,100)
for k=1:100
if theta(k)<=pi/2
y(k)=6*(2*theta(k)-0.5*sin(theta(k)))/pi;
elseif theta(k)<=2*pi/3
y(k)=6;
elseif theta(k)<=4*pi/3
y(k)=6-3*(1-0.5*cos(3*(theta(k)-2*pi/3)));
elseif theta(k)<=3*pi/2
y(k)=3;
elseif theta(k)<=7*pi/4
y(k)=3-1.5*((theta(k)-3*pi/2)/(pi/4))^2;
else
y(k)=0.75-0.75*(1-(theta(k)-7*pi/4)/(pi/4))^2;
end
end
plot(theta,y)
title('Cam Performance')
xlabel('Rotation Angle, rad')
ylabel('Follower Displacement, cm')

Dr./ Ahmed Nagib Elmekawy 22 of 23 Matlab Sheet 4 Solution

Faculty of Engineering Spring 2017 Mechanical Engineering Department

Figure Window:

Dr./ Ahmed Nagib Elmekawy 23 of 23 Matlab Sheet 4 Solution