96

PROBLEM SET 4.1

Problem 1

The problem is to …nd the zero of f (x) = x3 75. With f 0 (x) = 3x2 , Newton’s
formula is
x3 75
x x
3x2
Starting with x = 4, succesive iterations yield

43 75
x 4 = 4:229
3(4)2
4:2293 75
x 4:229 = 4:217
3(4:229)2
4:2173 75
x 4:217 = 4:217 J
3(4:217)2

Problem 2

f (x) = x3 3:23x2 5:54x + 9:84

We begin with a root search starting at x = 1 and launch bisection once the
root is bracketed.
x f (x) Interval
1:0 2:070
1:2 0:269
1:4 1:503 (1:2; 1:4)
(1:2 + 1:4)=2 = 1:3 0:624 (1:2; 1:3)
(1:2 + 1:3)=2 = 1:25 0:179 (1:2; 1:25)
(1:2 + 1:25)=2 = 1:225 0:045 (1:225; 1:25
(1:225 + 1:25)=2 = 1:2375 0:067 (1:225; 1:2375)
(1:225 + 1:2375)=2 = 1:2313 0:012 (1:225; 1:2313)
(1:225 + 1:2313)=2 = 1:2282 0:017 (1:2282; 1:2313)
(1:2282 + 1:2313)=2 = 1:2298 0:002 (1:2298; 1:2313)
(1:2298 + 1:2313)=2 = 1:2306 0:005 (1:2298; 1:2306
(1:2298 + 1:2306)=2 = 1:2302 0:002 (1:2298; 1:2302)

The root is x = 1:230 J

PROBLEM SET 4.1 97

Problem 3

f (x) = cosh x cos x 1

The starting points are

x1 = 4 x2 = 5

The …rst step is bisection, giving us the point

x3 = 4:5

Each subsequent step uses the quadratic interpolation formula

f2 f3 x1 (f2 f3 ) + f3 f1 x2 (f3 f1 ) + f1 f2 x3 (f1 f2 )
x=
(f1 f2 )(f2 f3 )(f3 f1 )
to compute the improved value of x, followed by reordering of data points
using the following scheme:

if x < x3 : x2 x3
if x > x3 : x1 x3
x3 x

Here are the results of the computations:

x1 x2 x3 f1 f2 f3 x f (x)
4:000 5:000 4:500 18:850 20:051 10:489 4:907 12:038
4:500 5:000 4:907 10:489 20:051 12:038 4:716 0:818
4:500 4:907 4:716 10:489 12:038 0:818 4:731 0:0582
4:716 4:907 4:731 0:818 12:038 0:0582 4:730

Hence the root is x = 4:730 J.

Problem 4

Newton’s formula is
f (x)
x x
f 0 (x)
where

f (x) = cosh x cos x 1

f 0 (x) = sinh x cos x cosh x sin x

98 PROBLEM SET 4.1

Starting with x = 4:5, successive applications of the formula yield

10:489
x 4:5 = 4:804
34:52
4:573
x 4:804 = 4:735
66:31
0:283
x 4:735 = 4:730
58:20
0:001
x 4:730 = 4:730 J
57:65

Problem 5

f (x) = tan x tanh x

x f (x) Interval
7:0 0:129
7:4 1:049 (7:0; 7:4)
(7:0 + 7:4)=2 = 7:2 0:305 (7:0; 7:2)
(7:0 + 7:2)=2 = 7:1 0:065 (7:0; 7:1)
(7:0 + 7:1)=2 = 7:05 0:036 (7:05; 7:1)
(7:05 + 7:1)=2 = 7:075 0:013 (7:05; 7:075)
(7:05 + 7:075)=2 = 7:063 0:011 (7:063; 7:075)
(7:063 + 7:075)=2 = 7:069 0:000

x = 7:069 J

Problem 6

f (x) = sin x + 3 cos x 2

f 0 (x) = cos x 3 sin x

f (x)
x x
f 0 (x)

PROBLEM 5 99
Starting with x = 2, successive applications of Newton’s iterative formula
yield

4:1577
x 2 = 0:2015
2:3117
0:7392
x 0:2015 = 0:6693
1:5801
0:2676
x 0:6693 = 0:5682
2:6456
0:0093
x 0:5682 = 0:5644
2:4571
0:0000
x 0:5644 = 0:5644 J
2:445
Starting with x = 2, we get

2:3391
x 2 = 1:2560
3:1440
0:1203
x 1:2560 = 1:2087
2:5430
0:0021
x 1:2087 = 1:2078
2:4512
0:0000
x 1:2078 = 1:2078 J
2:4495

Problem 7

f (x) = sin x + 3 cos x 2

xi xi 1
xi+1 = xi f (xi )
f (xi ) f (xi 1 )
Start with x0 = 2, x1 = 1:5:

1:5 ( 2)
x2 = 1:5 ( 2:7853) = 0:4852
2:7853 ( 4:1577)
0:4852 ( 1:5)
x3 = 0:4852 (0:1872) = 0:5491
0:1872 ( 2:7853)
0:5491 ( 0:4852)
x4 = 0:5491 (0:0369) = 0:5648
0:0369 0:1872
0:5648 ( 0:5491)
x5 = 0:5648 ( 0:0013) = 0:5643
0:0013 0:0369
0:5643 ( 0:5648)
x6 = 0:5643 (0:0000) = 0:5643 J
0:0000 ( 0:0013)

100 PROBLEM SET 4.1

Start with x0 = 2, x1 = 1:5:

1:5 2
x2 = 1:5 ( 0:7903) = 1:2449
0:7903 ( 2:3391)
1:2449 1:5)
x3 = 1:2449 ( 0:0921) = 1:2112
0:0921 ( 0:7903)
1:2112 1:2449
x4 = 1:2112 ( 0:0083) = 1:2079
0:0083 ( 0:0921)
1:2079 1:2112
x5 = 1:2079 ( 0:0001) = 1:2079 J
0:0001 ( 0:0083)

Problem 8

f (x) = cosh x cos x 1

(a)
We see from the plot that a root of f (x) = 0 is at approximately x = 4:75.

(b)
The …rst ”improved” value of x predicted by the Newton-Raphson formula is
at the intersection of the tangent at x = 4 and the x-axis. Since the tangent
is almost horizontal, the intersection point is o¤ the right end of the plot (in
x > 8). It is clear that subsequent iterations would keep x away from the root
near 4.75.

PROBLEM 8 101
We can con…rm our …ndings from the Newton-Raphson formula:
f (x) f (4) 18:85
x x =4 =4 = 10:66
f 0 (x) f 0 (4) 2:829
which is indeed ”o¤ the chart”.

Problem 9

f (x) = x3 1:2x2 8:19x + 13:23

f 0 (x) = 3x2 2:4x 8:19

In Example 4.7 it was suggested that if m is the multiplicity of the root,

convergence can be improved by using the modi…ed version
f (x)
x x m
f 0 (x)

of the Newton-Raphson formula (in our case m = 2).

Starting with x = 2, we get
0:05
x 2 2 = 2:1010
0:99
5:205 10 6
x 2:1010 2 = 2:1000
0:0103
0:000
x 2:1000 2 = 2:1000 J
1:02 10 6

Problem 10

The easiest way to handle this problem is to simply replace bisect with ridder
in Example 4.3. We chose a slightly di¤erent approach and wrote the function
rootsRidder loosely based on Example 4.3:

function rootsRidder(func,a,b,dx,tol)
% Computes all the roots of func(x) in the interval (a,b)
% with Ridder’s method.
% USAGE: roots(func,a,b,dx,tol)
% func = handle of function that returns f(x)
% dx = increment of x used in root search
% tol = error tolerance (default is 10.e-6)

102 PROBLEM SET 4.1

if nargin < 5; tol = 1.0e-6; end
fprintf(’Roots:\n’)
while 1
[x1,x2] = rootsearch(func,a,b,dx);
if isnan(x1)
fprintf(’Done’); break
else
a = x2;
x = ridder(func,x1,x2,tol);
if isnan(x); continue
else fprintf(’%16.6e\n’, x); end
end
end
The roots are now obtained with the program
% problem4_1_10
func = inline(’x*sin(x) + 3*cos(x) - x’,’x’);
rootsRidder(func,-6,6,0.5)

Roots:
-4.712389e+000
-3.208839e+000
1.570796e+000
Done
Note the use of MATLAB’s in-line function, which is passed to rootsRidder.
An in-line function can be evaluated by feval in the same manner as a function
stored in a M-…le. The advantage of an in-line function is that it does not
create a new M-…le.

Problem 11

The algorithm listed below is similar to rootsRidder in Problem 10.

function rootsNewton(func,dfunc,a,b,dx,tol)
% Computes all the roots of f(x) in the interval (a,b)
% with the Newton-Raphson method.
% USAGE: rootsNewton(func,dfunc,a,b,dx,tol)
% func = handle of function that returns f(x)
% dfunc = handle of function that returns f’(x)
% dx = increment of x used in root search
% tol = error tolerance (default is 10.e-6)

PROBLEM 11 103
if nargin < 6; tol = 1.0e-6; end
fprintf(’Roots:\n’)
while 1
[x1,x2] = rootsearch(func,a,b,dx);
if isnan(x1)
fprintf(’Done’); break
else
a = x2;
x = newtonRaphson(func,dfunc,x1,x2,tol);
if isnan(x); continue
else fprintf(’%16.6e\n’, x); end
end
end

% problem4_1_11
func = inline(’x*sin(x) + 3*cos(x) - x’,’x’);
dfunc = inline(’x*cos(x) - 2*sin(x) - 1’,’x’);
rootsNewton(func,dfunc,-6,6,0.5)

>> Roots:
-4.712389e+000
-3.208839e+000
1.570796e+000
Done

Problem 12

f (x) = x4 + 0:9x3 2:3x2 + 3:6x 25:2

f 0 (x) = 4x3 + 2:7x2 4:6x + 3:6

104 PROBLEM SET 4.1

Whenever possible, the function should plotted in order to gain information
about its behaviour and locate its zeros.

y 250

200

150

100

50

-4 -3 -2 -1 1 2 3 4
x

From the plot of f (x) we see that there are two roots, located in ( 3:2; 2:4). As
the derivative of the function is easily obtained, we use the Newton-Raphson
method due to its superior covergence. The program below calls rootsNewton
listed in Problem 11.

% problem4_1_12
func = inline(’x^4 + 0.9*x^3 - 2.3*x^2 + 3.6*x - 25.2’,’x’);
dfunc = inline(’4*x^3 + 2.7*x^2 - 4.6*x + 3.6’,’x’);
rootsNewton(func,dfunc,-3.2,2.4,2)

>> Roots:
-3.000000e+000
2.100000e+000
Done

Problem 13

f (x) = x4 + 2x3 7x2 + 3

f 0 (x) = 4x3 + 6x2 14x

PROBLEM 13 105
 y
6

0
0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0
x

By inspection of the plot of f (x) we see that the two positive roots are located
in (0:7; 1:7). Again we compute these roots with the function rootsNewton in
Problem 11.

% problem4_1_13
func = inline(’x^4 + 2*x^3 - 7*x^2 + 3’,’x’);
dfunc = inline(’4*x^3 + 6*x^2 - 14*x’,’x’);
rootsNewton(func,dfunc,0.7,1.7,0.5)

>> Roots:
7.912878e-001
1.618034e+000
Done

Problem 14

f (x) = sin x 0:1x

y
0.5

0.0
2 4 6 8 10 12 14
x
-0.5

-1.0

-1.5

-2.0

106 PROBLEM SET 4.1

The plot shows that all the positive nonzero roots are in the interval (2; 9).
Here we chose to compute the roots with rootsRidder listed in Problem 10.
% problem4_1_14
func = inline(’sin(x) - 0.1*x’,’x’);
rootsRidder(func,2,9,0.5)

>> Roots:
2.852342e+000
7.068174e+000
8.423204e+000
Done

Problem 15

f ( ) = cosh cos + 1
f 0 ( ) = sinh x cos x cosh x sin x

y 20

10

0
1 2 3 4 5
x
-10

The plot of f ( ) reveals that the roots lie in (1:8; 2:0) and (4:6; 4:8). The
following program uses newtonRaphson to compute these roots:
% problem4_1_15
func = inline(’cosh(x)*cos(x) + 1’,’x’);
dfunc = inline(’sinh(x)*cos(x) - cosh(x)*sin(x)’,’x’);
b = 0.025; h = 0.0025; rho = 7850; E = 200e9; L = 0.9;
I = b*h^3/12; m = rho*b*h*L;
C = sqrt(E*I/(m*L^3))/(2*pi);
bracket = [1.8 2; 4.6 4.8];
for i = 1:2
beta = newtonRaphson(func,dfunc,bracket(i,1),bracket(i,2));

PROBLEM 15 107
 freq = C*beta^2
end

>> freq =
2.5166
freq =
15.7713

Problem 16

1 s 1
f( ) = sinh = sinh 1:1
L
1 1
f 0( ) = cosh 2 sinh

y
0.05

0.00
0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0
x

-0.05

-0.10

According to the plot, the smallest positive root lies in (0:72; 0:80). We use
newtonRaphson to compute this root.
% problem4_1_16
func = inline(’sinh(x)/x - 1.1’,’x’);
dfunc = inline(’cosh(x)/x - sinh(x)/x^2’,’x’);
beta = newtonRaphson(func,dfunc,0.72,0.8)
gamma = 77000; L = 1000; s = 1100;
sigma0 = gamma*L/(2*beta);
max_stress = sigma0*cosh(beta)

>> beta =
0.7634
max_stress =
6.5855e+007

108 PROBLEM SET 4.1

Problem 17

We non-dimensionalize the secant formula by dividing both sides by E:

" r !#
ec L max
f = 1 + 2 sec
E E r 2r E E

Substituting

ec 85(170) L 7100
= = 0:7166 = = 25:0
r 2
(142)2 2r 2(142)
max 120 106 3
= = 1:6901 10
E 71 109
and using the notation u = =E, the secant formula is
p 3
f (u) = u 1 + 0:7166 sec 25 u 1:6901 10

0.02
y

0.01

0.00
0.001 0.002 0.003 0.004 0.005 0.006 0.007 0.008 0.009 0.010
x
-0.01

-0.02

The plot of f (u) shows that the smallest root is in the interval (0:0004; 0:0012).
We used Ridder’s method to compute this root:

% problem4_1_17
func = inline(’u*(1 + 0.7166/cos(25*sqrt(u))) - 1.6901e-3’,’u’);
u = ridder(func, 0.0004, 0.0012)

>> u =
8.6032e-004

P = A = AEu = (25 800 10 6 )(71 109 )(8:6032 10 4 )

= 1:576 106 N J

PROBLEM 17 109
Problem 18

Dividing both sides of the Bernoulli equation

Q2 Q2
+ h0 = +h+H
2gb2 h20 2gb2 h2
by h0 , we get
2
Q2 Q2 h0 h H
2 3
+ 1 = 2 3
+ +
2gb h0 2gb h0 h h0 h0
Introducing u = h=h0 and rearranging, this becomes
Q2 1 H
f (u) = 1 u+ 1 =0
2gb2 h30 u2 h0
Substituting
Q2 (1:2)2
= = 0:10487
2gb2 h30 2(9:81)(1:8)2 (0:6)3
H 0:075
1 = 1 = 0:875
h0 0:6
we obtain
1
f (u) = 0:10487 1 u + 0:875
u2
0.5
y

0.0
0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0
x

-0.5

-1.0

The plot of f (u) indicates that there are two roots. The smaller root, which is
in (0:4; 0:48), can be determined with the following program:
problem4_1_18
func = inline(’0.10487*(1 - 1/u^2) - u + 0.875’,’u’);
u = ridder(func, 0.4, 0.48)

>> u =
0.4412

110 PROBLEM SET 4.1

 h = uh0 = 0:4412(0:6) = 0:2647 m J
The larger root can be computed with the same program by changing the
brackets to (0:8; 0:84). It yields u = 0:8263, so that

h = 0:8263(0:6) = 0:4958 m J

Evidently the ‡uid ‡ow can exist in on of the two states under the given
conditions.

Problem 19
M0
v = u ln gt
M0 mt
_
We want the root of

1
f (t) = u ln gt vsound = 0
1 (m=M
_ 0 )t

where
m_ 13:3 103 1
= = 0:004 75 s
M0 2:8 106
Thus
1
f (t) = 2510 ln 9:81t 335
1 0:004 75t

300
y
200

100

0
10 20 30 40 50 60 70 80 90 100
x
-100

-200

-300

The plot of f (t) locates the root in (68 ; 76) s. The following program was used
for the computation of the root:

PROBLEM 19 111
% problem4_1_19
func = inline(’2510*log(1/(1 - 4.75e-3*t)) - 9.81*t - 335’,’t’);
t = ridder(func, 68, 76)

>> t =
70.8780

Problem 20

ln(T2 =T1 ) (1 T1 =T2 )

=
ln(T2 =T1 ) + (1 T1 =T2 )=( 1)

With = 0:3, 1 = 2=3 and the notation u = T2 =T1 , the equation becomes

ln u (1 1=u)
f (u) = 0:3 = 0
ln u + 1:5(1 1=u)

0.08
y
0.06
0.04
0.02
0.00
-0.02 3 4 5 6 7 8 9 10
x
-0.04
-0.06
-0.08
-0.10
-0.12
-0.14
-0.16

From the plot of f (u) we see that the root is in (5:2; 5:6). We found this root
with the following program:

% problem4_1_20
func = inline(’(log(u)-(1-1/u))/(log(u)+1.5*(1-1/u))-0.3’...
,’u’);
u = ridder(func,5.2,5.6)

>> u =
5.4125

112 PROBLEM SET 4.1

Problem 21

G= RT ln (T =T0 )5=2

5
f (T ) = G + RT ln(T =T0 )
2
Substituting
5 5
R = (8:314 41) = 20:7860
2 2
we get
f (T ) = 105 + 20:7860T ln(T =4:444 18)

x
y 100 200 300 400 500 600 700 800 900 1000
0

-20000

-40000

-60000

-80000

-1e+5

The plot of f (T ) shows a root in (880; 920). Here is the program that computes
this root:

% problem4_1_21
func = inline(’-1.0e5 + 20.7860*T*log(T/4.44418)’,’T’);
T = ridder(func,880,920)

>> T =
904.9435

Problem 22

(3 2 )2
f( ) = 249:2
(1 )3

PROBLEM 21 113
 y
400

200

0
0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0
x
-200

-400

The plot of f ( ) shows a root in (0:7; 0:9), which is computed with the program

% problem4_1_22
func = inline(’x*(3 - 2*x)^2/(1 - x)^3 - 249.2’,’x’);
xi = ridder(func,0.7,0.9)

>> xi =
0.8171

Problem 23

f1 (x; y) = (x 2)2 + y 2 4
f2 (x; y) = x2 + (y 3)2 4

3
2
x
2

The rough locations of the intersection points are (2; 2) and (0; 1). Letting
x = x1 and y = x2 , the following function de…ned the equations:

function y = p4_1_23(x)

114 PROBLEM SET 4.1

% Equations used in Problem 23, Problem Set 4.1
y = [(x(1) - 2)^2 + x(2)^2 - 4;
x(1)^2 + (x(2) - 3)^2 - 4];
The following command returns the coordinates of the …rst point:
>> newtonRaphson2(@p4_1_23,[2;2])
ans =
1.7206
1.9804
Changing the starting point to [0; 1], we obtain the coordinates of the second
point
>> newtonRaphson2(@p4_1_23,[0;1])
ans =
0.2794
1.0196

Problem 24

f1 (x; y) = sin x + 3 cos x 2

f2 (x; y) = cos x sin y + 0:2
The following function uses the notation x = x1 and y = x2 :
function y = p4_1_24(x)
% Equations used in Problem 24, Problem Set 4.1
y = [sin(x(1)) + 3*cos(x(2)) - 2;
cos(x(1)) - sin(x(2)) + 0.2];
The x and y-coordinates can now be obtained with the command
>> newtonRaphson2(@p4_1_24,[1;1])
ans =
0.7912
1.1267

Problem 25

tan x y = 1
cos x 3 sin y = 0

PROBLEM 24 115
It is not easy to search for the roots of simultaneous equations. Here we can
overcome the di¢ culty by solving the …rst equation for y and substituting
the result into the second equation. This gives us the single transcendental
equation
f (x) = cos x 3 sin(tan x 1) = 0

y
3

0
0.2 0.4 0.6 0.8 1.0 1.2 1.4
x
-1

-2

From the plot of f (x) we see that there are 5 roots in the interval (0; 1:5). The
…rst root is about x = 0:88; the other roots are closely spaced near the end
of the interval (the spacing of roots becomes in…nitesimal at x = =2). The
program listed below is based on rootsRidder in Problem 10. It searches for
the roots from x = 0:8 to 1:5 in increments of 0:025 (the increment has to be
small in order to catch all the roots).

% problem4_1_25
func = inline(’cos(x) - 3*sin(tan(x) - 1)’,’x’);
n = 5; a = 0.8; b = 1.5; dx = 0.025;
fprintf(’Roots:\n’)
while 1
[x1,x2] = rootsearch(func,a,b,dx);
if isnan(x1)
fprintf(’Done’); break
else
a = x2;
x = ridder(func,x1,x2);
if isnan(x); continue
else
y = tan(x) - 1;
fprintf(’%12.6f %12.6f\n’, x,y);
end
end
end

116 PROBLEM SET 4.1

>> Roots:
0.881593 0.213595
1.329402 3.061823
1.435176 6.328269
1.474872 9.392847
1.497350 12.590833
Done

Problem 26

(x a)2 + (y b)2 R2 = 0
Substituting the coordinates of the given points into the above equation, we
get
(8:21 a)2 + b2 R2 = 0
(0:34 a)2 + (6:62 b)2 R2 = 0
(5:96 a)2 + ( 1:12 b)2 R2 = 0
By plotting the points, we can estimate the parameters of the circle. It appears
that reasonable starting values are a = 5, b = 3 and R = 5.
y
6

4
Center of circle
2

0 x
0 2 4 6 8
-2
The following function uses the notation a = x1 , b = x2 , R = x3 :
function y = p4_1_26(x)
% Equations used in Problem 26, Problem Set 4.1
y = [(8.21 - x(1))^2 + x(2)^2 - x(3)^2;
(0.34 - x(1))^2 + (6.62 - x(2))^2 - x(3)^2;
(5.96 - x(1))^2 + (-1.12 - x(2))^2 - x(3)^2];

>> newtonRaphson2(@p4_1_26,[5;3;5])
ans =
4.8301
3.9699
5.2138

PROBLEM 26 117
Problem 27

C
R=0
1 + e sin( + )

After substituting the three sets of given data, we obtain the simulataneous
equations

C
6870 = 0
1 + e sin( =6 + )
C
6728 = 0
1 + e sin( )
C
6615 = 0
1 + e sin( =6 + )

The starting value C = 6800 seems reasonable, but e and are not easy to
guess. The orbit has some eccentricity, so that e = 0:5 should not be out of
line (e = 0 will not work because it results in a singular Jacobian matrix). We
also used = 0, which was a pure guess.
The minimum value of R is
C
Rmin =
1+e

occuring at

sin( + ) = 1 =
2

With the notation C = x1 , e = x2 and = x3 , we arrive at the following

program:

function y = p4_1_27(x)
% Equations used in Problem 27, Problem Set 4.1
y = [x(1)/(1 + x(2)*sin(-pi/6 + x(3))) - 6870;
x(1)/(1 + x(2)*sin(x(3))) - 6728;
x(1)/(1 + x(2)*sin(pi/6 + x(3))) - 6615];

>> format short e

>> newtonRaphson2(@p4_1_27,[6800;0.5;0])
ans =
6.8193e+003
4.0599e-002
3.4078e-001

118 PROBLEM SET 4.1

Problem 28

x = (v cos )t
1 2
y = gt + (v sin )t
2

We also need the expression for dy=dx:

dx dy
= v cos = gt + v sin
dt dt
dy gt + v sin
=
dx v cos
Letting denote the time of ‡ight, the speci…ed end conditions are

dy
x( ) = 300 m y( ) = 61 m = 1
dx

which yield the equations

(v cos ) 300 = 0
1 2
g + (v sin ) 61 = 0
2
g + v sin
+1 = 0
v cos
To estimate the initial values of the unknowns, we guess = 10 s and = 45 .
Then the …rst of the above equations yields v = 300= ( cos ) 57 m/s.
The following program uses the notation = x1 , v = x2 and = x3 :

function y = p4_1_28(x)
% Equations used in Problem 27, Problem Set 4.1
y = [x(1)*x(2)*cos(x(3)) - 300;
x(1)*x(2)*sin(x(3)) - 9.81*x(1)^2/2 - 61;
(-9.81*x(1) + x(2)*sin(x(3)))/(x(2)*cos(x(3))) + 1];

>> format short e

>> newtonRaphson2(@p4_1_28,[10;57;pi/4])
ans =
8.5789e+000
6.0353e+001
9.5279e-001

Thus = 8:579 s, v = 60:35 m/s and = 0:9528 rad = 54:6 deg. J

PROBLEM 28 119
Problem 29

150 cos 1+ 180 cos 2 200 cos 3 = 200

150 sin 1 + 180 sin 2 200 sin 3 = 0

mm
180
y θ2
y

mm
m

200
0m
15

θ1 200 mm θ3
x x

(a) (b)

Here is the function that re…nes de…nes the equations, given that 3 = 75 =
5 =12 rad:

function y = p4_1_29(x)
% Equations used in Problem 29, Problem Set 4.1
y = [150*cos(x(1)) + 180*cos(x(2)) - 200*cos(5*pi/12) - 200;
150*sin(x(1)) + 180*sin(x(2)) - 200*sin(5*pi/12)];

We estimate from the …gure that 1 = 60 and 2 = 20 in con…guration (a).

Using these as starting values, the solution is obtained with the command

>> newtonRaphson2(@p4_1_29,[pi/3;pi/9])
ans =
9.5960e-001
4.0148e-001

Therefore, 1 = 0:9596 rad = 55:0 , 2 = 0:4015 rad = 23:0 J

For con…guration (b) we estimate 1 = 20 and 2 = 60 . With these starting
values we have

>> newtonRaphson2(@p4_1_29,[pi/9;pi/3])

ans =
3.4939e-001
9.0752e-001

Here 1 = 0:3494 rad = 20:0 , 2 = 0:9075 rad = 52:0 J

120 PROBLEM SET 4.1

Problem 30
T
Letting x = 1 2 3 T , the equations to be solved are

function y = p4_1_30(x)
% Equations used in Prob. 30, Problem Set 4.1
y = [x(4)*(-tan(x(2)) + tan(x(1))) - 16;
x(4)*( tan(x(3)) + tan(x(2))) - 20;
-4*sin(x(1)) - 6*sin(x(2)) + 5*sin(x(3)) + 3;
4*cos(x(1)) + 6*cos(x(2)) + 5*cos(x(3)) - 12];

Rough estimates (starting values) of the variables are

1 = 0:8 rad 2 = 0:3 rad 3 = 0:4 rad T = 20 kN

The solution is obtained with the command

>> newtonRaphson2(@p4_1_30,[0.8;0.3;0.4;20])

ans =
0.9358
0.4334
0.5800
17.8884

Therefore, the solution is

1 = 0:9358 rad = 53:62 J

2 = 0:4334 rad = 24:83 J
3 = 0:5800 rad = 33:23 J
T = 17:89 kN

PROBLEM 30 121
122 PROBLEM SET 4.1
PROBLEM SET 4.2

Problem 1

P3 (x) = 3x3 + 7x2 36x + 20 r= 5

b 1 = a1 = 3
b2 = a2 + rb1 = 7 + ( 5)(3) = 8
b3 = a3 + rb2 = 36 + ( 5)( 8) = 4

) P2 = 3x2 8x + 4 J

