ParaMonte Fortran 2.0.0
Parallel Monte Carlo and Machine Learning Library
See the latest version documentation.
pm_mathRoot.F90
Go to the documentation of this file.
1!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3!!!! !!!!
4!!!! ParaMonte: Parallel Monte Carlo and Machine Learning Library. !!!!
5!!!! !!!!
6!!!! Copyright (C) 2012-present, The Computational Data Science Lab !!!!
7!!!! !!!!
8!!!! This file is part of the ParaMonte library. !!!!
9!!!! !!!!
10!!!! LICENSE !!!!
11!!!! !!!!
12!!!! https://github.com/cdslaborg/paramonte/blob/main/LICENSE.md !!!!
13!!!! !!!!
14!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
15!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
16
451
452!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
453
455
456 use pm_kind, only: SK, IK, LK
457
458 implicit none
459
460 character(*, SK), parameter :: MODULE_NAME = "@pm_mathRoot"
461
462!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
463
500 type, abstract :: method_type
501 end type
502
503!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
504
541 type, abstract, extends(method_type) :: bracket_type
542 end type
543
544!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
545
582 type, abstract, extends(method_type) :: iteration_type
583 end type
584
585!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
586
623 type, abstract, extends(method_type) :: hybrid_type
624 end type
625
626!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
627
670 type, extends(hybrid_type) :: brent_type
671 end type
672
709 type(brent_type), parameter :: brent = brent_type()
710#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
711 !DIR$ ATTRIBUTES DLLEXPORT :: brent
712#endif
713
714!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
715
758 type, extends(hybrid_type) :: toms748_type
759 end type
760
797 type(toms748_type), parameter :: toms748 = toms748_type()
798#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
799 !DIR$ ATTRIBUTES DLLEXPORT :: toms748
800#endif
801
802!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
803
846 type, extends(bracket_type) :: false_type
847 end type
848
885 type(false_type), parameter :: false = false_type()
886#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
887 !DIR$ ATTRIBUTES DLLEXPORT :: false
888#endif
889
890!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
891
934 type, extends(iteration_type) :: secant_type
935 end type
936
973 type(secant_type), parameter :: secant = secant_type()
974#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
975 !DIR$ ATTRIBUTES DLLEXPORT :: secant
976#endif
977
978!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
979
1022 type, extends(iteration_type) :: newton_type
1023 end type
1024
1061 type(newton_type), parameter :: newton = newton_type()
1062#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1063 !DIR$ ATTRIBUTES DLLEXPORT :: newton
1064#endif
1065
1066!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1067
1110 type, extends(iteration_type) :: halley_type
1111 end type
1112
1149 type(halley_type), parameter :: halley = halley_type()
1150#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1151 !DIR$ ATTRIBUTES DLLEXPORT :: halley
1152#endif
1153
1154!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1155
1199 end type
1200
1237 type(schroder_type), parameter :: schroder = schroder_type()
1238#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1239 !DIR$ ATTRIBUTES DLLEXPORT :: schroder
1240#endif
1241
1242!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1243
1286 type, extends(hybrid_type) :: ridders_type
1287 end type
1288
1325 type(ridders_type), parameter :: ridders = ridders_type()
1326#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1327 !DIR$ ATTRIBUTES DLLEXPORT :: ridders
1328#endif
1329
1330!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1331
1375 end type
1376
1414#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1415 !DIR$ ATTRIBUTES DLLEXPORT :: bisection
1416#endif
1417
1418!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1419
1586
1587 ! Default
1588
1589 interface getRoot
1590
1591 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1592
1593#if RK5_ENABLED
1594 recursive module function getRootDef_RK5(getFunc, lb, ub, abstol, neval, niter) result(root)
1595#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1596 !DEC$ ATTRIBUTES DLLEXPORT :: getRootDef_RK5
1597#endif
1598 use pm_kind, only: RKG => RK5
1599 procedure(real(RKG)) :: getFunc
1600 real(RKG) , intent(in) :: lb, ub
1601 real(RKG) , intent(in) , optional :: abstol
1602 integer(IK) , intent(in) , optional :: niter
1603 integer(IK) , intent(out) , optional :: neval
1604 real(RKG) :: root
1605 end function
1606#endif
1607
1608#if RK4_ENABLED
1609 recursive module function getRootDef_RK4(getFunc, lb, ub, abstol, neval, niter) result(root)
1610#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1611 !DEC$ ATTRIBUTES DLLEXPORT :: getRootDef_RK4
1612#endif
1613 use pm_kind, only: RKG => RK4
1614 procedure(real(RKG)) :: getFunc
1615 real(RKG) , intent(in) :: lb, ub
1616 real(RKG) , intent(in) , optional :: abstol
1617 integer(IK) , intent(in) , optional :: niter
1618 integer(IK) , intent(out) , optional :: neval
1619 real(RKG) :: root
1620 end function
1621#endif
1622
1623#if RK3_ENABLED
1624 recursive module function getRootDef_RK3(getFunc, lb, ub, abstol, neval, niter) result(root)
1625#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1626 !DEC$ ATTRIBUTES DLLEXPORT :: getRootDef_RK3
1627#endif
1628 use pm_kind, only: RKG => RK3
1629 procedure(real(RKG)) :: getFunc
1630 real(RKG) , intent(in) :: lb, ub
1631 real(RKG) , intent(in) , optional :: abstol
1632 integer(IK) , intent(in) , optional :: niter
1633 integer(IK) , intent(out) , optional :: neval
1634 real(RKG) :: root
1635 end function
1636#endif
1637
1638#if RK2_ENABLED
1639 recursive module function getRootDef_RK2(getFunc, lb, ub, abstol, neval, niter) result(root)
1640#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1641 !DEC$ ATTRIBUTES DLLEXPORT :: getRootDef_RK2
1642#endif
1643 use pm_kind, only: RKG => RK2
1644 procedure(real(RKG)) :: getFunc
1645 real(RKG) , intent(in) :: lb, ub
1646 real(RKG) , intent(in) , optional :: abstol
1647 integer(IK) , intent(in) , optional :: niter
1648 integer(IK) , intent(out) , optional :: neval
1649 real(RKG) :: root
1650 end function
1651#endif
1652
1653#if RK1_ENABLED
1654 recursive module function getRootDef_RK1(getFunc, lb, ub, abstol, neval, niter) result(root)
1655#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1656 !DEC$ ATTRIBUTES DLLEXPORT :: getRootDef_RK1
1657#endif
1658 use pm_kind, only: RKG => RK1
1659 procedure(real(RKG)) :: getFunc
1660 real(RKG) , intent(in) :: lb, ub
1661 real(RKG) , intent(in) , optional :: abstol
1662 integer(IK) , intent(in) , optional :: niter
1663 integer(IK) , intent(out) , optional :: neval
1664 real(RKG) :: root
1665 end function
1666#endif
1667
1668 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1669
1670 end interface
1671
1672 ! False
1673
1674 interface getRoot
1675
1676 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1677
1678#if RK5_ENABLED
1679 recursive module function getRootFalse_RK5(method, getFunc, lb, ub, abstol, neval, niter) result(root)
1680#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1681 !DEC$ ATTRIBUTES DLLEXPORT :: getRootFalse_RK5
1682#endif
1683 use pm_kind, only: RKG => RK5
1684 procedure(real(RKG)) :: getFunc
1685 type(false_type) , intent(in) :: method
1686 real(RKG) , intent(in) :: lb, ub
1687 real(RKG) , intent(in) , optional :: abstol
1688 integer(IK) , intent(in) , optional :: niter
1689 integer(IK) , intent(out) , optional :: neval
1690 real(RKG) :: root
1691 end function
1692#endif
1693
1694#if RK4_ENABLED
1695 recursive module function getRootFalse_RK4(method, getFunc, lb, ub, abstol, neval, niter) result(root)
1696#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1697 !DEC$ ATTRIBUTES DLLEXPORT :: getRootFalse_RK4
1698#endif
1699 use pm_kind, only: RKG => RK4
1700 procedure(real(RKG)) :: getFunc
1701 type(false_type) , intent(in) :: method
1702 real(RKG) , intent(in) :: lb, ub
1703 real(RKG) , intent(in) , optional :: abstol
1704 integer(IK) , intent(in) , optional :: niter
1705 integer(IK) , intent(out) , optional :: neval
1706 real(RKG) :: root
1707 end function
1708#endif
1709
1710#if RK3_ENABLED
1711 recursive module function getRootFalse_RK3(method, getFunc, lb, ub, abstol, neval, niter) result(root)
1712#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1713 !DEC$ ATTRIBUTES DLLEXPORT :: getRootFalse_RK3
1714#endif
1715 use pm_kind, only: RKG => RK3
1716 procedure(real(RKG)) :: getFunc
1717 type(false_type) , intent(in) :: method
1718 real(RKG) , intent(in) :: lb, ub
1719 real(RKG) , intent(in) , optional :: abstol
1720 integer(IK) , intent(in) , optional :: niter
1721 integer(IK) , intent(out) , optional :: neval
1722 real(RKG) :: root
1723 end function
1724#endif
1725
1726#if RK2_ENABLED
1727 recursive module function getRootFalse_RK2(method, getFunc, lb, ub, abstol, neval, niter) result(root)
1728#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1729 !DEC$ ATTRIBUTES DLLEXPORT :: getRootFalse_RK2
1730#endif
1731 use pm_kind, only: RKG => RK2
1732 procedure(real(RKG)) :: getFunc
1733 type(false_type) , intent(in) :: method
1734 real(RKG) , intent(in) :: lb, ub
1735 real(RKG) , intent(in) , optional :: abstol
1736 integer(IK) , intent(in) , optional :: niter
1737 integer(IK) , intent(out) , optional :: neval
1738 real(RKG) :: root
1739 end function
1740#endif
1741
1742#if RK1_ENABLED
1743 recursive module function getRootFalse_RK1(method, getFunc, lb, ub, abstol, neval, niter) result(root)
1744#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1745 !DEC$ ATTRIBUTES DLLEXPORT :: getRootFalse_RK1
1746#endif
1747 use pm_kind, only: RKG => RK1
1748 procedure(real(RKG)) :: getFunc
1749 type(false_type) , intent(in) :: method
1750 real(RKG) , intent(in) :: lb, ub
1751 real(RKG) , intent(in) , optional :: abstol
1752 integer(IK) , intent(in) , optional :: niter
1753 integer(IK) , intent(out) , optional :: neval
1754 real(RKG) :: root
1755 end function
1756#endif
1757
1758 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1759
1760 end interface
1761
1762 ! Bisection
1763
1764 interface getRoot
1765
1766 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1767
1768#if RK5_ENABLED
1769 recursive module function getRootBisection_RK5(method, getFunc, lb, ub, abstol, neval, niter) result(root)
1770#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1771 !DEC$ ATTRIBUTES DLLEXPORT :: getRootBisection_RK5
1772#endif
1773 use pm_kind, only: RKG => RK5
1774 procedure(real(RKG)) :: getFunc
1775 type(bisection_type) , intent(in) :: method
1776 real(RKG) , intent(in) :: lb, ub
1777 real(RKG) , intent(in) , optional :: abstol
1778 integer(IK) , intent(in) , optional :: niter
1779 integer(IK) , intent(out) , optional :: neval
1780 real(RKG) :: root
1781 end function
1782#endif
1783
1784#if RK4_ENABLED
1785 recursive module function getRootBisection_RK4(method, getFunc, lb, ub, abstol, neval, niter) result(root)
1786#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1787 !DEC$ ATTRIBUTES DLLEXPORT :: getRootBisection_RK4
1788#endif
1789 use pm_kind, only: RKG => RK4
1790 procedure(real(RKG)) :: getFunc
1791 type(bisection_type) , intent(in) :: method
1792 real(RKG) , intent(in) :: lb, ub
1793 real(RKG) , intent(in) , optional :: abstol
1794 integer(IK) , intent(in) , optional :: niter
1795 integer(IK) , intent(out) , optional :: neval
1796 real(RKG) :: root
1797 end function
1798#endif
1799
1800#if RK3_ENABLED
1801 recursive module function getRootBisection_RK3(method, getFunc, lb, ub, abstol, neval, niter) result(root)
1802#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1803 !DEC$ ATTRIBUTES DLLEXPORT :: getRootBisection_RK3
1804#endif
1805 use pm_kind, only: RKG => RK3
1806 procedure(real(RKG)) :: getFunc
1807 type(bisection_type) , intent(in) :: method
1808 real(RKG) , intent(in) :: lb, ub
1809 real(RKG) , intent(in) , optional :: abstol
1810 integer(IK) , intent(in) , optional :: niter
1811 integer(IK) , intent(out) , optional :: neval
1812 real(RKG) :: root
1813 end function
1814#endif
1815
1816#if RK2_ENABLED
1817 recursive module function getRootBisection_RK2(method, getFunc, lb, ub, abstol, neval, niter) result(root)
1818#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1819 !DEC$ ATTRIBUTES DLLEXPORT :: getRootBisection_RK2
1820#endif
1821 use pm_kind, only: RKG => RK2
1822 procedure(real(RKG)) :: getFunc
1823 type(bisection_type) , intent(in) :: method
1824 real(RKG) , intent(in) :: lb, ub
1825 real(RKG) , intent(in) , optional :: abstol
1826 integer(IK) , intent(in) , optional :: niter
1827 integer(IK) , intent(out) , optional :: neval
1828 real(RKG) :: root
1829 end function
1830#endif
1831
1832#if RK1_ENABLED
1833 recursive module function getRootBisection_RK1(method, getFunc, lb, ub, abstol, neval, niter) result(root)
1834#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1835 !DEC$ ATTRIBUTES DLLEXPORT :: getRootBisection_RK1
1836#endif
1837 use pm_kind, only: RKG => RK1
1838 procedure(real(RKG)) :: getFunc
1839 type(bisection_type) , intent(in) :: method
1840 real(RKG) , intent(in) :: lb, ub
1841 real(RKG) , intent(in) , optional :: abstol
1842 integer(IK) , intent(in) , optional :: niter
1843 integer(IK) , intent(out) , optional :: neval
1844 real(RKG) :: root
1845 end function
1846#endif
1847
1848 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1849
1850 end interface
1851
1852 ! Secant
1853
1854 interface getRoot
1855
1856 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1857
1858#if RK5_ENABLED
1859 recursive module function getRootSecant_RK5(method, getFunc, lb, ub, abstol, neval, niter) result(root)
1860#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1861 !DEC$ ATTRIBUTES DLLEXPORT :: getRootSecant_RK5
1862#endif
1863 use pm_kind, only: RKG => RK5
1864 procedure(real(RKG)) :: getFunc
1865 type(secant_type) , intent(in) :: method
1866 real(RKG) , intent(in) :: lb, ub
1867 real(RKG) , intent(in) , optional :: abstol
1868 integer(IK) , intent(in) , optional :: niter
1869 integer(IK) , intent(out) , optional :: neval
1870 real(RKG) :: root
1871 end function
1872#endif
1873
1874#if RK4_ENABLED
1875 recursive module function getRootSecant_RK4(method, getFunc, lb, ub, abstol, neval, niter) result(root)
1876#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1877 !DEC$ ATTRIBUTES DLLEXPORT :: getRootSecant_RK4
1878#endif
1879 use pm_kind, only: RKG => RK4
1880 procedure(real(RKG)) :: getFunc
1881 type(secant_type) , intent(in) :: method
1882 real(RKG) , intent(in) :: lb, ub
1883 real(RKG) , intent(in) , optional :: abstol
1884 integer(IK) , intent(in) , optional :: niter
1885 integer(IK) , intent(out) , optional :: neval
1886 real(RKG) :: root
1887 end function
1888#endif
1889
1890#if RK3_ENABLED
1891 recursive module function getRootSecant_RK3(method, getFunc, lb, ub, abstol, neval, niter) result(root)
1892#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1893 !DEC$ ATTRIBUTES DLLEXPORT :: getRootSecant_RK3
1894#endif
1895 use pm_kind, only: RKG => RK3
1896 procedure(real(RKG)) :: getFunc
1897 type(secant_type) , intent(in) :: method
1898 real(RKG) , intent(in) :: lb, ub
1899 real(RKG) , intent(in) , optional :: abstol
1900 integer(IK) , intent(in) , optional :: niter
1901 integer(IK) , intent(out) , optional :: neval
1902 real(RKG) :: root
1903 end function
1904#endif
1905
1906#if RK2_ENABLED
1907 recursive module function getRootSecant_RK2(method, getFunc, lb, ub, abstol, neval, niter) result(root)
1908#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1909 !DEC$ ATTRIBUTES DLLEXPORT :: getRootSecant_RK2
1910#endif
1911 use pm_kind, only: RKG => RK2
1912 procedure(real(RKG)) :: getFunc
1913 type(secant_type) , intent(in) :: method
1914 real(RKG) , intent(in) :: lb, ub
1915 real(RKG) , intent(in) , optional :: abstol
1916 integer(IK) , intent(in) , optional :: niter
1917 integer(IK) , intent(out) , optional :: neval
1918 real(RKG) :: root
1919 end function
1920#endif
1921
1922#if RK1_ENABLED
1923 recursive module function getRootSecant_RK1(method, getFunc, lb, ub, abstol, neval, niter) result(root)
1924#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1925 !DEC$ ATTRIBUTES DLLEXPORT :: getRootSecant_RK1
1926#endif
1927 use pm_kind, only: RKG => RK1
1928 procedure(real(RKG)) :: getFunc
1929 type(secant_type) , intent(in) :: method
1930 real(RKG) , intent(in) :: lb, ub
1931 real(RKG) , intent(in) , optional :: abstol
1932 integer(IK) , intent(in) , optional :: niter
1933 integer(IK) , intent(out) , optional :: neval
1934 real(RKG) :: root
1935 end function
1936#endif
1937
1938 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1939
1940 end interface
1941
1942 ! Brent
1943
1944 interface getRoot
1945
1946 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1947
1948#if RK5_ENABLED
1949 recursive module function getRootBrent_RK5(method, getFunc, lb, ub, abstol, neval, niter) result(root)
1950#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1951 !DEC$ ATTRIBUTES DLLEXPORT :: getRootBrent_RK5
1952#endif
1953 use pm_kind, only: RKG => RK5
1954 procedure(real(RKG)) :: getFunc
1955 type(brent_type) , intent(in) :: method
1956 real(RKG) , intent(in) :: lb, ub
1957 real(RKG) , intent(in) , optional :: abstol
1958 integer(IK) , intent(in) , optional :: niter
1959 integer(IK) , intent(out) , optional :: neval
1960 real(RKG) :: root
1961 end function
1962#endif
1963
1964#if RK4_ENABLED
1965 recursive module function getRootBrent_RK4(method, getFunc, lb, ub, abstol, neval, niter) result(root)
1966#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1967 !DEC$ ATTRIBUTES DLLEXPORT :: getRootBrent_RK4
1968#endif
1969 use pm_kind, only: RKG => RK4
1970 procedure(real(RKG)) :: getFunc
1971 type(brent_type) , intent(in) :: method
1972 real(RKG) , intent(in) :: lb, ub
1973 real(RKG) , intent(in) , optional :: abstol
1974 integer(IK) , intent(in) , optional :: niter
1975 integer(IK) , intent(out) , optional :: neval
1976 real(RKG) :: root
1977 end function
1978#endif
1979
1980#if RK3_ENABLED
1981 recursive module function getRootBrent_RK3(method, getFunc, lb, ub, abstol, neval, niter) result(root)
1982#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1983 !DEC$ ATTRIBUTES DLLEXPORT :: getRootBrent_RK3
1984#endif
1985 use pm_kind, only: RKG => RK3
1986 procedure(real(RKG)) :: getFunc
1987 type(brent_type) , intent(in) :: method
1988 real(RKG) , intent(in) :: lb, ub
1989 real(RKG) , intent(in) , optional :: abstol
1990 integer(IK) , intent(in) , optional :: niter
1991 integer(IK) , intent(out) , optional :: neval
1992 real(RKG) :: root
1993 end function
1994#endif
1995
1996#if RK2_ENABLED
1997 recursive module function getRootBrent_RK2(method, getFunc, lb, ub, abstol, neval, niter) result(root)
1998#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1999 !DEC$ ATTRIBUTES DLLEXPORT :: getRootBrent_RK2
2000#endif
2001 use pm_kind, only: RKG => RK2
2002 procedure(real(RKG)) :: getFunc
2003 type(brent_type) , intent(in) :: method
2004 real(RKG) , intent(in) :: lb, ub
2005 real(RKG) , intent(in) , optional :: abstol
2006 integer(IK) , intent(in) , optional :: niter
2007 integer(IK) , intent(out) , optional :: neval
2008 real(RKG) :: root
2009 end function
2010#endif
2011
2012#if RK1_ENABLED
2013 recursive module function getRootBrent_RK1(method, getFunc, lb, ub, abstol, neval, niter) result(root)
2014#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2015 !DEC$ ATTRIBUTES DLLEXPORT :: getRootBrent_RK1
2016#endif
2017 use pm_kind, only: RKG => RK1
2018 procedure(real(RKG)) :: getFunc
2019 type(brent_type) , intent(in) :: method
2020 real(RKG) , intent(in) :: lb, ub
2021 real(RKG) , intent(in) , optional :: abstol
2022 integer(IK) , intent(in) , optional :: niter
2023 integer(IK) , intent(out) , optional :: neval
2024 real(RKG) :: root
2025 end function
2026#endif
2027
2028 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2029
2030 end interface
2031
2032 ! Ridders
2033
2034 interface getRoot
2035
2036 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2037
2038#if RK5_ENABLED
2039 recursive module function getRootRidders_RK5(method, getFunc, lb, ub, abstol, neval, niter) result(root)
2040#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2041 !DEC$ ATTRIBUTES DLLEXPORT :: getRootRidders_RK5
2042#endif
2043 use pm_kind, only: RKG => RK5
2044 procedure(real(RKG)) :: getFunc
2045 type(ridders_type) , intent(in) :: method
2046 real(RKG) , intent(in) :: lb, ub
2047 real(RKG) , intent(in) , optional :: abstol
2048 integer(IK) , intent(in) , optional :: niter
2049 integer(IK) , intent(out) , optional :: neval
2050 real(RKG) :: root
2051 end function
2052#endif
2053
2054#if RK4_ENABLED
2055 recursive module function getRootRidders_RK4(method, getFunc, lb, ub, abstol, neval, niter) result(root)
2056#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2057 !DEC$ ATTRIBUTES DLLEXPORT :: getRootRidders_RK4
2058#endif
2059 use pm_kind, only: RKG => RK4
2060 procedure(real(RKG)) :: getFunc
2061 type(ridders_type) , intent(in) :: method
2062 real(RKG) , intent(in) :: lb, ub
2063 real(RKG) , intent(in) , optional :: abstol
2064 integer(IK) , intent(in) , optional :: niter
2065 integer(IK) , intent(out) , optional :: neval
2066 real(RKG) :: root
2067 end function
2068#endif
2069
2070#if RK3_ENABLED
2071 recursive module function getRootRidders_RK3(method, getFunc, lb, ub, abstol, neval, niter) result(root)
2072#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2073 !DEC$ ATTRIBUTES DLLEXPORT :: getRootRidders_RK3
2074#endif
2075 use pm_kind, only: RKG => RK3
2076 procedure(real(RKG)) :: getFunc
2077 type(ridders_type) , intent(in) :: method
2078 real(RKG) , intent(in) :: lb, ub
2079 real(RKG) , intent(in) , optional :: abstol
2080 integer(IK) , intent(in) , optional :: niter
2081 integer(IK) , intent(out) , optional :: neval
2082 real(RKG) :: root
2083 end function
2084#endif
2085
2086#if RK2_ENABLED
2087 recursive module function getRootRidders_RK2(method, getFunc, lb, ub, abstol, neval, niter) result(root)
2088#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2089 !DEC$ ATTRIBUTES DLLEXPORT :: getRootRidders_RK2
2090#endif
2091 use pm_kind, only: RKG => RK2
2092 procedure(real(RKG)) :: getFunc
2093 type(ridders_type) , intent(in) :: method
2094 real(RKG) , intent(in) :: lb, ub
2095 real(RKG) , intent(in) , optional :: abstol
2096 integer(IK) , intent(in) , optional :: niter
2097 integer(IK) , intent(out) , optional :: neval
2098 real(RKG) :: root
2099 end function
2100#endif
2101
2102#if RK1_ENABLED
2103 recursive module function getRootRidders_RK1(method, getFunc, lb, ub, abstol, neval, niter) result(root)
2104#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2105 !DEC$ ATTRIBUTES DLLEXPORT :: getRootRidders_RK1
2106#endif
2107 use pm_kind, only: RKG => RK1
2108 procedure(real(RKG)) :: getFunc
2109 type(ridders_type) , intent(in) :: method
2110 real(RKG) , intent(in) :: lb, ub
2111 real(RKG) , intent(in) , optional :: abstol
2112 integer(IK) , intent(in) , optional :: niter
2113 integer(IK) , intent(out) , optional :: neval
2114 real(RKG) :: root
2115 end function
2116#endif
2117
2118 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2119
2120 end interface
2121
2122 ! TOMS748
2123
2124 interface getRoot
2125
2126 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2127
2128#if RK5_ENABLED
2129 recursive module function getRootTOMS748_RK5(method, getFunc, lb, ub, abstol, neval, niter) result(root)
2130#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2131 !DEC$ ATTRIBUTES DLLEXPORT :: getRootTOMS748_RK5
2132#endif
2133 use pm_kind, only: RKG => RK5
2134 procedure(real(RKG)) :: getFunc
2135 type(toms748_type) , intent(in) :: method
2136 real(RKG) , intent(in) :: lb, ub
2137 real(RKG) , intent(in) , optional :: abstol
2138 integer(IK) , intent(in) , optional :: niter
2139 integer(IK) , intent(out) , optional :: neval
2140 real(RKG) :: root
2141 end function
2142#endif
2143
2144#if RK4_ENABLED
2145 recursive module function getRootTOMS748_RK4(method, getFunc, lb, ub, abstol, neval, niter) result(root)
2146#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2147 !DEC$ ATTRIBUTES DLLEXPORT :: getRootTOMS748_RK4
2148#endif
2149 use pm_kind, only: RKG => RK4
2150 procedure(real(RKG)) :: getFunc
2151 type(toms748_type) , intent(in) :: method
2152 real(RKG) , intent(in) :: lb, ub
2153 real(RKG) , intent(in) , optional :: abstol
2154 integer(IK) , intent(in) , optional :: niter
2155 integer(IK) , intent(out) , optional :: neval
2156 real(RKG) :: root
2157 end function
2158#endif
2159
2160#if RK3_ENABLED
2161 recursive module function getRootTOMS748_RK3(method, getFunc, lb, ub, abstol, neval, niter) result(root)
2162#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2163 !DEC$ ATTRIBUTES DLLEXPORT :: getRootTOMS748_RK3
2164#endif
2165 use pm_kind, only: RKG => RK3
2166 procedure(real(RKG)) :: getFunc
2167 type(toms748_type) , intent(in) :: method
2168 real(RKG) , intent(in) :: lb, ub
2169 real(RKG) , intent(in) , optional :: abstol
2170 integer(IK) , intent(in) , optional :: niter
2171 integer(IK) , intent(out) , optional :: neval
2172 real(RKG) :: root
2173 end function
2174#endif
2175
2176#if RK2_ENABLED
2177 recursive module function getRootTOMS748_RK2(method, getFunc, lb, ub, abstol, neval, niter) result(root)
2178#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2179 !DEC$ ATTRIBUTES DLLEXPORT :: getRootTOMS748_RK2
2180#endif
2181 use pm_kind, only: RKG => RK2
2182 procedure(real(RKG)) :: getFunc
2183 type(toms748_type) , intent(in) :: method
2184 real(RKG) , intent(in) :: lb, ub
2185 real(RKG) , intent(in) , optional :: abstol
2186 integer(IK) , intent(in) , optional :: niter
2187 integer(IK) , intent(out) , optional :: neval
2188 real(RKG) :: root
2189 end function
2190#endif
2191
2192#if RK1_ENABLED
2193 recursive module function getRootTOMS748_RK1(method, getFunc, lb, ub, abstol, neval, niter) result(root)
2194#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2195 !DEC$ ATTRIBUTES DLLEXPORT :: getRootTOMS748_RK1
2196#endif
2197 use pm_kind, only: RKG => RK1
2198 procedure(real(RKG)) :: getFunc
2199 type(toms748_type) , intent(in) :: method
2200 real(RKG) , intent(in) :: lb, ub
2201 real(RKG) , intent(in) , optional :: abstol
2202 integer(IK) , intent(in) , optional :: niter
2203 integer(IK) , intent(out) , optional :: neval
2204 real(RKG) :: root
2205 end function
2206#endif
2207
2208 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2209
2210 end interface
2211
2212 ! Newton
2213
2214 interface getRoot
2215
2216 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2217
2218#if RK5_ENABLED
2219 recursive module function getRootNewton_RK5(method, getFunc, lb, ub, abstol, neval, niter, init) result(root)
2220#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2221 !DEC$ ATTRIBUTES DLLEXPORT :: getRootNewton_RK5
2222#endif
2223 use pm_kind, only: RKG => RK5
2224 procedure(real(RKG)) :: getFunc
2225 type(newton_type) , intent(in) :: method
2226 real(RKG) , intent(in) :: lb, ub
2227 real(RKG) , intent(in) , optional :: abstol
2228 integer(IK) , intent(in) , optional :: niter
2229 integer(IK) , intent(out) , optional :: neval
2230 real(RKG) , intent(in) , optional :: init
2231 real(RKG) :: root
2232 end function
2233#endif
2234
2235#if RK4_ENABLED
2236 recursive module function getRootNewton_RK4(method, getFunc, lb, ub, abstol, neval, niter, init) result(root)
2237#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2238 !DEC$ ATTRIBUTES DLLEXPORT :: getRootNewton_RK4
2239#endif
2240 use pm_kind, only: RKG => RK4
2241 procedure(real(RKG)) :: getFunc
2242 type(newton_type) , intent(in) :: method
2243 real(RKG) , intent(in) :: lb, ub
2244 real(RKG) , intent(in) , optional :: abstol
2245 integer(IK) , intent(in) , optional :: niter
2246 integer(IK) , intent(out) , optional :: neval
2247 real(RKG) , intent(in) , optional :: init
2248 real(RKG) :: root
2249 end function
2250#endif
2251
2252#if RK3_ENABLED
2253 recursive module function getRootNewton_RK3(method, getFunc, lb, ub, abstol, neval, niter, init) result(root)
2254#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2255 !DEC$ ATTRIBUTES DLLEXPORT :: getRootNewton_RK3
2256#endif
2257 use pm_kind, only: RKG => RK3
2258 procedure(real(RKG)) :: getFunc
2259 type(newton_type) , intent(in) :: method
2260 real(RKG) , intent(in) :: lb, ub
2261 real(RKG) , intent(in) , optional :: abstol
2262 integer(IK) , intent(in) , optional :: niter
2263 integer(IK) , intent(out) , optional :: neval
2264 real(RKG) , intent(in) , optional :: init
2265 real(RKG) :: root
2266 end function
2267#endif
2268
2269#if RK2_ENABLED
2270 recursive module function getRootNewton_RK2(method, getFunc, lb, ub, abstol, neval, niter, init) result(root)
2271#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2272 !DEC$ ATTRIBUTES DLLEXPORT :: getRootNewton_RK2
2273#endif
2274 use pm_kind, only: RKG => RK2
2275 procedure(real(RKG)) :: getFunc
2276 type(newton_type) , intent(in) :: method
2277 real(RKG) , intent(in) :: lb, ub
2278 real(RKG) , intent(in) , optional :: abstol
2279 integer(IK) , intent(in) , optional :: niter
2280 integer(IK) , intent(out) , optional :: neval
2281 real(RKG) , intent(in) , optional :: init
2282 real(RKG) :: root
2283 end function
2284#endif
2285
2286#if RK1_ENABLED
2287 recursive module function getRootNewton_RK1(method, getFunc, lb, ub, abstol, neval, niter, init) result(root)
2288#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2289 !DEC$ ATTRIBUTES DLLEXPORT :: getRootNewton_RK1
2290#endif
2291 use pm_kind, only: RKG => RK1
2292 procedure(real(RKG)) :: getFunc
2293 type(newton_type) , intent(in) :: method
2294 real(RKG) , intent(in) :: lb, ub
2295 real(RKG) , intent(in) , optional :: abstol
2296 integer(IK) , intent(in) , optional :: niter
2297 integer(IK) , intent(out) , optional :: neval
2298 real(RKG) , intent(in) , optional :: init
2299 real(RKG) :: root
2300 end function
2301#endif
2302
2303 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2304
2305 end interface
2306
2307 ! Halley
2308
2309 interface getRoot
2310
2311 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2312
2313#if RK5_ENABLED
2314 recursive module function getRootHalley_RK5(method, getFunc, lb, ub, abstol, neval, niter, init) result(root)
2315#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2316 !DEC$ ATTRIBUTES DLLEXPORT :: getRootHalley_RK5
2317#endif
2318 use pm_kind, only: RKG => RK5
2319 procedure(real(RKG)) :: getFunc
2320 type(halley_type) , intent(in) :: method
2321 real(RKG) , intent(in) :: lb, ub
2322 real(RKG) , intent(in) , optional :: abstol
2323 integer(IK) , intent(in) , optional :: niter
2324 integer(IK) , intent(out) , optional :: neval
2325 real(RKG) , intent(in) , optional :: init
2326 real(RKG) :: root
2327 end function
2328#endif
2329
2330#if RK4_ENABLED
2331 recursive module function getRootHalley_RK4(method, getFunc, lb, ub, abstol, neval, niter, init) result(root)
2332#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2333 !DEC$ ATTRIBUTES DLLEXPORT :: getRootHalley_RK4
2334#endif
2335 use pm_kind, only: RKG => RK4
2336 procedure(real(RKG)) :: getFunc
2337 type(halley_type) , intent(in) :: method
2338 real(RKG) , intent(in) :: lb, ub
2339 real(RKG) , intent(in) , optional :: abstol
2340 integer(IK) , intent(in) , optional :: niter
2341 integer(IK) , intent(out) , optional :: neval
2342 real(RKG) , intent(in) , optional :: init
2343 real(RKG) :: root
2344 end function
2345#endif
2346
2347#if RK3_ENABLED
2348 recursive module function getRootHalley_RK3(method, getFunc, lb, ub, abstol, neval, niter, init) result(root)
2349#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2350 !DEC$ ATTRIBUTES DLLEXPORT :: getRootHalley_RK3
2351#endif
2352 use pm_kind, only: RKG => RK3
2353 procedure(real(RKG)) :: getFunc
2354 type(halley_type) , intent(in) :: method
2355 real(RKG) , intent(in) :: lb, ub
2356 real(RKG) , intent(in) , optional :: abstol
2357 integer(IK) , intent(in) , optional :: niter
2358 integer(IK) , intent(out) , optional :: neval
2359 real(RKG) , intent(in) , optional :: init
2360 real(RKG) :: root
2361 end function
2362#endif
2363
2364#if RK2_ENABLED
2365 recursive module function getRootHalley_RK2(method, getFunc, lb, ub, abstol, neval, niter, init) result(root)
2366#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2367 !DEC$ ATTRIBUTES DLLEXPORT :: getRootHalley_RK2
2368#endif
2369 use pm_kind, only: RKG => RK2
2370 procedure(real(RKG)) :: getFunc
2371 type(halley_type) , intent(in) :: method
2372 real(RKG) , intent(in) :: lb, ub
2373 real(RKG) , intent(in) , optional :: abstol
2374 integer(IK) , intent(in) , optional :: niter
2375 integer(IK) , intent(out) , optional :: neval
2376 real(RKG) , intent(in) , optional :: init
2377 real(RKG) :: root
2378 end function
2379#endif
2380
2381#if RK1_ENABLED
2382 recursive module function getRootHalley_RK1(method, getFunc, lb, ub, abstol, neval, niter, init) result(root)
2383#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2384 !DEC$ ATTRIBUTES DLLEXPORT :: getRootHalley_RK1
2385#endif
2386 use pm_kind, only: RKG => RK1
2387 procedure(real(RKG)) :: getFunc
2388 type(halley_type) , intent(in) :: method
2389 real(RKG) , intent(in) :: lb, ub
2390 real(RKG) , intent(in) , optional :: abstol
2391 integer(IK) , intent(in) , optional :: niter
2392 integer(IK) , intent(out) , optional :: neval
2393 real(RKG) , intent(in) , optional :: init
2394 real(RKG) :: root
2395 end function
2396#endif
2397
2398 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2399
2400 end interface
2401
2402 ! Schroder
2403
2404 interface getRoot
2405
2406 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2407
2408#if RK5_ENABLED
2409 recursive module function getRootSchroder_RK5(method, getFunc, lb, ub, abstol, neval, niter, init) result(root)
2410#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2411 !DEC$ ATTRIBUTES DLLEXPORT :: getRootSchroder_RK5
2412#endif
2413 use pm_kind, only: RKG => RK5
2414 procedure(real(RKG)) :: getFunc
2415 type(schroder_type) , intent(in) :: method
2416 real(RKG) , intent(in) :: lb, ub
2417 real(RKG) , intent(in) , optional :: abstol
2418 integer(IK) , intent(in) , optional :: niter
2419 integer(IK) , intent(out) , optional :: neval
2420 real(RKG) , intent(in) , optional :: init
2421 real(RKG) :: root
2422 end function
2423#endif
2424
2425#if RK4_ENABLED
2426 recursive module function getRootSchroder_RK4(method, getFunc, lb, ub, abstol, neval, niter, init) result(root)
2427#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2428 !DEC$ ATTRIBUTES DLLEXPORT :: getRootSchroder_RK4
2429#endif
2430 use pm_kind, only: RKG => RK4
2431 procedure(real(RKG)) :: getFunc
2432 type(schroder_type) , intent(in) :: method
2433 real(RKG) , intent(in) :: lb, ub
2434 real(RKG) , intent(in) , optional :: abstol
2435 integer(IK) , intent(in) , optional :: niter
2436 integer(IK) , intent(out) , optional :: neval
2437 real(RKG) , intent(in) , optional :: init
2438 real(RKG) :: root
2439 end function
2440#endif
2441
2442#if RK3_ENABLED
2443 recursive module function getRootSchroder_RK3(method, getFunc, lb, ub, abstol, neval, niter, init) result(root)
2444#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2445 !DEC$ ATTRIBUTES DLLEXPORT :: getRootSchroder_RK3
2446#endif
2447 use pm_kind, only: RKG => RK3
2448 procedure(real(RKG)) :: getFunc
2449 type(schroder_type) , intent(in) :: method
2450 real(RKG) , intent(in) :: lb, ub
2451 real(RKG) , intent(in) , optional :: abstol
2452 integer(IK) , intent(in) , optional :: niter
2453 integer(IK) , intent(out) , optional :: neval
2454 real(RKG) , intent(in) , optional :: init
2455 real(RKG) :: root
2456 end function
2457#endif
2458
2459#if RK2_ENABLED
2460 recursive module function getRootSchroder_RK2(method, getFunc, lb, ub, abstol, neval, niter, init) result(root)
2461#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2462 !DEC$ ATTRIBUTES DLLEXPORT :: getRootSchroder_RK2
2463#endif
2464 use pm_kind, only: RKG => RK2
2465 procedure(real(RKG)) :: getFunc
2466 type(schroder_type) , intent(in) :: method
2467 real(RKG) , intent(in) :: lb, ub
2468 real(RKG) , intent(in) , optional :: abstol
2469 integer(IK) , intent(in) , optional :: niter
2470 integer(IK) , intent(out) , optional :: neval
2471 real(RKG) , intent(in) , optional :: init
2472 real(RKG) :: root
2473 end function
2474#endif
2475
2476#if RK1_ENABLED
2477 recursive module function getRootSchroder_RK1(method, getFunc, lb, ub, abstol, neval, niter, init) result(root)
2478#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2479 !DEC$ ATTRIBUTES DLLEXPORT :: getRootSchroder_RK1
2480#endif
2481 use pm_kind, only: RKG => RK1
2482 procedure(real(RKG)) :: getFunc
2483 type(schroder_type) , intent(in) :: method
2484 real(RKG) , intent(in) :: lb, ub
2485 real(RKG) , intent(in) , optional :: abstol
2486 integer(IK) , intent(in) , optional :: niter
2487 integer(IK) , intent(out) , optional :: neval
2488 real(RKG) , intent(in) , optional :: init
2489 real(RKG) :: root
2490 end function
2491#endif
2492
2493 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2494
2495 end interface
2496
2497!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2498
2672
2673 ! False
2674
2675 interface setRoot
2676
2677 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2678
2679#if RK5_ENABLED
2680 recursive module subroutine setRootFalseFixed_RK5(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
2681#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2682 !DEC$ ATTRIBUTES DLLEXPORT :: setRootFalseFixed_RK5
2683#endif
2684 use pm_kind, only: RKG => RK5
2685 procedure(real(RKG)) :: getFunc
2686 type(false_type) , intent(in) :: method
2687 real(RKG) , intent(out) :: root
2688 real(RKG) , value :: lb, ub
2689 real(RKG) , value :: lf, uf
2690 real(RKG) , intent(in) :: abstol
2691 integer(IK) , intent(out) :: neval
2692 end subroutine
2693#endif
2694
2695#if RK4_ENABLED
2696 recursive module subroutine setRootFalseFixed_RK4(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
2697#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2698 !DEC$ ATTRIBUTES DLLEXPORT :: setRootFalseFixed_RK4
2699#endif
2700 use pm_kind, only: RKG => RK4
2701 procedure(real(RKG)) :: getFunc
2702 type(false_type) , intent(in) :: method
2703 real(RKG) , intent(out) :: root
2704 real(RKG) , value :: lb, ub
2705 real(RKG) , value :: lf, uf
2706 real(RKG) , intent(in) :: abstol
2707 integer(IK) , intent(out) :: neval
2708 end subroutine
2709#endif
2710
2711#if RK3_ENABLED
2712 recursive module subroutine setRootFalseFixed_RK3(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
2713#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2714 !DEC$ ATTRIBUTES DLLEXPORT :: setRootFalseFixed_RK3
2715#endif
2716 use pm_kind, only: RKG => RK3
2717 procedure(real(RKG)) :: getFunc
2718 type(false_type) , intent(in) :: method
2719 real(RKG) , intent(out) :: root
2720 real(RKG) , value :: lb, ub
2721 real(RKG) , value :: lf, uf
2722 real(RKG) , intent(in) :: abstol
2723 integer(IK) , intent(out) :: neval
2724 end subroutine
2725#endif
2726
2727#if RK2_ENABLED
2728 recursive module subroutine setRootFalseFixed_RK2(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
2729#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2730 !DEC$ ATTRIBUTES DLLEXPORT :: setRootFalseFixed_RK2
2731#endif
2732 use pm_kind, only: RKG => RK2
2733 procedure(real(RKG)) :: getFunc
2734 type(false_type) , intent(in) :: method
2735 real(RKG) , intent(out) :: root
2736 real(RKG) , value :: lb, ub
2737 real(RKG) , value :: lf, uf
2738 real(RKG) , intent(in) :: abstol
2739 integer(IK) , intent(out) :: neval
2740 end subroutine
2741#endif
2742
2743#if RK1_ENABLED
2744 recursive module subroutine setRootFalseFixed_RK1(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
2745#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2746 !DEC$ ATTRIBUTES DLLEXPORT :: setRootFalseFixed_RK1
2747#endif
2748 use pm_kind, only: RKG => RK1
2749 procedure(real(RKG)) :: getFunc
2750 type(false_type) , intent(in) :: method
2751 real(RKG) , intent(out) :: root
2752 real(RKG) , value :: lb, ub
2753 real(RKG) , value :: lf, uf
2754 real(RKG) , intent(in) :: abstol
2755 integer(IK) , intent(out) :: neval
2756 end subroutine
2757#endif
2758
2759 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2760
2761#if RK5_ENABLED
2762 recursive module subroutine setRootFalseNiter_RK5(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
2763#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2764 !DEC$ ATTRIBUTES DLLEXPORT :: setRootFalseNiter_RK5
2765#endif
2766 use pm_kind, only: RKG => RK5
2767 procedure(real(RKG)) :: getFunc
2768 type(false_type) , intent(in) :: method
2769 real(RKG) , intent(out) :: root
2770 real(RKG) , value :: lb, ub
2771 real(RKG) , value :: lf, uf
2772 real(RKG) , intent(in) :: abstol
2773 integer(IK) , intent(out) :: neval
2774 integer(IK) , intent(in) :: niter
2775 end subroutine
2776#endif
2777
2778#if RK4_ENABLED
2779 recursive module subroutine setRootFalseNiter_RK4(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
2780#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2781 !DEC$ ATTRIBUTES DLLEXPORT :: setRootFalseNiter_RK4
2782#endif
2783 use pm_kind, only: RKG => RK4
2784 procedure(real(RKG)) :: getFunc
2785 type(false_type) , intent(in) :: method
2786 real(RKG) , intent(out) :: root
2787 real(RKG) , value :: lb, ub
2788 real(RKG) , value :: lf, uf
2789 real(RKG) , intent(in) :: abstol
2790 integer(IK) , intent(out) :: neval
2791 integer(IK) , intent(in) :: niter
2792 end subroutine
2793#endif
2794
2795#if RK3_ENABLED
2796 recursive module subroutine setRootFalseNiter_RK3(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
2797#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2798 !DEC$ ATTRIBUTES DLLEXPORT :: setRootFalseNiter_RK3
2799#endif
2800 use pm_kind, only: RKG => RK3
2801 procedure(real(RKG)) :: getFunc
2802 type(false_type) , intent(in) :: method
2803 real(RKG) , intent(out) :: root
2804 real(RKG) , value :: lb, ub
2805 real(RKG) , value :: lf, uf
2806 real(RKG) , intent(in) :: abstol
2807 integer(IK) , intent(out) :: neval
2808 integer(IK) , intent(in) :: niter
2809 end subroutine
2810#endif
2811
2812#if RK2_ENABLED
2813 recursive module subroutine setRootFalseNiter_RK2(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
2814#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2815 !DEC$ ATTRIBUTES DLLEXPORT :: setRootFalseNiter_RK2
2816#endif
2817 use pm_kind, only: RKG => RK2
2818 procedure(real(RKG)) :: getFunc
2819 type(false_type) , intent(in) :: method
2820 real(RKG) , intent(out) :: root
2821 real(RKG) , value :: lb, ub
2822 real(RKG) , value :: lf, uf
2823 real(RKG) , intent(in) :: abstol
2824 integer(IK) , intent(out) :: neval
2825 integer(IK) , intent(in) :: niter
2826 end subroutine
2827#endif
2828
2829#if RK1_ENABLED
2830 recursive module subroutine setRootFalseNiter_RK1(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
2831#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2832 !DEC$ ATTRIBUTES DLLEXPORT :: setRootFalseNiter_RK1
2833#endif
2834 use pm_kind, only: RKG => RK1
2835 procedure(real(RKG)) :: getFunc
2836 type(false_type) , intent(in) :: method
2837 real(RKG) , intent(out) :: root
2838 real(RKG) , value :: lb, ub
2839 real(RKG) , value :: lf, uf
2840 real(RKG) , intent(in) :: abstol
2841 integer(IK) , intent(out) :: neval
2842 integer(IK) , intent(in) :: niter
2843 end subroutine
2844#endif
2845
2846 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2847
2848 end interface
2849
2850 ! Bisection
2851
2852 interface setRoot
2853
2854 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2855
2856#if RK5_ENABLED
2857 recursive module subroutine setRootBisectionFixed_RK5(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
2858#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2859 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBisectionFixed_RK5
2860#endif
2861 use pm_kind, only: RKG => RK5
2862 procedure(real(RKG)) :: getFunc
2863 type(bisection_type) , intent(in) :: method
2864 real(RKG) , intent(out) :: root
2865 real(RKG) , value :: lb, ub
2866 real(RKG) , value :: lf, uf
2867 real(RKG) , intent(in) :: abstol
2868 integer(IK) , intent(out) :: neval
2869 end subroutine
2870#endif
2871
2872#if RK4_ENABLED
2873 recursive module subroutine setRootBisectionFixed_RK4(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
2874#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2875 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBisectionFixed_RK4
2876#endif
2877 use pm_kind, only: RKG => RK4
2878 procedure(real(RKG)) :: getFunc
2879 type(bisection_type) , intent(in) :: method
2880 real(RKG) , intent(out) :: root
2881 real(RKG) , value :: lb, ub
2882 real(RKG) , value :: lf, uf
2883 real(RKG) , intent(in) :: abstol
2884 integer(IK) , intent(out) :: neval
2885 end subroutine
2886#endif
2887
2888#if RK3_ENABLED
2889 recursive module subroutine setRootBisectionFixed_RK3(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
2890#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2891 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBisectionFixed_RK3
2892#endif
2893 use pm_kind, only: RKG => RK3
2894 procedure(real(RKG)) :: getFunc
2895 type(bisection_type) , intent(in) :: method
2896 real(RKG) , intent(out) :: root
2897 real(RKG) , value :: lb, ub
2898 real(RKG) , value :: lf, uf
2899 real(RKG) , intent(in) :: abstol
2900 integer(IK) , intent(out) :: neval
2901 end subroutine
2902#endif
2903
2904#if RK2_ENABLED
2905 recursive module subroutine setRootBisectionFixed_RK2(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
2906#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2907 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBisectionFixed_RK2
2908#endif
2909 use pm_kind, only: RKG => RK2
2910 procedure(real(RKG)) :: getFunc
2911 type(bisection_type) , intent(in) :: method
2912 real(RKG) , intent(out) :: root
2913 real(RKG) , value :: lb, ub
2914 real(RKG) , value :: lf, uf
2915 real(RKG) , intent(in) :: abstol
2916 integer(IK) , intent(out) :: neval
2917 end subroutine
2918#endif
2919
2920#if RK1_ENABLED
2921 recursive module subroutine setRootBisectionFixed_RK1(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
2922#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2923 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBisectionFixed_RK1
2924#endif
2925 use pm_kind, only: RKG => RK1
2926 procedure(real(RKG)) :: getFunc
2927 type(bisection_type) , intent(in) :: method
2928 real(RKG) , intent(out) :: root
2929 real(RKG) , value :: lb, ub
2930 real(RKG) , value :: lf, uf
2931 real(RKG) , intent(in) :: abstol
2932 integer(IK) , intent(out) :: neval
2933 end subroutine
2934#endif
2935
2936 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2937
2938#if RK5_ENABLED
2939 recursive module subroutine setRootBisectionNiter_RK5(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
2940#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2941 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBisectionNiter_RK5
2942#endif
2943 use pm_kind, only: RKG => RK5
2944 procedure(real(RKG)) :: getFunc
2945 type(bisection_type) , intent(in) :: method
2946 real(RKG) , intent(out) :: root
2947 real(RKG) , value :: lb, ub
2948 real(RKG) , value :: lf, uf
2949 real(RKG) , intent(in) :: abstol
2950 integer(IK) , intent(out) :: neval
2951 integer(IK) , intent(in) :: niter
2952 end subroutine
2953#endif
2954
2955#if RK4_ENABLED
2956 recursive module subroutine setRootBisectionNiter_RK4(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
2957#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2958 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBisectionNiter_RK4
2959#endif
2960 use pm_kind, only: RKG => RK4
2961 procedure(real(RKG)) :: getFunc
2962 type(bisection_type) , intent(in) :: method
2963 real(RKG) , intent(out) :: root
2964 real(RKG) , value :: lb, ub
2965 real(RKG) , value :: lf, uf
2966 real(RKG) , intent(in) :: abstol
2967 integer(IK) , intent(out) :: neval
2968 integer(IK) , intent(in) :: niter
2969 end subroutine
2970#endif
2971
2972#if RK3_ENABLED
2973 recursive module subroutine setRootBisectionNiter_RK3(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
2974#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2975 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBisectionNiter_RK3
2976#endif
2977 use pm_kind, only: RKG => RK3
2978 procedure(real(RKG)) :: getFunc
2979 type(bisection_type) , intent(in) :: method
2980 real(RKG) , intent(out) :: root
2981 real(RKG) , value :: lb, ub
2982 real(RKG) , value :: lf, uf
2983 real(RKG) , intent(in) :: abstol
2984 integer(IK) , intent(out) :: neval
2985 integer(IK) , intent(in) :: niter
2986 end subroutine
2987#endif
2988
2989#if RK2_ENABLED
2990 recursive module subroutine setRootBisectionNiter_RK2(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
2991#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
2992 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBisectionNiter_RK2
2993#endif
2994 use pm_kind, only: RKG => RK2
2995 procedure(real(RKG)) :: getFunc
2996 type(bisection_type) , intent(in) :: method
2997 real(RKG) , intent(out) :: root
2998 real(RKG) , value :: lb, ub
2999 real(RKG) , value :: lf, uf
3000 real(RKG) , intent(in) :: abstol
3001 integer(IK) , intent(out) :: neval
3002 integer(IK) , intent(in) :: niter
3003 end subroutine
3004#endif
3005
3006#if RK1_ENABLED
3007 recursive module subroutine setRootBisectionNiter_RK1(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3008#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3009 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBisectionNiter_RK1
3010#endif
3011 use pm_kind, only: RKG => RK1
3012 procedure(real(RKG)) :: getFunc
3013 type(bisection_type) , intent(in) :: method
3014 real(RKG) , intent(out) :: root
3015 real(RKG) , value :: lb, ub
3016 real(RKG) , value :: lf, uf
3017 real(RKG) , intent(in) :: abstol
3018 integer(IK) , intent(out) :: neval
3019 integer(IK) , intent(in) :: niter
3020 end subroutine
3021#endif
3022
3023 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3024
3025 end interface
3026
3027 ! Secant
3028
3029 interface setRoot
3030
3031 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3032
3033#if RK5_ENABLED
3034 recursive module subroutine setRootSecantFixed_RK5(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3035#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3036 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSecantFixed_RK5
3037#endif
3038 use pm_kind, only: RKG => RK5
3039 procedure(real(RKG)) :: getFunc
3040 type(secant_type) , intent(in) :: method
3041 real(RKG) , intent(out) :: root
3042 real(RKG) , value :: lb, ub
3043 real(RKG) , value :: lf, uf
3044 real(RKG) , intent(in) :: abstol
3045 integer(IK) , intent(out) :: neval
3046 end subroutine
3047#endif
3048
3049#if RK4_ENABLED
3050 recursive module subroutine setRootSecantFixed_RK4(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3051#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3052 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSecantFixed_RK4
3053#endif
3054 use pm_kind, only: RKG => RK4
3055 procedure(real(RKG)) :: getFunc
3056 type(secant_type) , intent(in) :: method
3057 real(RKG) , intent(out) :: root
3058 real(RKG) , value :: lb, ub
3059 real(RKG) , value :: lf, uf
3060 real(RKG) , intent(in) :: abstol
3061 integer(IK) , intent(out) :: neval
3062 end subroutine
3063#endif
3064
3065#if RK3_ENABLED
3066 recursive module subroutine setRootSecantFixed_RK3(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3067#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3068 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSecantFixed_RK3
3069#endif
3070 use pm_kind, only: RKG => RK3
3071 procedure(real(RKG)) :: getFunc
3072 type(secant_type) , intent(in) :: method
3073 real(RKG) , intent(out) :: root
3074 real(RKG) , value :: lb, ub
3075 real(RKG) , value :: lf, uf
3076 real(RKG) , intent(in) :: abstol
3077 integer(IK) , intent(out) :: neval
3078 end subroutine
3079#endif
3080
3081#if RK2_ENABLED
3082 recursive module subroutine setRootSecantFixed_RK2(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3083#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3084 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSecantFixed_RK2
3085#endif
3086 use pm_kind, only: RKG => RK2
3087 procedure(real(RKG)) :: getFunc
3088 type(secant_type) , intent(in) :: method
3089 real(RKG) , intent(out) :: root
3090 real(RKG) , value :: lb, ub
3091 real(RKG) , value :: lf, uf
3092 real(RKG) , intent(in) :: abstol
3093 integer(IK) , intent(out) :: neval
3094 end subroutine
3095#endif
3096
3097#if RK1_ENABLED
3098 recursive module subroutine setRootSecantFixed_RK1(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3099#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3100 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSecantFixed_RK1
3101#endif
3102 use pm_kind, only: RKG => RK1
3103 procedure(real(RKG)) :: getFunc
3104 type(secant_type) , intent(in) :: method
3105 real(RKG) , intent(out) :: root
3106 real(RKG) , value :: lb, ub
3107 real(RKG) , value :: lf, uf
3108 real(RKG) , intent(in) :: abstol
3109 integer(IK) , intent(out) :: neval
3110 end subroutine
3111#endif
3112
3113 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3114
3115#if RK5_ENABLED
3116 recursive module subroutine setRootSecantNiter_RK5(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3117#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3118 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSecantNiter_RK5
3119#endif
3120 use pm_kind, only: RKG => RK5
3121 procedure(real(RKG)) :: getFunc
3122 type(secant_type) , intent(in) :: method
3123 real(RKG) , intent(out) :: root
3124 real(RKG) , value :: lb, ub
3125 real(RKG) , value :: lf, uf
3126 real(RKG) , intent(in) :: abstol
3127 integer(IK) , intent(out) :: neval
3128 integer(IK) , intent(in) :: niter
3129 end subroutine
3130#endif
3131
3132#if RK4_ENABLED
3133 recursive module subroutine setRootSecantNiter_RK4(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3134#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3135 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSecantNiter_RK4
3136#endif
3137 use pm_kind, only: RKG => RK4
3138 procedure(real(RKG)) :: getFunc
3139 type(secant_type) , intent(in) :: method
3140 real(RKG) , intent(out) :: root
3141 real(RKG) , value :: lb, ub
3142 real(RKG) , value :: lf, uf
3143 real(RKG) , intent(in) :: abstol
3144 integer(IK) , intent(out) :: neval
3145 integer(IK) , intent(in) :: niter
3146 end subroutine
3147#endif
3148
3149#if RK3_ENABLED
3150 recursive module subroutine setRootSecantNiter_RK3(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3151#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3152 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSecantNiter_RK3
3153#endif
3154 use pm_kind, only: RKG => RK3
3155 procedure(real(RKG)) :: getFunc
3156 type(secant_type) , intent(in) :: method
3157 real(RKG) , intent(out) :: root
3158 real(RKG) , value :: lb, ub
3159 real(RKG) , value :: lf, uf
3160 real(RKG) , intent(in) :: abstol
3161 integer(IK) , intent(out) :: neval
3162 integer(IK) , intent(in) :: niter
3163 end subroutine
3164#endif
3165
3166#if RK2_ENABLED
3167 recursive module subroutine setRootSecantNiter_RK2(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3168#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3169 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSecantNiter_RK2
3170#endif
3171 use pm_kind, only: RKG => RK2
3172 procedure(real(RKG)) :: getFunc
3173 type(secant_type) , intent(in) :: method
3174 real(RKG) , intent(out) :: root
3175 real(RKG) , value :: lb, ub
3176 real(RKG) , value :: lf, uf
3177 real(RKG) , intent(in) :: abstol
3178 integer(IK) , intent(out) :: neval
3179 integer(IK) , intent(in) :: niter
3180 end subroutine
3181#endif
3182
3183#if RK1_ENABLED
3184 recursive module subroutine setRootSecantNiter_RK1(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3185#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3186 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSecantNiter_RK1
3187#endif
3188 use pm_kind, only: RKG => RK1
3189 procedure(real(RKG)) :: getFunc
3190 type(secant_type) , intent(in) :: method
3191 real(RKG) , intent(out) :: root
3192 real(RKG) , value :: lb, ub
3193 real(RKG) , value :: lf, uf
3194 real(RKG) , intent(in) :: abstol
3195 integer(IK) , intent(out) :: neval
3196 integer(IK) , intent(in) :: niter
3197 end subroutine
3198#endif
3199
3200 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3201
3202 end interface
3203
3204 ! Brent
3205
3206 interface setRoot
3207
3208 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3209
3210#if RK5_ENABLED
3211 recursive module subroutine setRootBrentFixed_RK5(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3212#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3213 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBrentFixed_RK5
3214#endif
3215 use pm_kind, only: RKG => RK5
3216 procedure(real(RKG)) :: getFunc
3217 type(brent_type) , intent(in) :: method
3218 real(RKG) , intent(out) :: root
3219 real(RKG) , value :: lb, ub
3220 real(RKG) , value :: lf, uf
3221 real(RKG) , intent(in) :: abstol
3222 integer(IK) , intent(out) :: neval
3223 end subroutine
3224#endif
3225
3226#if RK4_ENABLED
3227 recursive module subroutine setRootBrentFixed_RK4(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3228#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3229 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBrentFixed_RK4
3230#endif
3231 use pm_kind, only: RKG => RK4
3232 procedure(real(RKG)) :: getFunc
3233 type(brent_type) , intent(in) :: method
3234 real(RKG) , intent(out) :: root
3235 real(RKG) , value :: lb, ub
3236 real(RKG) , value :: lf, uf
3237 real(RKG) , intent(in) :: abstol
3238 integer(IK) , intent(out) :: neval
3239 end subroutine
3240#endif
3241
3242#if RK3_ENABLED
3243 recursive module subroutine setRootBrentFixed_RK3(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3244#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3245 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBrentFixed_RK3
3246#endif
3247 use pm_kind, only: RKG => RK3
3248 procedure(real(RKG)) :: getFunc
3249 type(brent_type) , intent(in) :: method
3250 real(RKG) , intent(out) :: root
3251 real(RKG) , value :: lb, ub
3252 real(RKG) , value :: lf, uf
3253 real(RKG) , intent(in) :: abstol
3254 integer(IK) , intent(out) :: neval
3255 end subroutine
3256#endif
3257
3258#if RK2_ENABLED
3259 recursive module subroutine setRootBrentFixed_RK2(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3260#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3261 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBrentFixed_RK2
3262#endif
3263 use pm_kind, only: RKG => RK2
3264 procedure(real(RKG)) :: getFunc
3265 type(brent_type) , intent(in) :: method
3266 real(RKG) , intent(out) :: root
3267 real(RKG) , value :: lb, ub
3268 real(RKG) , value :: lf, uf
3269 real(RKG) , intent(in) :: abstol
3270 integer(IK) , intent(out) :: neval
3271 end subroutine
3272#endif
3273
3274#if RK1_ENABLED
3275 recursive module subroutine setRootBrentFixed_RK1(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3276#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3277 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBrentFixed_RK1
3278#endif
3279 use pm_kind, only: RKG => RK1
3280 procedure(real(RKG)) :: getFunc
3281 type(brent_type) , intent(in) :: method
3282 real(RKG) , intent(out) :: root
3283 real(RKG) , value :: lb, ub
3284 real(RKG) , value :: lf, uf
3285 real(RKG) , intent(in) :: abstol
3286 integer(IK) , intent(out) :: neval
3287 end subroutine
3288#endif
3289
3290 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3291
3292#if RK5_ENABLED
3293 recursive module subroutine setRootBrentNiter_RK5(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3294#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3295 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBrentNiter_RK5
3296#endif
3297 use pm_kind, only: RKG => RK5
3298 procedure(real(RKG)) :: getFunc
3299 type(brent_type) , intent(in) :: method
3300 real(RKG) , intent(out) :: root
3301 real(RKG) , value :: lb, ub
3302 real(RKG) , value :: lf, uf
3303 real(RKG) , intent(in) :: abstol
3304 integer(IK) , intent(out) :: neval
3305 integer(IK) , intent(in) :: niter
3306 end subroutine
3307#endif
3308
3309#if RK4_ENABLED
3310 recursive module subroutine setRootBrentNiter_RK4(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3311#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3312 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBrentNiter_RK4
3313#endif
3314 use pm_kind, only: RKG => RK4
3315 procedure(real(RKG)) :: getFunc
3316 type(brent_type) , intent(in) :: method
3317 real(RKG) , intent(out) :: root
3318 real(RKG) , value :: lb, ub
3319 real(RKG) , value :: lf, uf
3320 real(RKG) , intent(in) :: abstol
3321 integer(IK) , intent(out) :: neval
3322 integer(IK) , intent(in) :: niter
3323 end subroutine
3324#endif
3325
3326#if RK3_ENABLED
3327 recursive module subroutine setRootBrentNiter_RK3(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3328#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3329 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBrentNiter_RK3
3330#endif
3331 use pm_kind, only: RKG => RK3
3332 procedure(real(RKG)) :: getFunc
3333 type(brent_type) , intent(in) :: method
3334 real(RKG) , intent(out) :: root
3335 real(RKG) , value :: lb, ub
3336 real(RKG) , value :: lf, uf
3337 real(RKG) , intent(in) :: abstol
3338 integer(IK) , intent(out) :: neval
3339 integer(IK) , intent(in) :: niter
3340 end subroutine
3341#endif
3342
3343#if RK2_ENABLED
3344 recursive module subroutine setRootBrentNiter_RK2(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3345#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3346 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBrentNiter_RK2
3347#endif
3348 use pm_kind, only: RKG => RK2
3349 procedure(real(RKG)) :: getFunc
3350 type(brent_type) , intent(in) :: method
3351 real(RKG) , intent(out) :: root
3352 real(RKG) , value :: lb, ub
3353 real(RKG) , value :: lf, uf
3354 real(RKG) , intent(in) :: abstol
3355 integer(IK) , intent(out) :: neval
3356 integer(IK) , intent(in) :: niter
3357 end subroutine
3358#endif
3359
3360#if RK1_ENABLED
3361 recursive module subroutine setRootBrentNiter_RK1(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3362#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3363 !DEC$ ATTRIBUTES DLLEXPORT :: setRootBrentNiter_RK1
3364#endif
3365 use pm_kind, only: RKG => RK1
3366 procedure(real(RKG)) :: getFunc
3367 type(brent_type) , intent(in) :: method
3368 real(RKG) , intent(out) :: root
3369 real(RKG) , value :: lb, ub
3370 real(RKG) , value :: lf, uf
3371 real(RKG) , intent(in) :: abstol
3372 integer(IK) , intent(out) :: neval
3373 integer(IK) , intent(in) :: niter
3374 end subroutine
3375#endif
3376
3377 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3378
3379 end interface
3380
3381 ! Ridders
3382
3383 interface setRoot
3384
3385 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3386
3387#if RK5_ENABLED
3388 recursive module subroutine setRootRiddersFixed_RK5(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3389#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3390 !DEC$ ATTRIBUTES DLLEXPORT :: setRootRiddersFixed_RK5
3391#endif
3392 use pm_kind, only: RKG => RK5
3393 procedure(real(RKG)) :: getFunc
3394 type(ridders_type) , intent(in) :: method
3395 real(RKG) , intent(out) :: root
3396 real(RKG) , value :: lb, ub
3397 real(RKG) , value :: lf, uf
3398 real(RKG) , intent(in) :: abstol
3399 integer(IK) , intent(out) :: neval
3400 end subroutine
3401#endif
3402
3403#if RK4_ENABLED
3404 recursive module subroutine setRootRiddersFixed_RK4(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3405#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3406 !DEC$ ATTRIBUTES DLLEXPORT :: setRootRiddersFixed_RK4
3407#endif
3408 use pm_kind, only: RKG => RK4
3409 procedure(real(RKG)) :: getFunc
3410 type(ridders_type) , intent(in) :: method
3411 real(RKG) , intent(out) :: root
3412 real(RKG) , value :: lb, ub
3413 real(RKG) , value :: lf, uf
3414 real(RKG) , intent(in) :: abstol
3415 integer(IK) , intent(out) :: neval
3416 end subroutine
3417#endif
3418
3419#if RK3_ENABLED
3420 recursive module subroutine setRootRiddersFixed_RK3(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3421#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3422 !DEC$ ATTRIBUTES DLLEXPORT :: setRootRiddersFixed_RK3
3423#endif
3424 use pm_kind, only: RKG => RK3
3425 procedure(real(RKG)) :: getFunc
3426 type(ridders_type) , intent(in) :: method
3427 real(RKG) , intent(out) :: root
3428 real(RKG) , value :: lb, ub
3429 real(RKG) , value :: lf, uf
3430 real(RKG) , intent(in) :: abstol
3431 integer(IK) , intent(out) :: neval
3432 end subroutine
3433#endif
3434
3435#if RK2_ENABLED
3436 recursive module subroutine setRootRiddersFixed_RK2(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3437#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3438 !DEC$ ATTRIBUTES DLLEXPORT :: setRootRiddersFixed_RK2
3439#endif
3440 use pm_kind, only: RKG => RK2
3441 procedure(real(RKG)) :: getFunc
3442 type(ridders_type) , intent(in) :: method
3443 real(RKG) , intent(out) :: root
3444 real(RKG) , value :: lb, ub
3445 real(RKG) , value :: lf, uf
3446 real(RKG) , intent(in) :: abstol
3447 integer(IK) , intent(out) :: neval
3448 end subroutine
3449#endif
3450
3451#if RK1_ENABLED
3452 recursive module subroutine setRootRiddersFixed_RK1(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3453#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3454 !DEC$ ATTRIBUTES DLLEXPORT :: setRootRiddersFixed_RK1
3455#endif
3456 use pm_kind, only: RKG => RK1
3457 procedure(real(RKG)) :: getFunc
3458 type(ridders_type) , intent(in) :: method
3459 real(RKG) , intent(out) :: root
3460 real(RKG) , value :: lb, ub
3461 real(RKG) , value :: lf, uf
3462 real(RKG) , intent(in) :: abstol
3463 integer(IK) , intent(out) :: neval
3464 end subroutine
3465#endif
3466
3467 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3468
3469#if RK5_ENABLED
3470 recursive module subroutine setRootRiddersNiter_RK5(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3471#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3472 !DEC$ ATTRIBUTES DLLEXPORT :: setRootRiddersNiter_RK5
3473#endif
3474 use pm_kind, only: RKG => RK5
3475 procedure(real(RKG)) :: getFunc
3476 type(ridders_type) , intent(in) :: method
3477 real(RKG) , intent(out) :: root
3478 real(RKG) , value :: lb, ub
3479 real(RKG) , value :: lf, uf
3480 real(RKG) , intent(in) :: abstol
3481 integer(IK) , intent(out) :: neval
3482 integer(IK) , intent(in) :: niter
3483 end subroutine
3484#endif
3485
3486#if RK4_ENABLED
3487 recursive module subroutine setRootRiddersNiter_RK4(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3488#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3489 !DEC$ ATTRIBUTES DLLEXPORT :: setRootRiddersNiter_RK4
3490#endif
3491 use pm_kind, only: RKG => RK4
3492 procedure(real(RKG)) :: getFunc
3493 type(ridders_type) , intent(in) :: method
3494 real(RKG) , intent(out) :: root
3495 real(RKG) , value :: lb, ub
3496 real(RKG) , value :: lf, uf
3497 real(RKG) , intent(in) :: abstol
3498 integer(IK) , intent(out) :: neval
3499 integer(IK) , intent(in) :: niter
3500 end subroutine
3501#endif
3502
3503#if RK3_ENABLED
3504 recursive module subroutine setRootRiddersNiter_RK3(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3505#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3506 !DEC$ ATTRIBUTES DLLEXPORT :: setRootRiddersNiter_RK3
3507#endif
3508 use pm_kind, only: RKG => RK3
3509 procedure(real(RKG)) :: getFunc
3510 type(ridders_type) , intent(in) :: method
3511 real(RKG) , intent(out) :: root
3512 real(RKG) , value :: lb, ub
3513 real(RKG) , value :: lf, uf
3514 real(RKG) , intent(in) :: abstol
3515 integer(IK) , intent(out) :: neval
3516 integer(IK) , intent(in) :: niter
3517 end subroutine
3518#endif
3519
3520#if RK2_ENABLED
3521 recursive module subroutine setRootRiddersNiter_RK2(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3522#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3523 !DEC$ ATTRIBUTES DLLEXPORT :: setRootRiddersNiter_RK2
3524#endif
3525 use pm_kind, only: RKG => RK2
3526 procedure(real(RKG)) :: getFunc
3527 type(ridders_type) , intent(in) :: method
3528 real(RKG) , intent(out) :: root
3529 real(RKG) , value :: lb, ub
3530 real(RKG) , value :: lf, uf
3531 real(RKG) , intent(in) :: abstol
3532 integer(IK) , intent(out) :: neval
3533 integer(IK) , intent(in) :: niter
3534 end subroutine
3535#endif
3536
3537#if RK1_ENABLED
3538 recursive module subroutine setRootRiddersNiter_RK1(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3539#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3540 !DEC$ ATTRIBUTES DLLEXPORT :: setRootRiddersNiter_RK1
3541#endif
3542 use pm_kind, only: RKG => RK1
3543 procedure(real(RKG)) :: getFunc
3544 type(ridders_type) , intent(in) :: method
3545 real(RKG) , intent(out) :: root
3546 real(RKG) , value :: lb, ub
3547 real(RKG) , value :: lf, uf
3548 real(RKG) , intent(in) :: abstol
3549 integer(IK) , intent(out) :: neval
3550 integer(IK) , intent(in) :: niter
3551 end subroutine
3552#endif
3553
3554 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3555
3556 end interface
3557
3558 ! TOMS748
3559
3560 interface setRoot
3561
3562 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3563
3564#if RK5_ENABLED
3565 recursive module subroutine setRootTOMS748Fixed_RK5(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3566#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3567 !DEC$ ATTRIBUTES DLLEXPORT :: setRootTOMS748Fixed_RK5
3568#endif
3569 use pm_kind, only: RKG => RK5
3570 procedure(real(RKG)) :: getFunc
3571 type(toms748_type) , intent(in) :: method
3572 real(RKG) , intent(out) :: root
3573 real(RKG) , value :: lb, ub
3574 real(RKG) , value :: lf, uf
3575 real(RKG) , intent(in) :: abstol
3576 integer(IK) , intent(out) :: neval
3577 end subroutine
3578#endif
3579
3580#if RK4_ENABLED
3581 recursive module subroutine setRootTOMS748Fixed_RK4(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3582#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3583 !DEC$ ATTRIBUTES DLLEXPORT :: setRootTOMS748Fixed_RK4
3584#endif
3585 use pm_kind, only: RKG => RK4
3586 procedure(real(RKG)) :: getFunc
3587 type(toms748_type) , intent(in) :: method
3588 real(RKG) , intent(out) :: root
3589 real(RKG) , value :: lb, ub
3590 real(RKG) , value :: lf, uf
3591 real(RKG) , intent(in) :: abstol
3592 integer(IK) , intent(out) :: neval
3593 end subroutine
3594#endif
3595
3596#if RK3_ENABLED
3597 recursive module subroutine setRootTOMS748Fixed_RK3(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3598#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3599 !DEC$ ATTRIBUTES DLLEXPORT :: setRootTOMS748Fixed_RK3
3600#endif
3601 use pm_kind, only: RKG => RK3
3602 procedure(real(RKG)) :: getFunc
3603 type(toms748_type) , intent(in) :: method
3604 real(RKG) , intent(out) :: root
3605 real(RKG) , value :: lb, ub
3606 real(RKG) , value :: lf, uf
3607 real(RKG) , intent(in) :: abstol
3608 integer(IK) , intent(out) :: neval
3609 end subroutine
3610#endif
3611
3612#if RK2_ENABLED
3613 recursive module subroutine setRootTOMS748Fixed_RK2(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3614#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3615 !DEC$ ATTRIBUTES DLLEXPORT :: setRootTOMS748Fixed_RK2
3616#endif
3617 use pm_kind, only: RKG => RK2
3618 procedure(real(RKG)) :: getFunc
3619 type(toms748_type) , intent(in) :: method
3620 real(RKG) , intent(out) :: root
3621 real(RKG) , value :: lb, ub
3622 real(RKG) , value :: lf, uf
3623 real(RKG) , intent(in) :: abstol
3624 integer(IK) , intent(out) :: neval
3625 end subroutine
3626#endif
3627
3628#if RK1_ENABLED
3629 recursive module subroutine setRootTOMS748Fixed_RK1(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3630#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3631 !DEC$ ATTRIBUTES DLLEXPORT :: setRootTOMS748Fixed_RK1
3632#endif
3633 use pm_kind, only: RKG => RK1
3634 procedure(real(RKG)) :: getFunc
3635 type(toms748_type) , intent(in) :: method
3636 real(RKG) , intent(out) :: root
3637 real(RKG) , value :: lb, ub
3638 real(RKG) , value :: lf, uf
3639 real(RKG) , intent(in) :: abstol
3640 integer(IK) , intent(out) :: neval
3641 end subroutine
3642#endif
3643
3644 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3645
3646#if RK5_ENABLED
3647 recursive module subroutine setRootTOMS748Niter_RK5(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3648#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3649 !DEC$ ATTRIBUTES DLLEXPORT :: setRootTOMS748Niter_RK5
3650#endif
3651 use pm_kind, only: RKG => RK5
3652 procedure(real(RKG)) :: getFunc
3653 type(toms748_type) , intent(in) :: method
3654 real(RKG) , intent(out) :: root
3655 real(RKG) , value :: lb, ub
3656 real(RKG) , value :: lf, uf
3657 real(RKG) , intent(in) :: abstol
3658 integer(IK) , intent(out) :: neval
3659 integer(IK) , intent(in) :: niter
3660 end subroutine
3661#endif
3662
3663#if RK4_ENABLED
3664 recursive module subroutine setRootTOMS748Niter_RK4(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3665#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3666 !DEC$ ATTRIBUTES DLLEXPORT :: setRootTOMS748Niter_RK4
3667#endif
3668 use pm_kind, only: RKG => RK4
3669 procedure(real(RKG)) :: getFunc
3670 type(toms748_type) , intent(in) :: method
3671 real(RKG) , intent(out) :: root
3672 real(RKG) , value :: lb, ub
3673 real(RKG) , value :: lf, uf
3674 real(RKG) , intent(in) :: abstol
3675 integer(IK) , intent(out) :: neval
3676 integer(IK) , intent(in) :: niter
3677 end subroutine
3678#endif
3679
3680#if RK3_ENABLED
3681 recursive module subroutine setRootTOMS748Niter_RK3(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3682#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3683 !DEC$ ATTRIBUTES DLLEXPORT :: setRootTOMS748Niter_RK3
3684#endif
3685 use pm_kind, only: RKG => RK3
3686 procedure(real(RKG)) :: getFunc
3687 type(toms748_type) , intent(in) :: method
3688 real(RKG) , intent(out) :: root
3689 real(RKG) , value :: lb, ub
3690 real(RKG) , value :: lf, uf
3691 real(RKG) , intent(in) :: abstol
3692 integer(IK) , intent(out) :: neval
3693 integer(IK) , intent(in) :: niter
3694 end subroutine
3695#endif
3696
3697#if RK2_ENABLED
3698 recursive module subroutine setRootTOMS748Niter_RK2(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3699#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3700 !DEC$ ATTRIBUTES DLLEXPORT :: setRootTOMS748Niter_RK2
3701#endif
3702 use pm_kind, only: RKG => RK2
3703 procedure(real(RKG)) :: getFunc
3704 type(toms748_type) , intent(in) :: method
3705 real(RKG) , intent(out) :: root
3706 real(RKG) , value :: lb, ub
3707 real(RKG) , value :: lf, uf
3708 real(RKG) , intent(in) :: abstol
3709 integer(IK) , intent(out) :: neval
3710 integer(IK) , intent(in) :: niter
3711 end subroutine
3712#endif
3713
3714#if RK1_ENABLED
3715 recursive module subroutine setRootTOMS748Niter_RK1(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3716#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3717 !DEC$ ATTRIBUTES DLLEXPORT :: setRootTOMS748Niter_RK1
3718#endif
3719 use pm_kind, only: RKG => RK1
3720 procedure(real(RKG)) :: getFunc
3721 type(toms748_type) , intent(in) :: method
3722 real(RKG) , intent(out) :: root
3723 real(RKG) , value :: lb, ub
3724 real(RKG) , value :: lf, uf
3725 real(RKG) , intent(in) :: abstol
3726 integer(IK) , intent(out) :: neval
3727 integer(IK) , intent(in) :: niter
3728 end subroutine
3729#endif
3730
3731 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3732
3733 end interface
3734
3735 ! Newton
3736
3737 interface setRoot
3738
3739 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3740
3741#if RK5_ENABLED
3742 recursive module subroutine setRootNewtonFixed_RK5(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3743#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3744 !DEC$ ATTRIBUTES DLLEXPORT :: setRootNewtonFixed_RK5
3745#endif
3746 use pm_kind, only: RKG => RK5
3747 procedure(real(RKG)) :: getFunc
3748 type(newton_type) , intent(in) :: method
3749 real(RKG) , intent(inout) :: root
3750 real(RKG) , value :: lb, ub
3751 real(RKG) , value :: lf, uf
3752 real(RKG) , intent(in) :: abstol
3753 integer(IK) , intent(out) :: neval
3754 end subroutine
3755#endif
3756
3757#if RK4_ENABLED
3758 recursive module subroutine setRootNewtonFixed_RK4(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3759#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3760 !DEC$ ATTRIBUTES DLLEXPORT :: setRootNewtonFixed_RK4
3761#endif
3762 use pm_kind, only: RKG => RK4
3763 procedure(real(RKG)) :: getFunc
3764 type(newton_type) , intent(in) :: method
3765 real(RKG) , intent(inout) :: root
3766 real(RKG) , value :: lb, ub
3767 real(RKG) , value :: lf, uf
3768 real(RKG) , intent(in) :: abstol
3769 integer(IK) , intent(out) :: neval
3770 end subroutine
3771#endif
3772
3773#if RK3_ENABLED
3774 recursive module subroutine setRootNewtonFixed_RK3(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3775#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3776 !DEC$ ATTRIBUTES DLLEXPORT :: setRootNewtonFixed_RK3
3777#endif
3778 use pm_kind, only: RKG => RK3
3779 procedure(real(RKG)) :: getFunc
3780 type(newton_type) , intent(in) :: method
3781 real(RKG) , intent(inout) :: root
3782 real(RKG) , value :: lb, ub
3783 real(RKG) , value :: lf, uf
3784 real(RKG) , intent(in) :: abstol
3785 integer(IK) , intent(out) :: neval
3786 end subroutine
3787#endif
3788
3789#if RK2_ENABLED
3790 recursive module subroutine setRootNewtonFixed_RK2(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3791#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3792 !DEC$ ATTRIBUTES DLLEXPORT :: setRootNewtonFixed_RK2
3793#endif
3794 use pm_kind, only: RKG => RK2
3795 procedure(real(RKG)) :: getFunc
3796 type(newton_type) , intent(in) :: method
3797 real(RKG) , intent(inout) :: root
3798 real(RKG) , value :: lb, ub
3799 real(RKG) , value :: lf, uf
3800 real(RKG) , intent(in) :: abstol
3801 integer(IK) , intent(out) :: neval
3802 end subroutine
3803#endif
3804
3805#if RK1_ENABLED
3806 recursive module subroutine setRootNewtonFixed_RK1(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3807#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3808 !DEC$ ATTRIBUTES DLLEXPORT :: setRootNewtonFixed_RK1
3809#endif
3810 use pm_kind, only: RKG => RK1
3811 procedure(real(RKG)) :: getFunc
3812 type(newton_type) , intent(in) :: method
3813 real(RKG) , intent(inout) :: root
3814 real(RKG) , value :: lb, ub
3815 real(RKG) , value :: lf, uf
3816 real(RKG) , intent(in) :: abstol
3817 integer(IK) , intent(out) :: neval
3818 end subroutine
3819#endif
3820
3821 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3822
3823#if RK5_ENABLED
3824 recursive module subroutine setRootNewtonNiter_RK5(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3825#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3826 !DEC$ ATTRIBUTES DLLEXPORT :: setRootNewtonNiter_RK5
3827#endif
3828 use pm_kind, only: RKG => RK5
3829 procedure(real(RKG)) :: getFunc
3830 type(newton_type) , intent(in) :: method
3831 real(RKG) , intent(inout) :: root
3832 real(RKG) , value :: lb, ub
3833 real(RKG) , value :: lf, uf
3834 real(RKG) , intent(in) :: abstol
3835 integer(IK) , intent(out) :: neval
3836 integer(IK) , intent(in) :: niter
3837 end subroutine
3838#endif
3839
3840#if RK4_ENABLED
3841 recursive module subroutine setRootNewtonNiter_RK4(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3842#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3843 !DEC$ ATTRIBUTES DLLEXPORT :: setRootNewtonNiter_RK4
3844#endif
3845 use pm_kind, only: RKG => RK4
3846 procedure(real(RKG)) :: getFunc
3847 type(newton_type) , intent(in) :: method
3848 real(RKG) , intent(inout) :: root
3849 real(RKG) , value :: lb, ub
3850 real(RKG) , value :: lf, uf
3851 real(RKG) , intent(in) :: abstol
3852 integer(IK) , intent(out) :: neval
3853 integer(IK) , intent(in) :: niter
3854 end subroutine
3855#endif
3856
3857#if RK3_ENABLED
3858 recursive module subroutine setRootNewtonNiter_RK3(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3859#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3860 !DEC$ ATTRIBUTES DLLEXPORT :: setRootNewtonNiter_RK3
3861#endif
3862 use pm_kind, only: RKG => RK3
3863 procedure(real(RKG)) :: getFunc
3864 type(newton_type) , intent(in) :: method
3865 real(RKG) , intent(inout) :: root
3866 real(RKG) , value :: lb, ub
3867 real(RKG) , value :: lf, uf
3868 real(RKG) , intent(in) :: abstol
3869 integer(IK) , intent(out) :: neval
3870 integer(IK) , intent(in) :: niter
3871 end subroutine
3872#endif
3873
3874#if RK2_ENABLED
3875 recursive module subroutine setRootNewtonNiter_RK2(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3876#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3877 !DEC$ ATTRIBUTES DLLEXPORT :: setRootNewtonNiter_RK2
3878#endif
3879 use pm_kind, only: RKG => RK2
3880 procedure(real(RKG)) :: getFunc
3881 type(newton_type) , intent(in) :: method
3882 real(RKG) , intent(inout) :: root
3883 real(RKG) , value :: lb, ub
3884 real(RKG) , value :: lf, uf
3885 real(RKG) , intent(in) :: abstol
3886 integer(IK) , intent(out) :: neval
3887 integer(IK) , intent(in) :: niter
3888 end subroutine
3889#endif
3890
3891#if RK1_ENABLED
3892 recursive module subroutine setRootNewtonNiter_RK1(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
3893#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3894 !DEC$ ATTRIBUTES DLLEXPORT :: setRootNewtonNiter_RK1
3895#endif
3896 use pm_kind, only: RKG => RK1
3897 procedure(real(RKG)) :: getFunc
3898 type(newton_type) , intent(in) :: method
3899 real(RKG) , intent(inout) :: root
3900 real(RKG) , value :: lb, ub
3901 real(RKG) , value :: lf, uf
3902 real(RKG) , intent(in) :: abstol
3903 integer(IK) , intent(out) :: neval
3904 integer(IK) , intent(in) :: niter
3905 end subroutine
3906#endif
3907
3908 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3909
3910 end interface
3911
3912 ! Halley
3913
3914 interface setRoot
3915
3916 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3917
3918#if RK5_ENABLED
3919 recursive module subroutine setRootHalleyFixed_RK5(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3920#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3921 !DEC$ ATTRIBUTES DLLEXPORT :: setRootHalleyFixed_RK5
3922#endif
3923 use pm_kind, only: RKG => RK5
3924 procedure(real(RKG)) :: getFunc
3925 type(halley_type) , intent(in) :: method
3926 real(RKG) , intent(inout) :: root
3927 real(RKG) , value :: lb, ub
3928 real(RKG) , value :: lf, uf
3929 real(RKG) , intent(in) :: abstol
3930 integer(IK) , intent(out) :: neval
3931 end subroutine
3932#endif
3933
3934#if RK4_ENABLED
3935 recursive module subroutine setRootHalleyFixed_RK4(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3936#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3937 !DEC$ ATTRIBUTES DLLEXPORT :: setRootHalleyFixed_RK4
3938#endif
3939 use pm_kind, only: RKG => RK4
3940 procedure(real(RKG)) :: getFunc
3941 type(halley_type) , intent(in) :: method
3942 real(RKG) , intent(inout) :: root
3943 real(RKG) , value :: lb, ub
3944 real(RKG) , value :: lf, uf
3945 real(RKG) , intent(in) :: abstol
3946 integer(IK) , intent(out) :: neval
3947 end subroutine
3948#endif
3949
3950#if RK3_ENABLED
3951 recursive module subroutine setRootHalleyFixed_RK3(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3952#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3953 !DEC$ ATTRIBUTES DLLEXPORT :: setRootHalleyFixed_RK3
3954#endif
3955 use pm_kind, only: RKG => RK3
3956 procedure(real(RKG)) :: getFunc
3957 type(halley_type) , intent(in) :: method
3958 real(RKG) , intent(inout) :: root
3959 real(RKG) , value :: lb, ub
3960 real(RKG) , value :: lf, uf
3961 real(RKG) , intent(in) :: abstol
3962 integer(IK) , intent(out) :: neval
3963 end subroutine
3964#endif
3965
3966#if RK2_ENABLED
3967 recursive module subroutine setRootHalleyFixed_RK2(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3968#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3969 !DEC$ ATTRIBUTES DLLEXPORT :: setRootHalleyFixed_RK2
3970#endif
3971 use pm_kind, only: RKG => RK2
3972 procedure(real(RKG)) :: getFunc
3973 type(halley_type) , intent(in) :: method
3974 real(RKG) , intent(inout) :: root
3975 real(RKG) , value :: lb, ub
3976 real(RKG) , value :: lf, uf
3977 real(RKG) , intent(in) :: abstol
3978 integer(IK) , intent(out) :: neval
3979 end subroutine
3980#endif
3981
3982#if RK1_ENABLED
3983 recursive module subroutine setRootHalleyFixed_RK1(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
3984#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
3985 !DEC$ ATTRIBUTES DLLEXPORT :: setRootHalleyFixed_RK1
3986#endif
3987 use pm_kind, only: RKG => RK1
3988 procedure(real(RKG)) :: getFunc
3989 type(halley_type) , intent(in) :: method
3990 real(RKG) , intent(inout) :: root
3991 real(RKG) , value :: lb, ub
3992 real(RKG) , value :: lf, uf
3993 real(RKG) , intent(in) :: abstol
3994 integer(IK) , intent(out) :: neval
3995 end subroutine
3996#endif
3997
3998 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3999
4000#if RK5_ENABLED
4001 recursive module subroutine setRootHalleyNiter_RK5(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
4002#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
4003 !DEC$ ATTRIBUTES DLLEXPORT :: setRootHalleyNiter_RK5
4004#endif
4005 use pm_kind, only: RKG => RK5
4006 procedure(real(RKG)) :: getFunc
4007 type(halley_type) , intent(in) :: method
4008 real(RKG) , intent(inout) :: root
4009 real(RKG) , value :: lb, ub
4010 real(RKG) , value :: lf, uf
4011 real(RKG) , intent(in) :: abstol
4012 integer(IK) , intent(out) :: neval
4013 integer(IK) , intent(in) :: niter
4014 end subroutine
4015#endif
4016
4017#if RK4_ENABLED
4018 recursive module subroutine setRootHalleyNiter_RK4(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
4019#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
4020 !DEC$ ATTRIBUTES DLLEXPORT :: setRootHalleyNiter_RK4
4021#endif
4022 use pm_kind, only: RKG => RK4
4023 procedure(real(RKG)) :: getFunc
4024 type(halley_type) , intent(in) :: method
4025 real(RKG) , intent(inout) :: root
4026 real(RKG) , value :: lb, ub
4027 real(RKG) , value :: lf, uf
4028 real(RKG) , intent(in) :: abstol
4029 integer(IK) , intent(out) :: neval
4030 integer(IK) , intent(in) :: niter
4031 end subroutine
4032#endif
4033
4034#if RK3_ENABLED
4035 recursive module subroutine setRootHalleyNiter_RK3(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
4036#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
4037 !DEC$ ATTRIBUTES DLLEXPORT :: setRootHalleyNiter_RK3
4038#endif
4039 use pm_kind, only: RKG => RK3
4040 procedure(real(RKG)) :: getFunc
4041 type(halley_type) , intent(in) :: method
4042 real(RKG) , intent(inout) :: root
4043 real(RKG) , value :: lb, ub
4044 real(RKG) , value :: lf, uf
4045 real(RKG) , intent(in) :: abstol
4046 integer(IK) , intent(out) :: neval
4047 integer(IK) , intent(in) :: niter
4048 end subroutine
4049#endif
4050
4051#if RK2_ENABLED
4052 recursive module subroutine setRootHalleyNiter_RK2(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
4053#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
4054 !DEC$ ATTRIBUTES DLLEXPORT :: setRootHalleyNiter_RK2
4055#endif
4056 use pm_kind, only: RKG => RK2
4057 procedure(real(RKG)) :: getFunc
4058 type(halley_type) , intent(in) :: method
4059 real(RKG) , intent(inout) :: root
4060 real(RKG) , value :: lb, ub
4061 real(RKG) , value :: lf, uf
4062 real(RKG) , intent(in) :: abstol
4063 integer(IK) , intent(out) :: neval
4064 integer(IK) , intent(in) :: niter
4065 end subroutine
4066#endif
4067
4068#if RK1_ENABLED
4069 recursive module subroutine setRootHalleyNiter_RK1(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
4070#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
4071 !DEC$ ATTRIBUTES DLLEXPORT :: setRootHalleyNiter_RK1
4072#endif
4073 use pm_kind, only: RKG => RK1
4074 procedure(real(RKG)) :: getFunc
4075 type(halley_type) , intent(in) :: method
4076 real(RKG) , intent(inout) :: root
4077 real(RKG) , value :: lb, ub
4078 real(RKG) , value :: lf, uf
4079 real(RKG) , intent(in) :: abstol
4080 integer(IK) , intent(out) :: neval
4081 integer(IK) , intent(in) :: niter
4082 end subroutine
4083#endif
4084
4085 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
4086
4087 end interface
4088
4089 ! Schroder
4090
4091 interface setRoot
4092
4093 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
4094
4095#if RK5_ENABLED
4096 recursive module subroutine setRootSchroderFixed_RK5(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
4097#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
4098 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSchroderFixed_RK5
4099#endif
4100 use pm_kind, only: RKG => RK5
4101 procedure(real(RKG)) :: getFunc
4102 type(schroder_type) , intent(in) :: method
4103 real(RKG) , intent(inout) :: root
4104 real(RKG) , value :: lb, ub
4105 real(RKG) , value :: lf, uf
4106 real(RKG) , intent(in) :: abstol
4107 integer(IK) , intent(out) :: neval
4108 end subroutine
4109#endif
4110
4111#if RK4_ENABLED
4112 recursive module subroutine setRootSchroderFixed_RK4(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
4113#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
4114 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSchroderFixed_RK4
4115#endif
4116 use pm_kind, only: RKG => RK4
4117 procedure(real(RKG)) :: getFunc
4118 type(schroder_type) , intent(in) :: method
4119 real(RKG) , intent(inout) :: root
4120 real(RKG) , value :: lb, ub
4121 real(RKG) , value :: lf, uf
4122 real(RKG) , intent(in) :: abstol
4123 integer(IK) , intent(out) :: neval
4124 end subroutine
4125#endif
4126
4127#if RK3_ENABLED
4128 recursive module subroutine setRootSchroderFixed_RK3(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
4129#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
4130 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSchroderFixed_RK3
4131#endif
4132 use pm_kind, only: RKG => RK3
4133 procedure(real(RKG)) :: getFunc
4134 type(schroder_type) , intent(in) :: method
4135 real(RKG) , intent(inout) :: root
4136 real(RKG) , value :: lb, ub
4137 real(RKG) , value :: lf, uf
4138 real(RKG) , intent(in) :: abstol
4139 integer(IK) , intent(out) :: neval
4140 end subroutine
4141#endif
4142
4143#if RK2_ENABLED
4144 recursive module subroutine setRootSchroderFixed_RK2(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
4145#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
4146 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSchroderFixed_RK2
4147#endif
4148 use pm_kind, only: RKG => RK2
4149 procedure(real(RKG)) :: getFunc
4150 type(schroder_type) , intent(in) :: method
4151 real(RKG) , intent(inout) :: root
4152 real(RKG) , value :: lb, ub
4153 real(RKG) , value :: lf, uf
4154 real(RKG) , intent(in) :: abstol
4155 integer(IK) , intent(out) :: neval
4156 end subroutine
4157#endif
4158
4159#if RK1_ENABLED
4160 recursive module subroutine setRootSchroderFixed_RK1(method, getFunc, root, lb, ub, lf, uf, abstol, neval)
4161#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
4162 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSchroderFixed_RK1
4163#endif
4164 use pm_kind, only: RKG => RK1
4165 procedure(real(RKG)) :: getFunc
4166 type(schroder_type) , intent(in) :: method
4167 real(RKG) , intent(inout) :: root
4168 real(RKG) , value :: lb, ub
4169 real(RKG) , value :: lf, uf
4170 real(RKG) , intent(in) :: abstol
4171 integer(IK) , intent(out) :: neval
4172 end subroutine
4173#endif
4174
4175 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
4176
4177#if RK5_ENABLED
4178 recursive module subroutine setRootSchroderNiter_RK5(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
4179#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
4180 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSchroderNiter_RK5
4181#endif
4182 use pm_kind, only: RKG => RK5
4183 procedure(real(RKG)) :: getFunc
4184 type(schroder_type) , intent(in) :: method
4185 real(RKG) , intent(inout) :: root
4186 real(RKG) , value :: lb, ub
4187 real(RKG) , value :: lf, uf
4188 real(RKG) , intent(in) :: abstol
4189 integer(IK) , intent(out) :: neval
4190 integer(IK) , intent(in) :: niter
4191 end subroutine
4192#endif
4193
4194#if RK4_ENABLED
4195 recursive module subroutine setRootSchroderNiter_RK4(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
4196#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
4197 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSchroderNiter_RK4
4198#endif
4199 use pm_kind, only: RKG => RK4
4200 procedure(real(RKG)) :: getFunc
4201 type(schroder_type) , intent(in) :: method
4202 real(RKG) , intent(inout) :: root
4203 real(RKG) , value :: lb, ub
4204 real(RKG) , value :: lf, uf
4205 real(RKG) , intent(in) :: abstol
4206 integer(IK) , intent(out) :: neval
4207 integer(IK) , intent(in) :: niter
4208 end subroutine
4209#endif
4210
4211#if RK3_ENABLED
4212 recursive module subroutine setRootSchroderNiter_RK3(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
4213#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
4214 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSchroderNiter_RK3
4215#endif
4216 use pm_kind, only: RKG => RK3
4217 procedure(real(RKG)) :: getFunc
4218 type(schroder_type) , intent(in) :: method
4219 real(RKG) , intent(inout) :: root
4220 real(RKG) , value :: lb, ub
4221 real(RKG) , value :: lf, uf
4222 real(RKG) , intent(in) :: abstol
4223 integer(IK) , intent(out) :: neval
4224 integer(IK) , intent(in) :: niter
4225 end subroutine
4226#endif
4227
4228#if RK2_ENABLED
4229 recursive module subroutine setRootSchroderNiter_RK2(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
4230#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
4231 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSchroderNiter_RK2
4232#endif
4233 use pm_kind, only: RKG => RK2
4234 procedure(real(RKG)) :: getFunc
4235 type(schroder_type) , intent(in) :: method
4236 real(RKG) , intent(inout) :: root
4237 real(RKG) , value :: lb, ub
4238 real(RKG) , value :: lf, uf
4239 real(RKG) , intent(in) :: abstol
4240 integer(IK) , intent(out) :: neval
4241 integer(IK) , intent(in) :: niter
4242 end subroutine
4243#endif
4244
4245#if RK1_ENABLED
4246 recursive module subroutine setRootSchroderNiter_RK1(method, getFunc, root, lb, ub, lf, uf, abstol, neval, niter)
4247#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
4248 !DEC$ ATTRIBUTES DLLEXPORT :: setRootSchroderNiter_RK1
4249#endif
4250 use pm_kind, only: RKG => RK1
4251 procedure(real(RKG)) :: getFunc
4252 type(schroder_type) , intent(in) :: method
4253 real(RKG) , intent(inout) :: root
4254 real(RKG) , value :: lb, ub
4255 real(RKG) , value :: lf, uf
4256 real(RKG) , intent(in) :: abstol
4257 integer(IK) , intent(out) :: neval
4258 integer(IK) , intent(in) :: niter
4259 end subroutine
4260#endif
4261
4262 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
4263
4264 end interface
4265
4266!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
4267
4268end module pm_mathRoot
Generate and return a root of a specified continuous real-valued one-dimensional mathematical functio...
Return a root of a specified continuous real-valued one-dimensional mathematical function such that ...
This module defines the relevant Fortran kind type-parameters frequently used in the ParaMonte librar...
Definition: pm_kind.F90:268
integer, parameter RK5
Definition: pm_kind.F90:478
integer, parameter RK4
Definition: pm_kind.F90:489
integer, parameter RK2
Definition: pm_kind.F90:511
integer, parameter RK3
Definition: pm_kind.F90:500
integer, parameter LK
The default logical kind in the ParaMonte library: kind(.true.) in Fortran, kind(....
Definition: pm_kind.F90:541
integer, parameter IK
The default integer kind in the ParaMonte library: int32 in Fortran, c_int32_t in C-Fortran Interoper...
Definition: pm_kind.F90:540
integer, parameter SK
The default character kind in the ParaMonte library: kind("a") in Fortran, c_char in C-Fortran Intero...
Definition: pm_kind.F90:539
integer, parameter RK1
Definition: pm_kind.F90:522
This module contains classes and procedures for computing the roots of one-dimensional continuous mat...
type(brent_type), parameter brent
This is a scalar parameter object of type brent_type that is exclusively used to signify the use of B...
type(ridders_type), parameter ridders
This is a scalar parameter object of type ridders_type that is exclusively used to signify the use of...
type(toms748_type), parameter toms748
This is a scalar parameter object of type toms748_type that is exclusively used to signify the use of...
type(newton_type), parameter newton
This is a scalar parameter object of type newton_type that is exclusively used to signify the use of ...
type(bisection_type), parameter bisection
This is a scalar parameter object of type bisection_type that is exclusively used to signify the use ...
type(false_type), parameter false
This is a scalar parameter object of type false_type that is exclusively used to signify the use of F...
type(halley_type), parameter halley
This is a scalar parameter object of type halley_type that is exclusively used to signify the use of ...
type(schroder_type), parameter schroder
This is a scalar parameter object of type schroder_type that is exclusively used to signify the use o...
type(secant_type), parameter secant
This is a scalar parameter object of type secant_type that is exclusively used to signify the use of ...
character(*, SK), parameter MODULE_NAME
This is a concrete derived type whose instances are exclusively used to signify the use of the Bisect...
This is an abstract derived type for constructing concrete derived types to distinguish various proce...
This is a concrete derived type whose instances are exclusively used to signify the use of the Brent ...
This is a concrete derived type whose instances are exclusively used to signify the use of the False-...
This is a concrete derived type whose instances are exclusively used to signify the use of the Halley...
This is an abstract derived type for constructing concrete derived types to distinguish various proce...
This is an abstract derived type for constructing concrete derived types to distinguish various proce...
This is an abstract derived type for constructing concrete derived types to distinguish various proce...
This is a concrete derived type whose instances are exclusively used to signify the use of the Newton...
This is a concrete derived type whose instances are exclusively used to signify the use of the Ridder...
This is a concrete derived type whose instances are exclusively used to signify the use of the Schrod...
This is a concrete derived type whose instances are exclusively used to signify the use of the Secant...
This is a concrete derived type whose instances are exclusively used to signify the use of the TOMS74...