constr_precond_inner

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

0001 function expB = constr_precond_inner( A, X, mu )
0002 
0003     %   TTeMPS Toolbox.
0004     %   Michael Steinlechner, 2013-2016
0005     %   Questions and contact: michael.steinlechner@epfl.ch
0006     %   BSD 2-clause license, see LICENSE.txt
0007 
0008     n = size(A.L0, 1);
0009     sz = [X.rank(mu), X.size(mu), X.rank(mu+1)];
0010 
0011     B1 = zeros( X.rank(mu) );
0012     % calculate B1 part:
0013     for i = 1:mu-1
0014         % apply L to the i'th core
0015         tmp = X;
0016         tmp.U{i} = tensorprod_ttemps( tmp.U{i}, A.L0, 2 );
0017         B1 = B1 + innerprod( X, tmp, 'LR', mu-1);
0018     end
0019 
0020     % calculate B2 part:
0021     B2 = A.L0;
0022 
0023     B3 = zeros( X.rank(mu+1) );
0024     % calculate B3 part:
0025     for i = mu+1:A.order
0026         tmp = X;
0027         tmp.U{i} = tensorprod_ttemps( tmp.U{i}, A.L0, 2 );
0028         B3 = B3 + innerprod( X, tmp, 'RL', mu+1);
0029     end
0030     
0031     [V1,e1] = eig(B1);
0032     e1 = diag(e1);
0033     [V3,e3] = eig(B3);
0034     e3 = diag(e3);
0035 
0036     lmin = min(e1) + min(A.E_L) + min(e3);
0037     lmax = max(e1) + max(A.E_L) + max(e3);
0038 
0039     R = lmax/lmin;
0040     
0041     [omega, alpha] = load_coefficients( R );
0042 
0043     k = length(omega);
0044     omega = omega/lmin;
0045     alpha = alpha/lmin;
0046 
0047     expB = cell(3,k);
0048     
0049     for i = 1:k
0050         expB{1,i} = omega(i) * V1*diag( exp( -alpha(i)*e1 ))*V1';    % include omega in first part
0051         expB{2,i} = A.V_L*diag( exp( -alpha(i)*A.E_L ))*A.V_L';
0052         expB{3,i} = V3*diag( exp( -alpha(i)*e3 ))*V3';
0053     end
0054 end
0055