Problem 2

P4 (x) = x4 3x2 + 3x 1 r=1

b1 = a1 = 1
b2 = a2 + rb1 = 0 + 1(1) = 1
b3 = a3 + rb2 = 3 + 1(1) = 2
b4 = a4 + rb3 = 3 + 1( 2) = 1

) P3 = x 3 + x 2 2x + 1 J

Problem 3

P5 (x) = x5 30x4 + 361x3 2178x2 + 6588x 7992 r=6

b1 = a1 = 1
b2 = a2 + rb1 = 30 + 6(1) = 24
b3 = a3 + rb2 = 361 + 6( 24) = 217
b4 = a4 + rb3 = 2178 + 6(217) = 876
b5 = a5 + rb4 = 6588 + 6( 876) = 1332

P4 (x) = x4 24x3 + 217x2 876x + 1332

PROBLEM SET 4.2 123

Problem 4

P4 (x) = x4 5x3 2x2 20x 24 r = 2i

b1 = a1 = 1
b2 = a2 + rb1 = 5 + (2i)(1) = 5 + 2i
b3 = a3 + rb2 = 2 + (2i)( 5 + 2i) = 6 10i
b4 = a4 + rb3 = 20 + (2i)( 6 10i) = 12i

P3 (x) = x3 (5 2i)x2 (6 + 10i)x 12i J

