Added new files for Numerical Analysis

Bisection.m

@@ -0,0 +1,35 @@
+function [xList,errorList] = Bisection (f, lowerValue, upperValue, allowError)
+%
+%Uses the bisection method to find the possible answers to the root of the function f
+%
+
+  pkg load symbolic;
+  warning('off','OctSymPy:sym:rationalapprox');
+	%Setting necesary values for the function
+	cnt = 1;
+	maxItterations = 50;
+	errorValue = 1;
+	currentValue = 0.0;
+
+
+	%If the lower and upper bounds are mixed up swap them
+	if(double(subs(f,lowerValue)) > 0)
+		[lowerValue,upperValue] = swap(lowerValue, upperValue);
+	end
+
+	%Loop until the error is within bounds or the Maximum number of iterations is reached
+	while((abs(errorValue) > allowError) && (cnt < maxItterations))
+		currentValue = (lowerValue + upperValue)/ 2;
+		errorValue = double(subs(f,currentValue));
+		%Replace the correct value with the new value
+		%if error == 0 then the value has been found
+		if(errorValue < 0)
+			lowerValue = currentValue;
+		else
+			upperValue = currentValue;
+		end
+		xList(cnt) = currentValue;
+		errorList(cnt) = errorValue;
+		++cnt;
+	end
+end

FalsePosition.m

@@ -0,0 +1,31 @@
+function [xList, errorList] = FalsePosition(f, p0, p1, errorAllowed)
+%
+%FalsePosition(f, p0, p1, errorAllowed)
+%This function finds the root of a function using the method of False Position
+%
+
+	pkg load symbolic;
+	warning('off','OctSymPy:sym:rationalapprox');
+	maxIt = 50;
+	cnt = 2;
+	q0 = double(subs(f,p0));
+	q1 = double(subs(f,p1));
+	p = 0;
+	q = 0;
+	currentError = errorAllowed + 1;
+	while((cnt <= maxIt) && (currentError >= errorAllowed))
+		p = p1 - (q1 * (p1 - p0))/(q1 - q0);
+		currentError = abs(p - p1);
+		%Add Values to lists
+		xList(end+1) = p;
+		errorList(end+1) = currentError;
+		++cnt;
+		q = double(subs(f,p));
+		if((q * q1) < 0)
+			p0 = p1;
+			q0 = q1;
+		end
+		p1 = p;
+		q1 = q;
+	end
+end

Mullers.m

@@ -0,0 +1,44 @@
+function [xList, functionValueList] = Mullers(f, p0, p1, p2, errorAllowed)
+%
+%Mullers(f, p0, p1, p2, errorAllowed)
+%This function finds the root of a function using Muller's Method
+%
+
+	pkg load symbolic;
+	warning('off','OctSymPy:sym:rationalapprox');
+
+	h1 = p1 - p0;
+	h2 = p2 - p1;
+	s1 = (double(subs(f,p1)) - double(subs(f,p0))) / h1;
+	s2 = (double(subs(f,p2)) - double(subs(f,p2))) / h2;
+	d = (s2 - s1) / (h2 + h1);
+	cnt = 2;
+	maxIt = 50;
+
+	while(cnt < maxIt)
+		b = s2 + h2 * d;
+		D = (b^2 - (4 * double(subs(f,p2)) * d))^(1/2);
+		if(abs(b - D) < abs(b + D))
+			E = b + D;
+		else
+			E = b - D;
+		end
+		h = (02 * double(subs(f,p2)))/E;
+		p = p2 + h;
+		xList(end+1) = p;
+		functionValueList(end+1) = double(subs(f,p));
+		if(abs(h) < errorAllowed)
+			return;
+		else
+			p0 = p1;
+			p1 = p2;
+			p2 = p;
+			h1 = p1 - p0;
+			h2 = p2 - p1;
+			s1 = (double(subs(f,p1)) - double(subs(f,p0)))/h1;
+			s2 = (double(subs(f,p2)) - double(subs(f,p1)))/h2;
+			d = (s2 - s1)/(h2 + h1);
+			++cnt;
+		end
+	end
+end

NewNewton.m

@@ -0,0 +1,25 @@
+function [xList, errorList] = NewNewton (f, startingValue, errorAllowed)
+%
+%NewNewton (f, startingValue, errorAllowed)
+%This function computes the root of a function using the modified Newton's method
+%
+
+	pkg load symbolic;
+	warning('off','OctSymPy:sym:rationalapprox');
+	oldAnswer = startingValue;
+	newAnswer = 0;
+	currentError = errorAllowed + 1;
+	cnt = 1;
+	maxIt = 50;
+	fp = diff(f);
+	fpp = diff(fp);
+
+
+	while((currentError >= errorAllowed) && (cnt < maxIt))
+		newAnswer = oldAnswer - ((double(subs(f,oldAnswer)) * double(subs(fp,oldAnswer)))/(double(subs(fp,oldAnswer))^2 - (double(subs(f,oldAnswer)) * double(subs(fpp,oldAnswer)))));
+		currentError = abs(newAnswer - oldAnswer);
+		xList(end+1) = newAnswer;
+		errorList(end+1) = currentError;
+		oldAnswer = newAnswer;
+	end
+end

Newton.m

@@ -0,0 +1,24 @@
+function [xList,errorList] = Newton (f, startingValue, errorAllow)
+%
+%This function uses Newtons method to find a solution to the root of f
+%
+
+	pkg load symbolic;
+	warning('off','OctSymPy:sym:rationalapprox');
+	maxIt = 50;
+	fp = diff(f);
+	oldAnswer = startingValue;
+	newAnswer = 0;
+	cnt = 1;
+	currentError = errorAllow + 1;
+
+	%Loop until the error becomes small enough or the maximum number of itterations is met
+	while((cnt < maxIt) && (currentError > errorAllow))
+		newAnswer = oldAnswer - (double(subs(f,oldAnswer))/double(subs(fp,oldAnswer)));
+		currentError = abs(newAnswer - oldAnswer);
+		xList(end+1) = newAnswer;
+		errorList(end+1) = currentError;
+		oldAnswer = newAnswer;
+		++cnt;
+	end
+end

Secant.m

@@ -0,0 +1,30 @@
+function [xList, errorList] = Secant(f, p0, p1, errorAllowed)
+%
+%Secant(f, p0, p1, errorAllowed)
+%This function find the root of a function using the Secant Method
+%
+
+	pkg load symbolic;
+	warning('off','OctSymPy:sym:rationalapprox');
+	maxIt = 50;
+	cnt = 2;
+	q0 = double(subs(f,p0));
+	q1 = double(subs(f,p1));
+	currentError = errorAllowed + 1;
+	p = 0;
+
+
+	while((cnt <= maxIt) && (currentError >= errorAllowed))
+		p = p1 - (q1 * (p1 - p0))/(q1 - q0);
+		currentError = abs(p - p1);
+		%Add the x and error values to memory
+		xList(end+1) = p;
+		errorList(end+1) = currentError;
+		%Setup for the next run
+		++cnt;
+		p0 = p1;
+		q0 = q1;
+		p1 = p;
+		q1 = double(subs(f,p));
+	end
+end

swap.m

@@ -0,0 +1,8 @@
+function [x,y] = swap(first, second)
+%
+%Swap the first and second values
+%
+
+	x = second;
+	y = first;
+end