Newton Method

From CodeCodex

Newthon-Raphson Method In numerical analysis, Newton's method (also known as the Newton–Raphson method), named after Isaac Newton and Joseph Raphson, is a method for finding successively better approximations to the zeroes (or roots) of a real-valued function. The algorithm is first in the class of Householder's methods, succeeded by Halley's method. The Newton-Raphson method in one variable: Given a function ƒ(x) and its derivative ƒ '(x), we begin with a first guess x0. Provided the function is reasonably well-behaved a better approximation x1 is

Geometrically, x1 is the intersection point of the tangent line to the graph of f, with the x-axis. The process is repeated until a sufficiently accurate value is reached:

C[edit]

/*This program uses the Newton-Raphson method to determine roots of a function.
  The format is: newton <function name> <initial guess>
*/

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#define TOLERANCE 1E-6
#define MAX_ITER 12

double f(const char* funname, double x);
double fPrime(const char* funname, double x);
void printFunction();

/* Add global variables and helper functions here */

int main(int argc, char** argv)
{
	/* Add your implementation here */
	//Local Declaration
	double x, y, dy;
	int iterations;
	//Statements
	if(argc < 2) // Check if the input is in the correct format
	{
		printf("Usage: newton <poly1|sin|xsin|poly2|imaginary> <initial guess>\n");
		exit(1);
	}
	
	printFunction(argv[1]); // call print function

	x = atoi(argv[2]); // convert from string argument to double saves it as x
	y = f(argv[1], x); // find the function value of x 
	dy = fPrime(argv[1], x); // find the first derivative for x
	iterations = 0;
	do {
		printf("At iteration %d, x=%f, y=%f, y'=%f\n", iterations, x, y, dy);
		if(fabs(dy) <= 1e-6) // check if the absolute value of first derivative is less or equal to 0
		{
			printf("Error: at x=%f, f'(x)=%d\n", x, dy);
			exit(1);
		}
		if(iterations == 12)
		{
			printf("Error: after %d iterations, x=%f and f(x)=%f\n", iterations, x, y);
			exit(101);
		}
		x = x - (y / dy); // find the next x / xk+1
  		y = f(argv[1], x); // find the function value for xk+1
		dy = fPrime(argv[1], x); // find the first derivative for xk+1
  		iterations++;
	}while(fabs(y) > TOLERANCE && iterations <= MAX_ITER); //This process continues until x converges
	
	printf("At iteration %d, x=%f, y=%f, and y'=%f\n", iterations, x, y, dy);
	printf("Solution: iteration=%d x=%f y=%f\n", iterations, x, y);
		
	return 0;
}

/* Prints the function in readable form */
void printFunction(const char * funname)
{
        char * func;

	if(strcmp(funname, "poly1") == 0)
	  func = "y = x^2 - 4 = 0";
	else if(strcmp(funname, "sin") == 0)
	  func = "y = sin(x)-.5 = 0";
	else if(strcmp(funname, "xsin") == 0)
	  func = "y = x*sin(x)-10 = 0";
	else if(strcmp(funname, "poly2") == 0)
	  func = "y = x^3+3*x^2+4*x-1 = 0";
	else if(strcmp(funname, "imaginary") == 0)
	  func = "y = x^2+1 = 0";
	else
	{
	  printf("Error: %s is not a valid function name\n", funname);
		exit(1);
	}

        printf("Function: %s\n", func);
}

/* Evaluates f(x) */
double f(const char* funname, double x)
{
	if(strcmp(funname, "poly1") == 0)
		return (x*x - 4);
	else if(strcmp(funname, "sin") == 0)
	  return (sin(x)-.5);
	else if(strcmp(funname, "xsin") == 0)
	  return (x*sin(x)-10);
	else if(strcmp(funname, "poly2") == 0)
	  return (x*x*x+3*x*x+4*x-1);
	else if(strcmp(funname, "imaginary") == 0)
	  return (x*x+1);
	else
	{
	  printf("Error: %s is not a valid function name\n", funname);
		exit(1);
	}
}

/* Evaluates f'(x) */
double fPrime(const char* funname, double x)
{
	if(strcmp(funname, "poly1") == 0)
		return (2*x);
	else if(strcmp(funname, "sin") == 0)
	  return (cos(x));
	else if(strcmp(funname, "xsin") == 0)
	  return (sin(x)+x*cos(x));
	else if(strcmp(funname, "poly2") == 0)
	  return (3*x*x+6*x+4);
	else if(strcmp(funname, "imaginary") == 0)
	  return (2*x);
	else
	{
	  printf("Error: %s is not a valid function name\n", funname);
		exit(1);
	}
}

Tcl[edit]

# f  Procedure that returns the value the function at x
# f' Procedure that returns the derivative of the function at x 
package require math::calculus
proc f  x {return [expr {sin($x) - $x}]}
prof f' x {return [expr {cos($x) - 1}]}
set root [::math::calculus::newtonRaphson f f' 1.6]