Mercurial > octave
changeset 29221:56978d4b266c stable
pcg.m: Return correct FLAG and correct RELRES output (bug #59776)
* pcg.m: Calculate relative residual correctly by dividing by norm (B).
Add BIST test for bug #59776.
author | Rik <rik@octave.org> |
---|---|
date | Wed, 30 Dec 2020 10:52:51 -0800 |
parents | 0180eaf55fd0 |
children | c7b27a4a6b43 8acb139a923c |
files | scripts/sparse/pcg.m |
diffstat | 1 files changed, 11 insertions(+), 1 deletions(-) [+] |
line wrap: on
line diff
--- a/scripts/sparse/pcg.m Tue Dec 22 10:48:56 2020 -0800 +++ b/scripts/sparse/pcg.m Wed Dec 30 10:52:51 2020 -0800 @@ -411,7 +411,7 @@ elseif (resvec (1, 1) == 0) relres = 0; else - relres = resvec(iter_min+1, 1) ./ resvec(1, 1); + relres = resvec(iter_min+1, 1) ./ b_norm; endif iter -= 2; # compatibility @@ -672,3 +672,13 @@ %! [x, flag, relres] = pcg (A, b, [], 2, M); %! assert (norm (b - A * x) / norm (b), relres, 8 * eps); +%!test <*59776> +%! A = [ 1.00000000 -0.00054274 -0.00066848; +%! -0.00054274 1.00000000 -0.00060330; +%! -0.00066848 -0.00060330 1.00000000]; +%! b = [1 1 1]'; +%! [x, flag, relres, iter, resvec] = pcg (A, b, 1e-6, 4, [], [], [1; 1; 1]); +%! assert (flag, 0); +%! assert (relres, resvec(2) / norm (b)); +%! assert (iter, 1); +