%% Function sloshPar
function plotData = shallowPar() 


%% Create a data input objects only on processor 1
if labindex == 1
    
    input = initialize();
    
    
    P = input.P;
    m = input.m;
    n = input.n;
    
    pass(P).m = m;
    % constant data
    for i = 1:P
        pass(i).m = m;
        pass(i).n = n+2;
        pass(i).b = input.b;
        pass(i).a = input.a;
        pass(i).mu = input.mu;
        pass(i).dt = input.dt;
        pass(i).dx = input.dx;
        pass(i).runtime = input.Tmax;
        pass(i).P = input.P;
        pass(i).xc = input.xc;
        pass(i).yc = input.yc;
    end

    % vertex data
    pass(1).X = input.X(:,1:n+3);
    pass(1).Y = input.Y(:,1:n+3);
    pass(1).u = input.u(:,1:n+3);
    pass(1).v = input.v(:,1:n+3);
    pass(1).M = input.M(:,1:n+3);

    for k = 2:input.P-1
        pass(k).X = input.X(:,(k-1)*n-1:k*n+1);
        pass(k).Y = input.Y(:,(k-1)*n-1:k*n+1);
        pass(k).u = input.u(:,(k-1)*n-1:k*n+1);
        pass(k).v = input.v(:,(k-1)*n-1:k*n+1);
        pass(k).M = input.M(:,(k-1)*n-1:k*n+1);
    end

    pass(P).X = input.X(:,(P-1)*n-1:m+1);
    pass(P).Y = input.Y(:,(P-1)*n-1:m+1);
    pass(P).u = input.u(:,(P-1)*n-1:m+1);
    pass(P).v = input.v(:,(P-1)*n-1:m+1);
    pass(P).M = input.M(:,(P-1)*n-1:m+1);

    % triangle data
    pass(1).m1 = input.m1(:,1:n+2);
    pass(1).m2 = input.m2(:,1:n+2);

    for k = 2:P-1
        pass(k).m1 = input.m1(:,(k-1)*n-1:k*n);
        pass(k).m2 = input.m2(:,(k-1)*n-1:k*n);
    end

    pass(P).m1 = input.m1(:,(P-1)*n-1:m);
    pass(P).m2 = input.m2(:,(P-1)*n-1:m);
end

%% End setup, start distribution to labs

if labindex > 1
    data = labReceive(1);
else   % in this case labindex == 1
    for i = 2:P
        labSend(pass(i), i);
    end
    data = pass(1);
end

%% Set output object for constants
if labindex == 1
    consts.M = input.M;
    consts.MM = input.MM;
    consts.xc = input.xc;
    consts.yc = input.yc;
    consts.m = input.m;
    consts.n = input.n;
    consts.plotpoints = input.Tmax*10; %this is not yet dynamic
else
    consts = 0;
end

%% Begin run on each individual processor

runtime = data.runtime;
dt = data.dt;
P = data.P;

step=.1/dt; %used to determine when to plot
k = 0;      %index
t = 0;      %time elapsed

%% Set up plotting variables

plotData(runtime*10).empty = 0; %this is also not yet dynamic

c = 0; %plot count

%% Running Loop
while t <= runtime 
    
    k = k + 1;
    t = t + dt;
    
%% Computation 
    data = spOneStep(data);
    
%% Communication with neighbors

    % Share data with adjacent processors
    n = data.n;

    if labindex == 1
        labSend([data.X(:,n-3:n-2) data.Y(:,n-3:n-2) data.u(:,n-3:n-2) ...
            data.v(:,n-3:n-2)],2);
    else if labindex == P
            fromLeft = labReceive(P - 1);
        else
            fromLeft = labSendReceive(labindex + 1, labindex - 1, ...
                [data.X(:,n-1:n) data.Y(:,n-1:n) data.u(:,n-1:n) ...
                data.v(:,n-1:n)]);
        end
    end

    if labindex == 1
        fromRight = labReceive(labindex + 1);
    else if labindex == P
            labSend([data.X(:,3) data.Y(:,3) data.u(:,3) data.v(:,3)],...
                P-1);
        else if labindex == 2
                fromRight = labSendReceive(1, 3, ...
                    [data.X(:,3:5) data.Y(:,3:5) data.u(:,3:5) ...
                    data.v(:,3:5)]);
            else
                fromRight = labSendReceive(labindex - 1, labindex + 1, ...
                    [data.X(:,3) data.Y(:,3) data.u(:,3) data.v(:,3)]);
            end
        end
    end

    % Update matrices with data from adjacent processors

    if labindex ~= 1
        data.X(:, 1) = fromLeft(:,1);
        data.Y(:, 1) = fromLeft(:,3);
        data.u(:, 1) = fromLeft(:,5);
        data.v(:, 1) = fromLeft(:,7);
        data.X(:, 2) = fromLeft(:,2);
        data.Y(:, 2) = fromLeft(:,4);
        data.u(:, 2) = fromLeft(:,6);
        data.v(:, 2) = fromLeft(:,8);
    end
    if labindex ~= P
        if labindex == 1
            data.X(:, n-1) = fromRight(:,1);
            data.Y(:, n-1) = fromRight(:,4);
            data.u(:, n-1) = fromRight(:,7);
            data.v(:, n-1) = fromRight(:,10);
            data.X(:, n) = fromRight(:,2);
            data.Y(:, n) = fromRight(:,5);
            data.u(:, n) = fromRight(:,8);
            data.v(:, n) = fromRight(:,11);
            data.X(:, n+1) = fromRight(:,3);
            data.Y(:, n+1) = fromRight(:,6);
            data.u(:, n+1) = fromRight(:,9);
            data.v(:, n+1) = fromRight(:,12);
        else
            data.X(:, n+1) = fromRight(:,1);
            data.Y(:, n+1) = fromRight(:,2);
            data.u(:, n+1) = fromRight(:,3);
            data.v(:, n+1) = fromRight(:,4);
        end
    end
    
%% Plotting every 1/10 second 

    jj=k-step*floor(k/step);
    
    if jj == 0

        c = c + 1
        
        if labindex == 1

            plotData(c).X = data.X(:,1:end-3); 
            plotData(c).Y = data.Y(:,1:end-3); 
            plotData(c).u = data.u(:,1:end-3); 
            plotData(c).v = data.v(:,1:end-3); 

            plotData(c).A1 = data.A1(:,1:end-3);
            plotData(c).A2 = data.A2(:,1:end-3); 

        else if labindex == P

                plotData(c).X = data.X(:,3:end); 
                plotData(c).Y = data.Y(:,3:end); 
                plotData(c).u = data.u(:,3:end); 
                plotData(c).v = data.v(:,3:end); 

                plotData(c).A1 = data.A1(:,2:end);
                plotData(c).A2 = data.A2(:,2:end); 

            else

                plotData(c).X = data.X(:,3:end-1);
                plotData(c).Y = data.Y(:,3:end-1); 
                plotData(c).u = data.u(:,3:end-1); 
                plotData(c).v = data.v(:,3:end-1); 

                plotData(c).A1 = data.A1(:,2:end-1);
                plotData(c).A2 = data.A2(:,2:end-1); 

            end
        end
        
        plotData(c).t = t;
        
    end

%% Loop

%% Assign objects into 'base' workspace for easy transfer to client
    assignin('base', 'out', plotData);
    assignin('base', 'static', consts);
    

end