Problems with Matlab code for simulating spring motion

Join the discussion
Registration is free. Ask a follow-up in this thread, or start your own.
1 reply · 3K views
bcolson
Messages
3
Reaction score
0
The issue I am having is setting up the actual program for the problem. I am having trouble interpreting the code and, because of that, having issues getting it to work properly. The whole thing was written with Matlab and the goal was to write a function that could simulate the motion of a series of connected Hooke's law oscillators. Half of the system has stiffness 1 the other half has variable stiffness s.
Total size of the system should be 2*l.function [T,y] = mtxsolver(M,s,l)
%M == mass, s == stiffness, l == length

Size = 2*l;
size = Size-1;

%set up vectors
k = zeros((Size),1);
R = 1;

%set up k values in loops.
while R < (floor(Size/2))
k(R) = 1;
R = R + 1;
end

while R < (Size)
k(R) = s;
R = R + 1;
end

k(Size+1) = 0;

for c = 1:Size
for r = 1:Size
if r == c
%to diagonal
A(r,c) = (k(r)+k(r+1))/M;
elseif r-c == 1
%to left diagonal
A(r,c) = -k(r)/M;
elseif r-c == -1
%to right diagonal
A(r,c) = -k(r+1)/M;
else
%to all else
A(r,c) = 0;
end
end
end%set up initial conditions
initCond = zeros(2,floor((Size+1)/2));
initCond(1,1) = 1;

%ode function
fun = @(t,x)[x(2); - A*x(1)];
[T,y] = ode45(fun,[0,30],initCond);
endI can't wrap my head around the matrix size of the initial conditions. Why exactly does it have to be so much smaller than the size of the matrix describing the coupled stiffness? Is it something where ode45 is only using half of the matrix A for position calculations?
 
Physics news on Phys.org
bcolson said:
The issue I am having is setting up the actual program for the problem. I am having trouble interpreting the code and, because of that, having issues getting it to work properly. The whole thing was written with Matlab and the goal was to write a function that could simulate the motion of a series of connected Hooke's law oscillators.Half of the system has stiffness 1 the other half has variable stiffness s.
Total size of the system should be 2*l.

Matlab:
function [T,y] = mtxsolver(M,s,l)
%M == mass, s == stiffness, l == length

Size = 2*l;
size = Size-1;

%set up vectors
k = zeros((Size),1);
R = 1;

%set up k values in loops.
while R < (floor(Size/2))
    k(R) = 1;
    R = R + 1;
end

while R < (Size)
    k(R) = s;
    R = R + 1;
end

k(Size+1) = 0;

for c = 1:Size
  for r = 1:Size
      if r == c
           %to diagonal
           A(r,c) = (k(r)+k(r+1))/M;
       elseif r-c == 1
           %to left diagonal
           A(r,c) = -k(r)/M;
       elseif r-c == -1
           %to right diagonal
           A(r,c) = -k(r+1)/M;
       else
           %to all else
           A(r,c) = 0;
       end
   end
end%set up initial conditions
initCond = zeros(2,floor((Size+1)/2));
initCond(1,1) = 1;

%ode function
fun = @(t,x)[x(2); - A*x(1)];
[T,y] = ode45(fun,[0,30],initCond);
end

I can't wrap my head around the matrix size of the initial conditions. Why exactly does it have to be so much smaller than the size of the matrix describing the coupled stiffness? Is it something where ode45 is only using half of the matrix A for position calculations?
The ##[##code=matlab##]## and ##[##/code##]## make it look a lot better.

I know nothing about matlab, but you are start doubling the size. Why ? Could you adhere to the template and provide some more info on what is what ? For me size has the dimension of meters, but here it is the length of an array.
Note that you solve 'size' second order differential equations. Where are the velocities and their initial conditions ?