From: rusin@vesuvius.math.niu.edu (Dave Rusin) Newsgroups: sci.math Subject: Re: A Question Nobody Ever Answer Me??? Date: 8 Jul 1998 08:16:08 GMT Much to my surprise, this is a very interesting question! It began by someone asking >The question is, by knowing the coordinates of how many points on the >cylinder surface, I could tell the vector n and a, maybe as well as R. Two threads have now arisen. - . - . - . - . - . - . - . - . - . - . - . - . - . - . - . - . - . - On 7 Jul 1998 03:42:29 GMT, enoether@aol.com (ENoether) wrote: >The question is ill-posed. I can give a continuum of points on a cylinder that >do not uniquely define a cylinder (In fact, a continuum of points that lie on a >continuum of cylinders). On the other hand, consider the following: >1) 3 points (defining a circle, S) and one point x, with proper conditions on x >(i.e., loosely speaking it should be perpendicular to the circle) should be >enough. David Einstein wrote: > This does not seem to be enough. The three points define an >infinite number of ellipses, and x can fit in an infinite number of >ways. The second poster's answer gave me pause, too, but it's correct as given. For comparison, a circle can be reconstructed from any three of its points, and any three noncollinear points define a unique circle; yet one can specify a circle in the plane unambiguously by giving any antipodal pair, and any pair of points define a unique circle in this way. So one must indeed clarify the question: Do we mean, "what is the least number N for which almost all N-tuples of points on a cylinder lie on no other cylinder?" This is certainly the usual interpretation. But perhaps we mean, "what is the smallest number N for which a surjection (R^3)^N --> Cyl exists?" (where Cyl is the set of all (infinite) cylinders in R^3). Here the answer is at most N=3: simply give me an antipodal pair of points on one circular cross-section, and a point on the cylinder which also lies on the plane containing the first two points and the center line. (Phew! ASCII text is not the right medium for this.) From these three points, and the information that they are so arrayed, I can unambiguously reconstruct the original cylinder. (One could be particularly obstreperous and insist the answer is N=1 since there exist bijections (R^3)^1 --> (R^3)^N. ) - . - . - . - . - . - . - . - . - . - . - . - . - . - . - . - . - . - >>2) 5 points (x_1,..., x_5) with x_1, x_2 and x_3 colinear and every other >>combination of 3 of these points noncolinear. This should be sufficient (though >>not necessary) to uniquely define a cylinder. > This seems to be enough. Are 5 points in general sufficient? Some fast talking noise about dimensions of manifolds/varieties/parameter spaces/degrees of freedom are sufficient to show that N=5 is the right number. But now we have a very specific question: given 5 points in R^3, characterize the set of cylinders to which they belong. In particular, is there a unique cylinder for almost every N ? The answer seems to be ... No! I ran some experiments to answer this. Pick 5 points at random (I chose them to have integral coordinates in [0, 99]). We may characterize the cylinders to which they belong by finding the planes through the origin onto which they project to five points on a common circle. The situation may be made algebraic as follows. Introduce coordinates A,B,C to parameterize the planes through the origin as the sets of points satisfying A X + B Y + C Z = 0. To remove redundancy, we insist A^2+B^2+C^2=1, at which point the A,B,C corresponding to a given plane are unique up to a common sign. The projections P' of the five points P to this plane are easily computed as P' = P - (P . e) e . So now we ask, for what choices of A,B,C do these five P' lie on a common circle? Further introduce variables x0,y0,z0 to represent the centers of circles in that plane. These are then bound by the second equation Ax0+By0+Cz0=0. The condition that the five P' are co-circular is then the condition that the five values dist( P', [x0,y0,z0] ) be equal. Equating the first and second, second and third, and so on gives four more equations to be satisfied. So our enumeration of cylinders has become the question of counting points on an algebraic variety (defined over Q) given by 6 equations in 6 unknowns. Clearly there will be sets of five points for which these equations are degenerate, and so the variety has positive dimension; but in general we expect it to be 0-dimensional, that is, it's a finite set. OK, now guess: how many points are in this variety? One, in general? No -- we have already commented that if (A,B,C,x0,y0,z0) is a solution, then so is (-A,-B,-C,x0,y0,z0) (although they represent the same cylinder). Very well, is the cardinality two then? It appears after a few experiments that 1. There are _12_ points in this variety (defined over C ) 2. For some sets of five points, _none_ of these 12 is real 3. For other sets of five points, at least _four_ are real Indeed, after very limited experimentation (here it is late at night) I am inclined to think that four is the "usual" number of solutions, so that: 5 randomly selected points on a cylinder actually lie in precisely two cylinders! Here are the data from one such experiment, should you wish to try your hand at visualizing things. Given points [77, 3, 37], [49, 23, 57], [35, 95, 99], [61, 73, 94], and [48, 89, 28], find the planes onto which they project co-circularly. Give up? Here are the only two answers: plane: .009075354356*X-.8141492553*Y-.5805847294*Z circle center: [58.50139565, -.7220114753, 1.926929773] plane: .7627655842*X-.6293169729*Y-.1488247669*Z circle center: [60.86572683, 60.16482550, 57.54106658] I would appreciate hearing a geometric argument which explains why various sets of 5 points belong to 0, 1, 2, ..., 6 (or more?) cylinders. dave ============================================================================== From: rusin@vesuvius.math.niu.edu (Dave Rusin) Newsgroups: sci.math Subject: About those intersecting cylinders... (LONG) Date: 9 Jul 1998 02:37:49 GMT Yesterday I noticed that a question about determining cylinders with points led to this formulation: >But now we have a very specific question: given 5 points >in R^3, characterize the set of cylinders to which they belong. >In particular, is there a unique cylinder for almost every [such 5-tuple] ? Today I pursued this a little further, because I couldn't keep the darn cylinders from spinning around in my head, parrying and thrusting, fidgeting around their sets of five fixed points. Anyone similarly cursed by a fixation on geometric problems may appreciate a summary, or perhaps a challenge: Can you find five points in "general position" incident to six different cylinders? (Answer below). - . - . - . - . - . - . - . - . - . - . - . - . - . - . - . - . - . - . - . - Despite the geometric nature of the problem, most of my ideas have been algebraic. Sorry. The approach I suggested yesterday certainly seems to be adequate for any specific set of 5 points: consider all projections to planes and see how many cases leave the 5 images P' equidistant from some center point: >Introduce coordinates A,B,C to >parameterize the planes through the origin as the sets of points >satisfying A X + B Y + C Z = 0. To remove redundancy, we insist >A^2+B^2+C^2=1, at which point the A,B,C corresponding to a given >plane are unique up to a common sign. [...] >Further introduce variables x0,y0,z0 to represent the centers of >circles in that plane. These are then bound by the second equation >Ax0+By0+Cz0=0. The condition that the five P' are co-circular [...] >gives four more equations to be satisfied. With a "little" bit of algebra, these equations may be reduced to simple equations of the form x0 = f1(C^2), y0 = f2(C^2), z0 = f3(C^2), A = C f4(C^2), B = C f5(C^2), and f6( C^2 ) = 0 where all the f_i are polynomials in (C^2) (f6 is a sextic; the others are quintic) whose coefficients are rational functions of the coordinates of the initial points. (OK, that's really just a _conjecture_. Actually carrying out the computations on a couple of computer-algebra systems seems to swamp their resources, even though in principle this kind of reduction is nowadays straightforward. When the initial points' coordinates are _numerical_, this reduction is in fact remarkably fast with Magma; I was even able to allow one coordinate x1 to remain variable and get answers in a matter of minutes; the coefficients in f6 were always polynomials of degree 12 in x1.) So just solve f6=0, take a square root, and plug into the other f_i to find the cylinders. Thus the cylinders passing through the five points are in one-to-one correspondence with non-negative, real solutions of f6. It appears to be true that _all_ real roots of f6 are non-negative. I don't know why that's true. When the coefficients of f6 were allowed to be polynomial in one coordinate x1 of one point, the constant term of f6 was always the _square_ of an integral polynomial, so that even as x1 passed through a root of that polynomial, the sign of the polynomial -- and so, the sign of the product of the roots of f6 -- remained unchanged. Very well, then, as we change one coordinate at a time of the 5 points, how might the number of incident cylinders change? This requires a change in the number of real roots of f6. This, in turn, requires a sign change in the discriminant of f6. With specific values for all but one coordinate, I was able to compute that discriminant as a function of the remaining coordinate x1. It appears always to be a polynomial of degree 116 (!) in x1. If the examples I worked on are any indication, this discriminant always factors as F1^2*F2. Here F1 (of degree 47!) will not affect the sign of the discriminant, but it's of interest anyway: this polynomial shows up in the denominator of f1, f2, and f3, so that values of x1 which make it vanish appear to be those for which some of the "cylinders" degenerate to planes (i.e., the centers are at infinity). It is thus the vanishing of the other polynomial F2 (of degree 22) which ultimately affects the changes in the number of incident cylinders. In the examples I tried, F2 had only a few real roots, but that's enough: each meant a change in the number of incident cylinders. (Geometrically, the vanishing of the polynomial F2 appears to mean that two incident cylinders coalesce.) Assuming the examples I tried are representative, we draw the following conclusions: 1. Almost every 5-tuple of points is incident to either 0, 2, 4, or 6 cylinders 2. The set of 5-tuples incident to 6 cylinders is open (in R^15); so is the set of 5-tuples incident to 4 cylinders, and the sets of those incident to 2 or 0 cylinders. Of course we are now faced with the major question: which (if any) of these sets of 5-tuples are empty? In particular, is it really possible to find an _open set_ of 5-tuples (no degnerate configurations allowed!) each of which touches 6 distinct cylinders? After wasting zillions of machine cycles, I was led to an answer that's pretty obvious: yes, there are 5-tuples of points incident to 6 cylinders. Here's the easy construction: any sufficiently tall, regular, rectangular pyramid will do. Two cylinders pass (perpendicularly) atop the four points on the base, to skewer an opposite pair of top faces. Four more cylinders pass through a top face and out the bottom. Voila! To see that this configuration is really part of an open set, I distorted the pyramid into the corner of an irregular cube. Consider these five points (which are cube vertices + random noise): [0,0,0], [100,0,0], [1,99,0], [98,103,7], [x1,5,91]. For all x1 near 100, there are 6 cylinders incident to these points; two pass through a side of the cube and out the bottom, and four pass nearly straight across opposite sides. (There are two nearly-equal vertical cylinders, a result of having moved the top vertex of the pyramid so nearly above one of the base vertices.) As x1 varies past the four real roots of F2, the number of incident cylinders changes; here's how many there are on the intervals: x1= -oo -1360.44 12.7953 112.440 424.946 oo <--------------*-------------*------------*------------*-------------> # cyl: 2 4 6 4 2 Thus not only 6, but also 2 and 4 cylinders occur on open sets of 5-tuples. Other 5-tuples of points chosen at random within this cube indicate the 0-cylinder case is also a non-empty open set. (e.g. the points [0, 0, 0], [61, 0, 0], [59, 42, 0], [79, 86, 12], [x1, 57, 76], with 64.471984 < x1 < 82.499852. ) This kind of problem is sometimes known as "Enumerative Geometry", and although I know nothing of the particulars, there is quite a lot of machinery for tackling these problems. We even faced some before in this newsgroup, e.g. http://www.math.niu.edu/~rusin/known-math/96/skewer [URL updated 1999/01 -- djr] See also http://www.math.niu.edu/~rusin/known-math/index/14-XX.html dave [Cylinders expunged, vertices exorcised. Still, it's fun to compute the answers randomly and, using xmaple, watch the points line up to lie in a circle on the screen!] ============================================================================== From: rusin@vesuvius.math.niu.edu (Dave Rusin) Newsgroups: sci.math,sci.math.symbolic,sci.math.num-analysis Subject: Re: Solving nonlinear polynomial simultaneous equations?? Date: 30 Oct 1998 18:53:33 GMT Cai Lixin wrote: >Given the 3D coordinates of 5 random points on a cylinder surface >p1(x1,y1,z1), p2(x2,y2,z2), ... p5(x5,y5,z5). Asked to calculate the >parameters(center line and a point it passing through) of cylinders >passing through these 5 points. For those who came in late: the idea behind this method was developed in some posts earlier this year; see http://www.math.niu.edu/~rusin/known-math/98/5pt.cyl Here is some Maple code to carry this out; it works quite quickly. dave #The problem assumes 5 points: for i from 1 to 5 do P.i:=[x.i,y.i,z.i]:od: #By translations and rotations we can in principle always have x5:=0:y5:=0:z5:=0:z4:=0:y4:=0:z3:=0: #Without this, Maple is unable to solve "center:=..." below; I don't know why #In order to describe the cylinders, I'm looking for a plane through the origin #which has the property that the projections of these 5 points to this plane #lie on a circle in that plane. #The plane through the origin is described by giving a unit normal [A,B,C] #to it. The assumption that this is a unit vector gives one equation already: eq1:=A^2+B^2+C^2-1: #(Note too that the vectors [A,B,C] and [-A,-B,-C] give the same plane.) #The center of the hoped-for circle is a point in space ce:=[x0,y0,z0]: #which is to lie on that plane, giving another equation: eq2:=A*x0+B*y0+C*z0: #Now, what are these five projections? proj(v) = v - normal. #So define inner products ip:=proc(P) A*P[1]+B*P[2]+C*P[3]:end: #and then projection pro:=proc(P) local u: u:=ip(P): [P[1]-u*A,P[2]-u*B,P[3]-u*C]:end: #The five projections are co-circular if there is a "center" ce equidistant #from all 5. Here's the distance function (square) ... di:=proc(P) local i:[seq((P[i]-ce[i])^2,i=1..3)]:"[1]+"[2]+"[3]:end: #...applied to the projections: din:=proc(n) di(pro(P.n)) end: #We just want those five values equal. #So the plane specified by [A,B,C] is correct if there is an x0,y0,z0 #making all these equations hold: eq:=expand([eq1,eq2,din(1)-din(2),din(2)-din(3),din(3)-din(4),din(4)-din(5)]): #Now, the right thing to do here is to eliminate all but one variable all #at once. Magma can do this; apparently Maple cannot. I can help it along. #The last three equations are linear in x0,y0,z0 so we can easily eliminate #these variables (generically): center:=solve({eq[4],eq[5],eq[6]},{x0,y0,z0}): subs(",[eq[1],eq[2],eq[3]]): eqc:=numer("): #I have found that the equations are a wee bit redundant: eq3:=[eqc[1],eqc[2], numer(factor(eqc[3]/eqc[1])) ]: #Down to three unknowns A,B,C! I can go a step further and eliminate A: lessA:=simplify(eq3,{A^2=1-B^2-C^2}): getA:=solve("[2],A): #which we can use to simplify the other equations to justBC:=numer(subs(A=getA,[eq3[1],lessA[3]])): #We have to do some elimination readlib(eliminate): #but at this point the algebra is pretty horrible; I'd rather go numeric #the rest of the way. #The other coordinates I'll pick at random: randset:={seq(x.i=rand() mod 100,i=1..4),seq(y.i=rand() mod 100,i=1..3), seq(z.i=rand() mod 100,i=1..2)}: ex:=subs(randset,justBC): getB:=eliminate({op(ex)},{B}): justC:="[2][1]: #Now you have a polynomial in one variable C to solve [fsolve(justC)]: Cset:=[]: for c in "" do if c>0 then Cset:=[op(Cset),c]:fi:od: N:=nops(Cset): Bset:=[]:for i to N do Bset:=[op(Bset),subs(subs(C=Cset[i],getB[1]),B)]:od: Aset:=[]:for i to N do Aset:=[op(Aset),subs({B=Bset[i],C=Cset[i],op(randset)},getA)]:od: centerset:=[]: for i to N do centerset:=[op(centerset), subs({A=Aset[i],B=Bset[i],C=Cset[i],op(randset)},center)]:od: solutions:=[]: for i to N do solutions:=[op(solutions), [Aset[i],Bset[i],Cset[i],op(subs(centerset[i],ce))]]:od: #Here's the result of one run: subs(randset,[seq([x.i,y.i,z.i],i=1..5)]); # [[96, 57, 76], [79, 86, 12], [59, 42, 0], [61, 0, 0], [0, 0, 0]] for i to N do lprint(solutions[i]): od: #[.7029796131, .7075116912, .7243527737e-1, 39.68253968, -53.57142855, 169.0168651] #[.3178160211e-2, -.9961577220, .8751975693e-1, 22.94146826, .3083564478, 78.84072420] #[2.380952381, -1.916934473, .4959141161, -522.5127885, 134.8268235, 1636.524577] #I also use XMaple, which allows me to look at 3D objects from varying #perspectives. This routine works: # with(plots): which:=1: #which solution to view? 1..N aha:={A=solutions[which][1], B=solutions[which][2], C=solutions[which][3], x0=solutions[which][4], y0=solutions[which][5], z0=solutions[which][6]}: L1:=pointplot3d(subs(randset,[P1,P2,P3,P4,P5]),color=BLACK): L2:=pointplot3d(subs(aha,[[100*A,100*B,100*C],[10*A,10*B,10*C]]),color=BLUE): L3:=pointplot3d(subs(aha,[[x0,y0,z0]]),color=RED): display3d({L1,L2,L3}); #Click on the image to bring up menus; select BOX coords, NO perspective, #and CONSTRAINED view. Rotate box until the two blue dots align, and then the #red should be in the center of a circle of blacks. #Cool, huh? ============================================================================== From: spellucci@mathematik.tu-darmstadt.de (Peter Spellucci) Newsgroups: sci.math.num-analysis Subject: Re: HELP??? Date: 24 Nov 1998 11:04:01 GMT In article <22E71DAEC504D111B78100805FFE9DC714C832AB@PFS21>, Cai Lixin writes: |> Hello, All, |> |> I would appreciate if you could share your idea with me how to code the |> following 5 nonlinear equations. It is all the better if you are so |> interested that you would code it. |> |> |> [(1-a^2)*(xi-x0) - a*b*(yi-y0) - a*c*zi]^2 + [(1-b^2)*(yi-y0) - |> a*b*(xi-x0) - b*c*zi]^2 + [(1-c^2)*zi - a*c*(xi-x0) - b*c*(yi-y0)]^2 = |> R^2, |> |> a^2+b^2+c^2 = 1 |> |> Where, a, b, c, x0,y0 is the unknown parameters. |> xi, yi, zi, R is the known parameters. i runs from 1 to 4. One set of |> knowns I can give is (77, 3, 37, 32.3536), (49, 23, 57, 32.3536), (35, |> 95, 99, 32.3536), (61, 73, 94, 32.3536). |> subroutine fcn(n,x,f,ifail) implicit none integer n,ifail c***** these are dummys here double precision x(5),f(5) c**** local variables integer i,j double precision xi(4),yi(4),zi(4),R(4),a,b,c,x0,y0 data (xi(j),yi(j),zi(j),R(j),j=1,4)/ f 77, 3, 37, 32.3536 , 49, 23, 57, 32.3536 , f 35, 95, 99, 32.3536, 61, 73, 94, 32.3536 / c****** begin of computation ifail=0 a=x(1) b=x(2) c=x(3) x0=x(4) y0=x(5) do i=1,4 f(i)=(1.d0-a**2)*(xi(i)-x0) - a*b*(yi(i)-y0) - a*c*zi(i))**2 + f ((1.d0-b**2)*(yi(i)-y0) - f a*b*(xi(i)-x0) - b*c*zi)**2 + ((1.d0-c**2)*zi(i) - f a*c*(xi(i)-x0) - b*c*(yi(i)-y0))**2 - f R(i)**2 enddo f(5)=a**2+b**2+c**2-1.d0 return end to be used in cnjunction with nleq1.f from deuflhards group http://elib.zib.de (nleq1.f is to be found in the CODELIB hope this helps peter ============================================================================== Update: a paper has appeared on this and related topics: INSTITUT NATIONAL DE RECHERCHE EN INFORMATIQUE ET EN AUTOMATIQUE No. 4195: On circular cylinders by four or five points in space Olivier Devillers, Bernard Mourrain, Franco P. Preparata, Philippe Trebuchet Juin 2001 www-calfor.lip6.fr/~trebuche/RR-4195.ps.gz (Note: Yes, the English title has been copied accurately here.)