• Author / Uploaded
• Sidd265

Citation preview

Solution to Exercises Chapter 1 1.1 E1 = 1+5^3+3^3-153 E2 = 1+6^4+3^4+4^4-1634 E3 = 5^6+4^6+8^6+8^6+3^6+4^6-548834 E4 = 1+factorial(7)-71^2 E5 = factorial(1)+factorial(4)+factorial(5)+factorial(6)+factorial(7)+factorial(8)-215^2

1.2 pie = pi - (100-(2125^3+214^3+30^3+37^2)/82^5)^(1/4)

1.3 format long g I1 = 53453/log(53453) I2 = 613*exp(1)/37-35/991 format short

1.4 x = (49-27*sqrt(5)+3*sqrt(6)*sqrt(93-49*sqrt(5)))^(1/3) R = 0.5*sqrt((8*2^(2/3)-16*x+2^(1/3)*x^2)/(8*2^(2/3)-10*x+2^(1/3)*x^2)) Answers: x = 2.5375 + 1.9080i

R = 0.8161 - 0.0000i

1.5 E1 = cot(pi/5)-sqrt(25+10*sqrt(5))/5 E2 = sin(pi/15)-sqrt(7-sqrt(5)-sqrt(30-6*sqrt(5)))/4 E3 = pi-16*atan(1/5)+4*atan(1/239)

1.6 g = sqrt(17-sqrt(17)); b = sqrt(17+sqrt(17)); a = sqrt(34+6*sqrt(17)+sqrt(2)*(sqrt(17)-1)*g-8*b*sqrt(2)); E1 = sin(pi/17)-sqrt(2)/8*sqrt(g^2-sqrt(2)*(a+g))

1.7 a = 1; b = 2; E1 = exp(a)-sinh(a)-cosh(a) E2 = sinh(2*b)/(cosh(a+b)*cosh(a-b))-tanh(a+b)+tanh(a-b)

1.8 L = 1.5; h = 1; r = 1.6; V = L*(r^2*acos((r-h)/r)-(r-h)*sqrt(2*r*h-h^2)) © 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a1retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

Answer: V = 3.2209

1.9 g = 60*pi/180; a = 35*pi/180; n = 4/3; D = a-g+asin(n*sin(g-asin(sin(a)/n))) Answer: D = 0.4203

1.10 k = 7; n =12; Cnk = factorial(n)/factorial(k)/factorial(n-k) Answer: Cnk = 792

1.11 r = 2.5; I = (pi/8-8/9/pi)*r^4 Answer: I = 4.2874

1.12 c = 5; K = (4*c-1)/(4*c+4)+0.615/c Answer: K = 0.9147

1.13 B = 0.6; K = 3/(1-B)^3*(0.5-2*B+B*(1.5-log(B))) Answer: K =

23.7420

1.14 R = 30; r = 12; S = 50; theta = asin((R-r)/S); L = 2*S*cos(theta)+pi*(R+r)+2*theta*(R-r) Answer: L = 238.4998

1.15 Fn = 250; f = 0.35; r = 0.4; theta = 60*pi/180; T = (4*f*r*Fn*sin(theta/2))/(theta+sin(theta)) Answer: T = 36.5875 © 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a2retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

1.16 A = 1.7; B = 1.2; D = 1.265*((A*B)^3/(A+B))^(1/5) Answer: D = 1.5682

1.17 n = 6; M = 1/sin(pi/n); C = (1+M^2)/4/M; alpha = acos(sqrt(C^2+2)-C); aG = M*(1-M^2)*sin(alpha)/(1+M^2-2*M*cos(alpha))^2 Answer: aG = -1.3496

1.18 L = 3000; d = 45; V = 1600; Deltap = 0.03*L/d^1.24*(V/1000)^1.84 Answer: Deltap = 1.9048

