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

Support StaticArrays in X(t)_(inv)A_X(t) with ScalMat #181

Merged
merged 2 commits into from
Oct 2, 2023

Conversation

devmotion
Copy link
Member

Currently, ScalMat does not support X_A_Xt etc. with static arrays. For instance, running the tests in this PR on the master branch yields

StaticArrays: Error During Test at REPL[10]:42
  Test threw exception
  Expression: X_A_Xt(A, X) isa SMatrix{10, 10, Float64}
  setindex!(::SMatrix{10, 10, Float64, 100}, value, ::Int) is not defined.
   Hint: Use `MArray` or `SizedArray` to create a mutable static array
...
StaticArrays: Error During Test at REPL[10]:43
  Test threw exception
  Expression: X_A_Xt(A, X)  Matrix(X) * Matrix(A) * (Matrix(X))'
  setindex!(::SMatrix{10, 10, Float64, 100}, value, ::Int) is not defined.
   Hint: Use `MArray` or `SizedArray` to create a mutable static array
...
StaticArrays: Error During Test at REPL[10]:45
  Test threw exception
  Expression: X_invA_Xt(A, X) isa SMatrix{10, 10, Float64}
  setindex!(::SMatrix{10, 10, Float64, 100}, value, ::Int) is not defined.
   Hint: Use `MArray` or `SizedArray` to create a mutable static array
...
StaticArrays: Error During Test at REPL[10]:46
  Test threw exception
  Expression: X_invA_Xt(A, X)  Matrix(X) * (Matrix(A) \ (Matrix(X))')
  setindex!(::SMatrix{10, 10, Float64, 100}, value, ::Int) is not defined.
   Hint: Use `MArray` or `SizedArray` to create a mutable static array
...
StaticArrays: Error During Test at REPL[10]:48
  Test threw exception
  Expression: Xt_A_X(A, Y) isa SMatrix{10, 10, Float64}
  setindex!(::SMatrix{10, 10, Float64, 100}, value, ::Int) is not defined.
   Hint: Use `MArray` or `SizedArray` to create a mutable static array
...
StaticArrays: Error During Test at REPL[10]:49
  Test threw exception
  Expression: Xt_A_X(A, Y)  (Matrix(Y))' * Matrix(A) * Matrix(Y)
  setindex!(::SMatrix{10, 10, Float64, 100}, value, ::Int) is not defined.
   Hint: Use `MArray` or `SizedArray` to create a mutable static array
...
StaticArrays: Error During Test at REPL[10]:51
  Test threw exception
  Expression: Xt_invA_X(A, Y) isa SMatrix{10, 10, Float64}
  setindex!(::SMatrix{10, 10, Float64, 100}, value, ::Int) is not defined.
   Hint: Use `MArray` or `SizedArray` to create a mutable static array
...
StaticArrays: Error During Test at REPL[10]:52
  Test threw exception
  Expression: Xt_invA_X(A, Y)  (Matrix(Y))' * (Matrix(A) \ Matrix(Y))
  setindex!(::SMatrix{10, 10, Float64, 100}, value, ::Int) is not defined.
   Hint: Use `MArray` or `SizedArray` to create a mutable static array
...
Test Summary: | Pass  Error  Total  Time
StaticArrays  |   54      8     62  4.6s

The PR fixes this issue by changing the implementation to a non-mutating alternative, as done already e.g. for PDiagMat.

The disadvantage (also with the existing code for other psd matrix types) is that this leads to additional allocations when working with e.g. Array. I wonder if generally the non-mutating implementation should be the default (since it works with any matrix type) and additionally possibly internally mutating code for Array should be added (since it avoids allocations and this particular matrix type is so common).

@codecov-commenter
Copy link

codecov-commenter commented Sep 26, 2023

Codecov Report

All modified lines are covered by tests ✅

Files Coverage Δ
src/scalmat.jl 96.82% <100.00%> (ø)

... and 1 file with indirect coverage changes

📢 Thoughts on this report? Let us know!.

@devmotion
Copy link
Member Author

devmotion commented Sep 27, 2023

The PR is ready for review, @andreasnoack.

I wonder if generally the non-mutating implementation should be the default (since it works with any matrix type) and additionally possibly internally mutating code for Array should be added (since it avoids allocations and this particular matrix type is so common).

That's what I ended up doing for now. I only added specializations for Matrix - originally I thought StridedMatrix might be an alternative but AFAICT it is not guaranteed that a StridedMatrix is mutable.

@devmotion devmotion merged commit 77b5d59 into master Oct 2, 2023
11 checks passed
@devmotion devmotion deleted the dw/scalmat_staticarrays branch October 2, 2023 20:36
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