Quansight-Labs / numpy.net

A port of NumPy to .Net
BSD 3-Clause "New" or "Revised" License
131 stars 14 forks source link

np.linalg.eig Alternative #10

Closed AsiaMartini closed 2 years ago

AsiaMartini commented 3 years ago

Hi need to translate this peace of code:

 public static (NDarray, NDarray) _get_r_quat(NDarray d_matrix) {
    NDarray eigvals, eigvects;
    (eigvals, eigvects) = np.linalg.eig(d_matrix);
    NDarray beta_1 = np.argmax(eigvals);
    var r_quat = eigvects[":", beta_1];
    return (beta_1, r_quat);
}

I know the np.linalg.eig function is not available in this library, but I was wondering if you might help me to write it down in another way. I don't understand the maths behind it enough to break it down into little pieces.

Thanks!

KevinBaselinesw commented 3 years ago

The linear algebra module is a huge task and I don't have time to implement that right now. Also, I have no background knowledge of linear algebra so that makes it extra challenging.

I have done a little bit of research and I think this library might be able to offer some help. I was planning on experimenting with this library to see if I could reproduce the numpy linalg namespace outputs.

https://www.mathdotnet.com/ https://discuss.mathdotnet.com/t/evd-eigenvalues-ascending-eigenvectors-descending/52

I wonder if there is a linear algebra group somewhere on the internet that you could ask this question to. Someone with linear algebra expertise might be able to point us (you and I) towards a simple answer.

AsiaMartini commented 3 years ago

I'm giving a try with Math.NET. But comparing to Numpy.NET values, they seems different, and ordered in a different way (eigvals). This is causing a crash for me, cause later I look for the index of the maximum value and it's totally different (3 vs 0).

Here the results, same inputs (Math.NET to the left, Numpy.NET to the right): image

Here the code with Math.NET (and the error): image

Here the original Numpy.NET code:

NDarray eigvals, eigvects;
(eigvals, eigvects) = np.linalg.eig(d_matrix);
Console.WriteLine("-----------------");
Console.WriteLine($"\neigvals :\n{eigvals}");
Console.WriteLine($"\neigvects :\n{eigvects}");
NDarray beta_1 = np.argmax(eigvals);
Console.WriteLine($"\nbeta_1 :\n{beta_1}");
var r_quat = eigvects[":", beta_1];
return (beta_1, r_quat);
KevinBaselinesw commented 3 years ago

I think this is encouraging. The values appear to be almost exactly the same although in a different order.
Numpy.net is printing the values in exponential notation but seems to be the same value as Math.net.

I have no idea what the ordering issue is or how to fix it.

KevinBaselinesw commented 3 years ago

I wonder if the input array that you are passing to the Math.Net function is in the same order as Numpy. How did you get the array to pass to that function?

AsiaMartini commented 3 years ago

The input array is the c_matrix you can see in the logs screenshot above (as you can see, is absolutely equal in both procedures). The full code is in the repository:

https://github.com/AsiaMartini/simil-dotnet

The main function to call is always test(). The function in which you can find np.linalg.eig is: (ndarray, ndarray) _get_r_quat2(ndarray d_matrix).

The Simil.cs class contains the working code with Numpy.NET. If you want to give a try, I added a Simil-NumpyDotNet.cs class with the code I'm working on (using your library and Math.NET).

AsiaMartini commented 3 years ago

About the ordering, I can read in Math.NET documentation that EigenValues are sorted asc. The problem is I cannot find a way to avoid it, and also I cannot understand which kind of sort numpy is using...

AsiaMartini commented 3 years ago

The solution might be way more easier than we think!

This line returns an ndarray (in my case, containing a value of "3"): ndarray beta_1 = np.argmax(eigvals); {INTP 3}

I tried to replace this line: var r_quat = eigvects[":", beta_1];

With this: var r_quat = eigvects[":", 3];

This is the result: {DOUBLE { 0.4999999999999995, -0.5000000000000002, 0.4999999999999995, -0.5000000000000002 }}

This is the expected one: [ 0.5 -0.5 0.5 -0.5]

Great!

I only need your help to understand how can I get the "3" value out of the ndarray beta_1. Can you please help me?

KevinBaselinesw commented 3 years ago

np.argmax() returns the position in the array of the max value. Since the eigvals is different order, the position of the max value is going to be different.

AsiaMartini commented 3 years ago

I wrote that the order doesn't matter if the eigvects are ordered the same way. They only need to be synchronized.

In other words, both eigvects and eigvals are in another order, so the index reflects the same value.

AsiaMartini commented 3 years ago

Anyway, how can I extract the double value inside a ndarray? I feel dumb, I can't find a proper function.

KevinBaselinesw commented 3 years ago

ndarray a;

if single dimension array, you can use (int)a.GetItem(index) or (int)a[index];

if 2 dimensional array, you can use (int)a[x,y]; if 3 dimensional array, you can use (int)a[x,y,z];

AsiaMartini commented 3 years ago

Great, thanks! I already tried to use GetItem() but I thought it was 0-indexed. With GetItem(1) everything works fine.

Just one last (I promise :D) question... I finally get everything working but I would like to work with float insted of double type... But even if I initialize the np.array with dtype Float64, I get objects like this:

{DOUBLE 
{ { 0.0, 0.0, 0.0, 0.0 },
  { 0.0, 2.0, 2.0, 0.0 },
  { 2.0, 3.0, 1.0, 0.0 } }
}

And then errors like this: image

Am I doing something wrong?

KevinBaselinesw commented 3 years ago

When you create your array, specify Float32 instead. Float64 maps to System.Double. Float32 maps to System.Single (AKA float). These terms are from Numpy.

     var a = np.arange(2.5, 37.7, 2.2, dtype: np.Float32);

the dtype: variable is available on probably every ndarray creation routine.

If you have a double and you want to convert it to float, try something like this ndarray.AsType(np.Float32);

KevinBaselinesw commented 3 years ago

FYI, GetItem() IS 0 based. If you are getting the value you want with index 1, that suggests the returned array is not what you were expecting.