function [x,v,strikes,delv] = mover(x,v,npart, ...
xwall,mpv,vwall,tau);
% mover - Function to move particles by free flight
% Also handles collisions with walls
% Inputs
% x - Positions of the particles
% v - Velocities of the particles
% npart - Number of particles in the system
% xwall - Positions of walls
% mpv - Most probable velocity off the wall
% vwall - Wall velocities
% tau - Time step
% Output
% x,v - Updated positions and velocities
% strikes - Number of particles striking each wall
% delv - Change of y-velocity at each wall
%
strikes = [0 0]; delv = [0 0];
direction = [1 -1]; % Direction of particle leaving wall
x_old = x; % Remember original position
% Update x position of particle
x(:) = x_old(:) + v(:,1)*tau;
%x = rem(x+xwall(2),xwall(2)); % Periodic boundary conditions
%return;
for i=1:npart
if( x(i) <= xwall(1) )
flag=1; % Particle strikes left wall
elseif( x(i) >= xwall(2) )
flag=2; % Particle strikes right wall
else
flag=0; % Particle strikes neither wall
end
if( flag > 0 )
strikes(flag) = strikes(flag) + 1;
delv(flag) = delv(flag) - v(i,2);
% Reset velocities as biased Maxwellian
v(i,1) = direction(flag)*sqrt(-log(1-rand(1))) * mpv;
phi = 2*pi*rand(1);
v_parll = sqrt(-log(1-rand(1))) * mpv;
v(i,2) = v_parll*sin(phi) + vwall(flag);
v(i,3) = v_parll*cos(phi);
% Time of flight after leaving wall
dtr = tau*(x(i)-xwall(flag))/(x(i)-x_old(i));
% Position after leaving wall
x(i) = xwall(flag) + v(i,1)*dtr;
delv(flag) = delv(flag) + v(i,2);
end
end