diff --git a/src/mat.jl b/src/mat.jl index d8f1b32f..e05e0152 100644 --- a/src/mat.jl +++ b/src/mat.jl @@ -16,27 +16,11 @@ function destroy(M::AbstractMat{PetscLib}) where {PetscLib} return nothing end -""" - MatSeqAIJ{PetscLib, PetscScalar} - -PETSc sparse array using AIJ format (also known as a compressed sparse row or -CSR format). - -Memory allocation is handled by PETSc. - -# External Links -$(_doc_external("Mat/MatCreateSeqAIJ")) -""" -mutable struct MatSeqAIJ{PetscLib, PetscScalar} <: - AbstractMat{PetscLib, PetscScalar} - ptr::CMat -end - """ MatSeqAIJ(petsclib, num_rows, num_cols, nonzeros) -Create a PETSc sparse array using AIJ format (also known as a compressed sparse -row or CSR format) of size `num_rows X num_cols` with `nonzeros` per row +Create a PETSc serial sparse array using AIJ format (also known as a compressed +sparse row or CSR format) of size `num_rows X num_cols` with `nonzeros` per row If `nonzeros` is an `Integer` the same number of non-zeros will be used for each row, if `nonzeros` is a `Vector{PetscInt}` then one value must be specified for @@ -47,6 +31,11 @@ Memory allocation is handled by PETSc and garbage collection can be used. # External Links $(_doc_external("Mat/MatCreateSeqAIJ")) """ +mutable struct MatSeqAIJ{PetscLib, PetscScalar} <: + AbstractMat{PetscLib, PetscScalar} + ptr::CMat +end + function MatSeqAIJ( petsclib::PetscLib, num_rows::Integer, @@ -112,53 +101,54 @@ function MatSeqDense( end """ - assemble(M::AbstractMat[, t::MatAssemblyType = MAT_FINAL_ASSEMBLY) + assemble!(M::AbstractMat[, t::MatAssemblyType = MAT_FINAL_ASSEMBLY) Assemble the matrix `M` with assembly type `t`. -For overlapping assembly see [`assemblybegin`](@ref) and [`assemblyend`](@ref) +For overlapping assembly see [`assemblybegin!`](@ref) and [`assemblyend!`](@ref) # External Links $(_doc_external("Mat/MatAssemblyBegin")) $(_doc_external("Mat/MatAssemblyEnd")) """ -function assemble(M::AbstractMat, t::MatAssemblyType = MAT_FINAL_ASSEMBLY) - assemblybegin(M, t) - assemblyend(M, t) +function assemble!(M::AbstractMat, t::MatAssemblyType = MAT_FINAL_ASSEMBLY) + assemblybegin!(M, t) + assemblyend!(M, t) + return M end """ - assemblybegin(M::AbstractMat[, t::MatAssemblyType = MAT_FINAL_ASSEMBLY) + assemblybegin!(M::AbstractMat[, t::MatAssemblyType = MAT_FINAL_ASSEMBLY) Begin assembly of the matrix `M` with assembly type `t`; finished with -[`assemblyend`](@ref). +[`assemblyend!`](@ref). # External Links $(_doc_external("Mat/MatAssemblyBegin")) """ -function assemblybegin( +function assemblybegin!( M::AbstractMat{PetscLib}, t::MatAssemblyType = MAT_FINAL_ASSEMBLY, ) where {PetscLib} LibPETSc.MatAssemblyBegin(PetscLib, M, t) - return nothing + return M end """ - assemblyend(M::AbstractMat[, t::MatAssemblyType = MAT_FINAL_ASSEMBLY) + assemblyend!(M::AbstractMat[, t::MatAssemblyType = MAT_FINAL_ASSEMBLY) Finish assembly of the matrix `M` with assembly type `t`; start assembly with -[`assemblybegin`](@ref). +[`assemblybegin!`](@ref). # External Links $(_doc_external("Mat/MatAssemblyEnd")) """ -function assemblyend( +function assemblyend!( M::AbstractMat{PetscLib}, t::MatAssemblyType = MAT_FINAL_ASSEMBLY, ) where {PetscLib} LibPETSc.MatAssemblyEnd(PetscLib, M, t) - return nothing + return M end function Base.size(A::AbstractMat{PetscLib}) where {PetscLib} @@ -193,7 +183,7 @@ Base.show(io::IO, ::MIME"text/plain", mat::AbstractMat) = _show(io, mat) row0idxs::Vector{PetscInt}, col0idxs::Vector{PetscInt}, rowvals::Array{PetscScalar}, - insertmode::InsertMode; + insertmode::InsertMode = INSERT_VALUES; num_rows = length(row0idxs), num_cols = length(col0idxs) ) @@ -212,7 +202,7 @@ function setvalues!( row0idxs::Vector{PetscInt}, col0idxs::Vector{PetscInt}, rowvals::Array{PetscScalar}, - insertmode::InsertMode; + insertmode::InsertMode = INSERT_VALUES; num_rows = length(row0idxs), num_cols = length(col0idxs), ) where {PetscLib, PetscScalar, PetscInt} diff --git a/src/vec.jl b/src/vec.jl index 480ae2db..e35bb221 100644 --- a/src/vec.jl +++ b/src/vec.jl @@ -121,7 +121,8 @@ end An sequentially-stored MPI PETSc vector for `petsclib.PetscScalar` of local length `local_length` and global length `global_length` without ghost elements. -If `global_length isa Int` then `local_length` can be set to `PETSC_DECIDE`. +If `global_length isa Int` then `local_length` can be set to `PETSC_DECIDE` in +which case PETSc will decide the local_length. # External Links $(_doc_external("Vec/VecCreateMPI")) diff --git a/test/mat.jl b/test/mat.jl index ef0b14a2..c2254abe 100644 --- a/test/mat.jl +++ b/test/mat.jl @@ -14,14 +14,14 @@ using LinearAlgebra: norm, mul!, Adjoint, Transpose, issymmetric, ishermitian @test size(A) == (num_rows, num_cols) B = PETSc.MatSeqAIJ(petsclib, num_rows, num_cols, nz_int) - PETSc.assemble(A) - PETSc.assemble(B) + PETSc.assemble!(A) + PETSc.assemble!(B) # both empty @test A == B C = PETSc.MatSeqAIJ(petsclib, num_rows, num_cols, nz_vec) @test size(A) == (num_rows, num_cols) - PETSc.assemble(C) + PETSc.assemble!(C) # both empty @test A == C @@ -30,7 +30,7 @@ using LinearAlgebra: norm, mul!, Adjoint, Transpose, issymmetric, ishermitian D[1, [1, 2]] .= [1, 2] D[2, [3, 4]] .= [3, 4] D[5, [3, 4]] .= [3, 4] - PETSc.assemble(D) + PETSc.assemble!(D) DJ = zeros(PetscScalar, num_rows, num_cols) DJ[1, [1, 2]] .= [1, 2] @@ -41,7 +41,7 @@ using LinearAlgebra: norm, mul!, Adjoint, Transpose, issymmetric, ishermitian D[1, [1, 2]] .= [1, 2im] D[2, [3, 4]] .= [3, 4im] D[5, [3, 4]] .= [3, 4im] - PETSc.assemble(D) + PETSc.assemble!(D) DJ = zeros(PetscScalar, num_rows, num_cols) DJ[1, [1, 2]] .= [1, 2im] @@ -53,7 +53,7 @@ using LinearAlgebra: norm, mul!, Adjoint, Transpose, issymmetric, ishermitian E[2, [1, 3, 4]] .= [1, 3, 4] E[3, [3, 4]] .= [3, 4] E[4, [5]] .= [6] - PETSc.assemble(E) + PETSc.assemble!(E) @test A != E @test D != E @@ -87,13 +87,13 @@ using LinearAlgebra: norm, mul!, Adjoint, Transpose, issymmetric, ishermitian A[1, 1] = 1 A[2, 1] = -2 A[1, 2] = -2 - PETSc.assemble(A) + PETSc.assemble!(A) B = PETSc.MatSeqAIJ(petsclib, 5, 5, 2) B[1, 1] = 1 B[2, 1] = 2 B[1, 2] = -2 - PETSc.assemble(B) + PETSc.assemble!(B) @test issymmetric(A) @test ishermitian(A) @@ -104,13 +104,13 @@ using LinearAlgebra: norm, mul!, Adjoint, Transpose, issymmetric, ishermitian A[1, 1] = 1 A[2, 1] = -2 + im A[1, 2] = -2 - im - PETSc.assemble(A) + PETSc.assemble!(A) B = PETSc.MatSeqAIJ(petsclib, 5, 5, 2) B[1, 1] = 1 B[2, 1] = -2 + im B[1, 2] = -2 + im - PETSc.assemble(B) + PETSc.assemble!(B) @test !issymmetric(A) @test ishermitian(A) @test issymmetric(B)