Skip to content

Commit

Permalink
Merge pull request #971 from gridap/bugfix-block-assembly
Browse files Browse the repository at this point in the history
Bugfix: Incorrect block assembly when permuting variables
  • Loading branch information
JordiManyer authored Jan 11, 2024
2 parents 5c70c78 + 0a77dd0 commit f88e752
Show file tree
Hide file tree
Showing 2 changed files with 11 additions and 7 deletions.
1 change: 1 addition & 0 deletions NEWS.md
Original file line number Diff line number Diff line change
Expand Up @@ -18,6 +18,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
### Fixed

- Fixed issue where `FineToCoarseField` would not work for vector-valued fields. Since PR[#970](https://github.com/gridap/Gridap.jl/pull/970).
- Fixed issue where `BlockSparseMatrixAssembler` would not work correctly when permuting variables. Since PR[#971](https://github.com/gridap/Gridap.jl/pull/971).

## [0.17.21] - 2023-12-04

Expand Down
17 changes: 10 additions & 7 deletions src/MultiField/BlockSparseMatrixAssemblers.jl
Original file line number Diff line number Diff line change
Expand Up @@ -36,17 +36,20 @@ function FESpaces.get_vector_builder(a::BlockSparseMatrixAssembler{NB,NV}) where
return expand_blocks(a,builders)
end

expand_blocks(::BlockSparseMatrixAssembler{NB,NB},blocks::MatrixBlock) where NB = blocks
expand_blocks(::BlockSparseMatrixAssembler{NB,NB},blocks::VectorBlock) where NB = blocks

function expand_blocks(a::BlockSparseMatrixAssembler{NB,NV},blocks::MatrixBlock) where {NB,NV}
function expand_blocks(a::BlockSparseMatrixAssembler{NB,NV,SB,P},blocks::MatrixBlock) where {NB,NV,SB,P}
if (NB == NV) && all(x -> x[1] == x[2], enumerate(P))
return blocks
end
block_map = get_block_map(a)
ArrayBlockView(blocks,block_map)
return ArrayBlockView(blocks,block_map)
end

function expand_blocks(a::BlockSparseMatrixAssembler{NB,NV},blocks::VectorBlock) where {NB,NV}
function expand_blocks(a::BlockSparseMatrixAssembler{NB,NV,SB,P},blocks::VectorBlock) where {NB,NV,SB,P}
if (NB == NV) && all(x -> x[1] == x[2], enumerate(P))
return blocks
end
block_map = map(idx -> CartesianIndex(idx[1]), diag(get_block_map(a)))
ArrayBlockView(blocks,block_map)
return ArrayBlockView(blocks,block_map)
end

function get_block_ranges(NB::Integer,SB,P)
Expand Down

0 comments on commit f88e752

Please sign in to comment.