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

Faster construction of sparse block matrices #1326

Merged
merged 19 commits into from
Feb 13, 2025

Conversation

keileg
Copy link
Contributor

@keileg keileg commented Feb 10, 2025

Proposed changes

This PR introduces helper functions for fast construction of sparse block-diagonal matrices, of csr, csc and dia form. This allows us to avoid calling scipy.sparse.block_diag, which, though great, also takes account for all matter of edge cases, and is therefore slow. The new functions have been taken into use throughout the code.

Profiling (by line_profiler) indicates computational savings of ~90% on the construction of block diagonal discretization matrices under Ad parsing (which again consumed ~1/3 of all parsing time) on a 3d poromechanical case with a relatively coarse grid and not too many fractures. I do not know how this translates to other cases, but expect the speedup to be non-negligible in many cases.

I also renamed matrix_operations.csr_matrix_from_blocks -> matrix_operations.csr_matrix_from_dense_blocks, to enable more natural names for the new function csr_matrix_from_sparse_blocks (and similar for csc).

Tests are introduced, tests for csr_matrix_from_dense_blocks has been consolidated somewhat.

Resolves #1084.

Types of changes

What types of changes does this PR introduce to PorePy?
Put an x in the boxes that apply.

  • Minor change (e.g., dependency bumps, broken links).
  • Bugfix (non-breaking change which fixes an issue).
  • New feature (non-breaking change which adds functionality).
  • Breaking change (fix or feature that would cause existing functionality to not work as expected).
  • Testing (contribution related to testing of existing or new functionality).
  • Documentation (contribution related to adding, improving, or fixing documentation).
  • Maintenance (e.g., improve logic and performance, remove obsolete code).
  • Other:

Checklist

Put an x in the boxes that apply or explain briefly why the box is not relevant.

  • The documentation is up-to-date.
  • Static typing is included in the update.
  • This PR does not duplicate existing functionality.
  • The update is covered by the test suite (including tests added in the PR).
  • If new skipped tests have been introduced in this PR, pytest was run with the --run-skipped flag.

This is in preparation for the introduciton of a corresponding function for creating sparse matrices from sparse blocks
This is equivalent to, but faster than, sps.block_diag.

Also added a test
The method is faster than creating a csr-matrix or other alternatives

Also test
This should give notable computational savings
Mocked discretization matrices should be sparse to comply with the real framework to be tested
Copy link
Contributor

@jwboth jwboth left a comment

Choose a reason for hiding this comment

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

I did not review the entire code (yet). Had only a look so far at the main methods, which seemed good to me so far.

src/porepy/numerics/ad/_ad_utils.py Show resolved Hide resolved
Copy link
Contributor

@IvarStefansson IvarStefansson left a comment

Choose a reason for hiding this comment

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

Looks good in the main. Some questions and suggestions.

src/porepy/models/geometry.py Show resolved Hide resolved
src/porepy/numerics/ad/grid_operators.py Outdated Show resolved Hide resolved
src/porepy/numerics/fv/tpfa.py Show resolved Hide resolved
src/porepy/numerics/linalg/matrix_operations.py Outdated Show resolved Hide resolved
src/porepy/numerics/linalg/matrix_operations.py Outdated Show resolved Hide resolved
tests/numerics/ad/test_operators.py Show resolved Hide resolved
tests/utils/test_matrix_operations.py Show resolved Hide resolved
tests/utils/test_matrix_operations.py Show resolved Hide resolved
@Yuriyzabegaev
Copy link
Contributor

Great change! As far as I understood, this will fix #1084.

@keileg
Copy link
Contributor Author

keileg commented Feb 12, 2025

@IvarStefansson Thanks for the review. Most things should be taken care of now, though I left some comments for you to resolve.

From your comments, and my own reading of the code, it is clear that it is easy to misinterpret the different approaches to constructing sparse block matrices - it is actually necessary to read the documentation to understand what each function does. I tried to improve the documentation, without overemphasizing this point. Please have a look at this in particular.

EDIT: I upgraded the documentation of the functions with better descriptions, and inserted See also-references to the other methods. This should make users more aware that caution is needed when choosing a method.

PS: I'm also fighting mypy. Don't ask.

@keileg
Copy link
Contributor Author

keileg commented Feb 13, 2025

Great change! As far as I understood, this will fix #1084.

Yes, it should at least take away the main part of the cost.

@keileg keileg merged commit 628c26c into pmgbergen:develop Feb 13, 2025
6 checks passed
@keileg keileg deleted the sparse_diags branch February 13, 2025 16:05
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.

Performance bottleneck in MergedOperator.parse
4 participants