<div class="gmail_quote">On Sat, Dec 11, 2010 at 00:03, Luke Bloy <span dir="ltr">&lt;<a href="mailto:luke.bloy@gmail.com">luke.bloy@gmail.com</a>&gt;</span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;">
<div bgcolor="#ffffff" text="#000000">I&#39;m solving a diffusion problem. essentially I have 2,000,000
    possible states for my system to be in. The system evolves based on
    a markov matrix M, which describes the probability the system moves
    from one state to another. This matrix is extremely sparse on the
    &lt; 100,000,000 nonzero elements. The problem is to pump
    mass/energy into the system at certain states. What I&#39;m interested
    in is the steady state behavior of the system.<br>
    <br>
    basically the dynamics can be summarized as <br>
    <br>
    d_{t+1} = M d_{t} + d_i<br>
    <br>
    Where d_t is the state vector at time t and d_i shows the states I
    am pumping energy into. I want to find d_t as t goes to infinity.<br>
    <br>
    My current approach is to solve the following system.<br>
    <br>
    (I-M) d = d_i<br></div></blockquote><div><br></div><div>So you want to do this for some 500,000 d_i?  What problem are you really trying to solve?  Is it really to just brute-force compute states for all these inputs?  What are you doing with the resulting 500k states (all 8 terabytes of it)?  Are you, for example, looking for some d_i that would change the steady state d in a certain way?</div>
<div><br></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;"><div bgcolor="#ffffff" text="#000000"><div class="im"><blockquote type="cite"><div>2. If you can afford the memory, a direct solve probably
        makes sense.</div>
    </blockquote>
    <br></div>
    My understanding is the inverses would generally be dense. I
    certainly don&#39;t have any memory to hold a 2 million by 2 million
    dense matrix, I have about 40G to play with. So perhaps a
    decomposition might work? Which might you suggest?</div></blockquote></div><br><div>While inverses are almost always dense, sparse factorization is far from dense.  For PDE problems factored in an optimal ordering, the memory asymptotics are n*log n in 2D and n^{4/3} in 3D.  The time asymptotics are n^{3/2} and n^2 respectively.  Compare to n^2 memory, n^3 time for dense.</div>
<div><br></div><div>Jed</div>