Mercurial > matrix-functions
comparison mftoolbox/rootpm_newton.m @ 0:8f23314345f4 draft
Create local repository for matrix toolboxes. Step #0 done.
| author | Antonio Pino Robles <data.script93@gmail.com> |
|---|---|
| date | Wed, 06 May 2015 14:56:53 +0200 |
| parents | |
| children |
comparison
equal
deleted
inserted
replaced
| -1:000000000000 | 0:8f23314345f4 |
|---|---|
| 1 function [X,k] = rootpm_newton(A,p,c) | |
| 2 %ROOTPM_NEWTON Coupled Newton iteration for matrix pth root. | |
| 3 % [X,k] = ROOTPM_NEWTON(A,P,C) computes the principal | |
| 4 % Pth root of the matrix A. | |
| 5 % C (default 1) is a convergence parameter. | |
| 6 % k is the number of iterations. | |
| 7 | |
| 8 if nargin < 3, c = 1; end | |
| 9 | |
| 10 n = length(A); | |
| 11 M = A/c^p; | |
| 12 X = c*eye(n); | |
| 13 tol = mft_tolerance(A); | |
| 14 maxit = 20; | |
| 15 | |
| 16 relres = inf; | |
| 17 | |
| 18 for k=1:maxit | |
| 19 | |
| 20 X = ( ((p+1)*eye(n) - M)/p )\X; | |
| 21 M = power_binary( ((p+1)*eye(n) - M)/p, p) * M; | |
| 22 | |
| 23 relres_old = relres; | |
| 24 relres = norm(M-eye(n),inf); | |
| 25 | |
| 26 if relres <= tol || relres > relres_old/2, return, end | |
| 27 | |
| 28 end | |
| 29 error('Not converged after %2.0f iterations', maxit) |