Skip to content

Use internal sparse solver to replace EIGEN and MKL #77

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

Merged
merged 114 commits into from
Jul 15, 2022

Conversation

TonyXiang8787
Copy link
Member

@TonyXiang8787 TonyXiang8787 commented May 31, 2022

Introduction

This PR implements an internal sparse solver which replaces EIGEN and MKL solver.

Closes #75

Design choices

Adjustment outside sparse solver

  • In bus reordering of topology search, the fill-ins are pre-analyzed and recorded. So all the off-diagonal fill-ins are known in advance.
  • In the construction of Node Admittance Matrix, two sparse matrix structures are made
    • The node admittance matrix
    • The LU decomposition matrix with the fill-ins. This will have more entries. There is also a mapping between the linear index of node admittance matrix and LU matrix.
  • The structure of node admittance matrix is now shared between symmetric and asymmetric math solver.

Implementation of LU sparse solver

  • It uses standard right-looking Gaussian elimination for LU decomposition.
  • The buses are already pre-ordered, so no permutation or fill-in analysis is needed. At the beginning of the factorization, the LU matrix is pre-allocated including all the fill-ins.
  • Since power system matrices are always diagonal major, numeric pivoting across the sparse structure is not needed.
  • The entry in the matrix is either a scalar (real/complex number) or a block (2x2, 6x6)
    • For a scalar matrix it just uses the normal sparse LU decomposition.
    • For block matrix it factorizes the matrix per dense block. Inside the diagonal block a dense LU decomposition is used with full pivoting. The factorized diagonal LU block will be used to further factorize the off-diagonal elements. Because it is a full pivoting within the block, there will be row/column permutation within the block. This applies to the right hand side and unknown vectors.
  • The matrix is always structurally symmetric. The code uses this feature to quickly find some non-zero entries in the matrix: if position (i, j) is non-zero, position (j, i) is also non-zero.
    -The solver is a class template with types of block matrices. The compiler can then optimize (e.g. loop unrolling) for difference block sizes.

Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
@TonyXiang8787 TonyXiang8787 marked this pull request as ready for review June 17, 2022 11:39
Signed-off-by: Peter Salemink <[email protected]>
Signed-off-by: Peter Salemink <[email protected]>
Copy link
Contributor

@bramstoeller bramstoeller left a comment

Choose a reason for hiding this comment

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

I'm not done with my review yet, but here's some comments for the. sparse_lu_solver

@petersalemink95
Copy link
Member

For features this large: next time please only implement the feature and do refactoring/renaming afterwards in a new PR. This way it will be easier to review large changes like this.

@petersalemink95
Copy link
Member

In general, please add more documentation to classes/functions

bramstoeller
bramstoeller previously approved these changes Jul 13, 2022
Copy link
Contributor

@bramstoeller bramstoeller left a comment

Choose a reason for hiding this comment

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

Good work, very impressive!

In general, please try to avoid refactoring code in such a large pull request, because that makes it much harder to review. And in some cases it would be more readable to use a little less template magic (even at the cost of some code duplication, imho). For example, the is_block constexpr in SparseLUSolver basically splits most functions in two; one part for scalars and one for blocks. I would prefere to have two separate functions.

Signed-off-by: Tony Xiang <[email protected]>
Signed-off-by: Tony Xiang <[email protected]>
@TonyXiang8787 TonyXiang8787 merged commit f8e810b into release/1.4 Jul 15, 2022
@TonyXiang8787 TonyXiang8787 deleted the feature/internal-sparse-solver branch July 15, 2022 12:56
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.

[FEATURE] Use in-house implementation of sparse solver
3 participants