How to add AccuracyGoal and PrecisionGoal options in a numerical function?

Question

There are many built-in functions that contain AccuracyGoal, PrecisionGoal and WorkingPrecision.

For instance

NSolve

NSolve[x^3 + 2 x + 7 == 0, x, WorkingPrecision -> 40]

 {{x -> -1.568946403052382267352334751687751405502}, 
  {x -> 0.784473201526191133676167375843875702751 - 
     1.961171744579820388573529019015063895809 I}, 
  {x -> 0.784473201526191133676167375843875702751 + 
     1.961171744579820388573529019015063895809 I}}

FindRoot

FindRoot[Sin[x - 10] - x + 10, {x, 0}, AccuracyGoal -> 4, PrecisionGoal -> 4]

 {x -> 9.99852}

FindRoot[Cos[x^2] - x, {x, 1}, WorkingPrecision -> 40]

 {x -> 0.8010707652092183662168678540865655855422}

Here, I have a function Newton that using Newton method to find a numerical root of a equation.

Options[Newton] = {MaxIterations :> $RecursionLimit};

Newton::noconv = "Iteration did not converge in `1` steps.";

Newton::usage = "Newton[func, x0] finds a zero of the function func 
  using the initial guess x0 to start the iteration. "

Newton[func_, x0_, opts : OptionsPattern[]] :=
 Module[{res, maxIter, extraPre = 10},
  maxIter = OptionValue[MaxIterations];
  With[{prec = Precision[x0], fp = func'},
   Block[{$MaxPrecision = prec + extraPre},
    res = FixedPoint[(# - func[#]/fp[#]) &, x0, maxIter]
   ];
  If[! TrueQ[Abs@func[res] < 10^-prec],
    Message[Newton::noconv, maxIter]];
  res
 ]
]

Newton[expr_, x_, x0_, opts : OptionsPattern[]] :=
 Newton[Function[x, expr], x0, opts]

Newton Test

Newton[#^3 + 4 #^2 - 10 &, 1.2]

1.36523

Newton[x (x + 1), x, -1.2]

-1. 

Question

• How to add the options AccuracyGoal, PrecisionGoal and WorkingPrecisionto Newton. Namely,

  Options[Newton] = {MaxIterations :> $RecursionLimit, AccuracyGoal -> Automatic, PrecisionGoal -> Automatic, WorkingPrecision ->Automatic}

• Is there a general strategy for implementing the options AccuracyGoal, PrecisionGoal and WorkingPrecision in a numerical function?

However, I have no idea to achieve this.

In addition,

---

Update

As Jens said

This could also be done for precision instead. But right now I don't know how to write a general answer that would be applicable for all cases.

So is it possible to know how WRI implement these options(PrecisionGoal and AccuracyGoal ,and WorkingPrecision) in their built-ins?

Or how to add these options in my Newton function?

Answer 1

I'll mimic the precision/accuracy handling for FindRoot as indicated in its documentation:

• The default settings for AccuracyGoal and PrecisionGoal are WorkingPrecision/2.
• The setting for AccuracyGoal specifies the number of digits of accuracy to seek in both the value of the position of the root, and the value of the function at the root.
• The setting for PrecisionGoal specifies the number of digits of precision to seek in the value of the position of the root.

With both FindRoot and the OP's Newton there are two ways to measure the error, how close x is to the root and how close f[x] is to zero. Other functions, such as NIntegrate and NDSolve, have other notions of error.

The standard measure of error is indicated in the documentation for PrecisionGoal and AccuracyGoal.

With PrecisionGoal -> p and AccuracyGoal -> a, the Wolfram Language attempts to make the numerical error in a result of size x be less than 10^-a + |x| 10^-p.

With these hints we can construct a termination condition errortest for Newton. Note that since the size of the "result" func[x] is supposed to be zero, PrecisionGoal does not enter into its error measure; only AccuracyGoal is used, as indicated in the documentation for FindRoot.

I include a function checkopts to check the option settings and issue standard error messages. I also substituted the standard nonconvergence warning for the OP's Newton::noconv and reduced the default for MaxIterations.

ClearAll[Newton, checkopts];

Options[Newton] = {MaxIterations -> 100, AccuracyGoal -> Automatic, 
   PrecisionGoal -> Automatic, WorkingPrecision -> MachinePrecision};

Newton::usage = "Newton[func, x0] finds a zero of the function func 
    using the initial guess x0 to start the iteration. ";

checkopts[o : OptionsPattern[Newton]] :=
 With[
  {wp = OptionValue[WorkingPrecision],
   pg = OptionValue[PrecisionGoal],
   ag = OptionValue[AccuracyGoal],
   maxIter = OptionValue[MaxIterations]}, 
  If[! MatchQ[wp, Automatic | _?Positive],
   Message[Newton::wprec, wp]; Return[False]];
  If[! MatchQ[pg, Automatic | Infinity | _?Positive], 
   Message[Newton::precg, pg]; Return[False]];
  If[! MatchQ[ag, Automatic | Infinity | _?Positive], 
   Message[Newton::accg, ag]; Return[False]];
  If[! MatchQ[maxIter, Automatic | Infinity | _Integer?Positive], 
   Message[Newton::ioppfa, MaxIterations, maxIter]; Return[False]];
  True
]

Newton[func0_, x0_, opts : OptionsPattern[]] /; checkopts[opts] := 
  Module[{func, res, maxIter, ag, pg, wp, errortest, lasterror},

   (* initialization *)
   wp = OptionValue[WorkingPrecision] /. Automatic -> MachinePrecision;
   If[Precision[func] < wp,
    Message[Newton::precw, func, wp];
    func = SetPrecision[func0, wp],
    func = func0;
    ];
   pg = OptionValue[PrecisionGoal] /. Automatic -> wp/2;
   ag = OptionValue[AccuracyGoal] /. Automatic -> wp/2;
   maxIter = OptionValue[MaxIterations] /. Automatic -> 100;

   (* construct error test: False => goals are met *)
   errortest[x1_, x2_] := 
    lasterror = Abs[x1 - x2] >= 10^-ag +  x1 10^-pg || (* prec/acc of root *)
                Abs@func[x1] > 10^-ag;                 (* acc of function val *)

   (* iteration *)
   With[{fp = func'}, 
    res = NestWhile[
      (# - func[#]/fp[#]) &[
        SetPrecision[#, wp]] &,  (* reset precision of input at each step *)
      x0, 
      errortest, 2, maxIter]];

   (* check result and return *)
   If[lasterror, Message[Newton::cvmit, maxIter]];
   res
   ];

Examples:

FindRoot returns a solution with precision 20; this loses a slight amount of precision on the last calculation. I wasn't sure if this was worth fixing, or how exactly FindRoot fixes it, other than with SetPrecision. If you compare it to FindRoot[x^2 - 2, {x, 1.}, WorkingPrecision -> 20], Newton gives a more accurate answer.

Newton[#^2 - 2 &, 1., WorkingPrecision -> 20]

(*  1.4142135623730950488016896235025302436152453362821936500148`19.698970004335283  *)

Nonconvergence:

Newton[Exp, 1.]

Newton::cvmit: Failed to converge to the requested accuracy or precision within 100 iterations.

(*  -99.  *)

Bad option value:

Newton[Cos, 1., MaxIterations -> Pi]

Newton::ioppfa: Value of option MaxIterations -> π should be a positive integer, Infinity, or Automatic.

(*  Newton[Cos, 1., MaxIterations -> π]  *)

This cannot converge to the precision/accuracy goals with a working precision of only MachinePrecision (unless by accident).

Newton[Cos, 1., PrecisionGoal -> 17, AccuracyGoal -> 17, MaxIterations -> 1000]
Cos[%]

Newton::cvmit: Failed to converge to the requested accuracy or precision within 1000 iterations.

(*  1.5707963267948966`  *)    
(*  6.12323*10^-17  *)        (* <-- Cos[x] greater than accuracy goal 10^-17 *)