1.19 E1 = 206e9; E2 = 206e9; nu1 = 0.3; nu2 = 0.3; d1 = 0.038; d2 = 0.070; F = 450; z = 0.00025; a = (0.375*F*((1-nu1^2)/E1+(1-nu2^2)/E2)/(1/d1+1/d2))^(1/3); pmax = 1.5*F/pi/a^2; za = z/a; sxsy = -pmax*((1-za*atan(1/za))*(1-nu1)-0.5/(1+za^2)) sz = -pmax/(1+za^2) Answers: sxsy = 2.0779e+008 sz = -1.2421e+009

1.20 E1 = 206e9; E2 = 206e9; nu1 = 0.3; nu2 = 0.3; d1 = 0.038; d2 = 0.070; F = 450; z = 0.000025; L = 0.05; b = sqrt(2*F/pi/L*((1-nu1^2)/E1+(1-nu2^2)/E2 )/(1/d1+1/d2)); pmax = 2*F/pi/b/L; zb = z/b; sx = -pmax*2*nu2*(sqrt(1+zb^2)-zb) sy = -pmax*((2-1/(1+zb^2))*sqrt(1+zb^2)-2*zb) sz = -pmax/sqrt((1+zb^2)) Answers: sx = -5.0360e+007 sy = -3.5432e+007 sz = -1.3243e+008

1.21 © 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a3retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

e = 0.8; NL = pi*e*sqrt(pi^2*(1-e^2)+16*e^2)/(1-e^2)^2 Answer: NL = 72.0220

1.22 h = 0.03; d0 = 0.006; d1 = 0.016; E = 206e9; d2 = d1+h*tan(pi/6); k = pi*E*d0*tan(pi/6)/log((d2-d0)*(d1+d0)/(d2+d0)/(d1-d0)) Answer: k = 5.2831e+009

1.23 alpha = 2e-5; E = 206e9; nu = 0.3; Ta = 260; Tb = 150; a = 0.006; b = 0.012; r = 0.010; c1 = alpha*E*(Ta-Tb)/2/(1-nu)/log(b/a); c2 = a^2/(b^2-a^2)*log(b/a); br = b/r; sr = c1*(c2*(br^2-1)-log(br)) st = c1*(1-c2*(br^2+1)-log(br)) T = Tb+(Ta-Tb)*log(br)/log(b/a) Answers: sr = -3.7670e+007 st = 1.1859e+008 T = 178.9338

1.24 k = 1.4; pepo = 0.3; psi = sqrt(k/(k-1))*sqrt(pepo^(2/k)-pepo^((k+1)/k)) Answer: psi = 0.4271

1.25 x = 0.45; K = 1.2/x/(sqrt(16*x^2+1)+log(sqrt(16*x^2+1)+4*x)/4/x)^(2/3) Answer: K = 1.3394

1.26 c = sqrt(8)/9801; n = 0; oneopi0 = c*(1103+26390*n)*factorial(4*n)/factorial(n)^4/396^(4*n); n = 1; oneopi1 = c*(1103+26390*n)*factorial(4*n)/factorial(n)^4/396^(4*n); format long e pin0 = 1/oneopi0 diffpino = pin0-pi pin0n1 = 1/(oneopi0+oneopi1) diffpin0n1 = pin0n1-pi © 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a4retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

format short Answers: pin0 = 3.141592730013306e+000 diffpino = 7.642351240733092e-008 pin0n1 = 3.141592653589794e+000 diffpin0n1 = 4.440892098500626e-016

1.27 r = 10; rc = 3; k = 1.4; eta = 1-1/r^(k-1)*((rc^k-1)/k/(rc-1)) Answer: eta = 0.4803

1.28 M = 2; k = 1.4; AA = 1/M*(2/(k+1)*(1+(k-1)/2*M^2))^((k+1)/(2*(k-1))) Answer: AA = 1.6875

1.29 n = 5; ep = 0.1; w = 0.5 T = 1/sqrt(1+(ep*cos(n*acos(w))).^2) w=1 T = 1/sqrt(1+(ep*cos(n*acos(w))).^2) w = 1.5 T = 1/sqrt(1+(ep*cosh(n*acosh(w))).^2) Answers: w = 0.5000 T = 0.9988 w = 1 T = 0.9950 w = 1.5000 T = 0.1605

