[Math] Need matlab help related with the iteration method

MATLAB

I am reading an iteration method for computing the Moore- Penrose genralized inverse
of a given matrix $A$, which is given as follows:

$X_{k+1} = (1+\beta)X_{k} – \beta X_{k} A X_{k}$

where $X_{k}$, k = 0,1,… is a sequence of approximations for computing Moore- Penrose genralized inverse

$ X_{0} = \beta A' $ is the initial approximation , $0<\beta\leq 1$ and $A'$ is the transpose of matrix $A$

$d_{k} = \|X_{k+1} – X_{k}\|_{fro}$ is the error matrix norm (frobenius norm)

I have made following matlab program for computing Moore- Penrose genralized inverse by above mentioned method. But i am unable to make code for stopping criterion which says that.

perform the iteration untill

$|d_{k+1}/d_{k} – \beta -1|> 10^{-4}$

enter image description here

Please help me with this. I would be very much thankful to you.

Best Answer

The prep before your loop should stay the same. The appropriate script is

A = ...; % as you have given
beta = ...; % whatever you want
X0 = beta*A'; % calculate initial estimate

% (these initial values may need to be changed, I don't have a copy of
%    matlab in front of me)    
dklast = NaN; dk = NaN; % initialise to begin loop

iter = 0;
maxiter = 100;

while (abs(dk/dklast - beta - 1) > 1e-4) && (iter < maxiter) % loop until tolerance met
    iter = iter + 1; % keep count of iteration

    X1 = (1+beta)*X0 -beta*X0*A*X0; % calculate new iterate

    dklast = dk; % move old difference "new estimate to previous iterate"
    dk = norm(X1-X0,'fro'); % determine new difference 

    X0 = X1; % copy current iterate to "old" iterate for next iteration

end

I am wondering why you are using this convergence test at all. I would recommend using 

dk = norm(X1*X0-I,'fro');

which measures how close X1 is to the left inverse of $A$. Your termination criteria would then be

while dk > (some_tolerance) && iter < maxiter
    ....
end

As you currently have, you are measuring how much X1 changes from X0, which may be small, but still not an approximate inverse (or pseudoinverse) for $A$.

Related Question