Problem 5

P3 (x) = 3x3 19x2 + 45x 13 r=3 2i

b 1 = a1 = 3
b2 = a2 + rb1 = 19 + (3 2i)(3) = 10 6i
b3 = a3 + rb2 = 45 + (3 2i)( 10 6i) = 3 + 2i

P2 (x) = 3x2 (10 + 6i) x + (3 + 2i) J

Problem 6

P3 (x) = x3 + 1:8x2 9:01x 13:398 r= 3:3

b1 = a1 = 1
b2 = a2 + rb1 = 1:8 + ( 3:3)(1) = 1:5
b3 = a3 + rb2 = 9:01 + ( 3:3)( 1:5) = 4:06

P2 (x) = x2 1:5x 4:06

The roots are
p
1:5 1:52 + 4(4:06) 1:5 4:3 2: 9
r= = = J
2 2 1:4

124 PROBLEM SET 4.2

Problem 7

P3 (x) = x3 6:64x2 + 16:84x 8:32 r = 0:64

b 1 = a1 = 1
b2 = a2 + rb1 = 6:64 + 0:64(1) = 6
b3 = a3 + rb2 = 16:84 + 0:64( 6:0) = 13

P2 (x) = x2 6x + 13

The roots are

p
6 62 4(13) 6 4i 3 + 2i
r= = = J
2 2 3 2i

