On 08/16/2010 10:37 AM, Ether Jones wrote:
> Why does this not work ? There is a simple integer solution.
>
> (%i1) load("mnewton")$
>
> (%i2) mnewton([a*b+c+d-3, b*c+a+d-5, c*d+a+b-2, a*d+b+c-6], [a,b,c,d], [1,1,1,1]);
>
> Unable to compute the LU factorization
>
> #0: solve_by_lu(
>
> Maxima encountered a Lisp error:
>
> Error in MERROR [or a callee]: NIL is not of type CHARACTER.
>
> Automatically continuing.
>
> To enable the Lisp debugger set *debugger-hook* to nil.
>
> (%i3)
>
>
>
>
> _______________________________________________
> Maxima mailing list
> Maxima at math.utexas.edu
> http://www.math.utexas.edu/mailman/listinfo/maxima
>
>
On 08/16/2010 10:37 AM, Ether Jones wrote:
> Why does this not work ? There is a simple integer solution.
>
> (%i1) load("mnewton")$
>
> (%i2) mnewton([a*b+c+d-3, b*c+a+d-5, c*d+a+b-2, a*d+b+c-6], [a,b,c,d],
> [1,1,1,1]);
>
> Unable to compute the LU factorization
>
> #0: solve_by_lu(
>
> Maxima encountered a Lisp error:
>
> Error in MERROR [or a callee]: NIL is not of type CHARACTER.
>
> Automatically continuing.
>
> To enable the Lisp debugger set *debugger-hook* to nil.
>
> (%i3)
>
>
>
>
> _______________________________________________
> Maxima mailing list
> Maxima at math.utexas.edu
> http://www.math.utexas.edu/mailman/listinfo/maxima
>
The Jacobian determinant of your function f(a,b,c,d) is 0 at
[1,1,1,1]. Any newton method looking for a solution to f(x)=0, with
solution x_0 and initial guess x_1 requires that the derivative of the
function f(x) be invertible in a neighborhood U of x_0 and x_1 be in U.
The routine may work in some peculiar cases when the jacobian
determinant vanished, but it may break down as you observed.
You equations can be solved using 'solve'.
(%i45) solve([a*b+c+d-3, b*c+a+d-5, c*d+a+b-2, a*d+b+c-6], [a,b,c,d]);
(%o45) [[a = 2,b = 0,c = 0,d = 3]]
Curiously, for some guesses near this solution, mnewton does get close
to the solution.
(%i65) mnewton([a*b+c+d-3, b*c+a+d-5, c*d+a+b-2, a*d+b+c-6], [a,b,c,d],
[2,0,0,3]);
(%o65) [[a = 2.0,b = 0.0,c = 0.0,d = 3.0]]
(%i66) mnewton([a*b+c+d-3, b*c+a+d-5, c*d+a+b-2, a*d+b+c-6], [a,b,c,d],
[2,.2,.4,3]);
(%o66) [[a = 2.0,b = 1.182625636315788e-18,c = 3.549975962435521e-17,d =
3.0]]
-sen