%SYMKTENSOR Class for storing symmetric Kruskal tensor (decomposed).
%
%SYMKTENSOR Methods:
%   arrange     - Arranges the rank-1 components of a symktensor.
%   disp        - Command window display for a symktensor.
%   display     - Command window display for a symktensor.
%   double      - Convert a symktensor to a double array.
%   end         - Last index of indexing expression for symktensor.
%   entry       - Extract a single entry from a symktensor.
%   f_implicit  - Function F(M)=||X-M||^2 for two symktensors.
%   fg          - Master objective function for optimization of symmetric Kruskal model.
%   fg_explicit - Function and gradient of F(M)=||X-M||^2 for symktensor.
%   fg_implicit - Function and gradient of F(M)=||X-M||^2 for two symktensors.
%   fg_setup    - Setup for optimization of symmetric Kruskal model.
%   full        - Convert a symktensor to a symtensor.
%   g_implicit  - Gradient of F(M)=||X-M||^2 for two symktensors.
%   isequal     - True if each component of two symktensors is numerically equal.
%   isscalar    - False for symktensors.
%   issymmetric - Rhetorical function for a symktensor.
%   mtimes      - Implement A*B (scalar multiply) for symktensor.
%   ncomponents - Number of components for a symktensor.
%   ndims       - Number of modes for a symktensor.
%   norm        - Frobenius norm of a symktensor.
%   normalize   - Normalizes the columns of the factor matrix.
%   permute     - Permute dimensions of a symktensor.
%   randextract - Create new symktensor with random subset of components.
%   score       - Checks if two symktensors match except for permutation.
%   size        - Size of symktensor.
%   subsasgn    - Subscripted assignment for symktensor.
%   subsref     - Subscripted reference for a symktensor.
%   symktensor  - Tensor stored as a symmetric Kruskal operator (decomposed).
%   tovec       - Convert symktensor to vector representation.
%   ttsv        - Symktensor times same vector in multiple modes.
%   uminus      - Unary minus for symktensor. 
%   update      - Convert vector to symktensor.
%   uplus       - Unary plus for a symktensor. 
%
%   <a href="matlab:web(strcat('file://',fullfile(getfield(what('tensor_toolbox'),'path'),'doc','html','symktensor_doc.html')))">Documentation page for symmetric Kruskal tensor class</a>
%
%   See also TENSOR_TOOLBOX, SYMTENSOR
%
%Tensor Toolbox for MATLAB: <a href="https://www.tensortoolbox.org">www.tensortoolbox.org</a>

