1

Not sure if I should have posted this here or on the Maths stackexchange, so sorry if it's the wrong place. I'm very new to MATLAB and programming in general, and am having some troubles trying to solve an ODE problem using finite difference methods for an assignment.

My finite difference equation is: 

z(t+dt) = (dt^2*(γ^2*h*sin(γ*t)-β*z(t)) - z(t-dt)*(1-dt*α)+2*z(t))/(1 + dt*α)

Where t is a 51x1 array for the time increments. Basically I want to calculate z(t) for t values from 0 to 1 in increments of 0.02. I have the initial conditions z(0) = 0 and z(Δt) = 0.

My current code (not everything, but the bit that's giving me trouble:

    dt = 0.02

    t = [0:dt:T]';

    z(0) = 0
    z(dt)= 0

    for i = t

        z(i+dt) = (dt^2*(gamma^2.*h.*sin(gamma*t)-beta*z(i)) - z(i-dt)*(1-dt*alpha)+2*z(i))/(1 + dt*alpha)

    end

Alpha, beta and gamma are all constants in this case, they're defined earlier in the code.

I keep getting the error "Subscript indices must either be real positive integers or logicals." I understand that MATLAB arrays begin with element 1 and not 0, so trying to access element 0 will give an error.

I'm not sure if the error is with how I've entered my finite difference function, or the initial conditions. By setting i = t, am I running the for loop for those values of t, or for those elements in the matrix? E.g. when i = 0, is it trying to access the 0 element of the matrix, or is it setting the i variable in the equation to 0 like I want it to?

Any help would be greatly appreciated.

Thankyou!

Liam
  • 11
  • 2

1 Answers1

1

problem is with dt and I think dt is not integer value and cannot be used for indexing array. index of arrays are always integer values in Matlab. Please try below code and check if solves the problem

t = [0:dt:T]';

z(1) = 0
z(2)= 0

for i = 2 : length(t)

    z(i) = (dt^2*(gamma^2.*h.*sin(gamma*t(i))-beta*z(i)) - z(i-1)*(1-dt*alpha)+2*z(i))/(1 + dt*alpha)

end
User1551892
  • 3,236
  • 8
  • 32
  • 52
  • Thanks for the reply @User1551892. I tried your code and got another matrix error, "Index exceeds matrix dimensions." – Liam Sep 23 '16 at 14:00
  • So where you have for i = 2 : length(t) That's just running the for loop from z(Δt) through to z(1)? – Liam Sep 23 '16 at 14:07
  • Make sure that `z` is allocated first before running your code: `z = zeros(1, numel(t));` after the declaration of `t`. – rayryeng Sep 23 '16 at 14:22
  • Ah ok, so it was running for an undeclared output z matrix size? Thanks very much @rayryeng, it's working now, outputting a 1x51 matrix with my z values, which look pretty good to me. – Liam Sep 23 '16 at 14:29
  • @Liam Yes that's right. You were trying to index into a matrix / array that was not previously created. Also, you're welcome!... even though User1551892 addressed the more pressing issue regarding your code. – rayryeng Sep 23 '16 at 14:30
  • @User1551892: I have not upvoted your answer as it isn't entirely correct. If you correct the issue re: the allocation of `z`, you have my vote. – rayryeng Sep 23 '16 at 14:31
  • @rayryeng: Perhaps, i have some problem with understanding. This code works fine with prior declaration of kl. for i = 1:5 kl (i) = i ;end; and value 1 to 5 get stored in kl array – User1551892 Sep 23 '16 at 16:02
  • @User1551892 That is correct. The OP **did not declare** `z` which is why the code didn't work. Your code snippet in your comment assumes declaration of `kl`, which is why it works. The answer you provided does not assume so and ultimately why the OP couldn't get it to work initially. – rayryeng Sep 23 '16 at 16:05
  • @rayryeng : still I did not get. In my code snippet I am not declaring 'kl' prior to calling for loop and 'kl' is not available in work space either. Can you elaborate it why my code snippet is not giving error. – User1551892 Sep 23 '16 at 16:29
  • @User1551892 Ah right. Never mind. You are correct. What probably happened was that the OP already declared `z` to begin with but the size was smaller than the length of `t`. Without `z` being in the workspace, the variable should dynamically create itself and have the right space whenever you populate elements in. I apologize. I forget how MATLAB internals work sometimes. – rayryeng Sep 23 '16 at 16:31
  • @rayryeng: thanks. I think OP did not declare z(1)=0 and z(2)=0 as i did in my solution. – User1551892 Sep 23 '16 at 16:35
  • @User1551892 I see that now. OK, this is correct. You have my upvote. – rayryeng Sep 23 '16 at 16:42