21 views (last 30 days)

function za = dfdt(z, fs);

h = 1/fs;

N = length(z);

for i = 1:N

za(i) = (z(i+1) - z(i-1))/(2*h);

end

za(1) = za(2);

za(N) = za(N-1);

I know where the problem is: at the for-construction. It should be 2:N-1, but why?

The same goes for the 2 lines at the bottom. The problem is that my vectors are not the same length when they are not there, but can someone explain why this problem is solved by za(1) = za(2) and za(N) = za(N-1)?

John D'Errico
on 16 Jan 2021

Edited: John D'Errico
on 16 Jan 2021

The problem lies here:

for i = 1:N

za(i) = (z(i+1) - z(i-1))/(2*h);

end

When the loop starts out, what is i? i is 1, right?

Now, look at that line. What is z(i-1), thus z(0)?

Does MATLAB allow a zero index? (No.)

What was the error message?

Array indices must be positive integers or logical values

Do you see the problem? You will need to decide how to treat that first element.

Actually, you will also need to decide how to deal woth the LAST element too, since your vector has length N. But, when i==N, that line of code will try to acces element N+1. And that will cause another error.

You are taking the difference of the neighbors, then dividing by 2*h. So effectively, you are computing a second order approximation to the derivative at each point. There are two common options. A simple one is to use the common lower (first) order approximation to the derivative at the first and last points. You should know what that is.

Alternatively, you can use a non-central 3 point approximation. Either will work. Be careful, as this means you will use a forward difference at node 1, while you use a backwards difference at node N. Similarly, you will have subtly different 3 point approximations used at nodes 1 and N, if you choose that route.

Oh. I just saw you also asked this:

but can someone explain why this problem is solved by za(1) = za(2) and za(N) = za(N-1)?

And that is easy. Essentially this solution assumes the estimated derivative at node 1 is the same as that which was estimated at node 2. Since you cannot use the central difference approximation at node 1, the solution you saw just presumes that derivative did not change. In my opinion, this is a poor solution. Either of the approaches I suggest are better choices.

So your code might look like this:

za(1) = % YOUR CHOICE HERE!

for i = 2:N-1

za(i) = (z(i+1) - z(i-1))/(2*h);

end

za(N) = % also your choice here!

You need to fill in the holes I left behind, as that is why you have homework.

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!
## 0 Comments

Sign in to comment.