function ydot=langevin_ode(t,Y)
global kb T M N;
w=-.1;
for i=1:N
    v=[Y(2*N+2*i-1) Y(2*N+2*i)];
    
    
    F=-.02*M.*v;%+sqrt(kb*L*T*2*.02)*(-.5+rand());
    
    F(1)=F(1)+sqrt(kb*M*T*2*.02)*(-.5+rand());
    F(2)=F(2)+sqrt(kb*M*T*2*.02)*(-.5+rand());
    
    %F=-.02*L.*v+sqrt(kb*L*T*2*.02)*(-.5+rand());
    a=F/M;
%     ydot=[Y(2*N+2*i-1) Y(2*N+2*i) a(1) a(2)];
    ydot(2*i-1) = v(1);
    ydot(2*i)=v(2);
    ydot(2*N+2*i-1)=a(1);
    ydot(2*N+2*i)=a(2);
end
ydot=ydot';
end
