[genius] Sun Oct 20 23:44:25 2013 Jiri (George) Lebl <jirka 5z com>

[George Lebl <jirka src gnome org>]

commit 5319125087903f2199afd83e7aed57c91f6b9344
Author: Jiri (George) Lebl <jirka 5z com>
Date:   Sun Oct 20 23:44:28 2013 -0500

    Sun Oct 20 23:44:25 2013  Jiri (George) Lebl <jirka 5z com>
    
        * lib/equation_solving/newton.gel: Add HalleysMethod
    
        * src/geniustests.txt: do testing
    
        * help/C/genius.xml: Document Halleys method

 ChangeLog                       |    8 ++++++++
 NEWS                            |    2 +-
 help/C/genius.xml               |   26 ++++++++++++++++++++++++++
 lib/equation_solving/newton.gel |   21 ++++++++++++++++++++-
 lib/library-strings.c           |    1 +
 src/geniustests.txt             |    2 ++
 6 files changed, 58 insertions(+), 2 deletions(-)
---
diff --git a/ChangeLog b/ChangeLog
index 3ac0d93..7ac82f5 100644
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,11 @@
+Sun Oct 20 23:44:25 2013  Jiri (George) Lebl <jirka 5z com>
+
+       * lib/equation_solving/newton.gel: Add HalleysMethod
+
+       * src/geniustests.txt: do testing
+
+       * help/C/genius.xml: Document Halleys method
+
 Wed Oct 16 08:11:49 2013  Jiri (George) Lebl <jirka 5z com>
 
        * lib/functions/lambert.gel: slight improvement in speed.
diff --git a/NEWS b/NEWS
index 0c8f827..e16067a 100644
--- a/NEWS
+++ b/NEWS
@@ -1,7 +1,7 @@
 Changes to 1.0.18
 
 * New functions: KroneckerProduct (alias TensorProduct), NewtonsMethod,
-  LambertW, LambertWm1
+  HalleysMethod, LambertW, LambertWm1
 * New PlotCanvasFreeze/PlotCanvasThaw functions to improve flicker
   if doing animations with genius.
 * Allow setting color in LinePlotDrawLine with RGB vector.
diff --git a/help/C/genius.xml b/help/C/genius.xml
index e1cf9dc..aa1e4f8 100644
--- a/help/C/genius.xml
+++ b/help/C/genius.xml
@@ -7734,6 +7734,32 @@ and has period <userinput>b-a</userinput>.</para>
         </varlistentry>
 
         <varlistentry>
