Math-Numerical

 view release on metacpan or  search on metacpan

lib/Math/Numerical.pm  view on Meta::CPAN

The current implementation of this function is based on the Brent method
described in the
I<L<Numerical Recipes Third Edition|http://numerical.recipes/aboutNR3book.html>>
book, section 9.3.

The function supports the following parameters:

=over

=item C<max_iterations>

How many iterations of our algorithm will be applied at most while trying to
find a root for the given function. This gives an order of magnitude of the
number of times that C<$func> will be evaluated. Defaults to I<100>.

=item C<do_bracket>

Whether the C<L<bracket|/bracket>> function should be used to bracket a root of
the function before finding the root. If this is set to a false value, then the
passed C<$x1> and C<$x2> values must already form a bracket (that is, the
function must take values of opposite sign at these two points). Note that, when
they do, the C<L<bracket|/bracket>> function will immediately return these

lib/Math/Numerical.pm  view on Meta::CPAN

    $s->{b} += $s->{d};
  } else {
    $s->{b} += _sign($s->{tol1}, $s->{xm});
  }
  $s->{fb} = $f->($s->{b});
  return 0;
}

sub find_root ($func, $x1, $x2, %params) {
  my $do_bracket = $params{do_bracket} // 1;
  my $max_iter = $params{max_iterations} // $DEFAULT_MAX_ITERATIONS;
  my $f = _wrap_func($func);
  my ($xa, $xb, $fa, $fb);
  if ($do_bracket) {
    ($xa, $xb, $fa, $fb) = bracket($func, $x1, $x2, %params);
    croak 'Can’t bracket a root of the function' unless defined $xa;
  } else {
    ($xa, $xb) = ($x1, $x2);
    ($fa, $fb) = ($f->($xa), $f->($xb));
    croak 'A root must be bracketed in [\$x1; \$x2]'
        if ($fa > 0 && $fb > 0) || ($fa < 0 && $fb < 0);

lib/Math/Numerical.pm  view on Meta::CPAN

The current implementation is a mix of the inward and outward bracketing
approaches exposed in the
I<L<Numerical Recipes Third Edition|http://numerical.recipes/aboutNR3book.html>>
book, section 9.1.

The function supports the following parameters:

=over

=item C<max_iterations>

How many iterations of our algorithm will be applied at most while trying to
bracket the given function. This gives an order of magnitude of the number of
times that C<$func> will be evaluated. Defaults to I<100>.

=item C<do_outward>

Whether the function will try to bracket a root in an interval larger than the
one given by C<[$x1, $x2]>. Defaults to I<1>.

=item C<do_inward>

lib/Math/Numerical.pm  view on Meta::CPAN

    $s->{f2} = $f->($s->{x2});
  }
  return 0;
}

sub bracket ($func, $x1, $x2 = undef, %params) {
  if (!defined $x2 || $x1 == $x2) {
    Readonly my $LARGISH_FACTOR => 1000;
    $x2 += $LARGISH_FACTOR * $EPS;
  }
  my $max_iter = $params{max_iterations} // $DEFAULT_MAX_ITERATIONS;
  croak 'max_iterations must be positive' if $max_iter <= 0;

  my $f = _wrap_func($func);
  my $f1 = $f->($x1);

  my $inward_state;
  if ($params{do_inward} // 1) {
    $inward_state = _create_bracket_inward_state($x1, $x2, $f1, %params);
  }
  my $outward_state;
  if ($params{do_outward} // 1) {

t/bracket.t  view on Meta::CPAN

Readonly my $PI => 4 * atan2(1, 1);

is([bracket(\&CORE::cos, 0, 1)], [float_lt($PI / 2), float_gt($PI / 2), D(), D()]);

is([bracket(\&CORE::cos, 0, 0)], [float_lt($PI / 2), float_gt($PI / 2), D(), D()]);
is([bracket(\&CORE::cos, 0)], [float_lt($PI / 2), float_gt($PI / 2), D(), D()]);

like(scalar(eval { bracket(\&CORE::cos, 0, 0, do_inward => 0, do_outward => 0)}, $@),
     qr/One of do_outward and do_inward at least must be true/);

like(scalar(eval { bracket(\&CORE::cos, 0, 0, max_iterations => 0)}, $@),
     qr/max_iterations must be positive/);

like(scalar(eval { bracket(\&CORE::cos, 0, 0, outward_factor => 1)}, $@),
     qr/outward_factor must be larger than 1/);

like(scalar(eval { bracket(\&CORE::cos, 0, 0, inward_factor => 1.9)}, $@),
     qr/inward_factor must be at least 2/);

like(scalar(eval { bracket(\&CORE::cos, 0, 0, inward_split => 1.9)}, $@),
     qr/inward_split must be at least 2/);