diff --git a/MDCode/InitVStream.m b/MDCode/InitVStream.m
index 4290543..d2f214f 100644
--- a/MDCode/InitVStream.m
+++ b/MDCode/InitVStream.m
@@ -1,33 +1,51 @@
-doPlot = 1;
-dt = 5e-15;
-TStop = 3000 * dt;
-InitDist = 0.0;
-Method = 'VE'; % VE -- verlot; FD -- Forward Difference
-
-Mass0 = 14 * C.am; % Silicon
-Mass1 = 100 * C.am; % Argon
-
-AtomSpacing = 0.5430710e-9;
-LJSigma = AtomSpacing / 2^(1 / 6);
-LJEpsilon = 1e-21;
-
-PhiCutoff = 3 * AtomSpacing * 1.1;
-
-T = 30;
-
-AddRectAtomicArray(10, 10, 0, 0, 0, 0, 0, T, 0);
-% vy0 = -sqrt(0.02*Ep/Mass1);
-% AddRectAtomicArray(4,4,0,12*AtomSpacing,0,vy0,0,T,1);
-Ep = 2;
-AddParticleStream(5, 0.1, 10, -pi / 2, 1, Ep * C.q_0, 5);
-
-Size = 10*AtomSpacing;
-Limits = [-Size +Size -Size +Size]; % square is good
-PlDelt = 5 * dt;
-
-PlotFile = 'BlockSt.gif';
-doPlotImage = 1;
-PlotSize = [100, 100, 1049, 1049];
-
-ScaleV = .02e-11;
-ScaleF = 10;
+doPlot = 1;
+dt = 5e-15;
+TStop = 3000 * dt;
+InitDist = 0.0;
+Method = 'VE'; % VE -- verlot; FD -- Forward Difference
+
+Mass0 = 14 * C.am; % Silicon
+Mass1 = 100 * C.am; % Argon
+
+AtomSpacing = 0.5430710e-9;
+LJSigma = AtomSpacing / 2^(1 / 6);
+LJEpsilon = 1e-21;
+
+PhiCutoff = 3 * AtomSpacing * 1.1;
+
+T = 30;
+
+LAtoms=10;
+WAtoms=10;
+X0=0;
+Y0=0;
+VX0=0; 
+VY0=0;
+InitDist=0; 
+Temp=T;
+Type=0;
+AddRectAtomicArray(LAtoms, WAtoms, X0, Y0, VX0, VY0, InitDist, Temp, Type);
+
+% vy0 = -sqrt(0.02*Ep/Mass1);
+% AddRectAtomicArray(4,4,0,12*AtomSpacing,0,vy0,0,T,1);
+
+Ep = 2;
+numStream=5;
+x0=0.1;
+y0=10;
+PartAng=-pi/2;
+Type=1;
+Ep=Ep*C.q_0;
+Seper=5;
+AddParticleStream(numStream, x0, y0, PartAng, Type, Ep, Seper);
+
+Size = 10*AtomSpacing;
+Limits = [-Size +Size -Size +Size]; % square is good
+PlDelt = 5 * dt;
+
+PlotFile = 'BlockSt.gif';
+doPlotImage = 1;
+PlotSize = [100, 100, 1049, 1049];
+
+ScaleV = .02e-11;
+ScaleF = 10;
diff --git a/MDCode/PlotVars.m b/MDCode/PlotVars.m
index 76309f7..2ff4e4d 100644
--- a/MDCode/PlotVars.m
+++ b/MDCode/PlotVars.m
@@ -1,105 +1,111 @@
-function [ output_args ] = PlotVars(c, Limits)
-global Vx Vy L W x y Fx Fy C
-global Phi nAtoms time Mass0 Mass1 Pty0in Pty1in
-global LJEpsilon Phi0 PhiTot KETot MinX MaxX MinY MaxY
-global T T0 T1 ScaleV MarkerSize doPlotImage PlotCount
-global PlotFig map im PlotSize ScaleF
-
-if isempty(Limits)
-    Limits = Limits;
-end
-
-if doPlotImage == 1
-    if PlotCount == 0
-        close all;
-        set(0,'DefaultFigureWindowStyle','normal')
-        PlotFig = figure('Position', PlotSize);
-        PlotCount = 1;
-    else
-        PlotCount = PlotCount+1;
-    end
-end
-
-V2 = Vx.*Vx + Vy.*Vy;
-MaxV = max(sqrt(V2));
-dx = MaxX - MinX;
-dy =  MaxY - MinY;
-
-if c==0 && ~ScaleV
-    ScaleV = sqrt(dx*dx + dy*dy) / MaxV * 0.15;
-end
-
-subplot(3, 2, 1), plot(x(Pty0in), y(Pty0in), 'bo', 'markers',...
-    MarkerSize,'MarkerFaceColor', 'b');
-hold on
-subplot(3, 2, 1), plot(x(Pty1in), y(Pty1in), 'go', 'markers',...
-    MarkerSize,'MarkerFaceColor', 'g');
-subplot(3, 2, 1),quiver(x, y, Fx * ScaleF, Fy * ScaleF, 0, 'm', 'linewidth', 2);
-hold off
-axis(Limits);
-title('Atoms and Forces')
-xlabel('X')
-ylabel('Y')
-
-subplot(3, 2, 2), plot(time, KETot, 'linewidth', 2);
-hold on
-subplot(3, 2, 2), plot(time, PhiTot + LJEpsilon * nAtoms, 'g', 'linewidth',2);
-subplot(3, 2, 2), plot(time, PhiTot + KETot + LJEpsilon * nAtoms, 'k', 'linewidth',2);
-%     axis([0 TStop 0e-22 3e-23]);
-xlabel('time')
-ylabel('Energy')
-title('Total, Kinetic and Potential Energy')
-
-subplot(3,2,3), plot(x(Pty0in), y(Pty0in), 'bo', 'markers',...
-    MarkerSize,'MarkerFaceColor','b');
-hold on
-subplot(3,2,3), plot(x(Pty1in), y(Pty1in), 'go', 'markers',...
-    MarkerSize,'MarkerFaceColor','g');
-subplot(3,2,3),quiver(x, y, Vx * ScaleV, Vy * ScaleV, 0, 'r', 'linewidth', 2);
-hold off
-axis(Limits);
-title('Atoms and Velocities')
-xlabel('X')
-ylabel('Y')
-
-AE(Pty0in) = 1 / 2 * Mass0 * V2(Pty0in) + Phi(Pty0in) - Phi0;
-AE(Pty1in) = 1 / 2 * Mass1 * V2(Pty1in) + Phi(Pty1in) - Phi0;
-% if AddParticle, AE(end) = 0;end
-% AE = Phi;
-subplot(3, 2, 4),scatter3(x,y,AE,ones(1,nAtoms)*300,AE,'fill');
-axis(Limits);
-view(2);
-title('Atom Energy (KE + Pot)')
-
-subplot(3, 2, 5), plot(x(Pty0in), y(Pty0in), 'bo', 'markers',...
-    MarkerSize, 'MarkerFaceColor', 'b');
-hold on
-subplot(3, 2, 5), plot(x(Pty1in), y(Pty1in), 'go', 'markers',...
-    MarkerSize, 'MarkerFaceColor', 'g');
-hold off
-title('Atoms (All)')
-xlabel('X')
-ylabel('Y')
-
-subplot(3, 2, 6), plot(time, T, 'k', 'linewidth', 2);
-hold on
-subplot(3, 2, 6), plot(time, T0, 'b', 'linewidth', 2);
-subplot(3, 2, 6), plot(time, T1, 'g', 'linewidth', 2);
-%
-xlabel('time')
-ylabel('Temperature')
-title('Temperature: Each type and all')
-
-if doPlotImage
-    if PlotCount == 1
-        f = getframe(PlotFig);
-        [im, map] = rgb2ind(f.cdata, 256, 'nodither');
-        im(1, 1, 1, 2) = 0;
-
-    else
-        f = getframe(PlotFig);
-        im(:, :, 1, PlotCount) = rgb2ind(f.cdata, map, 'nodither');
-    end
-end
-
-end
+function [ output_args ] = PlotVars(c, Limits)
+global Vx Vy L W x y Fx Fy C
+global Phi nAtoms time Mass0 Mass1 Pty0in Pty1in
+global LJEpsilon Phi0 PhiTot KETot MinX MaxX MinY MaxY
+global T T0 T1 ScaleV MarkerSize doPlotImage PlotCount
+global PlotFig map im PlotSize ScaleF
+
+if isempty(Limits)
+    Limits = Limits;
+end
+
+if doPlotImage == 1
+    if PlotCount == 0
+        close all;
+        set(0,'DefaultFigureWindowStyle','normal')
+        PlotFig = figure('Position', PlotSize);
+        PlotCount = 1;
+    else
+        PlotCount = PlotCount+1;
+    end
+end
+
+V2 = Vx.*Vx + Vy.*Vy;
+MaxV = max(sqrt(V2));
+dx = MaxX - MinX;
+dy =  MaxY - MinY;
+
+if c==0 && ~ScaleV
+    ScaleV = sqrt(dx*dx + dy*dy) / MaxV * 0.15;
+end
+
+subplot(3, 2, 1), plot(x(Pty0in), y(Pty0in), 'bo', 'markers',...
+    MarkerSize,'MarkerFaceColor', 'b');
+hold on
+subplot(3, 2, 1), plot(x(Pty1in), y(Pty1in), 'go', 'markers',...
+    MarkerSize,'MarkerFaceColor', 'g');
+subplot(3, 2, 1),quiver(x, y, Fx * ScaleF, Fy * ScaleF, 0, 'm', 'linewidth', 2);
+hold off
+axis(Limits);
+title('Atoms and Forces')
+xlabel('X')
+ylabel('Y')
+
+subplot(3, 2, 2), plot(time, KETot, 'linewidth', 2);
+hold on
+subplot(3, 2, 2), plot(time, PhiTot + LJEpsilon * nAtoms, 'g', 'linewidth',2);
+subplot(3, 2, 2), plot(time, PhiTot + KETot + LJEpsilon * nAtoms, 'k', 'linewidth',2);
+%     axis([0 TStop 0e-22 3e-23]);
+xlabel('time')
+ylabel('Energy')
+title('Total, Kinetic and Potential Energy')
+
+subplot(3,2,3), plot(x(Pty0in), y(Pty0in), 'bo', 'markers',...
+    MarkerSize,'MarkerFaceColor','b');
+hold on
+subplot(3,2,3), plot(x(Pty1in), y(Pty1in), 'go', 'markers',...
+    MarkerSize,'MarkerFaceColor','g');
+subplot(3,2,3),quiver(x, y, Vx * ScaleV, Vy * ScaleV, 0, 'r', 'linewidth', 2);
+hold off
+axis(Limits);
+title('Atoms and Velocities')
+xlabel('X')
+ylabel('Y')
+
+AE0(Pty0in) = 1 / 2 * Mass0 * V2(Pty0in) + Phi(Pty0in) - Phi0;
+AE1(Pty1in) = 1 / 2 * Mass1 * V2(Pty1in) + Phi(Pty1in) - Phi0;
+AE1=AE1(AE1~=0);
+% if AddParticle, AE(end) = 0;end
+% AE = Phi;
+num0=length(AE0);
+num1=length(AE1);
+subplot(3, 2, 4),scatter3(x(1:num0),y(1:num0),AE0,ones(1,num0)*300,AE0,'fill');
+hold on
+subplot(3, 2, 4),scatter3(x(num0+1:end),y(num0+1:end),AE1,ones(1,num1)*100,AE1,'fill');
+axis(Limits);
+view(2);
+title('Atom Energy (KE + Pot)')
+hold off
+
+subplot(3, 2, 5), plot(x(Pty0in), y(Pty0in), 'bo', 'markers',...
+    MarkerSize, 'MarkerFaceColor', 'b');
+hold on
+subplot(3, 2, 5), plot(x(Pty1in), y(Pty1in), 'go', 'markers',...
+    MarkerSize, 'MarkerFaceColor', 'g');
+hold off
+title('Atoms (All)')
+xlabel('X')
+ylabel('Y')
+
+subplot(3, 2, 6), plot(time, T, 'k', 'linewidth', 2);
+hold on
+subplot(3, 2, 6), plot(time, T0, 'b', 'linewidth', 2);
+subplot(3, 2, 6), plot(time, T1, 'g', 'linewidth', 2);
+
+xlabel('time')
+ylabel('Temperature')
+title('Temperature: Each type and all')
+
+if doPlotImage == 1
+    if PlotCount == 1
+        f = getframe(PlotFig);
+        [im, map] = rgb2ind(f.cdata, 256, 'nodither');
+        im(1, 1, 1, 2) = 0;
+
+    else
+        f = getframe(PlotFig);
+        im(:, :, 1, PlotCount) = rgb2ind(f.cdata, map, 'nodither');
+    end
+end
+
+end
\ No newline at end of file