1.30 m = 1; n = 2; nu = 0.3; hR = 0.05; Rl = 0.1; Lm = pi/4*Rl*(4*m+1); Om = (1-nu^2)*Lm^4/(Lm^2+n^2+1.78*n^2*Lm^2); Om = Om+hR^2/12*(Lm^4+n^4+1.78*m^2*Lm^2) Answer: Om = 7.5159e-003

1.31 V = 2.9; D = 0.3; ep = 0.00025; nu = 1.012e-6; Re = V*D/nu; f = ((64/Re)^8+9.5*(log(ep/3.7/D+5.74/Re^0.9)-(2500/Re)^6)^-16)^(1/8) Answer: f = 0.0193

1.32 n = 3; th = pi/7; a = cos(th)+(n^2*cos(2*th)+sin(th)^4)/(n^2-sin(th)^2)^(3/2) © 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a5retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

Answer: a = 1.1168

1.33 [Note: These solutions are obtained with a large amount of user interaction, almost on a line by line basis.] (a) x = a*sin(t); y = a*cos(t)^2*(2+cos(t))/(3+sin(t)^2); xp = diff(x,t,1); xpp = diff(x,t,2); yp = diff(y,t,1); ypp = diff(y,t,2); [Nnum Nden] = numden(xp*ypp-yp*xpp); Nden = factor(Nden) Nnum = simple(Nnum) [Dnum Dden] = numden(simple(xp^2+yp^2)) Partial answers: Nden = (sin(t)^2 + 3)^3 Nnum = -a^2*cos(t)^2*(cos(t) + 2)^3*(9*cos(t) - 6) Dnum = a^2*(cos(2*t) + 1)*(9*cos(2*t) - 80*cos(t) + 73) Dden = 4*(cos(t) - 2)^4 At this point, we place these results in standard notation and make a few trigonometric substitutions as follows: 3'2

3

! =

a 2 cos 2 t (cos t + 2 ) (9cos(t ) 6 )(4(cos(t ) 2) 4 )

(sin

2

3

3/ 2

t + 3) (a 2 (cos(2t ) + 1)(9cos(2t ) 80cos(t ) + 73)) 3

=

6

24a 2 cos 2 t (cos t + 2 ) (3cos(t ) 2 )(cos(t ) 2 )

(4

3

3/ 2

cos 2 t ) (2a 2 cos 2 t )

3/ 2

(9cos(2t )

80cos(t ) + 73) 3

=

3

24a 2 cos 2 t #%(cos t + 2 )(cos(t ) 2 )& (3cos(t ) 2 )(cos(t ) 2 ) 3 3/ 2 23 / 2 a 3 cos3 t (4 cos 2 t ) (9cos(2t ) 80cos(t ) + 73) 3

3

=

(4 cos(t ) ) (3cos(t ) 2 )(cos(t ) 2 ) 6 2 a cos t (4 cos 2 t )3 (9cos(2t ) 80cos(t ) + 73)3 / 2

=

6 2 sec t (3cos(t ) 2 )(cos(t ) 2 ) a (9cos(2t ) 80cos(t ) + 73)3 / 2

3

(b) syms a t x = 3*t*a/(1+t^3); y = x*t; xp = diff(x,t,1); xpp = diff(x,t,2); yp = diff(y,t,1); ypp = diff(y,t,2); kn = simple(xp*ypp-yp*xpp); kd = simple(xp^2+yp^2); k = simple(kn/kd^(3/2)) Answer: k = (2*a^2)/(3*(t^3 + 1)^2*((a^2*(t^8 + 4*t^6 - 4*t^5 - 4*t^3 + 4*t^2 + 1))/(t^3 + 1)^4)^(3/2))

© 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a6retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

(c) syms a t x = a*(t-tanh(t)); y = a*sech(t); xp = diff(x,t,1); xpp = diff(x,t,2); yp = diff(y,t,1); ypp = diff(y,t,2); kn = simple(xp*ypp-yp*xpp); kn = subs(kn, cosh(t)^2 - 1, sinh(t)^2); kd = collect(simple(xp^2+yp^2), a); kd = subs(kd,1 - 1/cosh(t)^2, sinh(t)^2/cosh(t)^2); kd = subs(kd, (1/cosh(t)^2)-1, -sinh(t)^2/cosh(t)^2); k = simple(kn/kd^(3/2)) Answer: k = (a^2*sinh(t)^2)/(cosh(t)^3*((a^2*sinh(t)^2)/cosh(t)^2)^(3/2))

