Skip to content

Commit

Permalink
modified: solve/alstpz_solve.m;
Browse files Browse the repository at this point in the history 
	modified:   solve/amen_solve2.m;
	modified:   solve/dmrg_solve3.m:
Added automatic detection of absent solve3d_2, and turned MEX on by default;
More fair distribution of core sizes for MEX solver in dmrg_solve3
  • Loading branch information
Sergey Dolgov committed May 28, 2014
1 parent 62814cb commit f7927bb
Show file tree
Hide file tree
Showing 3 changed files with 78 additions and 27 deletions.
36 changes: 24 additions & 12 deletions solve/alstpz_solve.m
View file
Original file line number Diff line number Diff line change
Expand Up @@ -70,6 +70,7 @@
kickrank = 4;
x=[];
symm = false;
verb = 3;

for i=1:2:length(varargin)-1
switch lower(varargin{i})
Expand All @@ -95,12 +96,20 @@
resid_damp = varargin{i+1};
case 'symm'
symm=varargin{i+1};
case 'verb'
verb=varargin{i+1};

otherwise
error('Unrecognized option: %s\n',varargin{i});
end
end

% Disable MEX if it does not exist
if (ismex)&&(exist('solve3d_2', 'file')<2)
warning('MEX local solver is not found, disabled');
ismex = false;
end;

% Extract sizes and TT blocks
d = y.d;
n = A.n;
Expand All @@ -122,9 +131,9 @@

% Test output data -- e.g. for convergence tracking
testdata = cell(3,1);
testdata{2} = cell(nswp, d);
testdata{1} = zeros(nswp, d);
testdata{3} = zeros(nswp, 1);
testdata{2} = cell(d, nswp);
testdata{1} = zeros(d, nswp);
testdata{3} = zeros(1, nswp);
t_corr_amen = tic;

% ALS(t+z) sweeps
Expand Down Expand Up @@ -192,12 +201,12 @@
t_old = toc(t_corr_amen);
% Update via the ALS iteration
if (symm)
[x,td_als] = als_solve(A2, Ay, tol, 'max_full_size', max_full_size, 'x0', x, 'verb', 3, 'local_prec', local_prec, 'local_iters', local_iters, 'local_restart', local_restart, 'resid_damp', resid_damp, 'ismex', ismex);
[x,td_als] = als_solve(A2, Ay, tol, 'max_full_size', max_full_size, 'x0', x, 'verb', verb, 'local_prec', local_prec, 'local_iters', local_iters, 'local_restart', local_restart, 'resid_damp', resid_damp, 'ismex', ismex);
else
[x,td_als] = als_solve(A, y, tol, 'max_full_size', max_full_size, 'x0', x, 'verb', 3, 'local_prec', local_prec, 'local_iters', local_iters, 'local_restart', local_restart, 'resid_damp', resid_damp, 'ismex', ismex);
[x,td_als] = als_solve(A, y, tol, 'max_full_size', max_full_size, 'x0', x, 'verb', verb, 'local_prec', local_prec, 'local_iters', local_iters, 'local_restart', local_restart, 'resid_damp', resid_damp, 'ismex', ismex);
end;
testdata{1}(swp,:) = t_old+td_als{1}(:,1);
testdata{2}(swp,:) = td_als{2};
testdata{1}(:, swp) = t_old+td_als{1}(:,1);
testdata{2}(:, swp) = td_als{2};
testdata{3}(swp) = max(td_als{3});
end;

Expand Down Expand Up @@ -267,7 +276,7 @@
% This is some convergence output for test purposes
td = cell(3,1);
td{1} = zeros(d, 1);
td{2} = cell(1, d);
td{2} = cell(d, 1);
td{3} = zeros(d, 1);

t_amr_solve = tic;
Expand Down Expand Up @@ -396,11 +405,14 @@
sol = reshape(sol, rx(i), n(i), rx(i+1));
crx{i} = sol;
end;

if (verb>2)&&(i==d)
if (verb>2)
td{1}(i) = toc(t_amr_solve);
x = cell2core(x, crx); % for test
td{2}{i} = x;
if (verb>3)||(i==d) % each microstep is returned only if really asked for
% Otherwise the memory will blow up
x = cell2core(x, crx); % for test
td{2}{i} = x;
end;
end;

end;
Expand Down
17 changes: 13 additions & 4 deletions solve/amen_solve2.m
View file
Original file line number Diff line number Diff line change
Expand Up @@ -87,7 +87,7 @@

tol_exit = [];

ismex = false;
ismex = true;

% kicktype = 'svd';
kicktype = 'als';
Expand Down Expand Up @@ -203,6 +203,12 @@
if (strcmp(local_prec, 'rjacobi')); local_prec_char = 3; end;
% if (strcmp(trunc_norm, 'fro')); trunc_norm_char = 0; end;

% Disable MEX if it does not exist
if (ismex)&&(exist('solve3d_2', 'file')<2)
warning('MEX local solver is not found, disabled');
ismex = false;
end;

if (A.n~=A.m)
error(' AMEn does not know how to solve rectangular systems!\n Use amen_solve2(ctranspose(A)*A, ctranspose(A)*f, tol) instead.');
end;
Expand Down Expand Up @@ -756,10 +762,13 @@
crx{i} = sol;
end;

if (verb>2) % &&(i==d)
if (verb>2)
testdata{1}(i,swp) = toc(t_amen_solve);
x = cell2core(x, crx)*exp(sum(log(nrmsx))); % for test
testdata{2}{i,swp} = x;
if (verb>3)||(i==d) % each microstep is returned only if really asked for
% Otherwise the memory will blow up
x = cell2core(x, crx)*exp(sum(log(nrmsx))); % for test
testdata{2}{i,swp} = x;
end;
end;