Problem 8

P3 (x) = 2x3 13x2 + 32x 13 r=3 2i

b 1 = a1 = 2
b2 = a2 + rb1 = 13 + (3 2i) (2) = 7 4i
b3 = a3 + rb2 = 32 + (3 2i) ( 7 4i) = 3 + 2i

P2 (x) = 2x2 (7 + 4i)x + (3 + 2i)

Since complex roots come in conjugate pairs, we know that a zero of P2 (x) is

r = 3 + 2i J

b 1 = a1 = 2
b2 = a2 + rb1 = (7 + 4i) + (3 + 2i)(2) = 1

P1 (x) = 1 + 2x

The zero of P1 (x) is

r = 0:5 J

PROBLEM 7 125
Problem 9

P4 (x) = x4 3x3 + 10x2 6x 20 r = 1 + 3i

b1 = a1 = 1
b2 = a2 + rb1 = 3 + (1 + 3i)(1) = 2 + 3i
b3 = a3 + rb2 = 10 + (1 + 3i)( 2 + 3i) = 1 3i
b4 = a4 + rb3 = 6 + (1 + 3i)( 1 3i) = 2 6i

P3 (x) = x3 + ( 2 + 3i)x2 + ( 1 3i)x + (2 6i)

Another zero of P4 (x) is the conjugate of 1 + 3i, namely

