-------------------------------------------------------------------------------- # Notes on project 2: # # This problem involves ideas like those needed in project 2. # # Given a right triangle with legs a and b. Inscribe a circle c1, then # heading off towards one of the acute angles, inscribe an infinite # chain of circles cn each tangent to the preceeding circle and the sides # of the angle. Find the total area of the circles. # # Let rn be the radius of the nth circle. # # # First we get the radius of the 1st circle. # Let r be the radius. Then the center is at (r,r) # Let (x,y) be the point on the hypotenuse where the circles is # tangent. Then the slope of the segment from (r,r) to (x,y) is b/a, # and the slope of the segment from (b,0) to (x,y) is -a/b. We can use # this get x and y in terms of a,b, and r. # > restart;eq1 := (y-r)/(x-r)=b/a; y - r eq1 := ----- = b/a x - r > eq2 := (y-0)/(x-b) = -a/b; y eq2 := ----- = - a/b x - b > sol :=solve({eq1,eq2},{x,y}); 2 2 a (a r + b - b r) b (- a r + b r + a ) sol := {y = ------------------, x = --------------------} 2 2 2 2 a + b a + b # Now the length of the segment from (x,y) to (r,r) is r. # # So, we can get one equation in a,b, and r from this which we can solve # for r in terms of a and b. > > eq3 := r = simplify(subs(sol,sqrt((x-r)^2+(y-r)^2))); / 2\1/2 |(- b r + a b - a r) | eq3 := r = |--------------------| | 2 2 | \ a + b / # Solving for r, > sol2 := solve(eq3,r); 2 sol2 := RootOf(2 _Z + (- 2 a - 2 b) _Z + a b) # We can get all the roots of this with allvalues > sol2 := allvalues(sol2); 2 2 1/2 2 2 1/2 sol2 := 1/2 a + 1/2 b + 1/2 (a + b ) , 1/2 a + 1/2 b - 1/2 (a + b ) # By inspection, we can see that r must be the second root. Let's # make a function out of this. > rad := unapply(sol2[2],a,b); 2 2 1/2 rad := (a,b) -> 1/2 a + 1/2 b - 1/2 (a + b ) # Just to check if we're on the right track, let's draw the triangle # and the inscribed circle. We make these general so that a and b # can be varied. -------------------------------------------------------------------------------- > tri := (a,b) -> plot([[0,0],[b,0],[0,a],[0,0]],scaling=constrained); tri := (a,b) -> plot([[0, 0], [b, 0], [0, a], [0, 0]], scaling = constrained) -------------------------------------------------------------------------------- > circle := (c,r) -> plot([c[1]+r*cos(t),c[2]+r*sin(t),t=0..2*Pi]); circle := (c,r) -> plot([c[1] + r cos(t), c[2] + r sin(t), t = 0 .. 2 Pi]) -------------------------------------------------------------------------------- > with(plots); [animate, animate3d, conformal, contourplot, cylinderplot, densityplot, display, display3d, fieldplot, fieldplot3d, gradplot, gradplot3d, implicitplot, implicitplot3d, loglogplot, logplot, matrixplot, odeplot, pointplot, polarplot, polygonplot, polygonplot3d, polyhedraplot, replot, setoptions, setoptions3d, spacecurve, sparsematrixplot, sphereplot, surfdata, textplot, textplot3d, tubeplot] # -------------------------------------------------------------------------------- > inscribe := proc(a,b) \ local r;\ r := rad(a,b);\ plots[display]([tri(a,b),circle([r,r],r)]); end; inscribe := proc(a,b) local r; r := rad(a,b); plots[display]( [tri(a,b),circle([r,r],r)]) end -------------------------------------------------------------------------------- > inscribe(4,7); -------------------------------------------------------------------------------- # Now get the ratio of rn to r(n-1). Here we can use some trig. # # Look at the right triangle whose hypotenuse is the segment # connecting the centers of the n-1st and nth circles. The # hypotenuse is r(n-1)+rn, the angle alpha on the lower right # vertex is 1/2 arctan(a/b), and the side opposite this angle # is r(n-1) - rn, which is also (r(n-1)+rn) sin(alpha). This is eq4 > eq4 := r[n-1]-r[n]= (r[n-1]+r[n])*sin(arctan(a/b)/2); eq4 := r[n - 1] - r[n] = (r[n - 1] + r[n]) sin(1/2 arctan(a/b)) # Solve this equation for rn > sol := solve(eq4,r[n]); r[n - 1] - sin(1/2 arctan(a/b)) r[n - 1] sol := - ---------------------------------------- - 1 - sin(1/2 arctan(a/b)) > collect(sol,r[n-1]) ; (1 - sin(1/2 arctan(a/b))) r[n - 1] - ----------------------------------- - 1 - sin(1/2 arctan(a/b)) # Now we can see what the common ratio is: > rat := unapply( (1-sin(arctan(a/b)/2))/(1+sin(arctan(a/b)/2)) ,a,b); 1 - sin(1/2 arctan(a/b)) rat := (a,b) -> ------------------------ 1 + sin(1/2 arctan(a/b)) # Now we can write the total Area of the circles. -------------------------------------------------------------------------------- > Area := Sum(Pi*(rad(a,b)*rat(a,b)^i)^2,i=0..infinity); Area := infinity ----- \ 2 2 1/2 2 //1 - sin(1/2 arctan(a/b))\i\2 ) Pi (1/2 a + 1/2 b - 1/2 (a + b ) ) ||------------------------| | / \\1 + sin(1/2 arctan(a/b))/ / ----- i = 0 # This is a geometric series, which maple can sum. (We could too.) > form := value(Area ); 2 2 1/2 2 2 Pi (1/2 a + 1/2 b - 1/2 (a + b ) ) (1 + sin(1/2 arctan(a/b))) form := 1/4 ------------------------------------------------------------------ sin(1/2 arctan(a/b)) # Make a function out of this. > area := unapply(form,a,b); area := 2 2 1/2 2 2 Pi (1/2 a + 1/2 b - 1/2 (a + b ) ) (1 + sin(1/2 arctan(a/b))) (a,b) -> 1/4 ------------------------------------------------------------------ sin(1/2 arctan(a/b)) > area(4.,11.); 5.386806700 Pi # Now lets draw a picture of the triangle and the first few, say m, circles. > a := 4. : b := 7.: m := 7: -------------------------------------------------------------------------------- > for i from 2 to m do r.i := rat(a,b)^(i-1)*r1 od: -------------------------------------------------------------------------------- > x1 := evalf(r1); x1 := 1.468871126 # Get the x coordinate of the center of the ith circle i from 2 to m, # use some trig. -------------------------------------------------------------------------------- > for i from 2 to m do x.i := evalf(x.(i-1) + (r.i+r.(i-1))*cos(arctan(a/b)/2)) od: -------------------------------------------------------------------------------- > for i from 1 to m do cir.i := circle([x.i,r.i],r.i) od: -------------------------------------------------------------------------------- > display([tri(a,b),cir1,cir2,cir3,cir4,cir5 ,cir6,cir7]); # We can make a word out this. > drawem := proc(a,b,m)\ local i,r,x,cir,ar;\ r[1] := evalf(rad(a,b));\ for i from 2 to m do r[i]:= evalf(rat(a,b)^(i-1)*r[1] )od:\ x[1] := r[1]:\ for i from 2 to m do x[i] := \ evalf(x[i-1]+ (r[i]+r[i-1])*cos(arctan(a/b)/2)) od:\ for i from 1 to m do cir[i]:= circle([x[i],r[i]],r[i]) od:\ ar := evalf(area(a,b));\ plots[display]([tri(a,b),seq(cir[i],i=1..m)],title=cat(`area of circles is`,convert(ar,string))); \ end; drawem := proc(a,b,m) local i,r,x,cir,ar; r[1] := evalf(rad(a,b)); for i from 2 to m do r[i] := evalf(rat(a,b)^(i-1)*r[1]) od; x[1] := r[1]; for i from 2 to m do x[i] := evalf(x[i-1]+ (r[i]+r[i-1])* cos(1/2*arctan(a/b))) od; for i to m do cir[i] := circle([x[i],r[i]],r[i]) od; ar := evalf(area(a,b)); plots[display]([ tri(a,b),seq(cir[i],i = 1 .. m) ],title = cat( `area of circles is`, convert(ar,string))) end -------------------------------------------------------------------------------- > drawem(6,7,7); >