classdef symktensor

    properties
        lambda %< Vector of weights for each component.
        u      %< Factor matrix.
        m      %< Number of modes.
    end

    methods
        function t = symktensor(varargin)
            %SYMKTENSOR Tensor stored as a symmetric Kruskal operator (decomposed).
            %
            %   A symmetric Kruskal tensor is used to build a model of a M-way
            %   N-dimensional symmetric tensor. We have to specify the number of
            %   components, denoted by P. Each component comprises a scalar weight and
            %   the M-way symmetric outer product of a vector. We store the P weights
            %   in the length-P vector LAMBDA and the vectors in the NxP matrix X. The
            %   number of modes, M, is stored as a scalar.
            %
            %   S = SYMKTENSOR(LAMBDA,U,M) creates a symmetric M-way Kruskal tensor
            %   from its constituent parts. Here lambda is a K-vector and each U is a  
            %   matrix with K columns.
            %
            %   S = SYMKTENSOR(K) creates a symmetric Kruskal tensor by "symmetrizing"
            %   a ktensor, i.e., averaging the constituent factor matrices and taking
            %   care to get the signs aligned. 
            %
            %   S = SYMKTENSOR(S0) creates a symktensor by copying an existing
            %   symktensor. 
            %
            %   *** Below are specialized constructore for use in optimization. There
            %   are not recommended for general use. ***
            %
            %   S = SYMKTENSOR(P,A) creates a symmetric Krukal tensor with P components
            %   that is sized to match the tensor A. The LAMBDA is set to be all ones,
            %   and the X is initialized using to uniform random values in [0,1].
            %
            %   S = SYMKTENSOR(V,A) creates a symmetric Kruskal tensor from a
            %   vectorized verion V. The second argument is a tensor of the size that
            %   is being modeled.
            %
            %   S = SYMKTENSOR(V,A,NOLAMBDA) creates a symmetric Kruskal tensor from a
            %   vectorized version V. The second argument is a tensor of the size that
            %   is being modeled. The third argument indicates whether or not LAMBDA is
            %   stored in V. If NOLAMBDA=TRUE, then only X is stored (i.e., LAMBDA
            %   values are set to 1).  
            %
            %   S = SYMKTENSOR(V,M,P) creates an M-way symmetric Kruskal
            %   tensor with P components from a vectorized version V.  
            %
            %   S = SYMKTENSOR(V,M,P,NOLAMBDA) uses the NOLAMBDA parameter as described
            %   above. 
            %
            %   Examples
            %   LAMBDA = randn(2, 1); %<-- Lambda vector (should be a column)
            %   U = randn(4, 2); %<-- Factor matrix
            %   S = symktensor(LAMDBA, U, 3) %<--Declare symktensor 
            %
            %   See also KTENSOR, TOVEC
            %
            %Tensor Toolbox for MATLAB: <a href="https://www.tensortoolbox.org">www.tensortoolbox.org</a>

            % EMPTY CONSTRUCTOR
            if nargin == 0
                t.lambda = [];
                t.u = [];
                t.m = 0;
                return;
            end

            % Copy CONSTRUCTOR
            if (nargin == 1) && isa(varargin{1}, 'symktensor')
                t.lambda = varargin{1}.lambda;
                t.u = varargin{1}.u;
                t.m = varargin{1}.m;
                return;
            end

            % Symmetrize a KTENSOR
            if (nargin == 1) && isa(varargin{1}, 'ktensor')
                K = varargin{1};
                if ~issymmetric(K)
                     K = symmetrize(K);
                end
                t.lambda = K.lambda;
                t.u = K.U{1};
                t.m = length(K.U);
                return;
            end

            % Create a random SYMKTENSOR from the number of components and a tensor.
            if (nargin == 2) && isscalar(varargin{1})
                P = varargin{1};
                A = varargin{2};
                if  ~isa(A,'tensor') && ~isa(A,'symtensor')
                    error('Second argument should be a tensor when the first argument is a scalar');
                end
                M = ndims(A);
                N = size(A,1);
                t.lambda = ones(P,1); %Column vector
                t.u = rand(N,P);   %Robert development note: Document this better
                t.m = M;
                return;
            end

            % Convert from vector representation, given tensor to determine shape
            if nargin == 2
               v = varargin{1};
               A = varargin{2}; 
               if ~isvector(v) || ~(isa(varargin{2},'tensor') || isa(varargin{2},'symtensor'))
                   error('Wrong arguments');
               end
               n = size(A,1);
               p = length(v) / (n+1);
               if round(p) ~= p
                   error('Length of input vector is not evenly divisible by (n+1)');
               end
               t.lambda = v(1:p); 
               u_vec = v(p+1:end); % Renamed to avoid conflict with property t.u
               t.u = reshape(u_vec,n,p);
               t.m = ndims(A);
               return;
            end

            % Convert from vector representation, given tensor to determine shape
            if (nargin == 3) && (isa(varargin{2},'tensor') || isa(varargin{2},'symtensor'))
                
                tf = varargin{3}; %No lambda option
                
                if ~tf
                    % Call constructor again with different signature
                    t_temp = symktensor(varargin{1},varargin{2});
                    t.lambda = t_temp.lambda;
                    t.u = t_temp.u;
                    t.m = t_temp.m;
                    return;
                end
                
               v = varargin{1};
               A = varargin{2};
               if ~isvector(v)
                   error('Wrong arguments');
               end
               n = size(A,1);
               p = length(v) / n;
               if round(p) ~= p
                   error('Length of input vector is not evenly divisible by n');
               end
               t.lambda = ones(p,1); %Column vector
               u_vec = v; % Renamed to avoid conflict with property t.u
               t.u = reshape(u_vec,n,p);
               t.m = ndims(A);
               return;
            end

            % Convert from vector representation, given sizes
            if (nargin == 3) && isscalar(varargin{2})
                v = varargin{1};
                m_val = varargin{2}; % Renamed to avoid conflict with property t.m
                p = varargin{3};
                t.lambda = v(1:p); 
                u_vec = v(p+1:end); % Renamed to avoid conflict with property t.u
                t.u = reshape(u_vec,[],p);
                t.m = m_val;
                return;
            end

            % Convert from vector representation, given sizes and nolambda!
            if (nargin == 4) 
                tf = varargin{4};
                
                if ~tf
                    % Call constructor again with different signature
                    t_temp = symktensor(varargin{1},varargin{2},varargin{3});
                    t.lambda = t_temp.lambda;
                    t.u = t_temp.u;
                    t.m = t_temp.m;
                    return;
                end
                
                v = varargin{1};
                m_val = varargin{2}; % Renamed to avoid conflict with property t.m
                p = varargin{3};
                t.lambda = ones(p,1); %Column vector
                u_vec = v; % Renamed to avoid conflict with property t.u
                t.u = reshape(u_vec,[],p);
                t.m = m_val;
                return;
            end

            if nargin ~= 3
                error('Check arguments to create symktensor');
            end

            t.lambda = varargin{1};
            t.u = varargin{2};
            t.m = varargin{3};

            if ~isa(t.lambda,'numeric') || ~iscolumn(t.lambda)
                error('LAMBDA must be a column vector.'); 
            end
                
            % Check that each Um is indeed a matrix
            if ~ismatrix(t.u) 
                error(['Matrix U is not a matrix!']);
            end

            % Size error checking			     
            k = length(t.lambda); 
            if  size(t.u,2) ~= k
                   error(['Matrix U does not have ' int2str(k) ' columns.']);
            end
        end

        function b = saveobj(a)
            %SAVEOBJ Save a symktensor object.
            %
            %   B = SAVEOBJ(A) is called by SAVE when a symktensor object is
            %   saved. The result B is a struct containing the symktensor data.
            %
            %   See also SYMKTENSOR/LOADOBJ, SAVE, LOAD.
            %
            %Tensor Toolbox for MATLAB: <a href="https://www.tensortoolbox.org">www.tensortoolbox.org</a>
            
            b.lambda = a.lambda;
            b.u = a.u;
            b.m = a.m;
        end

    end % methods

    methods (Static)
        function t = loadobj(s)
            %LOADOBJ Load a symktensor object.
            %
            %   T = LOADOBJ(S) is called by LOAD when a symktensor object is
            %   loaded. S is a structure containing the symktensor data and T
            %   is the restored symktensor object. If S is an object, then the
            %   object is simply returned.
            %
            %   See also SYMKTENSOR/SAVEOBJ, SAVE, LOAD.
            %
            %Tensor Toolbox for MATLAB: <a href="https://www.tensortoolbox.org">www.tensortoolbox.org</a>
            
            if isa(s,'symktensor')
                t = s;
            else
                % s is a struct, create a new symktensor object
                % using the standard (lambda, u, m) constructor signature.
                t = symktensor(s.lambda, s.u, s.m);
            end
        end
    end % static methods
end
