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.

Sign up for GitHub

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/SparseMatricesCSR #118

Merged
merged 7 commits into from
Oct 29, 2019
Merged

Add/SparseMatricesCSR #118

Show file tree
Hide file tree
Changes from 3 commits
Commits
Show all changes
7 commits
Select commit Hold shift + click to select a range
52f0a23
Add SparseMatricesCSR
victorsndvg Oct 25, 2019
f10a408
Finish SparseMatricesCSR integration into Assemblers.jl
victorsndvg Oct 28, 2019
0c48dd0
MultiSparseMatrixAssembler produces SparseMatrices also from SparseMa…
victorsndvg Oct 28, 2019
53d940e
Solve review comments in #118
victorsndvg Oct 28, 2019
ee467fe
Solve review comments in #118
victorsndvg Oct 28, 2019
6ef127f
Add finalize_coo! call in assemble funcion for Assemblers and MultiAs…
victorsndvg Oct 28, 2019
08edf7e
Add new signatures for LinearFEOperator and NonLinearFEOperator funct…
victorsndvg Oct 29, 2019
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
6 changes: 6 additions & 0 deletions Manifest.toml
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -245,6 +245,12 @@ uuid = "6462fe0b-24de-5631-8697-dd941f90decc"
deps = ["LinearAlgebra", "Random"]
uuid = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"

[[SparseMatricesCSR]]
deps = ["LinearAlgebra", "SparseArrays"]
git-tree-sha1 = "7896f9c00e5413cbaecd758ade9fe67bf50ad32e"
uuid = "a0a7dd2c-ebf4-11e9-1f05-cf50bc540ca1"
version = "0.3.1"

[[SpecialFunctions]]
deps = ["BinDeps", "BinaryProvider", "Libdl"]
git-tree-sha1 = "3bdd374b6fd78faf0119b8c5d538788dbf910c6e"
Expand Down
1 change: 1 addition & 0 deletions Project.toml
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -13,6 +13,7 @@ NLsolve = "2774e3e8-f4cf-5e23-947b-6d7e65073b56"
QuadGK = "1fd47b50-473d-5c70-9696-f719f8f3bcdc"
Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
SparseMatricesCSR = "a0a7dd2c-ebf4-11e9-1f05-cf50bc540ca1"
StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
TensorPolynomialBases = "f77b161d-e03a-5d68-b365-d5306aabeaa6"
TensorValues = "31c64edf-cdeb-50e4-845d-ed2334360c20"
Expand Down
62 changes: 41 additions & 21 deletions src/FESpaces/Assemblers.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5,11 +5,14 @@ using Gridap
using Gridap.Helpers

using SparseArrays
using SparseMatricesCSR

export Assembler
export SparseMatrixAssembler
export assemble
export assemble!
export push_coo!
export sparse_from_coo

