# uberj/MTH351

lab3 initial commit

1 parent 8054162 commit e9665ebf00b0f7cdf1a0affadeec8dfd7e43307a committed Oct 22, 2012
4 labs/lab3/Makefile
 @@ -0,0 +1,4 @@ +all: + MATLABPATH=\$MATLABPATH:./ + export MATLABPATH + matlab -nosplash -nodisplay -nodesktop
1 labs/lab3/bisect.m
11 labs/lab3/lab3.m
 @@ -0,0 +1,11 @@ +%% Lab 3 + +the_root = 1 +six_digits_of_accuracy =.5E-6; +%% Problem 1 +%% + +k = 1; +while abs(eq1(reqx, k) - xtrue) > six_digits_of_accuracy + k = k + 1; +end
BIN labs/lab3/lab3.pdf
Binary file not shown.
9 labs/lab3/make_lib.m
 @@ -0,0 +1,9 @@ +function lib = make_lib() + lib.f = @the_function; + lib.the_root = 1; + lib.accuracy = 5E-6; +end + +function y = the_function(x) + (x ^ 5) - (x ^ 4) + x - 1; +end
1 labs/lab3/newton.m
75 labs/lab3/sample/bisect.m
 @@ -0,0 +1,75 @@ +function [it_count,root,xn]=bisect(fcn,a0,b0,ep,max_iterate) +% +% function bisect(fcn,a0,b0,ep,max_iterate) +% +% This is the bisection method for solving an equation f(x)=0. +% +% The function f=fcn is a string. The function f is +% to be continuous on the interval [a0,b0], and it is to be of +% opposite signs at a0 and b0. The quantity ep is the error +% tolerance. The routine guarantees this as an error bound +% provided: (1) the restrictions on the initial interval are +% correct, and (2) ep is not too small when the machine epsilon +% is taken into account. Most of these conditions are not +% checked in the program! The parameter max_iterate is an upper +% limit on the number of iterates to be computed. +% +% For the given function f(x), an example of a calling sequence +% might be the following: +% root = bisect('x^2-1',1,1.5,1.0E-6,10) +% The parameter index_f specifies the function to be used. +% +% The following will optionally print out for each iteration the values +% of count, a, b, c, f(c), (b-a)/2 +% with c the current iterate and (b-a)/2 the error bound for c. +% The variable count is the index of the current interate. + +it_count=0; +xn=a0; +if a0 >= b0 + disp('a0 < b0 is not true. Stop!') + root=NaN; + return +end + +a = a0; b = b0; +f=inline(fcn); +fa = feval(f,a); fb = feval(f,b); + +if sign(fa)*sign(fb) > 0 + disp('f(a0) and f(b0) are of the same sign. Stop!') + root=NaN; + return +end + +c = (a+b)/2; +xsave=zeros(max_iterate,1); % for saving each iteration +xsave(it_count+1)=c; + +while b-c > ep & it_count < max_iterate + it_count = it_count + 1; + fc = feval(f,c); + %Internal print of bisection method: + %format short e + %iteration = [it_count a b c fc b-c] + if sign(fb)*sign(fc) <= 0 + a = c; + fa = fc; + else + b = c; + fb = fc; + end + c = (a+b)/2; + xsave(it_count+1)=c; +end + +root = c; +xn=xsave(1:it_count+1); + +%format long +%root +%format short e +%error_bound = b-c +%format short +%it_count +
66 labs/lab3/sample/newton.m
 @@ -0,0 +1,66 @@ +function [it_count,root,xn] = newton(fcn,dfcn,x0,error_bd,max_iterate) +% +% function newton(fcn,x0,error_bd,max_iterate) +% +% This is Newton's method for solving an equation f(x) = 0. +% +% The functions f(x)=fcn and deriv_f(x)=dfcn are strings. +% The parameter error_bd is used in the error test for the +% accuracy of each iterate. The parameter max_iterate +% is an upper limit on the number of iterates to be +% computed. An initial guess x0 must also be given. +% +% For the given function f(x), an example of a calling sequence +% might be the following: +% root = newton('x^2-1','2*x',.1,1.0E-12,10) +% +% The program _optionally_ prints the iteration values +% iterate_number, x, f(x), deriv_f(x), error +% The value of x is the most current initial guess, called +% previous_iterate here, and it is updated with each iteration. +% The value of error is +% error = newly_computed_iterate - previous_iterate +% and it is an estimated error for previous_iterate. + + +f=inline(fcn); +df=inline(dfcn); +error = 1; +it_count = 0; +xsave=zeros(max_iterate,1); % for saving each iteration +xsave(1)=x0; + +while abs(error) > error_bd & it_count < max_iterate + fx = feval(f,x0); + dfx = feval(df,x0); + if dfx == 0 + disp('The derivative is zero. Stop') + root=NaN; + return + end + it_count = it_count + 1; + x1 = x0 - fx/dfx; + xsave(it_count+1)=x1; + error = x1 - x0; + %Internal print of newton method: + %format short e + %iteration = [it_count x0 fx dfx error] + x0 = x1; +end + +root = x1; +xn=xsave(1:it_count+1); + +if it_count > max_iterate + disp('The number of iterates calculated exceeded') + disp('max_iterate. An accurate root was not') + disp('calculated.') +end + +% format long +% root +% format short e +% error +% format short +% it_count +
