cgm = conjugate_gradient_method(A, b, x0, 1e-6, 1000);
pcgm = preconditioned_conjugate_gradient_method(A, M, b, x0, 1e-6, 1000);
diff = norm(cgm - pcgm);
cgm = conjugate_gradient_method(A, b, x0, 1e-6, 1000);
pcgm = preconditioned_conjugate_gradient_method(A, M, b, x0, 1e-6, 1000);
diff = norm(cgm - pcgm);