"""
Abstract assembly operator
Expand Down Expand Up @@ -80,21 +83,28 @@ function assemble!(r,a::Assembler,cv::CellMatrix)
end

"""
Assembler that produces SparseMatrices from the SparseArrays package
Assembler that produces SparseMatrices from the SparseArrays and SparseMatricesCSR packages
"""
struct SparseMatrixAssembler{E} <: Assembler{SparseMatrixCSC{E,Int},Vector{E}}
struct SparseMatrixAssembler{E,M} <: Assembler{AbstractSparseMatrix{E,Int},Vector{E}}
testfesp::FESpace
trialfesp::FESpace

function SparseMatrixAssembler(
::Type{M},
testfesp::FESpace{D,Z,T},
trialfesp::FESpace{D,Z,T}) where {M<:AbstractSparseMatrix,D,Z,T}
E = eltype(T)
new{E,M}(testfesp,trialfesp)
end
end

function SparseMatrixAssembler(test::FESpace{D,Z,T}, trial::FESpace{D,Z,T}) where {D,Z,T}
E = eltype(T)
SparseMatrixAssembler{E}(trial,test)
SparseMatrixAssembler(SparseMatrixCSC, trial,test)
end

function assemble(
this::SparseMatrixAssembler{E},
vals::Vararg{Tuple{<:CellVector,<:CellNumber}}) where E
this::SparseMatrixAssembler{E,M},
vals::Vararg{Tuple{<:CellVector,<:CellNumber}}) where {E,M}

n = num_free_dofs(this.testfesp)
vec = zeros(E,n)
Expand All @@ -104,8 +114,8 @@ end

function assemble!(
vec::Vector{E},
this::SparseMatrixAssembler{E},
allvals::Vararg{Tuple{<:CellVector,<:CellNumber}}) where E
this::SparseMatrixAssembler{E,M},
allvals::Vararg{Tuple{<:CellVector,<:CellNumber}}) where {E,M}

vec .= zero(E)
rows = celldofids(this.testfesp)
Expand All @@ -127,35 +137,32 @@ function _assemble_vector!(vec,vals,rows)
end

function assemble(
this::SparseMatrixAssembler{E},
allvals::Vararg{Tuple{<:CellMatrix,<:CellNumber,<:CellNumber}}) where E
this::SparseMatrixAssembler{E,M},
allvals::Vararg{Tuple{<:CellMatrix,<:CellNumber,<:CellNumber}}) where {E,M}

I = Int
aux_row = I[]; aux_col = I[]; aux_val = E[]

_rows_m = celldofids(this.testfesp)
_cols_m = celldofids(this.trialfesp)

for (vals, cellids_row, cellids_col) in allvals
_vals = apply_constraints_rows(this.testfesp, vals, cellids_row)
rows_m = reindex(_rows_m, cellids_row)
vals_m = apply_constraints_cols(this.trialfesp, _vals, cellids_col)
cols_m = reindex(_cols_m, cellids_col)
_assemble_sparse_matrix_values!(
_assemble_sparse_matrix_values!(M,
aux_row,aux_col,aux_val,vals_m,rows_m,cols_m)
end
sparse(aux_row,aux_col,aux_val)
sparse_from_coo(M,aux_row,aux_col,aux_val)
end

function _assemble_sparse_matrix_values!(aux_row,aux_col,aux_val,vals,rows,cols)
function _assemble_sparse_matrix_values!(::Type{M},aux_row,aux_col,aux_val,vals,rows,cols) where {M}
for (rows_c, cols_c, vals_c) in zip(rows,cols,vals)
for (j,gidcol) in enumerate(cols_c)
if gidcol > 0
for (i,gidrow) in enumerate(rows_c)
if gidrow > 0
push!(aux_row, gidrow)
push!(aux_col, gidcol)
push!(aux_val, vals_c[i,j])
push_coo!(M,aux_row, aux_col, aux_val, gidrow, gidcol, vals_c[i,j])
end
end
end
Expand All @@ -164,13 +171,26 @@ function _assemble_sparse_matrix_values!(aux_row,aux_col,aux_val,vals,rows,cols)
end

function assemble!(
mat::SparseMatrixCSC{E},
this::SparseMatrixAssembler{E},
vals::Vararg{Tuple{<:CellMatrix,<:CellNumber,<:CellNumber}}) where E
mat::AbstractSparseMatrix{E,Int},
this::SparseMatrixAssembler{E,M},
vals::Vararg{Tuple{<:CellMatrix,<:CellNumber,<:CellNumber}}) where {E,M}
# This routine can be optimized a lot taking into a count the sparsity graph of mat
# For the moment we create an intermediate matrix and then transfer the nz values
m = assemble(this,vals...)
mat.nzval .= m.nzval
nonzeros(mat) .= nonzeros(m)
end

function sparse_from_coo(::Type{<:SparseMatrixCSC}, args...)
sparse(args...)
end

function sparse_from_coo(::Type{<:SparseMatrixCSR}, args...)
sparsecsr(args...)
end

function sparse_from_coo(::Type{<:SymSparseMatrixCSR}, args...)
symsparsecsr(args...)
end


end # module
80 changes: 59 additions & 21 deletions src/MultiField/MultiAssemblers.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -12,6 +12,9 @@ export MultiSparseMatrixAssembler
import Gridap.Assemblers: assemble
import Gridap.Assemblers: assemble!
import Gridap.Assemblers: SparseMatrixAssembler
import Gridap.Assemblers: SparseMatrixAssembler
import Gridap.Assemblers: push_coo!
import Gridap.Assemblers: sparse_from_coo

abstract type MultiAssembler{M<:AbstractMatrix,V<:AbstractVector} end

Expand Down Expand Up @@ -70,34 +73,70 @@ end
"""
MultiAssembler that produces SparseMatrices from the SparseArrays package
"""
struct MultiSparseMatrixAssembler{E} <: MultiAssembler{SparseMatrixCSC{E,Int},Vector{E}}
struct MultiSparseMatrixAssembler{E,M} <: MultiAssembler{AbstractSparseMatrix{E,Int},Vector{E}}
testfesps::MultiFESpace{E}
trialfesps::MultiFESpace{E}

function MultiSparseMatrixAssembler(
::Type{M},
testfesp::MultiFESpace{E},
trialfesp::MultiFESpace{E}) where {M<:AbstractSparseMatrix,E}
new{E,M}(testfesp,trialfesp)
end
end

function MultiSparseMatrixAssembler(
testfesp::MultiFESpace{E},
trialfesp::MultiFESpace{E}) where {E}
MultiSparseMatrixAssembler(SparseMatrixCSC, testfesps, trialfesp)
end

function MultiSparseMatrixAssembler(
::Type{M},
testfesps::Vector{<:FESpaceWithDirichletData},
trialfesps::Vector{<:FESpaceWithDirichletData})
trialfesps::Vector{<:FESpaceWithDirichletData}) where {M<:AbstractSparseMatrix}
V = MultiFESpace(testfesps)
U = MultiFESpace(trialfesps)
MultiSparseMatrixAssembler(V,U)
MultiSparseMatrixAssembler(M,V,U)
end

function MultiSparseMatrixAssembler(
testfesps::Vector{<:FESpaceWithDirichletData},
trialfesps::Vector{<:FESpaceWithDirichletData})
MultiSparseMatrixAssembler(SparseMatrixCSC,V,U)
end

function SparseMatrixAssembler(
::Type{M},
testfesps::Vector{<:FESpaceWithDirichletData},
trialfesps::Vector{<:FESpaceWithDirichletData}) where {M<:AbstractSparseMatrix}
MultiSparseMatrixAssembler(M,testfesps,trialfesps)
end

function SparseMatrixAssembler(
testfesps::Vector{<:FESpaceWithDirichletData},
trialfesps::Vector{<:FESpaceWithDirichletData})
MultiSparseMatrixAssembler(testfesps,trialfesps)
MultiSparseMatrixAssembler(SparseMatrixCSC,testfesps,trialfesps)
end

function MultiSparseMatrixAssembler(
testfesps::MultiFESpace{E,S}, trialfesps::MultiFESpace{E,S}) where {E,S}
::Type{M},
testfesps::MultiFESpace{E,S},
trialfesps::MultiFESpace{E,S}) where {M<:AbstractSparseMatrix,E,S}
@notimplementedif S != ConsequtiveMultiFieldStyle
MultiSparseMatrixAssembler{E}(testfesps,trialfesps)
MultiSparseMatrixAssembler{E,M}(M,testfesps,trialfesps)
end

function MultiSparseMatrixAssembler(
testfesps::MultiFESpace{E,S},
trialfesps::MultiFESpace{E,S}) where {E,S}
@notimplementedif S != ConsequtiveMultiFieldStyle
MultiSparseMatrixAssembler(SparseMatrixCSC,testfesps,trialfesps)
end

function assemble(
this::MultiSparseMatrixAssembler{E},
allvals::Vararg{Tuple{<:MultiCellVector,<:CellNumber}}) where E
this::MultiSparseMatrixAssembler{E,M},
allvals::Vararg{Tuple{<:MultiCellVector,<:CellNumber}}) where {E,M}

n = num_free_dofs(this.testfesps)
vec = zeros(E,n)
Expand All @@ -107,8 +146,8 @@ end

function assemble!(
vec::Vector{E},
this::MultiSparseMatrixAssembler{E},
allmcv::Vararg{Tuple{<:MultiCellVector,<:CellNumber}}) where E
this::MultiSparseMatrixAssembler{E,M},
allmcv::Vararg{Tuple{<:MultiCellVector,<:CellNumber}}) where {E,M}

vec .= zero(E)
V = this.testfesps
Expand All @@ -123,8 +162,8 @@ function assemble!(
end

function assemble(
this::MultiSparseMatrixAssembler{E},
allmcm::Vararg{Tuple{<:MultiCellMatrix,<:CellNumber,<:CellNumber}}) where {E}
this::MultiSparseMatrixAssembler{E,M},
allmcm::Vararg{Tuple{<:MultiCellMatrix,<:CellNumber,<:CellNumber}}) where {E,M}

V = this.testfesps
U = this.trialfesps
Expand All @@ -140,21 +179,21 @@ function assemble(
mf_vals = apply_constraints_cols(U,mf_vals,cellids_col)
mf_cols = reindex(_mf_cols, cellids_col)
i_to_fieldid = mf_vals.fieldids
_assemble_mat!(
_assemble_mat!(M,
aux_row,aux_col,aux_val,mf_rows,mf_cols,mf_vals,
i_to_fieldid,offsets_row,offsets_col)
end
sparse(aux_row, aux_col, aux_val)
sparse_from_coo(M,aux_row, aux_col, aux_val)
end

function assemble!(
mat::SparseMatrixCSC{E},
mat::AbstractSparseMatrix{E},
this::MultiSparseMatrixAssembler{E},
allvals::Vararg{Tuple{<:MultiCellMatrix,<:CellNumber,<:CellNumber}}) where E
# This routine can be optimized a lot taking into a count the sparsity graph of mat
# For the moment we create an intermediate matrix and then transfer the nz values
m = assemble(this,allvals...)
mat.nzval .= m.nzval
nonzeros(mat) .= nonzeros(m)
end

function _assemble_vec!(vec,mf_rows,mf_vals,i_to_fieldid,offsets)
Expand All @@ -172,9 +211,9 @@ function _assemble_vec!(vec,mf_rows,mf_vals,i_to_fieldid,offsets)
end
end

function _assemble_mat!(
function _assemble_mat!(::Type{M},
aux_row,aux_col,aux_val,mf_rows,mf_cols,mf_vals,
i_to_fieldid,offsets_row,offsets_col)
i_to_fieldid,offsets_row,offsets_col) where {M<:AbstractSparseMatrix}

for (mf_rows_c,mf_cols_c,mf_vals_c) in zip(mf_rows,mf_cols,mf_vals)
for (i,vals_c) in enumerate(mf_vals_c)
Expand All @@ -187,9 +226,8 @@ function _assemble_mat!(
if gidcol > 0
for (lidrow,gidrow) in enumerate(rows_c)
if gidrow > 0
push!(aux_row, gidrow+offset_row)
push!(aux_col, gidcol+offset_col)
push!(aux_val, vals_c[lidrow,lidcol])
push_coo!(M,aux_row, aux_col, aux_val,
gidrow+offset_row, gidcol+offset_col, vals_c[lidrow,lidcol])
end
end
end
Expand Down
Loading