function varargout = scalaEig(A, isvec)
% SCALAEIG  ScaLAPACK symmetric eigenvalue solver

%   Copyright 2006-2020 The MathWorks, Inc.

nargoutchk(0, 3);

% When there are no output arguments, the function should still return
% the eigenvalues.
numOutputs = max(1, nargout);

% Redistribute input according to the default 2dbc distribution scheme.
if ~isa(getCodistributor(A), 'codistributor2dbc')
    A = redistribute(A, codistributor2dbc());
end

aDist = getCodistributor(A);
lbgrid = aDist.LabGrid;

localA = hGetLocalPartWithCorrectComplexity(A);

% Create array descriptor required by ScaLAPACK
descA = arraydescriptor(aDist);

% Output from ScalaPack is ordered with eigenvalues first, then
% eigenvectors. Third output, if requested, is created later on in this
% file and we only request a max of 2 outputs from ScaLAPACK.
[eigOut{1:min(2, numOutputs)}] = parallel.internal.codistextern.denseEig(localA, ...
    descA, lbgrid(1), lbgrid(2), aDist.Orientation, isreal(A));

% We always need the values.
% The ScaLAPACK output is replicated but we want codistributed output.
eigenVals = codistributed.pConstructFromReplicated(eigOut{1}, codistributor1d(1));
if ~isvec
    eigenVals = diag(eigenVals);
end

if (numOutputs == 1)
    varargout{1} = eigenVals;
else
    % When there are multiple output arguments, they need to be rearranged.  The
    % ScaLAPACK ordering is [evals, evecs] while MATLAB expects [evecs, evals]
    codistr = codistributor2dbc(lbgrid, aDist.BlockSize, aDist.Orientation, size(A));
    A = codistributed.pDoBuildFromLocalPart(eigOut{2}, codistr);
    varargout{1} = redistribute(A, codistributor1d());
    varargout{2} = eigenVals;
    % The 3rd output (left eigenvectors) in the symmetric case is identical
    % to the 1st output (right eigenvectors).
    if nargout==3
        varargout{3} = varargout{1};
    end
end
end
