jacobwilliams / roots-fortran

A modern Fortran library for finding the roots of continuous scalar functions of a single real variable, using derivative-free methods.
Other
37 stars 6 forks source link
bisection brent-dekker fortran fortran-package-manager muller muller-s-method regula-falsi root-finding zeroin

roots-fortran

roots-fortran: root solvers for modern Fortran

Language GitHub release CI Status codecov last-commit

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:

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

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

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

See also