14 labs/lab3/sample/samp_prob1.m
 @@ -0,0 +1,14 @@ +tol=1e-6; +max_its=100; +fcn='x^6-x-1'; +trueroot=[1.13472414]; +intervals=[0.9, 2.7; 1, 2; 0.75, 1.25]; + +disp(sprintf('\nBisection estimate for root of %s with accuracy of %g:',fcn,tol)); +disp(sprintf('_interval_ \t _estimate_ \t _error_ \t _iterations_')) + +for i=1:size(intervals,1), + [itiB(i),rootiB(i)]=bisect(fcn,intervals(i,1),intervals(i,2),tol,max_its); + disp(sprintf('[%g,%g] \t %0.8f \t %0.5e \t %d',... + intervals(i,:),rootiB(i),abs(trueroot-rootiB(i)),itiB(i))); +end
16 labs/lab3/sample/samp_prob2.m
 @@ -0,0 +1,16 @@ +tol=1e-6; +max_its=100; +fcn='x^6-x-1'; +dfcn='6*x^5-1'; +trueroot=[1.13472414]; +inits=[1.1,1.2,1.3,1.4,1.5]; + +disp(sprintf('\nNewton estimate of root of %s with accuracy of %g',fcn,tol)); +disp(sprintf('_initial_ \t _estimate_ \t _error_ \t _iterations_')) +for i=1:length(inits), + + [itiN(i),rootiN(i)]=newton(fcn,dfcn,inits(i),tol,max_its); + disp(sprintf('%g\t\t %0.8f \t %0.5e \t %d',... + inits(i),rootiN(i),abs(trueroot-rootiN(i)),itiN(i))); +end +
70 labs/lab3/sample/secant.m
 @@ -0,0 +1,70 @@ +function [it_count,root,xn] = secant(fcn,x0,x1,error_bd,max_iterate) +% +% function secant(fcn,x0,x1,error_bd,max_iterate) +% +% This implements the secant method for solving an +% equation f(x) = 0. The function f=fcn is a string. +% +% The parameter error_bd is used in the error test for the +% accuracy of each iterate. The parameter max_iterate is an +% upper limit on the number of iterates to be computed. Two +% initial guesses, x0 and x1, must also be given. +% +% For the given function f(x), an example of a calling sequence +% might be the following: +% root = secant('x^2-1',x0,x1,1.0E-12,10) +% +% The program optionally prints the iteration values +% iterate_number, x, f(x), error +% The value of x is the most current initial guess, called +% previous_iterate here, and it is updated with each iteration. +% The value of error is +% error = newly_computed_iterate - previous_iterate +% and it is an estimated error for previous_iterate. + +error = 1; +f=inline(fcn); +fx0 = feval(f,x0); +it_count = 0; +%format short e +%iteration = [it_count x0 fx0] +xsave=zeros(max_iterate,1); xn=x0; % for saving each iteration +xsave(1)=x0; +xsave(2)=x1; + +while abs(error) > error_bd & it_count <= max_iterate + fx1 = feval(f,x1); + if fx1 - fx0 == 0 + disp('f(x1) = f(x0); Division by zero; Stop') + root=NaN; + return + end + it_count = it_count + 1; + x2 = x1 - fx1*(x1-x0)/(fx1-fx0); + xsave(it_count+2)=x2; + error = x2 - x1; + %Internal print of secant method: + %format short e + %iteration = [it_count x1 fx1 error] + x0 = x1; + x1 = x2; + fx0 = fx1; +end + +root = x2; +xn=xsave(1:it_count+2); + +if it_count > max_iterate + disp('The number of iterates calculated exceeded') + disp('max_iterate. An accurate root was not') + disp('calculated.') +end + +%format long +%root +%format short e +%error +%format short +%it_count + +
1 labs/lab3/secant.m
3 labs/lab3/the_function.m
 @@ -0,0 +1,3 @@ +function y = the_function(x) + (x ^ 5) - (x ^ 4) + x - 1; +end