Nov 15, 2013 at 12:57 PM
Edited Nov 15, 2013 at 2:49 PM

Hi to all,
I am having a problem with the solvers provided in MathNet. Let me explain what I am trying to do, I have a very big nonsquare sparse matrix (five nonzeroes per row), of a size of 20,000x12,000, and I would like to solve it with the least squares method,
but I am stucked.
I could prepare the matrix to solve it with a direct solver like LUby setting some fixed values and getting a square matrix and some other things (it is the way I was trying to implement). But, how can I solve this system by least squares method? Can I do it
without any preprocessing directly over the Mx = b system?
I have this variables:
MathNet.Numerics.LinearAlgebra.Double.SparseMatrix M = new MathNet.Numerics.LinearAlgebra.Double.SparseMatrix(m,n, 0.0);
//I populate the matrix with the .SetRow method
MathNet.Numerics.LinearAlgebra.Double.DenseVector b = new MathNet.Numerics.LinearAlgebra.Double.DenseVector(m); //b is all zeroes.
But when I try to solve it (directly, without preparing it, just to see what happens):
MathNet.Numerics.LinearAlgebra.Double.DenseVector result = M.LU().Solve(b);
It appears a message saying that SparseMatrix does not contain a definition called LU or an extension method called like that...
I have tried to define the matrixv M as a DenseMatrix, but I get an OutOfMemory exception, but defining a smaller
DenseMatrix I get the same message...
What is happening? I do not know to use the solvers, any help?
Thanks in advanced!


Coordinator
Nov 15, 2013 at 10:49 PM
Edited Nov 15, 2013 at 11:11 PM

Hi,
While it is technically possible to go for a direct solver/factorization, for a matrix of this size (~1.8 GB) this is currently only feasible in practice (if at all) when using our native MKL provider and an actual dense matrix (not sparse, since it will perform
badly  this is something we hope to adress for the v3 release though).
I assume the "message" you see is actually a compilation error? In V2, LU() was indeed defined in an extension method so you may be missing its namespace. However, in V3 this is a proper method.
However, you may want to give the iterative sparse solvers a trial, especially with the new MILU0 preconditioner in v3.0.0alpha5:
NuGet packages: MathNet.Numerics v3.0.0alpha5, MathNet.Numerics.Data.Text v3.0.0alpha5
Namespaces: MathNet.Numerics.Data.Text, MathNet.Numerics.LinearAlgebra, MathNet.Numerics.LinearAlgebra.Double, MathNet.Numerics.LinearAlgebra.Double.Solvers, MathNet.Numerics.LinearAlgebra.Solvers
I'm using a matrix from
NIST MatrixMarket here, which is 1633x1633, 1.75% nonzero, and compute a b vector (that's why I'm also using the .Data.Text package here, you won't need it in your case). But it works also well e.g. with
UTM5940 which is a bit larger, 5940x5940, 0.24% nonzero. Both should converge in less than 1020ms.
// m is Matrix<double> but its value is actually of type Double.SparseMatrix
var m = MatrixMarketReader.ReadMatrix<double>(@"fidap007.mtx");
var b = m * DenseVector.Create(m.RowCount, 1.0);
Then:
var solver = new BiCgStab();
var preconditioner = new MILU0Preconditioner();
var iterator = new Iterator<double>(
new ResidualStopCriterium(1.0e8),
new IterationCountStopCriterium<double>(1000));
var x = m.SolveIterative(b, solver, iterator, preconditioner);
var status = iterator.Status;
However, how well this works depends on your data. You may need to tweak the parameters a bit, the stop criteria, preconditioner, or try a solver other than BiCgStab.
Thanks,
Christoph

