Help := proc();

print(`Newton, Secant, ModifiedNewton, Steffensen, LaInter`);

end:

# Section 2.3 Newton's Method
# Program from the book 64.

# Input: function f, the initial point start,
# tolerance TOL, number of iteration N0;
# Output: the estimation of the root of f with 
# error less than TOL or fail message

# Try: problem 1 section 2.3 in the book 
# Newton(x->x^2-6,1,10^(-5),2); 

# Another example, problem 22;
# Digits := 102;
# Newton(x->cos(x)-x,Pi/4,10^(-100),10);

Newton := proc(f,start,TOL,N0) local p0,p,i,fprime; 

	fprime := diff(f(x),x); 
	
	p0 := evalf(start);
	
	print(0,p0);
	i := 1;

	while i <= N0 do
		
		p := evalf(p0 - (f(p0)/subs(x = p0,fprime)));
		print(i,p);

		if abs(p-p0) < TOL then
			return(p);
		fi:
		i := i+1;
		p0 := p;
	od:

	print( `method fail after `, N0);
	return(FAIL);

end:

# Section 2.3 Secant Method
# Program from the book 68.

# Input: function f, two initial point ps0 and ps1,
# tolerance TOL, number of iteration N0;
# Output: the estimation of the root of f 
# with error less than TOL or fail message

# Try: Secant(x->x^2-6,1,2,10^(-5),10);

Secant := proc(f,ps0,ps1,TOL,N0) local p,p0,p1,q0,q1,i; 
	
	p0 := evalf(ps0);
	p1 := evalf(ps1);
	i := 2;
	q0 := f(p0);
	q1 := f(p1);

	print(0,p0);
	print(1,p1);
	
	while i <= N0 do
		
		p := evalf(p1 - q1*(p1-p0)/(q1-q0));
		print(i,p);

		if abs(p-p1) < TOL then
			return(p);
		fi:

		i := i+1;

		p0 := p1;
		q0 := q1;
		p1 := p;
		q1 := f(p);

	od:

	print( `method fail after `, N0);
	return(FAIL);

end:

#section 2.4 Modified Newton's Method

# Input: function f, two initial point start, 
# tolerance TOL, number of iteration N0;
# Output: the estimation of the root of f with 
# error less than TOL or fail message

# Try: problem 3a) section 2.4
# ModifiedNewton(x->x^2-2*x*exp(-x)+exp(-2*x),0.5,10^(-5),10);

ModifiedNewton := proc(f,start,TOL,N0) 

local p0,p,i,fprime,fprimeprime; 

	fprime := diff(f(x),x); 
	fprimeprime := diff(fprime,x);
	
	p0 := evalf(start);
	i := 1;
	print(i,p0);
	while i <= N0 do

		p := evalf(p0 - f(p0)*subs(x = p0,fprime)/
		           (subs(x = p0,fprime)^2-f(p0)*subs(x = 
p0,fprimeprime)));
		print(i,p);
		if abs(p-p0) < TOL then
			return(p);
		fi:
		i := i+1;
		p0 := p;
	od:

	print( `method fail after `, N0);
	return(FAIL);

end:

# Section 2.5 Steffensen's method

# Input function g , initial point start, 
# tolerance TOL, number of iteration N0;
# Output p_0_(n) if the error of p_0_(n) 
# and fixed point p less than TOL or fail 
# message

# Try problem 3 section 2.5 
# and check the answer with the back of the book.
# Steffensen(x-> cos(x-1),2,10^(-5),1);

Steffensen := proc(g,start,TOL,N0) local i, P, p;

i := 1;

P[1] := start;

while i <= N0 do

	P[2] := evalf(g(P[1]));
	P[3] := evalf(g(P[2]));

	p := P[1]-(P[2]-P[1])^2/(P[3]-2*P[2]+P[1]);
	
	print(i,p);

	if abs(p-P[1]) < TOL then
		return(p);
	fi:

	i := i+1;
	P[1] := p;
od:

print(`method fail after`, N0);
return(FAIL);

end:


#Section 3.1 Lagrange interpolating polynomial

# Input: the list of point x , the list of 
# point y and variable x;
# Output: Lagrange interpolating polynomial 
# of degree nops(x)-1

# Try: problem 1a) (with 3 points) section 3.1 
#      and check the answer with the back of the book
# LaInter([0,0.6,0.9],[cos(0),cos(0.6),cos(0.9)],x);

# Try: problem 5a) (with 4 points) section 3.1 
#      and check the answer with the back of the book
# LaInter([8.1,8.3,8.6,8.7],
#         [16.94410,17.56492,18.50515,18.82091],x);

# if you want to evaluate polynomial 
# at point x_0 (in problem 6)
# , use command subs
# Try: subs(x=2, x+5);

LaInter := proc(X,Y,x) local temp, i,j, L;

L := [0$nops(X)];

for i from 1 to nops(X) do
	temp := mul(x-X[j], j=1..i-1)*mul(x-X[j], j=i+1..nops(X));
	L[i] := temp/subs(x=X[i],temp);
od:

simplify(sum(Y[j]*L[j], j = 1..nops(Y)));

end: