Services & Resources / Wolfram Forums
-----
 /
MathGroup Archive
1992
*January
*February
*March
*April
*May
*June
*July
*August
*September
*October
*November
*December
*Archive Index
*Ask about this page
*Print this page
*Give us feedback
*Sign up for the Wolfram Insider

MathGroup Archive 1992

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

Search the Archive

ComplexRoots

  • To: mathgroup at yoda.physics.unc.edu
  • Subject: ComplexRoots
  • From: akyildiz at maths.ox.ac.uk ( Yilmaz Akyildiz tel 2-73558)
  • Date: Thu, 19 Nov 92 21:04:08 GMT

(* 
This is a Mathematics V.2.0  programme allowing  the user to
numerically compute all the common roots of the two real functions 
f1 and f2 of two variables x1 and x2 on a given rectangle in the 
(x1, x2)-plane. It can easily be generalized to higher dimensions.

	comroots[{f1, f2}, {x1, a, b}, {x2, c, d}, tolerance] 

will find all the common roots of the functions f1(x1, x2) and 
f2(x1, x2) with accuracy correct to the tolerance. It uses a version 
of Interval Bisection Method. 
For example: 
	comroots[{x1 Sin[x2] -1, x2 Cos[x1] + 1}, 
				{x1, 0, 5}, {x2, 0, 9}] 

finds all the roots of  {x1 Sin[x2] -1 = 0, x2 Cos[x1] + 1 = 0} 
in the intervals 0 <= x1 <= 5 and 0 <= x2 <= 9 with accuracy correct to 
(the default) 3 decimal digits. 

In case f1 and f2 are real and imaginary parts of a single
function f, the result will be all the (real and complex) roots of the 
function f(x) in the specified rectangle. 
For example,
	comroots[{Exp[2x1] Cos[2x2] + 3, 
		  Exp[2x1] Sin[2x2]}, {x1, 0, 1}, {x2, 0, 30}]

finds all the (complex) roots of Tanh(x) = 2 on the rectangle 
[0, 1] x [0, 30]. 

In case f(x) is a POLYNOMIAL with real or complex coefficients, after 
running the accompanying pacakge reim. m (from Maeder), 

	cxroots[f(x), {x1, a, b}, {x2, c, d}, tol]  
	
will give all the real and complex roots of the polynomial f in the 
rectangle [a, b] x [c, d].

For example, try 
	cxroots[x^3 +1, {x1, -2, 1}, {x2, -1, 1}]
(Sorry for the repetitions in the output...).

The output will be displayed in terms of the Mathematica built-in 
function RealInterval.  {RealInterval[{x1,x2}], RealInterval[{y1,y2}]
will stand for the root which is contained in the closed
rectangle [x1, x2] x [y1, y2].  This output format is preferable
since Mathematica can do arithmetic with RealInterval objects.

This is a rather slow programme. Faster versions, using 
Secant and Newton's methods, are also available. 

The built-in function RealInterval is not yet bug-free.  
I would be interested to know if you come across any anomolies
related to this package.

Yilmaz Akyildiz 
akyildiz at maths.oxford.ac.uk
*)

BeginPackage["ReIm`"];

Begin["`Private`"];

protected = Unprotect[Re, Im]

ri[x_,y_] := RealInterval[{x, y}]

Re[x_] := x /; Im[x] == 0

Re[x_ + y_] := Re[x] + Re[y]
Im[x_ + y_] := Im[x] + Im[y]

Re[x_ y_] := Re[x] Re[y] - Im[x] Im[y]
Im[x_ y_] := Re[x] Im[y] + Im[x] Re[y]

reim[f_] := {Re[Expand[f/.x->x+I y]], Im[Expand[f/.x->x+I y]]} /.
                {Im[x]->0,Im[y]->0,Re[x]->x,Re[y]->y,
                Im[x^n_]->0,Im[y^n_]->0,Re[x^n_]->x^n,Re[y^n_]->y^n}

Protect[ Evaluate[protected] ]

EndPackage[ ]

ri[x_, y_] := RealInterval[{x, y}]

leftup[{ri[x1_, x2_], ri[y1_, y2_]}] :=
	{ri[x1, (x1+x2)/2.], ri[(y1+y2)/2., y2]} // N

leftdown[{ri[x1_, x2_], ri[y1_, y2_]}] :=
	{ri[x1, (x1+x2)/2.], ri[y1, (y1+y2)/2.]} // N

rightup[{ri[x1_, x2_], ri[y1_, y2_]}] :=
	{ri[(x1+x2)/2., x2], ri[(y1+y2)/2., y2]}  // N

rightdown[{ri[x1_, x2_], ri[y1_, y2_]}] :=
	{ri[(x1+x2)/2., x2], ri[y1, (y1+y2)/2.]}  // N

quadrasect[rect_] := {leftup[rect], leftdown[rect], 
			rightup[rect], rightdown[rect]}

containszero[{ri[x1_, x2_], ri[y1_, y2_]}] :=  
			(x1 <= 0 <= x2) && (y1 <= 0 <= y2)

F[{f1_,f2_}, {a_RealInterval,b_RealInterval}] := 
			{f1, f2} /. {x1->a, x2->b}

comroots[{f1_, f2_}, {x1_, a_, b_}, {x2_, c_, d_}, tol_:0.0001] :=
  With[{step = Select[Flatten[(quadrasect/@#),1], 
			containszero[N[F[{f1,f2},#]]]&]& },
	Nest[step, N[{{ri[a, b], ri[c, d]}}], 
			Ceiling[Log[2, N[Max[b-a,d-c]/tol]]]]
	]
  
cxroots[f_, {x1_, a_, b_}, {x2_, c_, d_}, tol_:0.00001] :=
		comroots[reim[f]/.{x->x1, y->x2}, 
			{x1, a, b}, {x2, c, d}, tol]








  • Prev by Date: Re: HypergeometricPFQ
  • Next by Date: Math typesetting in plot labels
  • Previous by thread: Save and big arrays
  • Next by thread: ComplexRoots