shift

   TTeMPS Toolbox. 
   Michael Steinlechner, 2013-2016
   Questions and contact: michael.steinlechner@epfl.ch
   BSD 2-clause license, see LICENSE.txt

0001 %   TTeMPS Toolbox.
0002 %   Michael Steinlechner, 2013-2016
0003 %   Questions and contact: michael.steinlechner@epfl.ch
0004 %   BSD 2-clause license, see LICENSE.txt
0005 function x = shift( x, nu, tol, maxrank )
0006     
0007     r = x.rank;
0008     n = x.size;
0009     p = x.p;
0010     mu = x.mu;
0011 
0012     if ~exist('tol', 'var')
0013         tol = eps;
0014     end
0015     if ~exist('maxrank', 'var')
0016         maxrank = inf;%r(mu+1);
0017     end
0018 
0019 
0020     if mu == nu-1   
0021         % shift block one to the right
0022         U = permute( x.U{mu}, [1, 2, 4, 3] );
0023         U = reshape( U, [r(mu)*n(mu), p*r(mu+1)] );
0024 
0025         [U,S,V] = svd( U, 'econ' );
0026         if p == 1 
0027             s = length(diag(S));
0028         else
0029             s = trunc_singular( diag(S), tol );
0030         end
0031         if length(diag(S)) >= s+1
0032             disp(['cut singular value of rel. magnitude (s_{i+1}/s_1): ', ...
0033                         num2str(S(s+1,s+1)/S(1,1))])
0034         end
0035         U = U(:,1:s);
0036         x.U{mu} = reshape( U, [r(mu), n(mu), s] );
0037         W = S(1:s,1:s)*V(:,1:s)';
0038         W = reshape( W, [s, p, r(mu+1)]);
0039         W = permute( W, [1, 3, 2]);
0040         
0041         C = zeros( [s, n(nu), r(nu+1), p] ); 
0042         for k = 1:p
0043             C(:,:,:,k) = tensorprod_ttemps( x.U{nu}, W(:,:,k), 1);
0044         end
0045         
0046         x.U{nu} = C;
0047         x.mu = nu;
0048 
0049     elseif x.mu == nu+1   
0050         % shift block one to the left
0051         V = permute( x.U{mu}, [1, 4, 2, 3] );
0052         V = reshape( V, [r(mu)*p, n(mu)*r(mu+1)] );
0053 
0054         [U,S,V] = svd( V, 'econ' );
0055         if p == 1
0056             s = length(diag(S));
0057         else
0058             s = trunc_singular( diag(S), tol );
0059         end
0060         if length(diag(S)) >= s+1
0061             disp(['cut singular value of rel. magnitude (s_{i+1}/s_1): ', ...
0062                         num2str(S(s+1,s+1)/S(1,1))])
0063         end
0064         V = V(:,1:s)';
0065         x.U{mu} = reshape( V, [s, n(mu), r(mu+1)] );
0066 
0067         W = U(:,1:s)*S(1:s,1:s);
0068         W = reshape( W, [r(mu), p, s]);
0069         W = permute( W, [1, 3, 2]);
0070 
0071         C = zeros( [r(nu), n(nu), s, p] ); 
0072         for k = 1:p
0073             C(:,:,:,k) = tensorprod_ttemps( x.U{nu}, W(:,:,k)', 3);
0074         end
0075         
0076         x.U{nu} = C;
0077         x.mu = nu;
0078     else
0079         error('Can only shift the superblock one core left or right')
0080     end
0081 
0082 
0083 end