Merge pull request #680 from gridap/replace_fill_entries_scale_entrie…

…s_by_stlib_functions Replacing fill_stored! by fill! or LinearAlgebra.fillstored! for arrays
gridap · Oct 20, 2021 · d24f5dd · d24f5dd
2 parents 44af28b + a743312
commit d24f5dd
Show file tree

Hide file tree

Showing 15 changed files with 29 additions and 67 deletions.
diff --git a/Project.toml b/Project.toml
@@ -43,8 +43,8 @@ LineSearches = "7.0.1"
 NLsolve = "4.3.0"
 NearestNeighbors = "0.4.8"
 QuadGK = "2.3.1, 2.4"
-SparseMatricesCSR = "0.6"
 StaticArrays = "0.12.1, 1.0"
+SparseMatricesCSR = "0.6.4"
 WriteVTK = "1.7, 1.8"
 julia = "1.3"
 

diff --git a/src/Algebra/Algebra.jl b/src/Algebra/Algebra.jl
@@ -26,7 +26,6 @@ export allocate_matrix_and_vector
 export allocate_in_domain
 export allocate_in_range
 export add_entries!
-export scale_entries!
 export muladd!
 export nz_counter
 export nz_allocation
@@ -45,7 +44,6 @@ export is_entry_stored
 export finalize_coo!
 export sparse_from_coo
 export add_entry!
-export fill_entries!
 export copy_entries!
 export allocate_coo_vectors
 export push_coo!

diff --git a/src/Algebra/AlgebraInterfaces.jl b/src/Algebra/AlgebraInterfaces.jl
@@ -61,17 +61,6 @@ function allocate_in_domain(::Type{V},matrix) where V
   allocate_vector(V,n)
 end
 
-"""
-    fill_entries!(a,v)
-
-Fill the entries of array `a` with the value `v`. Returns `a`.
-For sparse matrices it only fills the non-zero entries.
-"""
-function fill_entries!(a,v)
-  fill!(a,v)
-  a
-end
-
 """
     copy_entries!(a,b)
 
@@ -196,19 +185,6 @@ end
   A
 end
 
-"""
-    scale_entries!(a,v)
-
-Scale the entries of array `a` with the value `v`. Returns `a`.
-"""
-function scale_entries!(a,b)
-  @inbounds for i in eachindex(a)
-    a[i] = b*a[i]
-  end
-  a
-end
-
-# Base.mul!
 
 """
     muladd!(c,a,b)
@@ -280,7 +256,7 @@ end
 # We can also do a loop and update
 # the entries of c
 #
-#    fill_entries!(c,0)
+#    fill!(c,0) or LinearAlgebra.fillstored!(c,0)
 #    add_entry!(c,v,i,j)
 #    add_entries!(c,vs,is,js)
 #
@@ -466,13 +442,7 @@ function copy_entries!(a::T,b::T) where T<:AbstractSparseMatrix
   end
 end
 
-function fill_entries!(A::AbstractSparseMatrix,v)
-  nonzeros(A) .= v
-  A
-end
-
 function allocate_coo_vectors(
    ::Type{<:AbstractSparseMatrix{Tv,Ti}},n::Integer) where {Tv,Ti}
   (zeros(Ti,n), zeros(Ti,n), zeros(Tv,n))
 end
-
diff --git a/src/Algebra/LinearSolvers.jl b/src/Algebra/LinearSolvers.jl
@@ -43,7 +43,7 @@ end
 
 function zero_initial_guess(op::AffineOperator)
   x = allocate_in_domain(typeof(op.vector),op.matrix)
-  fill_entries!(x,zero(eltype(x)))
+  fill!(x,zero(eltype(x)))
   x
 end
 
@@ -174,12 +174,12 @@ function solve!(x::AbstractVector,
                 ls::LinearSolver,
                 op::NonlinearOperator,
                 cache::Nothing)
-  fill_entries!(x,zero(eltype(x)))
+  fill!(x,zero(eltype(x)))
   b = residual(op, x)
   A = jacobian(op, x)
   ss = symbolic_setup(ls, A)
   ns = numerical_setup(ss,A)
-  scale_entries!(b,-1)
+  rmul!(b,-1)
   solve!(x,ns,b)
   LinearSolverCache(A,b,ns)
 end
@@ -188,13 +188,13 @@ function solve!(x::AbstractVector,
                 ls::LinearSolver,
                 op::NonlinearOperator,
                 cache)
-  fill_entries!(x,zero(eltype(x)))
+  fill!(x,zero(eltype(x)))
   b = cache.b
   A = cache.A
   ns = cache.ns
   residual!(b, op, x)
   numerical_setup!(ns,A)
-  scale_entries!(b,-1)
+  ruml!(b,-1)
   solve!(x,ns,b)
   cache
 end

diff --git a/src/Algebra/NLSolvers.jl b/src/Algebra/NLSolvers.jl
@@ -54,7 +54,7 @@ function _solve_nr!(x,A,b,dx,ns,nls,op)
   for nliter in 1:nls.max_nliters
 
     # Solve linearized problem
