Back
C
/******************************************************************************
* LICENSE *
******************************************************************************
* This file is part of mitx_mathematics_programming_examples. *
* *
* mitx_mathematics_programming_examples is free software: you can *
* redistribute it and/or modify it under the terms of the GNU General *
* Public License as published by the Free Software Foundation, either *
* version 3 of the License, or (at your option) any later version. *
* *
* mitx_mathematics_programming_examples is distributed in the hope that *
* it will be useful but WITHOUT ANY WARRANTY; without even the implied *
* warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. *
* See the GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with mitx_mathematics_programming_examples. If not, see *
* <https://www.gnu.org/licenses/>. *
******************************************************************************
* Purpose: *
* Calculates the root of a function using Steffensen's method. *
******************************************************************************
* Author: Ryan Maguire *
* Date: 2025/04/17 *
******************************************************************************/
/* stdio.h provides the "printf" function, used for printing text. */
#include <stdio.h>
/* Floating-point absolute value function, fabs, provided here. */
#include <math.h>
/* Function pointer notation is a little confusing. Create a typedef for it *
* so we do not need to explicitly use it later. */
typedef double (*function)(double);
/* Computes the root of a function using Steffensen's method. */
static double steffensens_method(function f, double x)
{
/* Steffensen's method is iterative and converges very quickly. *
* Because of this we may exit the function after a few iterations. */
const unsigned int maximum_number_of_iterations = 16U;
/* The maximum allowed error. This is 4x double precision epsilon. */
const double epsilon = 8.881784197001252E-16;
/* Variable for keeping track of how many iterations we have performed. */
unsigned int iters;
/* The method starts at the provided guess point and updates iteratively.*/
double xn = x;
/* Iteratively apply Steffensen's method to find the root. */
for (iters = 0; iters < maximum_number_of_iterations; ++iters)
{
/* Steffensen's method needs the evaluations f(x) and f(x + f(x)), *
* in particular the denominator is f(x + f(x)) / f(x) - 1. Compute. */
const double f_xn = f(xn);
const double g_xn = f(xn + f_xn) / f_xn - 1.0;
/* Like Newton's method, the new point is obtained by subtracting *
* the ratio. g(x) = f(x + f(x)) / f(x) - 1 acts as the derivative *
* of f, but we do not explicitly need to calculate f'(x). */
xn = xn - f_xn / g_xn;
/* If f(x) is very small, we are close to a root and can break out *
* of this for loop. Check for this. */
if (fabs(f_xn) < epsilon)
break;
}
/* Like Newton's method, and like Heron's method, the convergence is *
* quadratic. After a few iterations we will be very to close a root. */
return xn;
}
/* End of steffensens_method. */
/* sqrt(2) is a root to the function f(x) = 2 - x^2. Provide this. */
static double func(double x)
{
return 2.0 - x*x;
}
/* End of func. */
/* Main routine used for testing our implementation of Steffensen's method. */
int main(void)
{
/* The initial guess point for Steffensen's method. */
const double x = 2.0;
/* Calculate the square root and print it to the screen. If we have *
* written things correctly, we should get 1.414..., which is sqrt(2). */
const double sqrt_x = steffensens_method(func, x);
printf("sqrt(%.1f) = %.16f\n", x, sqrt_x);
return 0;
}
C3
/******************************************************************************
* LICENSE *
******************************************************************************
* This file is part of mitx_mathematics_programming_examples. *
* *
* mitx_mathematics_programming_examples is free software: you can *
* redistribute it and/or modify it under the terms of the GNU General *
* Public License as published by the Free Software Foundation, either *
* version 3 of the License, or (at your option) any later version. *
* *
* mitx_mathematics_programming_examples is distributed in the hope that *
* it will be useful but WITHOUT ANY WARRANTY; without even the implied *
* warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. *
* See the GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with mitx_mathematics_programming_examples. If not, see *
* <https://www.gnu.org/licenses/>. *
******************************************************************************
* Purpose: *
* Calculates the root of a function using Steffensen's method. *
******************************************************************************
* Author: Ryan Maguire *
* Date: 2025/05/25 *
******************************************************************************/
/* std::io provides the "printfn" function, used for printing text. */
import std::io;
/* Floating-point absolute value function, abs, provided here. */
import std::math;
/* Create an alias for functions of the form f: R -> R. */
alias RealFunc = fn double(double);
/* Computes the root of a function using Steffensen's method. */
fn double steffensens_method(RealFunc f, double x)
{
/* Steffensen's method is iterative and converges very quickly. *
* Because of this we may exit the function after a few iterations. */
const uint MAXIMUM_NUMBER_OF_ITERATIONS = 16U;
/* The maximum allowed error. This is 4x double precision epsilon. */
const double EPSILON = 8.881784197001252E-16;
/* Variable for keeping track of how many iterations we have performed. */
uint iters;
/* The method starts at the provided guess point and updates iteratively.*/
double xn = x;
/* Iteratively apply Steffensen's method to find the root. */
for (iters = 0; iters < MAXIMUM_NUMBER_OF_ITERATIONS; ++iters)
{
/* Steffensen's method needs the evaluations f(x) and f(x + f(x)), *
* in particular the denominator is f(x + f(x)) / f(x) - 1. Compute. */
double f_xn = f(xn);
double g_xn = f(xn + f_xn) / f_xn - 1.0;
/* Like Newton's method, the new point is obtained by subtracting *
* the ratio. g(x) = f(x + f(x)) / f(x) - 1 acts as the derivative *
* of f, but we do not explicitly need to calculate f'(x). */
xn = xn - f_xn / g_xn;
/* If f(x) is very small, we are close to a root and can break out *
* of this for loop. Check for this. */
if (math::abs(f_xn) < EPSILON) {
break;
}
}
/* Like Newton's method, and like Heron's method, the convergence is *
* quadratic. After a few iterations we will be very to close a root. */
return xn;
}
/* End of steffensens_method. */
/* sqrt(2) is a root to the function f(x) = 2 - x^2. Provide this. */
fn double func(double x)
{
return 2.0 - x*x;
}
/* End of func. */
/* Main routine used for testing our implementation of Steffensen's method. */
fn void main()
{
/* The initial guess point for Steffensen's method. */
const double X = 2.0;
/* Calculate the square root and print it to the screen. If we have *
* written things correctly, we should get 1.414..., which is sqrt(2). */
double sqrt_x = steffensens_method(&func, X);
io::printfn("sqrt(%.1f) = %.16f", X, sqrt_x);
}
C++
/******************************************************************************
* LICENSE *
******************************************************************************
* This file is part of mitx_mathematics_programming_examples. *
* *
* mitx_mathematics_programming_examples is free software: you can *
* redistribute it and/or modify it under the terms of the GNU General *
* Public License as published by the Free Software Foundation, either *
* version 3 of the License, or (at your option) any later version. *
* *
* mitx_mathematics_programming_examples is distributed in the hope that *
* it will be useful but WITHOUT ANY WARRANTY; without even the implied *
* warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. *
* See the GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with mitx_mathematics_programming_examples. If not, see *
* <https://www.gnu.org/licenses/>. *
******************************************************************************
* Purpose: *
* Calculates the root of a function using Steffensen's method. *
******************************************************************************
* Author: Ryan Maguire *
* Date: 2025/05/22 *
******************************************************************************/
/* cstdio provides the "printf" function, used for printing text. */
#include <cstdio>
/* Floating-point absolute value function, fabs, provided here. */
#include <cmath>
/* Function pointer notation is a little confusing. Create a typedef for it *
* so we do not need to explicitly use it later. */
typedef double (*function)(double);
/* Computes the root of a function using Steffensen's method. */
class Steffensen {
/* Steffensen's method is iterative and converges very quickly. *
* Because of this we may exit the function after a few iterations. */
static const unsigned int maximum_number_of_iterations = 16U;
/* We want the function visible outside the class. Declare it public. */
public:
/* Computes the root of a function using Steffensen's method. */
static double root(function f, double x)
{
/* Maximum allowed error. This is 4x double precision epsilon. */
const double epsilon = 8.881784197001252E-16;
/* Variable keeping track of how many iterations we perform. */
unsigned int iters;
/* The method starts at the guess point and updates iteratively. */
double xn = x;
/* Iteratively apply Steffensen's method to find the root. */
for (iters = 0; iters < maximum_number_of_iterations; ++iters)
{
/* Steffensen's method needs both f(x) and f(x + f(x)), *
* in particular the denominator is f(x + f(x)) / f(x) - 1. */
double f_xn = f(xn);
double g_xn = f(xn + f_xn) / f_xn - 1.0;
/* Like Newton's method the new point is obtained by *
* subtracting the ratio. g(x) = f(x + f(x))/f(x) - 1 acts *
* as the derivative of f, but we do not explicitly need to *
* calculate f'(x). */
xn = xn - f_xn / g_xn;
/* If f(x) is very small, we are close to a root and can *
* break out of this for loop. Check for this. */
if (std::fabs(f_xn) < epsilon)
break;
}
/* Like Newton's method and Heron's method, the convergence is *
* quadratic. After a few iterations we will be close a root. */
return xn;
}
/* End of root. */
};
/* sqrt(2) is a root to the function f(x) = 2 - x^2. Provide this. */
static double func(double x)
{
return 2.0 - x*x;
}
/* End of func. */
/* Main routine used for testing our implementation of Steffensen's method. */
int main(void)
{
/* The initial guess point for Steffensen's method. */
const double x = 2.0;
/* Calculate the square root and print it to the screen. If we have *
* written things correctly, we should get 1.414..., which is sqrt(2). */
const double sqrt_x = Steffensen::root(func, x);
std::printf("sqrt(%.1f) = %.16f\n", x, sqrt_x);
return 0;
}
C#
/******************************************************************************
* LICENSE *
******************************************************************************
* This file is part of mitx_mathematics_programming_examples. *
* *
* mitx_mathematics_programming_examples is free software: you can *
* redistribute it and/or modify it under the terms of the GNU General *
* Public License as published by the Free Software Foundation, either *
* version 3 of the License, or (at your option) any later version. *
* *
* mitx_mathematics_programming_examples is distributed in the hope that *
* it will be useful but WITHOUT ANY WARRANTY; without even the implied *
* warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. *
* See the GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with mitx_mathematics_programming_examples. If not, see *
* <https://www.gnu.org/licenses/>. *
******************************************************************************
* Purpose: *
* Calculates the root of a function using Steffensen's method. *
******************************************************************************
* Author: Ryan Maguire *
* Date: 2025/05/21 *
******************************************************************************/
/* Console.WriteLine and Math.Abs are both provided here. */
using System;
/* Create an alias for functions of the form f: R -> R. */
using function = System.Func<double, double>;
/* Class providing an implementation of Steffensen's method. */
class Steffensen {
/* Steffensen's method is iterative and converges very quickly. *
* Because of this we may exit the function after a few iterations. */
const uint maximumNumberOfIterations = 16U;
/* The maximum allowed error. This is 4x double precision epsilon. */
const double epsilon = 8.881784197001252E-16;
/* Computes the root of a function using Steffensen's method. */
static double Root(function f, double x)
{
/* Variable for keeping track of how many iterations we perform. */
uint iters;
/* The method starts at the guess point and updates iteratively. */
double xn = x;
/* Iteratively apply Steffensen's method to find the root. */
for (iters = 0; iters < maximumNumberOfIterations; ++iters)
{
/* Steffensen's method needs the evaluations f(x) and f(x+f(x)), *
* in particular the denominator is f(x + f(x)) / f(x) - 1. */
double f_xn = f(xn);
double g_xn = f(xn + f_xn) / f_xn - 1.0;
/* Like Newton's method the new point is obtained by subtracting *
* the ratio. g(x) = f(x + f(x))/f(x) - 1 acts as the derivative *
* of f, but we do not explicitly need to calculate f'(x). */
xn = xn - f_xn / g_xn;
/* If f(x) is very small, we are close to a root and can break *
* out of this for loop. Check for this. */
if (Math.Abs(f_xn) < epsilon)
break;
}
/* Like Newton's method, and like Heron's method, the convergence is *
* quadratic. After a few iterations we will be very to close a root.*/
return xn;
}
/* End of Root. */
/* The function we'll test. The roots are +/- sqrt(2). */
static double func(double x)
{
return 2.0 - x*x;
}
/* Main routine used for testing our Steffensen implementation. */
static void Main()
{
/* The initial guess point for Steffensen's method. */
const double x = 2.0;
/* Calculate the square root and print it to the screen. If we have *
* written things correctly we should get 1.414..., which is sqrt(2).*/
double sqrt_x = Root(func, x);
Console.WriteLine($"sqrt({x}) = {sqrt_x}");
}
}
D
/******************************************************************************
* LICENSE *
******************************************************************************
* This file is part of mitx_mathematics_programming_examples. *
* *
* mitx_mathematics_programming_examples is free software: you can *
* redistribute it and/or modify it under the terms of the GNU General *
* Public License as published by the Free Software Foundation, either *
* version 3 of the License, or (at your option) any later version. *
* *
* mitx_mathematics_programming_examples is distributed in the hope that *
* it will be useful but WITHOUT ANY WARRANTY; without even the implied *
* warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. *
* See the GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with mitx_mathematics_programming_examples. If not, see *
* <https://www.gnu.org/licenses/>. *
******************************************************************************
* Purpose: *
* Calculates the root of a function using Steffensen's method. *
******************************************************************************
* Author: Ryan Maguire *
* Date: 2025/05/22 *
******************************************************************************/
/* Create an alias for real-valued functions f: R -> R. */
alias realfunc = double function(double) pure nothrow @nogc @safe;
/* Computes the root of a function using Steffensen's method. */
static double steffensens_method(realfunc f, double x)
{
/* Floating-point absolute value (fabs) provided here. */
import std.math : fabs;
/* Steffensen's method is iterative and converges very quickly. *
* Because of this we may exit the function after a few iterations. */
const uint maximum_number_of_iterations = 16U;
/* The maximum allowed error. This is 4x double precision epsilon. */
const double epsilon = 8.881784197001252E-16;
/* Variable for keeping track of how many iterations we have performed. */
uint iters;
/* The method starts at the provided guess point and updates iteratively.*/
double xn = x;
/* Iteratively apply Steffensen's method to find the root. */
for (iters = 0; iters < maximum_number_of_iterations; ++iters)
{
/* Steffensen's method needs the evaluations f(x) and f(x + f(x)), *
* in particular the denominator is f(x + f(x)) / f(x) - 1. Compute. */
const double f_xn = f(xn);
const double g_xn = f(xn + f_xn) / f_xn - 1.0;
/* Like Newton's method, the new point is obtained by subtracting *
* the ratio. g(x) = f(x + f(x)) / f(x) - 1 acts as the derivative *
* of f, but we do not explicitly need to calculate f'(x). */
xn = xn - f_xn / g_xn;
/* If f(x) is very small, we are close to a root and can break out *
* of this for loop. Check for this. */
if (fabs(f_xn) <= epsilon)
break;
}
/* Like Newton's method, and like Heron's method, the convergence is *
* quadratic. After a few iterations we will be very to close a root. */
return xn;
}
/* End of steffensens_method. */
/* sqrt(2) is a root to the function f(x) = 2 - x^2. Provide this. */
static double func(double x)
pure nothrow @nogc @safe
{
return 2.0 - x*x;
}
/* End of func. */
/* Main routine used for testing our implementation of Steffensen's method. */
int main()
{
/* stdio contains "printf", used for printing text. */
import std.stdio : printf;
/* The initial guess point for Steffensen's method. */
const double x = 2.0;
/* Calculate the square root and print it to the screen. If we have *
* written things correctly, we should get 1.414..., which is sqrt(2). */
const double sqrt_x = steffensens_method(&func, x);
printf("sqrt(%.1f) = %.16f\n", x, sqrt_x);
return 0;
}
Fortran
!------------------------------------------------------------------------------!
! LICENSE !
!------------------------------------------------------------------------------!
! This file is part of mitx_mathematics_programming_examples. !
! !
! mitx_mathematics_programming_examples is free software: you can !
! redistribute it and/or modify it under the terms of the GNU General Public !
! License as published by the Free Software Foundation, either version 3 of !
! the License, or (at your option) any later version. !
! !
! mitx_mathematics_programming_examples is distributed in the hope that it !
! will be useful but WITHOUT ANY WARRANTY; without even the implied warranty !
! of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the !
! GNU General Public License for more details. !
! !
! You should have received a copy of the GNU General Public License !
! along with mitx_mathematics_programming_examples. If not, see !
! <https://www.gnu.org/licenses/>. !
!------------------------------------------------------------------------------!
! Purpose: !
! Calculates roots using Steffensen's method. !
!------------------------------------------------------------------------------!
! Author: Ryan Maguire !
! Date: 2025/05/25 !
!------------------------------------------------------------------------------!
MODULE STEFFENSEN
IMPLICIT NONE
! Steffensen's method is iterative. The convergence is quadratic, meaning
! the number of accurate decimals doubles with each iteration. Because of
! this we can halt the algorithm after a few steps.
INTEGER :: MAXIMUM_NUMBER_OF_ITERATIONS = 16
! Maximum allowed error. This is 4x double precision epsilon.
REAL :: EPSILON = 8.881784197001252E-16
CONTAINS
!--------------------------------------------------------------------------!
! Function: !
! STEFFENSENS_METHOD !
! Purpose: !
! Computes roots using Steffensen's method. !
! Arguments: !
! FUNC (FUNCTION): !
! A function f: R -> R. !
! X (REAL): !
! The initial guess for the root of func. !
! OUTPUT: !
! ROOT (REAL): !
! A root for FUNC. !
!--------------------------------------------------------------------------!
FUNCTION STEFFENSENS_METHOD(FUNC, X)
IMPLICIT NONE
! FUNC is a function f: R -> R. Define this.
INTERFACE
FUNCTION FUNC(X)
IMPLICIT NONE
REAL, INTENT(IN) :: X
REAL :: FUNC
END FUNCTION FUNC
END INTERFACE
! The input is a real number, the guess for Steffensen's method.
REAL, INTENT(IN) :: X
! The output is also a positive real number, the root of FUNC.
REAL :: STEFFENSENS_METHOD
! Dummy variable for keeping track of how many iterations we've done.
INTEGER :: ITERS
! Steffensen's method is similar to Newton's method. At each step we
! Calculate the function F evaluated at the new approximation, and
! also compute a function G that acts as the derivative of F. Provide
! variables for these two evaluations.
REAL :: F_XN, G_XN
! Variable for the updated approximations.
REAL :: XN
! Steffensen's method starts the process at the guess point.
XN = X
! Iteratively perform Steffensen's method.
DO ITERS = 0, MAXIMUM_NUMBER_OF_ITERATIONS
! Steffensen's method needs the evaluations f(x) and f(x + f(x)),
! in particular the denominator is f(x + f(x)) / f(x) - 1.
F_XN = FUNC(XN)
G_XN = FUNC(XN + F_XN) / F_XN - 1.0
! Like Newton's method the new point is obtained by subtracting
! the ratio. g(x) = f(x + f(x))/f(x) - 1 acts as the derivative
! of f, but we do not explicitly need to calculate f'(x).
XN = XN - F_XN / G_XN
! If f(x) is very small, we are close to a root and can break
! out of this loop. Check for this.
IF (ABS(F_XN) .LE. EPSILON) THEN
EXIT
END IF
END DO
! Like Newton's method, and like Heron's method, the convergence is
! quadratic. After a few iterations we will be very to close a root.
STEFFENSENS_METHOD = XN
END FUNCTION STEFFENSENS_METHOD
! Provide the function f(x) = 2 - x^2. sqrt(2) is a root.
FUNCTION FUNC(X)
IMPLICIT NONE
REAL, INTENT(IN) :: X
REAL :: FUNC
FUNC = 2.0 - X * X
END FUNCTION FUNC
END MODULE STEFFENSEN
! Program for testing our implementation of Steffensen's method.
PROGRAM MAIN
USE STEFFENSEN
IMPLICIT NONE
! 2 is close enough to sqrt(2) for Steffensen's method to converge.
REAL :: X = 2.0
! Variable for the output.
REAL :: SQRT_2
! Run the routine, compute sqrt(2), and print the result.
SQRT_2 = STEFFENSENS_METHOD(FUNC, X)
PRINT "(A,F18.16)", "sqrt(2) = ", SQRT_2
END PROGRAM MAIN
Go
/******************************************************************************
* LICENSE *
******************************************************************************
* This file is part of mitx_mathematics_programming_examples. *
* *
* mitx_mathematics_programming_examples is free software: you can *
* redistribute it and/or modify it under the terms of the GNU General *
* Public License as published by the Free Software Foundation, either *
* version 3 of the License, or (at your option) any later version. *
* *
* mitx_mathematics_programming_examples is distributed in the hope that *
* it will be useful but WITHOUT ANY WARRANTY; without even the implied *
* warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. *
* See the GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with mitx_mathematics_programming_examples. If not, see *
* <https://www.gnu.org/licenses/>. *
******************************************************************************
* Purpose: *
* Calculates the root of a function using Steffensen's method. *
******************************************************************************
* Author: Ryan Maguire *
* Date: 2025/05/22 *
******************************************************************************/
package main
/* Only standard library imports are needed. */
import (
"fmt" /* Printf provided here, used for printing text to the screen. */
"math" /* Abs, floating-point absolute value function, found here. */
)
/* Type for a function of the form f: R -> R. */
type realfunc func(x float64) float64
/* Computes the root of a function using Steffensen's method. */
func steffensens_method(f realfunc, x float64) float64 {
/* Steffensen's method is iterative and converges very quickly. *
* Because of this we may exit the function after a few iterations. */
const maximum_number_of_iterations uint32 = 16
/* The maximum allowed error. This is 4x double precision epsilon. */
const epsilon float64 = 8.881784197001252E-16
/* Variable for keeping track of how many iterations we have performed. */
var iters uint32
/* The method starts at the provided guess point and updates iteratively.*/
var xn float64 = x
/* Iteratively apply Steffensen's method to find the root. */
for iters = 0; iters < maximum_number_of_iterations; iters += 1 {
/* Steffensen's method needs the evaluations f(x) and f(x + f(x)), *
* in particular the denominator is f(x + f(x)) / f(x) - 1. Compute. */
var f_xn float64 = f(xn)
var g_xn float64 = f(xn + f_xn) / f_xn - 1.0
/* Like Newton's method, the new point is obtained by subtracting *
* the ratio. g(x) = f(x + f(x)) / f(x) - 1 acts as the derivative *
* of f, but we do not explicitly need to calculate f'(x). */
xn = xn - f_xn / g_xn
/* If f(x) is very small, we are close to a root and can break out *
* of this for loop. Check for this. */
if math.Abs(f_xn) < epsilon {
break
}
}
/* Like Newton's method, and like Heron's method, the convergence is *
* quadratic. After a few iterations we will be very to close a root. */
return xn
}
/* End of steffensens_method. */
/* sqrt(2) is a root to the function f(x) = 2 - x^2. Provide this. */
func f(x float64) float64 {
return 2.0 - x*x
}
/* End of func. */
/* Main routine used for testing our implementation of Steffensen's method. */
func main() {
/* The initial guess point for Steffensen's method. */
const x float64 = 2.0
/* Calculate the square root and print it to the screen. If we have *
* written things correctly, we should get 1.414..., which is sqrt(2). */
var sqrt_x float64 = steffensens_method(f, x)
fmt.Printf("sqrt(%.1f) = %.16f\n", x, sqrt_x)
}
Java
/******************************************************************************
* LICENSE *
******************************************************************************
* This file is part of mitx_mathematics_programming_examples. *
* *
* mitx_mathematics_programming_examples is free software: you can *
* redistribute it and/or modify it under the terms of the GNU General *
* Public License as published by the Free Software Foundation, either *
* version 3 of the License, or (at your option) any later version. *
* *
* mitx_mathematics_programming_examples is distributed in the hope that *
* it will be useful but WITHOUT ANY WARRANTY; without even the implied *
* warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. *
* See the GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with mitx_mathematics_programming_examples. If not, see *
* <https://www.gnu.org/licenses/>. *
******************************************************************************
* Purpose: *
* Calculates the root of a function using Steffensen's method. *
******************************************************************************
* Author: Ryan Maguire *
* Date: 2025/05/25 *
******************************************************************************/
/* Class providing an implementation of Steffensen's method. */
final public class Steffensen {
/* Interface used for defining functions of the type f: R -> R. */
interface Function
{
double evaluate(double x);
}
/* Steffensen's method is iterative and converges very quickly. *
* Because of this we may exit the function after a few iterations. */
static public final int maximumNumberOfIterations = 16;
/* The maximum allowed error. This is 4x double precision epsilon. */
static public final double epsilon = 8.881784197001252E-16;
/* Computes the root of a function using Steffensen's method. */
static double root(Function f, double x)
{
/* Variable for keeping track of how many iterations we perform. */
int iters;
/* The method starts at the guess point and updates iteratively. */
double xn = x;
/* Iteratively apply Steffensen's method to find the root. */
for (iters = 0; iters < maximumNumberOfIterations; ++iters)
{
/* Steffensen's method needs the evaluations f(x) and f(x+f(x)), *
* in particular the denominator is f(x + f(x)) / f(x) - 1. */
double f_xn = f.evaluate(xn);
double g_xn = f.evaluate(xn + f_xn) / f_xn - 1.0;
/* Like Newton's method the new point is obtained by subtracting *
* the ratio. g(x) = f(x + f(x))/f(x) - 1 acts as the derivative *
* of f, but we do not explicitly need to calculate f'(x). */
xn = xn - f_xn / g_xn;
/* If f(x) is very small, we are close to a root and can break *
* out of this for loop. Check for this. */
if (Math.abs(f_xn) < epsilon)
break;
}
/* Like Newton's method, and like Heron's method, the convergence is *
* quadratic. After a few iterations we will be very to close a root.*/
return xn;
}
/* End of root. */
/* Main routine used for testing our Steffensen implementation. */
static public void main(String[] args)
{
/* Our guess point for sqrt(2). 2 is close enough for Steffensen's *
* method to work. */
final double x = 2.0;
/* Create a function that evaluates 2 - x^2. sqrt(2) is a root. */
Function func = new Function()
{
public double evaluate(double x)
{
return 2.0 - x*x;
}
};
/* Calculate the square root and print it to the screen. If we have *
* written things correctly we should get 1.414..., which is sqrt(2).*/
double sqrt_2 = Steffensen.root(func, x);
System.out.printf("pi = %.16f\n", sqrt_2);
}
}
Julia
################################################################################
# LICENSE #
################################################################################
# This file is part of mitx_mathematics_programming_examples. #
# #
# mitx_mathematics_programming_examples is free software: you can #
# redistribute it and/or modify it under the terms of the GNU General Public #
# License as published by the Free Software Foundation, either version 3 of #
# the License, or (at your option) any later version. #
# #
# mitx_mathematics_programming_examples is distributed in the hope that it #
# will be useful, but WITHOUT ANY WARRANTY; without even the implied #
# warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the #
# GNU General Public License for more details. #
# #
# You should have received a copy of the GNU General Public License #
# along with mitx_mathematics_programming_examples. If not, see #
# <https://www.gnu.org/licenses/>. #
################################################################################
# Purpose: #
# Calculates roots of a function using Steffensen's method. #
################################################################################
# Author: Ryan Maguire #
# Date: May 25, 2025. #
################################################################################
################################################################################
# Function: #
# steffensen #
# Purpose: #
# Computes the root of a function f: R -> R. #
# Arguments: #
# f (function): #
# A function f: R -> R. The root of f is computed. #
# x (float): #
# The guess point for the root of f. #
# Output: #
# root (float): #
# A root of the function f near x. #
################################################################################
function steffensen(f, x)
# Steffensen's method is iterative. The convergence is quadratic, meaning
# the number of accurate decimals doubles with each iteration. Because of
# this we can halt the algorithm after a few steps.
maximum_number_of_iterations = 16
# Maximum allowed error. This is 4x double precision epsilon.
epsilon = 8.881784197001252E-16
# Steffensen's method starts the process at the guess point.
xn = x
# Iteratively perform Steffensen's method.
for _ = 1:maximum_number_of_iterations
# Steffensen's method needs the evaluations f(x) and f(x + f(x)),
# in particular the denominator is f(x + f(x)) / f(x) - 1.
f_xn = f(xn)
g_xn = f(xn + f_xn) / f_xn - 1.0
# Like Newton's method the new point is obtained by subtracting
# the ratio. g(x) = f(x + f(x))/f(x) - 1 acts as the derivative
# of f, but we do not explicitly need to calculate f'(x).
xn = xn - f_xn / g_xn
# If f(x) is very small, we are close to a root and can break
# out of this for loop. Check for this.
if abs(f_xn) < epsilon
break
end
end
# Like Newton's method, and like Heron's method, the convergence is
# quadratic. After a few iterations we will be very to close a root.
return xn
end
# Input function for Steffensen's method.
function func(x)
"""
Computes f(x) = 2 - x^2. sqrt(2) is a root of this, we can compute
this value by applying Steffensen's method to this function.
"""
return 2.0 - x*x
end
# 2 is close enough to sqrt(2) for Steffensen's method to work.
X = 2.0
# Compute sqrt(2) using Steffensen's method and print the result.
sqrt_2 = steffensen(func, X)
println("sqrt(2) = ", sqrt_2)
JavaScript
/******************************************************************************
* LICENSE *
******************************************************************************
* This file is part of mitx_mathematics_programming_examples. *
* *
* mitx_mathematics_programming_examples is free software: you can *
* redistribute it and/or modify it under the terms of the GNU General *
* Public License as published by the Free Software Foundation, either *
* version 3 of the License, or (at your option) any later version. *
* *
* mitx_mathematics_programming_examples is distributed in the hope that *
* it will be useful but WITHOUT ANY WARRANTY; without even the implied *
* warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. *
* See the GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with mitx_mathematics_programming_examples. If not, see *
* <https://www.gnu.org/licenses/>. *
******************************************************************************
* Purpose: *
* Calculates the root of a function using Steffensen's method. *
******************************************************************************
* Author: Ryan Maguire *
* Date: 2025/05/22 *
******************************************************************************/
/* Computes the root of a function using Steffensen's method. */
function steffensensMethod(f, x) {
/* Steffensen's method is iterative and converges very quickly. *
* Because of this we may exit the function after a few iterations. */
const MAXIMUM_NUMBER_OF_ITERATIONS = 16;
/* The maximum allowed error. This is 4x double precision epsilon. */
const EPSILON = 8.881784197001252E-16;
/* The method starts at the provided guess point and updates iteratively.*/
let xn = x;
/* Variable for keeping track of how many iterations we perform. */
let iters;
/* Iteratively apply Steffensen's method to find the root. */
for (iters = 0; iters < MAXIMUM_NUMBER_OF_ITERATIONS; ++iters) {
/* Steffensen's method needs the evaluations f(x) and f(x + f(x)), *
* in particular the denominator is f(x + f(x)) / f(x) - 1. Compute. */
const fXn = f(xn);
const gXn = f(xn + fXn) / fXn - 1.0;
/* Like Newton's method, the new point is obtained by subtracting *
* the ratio. g(x) = f(x + f(x)) / f(x) - 1 acts as the derivative *
* of f, but we do not explicitly need to calculate f'(x). */
xn = xn - fXn / gXn;
/* If f(x) is very small, we are close to a root and can break out *
* of this for loop. Check for this. */
if (Math.abs(fXn) < EPSILON) {
break;
}
}
/* Like Newton's method, and like Heron's method, the convergence is *
* quadratic. After a few iterations we will be very to close a root. */
return xn;
}
/* End of steffensensMethod. */
/* sqrt(2) is a root to the function f(x) = 2 - x^2. Provide this. */
function f(x) {
return 2.0 - x*x;
}
/* The initial guess point for Steffensen's method. */
let x = 2.0;
/* Calculate the square root and print it to the screen. If we have *
* written things correctly, we should get 1.414..., which is sqrt(2). */
let sqrtX = steffensensMethod(f, x);
console.log("sqrt(" + x.toFixed(1) + ") = " + sqrtX.toFixed(16));
Kotlin
/******************************************************************************
* LICENSE *
******************************************************************************
* This file is part of mitx_mathematics_programming_examples. *
* *
* mitx_mathematics_programming_examples is free software: you can *
* redistribute it and/or modify it under the terms of the GNU General *
* Public License as published by the Free Software Foundation, either *
* version 3 of the License, or (at your option) any later version. *
* *
* mitx_mathematics_programming_examples is distributed in the hope that *
* it will be useful but WITHOUT ANY WARRANTY; without even the implied *
* warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. *
* See the GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with mitx_mathematics_programming_examples. If not, see *
* <https://www.gnu.org/licenses/>. *
******************************************************************************
* Purpose: *
* Calculates the root of a function using Steffensen's method. *
******************************************************************************
* Author: Ryan Maguire *
* Date: 2025/05/25 *
******************************************************************************/
import kotlin.math.abs
import kotlin.math.sin
/* Type for a function of the form f: R -> R. */
typealias RealFunc = (Double) -> Double
/* Computes the root of a function using Steffensen's method. */
fun steffensensMethod(f: RealFunc, x: Double): Double {
/* Steffensen's method is iterative and converges very quickly. *
* Because of this we may exit the function after a few iterations. */
val maximum_number_of_iterations: Int = 16
/* The maximum allowed error. This is 4x double precision epsilon. */
val epsilon: Double = 8.881784197001252E-16
/* The method starts at the provided guess point and updates iteratively.*/
var xn: Double = x
/* Iteratively apply Steffensen's method to find the root. */
for (iters in 0 .. maximum_number_of_iterations) {
/* Steffensen's method needs the evaluations f(x) and f(x + f(x)), *
* in particular the denominator is f(x + f(x)) / f(x) - 1. Compute. */
val f_xn: Double = f(xn)
val g_xn: Double = f(xn + f_xn) / f_xn - 1.0
/* Like Newton's method, the new point is obtained by subtracting *
* the ratio. g(x) = f(x + f(x)) / f(x) - 1 acts as the derivative *
* of f, but we do not explicitly need to calculate f'(x). */
xn = xn - f_xn / g_xn
/* If f(x) is very small, we are close to a root and can break out *
* of this for loop. Check for this. */
if (abs(f_xn) < epsilon) {
break
}
}
/* Like Newton's method, and like Heron's method, the convergence is *
* quadratic. After a few iterations we will be very to close a root. */
return xn
}
/* End of steffensensMethod. */
/* sqrt(2) is a root to the function f(x) = 2 - x^2. Provide this. */
fun f(x: Double): Double {
return 2.0 - x*x
}
fun main() {
/* The initial guess point for Steffensen's method. */
val x: Double = 2.0
/* Calculate the square root and print it to the screen. If we have *
* written things correctly, we should get 1.414..., which is sqrt(2). */
val sqrt_x: Double = steffensensMethod(::f, x)
println("sqrt(2) = $sqrt_x")
}
MATLAB
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% LICENSE %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% This file is part of mitx_mathematics_programming_examples. %
% %
% mitx_mathematics_programming_examples is free software: you can %
% redistribute it and/or modify it under the terms of the GNU General Public %
% License as published by the Free Software Foundation, either version 3 of %
% the License, or (at your option) any later version. %
% %
% mitx_mathematics_programming_examples is distributed in the hope that it %
% will be useful, but WITHOUT ANY WARRANTY; without even the implied %
% warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the %
% GNU General Public License for more details. %
% %
% You should have received a copy of the GNU General Public License %
% along with mitx_mathematics_programming_examples. If not, see %
% <https://www.gnu.org/licenses/>. %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Purpose: %
% Calculates roots of a function using Steffensen's method. %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Author: Ryan Maguire %
% Date: May 25, 2025. %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Function: %
% steffensen %
% Purpose: %
% Computes the root of a function f: R -> R. %
% Arguments: %
% f (function): %
% A function f: R -> R. The root of f is computed. %
% x (float): %
% The guess point for the root of f. %
% Output: %
% root (float): %
% A root of the function f near x. %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function root = steffensen(f, x)
% Steffensen's method is iterative. The convergence is quadratic, meaning
% the number of accurate decimals doubles with each iteration. Because of
% this we can halt the algorithm after a few steps.
maximum_number_of_iterations = 16;
% Maximum allowed error. This is 4x double precision epsilon.
epsilon = 8.881784197001252E-16;
% Steffensen's method starts the process at the guess point.
xn = x;
% Iteratively perform Steffensen's method.
for _ = 1:maximum_number_of_iterations
% Steffensen's method needs the evaluations f(x) and f(x + f(x)),
% in particular the denominator is f(x + f(x)) / f(x) - 1.
f_xn = f(xn);
g_xn = f(xn + f_xn) / f_xn - 1.0;
% Like Newton's method the new point is obtained by subtracting
% the ratio. g(x) = f(x + f(x))/f(x) - 1 acts as the derivative
% of f, but we do not explicitly need to calculate f'(x).
xn = xn - f_xn / g_xn;
% If f(x) is very small, we are close to a root and can break
% out of this for loop. Check for this.
if abs(f_xn) < epsilon;
root = xn;
return;
end
end
% Like Newton's method, and like Heron's method, the convergence is
% quadratic. After a few iterations we will be very to close a root.
root = xn;
end
% Input function for Steffensen's method.
function out = func(x)
out = 2.0 - x*x;
end
% 2 is close enough to sqrt(2) for Steffensen's method to work.
x = 2.0;
% Compute sqrt(2) using Steffensen's method.
sqrt_2 = steffensen(@func, x);
printf("sqrt(2) = %.16f\n", sqrt_2);
Objective-C
/******************************************************************************
* LICENSE *
******************************************************************************
* This file is part of mitx_mathematics_programming_examples. *
* *
* mitx_mathematics_programming_examples is free software: you can *
* redistribute it and/or modify it under the terms of the GNU General *
* Public License as published by the Free Software Foundation, either *
* version 3 of the License, or (at your option) any later version. *
* *
* mitx_mathematics_programming_examples is distributed in the hope that *
* it will be useful but WITHOUT ANY WARRANTY; without even the implied *
* warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. *
* See the GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with mitx_mathematics_programming_examples. If not, see *
* <https://www.gnu.org/licenses/>. *
******************************************************************************
* Purpose: *
* Calculates the root of a function using Steffensen's method. *
******************************************************************************
* Author: Ryan Maguire *
* Date: 2025/05/25 *
******************************************************************************/
/* NSObject, base class for all Objective-C classes, found here. math.h and *
* stdio.h are also included in this header file. */
#import <Foundation/Foundation.h>
/* Function pointer notation is a little confusing. Create a typedef for it *
* so we do not need to explicitly use it later. */
typedef double (*function)(double);
/* Class providing an implementation of Steffensen's method. */
@interface Steffensen: NSObject
{
/* How many iterations of Steffensen's method we do before exiting. */
const unsigned int maximumNumberOfIterations;
/* Allowed tolerance for the root-finding method. */
const double epsilon;
}
/* Given a function f: R -> R and a guess point x, this finds a root of *
* f that is near the guess point. */
+ (double) root: (function)f withInitialGuess: (double)x;
@end
@implementation Steffensen
/* Steffensen's method is iterative and converges very quickly. *
* Because of this we may exit the function after a few iterations. */
static const unsigned int maximum_number_of_iterations = 16U;
/* The maximum allowed error. This is 4x double precision epsilon. */
static const double epsilon = 8.881784197001252E-16;
/* Computes the root of a function using Steffensen's method. */
+ (double) root: (function) f withInitialGuess: (double) x
{
/* Variable for keeping track of how many iterations we perform. */
unsigned int iters;
/* The method starts at the guess point and updates iteratively. */
double xn = x;
/* Iteratively apply Steffensen's method to find the root. */
for (iters = 0; iters < maximum_number_of_iterations; ++iters)
{
/* Steffensen's method needs the evaluations f(x) and f(x+f(x)), *
* in particular the denominator is f(x + f(x)) / f(x) - 1. */
const double f_xn = f(xn);
const double g_xn = f(xn + f_xn) / f_xn - 1.0;
/* Like Newton's method the new point is obtained by subtracting *
* the ratio. g(x) = f(x + f(x))/f(x) - 1 acts as the derivative *
* of f, but we do not explicitly need to calculate f'(x). */
xn = xn - f_xn / g_xn;
/* If f(x) is very small, we are close to a root and can break *
* out of this for loop. Check for this. */
if (fabs(f_xn) < epsilon)
break;
}
/* Like Newton's method, and like Heron's method, the convergence is *
* quadratic. After a few iterations we will be very to close a root.*/
return xn;
}
/* End of root. */
@end
/* End of Steffensen implementation. */
/* sqrt(2) is a root to the function f(x) = 2 - x^2. Provide this. */
static double func(double x)
{
return 2.0 - x*x;
}
/* End of func. */
/* Main routine used for testing our implementation of Steffensen's method. */
int main(void)
{
/* The initial guess point for Steffensen's method. */
const double x = 2.0;
/* Calculate the square root and print it to the screen. If we have *
* written things correctly, we should get 1.414..., which is sqrt(2). */
const double sqrt_x = [Steffensen root: func withInitialGuess: x];
printf("sqrt(%.1f) = %.16f\n", x, sqrt_x);
return 0;
}
Pascal
(******************************************************************************
* LICENSE *
******************************************************************************
* This file is part of mitx_mathematics_programming_examples. *
* *
* mitx_mathematics_programming_examples is free software: you can *
* redistribute it and/or modify it under the terms of the GNU General Public*
* License as published by the Free Software Foundation, either version 3 of *
* the License, or (at your option) any later version. *
* *
* mitx_mathematics_programming_examples is distributed in the hope that it *
* will be useful but WITHOUT ANY WARRANTY; without even the implied warranty*
* of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with mitx_mathematics_programming_examples. If not, see *
* <https://www.gnu.org/licenses/>. *
******************************************************************************
* Purpose: *
* Computes roots using Steffensen's method. *
******************************************************************************
* Author: Ryan Maguire *
* Date: 2025/05/25 *
******************************************************************************)
PROGRAM Steffensen;
(* We define two constants: The allowed tolerance (or epsilon value), *
* and the maximum number of iterations we allow for Steffensen's method. *)
CONST
(* Steffensen's method is iterative and converges very quickly. *
* Because of this we may exit the function after a few iterations. *)
MaximumNumberOfIterations: Integer = 16;
(* The maximum allowed error. This is 4x double precision epsilon. *)
Epsilon: Real = 8.881784197001252E-16;
(* 2 is close enough to sqrt(2) for Steffensen's method to work. *)
Guess: Real = 2.0;
(* Create an alias for real functions f: R -> R. This makes the syntax of *
* the parameters for SteffensensMethod function a little more readable. *)
TYPE
RealFunc = Function(Const X: Real): Real;
(* The program has one variable, sqrt(2) computed via Steffensen's method. *)
VAR
Sqrt2: Real;
(* f(x) = 2 - x^2 has sqrt(2) as a root. This is the input function. *)
Function Func(const X: Real) : Real;
BEGIN
Func := 2.0 - X*X;
END;
(******************************************************************************
* Function: *
* SteffensensMethod *
* Purpose: *
* Computes roots using Steffensen's method. *
* Arguments: *
* F (RealFunc): *
* A continuous function f: R -> R. *
* X (Real): *
* The initial guess point for Steffensen's method. *
* OUTPUT: *
* Root (Real): *
* A root for F near X. *
******************************************************************************)
Function SteffensensMethod(Const F: RealFunc; Const X: Real) : Real;
VAR
(* Steffensen's method is iterative, like Newton's method, but avoids *
* explicit use of the derivative of F. Instead it uses the function *
* G(x) = F(x + F(x)) / F(X) - 1, which acts like F'(x). Provide a *
* variable for each of these expressions, Xn is the nth approximation. *)
Xn, Fxn, Gxn: Real;
(* Dummy variable for tracking how many iterations we've performed. *)
Iters: Integer;
BEGIN
(* Steffensen's method starts the iterative process at the guess point. *)
Xn := X;
(* Iteratively perform Steffensen's method. *)
FOR Iters := 0 TO MaximumNumberOfIterations - 1 DO
BEGIN
(* Steffensen's method needs the evaluations f(x) and f(x + f(x)), *
* in particular the denominator is f(x + f(x)) / f(x) - 1. Compute. *)
Fxn := F(Xn);
Gxn := F(Xn + Fxn) / Fxn - 1.0;
(* Like Newton's method, the new point is obtained by subtracting *
* the ratio. g(x) = f(x + f(x)) / f(x) - 1 acts as the derivative *
* of f, but we do not explicitly need to calculate f'(x). *)
Xn := Xn - Fxn / Gxn;
(* If f(x) is very small, we are close to a root and can break out *
* of this for loop. Check for this. *)
IF (ABS(Fxn) < Epsilon) THEN BREAK;
END;
(* Like Newton's method, and like Heron's method, the convergence is *
* quadratic. After a few iterations we will be very to close a root. *)
SteffensensMethod := Xn;
END;
(* Program for testing our implementation of Steffensen's method. *)
BEGIN
(* Compute sqrt(2) using Steffensen's method and print the result. *)
Sqrt2 := SteffensensMethod(@Func, Guess);
WriteLn('sqrt(2) = ', Sqrt2:0:16);
END.
IDL
;------------------------------------------------------------------------------;
; LICENSE ;
;------------------------------------------------------------------------------;
; This file is part of mitx_mathematics_programming_examples. ;
; ;
; mitx_mathematics_programming_examples is free software: you can ;
; redistribute it and/or modify it under the terms of the GNU General Public ;
; License as published by the Free Software Foundation, either version 3 of ;
; the License, or (at your option) any later version. ;
; ;
; mitx_mathematics_programming_examples is distributed in the hope that it ;
; will be useful but WITHOUT ANY WARRANTY; without even the implied warranty ;
; of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ;
; GNU General Public License for more details. ;
; ;
; You should have received a copy of the GNU General Public License ;
; along with mitx_mathematics_programming_examples. If not, see ;
; <https://www.gnu.org/licenses/>. ;
;------------------------------------------------------------------------------;
; Purpose: ;
; Calculates roots using Steffensen's method. ;
;------------------------------------------------------------------------------;
; Author: Ryan Maguire ;
; Date: 2025/05/22 ;
;------------------------------------------------------------------------------;
FUNCTION STEFFENSEN, F, X
; Tells the compiler that integers should be 32 bits, not 16.
COMPILE_OPT IDL2
; Error checking code.
ON_ERROR, 2
; Steffensen's method is iterative. The convergence is quadratic, meaning
; the number of accurate decimals doubles with each iteration. Because of
; this we can halt the algorithm after a few steps.
MAXIMUM_NUMBER_OF_ITERATIONS = 16
; Maximum allowed error. This is 4x double precision epsilon.
EPSILON = 8.881784197001252E-16
; Variable used for Steffensen's approximation. We will iteratively update
; this value with better approximations using Steffensen's method.
XN = X
; Iteratively perform Steffensen's method.
FOR ITERS = 1, MAXIMUM_NUMBER_OF_ITERATIONS DO BEGIN
; Steffensen's method needs the evaluations f(x) and f(x + f(x)),
; in particular the denominator is f(x + f(x)) / f(x) - 1. Compute.
F_XN = CALL_FUNCTION(F, XN)
G_XN = CALL_FUNCTION(F, XN + F_XN) / F_XN - 1.0
; Like Newton's method, the new point is obtained by subtracting
; the ratio. g(x) = f(x + f(x)) / f(x) - 1 acts as the derivative
; of f, but we do not explicitly need to calculate f'(x).
XN = XN - F_XN / G_XN
; If f(x) is very small, we are close to a root and can break out
; of this for loop. Check for this.
IF (ABS(F_XN) LE EPSILON) THEN BREAK
END
; Like Newton's method, and like Heron's method, the convergence is
; quadratic. After a few iterations we will be very to close a root.
RETURN, XN
END
; sqrt(2) is the root of the function 2 - x^2. Define this.
FUNCTION FUNC, X
COMPILE_OPT IDL2
RETURN, 2.0 - X*X
END
; Program for testing our implementation of Steffensen's method.
PRO MAIN
COMPILE_OPT IDL2
; The input for the method. We'll compute sqrt(2).
X = DOUBLE(2.0)
; Run the routine, computing sqrt(x), and print the result.
SQRT_2 = STEFFENSEN("FUNC", X)
PRINT, X, SQRT_2, FORMAT = 'SQRT(%3.1f) = %18.16f'
END
Python
"""
################################################################################
# LICENSE #
################################################################################
# This file is part of mitx_mathematics_programming_examples. #
# #
# mitx_mathematics_programming_examples is free software: you can #
# redistribute it and/or modify it under the terms of the GNU General Public #
# License as published by the Free Software Foundation, either version 3 of #
# the License, or (at your option) any later version. #
# #
# mitx_mathematics_programming_examples is distributed in the hope that it #
# will be useful, but WITHOUT ANY WARRANTY; without even the implied #
# warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the #
# GNU General Public License for more details. #
# #
# You should have received a copy of the GNU General Public License #
# along with mitx_mathematics_programming_examples. If not, see #
# <https://www.gnu.org/licenses/>. #
################################################################################
# Purpose: #
# Calculates roots of a function using Steffensen's method. #
################################################################################
# Author: Ryan Maguire #
# Date: May 25, 2025. #
################################################################################
"""
# Pylint doesn't like "x" as a variable name. Disable this warning.
# pylint: disable = invalid-name
# Computes the root of a function via Steffensen's method.
def steffensen(f, x):
"""
Function:
steffensen
Purpose:
Computes the root of a function f: R -> R.
Arguments:
f (function):
A function f: R -> R. The root of f is computed.
x (float):
The guess point for the root of f.
Output:
root (float):
A root of the function f near x.
"""
# Steffensen's method is iterative. The convergence is quadratic, meaning
# the number of accurate decimals doubles with each iteration. Because of
# this we can halt the algorithm after a few steps.
maximum_number_of_iterations = 16
# Maximum allowed error. This is 4x double precision epsilon.
epsilon = 8.881784197001252E-16
# Steffensen's method starts the process at the guess point.
xn = x
# Iteratively perform Steffensen's method.
for _ in range(maximum_number_of_iterations):
# Steffensen's method needs the evaluations f(x) and f(x + f(x)),
# in particular the denominator is f(x + f(x)) / f(x) - 1.
f_xn = f(xn)
g_xn = f(xn + f_xn) / f_xn - 1.0
# Like Newton's method the new point is obtained by subtracting
# the ratio. g(x) = f(x + f(x))/f(x) - 1 acts as the derivative
# of f, but we do not explicitly need to calculate f'(x).
xn = xn - f_xn / g_xn
# If f(x) is very small, we are close to a root and can break
# out of this for loop. Check for this.
if abs(f_xn) < epsilon:
break
# Like Newton's method, and like Heron's method, the convergence is
# quadratic. After a few iterations we will be very to close a root.
return xn
# Input function for Steffensen's method.
def func(x):
"""
Computes f(x) = 2 - x^2. sqrt(2) is a root of this, we can compute
this value by applying Steffensen's method to this function.
"""
return 2.0 - x*x
# 2 is close enough to sqrt(2) for Steffensen's method to work.
X = 2.0
# Compute sqrt(2) using Steffensen's method and print the result.
sqrt_2 = steffensen(func, X)
print(f'sqrt(2) = {sqrt_2}')
R
################################################################################
# LICENSE #
################################################################################
# This file is part of mitx_mathematics_programming_examples. #
# #
# mitx_mathematics_programming_examples is free software: you can #
# redistribute it and/or modify it under the terms of the GNU General Public #
# License as published by the Free Software Foundation, either version 3 of #
# the License, or (at your option) any later version. #
# #
# mitx_mathematics_programming_examples is distributed in the hope that it #
# will be useful, but WITHOUT ANY WARRANTY; without even the implied #
# warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the #
# GNU General Public License for more details. #
# #
# You should have received a copy of the GNU General Public License #
# along with mitx_mathematics_programming_examples. If not, see #
# <https://www.gnu.org/licenses/>. #
################################################################################
# Purpose: #
# Calculates roots of a function using Steffensen's method. #
################################################################################
# Author: Ryan Maguire #
# Date: May 25, 2025. #
################################################################################
# Computes the root of a function using Steffensen's method.
steffensen <- function(f, x) {
# Steffensen's method is iterative and converges very quickly.
# Because of this we may exit the function after a few iterations.
maximum_number_of_iterations <- 64
# The maximum allowed error. This is 4x double precision epsilon.
epsilon <- 8.881784197001252E-16
# The method starts at the provided guess point and updates iteratively.
xn <- x
# Iteratively apply Steffensen's method to find the root.
for (iters in 1:maximum_number_of_iterations) {
# Steffensen's method needs the evaluations f(x) and f(x + f(x)),
# in particular the denominator is f(x + f(x)) / f(x) - 1. Compute.
f_xn <- f(xn)
g_xn <- f(xn + f_xn) / f_xn - 1.0
# Like Newton's method, the new point is obtained by subtracting
# the ratio. g(x) = f(x + f(x)) / f(x) - 1 acts as the derivative
# of f, but we do not explicitly need to calculate f'(x).
xn <- xn - f_xn / g_xn
# If f(x) is very small, we are close to a root and can break out
# of this for loop. Check for this.
if (abs(f_xn) < epsilon) {
break
}
}
# Like Newton's method, and like Heron's method, the convergence is
# quadratic. After a few iterations we will be very to close a root.
return(xn)
}
# End of steffensen.
# sqrt(2) is a root to the function f(x) = 2 - x^2.
f <- function(x) {
return(2.0 - x*x)
}
# The initial guess point for Steffensen's method.
x <- 2.0
# Calculate the square root and print it to the screen. If we have
# written things correctly, we should get 1.414..., which is sqrt(2).
sqrt_x <- steffensen(f, x)
output = sprintf("sqrt(2) = %.16f", sqrt_x)
message(output)
Rust
/******************************************************************************
* LICENSE *
******************************************************************************
* This file is part of mitx_mathematics_programming_examples. *
* *
* mitx_mathematics_programming_examples is free software: you can *
* redistribute it and/or modify it under the terms of the GNU General *
* Public License as published by the Free Software Foundation, either *
* version 3 of the License, or (at your option) any later version. *
* *
* mitx_mathematics_programming_examples is distributed in the hope that *
* it will be useful but WITHOUT ANY WARRANTY; without even the implied *
* warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. *
* See the GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with mitx_mathematics_programming_examples. If not, see *
* <https://www.gnu.org/licenses/>. *
******************************************************************************
* Purpose: *
* Calculates the root of a function using Steffensen's method. *
******************************************************************************
* Author: Ryan Maguire *
* Date: 2025/05/22 *
******************************************************************************/
/* Type for a function of the form f: R -> R. */
type RealFunc = fn(f64) -> f64;
/* Computes the root of a function using Steffensen's method. */
fn steffensens_method(f: RealFunc, x: f64) -> f64 {
/* Steffensen's method is iterative and converges very quickly. *
* Because of this we may exit the function after a few iterations. */
const MAXIMUM_NUMBER_OF_ITERATIONS: u32 = 16;
/* The maximum allowed error. This is 4x double precision epsilon. */
const EPSILON: f64 = 8.881784197001252E-16;
/* The method starts at the provided guess point and updates iteratively.*/
let mut xn: f64 = x;
/* Iteratively apply Steffensen's method to find the root. */
for _ in 0 .. MAXIMUM_NUMBER_OF_ITERATIONS {
/* Steffensen's method needs the evaluations f(x) and f(x + f(x)), *
* in particular the denominator is f(x + f(x)) / f(x) - 1. Compute. */
let f_xn: f64 = f(xn);
let g_xn: f64 = f(xn + f_xn) / f_xn - 1.0;
/* Like Newton's method, the new point is obtained by subtracting *
* the ratio. g(x) = f(x + f(x)) / f(x) - 1 acts as the derivative *
* of f, but we do not explicitly need to calculate f'(x). */
xn = xn - f_xn / g_xn;
/* If f(x) is very small, we are close to a root and can break out *
* of this for loop. Check for this. */
if f_xn.abs() < EPSILON {
break;
}
}
/* Like Newton's method, and like Heron's method, the convergence is *
* quadratic. After a few iterations we will be very to close a root. */
return xn;
}
/* End of steffensens_method. */
/* sqrt(2) is a root to the function f(x) = 2 - x^2. Provide this. */
fn f(x: f64) -> f64 {
return 2.0 - x*x;
}
/* Main routine used for testing our implementation of Steffensen's method. */
fn main() {
/* The initial guess point for Steffensen's method. */
const X: f64 = 2.0;
/* Calculate the square root and print it to the screen. If we have *
* written things correctly, we should get 1.414..., which is sqrt(2). */
let sqrt_x: f64 = steffensens_method(f, X);
println!("sqrt({}) = {}", X, sqrt_x);
}
Swift
/******************************************************************************
* LICENSE *
******************************************************************************
* This file is part of mitx_mathematics_programming_examples. *
* *
* mitx_mathematics_programming_examples is free software: you can *
* redistribute it and/or modify it under the terms of the GNU General *
* Public License as published by the Free Software Foundation, either *
* version 3 of the License, or (at your option) any later version. *
* *
* mitx_mathematics_programming_examples is distributed in the hope that *
* it will be useful but WITHOUT ANY WARRANTY; without even the implied *
* warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. *
* See the GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with mitx_mathematics_programming_examples. If not, see *
* <https://www.gnu.org/licenses/>. *
******************************************************************************
* Purpose: *
* Calculates the root of a function using Steffensen's method. *
******************************************************************************
* Author: Ryan Maguire *
* Date: 2025/05/22 *
******************************************************************************/
/* fabs (floating-point absolute value) found here. */
import Foundation
/* Type for a function of the form f: R -> R. */
typealias RealFunc = (Double) -> Double
/* Computes the root of a function using Steffensen's method. */
func steffensensMethod(f: RealFunc, x: Double) -> Double {
/* Steffensen's method is iterative and converges very quickly. *
* Because of this we may exit the function after a few iterations. */
let maximum_number_of_iterations: UInt32 = 16
/* The maximum allowed error. This is 4x double precision epsilon. */
let epsilon: Double = 8.881784197001252E-16
/* The method starts at the provided guess point and updates iteratively.*/
var xn: Double = x
/* Iteratively apply Steffensen's method to find the root. */
for _ in 0 ..< maximum_number_of_iterations {
/* Steffensen's method needs the evaluations f(x) and f(x + f(x)), *
* in particular the denominator is f(x + f(x)) / f(x) - 1. Compute. */
let f_xn: Double = f(xn)
let g_xn: Double = f(xn + f_xn) / f_xn - 1.0
/* Like Newton's method, the new point is obtained by subtracting *
* the ratio. g(x) = f(x + f(x)) / f(x) - 1 acts as the derivative *
* of f, but we do not explicitly need to calculate f'(x). */
xn = xn - f_xn / g_xn
/* If f(x) is very small, we are close to a root and can break out *
* of this for loop. Check for this. */
if fabs(f_xn) < epsilon {
break
}
}
/* Like Newton's method, and like Heron's method, the convergence is *
* quadratic. After a few iterations we will be very to close a root. */
return xn
}
/* End of steffensensMethod. */
/* sqrt(2) is a root to the function f(x) = 2 - x^2. Provide this. */
func f(x: Double) -> Double {
return 2.0 - x*x
}
/* The initial guess point for Steffensen's method. */
let x: Double = 2.0
/* Calculate the square root and print it to the screen. If we have *
* written things correctly, we should get 1.414..., which is sqrt(2). */
let sqrt_x: Double = steffensensMethod(f: f, x: x)
print("sqrt(\(x)) = \(sqrt_x)")
TypeScript
/******************************************************************************
* LICENSE *
******************************************************************************
* This file is part of mitx_mathematics_programming_examples. *
* *
* mitx_mathematics_programming_examples is free software: you can *
* redistribute it and/or modify it under the terms of the GNU General *
* Public License as published by the Free Software Foundation, either *
* version 3 of the License, or (at your option) any later version. *
* *
* mitx_mathematics_programming_examples is distributed in the hope that *
* it will be useful but WITHOUT ANY WARRANTY; without even the implied *
* warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. *
* See the GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with mitx_mathematics_programming_examples. If not, see *
* <https://www.gnu.org/licenses/>. *
******************************************************************************
* Purpose: *
* Calculates the root of a function using Steffensen's method. *
******************************************************************************
* Author: Ryan Maguire *
* Date: 2025/05/22 *
******************************************************************************/
/* Type for a function of the form f: R -> R. */
type realfunc = (x: number) => number;
/* Computes the root of a function using Steffensen's method. */
function steffensensMethod(f: realfunc, x: number): number {
/* Steffensen's method is iterative and converges very quickly. *
* Because of this we may exit the function after a few iterations. */
const MAXIMUM_NUMBER_OF_ITERATIONS: number = 16;
/* The maximum allowed error. This is 4x double precision epsilon. */
const EPSILON: number = 8.881784197001252E-16;
/* The method starts at the provided guess point and updates iteratively.*/
let xn: number = x;
/* Variable for keeping track of how many iterations we perform. */
let iters: number;
/* Iteratively apply Steffensen's method to find the root. */
for (iters = 0; iters < MAXIMUM_NUMBER_OF_ITERATIONS; ++iters) {
/* Steffensen's method needs the evaluations f(x) and f(x + f(x)), *
* in particular the denominator is f(x + f(x)) / f(x) - 1. Compute. */
const fXn: number = f(xn);
const gXn: number = f(xn + fXn) / fXn - 1.0;
/* Like Newton's method, the new point is obtained by subtracting *
* the ratio. g(x) = f(x + f(x)) / f(x) - 1 acts as the derivative *
* of f, but we do not explicitly need to calculate f'(x). */
xn = xn - fXn / gXn;
/* If f(x) is very small, we are close to a root and can break out *
* of this for loop. Check for this. */
if (Math.abs(fXn) < EPSILON) {
break;
}
}
/* Like Newton's method, and like Heron's method, the convergence is *
* quadratic. After a few iterations we will be very to close a root. */
return xn;
}
/* End of steffensensMethod. */
/* sqrt(2) is a root to the function f(x) = 2 - x^2. Provide this. */
function f(x: number): number {
return 2.0 - x*x;
}
/* The initial guess point for Steffensen's method. */
let x: number = 2.0;
/* Calculate the square root and print it to the screen. If we have *
* written things correctly, we should get 1.414..., which is sqrt(2). */
let sqrtX: number = steffensensMethod(f, x);
console.log("sqrt(" + x.toFixed(1) + ") = " + sqrtX.toFixed(16));