r=1 3i

b1 = a1 = 1
b2 = a2 + rb1 = ( 2 + 3i) + (1 3i)(1) = 1
b3 = a3 + rb2 = ( 1 3i) + (1 3i)( 1) = 2

P2 (x) = x2 x+ 2

The roots of the quadratic are

p
1 12 + 4(2) 2
r= =
2 1

Thus the roots of P4 (x) are 1 3i, 2 and 1. J

Problem 10

>> polyroots([1 2.1 -2.52 2.1 -3.52])

ans =
-0.0000 - 1.0000i
1.1000
0.0000 + 1.0000i
-3.2000

126 PROBLEM SET 4.2

Problem 11

>> polyroots([1 -156 -5 780 4 -624])

ans =
1.0000
-1.0000
-2.0000
2.0000
156.0000

Problem 12

>> polyroots([1 4 -8 -34 57 130 -150])

ans =
2.0000 - 1.0000i
1.0000
-3.0000
2.0000 + 1.0000i
-3.0000 - 1.0000i
-3.0000 + 1.0000i

Problem 13

>> polyroots([8 28 34 -13 -124 19 220 -100])

ans =
-2.0000
0.5000
-2.0000
1.0000 - 0.5000i
1.0000 + 0.5000i
-1.0000 - 2.0000i
-1.0000 + 2.0000i

