libmir/numir

Create diagonal matrix from slice

Closed this issue · 6 comments

I was trying to create a diagonal matrix and found it more painful than I had been expecting. In the code below, I had tried the commented line first, but that did not work. I then tried to use a foreach (which took some time to realized I needed to enumerate the temp variable), but that also resulted in a number of errors.

I'm sure I can figure out some other way to do this. However, I would think that the diag function would do this.

import mir.ndslice.algorithm : sliced;
import numir.core;
import std.range : enumerate;

void main()
{
	//auto x = linspace(2, 3, 5).diag;
	auto temp = linspace(2, 3, 5);
	auto x = zeros(5);
	foreach(i, t; temp.enumerate)
	{
		x[i, i] = t;
	}
}
9il commented
import mir.ndslice: diagonal;
import numir.core;

auto matrix = zeros(5);
matrix.diagonal[] = linspace(2, 3, 5);

I suppose my point is more that it needs to be modified if they want the same sort of approach as numpy. In particular, I was thinking something like below. I tried to get it to work, adding in some isMatrix/isVector code, but linking was taking forever.

I recommend adding isMatrix/isVector to mir-algorithm.


auto diag(S)(S s, long k=0) pure
     if (isMatrix!S)
{
     import mir.ndslice : diagonal;

     auto sk = k >= 0 ?  s[0 .. $, k .. $] : s[-k .. $, 0 .. $];
     return sk.diagonal;
}

auto diag(S)(S s) pure
     if (isVector!S)
{
     import mir.ndslice : diagonal, zeros;

     auto result = zeros(s.length, s.length);
     result.diagonal[] = s;
     return result;
}

EDIT: slight correction to zeros line.

@9il Just as an FYI, the following code does not compile because diagonal cannot be an lvalue

auto temp = [10, 10, 10].sliced;
auto result = iota(3, 3);
result.diagonal[] = temp;

EDIT: This works though:

	auto temp = [10, 10, 10].sliced;
	auto result = iota(3, 3);
	size_t i = 0;
	foreach (ref e; result.diagonal)
	{
		e = temp[i];
		i++;
	}

EDIT: Arg, I spoke too soon. It compiles, but the actual value doesn't change. The example in the diagonal documentation in mir-algorithm uses mir.ndslice.allocation.slice to allocate it. The same approach didn't work with iota.

@9il thank you for answering!

@jmh530 So is your point is diag to do this?

>>> numpy.diag(numpy.linspace(2, 3, 5))
array([[ 2.  ,  0.  ,  0.  ,  0.  ,  0.  ],
       [ 0.  ,  2.25,  0.  ,  0.  ,  0.  ],
       [ 0.  ,  0.  ,  2.5 ,  0.  ,  0.  ],
       [ 0.  ,  0.  ,  0.  ,  2.75,  0.  ],
       [ 0.  ,  0.  ,  0.  ,  0.  ,  3.  ]])

I think your proposal is useful too. Why don't you PR?

BTW in your first example

auto result = iota(3, 3);

This should be auto result = iota(3, 3).slice; as you mentioned mir.ndslice.allocation.slice works right (I confirmed that result == [[10, 1, 2], [3, 10, 5], [6, 7, 10]]). Is it still not OK?

@ShigekiKarita I have submitted a pull request to address this issue.