MATLAB: How parallelize the solution of sparse matrices using mldivide

MATLABmldivideParallel Computing Toolboxparallelizationsparse matrix

I am trying to parallelize the solution of
x= A\B (mldivide.
My variables are: x = V(I*J*K), A = A(I*J*K x I*J*K) sparse matrix, vec = u (I,J,K) +V/constant +Bswitch (I*J*K x I*J*K)*V
To do this without parallelization, my code currently does this:
V_stacked = reshape(V,I*J*L,1);
vec = u_stacked + V_stacked/Delta + Bswitch*V_stacked;
V_stacked = A\vec;
To parallelize I have tried
u_stacked = reshape(u,I*J,L);
    V_stacked = reshape(V,I*J,L);
    BswitchTimesVstacked = Bswitch*reshape(V,I*J*L,1);
    BswitchTimesVstacked = reshape(BswitchTimesVstacked,I*J,L);
    vec = u_stacked + V_stacked/Delta + BswitchTimesVstacked;
    tic
    parfor l = 1:L 
        V_stacked(:,l) = A(:,:,l)\vec(:,l);
    end
But as A is still I*J*L times I*J*L, it wont work. I am not sure if 1. what I am doing so far is correct and 2. how to reshape B appropriately.
Mathematical info is here: http://www.princeton.edu/~moll/HACTproject/two_asset_kinked.pdf (section 4)
Any help is highly appreciated 🙂

Best Answer

You can make a 3D stack of sparse matrices A(:,:,l) by converting them to ndSparse type. Then, the parfor construct will work. Alternatively, you can make a block diagonal matrix where all the A(:,:,l) form the diagonal blocks. Then you can solve all systems simultaneously, which would take advantage of Matlab's internal parallelization. I don't know which would be faster.
Related Question