Merge pull request #56 from fverdugo/example004

Example004
GalerkinToolkit · Feb 9, 2024 · 174cd60 · 174cd60
2 parents 075bc4f + 2430b40
commit 174cd60
Show file tree

Hide file tree

Showing 8 changed files with 859 additions and 76 deletions.
diff --git a/extensions/GalerkinToolkitExamples/src/GalerkinToolkitExamples.jl b/extensions/GalerkinToolkitExamples/src/GalerkinToolkitExamples.jl
@@ -3,5 +3,6 @@ module GalerkinToolkitExamples
 include("example001.jl")
 include("example002.jl")
 include("example003.jl")
+include("example004.jl")
 
 end # module
diff --git a/extensions/GalerkinToolkitExamples/src/example003.jl b/extensions/GalerkinToolkitExamples/src/example003.jl
@@ -6,7 +6,7 @@ using StaticArrays
 using LinearAlgebra
 using SparseArrays
 using WriteVTK
-using PartitionedArrays: JaggedArray, val_parameter, nzindex
+using PartitionedArrays: JaggedArray, sparse_matrix, sparse_matrix!
 using TimerOutputs
 using NLsolve
 using Random
@@ -65,7 +65,7 @@ function default_params()
     params[:neumann_tags] = String[]
     params[:integration_degree] = 2
     params[:p] = 2
-    params[:solver] = nlsolve_solver()
+    params[:solver] = nlsolve_solver(;method=:newton)
     params[:example_path] = joinpath(outdir,"example003")
     params[:export_vtu] = true
     params[:autodiff] = :hand
@@ -696,77 +696,4 @@ function export_results(uh,params,state)
     nothing
 end
 
-# TODO move this to partitioned arrays
-struct FilteredCooVector{F,A,B,C,T} <: AbstractVector{T}
-    f::F
-    I::A
-    J::B
-    V::C
-    function FilteredCooVector(f::F,I::A,J::B,V::C) where {F,A,B,C}
-        T = eltype(C)
-        new{F,A,B,C,T}(f,I,J,V)
-    end
-end
-Base.size(a::FilteredCooVector) = size(a.V)
-Base.IndexStyle(::Type{<:FilteredCooVector}) = IndexLinear()
-Base.@propagate_inbounds function Base.getindex(a::FilteredCooVector,k::Int)
-    i = a.I[k]
-    j = a.J[k]
-    v = a.V[k]
-    if i < 1 || j < 1
-        return a.f(v)
-    end
-    v
-end
-
-function sparse_matrix(I,J,V,m,n;kwargs...)
-    sparse_matrix(sparse,I,J,V,m,n;kwargs...)
-end
-function sparse_matrix(f,I,J,V,m,n;reuse=Val(false),skip_out_of_bounds=true)
-    if !skip_out_of_bounds
-        I2 = I
-        J2 = J
-        V2 = V
-    elseif m*n == 0
-        Ti = eltype(I)
-        T = eltype(V)
-        I2 = Ti[]
-        J2 = Ti[]
-        V2 = Tv[]
-    else
-        I2 = FilteredCooVector(one,I,J,I)
-        J2 = FilteredCooVector(one,I,J,J)
-        V2 = FilteredCooVector(zero,I,J,V)
-    end
-    A = f(I2,J2,V2,m,n)
-    if val_parameter(reuse)
-        K = precompute_nzindex(A,I,J)
-        return A,K
-    end
-    A
-end
-
-function precompute_nzindex(A,I,J)
-    K = zeros(Int32,length(I))
-    for (p,(i,j)) in enumerate(zip(I,J))
-        if i < 1 || j < 1
-            continue
-        end
-        K[p] = nzindex(A,i,j)
-    end
-    K
-end
-
-function sparse_matrix!(A,V,K)
-    LinearAlgebra.fillstored!(A,0)
-    A_nz = nonzeros(A)
-    for (k,v) in zip(K,V)
-        if k < 1
-            continue
-        end
-        A_nz[k] += v
-    end
-    A
-end
-
 end # module