For the one-dimensional case it is simple enough to check for zero and give
an appropriate message.
I suspect that mnewton could be fixed as well. At least to give a better
message than divide by zero.
Here is 1-d newton code..
newton(exp,var,x0,eps):=block([xn,s,numer],numer:true,
s:diff(exp,var),xn:x0,
loop,if abs(subst(xn,var,exp)) < eps then return(xn),
xn:xn-subst(xn,var,exp)/subst(xn,var,s),
go(loop))$
a possible fix could be
newton(exp,var,x0,eps):=block([xn,s,numer,dv],numer:true, /*changed */
s:diff(exp,var),
xn:x0,
loop,if abs(subst(xn,var,exp)) < eps then return(xn),
dv:subst(xn,var,s),
if equal(dv,0) then error ("derivative is zero at ", " newton failed"),
/* changed */
xn:xn-subst(xn,var,exp)/ dv,
go(loop))$
> -----Original Message-----
> From: maxima-bounces at math.utexas.edu
> [mailto:maxima-bounces at math.utexas.edu] On Behalf Of reyssat
> Sent: Friday, April 04, 2008 5:15 AM
> To: Edwin Woollett
> Cc: maxima mailing list
> Subject: Re: [Maxima] mnewton divide by zero
>
> Edwin Woollett a ?crit :
> > lisp error message with mnewton:
> >
> > (%i18) eqns : [ x + y = 3, x^2 + y^2 = 9 ]$
> > (%i19) mn(x1, x2) := mnewton( eqns, [x, y], [x1, x2] )$
> > (%i20) solve( eqns, [x, y] );
> > (%o20) [[x = 3, y = 0], [x = 0, y = 3]]
> > (%i21) mn(1, 2);
> > (%o21) [ [x = 6.37714165E-17, y = 3.0] ]
> > ......
> > (%i24) mn(0, 0);
> > Maxima encountered a Lisp error:
> >
> > Error in FUNCALL [or a callee]: Zero divisor.
> >
> > Automatically continuing.
> > To reenable the Lisp debugger set *debugger-hook* to nil.
> > (%i25) mn(0, 1);
> > (%o25) [ [x = 1.13978659E-16, y = 3.0] ]
> >
> > Is this a bug?
> >
> In my opinion, this is not a bug but a feature of Newton's
> method. The
> same thing happens for 1-dimensional case :
>
> load(newton1);
> newton(cos(x),x,0,10^-7); --------> Division by 0
>
> The principle of Newton's method is to take the intersection point of
> the x-axis with the tangent at starting point and iterate. But if the
> starting point has an horizontal tangent, a division by zero
> occurs. So
> it is your responsability to start from a point where the
> derivative (or
> partial derivatives) is not zero. You could modify your
> starting point
> by a small random number for instance.
> But there is no guarantee of convergence anyway in general.
>
> Eric Reyssat
>
>
> > Ted Woollett
> > _______________________________________________
> > Maxima mailing list
> > Maxima at math.utexas.edu
> > http://www.math.utexas.edu/mailman/listinfo/maxima
> >
> >
> _______________________________________________
> Maxima mailing list
> Maxima at math.utexas.edu
> http://www.math.utexas.edu/mailman/listinfo/maxima
>