i am writing in matlab a code to simulate a mass spring system, here there is the code:

test.m

- Code: Select all
`x0=[-1 0 1];`

y0=[0 1/2 0];

z0=[0 0 0];

Vx0=[0 0 0];

Vy0=[0 1 0];

Vz0=[0 0 0];

A=[0 1 0; 1 0 1; 0 1 0];

l = A ;

k = A ;

m = [1 1 1];

u0 = [x0';y0';z0';Vx0';Vy0';Vz0'];

tspan=[0 10];

[t,u] = ode45(@molle,tspan,u0,[],A,l,k,m);

x=u(:,1:3:3*size(A,1));

y=u(:,2:3:3*size(A,1));

z=u(:,3:3:3*size(A,1));

vx=u(:,3*size(A,1)+1:3:6*size(A,1));

vy=u(:,3*size(A,1)+2:3:6*size(A,1));

vz=u(:,3*size(A,1)+3:3:6*size(A,1));

maxX=max(max(x));

maxY=max(max(y));

maxZ=max(max(z));

minX=min(min(x));

minY=min(min(y));

minZ=min(min(z));

figure,

for k = 1 : length(t)

[aa bb]=view;

plot3(x(k,:),y(k,:),z(k,:),'*')

view(aa,bb)

axis([minX maxX minY maxY minZ maxZ])

grid on

pause(0.5)

end

molle.m

- Code: Select all
`function u = molle(tspan,u0,A,l,k,m)`

N=size(A,1)

x=u0(1:3:3*N);

y=u0(2:3:3*N);

z=u0(3:3:3*N);

vx=u0(3*N+1:3:6*N);

vy=u0(3*N+2:3:6*N);

vz=u0(3*N+3:3:6*N);

for i = 1 : N

Fx=0;

Fy=0;

Fz=0;

for j = 1 : N

if A(i,j) ~= 0

d=((x(i)-x(j))^2+(y(i)-y(j))^2+(z(i)-z(j))^2)^0.5;

Fx=Fx+k(i,j)*(abs(x(j)-x(i))-l(i,j))*(x(j)-x(i))*A(i,j)/d;

Fy=Fy+k(i,j)*(abs(y(j)-y(i))-l(i,j))*(y(j)-y(i))*A(i,j)/d;

Fz=Fz+k(i,j)*(abs(z(j)-z(i))-l(i,j))*(z(j)-z(i))*A(i,j)/d;

end

end

u(1+3*(i-1))=u0(3*N+1+3*(i-1)); % dx

u(2+3*(i-1))=u0(3*N+2+3*(i-1)); % dy

u(3+3*(i-1))=u0(3*N+3+3*(i-1)); % dz

u(3*N+1+3*(i-1))=Fx/m(i) ; % ax=F/m

u(3*N+2+3*(i-1))=Fy/m(i) ; % ay=F/m

u(3*N+3+3*(i-1))=Fz/m(i) ; % az=F/m

end

u=u';

i don't know how to write in python ode45 function....could you help me please????

Thanks in advance...

Francesco