# Homework Help: 2nd order initial value problem in matlab

1. Jan 29, 2013

### epic325

1. The problem statement, all variables and given/known data

1.)I want to write a function in matlab that contains the 2nd order function:

20*d$^{2}$x;(dt$^{2}$)+5*dx/dt + 20*x=0 (dampened spring)

-The function should have 2 inputs (time,[initial values]) initial values should be a vector of 2 values
-The function should output a vector of the time derivatives

2.)I want to write a script that that uses the function to plot the function using euler's method with a for loop.

2. Relevant equations
I solved for the time derivatives

[w(2) ; -.25w(2)-w(1)] where my initial values are [1; 0] x(0)=1 and dx/dt(0)=0

3. The attempt at a solution
I think I have figured the code out but it's not plotting the right thing. it should be a dampened oscillator, but its plotting a line. I just need some help with the debugging. I will paste the function code followed by the separate script:

%This is the function where the input will be t and wVec=[initial values]

function wDerVec = odeFunc(t,wVec)
wDerVec=[wVec(2) -.25*wVec(2)-wVec(1)];
end

clear
clc
h=0.1; %h is the time step (.1s)
t=0:h:15; %time interval (0 to 15 s)
wVec=[1 0]; %initial condition

for i=2:length(t),
odeFunc(t,wVec);
wVec(i,:)=wVec(i-1,:)+ans*h;
end

plot(t,wVec(:,1)); %plots the distance with respect to time
title('Euler Approximation, damped spring');
xlabel('Time');
ylabel('distance, x');

I'm not very good at matlab so any help would be greatly appreciated! Thanks!

2. Jan 31, 2013

### jfgobin

Good Afternoon epic325,

I tried your equation in Octave and had no problem. Here is the sequence of commands I used.

Code (Text):

function Y=odeFCT(x,t)
Y=zeros(2,1);
Y(1)=x(2);
Y(2)=(-1/20)*(5*x(2)+20*x(1));
endfunction
t=linspace(0,50,1001);
y=lsode(@odeFCT,[1 0]',t);
plot(t,y(:,1))

Here is the graph obtained (plotted with gnuplot instead of Octave.)

Let me dive into your code.

#### Attached Files:

• ###### dampened-spring.png
File size:
18.3 KB
Views:
102
3. Jan 31, 2013

### jfgobin

Quick question -

In your code, you have

Code (Text):

for i=2:length(t),
odeFunc(t,wVec);
wVec(i,:)=wVec(i-1,:)+ans*h;
end

Shouldn't it be

Code (Text):

for i=2:length(t),
odeFunc(t,wVec(i-1,:));
wVec(i,:)=wVec(i-1,:)+ans*h;
end

?

4. Feb 3, 2013

### epic325

Yes, jfgobin. That was the problem. Thank you very much for your help!

Share this great discussion with others via Reddit, Google+, Twitter, or Facebook