Question

Write a single MATLAB function called numdiff, saved in a file called numdiff.m. Submit that single file here in Canvas. To make grading easier, please do the following: • Close your function with "end". • Don't have your function print anything to the console; just return the output. To accomplish this: o Put a semicolon ";" at the end of each statement. o Put a percent "%" before each disp or fprintf statement. The function numdiff should take the following inputs, in order: • A real-valued function f. • A positive number h. • A real number a. The function numdiff should return as output a 3 x 9 matrix as follows: The first row should contain, in order, the following approximations for f'(a) using step size h: o The actual value of f'(a), computed symbolically. o The forward-difference formula. o The backward-difference formula. o The three-point midpoint formula o The three-point endpoint formula, assuming a is the left endpoint. o The three-point endpoint formula, assuming a is the right endpoint. . The five-point midpoint formula o The five-point endpoint formula, assuming a is the left endpoint. The five-point endpoint formula, assuming a is the right endpoint. • The second row should contain the absolute error in each approximation. • The third row should contain the relative error in each approximation.

You may assume that f(a), f(a+h), f(a + 2h), and f'(a) are all defined. However, you should NOT assume that f'(a) is necessarily nonzero; if f'(a) = 0, your function should still return each approximation and its absolute error. Use the MATLAB command "diff" to calculate the derivative of f(a). As many of you discovered in project 2, one way to calculate f'(a) is the following (assuming that the function handle f and the number a have already been defined): syms x df = diff( f(x) ) Df = eval( subs( df, x, a ) )

Answer 1

function dMat = numdiff(f,h,a)
syms x
df = diff(f(x));
trueValue = eval(subs(df,x,a));
dMat=zeros(3,9);
dMat(1,1)=trueValue;
dMat(1,2)=(f(h+a)-f(a))/h;
dMat(1,3)=(f(a)-f(a-h))/h;
dMat(1,4)=(f(h+a)-f(a-h))/(2*h);
dMat(1,5)=(-f(a+2*h)+4*f(a+h)-3*f(a))/(2*h);
dMat(1,6)=(f(a-2*h)-4*f(a-h)+3*f(a))/(2*h);
dMat(1,7)=(f(a-2*h)-8*f(a-h)+8*f(a+h)-f(a+2*h))/(12*h);
dMat(1,8)=(-25*f(a)+48*f(a+h)-36*f(a+2*h)+16*f(a+3*h)-3*f(a+4*h))/(12*h);
dMat(1,9)=(25*f(a)-48*f(a-h)+36*f(a-2*h)-16*f(a-3*h)+3*f(a-4*h))/(12*h);
dMat(2,:)=abs(trueValue-dMat(1,:));
if trueValue~=0%relative error will be calculated only if derivativeis not 0
    dMat(3,:)=dMat(2,:)/trueValue;
end