+         <term><anchor id="gel-function-HalleysMethod"/>HalleysMethod</term>
+         <listitem>
+          <synopsis>HalleysMethod (f,df,ddf,guess,epsilon,maxn)</synopsis>
+         <para>Find zeros using Halleys's method.  <varname>f</varname> is
+                 the function, <varname>df</varname> is the derivative of
+                 <varname>f</varname>, and <varname>ddf</varname> is the second derivative of
+                 <varname>f</varname>.  <varname>guess</varname> is the initial
+                 guess.  The function returns after two successive values are
+                 within <varname>epsilon</varname> of each other, or after <varname>maxn</varname> tries, in 
which case the function returns <constant>null</constant> indicating failure.
+         </para>
+         <para>
+         See also <link linkend="gel-function-NewtonsMethod"><function>NewtonsMethod</function></link> and 
<link linkend="gel-function-SymbolicDerivative"><function>SymbolicDerivative</function></link>.
+         </para>
+         <para>
+           Example to find the square root of 10:
+          <screen><prompt>genius></prompt> 
<userinput>HalleysMethod(`(x)=x^2-10,`(x)=2*x,`(x)=2,3,10^-10,100)</userinput>
+</screen>
+         </para>
+          <para>
+           See
+           <ulink url="http://en.wikipedia.org/wiki/Halley%27s_method";>Wikipedia</ulink> for more 
information.
+         </para>
+         </listitem>
+        </varlistentry>
+
+        <varlistentry>
          <term><anchor id="gel-function-NewtonsMethod"/>NewtonsMethod</term>
          <listitem>
           <synopsis>NewtonsMethod (f,df,guess,epsilon,maxn)</synopsis>
diff --git a/lib/equation_solving/newton.gel b/lib/equation_solving/newton.gel
index 34449fa..51f9d88 100644
--- a/lib/equation_solving/newton.gel
+++ b/lib/equation_solving/newton.gel
@@ -1,4 +1,4 @@
-# Newton's method
+# Newton's method and related
 
 SetHelp("NewtonsMethod","equation_solving","Attempt to find a zero of a functionf with derivative df using 
Newton's method, returning after two successive values are within epsilon or after maxn tries (then returns 
null)")
 function NewtonsMethod(f,df,guess,epsilon,maxn) = (
@@ -17,4 +17,23 @@ function NewtonsMethod(f,df,guess,epsilon,maxn) = (
        null
 )
 
+SetHelp("HalleysMethod","equation_solving","Attempt to find a zero of a functionf with derivative df and 
second derivative ddf using Halleys's method, returning after two successive values are within epsilon or 
after maxn tries (then returns null)")
+function HalleysMethod(f,df,ddf,guess,epsilon,maxn) = (
+       guess := float(guess);
+       for n=1 to maxn do (
+               fg = f(guess);
+               dfg = df(guess);
+               denom := 2*dfg^2 - fg*ddf(guess);
+               if denom == 0.0 then (
+                       error ("HalleysMethod: division by zero");
+                       bailout
+               );
+               guessn := guess - (2*fg*dfg)/denom;
+               if |guessn-guess| <= epsilon then
+                       return guessn;
+               guess := guessn
+       );
+       null
+)
+
 
diff --git a/lib/library-strings.c b/lib/library-strings.c
index 103df5f..7b962c7 100644
--- a/lib/library-strings.c
+++ b/lib/library-strings.c
@@ -237,6 +237,7 @@ char *fake = N_("Find root of a function using the bisection method to within TO
 char *fake = N_("Find root of a function using the method of false position to within TOL tolerance in up to 
N iterations.  f(a) and f(b) must have opposite signs.");
 char *fake = N_("Find root of a function using the Muller's method");
 char *fake = N_("Find root of a function using the secant method to within TOL tolerance in up to N 
iterations.  f(a) and f(b) must have opposite signs.");
+char *fake = N_("Attempt to find a zero of a functionf with derivative df and second derivative ddf using 
Halleys's method, returning after two successive values are within epsilon or after maxn tries (then returns 
null)");
 char *fake = N_("Attempt to find a zero of a functionf with derivative df using Newton's method, returning 
after two successive values are within epsilon or after maxn tries (then returns null)");
 char *fake = N_("Find roots of a polynomial (given as vector of coefficients)");
 char *fake = N_("Find roots of a quartic polynomial (given as vector of coefficients)");
diff --git a/src/geniustests.txt b/src/geniustests.txt
index 9d49d36..863536a 100644
--- a/src/geniustests.txt
+++ b/src/geniustests.txt
@@ -1137,10 +1137,12 @@ KroneckerProduct(null,[1,2;3,4])+0                              ((null)+0)
 KroneckerProduct([0],null)+0                                   ((null)+0)
 |NewtonsMethod (`(x)=x^2-10,`(x)=2*x,1.0,0.0001,100)-sqrt(10)|<0.0001  true
 |NewtonsMethod (`(x)=x^2-2,`(x)=2*x,1.0,0.00001,100)-sqrt(2)|<0.00001  true
+|HalleysMethod (`(x)=x^2-2,`(x)=2*x,`(x)=2,1.0,0.00001,100)-sqrt(2)|<0.00001   true
 |NewtonsMethod (`(x)=x*e^x-0.5,`(x)=x*e^x+e^x,1.0,0.0001,100) - LambertW(0.5)|<0.0001  true
 |NewtonsMethod (`(x)=x*e^x+0.1,`(x)=x*e^x+e^x,1.0,0.0001,100) - LambertW(-0.1)|<0.0001 true
 |NewtonsMethod (`(x)=x*e^x-10,`(x)=x*e^x+e^x,1.0,0.0001,100) - LambertW(10)|<0.0001    true
 |NewtonsMethod (`(x)=x*e^x+0.1,`(x)=x*e^x+e^x,-2.0,0.0001,100) - LambertWm1(-0.1)|<0.0001      true
+|HalleysMethod (`(x)=x*e^x-0.5,`(x)=x*e^x+e^x,`(x)=x*e^x+2*e^x,1.0,0.0001,100) - LambertW(0.5)|<0.0001 true
 load "nullspacetest.gel"                                       true
 load "longtest.gel"                                            true
 load "testprec.gel"                                            true