Why is my custom Conjugate gradient solver and MATLAB pcg() do not have the same result?

6 views (last 30 days)
Hello
Based on some literature I tried to implement a C code based on this Matlab code for conjugate gradient solver. However, When I tried to compare the code published on Wikipedia with MATLAB pcg() function as follows:
rng default
A = sprand(400,400,.5);
A = A'*A;
b = sum(A,2);
x1 = pcg(A, b, 1e-4, 150);
x2 = conj_grad(A, b, 1e-4, 150);
cmp = abs(x1-x2) < 1e6*eps(min(abs(x1),abs(x2))); % to compare
function x_i = conj_grad(m_A, v_b, tol, iter)
x_i = v_b;
r = v_b - m_A * x_i;
p = r;
rsold = r' * r;
for i = 1:iter
Ap = m_A * p;
alpha = rsold / (p' * Ap);
x_i = x_i + alpha * p;
r = r - alpha * Ap;
rsnew = r' * r;
if sqrt(rsnew) < tol
fprintf('cg converged at iteration %d to a solution with relative residual %f', i, sqrt(rsnew));
break;
end
p = r + (rsnew / rsold) * p;
rsold = rsnew;
end
end
I discovered that the custom made function does not converge as pcg() and the results are different for example:
pcg() conj_grad()
---------------------------------------------
0.9906 0.9959
1.0118 1.0056
0.9989 0.9999
1.0201 0.9970
1.0049 1.0139
......
Could someone give some hints please?
Thank you

Answers (0)

Products


Release

R2021a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!