Reputation: 1
I'm trying to implement TN with preconditioning on MATLAB, and evaluate it for different functions, with different starting points. I'm also trying to confront it with TN without preconditioning.
What happens is that for one function (Generalized Broyden tridiagonal function) and some initial points, TN with pre is better than TN without pre in terms of execution time and number of iterations (with the same result). But with other two functions (a Variably dimensioned function and a Penalty function), TN with pre is worse than TN without pre. I suppose I made some mistakes in my preconditioning code.
function xk = cg_curvtrun_newt_pre( A, b, x0 , tol)
% INPUTS :
% A = matrix of the coefficients of the lienar system
% b = vector of constants
% x0 = n- dimensional column vector ;.
% tol = tolerance on the norm of the relative residual
% OUTPUT :
% xk = the last x computed by the function ;
% preconditioning operations
%either use the ichol when possible or diag ( diag (A)) as preconditioner
try
L = ichol (sparse(A));
catch exception
L = diag ( diag(A));
end
%L = diag ( diag (A));
L =L'*L;
b=b';
L_i = L \eye( size (L ,1) ); % compute the inverse
A = L_i * A; %new matrix
b = L_i * b; % new vector
% initilizations
xk = x0;
rk = b - A * xk;
pk = rk;
norm_b = norm (b);
relres = norm (rk) / norm_b ; % relative residual
kmax = 1000;
k = 0;
% compute pcg iterations until the relative residual gets smaller
% than the tolerance or you reach the max number of iterations
while relres > tol && k < kmax
if pk' * A * pk <= 0 %if the negative curvature condition is satisfied
if k == 0 % if its the first iteration , compute the new xk and stop
zk = A * pk;
alphak = (rk' * pk )/ (pk' * zk);
xk = xk + alphak * pk;
end
break
% stop computing xk and return it
end
zk = A * pk;
alphak = (rk' * pk )/ (pk' * zk);
xk = xk + alphak * pk;
rk = rk - alphak * zk;
betak = -(rk' * zk) / (rk' * zk);
pk = rk + betak * pk;
relres = norm (rk) / norm_b ;
k = k + 1;
end
end
Upvotes: 0
Views: 31