1.34 (a) z = vpa('cos(pi*cos(pi*cos(log(pi+20))))+1', 40) Answer: 0.000000000000000000000000000000000039321609261272240194541

(b) z = vpa('exp(pi*sqrt(163))', 35) Answer: z = 262537412640768743.99999999999925007

1.35 syms x tay1 = taylor(x*(7*x+1)/(1-x)^3, 5, x, 0) tay2 = taylor((1+x-sqrt(1-6*x+x^2))/4,5,x,0) tay3 = taylor(2*x/(1-x)^3,5,x,0) tay4 = taylor(x*(x^2+4*x+1)/(1-x)^3,5,x,0) Answers: tay1 = 52*x^4 + 27*x^3 + 10*x^2 + x tay2 = 11*x^4 + 3*x^3 + x^2 + x tay3 = 20*x^4 + 12*x^3 + 6*x^2 + 2*x tay4 = 37*x^4 + 19*x^3 + 7*x^2 + x

1.36 syms x n e a t L1 = limit((x^e-1)/e, e, 0) L2 = limit((1-sin(2*x))^(1/x), x, 0) L3 = limit(((exp(t)-1)/t)^(-a), t, 0) L4 = limit(log(x^n)/(1-x^2), x, 1) Answers: L1 = log(x)

L2 = 1/exp(2)

L3 = 1 L4 = -n/2

© 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a7retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

1.37 syms s t rD = roots(sym2poly(s^4+0.282*s^3+4.573*s^2+0.4792*s+2.889)); rN = roots(sym2poly(0.1*s^3+0.0282*s^2-0.0427*s+0.0076)); Xt = vpa(ilaplace(((s-rN(1))*(s-rN(2))*(s-rN(3)))/((s-rD(1))*(s-rD(2))*(s-rD(3))*(s-rD(4))), s, t),5) Answer: Xt = (1.3974*(cos(1.9456*t)+0.050033*sin(1.9456*t)))/exp(0.097401*t)(0.39736*(cos(0.87145*t)+0.049799*sin(0.87145*t)))/exp(0.043599*t)

© 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a8retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

Solution to Exercises Chapter 2 2.1 n = 0:7; a = 2*n-1; b = 2*n+1; aplusb = a+b aminusb = a-b atimesb = a'*b detatimesb = det(atimesb) atimesbtranspose = a*b' Answers: aplusb = 0 4 8 12 16 20 24 28 aminusb = -2 -2 -2 -2 -2 -2 -2 -2 atimesb = -1 -3 -5 -7 -9 -11 -13 -15 1 3 5 7 9 11 13 15 3 9 15 21 27 33 39 45 5 15 25 35 45 55 65 75 7 21 35 49 63 77 91 105 9 27 45 63 81 99 117 135 11 33 55 77 99 121 143 165 13 39 65 91 117 143 169 195 detatimesb = 0 atimesbtranspose = 552

2.2 x = sort( [17 -3 -47 5 29 -37 51 -7 19], 'descend'); xs = [x(find(sort(x, 'descend')0))] Answer: w = -3

-7 -37 -47

51

29

19

17

5