end;
Expand Down
52 changes: 41 additions & 11 deletions solve/dmrg_solve3.m
View file
Original file line number Diff line number Diff line change
Expand Up @@ -65,7 +65,7 @@
local_solver = 'gmres';
% local_solver = 'pcg';

ismex = false;
ismex = true;

% dirfilter = -1; % only backward
% dirfilter = 1; % only forward
Expand Down Expand Up @@ -137,6 +137,12 @@
error(' DMRG does not know how to solve rectangular systems!\n Use dmrg_solve3(ctranspose(A)*A, ctranspose(A)*f, tol) instead.');
end;

% Disable MEX if it does not exist
if (ismex)&&(exist('solve3d_2', 'file')<2)
warning('MEX local solver is not found, disabled');
ismex = false;
end;

d = y.d;
n = A.n;
if (isempty(x))
Expand Down Expand Up @@ -171,7 +177,7 @@

somedata = cell(3,1);
somedata{1}=zeros(d,nswp*2);
somedata{2}=cell(nswp*2, d);
somedata{2}=cell(d, nswp*2);
somedata{3}=zeros(d,nswp*2);

t_dmrg_solve = tic;
Expand Down Expand Up @@ -407,12 +413,33 @@
sol = sol_prev + dsol;

else % ismex
Amex = reshape(A1, ra(i)*n(i)*n(i), ra(i+1));
Amex = Amex*reshape(A2, ra(i+1), n(i+1)*n(i+1)*ra(i+2));
Amex = reshape(Amex, ra(i), n(i), n(i), n(i+1), n(i+1), ra(i+2));
Amex = permute(Amex, [1, 2, 4, 3, 5, 6]);
Amex = reshape(Amex, ra(i), n(i)*n(i+1), n(i)*n(i+1), ra(i+2));
sol = solve3d_2(permute(Phi1,[1,3,2]), Amex, permute(Phi2, [3,2,1]), rhs, real_tol, 1, sol_prev, local_prec_char, local_restart, local_iters, max(verb-1,0));
% Since we only have MEX for 3 blocks, we have to merge some
if (rx(i)*n(i)<=min(n(i)*n(i+1), n(i+1)*rx(i+2)))
Phi1mex = reshape(permute(Phi1, [1, 3, 2]), rx(i)*rx(i), ra(i));
Phi1mex = Phi1mex*reshape(A1, ra(i), n(i)*n(i)*ra(i+1));
Phi1mex = reshape(Phi1mex, rx(i), rx(i), n(i), n(i), ra(i+1));
Phi1mex = permute(Phi1mex, [1, 3, 2, 4, 5]);
Phi1mex = reshape(Phi1mex, rx(i)*n(i), rx(i)*n(i), ra(i+1));
Amex = reshape(A2, ra(i+1), n(i+1), n(i+1), ra(i+2));
Phi2mex = permute(Phi2, [3,2,1]);
elseif (n(i)*n(i+1)<=min(rx(i)*n(i), n(i+1)*rx(i+2)))
Amex = reshape(A1, ra(i)*n(i)*n(i), ra(i+1));
Amex = Amex*reshape(A2, ra(i+1), n(i+1)*n(i+1)*ra(i+2));
Amex = reshape(Amex, ra(i), n(i), n(i), n(i+1), n(i+1), ra(i+2));
Amex = permute(Amex, [1, 2, 4, 3, 5, 6]);
Amex = reshape(Amex, ra(i), n(i)*n(i+1), n(i)*n(i+1), ra(i+2));
Phi1mex = permute(Phi1, [1,3,2]);
Phi2mex = permute(Phi2, [3,2,1]); % rx', ra, rx -> rx,ra,rx'
else
Phi1mex = permute(Phi1, [1,3,2]);
Amex = reshape(A1, ra(i), n(i), n(i), ra(i+1));
Phi2mex = reshape(A2, ra(i+1)*n(i+1)*n(i+1), ra(i+2));
Phi2mex = Phi2mex*reshape(permute(Phi2, [2, 1, 3]), ra(i+2), rx(i+2)*rx(i+2));
Phi2mex = reshape(Phi2mex, ra(i+1), n(i+1), n(i+1), rx(i+2), rx(i+2));
Phi2mex = permute(Phi2mex, [3,5, 1, 2,4]);
Phi2mex = reshape(Phi2mex, n(i+1)*rx(i+2), ra(i+1), n(i+1)*rx(i+2));
end;
sol = solve3d_2(Phi1mex, Amex, Phi2mex, rhs, real_tol, 1, sol_prev, local_prec_char, local_restart, local_iters, 0);
flg = 0;
iter = 1;
end;
Expand Down Expand Up @@ -590,10 +617,13 @@
crx{i} = u;
crx{i+1} = v;

if (verb>2) % &&(i==d-1)
if (verb>2)
somedata{1}(i,(swp-1)*2+1.5-dir/2) = toc(t_dmrg_solve);
x = cell2core(tt_tensor, crx);
somedata{2}{(swp-1)*2+1.5-dir/2, i}=x;
if (verb>3)||(i==d-1) % each microstep is returned only if really asked for
% Otherwise the memory will blow up
x = cell2core(tt_tensor, crx);
somedata{2}{i, (swp-1)*2+1.5-dir/2}=x;
end;
end;

i = i+dir;
Expand Down

0 comments on commit f7927bb

Please sign in to comment.