Use of a time integrator to find a matrix inverse

35 Views Asked by At

Under what conditions can we use a time integrator (ODE solver) to find matrix inverse? Specifically, say we have a linear system $Ax=b$. Under some conditions on $A$, we can rewrite this equation as $$ \dot{x} + Ax = b $$ and then apply a numerical method to solve for $x$. For example, we can use a forward Euler scheme to generate the iteration $$ x_{k+1} = x_k + \Delta t(-Ax_k + b) $$ In terms of code, the following MATLAB/Octave script demonstrates that this works

% Set the size                                                                  
m = 5;                                                                          

% Set the number of iterations
niter = 100000;

% Set the timestep
dt = 0.001;

% Find the operator and rhs                                                     
A = randn(m);                                                                   
b = randn(m,1);   
A = triu(A) + triu(A,1)';
[v d]=eig(A);
d = rand(m,1)*100;
A = v*diag(d)*v';

% Find the initial guess                                                        
x = zeros(m,1);                                                                 

% Run the iteration                                                             
for i=1:niter                                                                   
    x = x + dt*(-A*x+b);                                                        
end                                                                             

% Find the residual                                                             
xtrue = A\b;                                                                    
fprintf('|| x - xtrue ||: %e\n',norm(x-xtrue)); 

> test01
|| x - xtrue ||: 1.551632e-16

Now, of course it doesn't always work. From what I can tell, it works for symmetric matrices where the eigenvalues are positive. That said, I don't have a proof nor better conditions when it should or not occur. As such, under what conditions for $A$ is this valid? Further, why does this work under some situations?