Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add procedures to create and handle diagonal matrices #597

Merged
merged 4 commits into from
Oct 18, 2023

Conversation

AngelEzquerra
Copy link
Contributor

Add procedures to create diagonal matrices as well as to set and get diagonals from matrices.

Copy link
Collaborator

@Vindaar Vindaar left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Generally great work! A few minor things as noted in the comments.

And another one: Shouldn't we implement diag using set_diagonal as well? Then diag can also support anti and we remove some duplicated code.

@@ -47,3 +47,137 @@ proc vandermonde*[T](x: Tensor[T], order: int): Tensor[float] =
let orders = arange(order.float)
for i, ax in enumerateAxis(result, axis = 1):
result[_, i] = (x ^. orders[i]).unsqueeze(axis = 1)

proc diag*[T](d: Tensor[T], k = 0): Tensor[T] {.noInit.}=
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
proc diag*[T](d: Tensor[T], k = 0): Tensor[T] {.noInit.}=
proc diag*[T](d: Tensor[T], k = 0): Tensor[T] {.noInit.} =

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You are right, implementing the diag procedure based on set_diagonal makes more sense as it reduces code duplication and let's us easily add support for the anti flag as well.
I will fix the missing spaces as well.

for i in 0 ..< d.size:
result[i-k, i] = d[i]

proc diagonal*[T](a: Tensor[T], k = 0, anti = false): Tensor[T] {.noInit.}=
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
proc diagonal*[T](a: Tensor[T], k = 0, anti = false): Tensor[T] {.noInit.}=
proc diagonal*[T](a: Tensor[T], k = 0, anti = false): Tensor[T] {.noInit.} =

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I will fix this.

Comment on lines 125 to 127
when compileOption("boundChecks"):
let size = min(a.shape[0] - abs(k), a.shape[1])
assert size == d.size, &"Input size ({d.size}) does not match the {k}-th upper anti-diagonal size ({size})"
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'm not sure if this is what we would "intend". For example if the user wishes bound checks on a danger build the assertion is still going to be ignored. If we explicitly activate the assertion only for bounds checks being active, then maybe we should use doAssert to also have it work on danger builds or use a regular exception? Same for the other when below of course.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I will replace the calls to assert with doAssert as you suggest.

src/arraymancer/linear_algebra/special_matrices.nim Outdated Show resolved Hide resolved
src/arraymancer/linear_algebra/special_matrices.nim Outdated Show resolved Hide resolved
src/arraymancer/linear_algebra/special_matrices.nim Outdated Show resolved Hide resolved
Comment on lines 174 to 179
## Return a 2-D tensor with ones on the diagonal and zeros elsewhere
##
## Input:
## - Rank-1 tensor containg the elements of the diagonal
## - The index of the diagonal that will be set. The default is 0.
## Use k>0 for diagonals above the main diagonal, and k<0 for diagonals below the main diagonal.
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The doc string seems wrong. And do I understand correctly that the difference between eye and identity is just that it allows to construct rectangular matrices? What happens with shape.len > 2?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You are right, the comment is wrong. I will fix it and also add an assertion when the input is not a matrix.

…the procedure is not a matrix.

Also add a comment explaining the difference between eye() and identity().
This reduces code duplication and adds support for the "anti" argument to diag().
Add some missing assertion messages and replace some assertions with doAssert to trigger them when using -d:danger but --boundChecks:on.
@AngelEzquerra
Copy link
Contributor Author

The new changes should hopefully address all the existing comments.

Copy link
Collaborator

@Vindaar Vindaar left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Looks great now!

@Vindaar Vindaar merged commit c78981e into mratsim:master Oct 18, 2023
6 checks passed
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

3 participants