[sage-support] Re: find_root reports zeroes of denominator as roots

2009-03-16 Thread jpc

Using symbolic computation:

sage: x=var('x')
sage: solve(1/x==0)
[ ]
sage:

However, and like matlab fzero command, the numerical answer of the
equation is $x \approx 0$ because of the change of sign of $1/x \in
[-1,1]$.

Maybe the find_root documentation needs to mention that continuity
of $f$ is necessary.

jpc



On Mar 15, 3:54 am, yacob imya...@gmail.com wrote:
 f(x)=1/x
 r1=find_root(f,-1,1);r1
 -1.8189894035440371e-12

--~--~-~--~~~---~--~~
To post to this group, send email to sage-support@googlegroups.com
To unsubscribe from this group, send email to 
sage-support-unsubscr...@googlegroups.com
For more options, visit this group at 
http://groups.google.com/group/sage-support
URLs: http://www.sagemath.org
-~--~~~~--~~--~--~---



[sage-support] Re: find_root reports zeroes of denominator as roots

2009-03-16 Thread mark mcclure

On Mar 16, 6:14 am, jpc pedrocruzave...@gmail.com wrote:
 like matlab fzero command, the numerical answer of the
 equation is $x \approx 0$ because of the change of sign of
 $1/x \in [-1,1]$.


This is not really a complete description of Matlab's behavior
in this situation.  If you simply define f(x) = 1/x and apply
fzero, without requesting any additional information, you get:

 f = @(x) 1/x;
 fzero(f,-1,1)
ans =
  -2.6773e-16

However, you can ask Matlab to provide more information via:

 [x,fval,exitflag,output] = fzero(f,-1,1);

Now, fval is a huge number and output.message indicates that
the function might be singular.  Similarly, Mathematica
yields a numerical zero but issues a warning.


Now, you can find that Sage uses SciPy's implementation
of Brent's algorithm when tackling this problem by typing
find_root?? and SciPy provides similar functionality as
Matlab, in this case.  What I find particularly puzzling,
though, is that scipy.opitimize.brentq, raises an error
for this problem, at least when applied to the interval
[-1.0, 1.0].  On the other hand, brentq runs without
reporting a complaint when applied to the interval
[1.0, pi].  See below.

Mark McClure


from scipy.optimize import brentq
def f(x): return 1/x
[x,info] = brentq(f, -1.0, pi, full_output=True)
[x,info.converged]

   [-1.155868814851219e-12, True]



brentq(f, -1.0, 1.0, full_output=True)


Traceback (most recent call last):
  File stdin, line 1, in module
  File /Users/markmcclure/.sage/sage_notebook/worksheets/admin/31/
code/16.py, line 9, in module
[a,b] = brentq(f, -_sage_const_1p0 , _sage_const_1p0 ,
full_output=True)
  File /Applications/Sage.app/Contents/Resources/sage/local/lib/
python2.5/site-packages/SQLAlchemy-0.4.6-py2.5.egg/, line 1, in
module

  File /Applications/Sage.app/Contents/Resources/sage/local/lib/
python2.5/site-packages/scipy/optimize/zeros.py, line 223, in brentq
r = _zeros._brentq(f,a,b,xtol,maxiter,args,full_output,disp)
  File /Users/markmcclure/.sage/sage_notebook/worksheets/admin/31/
code/16.py, line 8, in f
def f(x): return _sage_const_1 /x
  File element.pyx, line 1201, in
sage.structure.element.RingElement.__div__ (sage/structure/element.c:
9116)
  File coerce.pyx, line 672, in
sage.structure.coerce.CoercionModel_cache_maps.bin_op (sage/structure/
coerce.c:5437)
ZeroDivisionError: float division


--~--~-~--~~~---~--~~
To post to this group, send email to sage-support@googlegroups.com
To unsubscribe from this group, send email to 
sage-support-unsubscr...@googlegroups.com
For more options, visit this group at 
http://groups.google.com/group/sage-support
URLs: http://www.sagemath.org
-~--~~~~--~~--~--~---