print(`pinsky: a small Maple package for answering`): print(`a Calculus exercise posed by Mark Pinsky`): print(`Version of April 7, 1997`): print(`written by Doron Zeilberger(zeilberg@math.temple.edu).`): print(`The most current version of the package`): print(` is available from`): print(`http://www.math.temple.edu/~zeilberg`): print(`For general help, and a list of the procedures`): print(` type "ezra():" (no quotes), for help with`): print(`a specific procedure, type "ezra(procedure_name):" (no quotes).`): print(``): ezra:=proc() if args=NULL then print(`This is pinsky, a small Maple package that answers `): print(`the following question posed by Mark Pinsky.`): lprint(``): print(`Let (x0,y0,z0) be a point in the positive octant of`): print(` Euclidean 3D space. Every plane passing through it`): print(` meets the coordinate planes in 3 edges that form a triangle.`): print(`Find the minimum area of the triangle formed, and the plane`): print(`that achieves that minimum.`): lprint(``): print(`It turns out that the answer can be expressed in terms of the`): print(`root of an algebraic equation of degree 7 with coefficients`): print(`that are rational functions in the coordinates of the point`): print(``): print(`Contains the following procedures:`): print(`pinsky, pinskyG , pinskyGf, MinArea, MinArea4Sq `): fi: if nops([args])=1 and op(1,[args])=`MinArea4Sq` then print(`MinArea4Sq(x0,y0,z0): inputs the coordinates of a point`): print(`(either numerically or symbolically)`): print(`in the first octant, and outputs 4x the (area)^2 of the triangle`): print(` of least possible area, that can be formed by the edges`): print(`formed from the intersection of a plane`): print(`passing through the given point, with the coordinate axes`): fi: if nops([args])=1 and op(1,[args])=`MinArea` then print(`MinArea(x0,y0,z0): inputs the coordinates of a point`): print(`(either numerically or symbolically)`): print(`in the first octant, and outputs the area of the triangle`): print(` of least possible area, that can be formed by the edges`): print(`formed from the intersection of a plane`): print(`passing through the given point, with the coordinate axes`): fi: if nops([args])=1 and op(1,[args])=`pinsky` then print(`pinsky(x0,y0,z0) inputs 3 specific positive numbers`): print(`and gives the answer in floating point`): print(`For example, try pinsky(1,4,6):, and you would get`): print(`the eq. of the plane that minimizes the area of the triangle`): print(`formed by intersection-edges with the coordinate axes`): print(`and the corresponding area`): fi: if nops([args])=1 and op(1,[args])=`pinskyG` then print(`pinskyG(x0,y0,z0) inputs 3 specific positive numbers`): print(` or symbolic numbers (or any combination)`): print(`and gives the answer in symbolic form.`): lprint(``): print(`The inputs (i.e. coordinates) should be all distinct.`): lprint(``): print(`Without loss of generality, the eq. of any plane passing`): print(`through the point (x0,y0,z0) may be written`): print(` (x-x0)+B(y-y0)+C(z-z0)=0.`): print(` The procedure gives the minimal polynomial that B `): print(` is a root of, with coefficients that are rational functions`): print(`of the coordinates (x0,y0,z0).`): lprint(``): print(`To get the complete answer to Mark Pinsky's question`): print(`you should type pinskyG(x0,y0,z0):`): fi: if nops([args])=1 and op(1,[args])=`pinskyGf` then print(`Does what pinsky does, but via the method of pinskyG`): print(`and then converting the output to floats`): print(`It verifies that pinskyG is valid`): fi: end: pinsky:=proc(x0,y0,z0) local A,B,C,eq,Pl,f,lu,x,y,z,x1,y1,z1,bas,shoreshB,shoreshC,Bpol,eq1: with(grobner): if x0<=0 or y0<=0 or z0<=0 then ERROR(`All arguments must be positive numbers`): fi: A:=1; Pl:=A*(x-x0)+B*(y-y0)+C*(z-z0): x1:=solve(subs({y=0,z=0},Pl),x): y1:=solve(subs({x=0,z=0},Pl),y): z1:=solve(subs({x=0,y=0},Pl),z): f:=(x1*z1)^2+(x1*y1)^2+(y1*z1)^2: lu:=(x0+B*y0+C*z0)^3: eq:={ factor(numer(normal(diff(f,B)))), factor(numer(normal(diff(f,C)))) }: eq:={normal(op(1,eq)/lu),normal(op(2,eq)/lu)}: print(`The best plane is (x-x0)+B(y-y0)+C(z-z0)=0, where B is a root of`): bas:=gbasis(eq,[C,B],plex): Bpol:=bas[nops(bas)]: print(Bpol): print(`its real roots are`): shoreshB:=fsolve(Bpol,B): print(shoreshB): print(`C is a solution of`): print(bas): print(`whose corresponding value, to the above value of B is`): eq1:=subs(B=shoreshB, convert([op(1..nops(bas)-1,bas)],set)) minus {0}: shoreshC:=subs(fsolve(eq1,C),C): print( shoreshC): lu:=subs({B=shoreshB, C=shoreshC},f): lu:=sqrt(lu)/2: print(`The minimal area of a triangle made by intersection of the`): print(`coordinate planes, passing through the point`): print([x0,y0,z0]): print(` happens to be `,lu): print(`it is achieved by the plane`): print((x-x0)+(y-y0)*shoreshB+(z-z0)*shoreshC=0): lu,shoreshB,shoreshC: end: #Given an algebraic number that is a solution of POL(x)=0 #and polynomials A and B in x, expresses A/B as a polynomial #in x over the ground field with(grobner): khalek:=proc(A,B,POL,x) local i,A1,B1,eq,var,c,deg,C1,gu: deg:=degree(POL,x)-1: A1:=normalf(A,[POL],[x]): B1:=normalf(B,[POL],[x]): eq:={}: C1:=0: var:={}: for i from 0 to deg do C1:=C1+c[i]*x^i: var:=var union {c[i]}: od: gu:=B1*C1-A1: gu:=expand(gu): gu:=normalf(gu,[POL],[x]): gu:=expand(gu): eq:={}: for i from 0 to deg do eq:=eq union {coeff(gu,x,i)}: od: var:=solve(eq,var): C1:=subs(var,C1): C1:=expand(C1): C1: end: pinskyG:=proc(x0,y0,z0) local i,A,B,C,eq,Pl,f,lu,x,y,z,x1,y1,z1,bas,Bpol,topf,botf,C1,Bpol1,Y,f1: with(grobner): if x0=y0 or y0=z0 then ERROR(`the the 3 coordinates should be all different`): fi: A:=1; Pl:=A*(x-x0)+B*(y-y0)+C*(z-z0): x1:=solve(subs({y=0,z=0},Pl),x): y1:=solve(subs({x=0,z=0},Pl),y): z1:=solve(subs({x=0,y=0},Pl),z): f:=(x1*z1)^2+(x1*y1)^2+(y1*z1)^2: lu:=(x0+B*y0+C*z0)^3: eq:={ factor(numer(normal(diff(f,B)))), factor(numer(normal(diff(f,C)))) }: eq:={normal(op(1,eq)/lu),normal(op(2,eq)/lu)}: bas:=gbasis(eq,[C,B],plex): Bpol:=bas[nops(bas)]: C1:=solve(bas[1],C): Bpol1:=factor(Bpol): for i from 1 to nops(Bpol1) do if degree(op(i,Bpol1),B)=7 then Bpol1:=op(i,Bpol1): break: fi: od: C1:=normalf(C1,[Bpol1],[B]): f:=normal(f); f:=subs(C=C1,f): topf:=numer(f): botf:= denom(f): topf:=normalf(topf,[Bpol1],[B]): botf:=normalf(botf,[Bpol1],[B]): print(`The minimal area is (1/2) times the square root of`): f:=khalek(topf,botf,Bpol1,B): print(f): print(`where B is the real root of the polynomial`): print(Bpol1): print(`Also, the square of the area times 2 is the real root`): print(` of the polynomial`): f1:=root1(f,Bpol1,B,Y): print(f1): f,Bpol1,f1: end: MinArea:=proc(x0,y0,z0) local i,A,B,C,eq,Pl,f,lu,x,y,z,x1,y1,z1,bas,Bpol,topf,botf,C1,Bpol1,Y,f1: with(grobner): if x0=y0 or y0=z0 then ERROR(`the the 3 coordinates should be all different`): fi: A:=1; Pl:=A*(x-x0)+B*(y-y0)+C*(z-z0): x1:=solve(subs({y=0,z=0},Pl),x): y1:=solve(subs({x=0,z=0},Pl),y): z1:=solve(subs({x=0,y=0},Pl),z): f:=(x1*z1)^2+(x1*y1)^2+(y1*z1)^2: lu:=(x0+B*y0+C*z0)^3: eq:={ factor(numer(normal(diff(f,B)))), factor(numer(normal(diff(f,C)))) }: eq:={normal(op(1,eq)/lu),normal(op(2,eq)/lu)}: bas:=gbasis(eq,[C,B],plex): Bpol:=bas[nops(bas)]: C1:=solve(bas[1],C): Bpol1:=factor(Bpol): for i from 1 to nops(Bpol1) do if degree(op(i,Bpol1),B)=7 then Bpol1:=op(i,Bpol1): break: fi: od: C1:=normalf(C1,[Bpol1],[B]): f:=normal(f); f:=subs(C=C1,f): topf:=numer(f): botf:= denom(f): topf:=normalf(topf,[Bpol1],[B]): botf:=normalf(botf,[Bpol1],[B]): f:=khalek(topf,botf,Bpol1,B): f1:=root1(f,Bpol1,B,Y): (1/2)*sqrt(solve(f1)): end: MinArea4Sq:=proc(x0,y0,z0) local i,A,B,C,eq,Pl,f,lu,x,y,z,x1,y1,z1,bas,Bpol,topf,botf,C1,Bpol1,Y,f1: with(grobner): if x0=y0 or y0=z0 then ERROR(`the the 3 coordinates should be all different`): fi: A:=1; Pl:=A*(x-x0)+B*(y-y0)+C*(z-z0): x1:=solve(subs({y=0,z=0},Pl),x): y1:=solve(subs({x=0,z=0},Pl),y): z1:=solve(subs({x=0,y=0},Pl),z): f:=(x1*z1)^2+(x1*y1)^2+(y1*z1)^2: lu:=(x0+B*y0+C*z0)^3: eq:={ factor(numer(normal(diff(f,B)))), factor(numer(normal(diff(f,C)))) }: eq:={normal(op(1,eq)/lu),normal(op(2,eq)/lu)}: bas:=gbasis(eq,[C,B],plex): Bpol:=bas[nops(bas)]: C1:=solve(bas[1],C): Bpol1:=factor(Bpol): for i from 1 to nops(Bpol1) do if degree(op(i,Bpol1),B)=7 then Bpol1:=op(i,Bpol1): break: fi: od: C1:=normalf(C1,[Bpol1],[B]): f:=normal(f); f:=subs(C=C1,f): topf:=numer(f): botf:= denom(f): topf:=normalf(topf,[Bpol1],[B]): botf:=normalf(botf,[Bpol1],[B]): f:=khalek(topf,botf,Bpol1,B): f1:=root1(f,Bpol1,B,Y): end: pinskyGf:=proc(x0,y0,z0) local A,B,C,eq,Pl,f,lu,x,y,z,x1,y1,z1,bas,Bpol,topf,botf,C1,Bpol1, Bshoresh: with(grobner): A:=1; Pl:=A*(x-x0)+B*(y-y0)+C*(z-z0): x1:=solve(subs({y=0,z=0},Pl),x): y1:=solve(subs({x=0,z=0},Pl),y): z1:=solve(subs({x=0,y=0},Pl),z): f:=(x1*z1)^2+(x1*y1)^2+(y1*z1)^2: lu:=(x0+B*y0+C*z0)^3: eq:={ factor(numer(normal(diff(f,B)))), factor(numer(normal(diff(f,C)))) }: eq:={normal(op(1,eq)/lu),normal(op(2,eq)/lu)}: bas:=gbasis(eq,[C,B],plex): Bpol:=bas[nops(bas)]: C1:=solve(bas[1],C): Bpol1:=factor(Bpol): if degree(op(1,Bpol1),B)=7 then Bpol1:=op(1,Bpol1): else Bpol1:=op(2,Bpol1): fi: C1:=normalf(C1,[Bpol1],[B]): f:=normal(f); f:=subs(C=C1,f): topf:=numer(f): botf:= denom(f): topf:=normalf(topf,[Bpol1],[B]): botf:=normalf(botf,[Bpol1],[B]): print(`The minimal area is is (1/2) times the square root of`): f:=khalek(topf,botf,Bpol1,B): print(f): print(`where B is the real root of the polynomial`): print(Bpol1): print(`This real root, in floating point, happens to be`): Bshoresh:=fsolve(Bpol1,B): f:=subs(B=Bshoresh,f): f:=sqrt(f)/2: print(`So in floating point, the minimal area is`): print(f): C1:=subs(B=Bshoresh,C1): print(`And the plane that achieves it is `): Pl:=(x-x0)+(y-y0)*Bshoresh+(z-z0)*C1: print(Pl): f,Pl: end: #root1(P,POL,x,y): Given a polynomial P in x and a pol POL #such that POL(x)=0 finds the polynomial Q(y) such that Q(P(x))=0 root1:=proc(P,POL,x,y) local bas,P1,var,eq: with(grobner): P1:=normalf(P,[POL],[x]): eq:={y-P1,POL}: var:=[x,y]: bas:=gbasis(eq,var,plex): RETURN(bas[nops(bas)]): end: