0
$\begingroup$

My objective is to build a function using modules that runs the bisection method. I have a working bisection method but when I try to implement it as a module I keep getting this error and I can't figure out why.

Small thing to mention:
I'm aware the computations of the code are not giving me the bisection method EXACTLY. but the code itself works outside of the module it's own loop and does what its meant to do which is the important part.

F[z_]:=3 + z + z^2 * Cos[z]

bisectionMethod[f_, a_, b_, n_] := 
  Module[{},
    For[i = 0, i < n, i++,
      Print[StringForm["Loop #: ``", i + 1]];
      c = (a + b)/2;
      y = f[c];
      If[y < 0, a = c];
      If[y > 0, b = c];
      Print[StringForm["a = ``", N[a]]];
      Print[StringForm["b  = ``", N[b]]];
      c = (a + b)/2;
      Print[StringForm["c = ``", N[c]]];
      Print[N[f[c]]];]]

bisectionMethod[F, -4 , 2 , 10]

After being run

Loop #: 1
Set::setraw: Cannot assign to raw object 2.
a = -4.
b = 2.`
c = -1.
2.5403

Loop #: 2
Set::setraw: Cannot assign to raw object 2.
a = -4.
b = 2.`
c = -1.
2.5403

The issue appears to be specific to the variable b, but I have no idea how b is any different from the others or how to tackle fixing this.

I did a little bit of experimentation and it seems like the line

If[y > 0, b = c];

is problematic

$\endgroup$
1
  • $\begingroup$ In the code for bisectionMethod[F, -4 , 2 , 10], a represents, not a variable, but the raw object -4, and b represents, not a variable, but the raw object 2. Declare some variables in your Module and initialize them to a and b. $\endgroup$
    – Michael E2
    Oct 24 '20 at 1:24
1
$\begingroup$

Mathematica functions without attributes that give them non-standard evaluation pass in their arguments by value. Hence, a, and b in your function call are passed in as -4 and 2, which are variables and which can not be assigned to. The solution is to define local variables in your module and initialize them with the values passed in. Like so;

bisectionMethod[f_, a_, b_, n_] :=
  Module[{aa, bb, c, y},
    aa = a; bb = b;
    c = (aa + bb)/2;
    For[i = 0, i < n, i++,
    Print[StringForm["Loop #: ``", i + 1]];
    Print[StringForm["a = ``", N[aa]]];
    Print[StringForm["b  = ``", N[bb]]];
    y = f[c];
    If[y < 0, aa = c, If[y > 0, bb = c, Break[]]];
    c = (aa + bb)/2;
    Print[StringForm["c = ``", N[c]]];
    Print[N[f[c]]]]]

F[z_] := 3 + z + z^2 Cos[z]

bisectionMethod[F, -4, 2, 15]

first

...

last

But here is better (but still no really good) way to implement bisectionMethod:

Clear[bisectionMethod]
bisectionMethod[f_, a_, b_, ϵ_ : 10.^-6] :=
  Module[{aa, bb, c, i, y},
    aa = N[a]; bb = N[b];
    y = f[c = (aa + bb)/2];
    i = 0;
    While[Abs[y] > ϵ,
      If[y < 0, aa = c, If[y >= 0, bb = c]];
      y = f[c = (aa + bb)/2];
      Print @ 
        Row[{++i, ": a = ", aa, ", b = ", bb, ", c = ", c, ", f[c] = ", y}]]]

bisectionMethod[F, -4, 2]
1: a = -4., b = -1., c = -2.5, f[c] = -4.50715
2: a = -2.5, b = -1., c = -1.75, f[c] = 0.704121
...   
21: a = -1.88803, b = -1.88803, c = -1.88803, f[c] = 4.46423×10^-6
22: a = -1.88803, b = -1.88803, c = -1.88803, f[c] = 4.84041×10^-7

Note that in this version the function completes when a convergence test is passed. This is a better way to determine completion than using an arbitrary loop count.

Update

Here is an implementation of bisectionMethod that follows better Mathematica practice. It is robust in the sense that it would be easier to maintain or repurpose then the preceding implementations. It is also more efficient (runs faster).

Clear[bisectionMethod]
bisectionMethod[f_, a_, b_, ϵ_ : 10.^-6] :=
  Block[{h, pts},
    h =
      Module[{aa, bb, c},
        aa = #[[1]]; bb = #[[2]]; c = #[[3]];
        With[{y = f[c]},
          If[y < 0, {c, bb, (c + bb)/2}, If[y >= 0, {aa, c, (aa + c)/2}]]]]&;
    pts =
      FixedPointList[h, {a, b, (a + b)/2} // N, 
        SameTest -> (Abs[#1[[3]] - #2[[3]]] < ϵ &)];
    Column @
      MapIndexed[
        Row[
          {#2[[1]], ": a = ", #1[[1]], ", b = ", #1[[2]], ", c = ", #1[[3]], 
           ", f[c] = ", f[#1[[3]]]}] &, pts]]
bisectionMethod[F, -4, 2]

output

To convince you of its robustness, see how easy it is to convert it from a program for studying the behavior of the bisection method to a real solver.

bisectionSlover[f_, a_, b_, ϵ_ : 10.^-6] := 
  Block[{h},
    h =
      Module[{aa, bb, c},
        aa = #[[1]]; bb = #[[2]]; c = #[[3]];
        With[{y = f[c]},
          If[y < 0, {c, bb, (c + bb)/2}, If[y >= 0, {aa, c, (aa + c)/2}]]]]&;
    FixedPoint[h, {a, b, (a + b)/2} // N, 
      SameTest -> (Abs[#1[[3]] - #2[[3]]] < ϵ &)][[3]]]
With[{x = bisectionSlover[F, -4, 2]}, Row[{"F[x] = ", F[x], "at x = ", x}]]

solver

$\endgroup$

Not the answer you're looking for? Browse other questions tagged or ask your own question.