MathGroup Archive 2002

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Re: Intersection Ellipse & Circle


A few simple replacements can factor out most of the remaining complexity:

SetOptions[Roots, Cubics->False, Quartics->False];
Solve[{(x - c)^2/b^2 + (y - d)^2/a^2 == 1, x^2 + y^2 == 1}, {x, y}] /.
 a^4 - 2*a^4*b^2 + .. <snip> .. + b^4*#1^4 -> P

where P is the first argument of Root[...] in the output. The result is:
x -> (b^2*(d - y)^2 - a^2*(-1 + b^2 - c^2 + y^2))/(2*a^2*c)
y->Table[Root[P&,n],{n,1,4}]

Note: You probably want to filter out complex solutions.


> One thing to realize is you are getting solutions invloving symbolic
> roots of quartics. These can be quite large and moreover are not
> necessarily useful for plugging in numeric values after the fact (due to
> instability). You can reduce to a more manageable size by forcing the
> underlying Roots code to avoid the Cardano-Tartaglia formulas. The code
> below indicates that most of the time is actually spent in messing with
> the quartic formula so disabling it will save time as well as memory in
> this case.
>
> In[2]:= Timing[s1 = Solve[{(x-c)^2/b^2 + (y-d)^2/a^2 == 1, x^2 + y^2 ==
> 1}, {x,y}];]
> Out[2]= {2.48 Second, Null}
>
> In[3]:= LeafCount[s1]
> Out[3]= 286345
>
> In[4]:= SetOptions[Roots, Cubics->False, Quartics->False];
>
> In[5]:= Timing[s2 = Solve[{(x-c)^2/b^2 + (y-d)^2/a^2 == 1, x^2 + y^2 ==
> 1}, {x,y}];]
> Out[5]= {0.25 Second, Null}
>
> In[6]:= LeafCount[s2]
> Out[6]= 4809
>
>
> Daniel Lichtblau
> Wolfram Research
>
>



  • Prev by Date: Re: Unloading packages
  • Next by Date: Re: principle root? problem
  • Previous by thread: Re: Intersection Ellipse & Circle
  • Next by thread: Re: Intersection Ellipse & Circle