function [ x_res, vec_rhos, sqrt_vec_rhos ] = mine_pcg( A, b, L, LT, stop_threshold, num_iter )

n = length(A);
x = zeros(n,1);
r = b;

stop_rho = stop_threshold * stop_threshold * (b' *b);
z = LT\(L\r);
rho = r'* z;

vec_rhos = [];sqrt_vec_rhos = [];
vec_rhos = [vec_rhos, rho];
sqrt_vec_rhos = [sqrt_vec_rhos, sqrt(r'*z)];
p = 0;
z = 0;
for k = 1:num_iter
    % whether to stop
    if sqrt(rho) < sqrt(stop_rho)
        break;
    end
    
    if k == 1
        p = LT\(L\r);
    else
        p = z + (vec_rhos(k)/vec_rhos(k-1))*p;
    end
    
    w = A * p;
    alpha = vec_rhos(k)/(p'*w);
    x = x + alpha * p;
    r = r - alpha * w;
    z = LT\(L\r);
    rho = r'* z;
    vec_rhos = [vec_rhos, rho];
    sqrt_vec_rhos = [sqrt_vec_rhos, sqrt(r'*r)];
end

x_res = x;

end