Iterative solver for matrix division

Question

I have two huge matrices $A$ and $B$. I am trying to find some iterative solvers like bcg or lsqr in Matlab.

I mean if my matrix $A$ is sparse and $n\times n$ and $B$ is a column vector of size $n\times 1$, I can use lsqr and bcg to solve the equation $Ax=B$.

Now suppose I need to solve $XD=C$, so I need to calculate $CD^{-1}$ where both $C$ and $D$ are huge matrices. If I use matlab's C/D operation directly, it consumes lots of memory and crashes. Are there any iterative solvers for this operation instead of using the forward slash operator directly?

Answer 1

The function GMRES offers the best speed, though I think QMR uses less memory. Otherwise the lu function allows you to recompose the matrix into an upper and a lower matrix like so: [L,U,P] = lu(A); x = U(L(P*b)); Where A*x = b

Answer 2

Let's frame your problem as a Sylvester equation:

$$ AX+XD - C = 0$$

with $A = 0$. We can solve Sylvester equations using lyap.

As per my comment, the following will work:

n = length(D);
lyap(zeros(n),D,-C);

Example:

n = 4;
B = rand(n);
C = rand(n);
C*inv(B)-lyap(zeros(n),B,-C)

ans =

  1.0e-014 *

   -0.2220   -0.2220    0.1221    0.3997
   -0.0666   -0.1582   -0.0278   -0.0611
    0.0666    0.0555    0.0444    0.5107
    0.0222   -0.0222    0.1332    0.7883

This has the benefit of the fact that the code that ships with lyap calls directly into LAPACK routines that should be capable of calling the appropriate sparse solvers/etc.

Answer 3

How about this: rewrite $XD = C$ as $D^T X^T = C^T$. Now solve for $X^T$ column by column.

For example, to get the first column of $X^T$, you would solve $D^T x = c_1$, where $c_1$ is the first column of $C^T$.

You will have to solve many systems of equations (one for each column of $X^T$), but that seems unavoidable, and at least each system is tractable.

If you use a method that first factors $D^T$, then (after the expensive factorization step, which is a one time cost) each system can be solved relatively cheaply using this factorization. (But I'm not sure if a factorization method is appropriate for your large problem.)