-    scale_entries!(b,-1)
+    rmul!(b,-1)
     solve!(dx,ns,b)
     x .+= dx
 
@@ -183,5 +183,3 @@ function _update_nlsolve_cache!(cache,x0,op)
   numerical_setup!(ns,j0)
   NLSolversCache(f0,j0,df,ns,nothing)
 end
-
-
diff --git a/src/Arrays/Arrays.jl b/src/Arrays/Arrays.jl
@@ -26,7 +26,6 @@ import Base: similar
 import Base: IndexStyle
 
 import Gridap.Algebra: scale_entries!
-import Gridap.Algebra: fill_entries!
 
 # CachedArray
 

diff --git a/src/CellData/AttachDirichlet.jl b/src/CellData/AttachDirichlet.jl
@@ -83,11 +83,9 @@ end
   cm, cv = cache
   if mask
     vec_with_bcs = evaluate!(cm,*,mat,vals)
-    scale_entries!(vec_with_bcs,-1)
+    rmul!(vec_with_bcs,-1)
   else
     vec_with_bcs = evaluate!(cv,ZeroVectorMap(),mat)
   end
   (mat, vec_with_bcs)
 end
-
-
diff --git a/src/CellData/CellData.jl b/src/CellData/CellData.jl
@@ -32,7 +32,7 @@ import Gridap.Geometry: num_cells
 import Gridap.Geometry: get_triangulation
 
 import Gridap.TensorValues: inner, outer, double_contraction, symmetric_part
-import LinearAlgebra: det, tr, cross, dot, ⋅
+import LinearAlgebra: det, tr, cross, dot, ⋅, rmul!
 import Base: inv, abs, abs2, *, +, -, /, adjoint, transpose, real, imag, conj
 
 export gradient, ∇

diff --git a/src/FESpaces/SparseMatrixAssemblers.jl b/src/FESpaces/SparseMatrixAssemblers.jl
@@ -63,7 +63,7 @@ function allocate_vector(a::SparseMatrixAssembler,vecdata)
 end
 
 function assemble_vector!(b,a::SparseMatrixAssembler,vecdata)
-  fill_entries!(b,zero(eltype(b)))
+  fill!(b,zero(eltype(b)))
   assemble_vector_add!(b,a,vecdata)
 end
 
@@ -91,7 +91,7 @@ function allocate_matrix(a::SparseMatrixAssembler,matdata)
 end
 
 function assemble_matrix!(mat,a::SparseMatrixAssembler,matdata)
-  fill_entries!(mat,zero(eltype(mat)))
+  LinearAlgebra.fillstored!(mat,zero(eltype(mat)))
   assemble_matrix_add!(mat,a,matdata)
 end
 
@@ -122,8 +122,8 @@ function allocate_matrix_and_vector(a::SparseMatrixAssembler,data)
 end
 
 function assemble_matrix_and_vector!(A,b,a::SparseMatrixAssembler, data)
-  fill_entries!(A,zero(eltype(A)))
-  fill_entries!(b,zero(eltype(b)))
+  LinearAlgebra.fillstored!(A,zero(eltype(A)))
+  fill!(b,zero(eltype(b)))
   assemble_matrix_and_vector_add!(A,b,a,data)
 end
 
@@ -347,4 +347,3 @@ end
     evaluate!(add_vec_cache,add!,b,vecvals,rows)
   end
 end
-
diff --git a/src/Fields/ArrayBlocks.jl b/src/Fields/ArrayBlocks.jl
@@ -1022,10 +1022,10 @@ function Base.:*(a::ArrayBlock{A,2},b::ArrayBlock{B,2}) where {A,B}
   ArrayBlock(array,touched)
 end
 
-function Algebra.scale_entries!(a::ArrayBlock,β)
+function LinearAlgebra.rmul!(a::ArrayBlock,β)
   for i in eachindex(a.touched)
     if a.touched[i]
-      scale_entries!(a.array[i],β)
+      rmul!(a.array[i],β)
     end
   end
 end