PROBLEM 11 127
Problem 14

>> polyroots([1 -7 7 25 24 -98 -472 440 800])

ans =
2.0000
-1.0000
-2.0000
3.0000 - 1.0000i
-1.0000 + 2.0000i
3.0000 + 1.0000i
-1.0000 - 2.0000i
4.0000

Problem 15

>> polyroots([1 5+1i -8+5i 30-14i -84])

ans =
2.0000
-0.0000 + 2.0000i
0.0000 - 3.0000i
-7.0000

Note that the complex roots do not appear in conjugate pairs if the coe¢ cients
of the polynomial are not real.

Problem 16
2
c 3 k c k k
!4 + 2 ! + 3 !2 + !+ =0
m m mm m
1
With c=m = 12 s and k=m = 1500 s 2 , we get

! 4 + 24! 3 + 4500! 2 + 18 103 ! + 2:25 106 = 0

> polyroots([1 24 4500 18e3 2.25e6])

ans =
-0.6230 -24.0302i
-0.6230 +24.0302i
-11.3770 +61.3545i
-11.3770 -61.3545i

128 PROBLEM SET 4.2

The two combinations of (! r ; ! i ) are

0:0623 s 1 ; 24:03 s 1
and ( 11:38 s 1 ; 61:35 s 1 ) J

Problem 17

The slope of the beam is

w0
y0 = (5x4 9L2 x2 + 6L3 x)
120EI
w0 L4
= (5 4 9 2
+6 )
120EI
where = x=L. Since y 0 = 0 at the point of maximum displacement, the value
of that we are looking for is a root of
4 2
P4 ( ) = 5 9 +6

We could …nd the our roots of this equation with the function polyroots, but
this is unnecessary. Because the slope of the beam is zero at supports, we know
that two of the roots are = 0 and = 1. Factoring out these roots, we have
2
P4 ( ) = ( 1)(b1 + b2 + b3 )

The b’s are obtained by Horner’s algorithm:

b 1 = a1 = 5
b2 = a2 + rb1 = 0 + 1(5) = 5
b3 = a3 + rb2 = 9 + 1(5) = 4

We have now reduced the problem to …nding the roots of the quadratic equation
2
5 +5 4=0

The positive root is

p
5+ 52 4(5)( 4)
= = 0:5247 J
10

PROBLEM 17 129