2.3 y = [0, -0.2, 0.4, -0.6, 0.8, -1.0, -1.2, -1.4, 1.6]; sy = sin(y); % (a) mx = max(sy(find(sy1); d1 = On(1) tw = find(dn>2); d2 = tw(1) th = find(dn>3); d3 = th(1) Answer: d1 = 4 d2 = 31 d3 = 227

2.20 N = 10; k = 0:N; e = sum(1./factorial(k))-exp(1) Jo = sum((-1).^k./factorial(k).^2)-besselj(0,2) co = sum((-1).^k./factorial(2*k))-cos(1) Answers: e = -2.7313e-008 Jo = 6.3838e-016 co = -1.1102e-016

© 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a6retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

2.21 K = 41; k = 1:2:K; y = pi/3; x = 0.75; SKdiff = sum(exp(-k*x).*sin(k*y)./k)-0.5*atan(sin(y)/sinh(x)) Answer: SKdiff = -2.2204e-016

2.22 K = 14; k = 0:K; a = 3; x = -2; axdiff = sum((x*log(a)).^k./factorial(k))-a^x Answer: axdiff = 9.0230e-008

2.23 N = 500; k = 1:N; z = 10; c = pi*coth(pi*z)/(2*z)-1/(2*z^2) S1 = cumsum(1./(k.^2+z^2)); K = find(abs(S1-c)1e-8 xnew = xold-(cos(xold)+sin(xold))/(-sin(xold)+cos(xold)); fx = abs(cos(xnew)+sin(xnew)); cnt = cnt+1; xold = xnew; end disp(['At x = ' num2str(xold,8) ', f(x) = ' num2str(fx) ' after ' int2str(cnt) ' iterations']) Answer: At x = 2.3561945, f(x) = 1.1102e-016 after 4 iterations

4.5 ir = 0.045; P = 250000; A = 25000; disp('Principal Annuity Interest No. Remain') disp(' ($) ($/yr) (%) years ($)') R = P-A; n = 1; while R>A n = n+1; R = R*(1+ir)-A; end disp([num2str(P,6) blanks(6) num2str(A,5) blanks(8) num2str(ir*100,3) blanks(9) num2str(n,3) blanks(6) num2str(R,5)]) Answer: Principal Annuity Interest

No.

Remain

© 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a2retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

($) 250000

($/yr) 25000

(%) 4.5

years 12

($) 19112

4.6 alp = pi/4; tol = 1e-5; ao = 1; bo = cos(alp); co = sin(alp); while abs(co)>tol an = (ao+bo)/2; bn = sqrt(ao*bo); co = (ao-bo)/2; ao = an; bo = bn; end format long e Kalp = pi/ao/2 K = ellipke(sin(alp)^2) format short Answers: Kalp = 1.854074677301372e+000 K = 1.854074677301372e+000

4.7 h = [1 3 6 -7 -45 12 17 9]; a = 3; b=13; n = zeros(length(h),1); for k = 1:length(h) if h(k)=b n(k) = 0; else n(k) = 1; end end disp(['n = ' int2str(n')]) Answer: n = 0 0 1 0 0 1 0 1

4.8 N = input('Enter a positive integer < 10: '); h = zeros(N,N); for m = 1:N for n = 1:N if (n+m-1) 19)||(rem(num,1)~=0) num = input(' Enter any integer from 1 to 19: '); end

4.10 x = zeros(1000,1); x(1) = 3; cnt = 1; while x(cnt)~=1 cnt = cnt+1; if rem(x(cnt-1), 2) > 0 x(cnt) = 3*x(cnt-1)+1; else x(cnt) = x(cnt-1)/2; end end disp(['x = ' num2str(x(1:cnt)')]) Answer: x = 3 10 5 16 8 4 2 1

4.11 x(1) = 0; n = 201; x = zeros(n,1); for k = 2:n x(k) = x(k-1)^2+0.25; end plot(0:5:n-1, x(1:5:n), 'ks') m = 1; xw(1) = 0; xw = zeros(201, 1); while m1e-6 xnew = 0.5*(xold+a/xold); test = abs(xnew-xold); xold = xnew; cnt = cnt+1; end disp(['For xo = ' int2str(x(k)) ' the number of iterations = ' num2str(cnt)]) end Answers: For xo = 3 the number of iterations = 4 For xo = 100 the number of iterations = 10

4.13 for k = 1:3 switch k case 1 p=[1 2 3 4]; s = [10 20 30 40]; case 2 p = [11 12 13 14]; s = [101 102]; case 3 p = [43 54 55]; s = [77 66 88 44 33]; end n = length(p); m = length(s); d = n-m; if d==0 h = p+s elseif d1e-4) cnt = cnt+1; xnew = xold*a(k)*(1-xold); test = abs((xnew-xold)/xnew); xold = xnew; if cnt>N break end end disp(['For a = ' num2str(a(k)) ', x(n=' int2str(cnt) ') = ' num2str(xold)]) end Answers: For a = 1.45, x(n=19) = 0.31031 For a = 2.75, x(n=31) = 0.63634 For a = 3.2, x(n=201) = 0.51304 For a = 4, x(n=201) = 0.08737

4.15 Nu = 10; for k = 1:2 test = -1; while test < 0 NM(k) = input(['Enter a positive integer2): '); rb = 1.5; phi = linspace(0, 2*pi); rs = rb*sin(pi/n)/(1-sin(pi/n)); th = (1:n)*2*pi/n; plot(rb*cos(phi), rb*sin(phi), 'k', 0, 0, 'k+', (rb+rs)*cos(th), (rb+rs)*sin(th), 'k+') axis equal [x, phd] = meshgrid((rb+rs)*cos(th), phi); y = meshgrid((rb+rs)*sin(th), phi); hold on axis off plot(x+rs.*cos(phd), y+rs.*sin(phd), 'k-') Answer:

6.5 N = 15000; x = zeros(N,1); y=x; x(1) = 1; y(1) = 3.65; for n = 2:N x(n) = 1-y(n-1)+abs(x(n-1)); © 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a7retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

y(n) = x(n-1); end plot(x, y, 'k.') axis equal off Answer:

6.6 eta = linspace(0, 1, 10); phi = linspace(0, 2*pi, 200); c = cos(phi); s = sin(phi); a = 1.5; N = 4; for k = 1:length(eta) D = cosh(eta(k))-c; x = sinh(eta(k))./D; y = s./D; plot(a*x, a*y, 'k-', -a*x, -a*y, 'k-') hold on end phi = linspace(0, 2*pi, 30); for k = 1:length(phi) D = cosh(eta)-cos(phi(k)); x = sinh(eta)./D; y = sin(phi(k))./D; plot(a*x, a*y, 'k-', -a*x, -a*y, 'k-') hold on end axis equal off axis([-N N -N N]) Answer:

© 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a8retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

6.7 A = 2; N = 20; t = linspace(-pi, pi, 400); T = t(end)-t(1); signal = A*sin(t); dt = t(2)-t(1); x = linspace(-A, A, N); dx = x(2)-x(1); L = zeros(N-1,1); for k = 1:N-1 Index = find((signal>=x(k))&(signal

R(
)

© 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a22 retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

Solution to Exercises Chapter 7 7.1 % Spherical helix t = linspace(0, 10*pi, 300); x = t/(2*5); plot3(sin(x).*cos(t), sin(x).*sin(t), cos(x), 'k-') axis equal Answer: 1

0.5

0

-0.5

-1

0.5 0.5

0 0 -0.5

-0.5

% Toroidal spiral a = 0.2; b = 0.8; c = 20; t = linspace(0, 2*pi, 400); plot3((b+a*sin(c*t)).*cos(t), (b+a*sin(c*t)).*sin(t), a*cos(c*t), 'k-') axis equal Answer: 0.2 0.1 0 -0.1

0.8 0.6 0.4

0.8 0.2

0.6 0.4

0

0.2

-0.2

0

-0.4

-0.2 -0.4

-0.6 -0.8

-0.6 -0.8

% Sine wave on sphere a = 10; b = 1; c = 0.3; t = linspace(0, 2*pi, 200); z = sqrt(b^2-c^2*cos(a*t).^2); plot3(z.*cos(t), z.*sin(t), c*cos(a*t), 'k-') axis equal Answer:

© 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a1retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

0.2 0.1 0 -0.1 -0.2

0.8 0.6 0.4

0.8 0.2

0.6 0.4

0

0.2

-0.2

0 -0.4

-0.2 -0.4

-0.6

-0.6

-0.8

-0.8

% Concho-spiral a = 1; b = 1.05; c = 3; u = linspace(0, 12*pi, 200); bu = b.^u; plot3(a*bu.*cos(u), a*bu.*sin(u), c*bu, 'k-') axis equal Answer: 18 16 14 12 10 8 6 4

4 2

5 0 -2

0 -4 -5

% Intersection of two cylinders a = 1; b = 1.3; phi = linspace(0, 2*pi, 200); plot3(a*cos(phi), a*sin(phi), sqrt(b^2-a^2*sin(phi).^2), 'k-') axis equal view([-15, 20]) Answer: 1.2 1.1 1 0.9

0.5 0 -0.5

-0.8

-0.6

-0.4

-0.2

0

0.2

0.4

0.6

0.8

1

% Baseball seam a = 0.4; t = linspace(0, 4*pi, 200); g = pi/2-(pi/2-a)*cos(t); h = t/2+a*sin(2*t); plot3(sin(g).*cos(h), sin(g).*sin(h), cos(g), 'k-') © 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a2retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

axis equal Answer: 0.8 0.6 0.4 0.2 0 -0.2 -0.4 -0.6 -0.8

0.5 0.5 0 0 -0.5

-0.5

% Spherical spiral a = 0.08; phi = linspace(-12*pi, 12*pi, 400); p = atan(a*phi); plot3(cos(phi).*cos(p), sin(phi).*cos(p), -sin(p), 'k-') axis equal Answer: 0.8 0.6 0.4 0.2 0 -0.2 -0.4 -0.6 -0.8

0.5 0.5

0 0 -0.5

-0.5

7.2 % Seashell vv = linspace(0, 2*pi, 25); uu = linspace(0, 6*pi, 45); [u, v] = meshgrid(uu, vv); E = exp(u/6/pi); C2 = cos(0.5*v).^2; x = 2*(1-E).*cos(u).*C2; y = 2*(-1+E).*sin(u).*C2; z = 1-exp(u/3/pi)-sin(v)+E.*sin(v); surf(x, y, z) axis vis3d equal Answer:

© 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a3retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

0 -1 -2 -3 -4 -5 -6 -7 -8 2 2

1 1

0

0

-1

-1 -2

-2 -3

% Figure eight torus c = 1; vv = linspace(-pi, pi, 25); uu = linspace(-pi, pi, 45); [u, v] = meshgrid(uu, vv); C = c+sin(v).*cos(u)-sin(2*v).*sin(u)/2; x = cos(u).*C; y = sin(u).*C; z = sin(u).*sin(v)+cos(u).*sin(2*v)/2; surf(x, y, z) axis vis3d equal Answer:

1 0.5

0 -0.5

-1 1.5

2 1 1

0.5 0

0 -0.5 -1

-1 -1.5

-2

% Helical spring r1 = 0.25; r2 = 0.25; T = 2.0; n = 6; [u, v] = meshgrid(linspace(0, 2*n*pi, 140), linspace(0, 2*pi, 25)); x = [1-r1.*cos(v)].*cos(u); y = [1-r1*cos(v)].*sin(u); z = r2*[sin(v)+T*u/pi]; surf(x, y, z) view([-64 0]) surf(x, y, z) axis vis3d equal off Answer:

© 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a4retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

% Cornucopia a = 0.3; b = 0.5; [u, v] = meshgrid(linspace(0, 2*pi, 22), linspace(-3, 3, 20)); x = exp(b*v).*cos(v)+exp(a*v).*cos(v).*cos(u); y = exp(b*v).*sin(v)+exp(a*v).*sin(v).*cos(u); z = exp(a*v).*sin(u); surf(x, y, z) Answer:

% Astroidal ellipsoid [u, v] = meshgrid(linspace(-pi/2, pi/2, 40), linspace(-pi, pi, 40)); a = 1; b = 1; c = 1; x = (a*cos(u).*cos(v)).^3; y = (b*sin(u).*cos(v)).^3; z = (c*sin(v)).^3; surf(x, y, z) axis vis3d equal off Answer:

© 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a5retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

% Mobius strip [s, t] = meshgrid(linspace(0, 2*pi, 30), linspace(-0.4, 0.4, 15)); x = cos(s)+t.*cos(s/2).*cos(s); y = sin(s)+t.*cos(s/2).*sin(s); z = t.*sin(s/2); surf(x, y, z) axis vis3d equal off Answer:

% Bow curve T = 0.7; [u, v] = meshgrid(linspace(0, 2*pi, 20), linspace(0, 2*pi, 40)); x = (2+T*sin(u)).*sin(2*v); y = (2+T*sin(u)).*cos(2*v); z = T*cos(u)+3*cos(v); surf(x, y, z) axis off equal vis3d shading interp

Answer:

© 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a6retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

% Hyperbolic helicoid tau = 7; [u, v] = meshgrid(linspace(-pi, pi, 100), linspace(0, 0.5, 8)); D = 1+cosh(u).*cosh(v); x = sinh(v).*cos(tau*u)./D; y = sinh(v).*sin(tau*u)./D; z = sinh(u).*cosh(v)./D; surf(x, y, z) axis vis3d equal off Answer:

% Apple surface [u, v]= meshgrid(linspace(0, 2*pi, 30), linspace(-pi, pi, 50)); x = cos(u).*(4 + 3.8*cos(v)); y = sin(u).*(4 + 3.8*cos(v)); z = (cos(v) + sin(v) - 1).*(1 + sin(v)).* log(1 - pi*v / 10) + 7.5*sin(v); h = surf(x, y, z); set(h, 'FaceAlpha', 0.4) axis off equal vis3d shading interp Answer:

© 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a7retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

7.3 function Exercise7_3 the = linspace(0, 2*pi, 30); rad = linspace(0, 1, 15)'; xc = rad*cos(the); yc = rad*sin(the); options = optimset('display','off'); nm = 2; nk = 3; n = 0; for m = 0:1:nm vs = m+2.8; for k = 1:nk v = fzero(@circularplate, [vs vs+3], options, m); mode = circularplatemode(m, v, rad, the); vs = 1.2*v; n = n+1; subplot(nm+1, nk, n) meshc(yc, xc, mode) colormap([0 0 1]) title([num2str(v) ' m=' num2str(m) ' n=' num2str(k)]) axis off end end function d = circularplate(x, n) d = besselj(n, x)*besseli(n+1, x)+besseli(n, x)*besselj(n+1, x); function [mod]= circularplatemode(m, omn, rad, the) c1mn = -besseli(m, omn)/besselj(m, omn); mod = (c1mn*besselj(m, omn*rad)+ besseli(m, omn*rad))*cos(the*m); mod = mod/max(max(abs(mod))); Answer:

© 2011 Pearson Education, Inc., Upper Saddle River, NJ. All rights reserved. This publication is protected by Copyright and written permission should be obtained from the publisher prior to any prohibited reproduction, storage in a8retrieval system, or transmission in any form or by any means, electronic, mechanical, photocopying, recording, or likewise. For information regarding permission(s), write to: Rights and Permissions Department, Pearson Education, Inc., Upper Saddle River, NJ 07458.

3.1962 m=0 n=1

6.3064 m=0 n=2

9.4395 m=0 n=3

4.6109 m=1 n=1

7.7993 m=1 n=2

10.9581 m=1 n=3

5.9057 m=2 n=1

9.1969 m=2 n=2

12.4022 m=2 n=3

7.4 function Exercise7_4 heatslab = inline('cos(x)-x.*sin(x)/b', 'x', 'b'); tau = [linspace(0, 0.5, 9) linspace(0.6, 2, 9)]; eta = linspace(0, 1, 11); bi = 0.7; Nroot = 20; %x = linspace(0, 20*pi, 200); r = FindZeros(heatslab, Nroot, linspace(0, 20*pi, 200), bi); s = meshgrid(sin(r)./(r+sin(r).*cos(r)), eta); th = (2*cos(r*eta).*s')'*exp(-r.^2*tau); surf(tau, eta, th) xlabel('\tau') ylabel('\eta') zlabel('\theta/\theta_i') title(['Boit number = ',num2str(bi)]) axis vis3d view(-32.5, 36) function Rt = FindZeros(FunName, Nroot, x, w) f = feval(FunName, x, w); indx = find(f(1:end-1).*f(2:end)