@@ -1057,7 +1057,7 @@ function LinearAlgebra.mul!(
   @check nj == size(b.array,1)
   for i in 1:ni
     if β!=1 && c.touched[i]
-      scale_entries!(c.array[i],β)
+      rmul!(c.array[i],β)
     end
     for j in 1:nj
       if a.touched[i,j] && b.touched[j]
@@ -1085,7 +1085,7 @@ function LinearAlgebra.mul!(
   for i in 1:ni
     for j in 1:nj
       if β!=1 && c.touched[i,j]
-        scale_entries!(c.array[i,j],β)
+        rmul!(c.array[i,j],β)
       end
       for k in 1:nk
         if a.touched[i,k] && b.touched[k,j]

diff --git a/src/Fields/Fields.jl b/src/Fields/Fields.jl
@@ -10,7 +10,6 @@ using Gridap.Helpers: @abstractmethod, @notimplemented
 using Gridap.Helpers: @notimplementedif, @unreachable, @check
 
 using Gridap.Algebra: mul!
-using Gridap.Algebra: fill_entries!
 
 using Gridap.TensorValues
 using Gridap.Algebra

diff --git a/test/AlgebraTests/AlgebraInterfacesTests.jl b/test/AlgebraTests/AlgebraInterfacesTests.jl
@@ -1,6 +1,7 @@
 module AlgebraInterfacesTests
 
 using Gridap.Algebra
+using LinearAlgebra
 using Test
 
 a = allocate_vector(Vector{Int},10)
@@ -19,7 +20,7 @@ b = allocate_in_range(Vector{Int},a)
 b = allocate_in_domain(Vector{Int},a)
 @test length(b) == 6
 
-fill_entries!(b,4)
+fill!(b,4)
 @test all( b .== 4 )
 
 a = rand(6)
@@ -44,7 +45,7 @@ a .-= b
 
 a = rand(6)
 c = copy(a)
-scale_entries!(a,10)
+rmul!(a,10)
 @test all( a .== 10*c)
 
 a = rand(4,6)
@@ -139,7 +140,7 @@ for A in (SparseMatrixCSC{Float64,Int},SparseMatrixCSC{Float64,Int32})
   add_entry!(c,nothing,1,1)
   add_entry!(c,nothing,3,1)
   add_entry!(c,3.0,4,9)
-  
+
   a = nz_counter(builder,(rows,cols))
   add_entries!(a,[1.0 -1.0; -1.0 1.0],[1,-1],[-1,1])
   add_entries!(a,nothing,[1,1],[1,-1])
@@ -188,7 +189,7 @@ for A in (SparseMatrixCSC{Float64,Int},SparseMatrixCSC{Float64,Int32})
   add_entry!(c,nothing,1,1)
   add_entry!(c,nothing,3,1)
   add_entry!(c,3.0,4,9)
-  
+
   a = nz_counter(builder,(rows,cols))
   add_entries!(a,[1.0 -1.0; -1.0 1.0],[1,-1],[-1,1])
   add_entries!(a,nothing,[1,1],[1,-1])
@@ -256,7 +257,7 @@ for A in (SparseMatrixCSR{1,Float64,Int},SparseMatrixCSR{0,Float64,Int32})
   add_entry!(c,nothing,1,1)
   add_entry!(c,nothing,1,3)
   add_entry!(c,3.0,9,4)
-  
+
   a = nz_counter(builder,(rows,cols))
   add_entries!(a,[1.0 -1.0; -1.0 1.0],[1,-1],[-1,1])
   add_entries!(a,nothing,[1,1],[1,-1])
@@ -305,7 +306,7 @@ for A in (SparseMatrixCSR{0,Float64,Int},SparseMatrixCSR{1,Float64,Int32})
   add_entry!(c,nothing,1,1)
   add_entry!(c,nothing,3,1)
   add_entry!(c,3.0,9,4)
-  
+
   a = nz_counter(builder,(rows,cols))
   add_entries!(a,[1.0 -1.0; -1.0 1.0],[1,-1],[-1,1])
   add_entries!(a,nothing,[1,1],[1,-1])
@@ -362,7 +363,7 @@ for A in (
   add_entry!(c,3.0,2,2)
   add_entry!(c,3.0,2,6)
   add_entry!(c,3.0,6,2)
-  
+
   a = nz_counter(builder,(rows,cols))
   add_entries!(a,[1.0 -1.0; -1.0 1.0],[1,-1],[-1,1])
   add_entries!(a,nothing,[1,1],[1,-1])

diff --git a/test/AlgebraTests/SparseMatrixCSCTests.jl b/test/AlgebraTests/SparseMatrixCSCTests.jl
@@ -56,7 +56,7 @@ for T in (Int32,Int,Float32,Float64)
   add_entry!(+,CSC,1,maxrows,maxcols)
   @test getindex(CSC,maxrows,maxcols) == vold+1
 
-  fill_entries!(CSC,0)
+  LinearAlgebra.fillstored!(CSC,0)
   @test all(x->x==0, nonzeros(CSC))
 
 end

diff --git a/test/AlgebraTests/SparseMatrixCSRTests.jl b/test/AlgebraTests/SparseMatrixCSRTests.jl
@@ -59,7 +59,7 @@ for Ti in int_types
       add_entry!(+,CSR,1,maxrows,maxcols)
       @test getindex(CSR,maxrows,maxcols) == vold+1
 
-      fill_entries!(CSR,0)
+      LinearAlgebra.fillstored!(CSR,0)
       @test all(x->x==0, nonzeros(CSR))
 
     end

diff --git a/test/AlgebraTests/SymSparseMatrixCSRTests.jl b/test/AlgebraTests/SymSparseMatrixCSRTests.jl
@@ -57,7 +57,7 @@ for Ti in int_types
       add_entry!(+,SYMCSR,1,maxrows,maxcols)
       @test getindex(SYMCSR,maxrows,maxcols) == vold+1
 
-      fill_entries!(SYMCSR,0)
+      LinearAlgebra.fillstored!(SYMCSR,0)
       @test all(x->x==0, nonzeros(SYMCSR))
 
     end