COST (GBP)

0.00

Maths › Rootfinding ›

Muller

Calculates the zeros of a function using Muller's method.

Controller: CodeCogs

Contents

Interface
Muller
Page Comments

Interface

C++

Muller

doublemuller(	double	(`f`)(double)[function pointer]*
	double	`x0` `= 0`
	double	`d` `= 3`
	double	`eps` `= 1E-10`
	int	`maxit` `= 1000`	)

This method extends the idea of the secant method which works with a linear polynomial, to a quadratic polynomial. Given three previous estimates $\inline z^{(k - 2)}$ , $\inline z^{(k - 1)}$ and $\inline z^{(k)}$ , for an unknown root, a new value is computed by

$z^{(k + 1)} = z^{(k)} + h_k q_k$

(2)

where

$h_k = z^{(k)} - z^{(k - 1)} \qquad r_k = h_k / h_{k - 1} \qquad q_k = \frac { -2 C_k } { B_k \pm \sqrt{ B_k^2 - 4 A_k C_k } }$

(3)

and

$A_k = r_k P \left( z^{(k)} \right) - r_k (1 + r_k) P \left( z^{(k - 1)} \right) + r_k^2 P \left( z^{(k - 2)} \right)$

(4)

$B_k = (2 r_k + 1) P \left( z^{(k)} \right) - (1 + r_k)^2 P \left( z^{(k - 1)} \right) + r_k^2 P \left( z^{(k - 2)} \right)$

(5)

$C_k = (1 + r_k) P \left( z^{(k)} \right)$

(6)

The values $\inline q_k$ may yield too large changes for $\inline z^{(k)}$ which possibly leads to another root and causes slow convergence. This can be circumvented by allowing a fixed maximum increase of $\inline |q_k|$ from one iteration to another. Care must also be taken when computing $\inline P \left( z^{(k)} \right)$ which is necessary to compute $\inline A_k$ , $\inline B_k$ and $\inline C_k$ . If an estimate of $\inline \left| P \left( z^{(k)} \rigth) \right|$ indicates a value greater than the maximum possible number, we choose

$z^{(k + 1)} = z^{(k)} + \frac {h_k q_k} {2}$

(7)

in place of the original relation, and repeat this until no overflow occurs. The algorithm stops whenever the actual value $\inline \left| P \left( z^{(k)} \rigth) \right|$ is smaller than the smallest value $\inline \left| P(z_{min}) \right|$ until now and

$\left| \frac { z_{min} - z^{(k + 1)} } { z_{min} } \right| < \epsilon$

(8)

holds where $\inline \epsilon$ is some small number depending on the computer accuracy. To avoid a lot of iterations where the above condition fails, we allow only a fixed maximum number of iterations.

Convergence for this method is superlinear. However, it is one of those methods that will converge to both real and complex roots from a real initial approximation.

This algorithm finds the roots of the user-defined function f starting with an initial guess x0 and iterating the sequence above until either the accuracy eps is achieved or the maximum number of iterations maxit is exceeded. Another required parameter is the bound on the error of the initial guess d.

References:

Jean-Pierre Moreau's Home Page, http://perso.wanadoo.fr/jean-pierre.moreau/
F.R. Ruckdeschel, "BASIC Scientific Subroutines", Vol. II, BYTE/McGRAWW-HILL, 1981
Wankere R. Mekwi, "Iterative Methods for Roots of Polynomials", Exeter College, University of Oxford

Example 1

#include <codecogs/maths/rootfinding/muller.h>
 
#include <iostream>
#include <iomanip>
 
// user-defined function
double f(double x) {
    return sin(x);
}
 
int main() 
{
  double x = Maths::RootFinding::muller(f, 4);
 
  std::cout << "The calculated zero is X = " << std::setprecision(15) << x << std::endl;
  std::cout << "The associated ordinate value is Y = " << f(x) << std::endl;
  return 0;
}

Output:

The calculated zero is X = 3.14159265358979
The associated ordinate value is Y = 1.22460635382238e-016

Parameters

f	the user-defined function
x0	Default value = 0
d	Default value = 3
eps	Default value = 1E-10
maxit	Default value = 1000

Authors

Lucian Bentea (August 2005)

Source Code

Source code is available when you buy a Commercial licence.

Not a member, then Register with CodeCogs. Already a Member, then Login.