2

I am trying to solve a Matrix in Math.Net when one of the actual solutions to the matrix is 0, but I am getting -NaN- as results.

Here is an example matrix which has already been reduced for simplicity.

1 0  1 | 10000 
0 1 -1 | 1000
0 0  0 | 0

Code example:

public void DoExample()
{
    Matrix<double> A = Matrix<double>.Build.DenseOfArray(new double[,] {
        { 1, 0, 1 }, 
        { 0, 1, -1 }, 
        { 0, 0, 0 }, 
    });

    Vector<double> B = Vector<double>.Build.Dense(new double[] { 10000, 1000, 0 });

    var result = A.Solve(B);
}

The solution I am hoping to get to is [ 10000, 1000, 0 ].

As you can see, the result I want is already the augment vector. This is because I simplified the matrix to reduced row echelon form (RREF) by hand using Gauss-Jordan for this example. If I could somehow use a Gauss-Jordan operations within Math.Net to do this, I could check for the scenario where an all 0 row exists in the RREF matrix. Can this be done?

Otherwise, is there any way I can recognize when 0 is the only possible solution for one of the variables using the existing Math.Net linear algebra solver operations?

Thanks!

foreach_potato
  • 23
  • 5
  • Try an SVD solver. It'll do the best it can with your degenerate matrix. As you can see, you have three unknowns but only two equations. The best you can do is solve for one of the three in terms of the other two. – duffymo May 11 '17 at 09:07
  • duffymo, Thank you. I was so hung up on the problem being that one of my expected values was 0 that I missed this. Perhaps I can start substituting in previously known "good" values with my application in order to cheat towards a reasonable solution in these situations. Is there a way to detect which variable is free (or lacking sufficient clarity)? – foreach_potato May 11 '17 at 15:14
  • Any one of the three variables can be considered "free". In this case, the third one is natural b/c its row is all zeros. – duffymo May 11 '17 at 15:17
  • Oh, and the SVD solver also fails when using the solve method. I don't know if using the iterative solver methods will help. (I am utterly lost on those as the online documentation is intimidating and I can find no good examples.) – foreach_potato May 11 '17 at 15:19
  • Search my name for SVD. Went through it in detail a while back. http://stackoverflow.com/questions/19763698/solving-non-square-linear-system-with-r/19767525#19767525 – duffymo May 11 '17 at 15:20
  • Thanks duffymo, I'm appreciative to how in-depth that explanation was. Unfortuantely, I am looking for an answer specific to the Math.Net libraries (so I didn't have to implement your solution myself) which doesn't seem to exist just yet. – foreach_potato May 11 '17 at 16:16
  • As it turns out, I can "cheat" by substituting one of my expected values back into the equation to get the other two. It doesn't matter which of the values it is either, so long as my substitution brings me to where my rank = my row count. I may implement some type of consistency check and and raise a warning if any of the suspected substitutions lead to a different result. Obviously, I'd rather have the math answer it for me, but I'm asking for too much here. – foreach_potato May 11 '17 at 16:17

2 Answers2

2

This is degenerate matrix with rank 2, and you cannot expect to get true solution (there are infinity number of solutions)

MBo
  • 77,366
  • 5
  • 53
  • 86
  • @foreach_potato, Just to give an example, `[9999, 1001, 1]` is as good solution to this system of equations as the one you were hoping for – SergGr May 11 '17 at 14:22
  • Thanks MBo and SergGr for illustrating this for me. I lacked the imagination to come up with any other solutions this morning. – foreach_potato May 11 '17 at 15:15
1

The iterative solver can actually handle this, for example

using MathNet.Numerics.LinearAlgebra.Double.Solvers;
A.SolveIterative(B, new MlkBiCgStab());

returns

[10000, 1000, 0]

Interestingly, with the MKL Native Provider this also works with the normal Solve routine, but not with the managed provider (as you have found out) nor with e.g. the OpenBLAS native provider.

Christoph Rüegg
  • 4,626
  • 1
  • 20
  • 34
  • Wow, that is even more amazing than I could have hoped for! Rather than substituting in suspected variables (as I mentioned above), I'll try this and simply check the result. Thanks, Christoph! – foreach_potato May 11 '17 at 17:18
  • That solver tends to get the solution, only not always the one I am expecting. For example, running the same matrix through with row 1 moved to the front gives me a different solution of [0, 11000, 10000] So same matrix, yet different result based on whatever the algorithm is doing being the scenes. I cant really rely on the result without providing it some input. Can I use the preconditioner to provide acceptable values to try to adhere to? I'm sorry, but I'm finding no examples of how the preconditioner works. – foreach_potato May 11 '17 at 20:45
  • I'm not sure I follow what you mean by moving row 1 to the front - if I switch rows 0 and 1 of both A and B I get exactly the same result. Is what you get in your case still a correct solution? – Christoph Rüegg May 11 '17 at 21:59
  • Sorry, I mistyped. If you move row 3 to the top of the matrix (so that it is now row 1), you'll get the other solution I mentioned. In other words: A: { 0, 0, 0 }, { 1, 0, 1 }, { 0, 1, -1 } B: {0, 10000, 1000} Result: {0, 11000, 10000} – foreach_potato May 11 '17 at 22:31