roots-fortran: root solvers for modern Fortran

Description

A modern Fortran library for finding the roots of continuous scalar functions of a single real variable.

Compiling

A fpm.toml file is provided for compiling roots-fortran with the Fortran Package Manager. For example, to build:

fpm build --profile release

By default, the library is built with double precision (real64) real values. Explicitly specifying the real kind can be done using the following processor flags:

Preprocessor flag	Kind	Number of bytes
REAL32	real(kind=real32)	4
REAL64	real(kind=real64)	8
REAL128	real(kind=real128)	16

For example, to build a single precision version of the library, use:

fpm build --profile release --flag "-DREAL32"

To run the unit tests:

fpm test

To use roots-fortran within your fpm project, add the following to your fpm.toml file:

[dependencies]
roots-fortran = { git="https://github.com/jacobwilliams/roots-fortran.git" }

or, to use a specific version:

[dependencies]
roots-fortran = { git="https://github.com/jacobwilliams/roots-fortran.git", tag = "1.0.0"  }

To generate the documentation using ford, run: ford roots-fortran.md

Usage

Methods

The module contains the following methods:

Procedure	Description	Reference
bisection	Classic bisection method	Bolzano (1817)
regula_falsi	Classic regula falsi method	?
muller	Improved Muller method (for real roots only)	Muller (1956)
brent	Classic Brent's method (a.k.a. Zeroin)	Brent (1971)
brenth	SciPy variant of brent	SciPy
brentq	SciPy variant of brent	SciPy
illinois	Illinois method	Dowell & Jarratt (1971)
pegasus	Pegasus method	Dowell & Jarratt (1972)
anderson_bjorck	Anderson-Bjorck method	King (1973)
anderson_bjorck_king	a variant of anderson_bjorck	King (1973)
ridders	Classic Ridders method	Ridders (1979)
toms748	Algorithm 748	Alefeld, Potra, Shi (1995)
chandrupatla	Hybrid quadratic/bisection algorithm	Chandrupatla (1997)
bdqrf	Bisected Direct Quadratic Regula Falsi	Gottlieb & Thompson (2010)
zhang	Zhang's method (with corrections from Stage)	Zhang (2011)
itp	Interpolate Truncate and Project method	Oliveira & Takahashi (2020)
barycentric	Barycentric interpolation method	Mendez & Castillo (2021)
blendtf	Blended method of trisection and false position	Badr, Almotairi, Ghamry (2021)

In general, all the methods are guaranteed to converge. Some will be more efficient (in terms of number of function evaluations) than others for various problems. The methods can be broadly classified into three groups:

• Simple classical methods (bisection, regula_falsi, illinois, ridders).
• Newfangled methods (zhang, barycentric, blendtf, bdqrf, anderson_bjorck_king). These rarely or ever seem to be better than the best methods.
• Best methods (anderson_bjorck, muller, pegasus, toms748, brent, brentq, brenth, chandrupatla, itp). Generally, one of these will be the most efficient method.

Note that some of the implementations in this library contain additional checks for robustness, and so may behave better than naive implementations of the same algorithms. In addition, all methods have an option to fall back to bisection if the method fails to converge.

Functional Interface Example

program main

  use root_module, wp => root_module_rk

  implicit none

  real(wp) :: x, f
  integer :: iflag

  call root_scalar('bisection',func,-9.0_wp,31.0_wp,x,f,iflag)

  write(*,*) 'f(',x,') = ', f
  write(*,*) 'iflag    = ', iflag

contains

  function func(x) result(f)

  implicit none

  real(wp),intent(in) :: x
  real(wp) :: f

  f = -200.0_wp * x * exp(-3.0_wp*x)

  end function func

end program main

Object Oriented Interface Example

program main

  use root_module, wp => root_module_rk

  implicit none

  type,extends(bisection_solver) :: my_solver
  end type my_solver

  real(wp) :: x, f
  integer :: iflag
  type(my_solver) :: solver

  call solver%initialize(func)
  call solver%solve(-9.0_wp,31.0_wp,x,f,iflag)

  write(*,*) 'f(',x,') = ', f
  write(*,*) 'iflag    = ', iflag

contains

  function func(me,x)

    class(root_solver),intent(inout) :: me
    real(wp),intent(in) :: x
    real(wp) :: f

    f = -200.0_wp * x * exp(-3.0_wp*x)

  end function func

end program main

Result

 f( -2.273736754432321E-013 ) =   4.547473508867743E-011
 iflag    =            0

Notes

• Originally this was an idea for the Fortran stdlib. See: #87: Optimization, Root finding, and Equation Solvers. Eventually, it may be merged into the other one in one form or another.

Documentation

The latest API documentation for the master branch can be found here. This was generated from the source code using FORD.

License

The roots-fortran source code and related files and documentation are distributed under a permissive free software license (BSD-style).

Related Fortran libraries

• polyroots-fortran -- For finding all the roots of polynomials with real or complex coefficients.

Similar libraries in other programming languages

Language	Library
C	GSL
C++	Boost Math Toolkit
Julia	Roots.jl
R	uniroot
Rust	roots
MATLAB	fzero
Python	scipy.optimize.root_scalar

References

• D. E. Muller, "A Method for Solving Algebraic Equations Using an Automatic Computer", Mathematical Tables and Other Aids to Computation, 10 (1956), 208-215.
• M. Dowell, P. Jarratt, "A modified regula falsi method for computing the root of an equation", BIT 11 (1971), 168-174.
• R. P. Brent, "An algorithm with guaranteed convergence for finding a zero of a function", The Computer Journal, Vol 14, No. 4., 1971.
• R. P. Brent, "Algorithms for minimization without derivatives", Prentice-Hall, Inc., 1973.
• Ridders, C., "A new algorithm for computing a single root of a real continuous function", IEEE Trans. on Circuits and Systems, Vol 26, Issue 11, Nov 1979.
• G. E. Alefeld, F. A. Potra and Yixun Shi, "Algorithm 748: Enclosing Zeros of Continuous Functions", ACM Transactions on Mathematical Software, Vol. 21. No. 3. September 1995. Pages 327-344.
• G.E. Mullges & F. Uhlig, "Numerical Algorithms with Fortran", Springer, 1996. Section 2.8.1, p 32-34.
• T.R. Chandrupatla, "A new hybrid quadratic/bisection algorithm for finding the zero of a nonlinear function without derivatives", Advances in Engineering Software, Vol 28, 1997, pp. 145-149.
• R. G. Gottlieb, B. F. Thompson, "Bisected Direct Quadratic Regula Falsi", Applied Mathematical Sciences, Vol. 4, 2010, no. 15, 709-718.
• A. Zhang, "An Improvement to the Brent's Method", International Journal of Experimental Algorithms (IJEA), Volume (2) : Issue (1) : 2011.
• E Badr, S Almotairi, A El Ghamry, "A Comparative Study among New Hybrid Root Finding Algorithms and Traditional Methods", Mathematics 2021, 9, 1306.

See also

• Code coverage statistics [codecov.io]