MATLAB question re interpolation and approximating derivative.

[peripatein]

Hi,
I'd sincerely appreciate it if someone were willing to review the few lines of MATLAB code below and indicate why they don't quite yield the expected output.

Homework Statement

I am asked to generate using MATLAB approximated values of f(x)=cos(x) at nodes x+h,x-h with random errors <=5*10^-6 (using rand) for h=10^-8,10^-7,...,10^-1. Hence,
f^~(x+h)=f(x+h)+e(x+h)
f^~(x-h)=f(x-h)+e(x-h)
where |e(x)|<=5*10^-6
in order to then find approximation for f'(1.2) by using the approximation:
f'(x)=[f(x+h)-f(x-h)]/2h
I am finally asked to plot the error with respect to the value of h.

Below is my code. I am not really sure why it yields one line across the y-axis and another across the x axis.

Homework Equations

The Attempt at a Solution

h=(10^-1).^[1:8];
x=1.2;
fminush=cos(x-h)+(5e-6)*rand(1,1);
fplush=cos(x+h)+(5e-6)*rand(1,1);
fder=(fplush-fminush)./(2*h);
plot(h,abs(-sin(x)-fder))

 

[mfb]

What do you get for fder (typical values)?
Do you use a logarithmic scale for h?

 

[peripatein]

h was defined as a vector whose values are:
10^-8,10^-7,10^-6,10^-5,10^-4,10^-3,10^-2,10^-1
What do you mean by typical values? Are you asking which values I obtain for that parameter (viz. fder) as I execute the code?

 

[mfb]

h was defined as a vector whose values are:
10^-8,10^-7,10^-6,10^-5,10^-4,10^-3,10^-2,10^-1

Sure.

What do you mean by typical values? Are you asking which values I obtain for that parameter (viz. fder) as I execute the code?

Yes. Please post them, so it is easier to understand what went wrong.
I said "typical" as your code has random numbers in it, so the result won't be the same every time you run the code.

 

[peripatein]

fder:
2.27670626629095e-06
2.27670626684606e-05
0.000227670626684606
0.00227670626684606
0.0227670626684606
0.227670626684606
2.27670626684606
22.7670626684606

 

[mfb]

Hmm, you are adding the same random value every time. That's probably not the intended way.

Apart from that: okay. Use a logarithmic scale for both h and fder.

 

[peripatein]

I am not sure how to go about these two. Could you kindly guide me through? I am not asking you to state precisely how it should be carried out, but elaborating a bit more would prove very helpful.
For instance, need I to reset the rand seed to solve problem #1?
For #2, did you mean I should log both sides of the equation?
What about the plotting itself? How may I correct that?

 

[mfb]

For instance, need I to reset the rand seed to solve problem #1?

No, it is sufficient to let MATLAB generate 16 random numbers. You can generate them as array, or try .+ or .* in the lines where you use them.

For #2, did you mean I should log both sides of the equation?

The scale in the plot should be logarithmic.

 

[peripatein]

How may I change the scale in the plot to be logarithmic? Did you by any chance mean something like this?:
fminush=cos(x-h).+(5e-6)*rand(1,1);
fplush=cos(x+h).+(5e-6)*rand(1,1);
plot(log(h),log(abs(-sin(x)-fder)))

 

[mfb]

That could work. Just test it? I don't have access to MATLAB currently.
If it does not work, the manual will explain it, together with many websites.

 

[peripatein]

Well, in that case I have a few questions, if I may, concerning what you wrote, but also regarding the task in general.
(1) Does e(x+h) differ in any way from e(x)? I mean, I understand that e(x+h) designate the error for f(x+h) whereas e(x) designate that for f(x), but do they differ in value and, mainly, in function/expression?
(2) Should e(x) be defined as a function? For instance, the algebraic expression for the remainder in a Taylor series could be used, wherein values of h should/could be substituted, or, at least, which could be limited to random number which <= 5e-6. What I am asking is basically this - should e(x) (and e(x+h)) be evaluated algebraically or does the simple, most straightforward use of rand() in this case sufficient?
(3) How could f^~(x+h) be only an approximation for f(x+h) whereas the first term is f(x+h) itself?
(4) Why did you advise using a logarithmic scale and why wouldn't the plot be visible otherwise? Using a logarithmic scale was NOT part of the instructions.
(5) Below is my revised code:

h=(10^-1).^[1:8];
x=pi;
r=5e-6.*rand(16,1);
for i=1:8
fminush(i)=cos(x-h(i))+r(i);
fplush(i)=cos(x+h(i))+r(i+8);
end
fder=(fplush-fminush)./(2*h);
plot(log(h),log(abs(-sin(x)-fder)))

is there a way to avoid the loop? I have tried redefining r as (8,1) but I keep getting the error message: "Matrix dimentions must fit" or something of the sort.
(6) Here are the values I got for fder, upon executing the revised code:
4.617329096845424e-06
0.000146689415891910
-0.000499959502409109
0.00630868582096600
0.0195747882347774
-0.798172721805379
-6.79296949412933
107.411448829753

and attached is the plot. Is it actually supposed to look like that?

 

[mfb]

Well, in that case I have a few questions, if I may, concerning what you wrote, but also regarding the task in general.
(1) Does e(x+h) differ in any way from e(x)? I mean, I understand that e(x+h) designate the error for f(x+h) whereas e(x) designate that for f(x), but do they differ in value and, mainly, in function/expression?

The expression |e(x)|<c means "the function is smaller than c for every x" - it does not matter what you insert for x (as long as it is a real number). |e(x+h)|<c would have exactly the same meaning, just written in a more confusing way.

(2) Should e(x) be defined as a function?

Technically it is, but please don't go that way, this will just lead to confusion. It is some numerical error you introduce manually with the rand() (because MATLAB is more accurate than the calculation you simulate), done.

(3) How could f^~(x+h) be only an approximation for f(x+h) whereas the first term is f(x+h) itself?

You simulate a calculation with a larger error, so this calculation won't be exactly f(x+h).

(4) Why did you advise using a logarithmic scale and why wouldn't the plot be visible otherwise? Using a logarithmic scale was NOT part of the instructions.

Both values vary by 8 orders of magnitude. A logarithmic plot is the only reasonable way to make them visible at the same time.

Using a logarithmic scale was NOT part of the instructions.

At some point you'll have to learn how to plot things on your own,

and attached is the plot. Is it actually supposed to look like that?

For the function you plotted and x=pi, yes.
I think the plot gets more interesting for other values of x.

 

[peripatein]

Thank you very much! Is there a way to avoid using a for loop? I couldn't simply do cos(x-h)+r or even cos(x-h).r. Matlab won't allow that, even if I redefined r to be of size 8 instead of 16.

 

[mfb]

I would have expected that .+ works with 8 elements for r. If it does not, I don't know.