I have a matrix which is complex symmetric. It is around 50,000 elements per side. It is a Method of Moments matrix. Is it feasible to use a standard direct solver such as Lapack to do a matrix decomposition such as $LU$ or $LDL^T$. I need to solve approximately 7200 right hand side vectors with this matrix, so if I were to take an iterative approach I would need to run it 7200 times.
Tell me more
×
Computational Science Stack Exchange is a question and answer site for
scientists using computers to solve scientific problems. It's 100% free, no registration required.
|
|
[Moving the comment thread down here] A direct solve will be pretty expensive, likely on the order of an hour for the sort of machine you describe. If you want invest time in something faster using more specialized methods, I suggest using FMM and a block Krylov method. FMM applies the operator in $\mathcal O(N)$ memory and time. FMM has a big constant, but less than $N=50000$, so it should be faster than dense matrix multiply. I don't know if you'll find a library implementation that applies to multiple vectors at once. Block Krylov methods apply to $k$ right hand sides, converging as much as $k$ times faster than with a single right hand side. The cost per iteration is $\mathcal O(N k)$ to apply the operator using FMM plus $\mathcal O(N k^2)$ to orthogonalize (with an additional factor for the restart length if using GMRES) and $\mathcal O(k^3)$ "scalar" work. You can choose $k=7200$ or some smaller number, to reduce the $Nk^2$ and $k^3$ cost at the expense of (slightly) slower convergence. Given its source in integral equations, I assume that the continuum operator associated with your problem is of the form $I + K$ where $K$ is a compact operator. Being compact, $K$ has essentially finite rank, therefore the discretized system has an essentially constant number of eigenvalues that are not "near" $1$, thus we expect the block Krylov method to converge very fast. |
|||
|
|
It will be slow but feasible. (You can find out the likely time it takes by solving first for sevaral $k$ reduced problems with dimension $50000/k$ and $7200/k$ right hand sides, and multiply the resulting time with $k^3$.) You need to do pivoting for accuracy (as $LDL^T$ may be unstable http://eprints.ma.man.ac.uk/344/), and you should process all right hand sides simultaneously (i.e., solve $AX=B$ with a matrix $B$ composed of all right hand sides) to take maximal advantage of the BLAS3 kernel. But it is still a brute force approach, and a multipole method is the more elegant choice. Moreover, you can combine the latter with the method in http://scicomp.stackexchange.com/a/2422/1117 to avoid having to solve 7200 independent problems. |
|||
|
|
