Skip to content

Commit

Permalink
DofHandler with any <: AbstractGrid (#379)
Browse files Browse the repository at this point in the history
  • Loading branch information
koehlerson authored Apr 22, 2022
1 parent 9d40af5 commit b46853e
Show file tree
Hide file tree
Showing 4 changed files with 104 additions and 11 deletions.
17 changes: 9 additions & 8 deletions src/Dofs/DofHandler.jl
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@ Operates slightly faster than [`MixedDofHandler`](@docs). Supports:
- `Grid`s with a single concrete cell type.
- One or several fields on the whole domaine.
"""
struct DofHandler{dim,C,T} <: AbstractDofHandler
struct DofHandler{dim,T,G<:AbstractGrid{dim}} <: AbstractDofHandler
field_names::Vector{Symbol}
field_dims::Vector{Int}
# TODO: field_interpolations can probably be better typed: We should at least require
Expand All @@ -19,11 +19,11 @@ struct DofHandler{dim,C,T} <: AbstractDofHandler
cell_dofs::Vector{Int}
cell_dofs_offset::Vector{Int}
closed::ScalarWrapper{Bool}
grid::Grid{dim,C,T}
grid::G
ndofs::ScalarWrapper{Int}
end

function DofHandler(grid::Grid)
function DofHandler(grid::AbstractGrid)
isconcretetype(getcelltype(grid)) || error("Grid includes different celltypes. Use MixedDofHandler instead of DofHandler")
DofHandler(Symbol[], Int[], Interpolation[], BCValues{Float64}[], Int[], Int[], ScalarWrapper(false), grid, Ferrite.ScalarWrapper(-1))
end
Expand Down Expand Up @@ -302,8 +302,8 @@ function celldofs!(global_dofs::Vector{Int}, dh::DofHandler, i::Int)
return global_dofs
end

function cellnodes!(global_nodes::Vector{Int}, grid::Grid{dim,C}, i::Int) where {dim,C}
nodes = grid.cells[i].nodes
function cellnodes!(global_nodes::Vector{Int}, grid::AbstractGrid{dim}, i::Int) where {dim,C}
nodes = getcells(grid,i).nodes
N = length(nodes)
@assert length(global_nodes) == N
for j in 1:N
Expand All @@ -312,15 +312,16 @@ function cellnodes!(global_nodes::Vector{Int}, grid::Grid{dim,C}, i::Int) where
return global_nodes
end

function cellcoords!(global_coords::Vector{Vec{dim,T}}, grid::Grid{dim,C}, i::Int) where {dim,C,T}
nodes = grid.cells[i].nodes
function cellcoords!(global_coords::Vector{Vec{dim,T}}, grid::AbstractGrid{dim}, i::Int) where {dim,C,T}
nodes = getcells(grid,i).nodes
N = length(nodes)
@assert length(global_coords) == N
for j in 1:N
global_coords[j] = grid.nodes[nodes[j]].x
global_coords[j] = getcoordinates(getnodes(grid,nodes[j]))
end
return global_coords
end

cellcoords!(global_coords::Vector{<:Vec}, dh::DofHandler, i::Int) = cellcoords!(global_coords, dh.grid, i)

function celldofs(dh::DofHandler, i::Int)
Expand Down
6 changes: 3 additions & 3 deletions src/iterators.jl
Original file line number Diff line number Diff line change
Expand Up @@ -41,7 +41,7 @@ struct CellIterator{dim,C,T,DH<:Union{AbstractDofHandler,Nothing}}
dh::Union{DH,Nothing}
celldofs::Vector{Int}

function CellIterator{dim,C,T}(dh::Union{DofHandler{dim,C,T},MixedDofHandler{dim,T,G},Nothing}, cellset::Union{AbstractVector{Int},Nothing}, flags::UpdateFlags) where {dim,C,T,G}
function CellIterator{dim,C,T}(dh::Union{DofHandler{dim,T,G},MixedDofHandler{dim,T,G},Nothing}, cellset::Union{AbstractVector{Int},Nothing}, flags::UpdateFlags) where {dim,C,T,G}
isconcretetype(C) || _check_same_celltype(dh.grid, cellset)
N = nnodes_per_cell(dh.grid, cellset === nothing ? 1 : first(cellset))
cell = ScalarWrapper(0)
Expand All @@ -64,8 +64,8 @@ end

CellIterator(grid::Grid{dim,C,T}, cellset::Union{AbstractVector{Int},Nothing}=nothing, flags::UpdateFlags=UpdateFlags()) where {dim,C,T} =
CellIterator{dim,C,T}(grid, cellset, flags)
CellIterator(dh::DofHandler{dim,C,T}, cellset::Union{AbstractVector{Int},Nothing}=nothing, flags::UpdateFlags=UpdateFlags()) where {dim,C,T} =
CellIterator{dim,C,T}(dh, cellset, flags)
CellIterator(dh::DofHandler{dim,T}, cellset::Union{AbstractVector{Int},Nothing}=nothing, flags::UpdateFlags=UpdateFlags()) where {dim,C,T} =
CellIterator{dim,getcelltype(dh.grid),T}(dh, cellset, flags)
CellIterator(dh::MixedDofHandler{dim,T}, cellset::Union{AbstractVector{Int},Nothing}=nothing, flags::UpdateFlags=UpdateFlags()) where {dim,T} =
CellIterator{dim,getcelltype(dh.grid),T}(dh, cellset, flags)

Expand Down
1 change: 1 addition & 0 deletions test/runtests.jl
Original file line number Diff line number Diff line change
Expand Up @@ -16,6 +16,7 @@ include("test_assemble.jl")
include("test_dofs.jl")
include("test_constraints.jl")
include("test_grid_dofhandler_vtk.jl")
include("test_abstractgrid.jl")
include("test_mixeddofhandler.jl")
include("test_l2_projection.jl")
include("test_pointevaluation.jl")
Expand Down
91 changes: 91 additions & 0 deletions test/test_abstractgrid.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,91 @@
@testset "AbstractGrid" begin

struct SmallGrid{dim,N,C<:Ferrite.AbstractCell} <: Ferrite.AbstractGrid{dim}
nodes_test::Vector{NTuple{dim,Float64}}
cells_test::NTuple{N,C}
end

Ferrite.getcells(grid::SmallGrid) = grid.cells_test
Ferrite.getcells(grid::SmallGrid, v::Union{Int, Vector{Int}}) = grid.cells_test[v]
Ferrite.getncells(grid::SmallGrid{dim,N}) where {dim,N} = N
Ferrite.getcelltype(grid::SmallGrid) = eltype(grid.cells_test)
Ferrite.getcelltype(grid::SmallGrid, i::Int) = typeof(grid.cells_test[i])
Ferrite.getcoordinates(x::NTuple{dim,Float64}) where dim = Vec{dim,Float64}(x)

Ferrite.getnodes(grid::SmallGrid) = grid.nodes_test
Ferrite.getnodes(grid::SmallGrid, v::Union{Int, Vector{Int}}) = grid.nodes_test[v]
Ferrite.getnnodes(grid::SmallGrid) = length(grid.nodes_test)
Ferrite.nnodes_per_cell(grid::SmallGrid, i::Int=1) = Ferrite.nnodes(grid.cells_test[i])
Ferrite.n_faces_per_cell(grid::SmallGrid) = nfaces(eltype(grid.cells_test))
function Ferrite.getcoordinates!(x::Vector{Vec{dim,T}}, grid::SmallGrid, cell::Int) where {dim,T}
for i in 1:length(x)
x[i] = Vec{dim,T}(grid.nodes_test[grid.cells_test[cell].nodes[i]])
end
end
function Ferrite.getcoordinates(grid::SmallGrid{dim}, cell::Int) where dim
nodeidx = grid.cells_test[cell].nodes
return [Vec{dim,Float64}(grid.nodes_test[i]) for i in nodeidx]::Vector{Vec{dim,Float64}}
end

nodes = [(-1.0,-1.0); (0.0,-1.0); (1.0,-1.0); (-1.0,0.0); (0.0,0.0); (1.0,0.0); (-1.0,1.0); (0.0,1.0); (1.0,1.0)]
cells = (Quadrilateral((1,2,5,4)), Quadrilateral((2,3,6,5)), Quadrilateral((4,5,8,7)), Quadrilateral((5,6,9,8)))
subtype_grid = SmallGrid(nodes,cells)
reference_grid = generate_grid(Quadrilateral, (2,2))

ip = Lagrange{2, RefCube, 1}()
qr = QuadratureRule{2, RefCube}(2)
cellvalues = CellScalarValues(qr, ip);

dhs = [DofHandler(grid) for grid in (subtype_grid, reference_grid)]
u1 = Vector{Float64}(undef, 9)
u2 = Vector{Float64}(undef, 9)
∂Ω = union(getfaceset.((reference_grid, ), ["left", "right", "top", "bottom"])...)
dbc = Dirichlet(:u, ∂Ω, (x, t) -> 0)

function doassemble!(cellvalues::CellScalarValues{dim}, K::SparseMatrixCSC, dh::DofHandler) where {dim}
n_basefuncs = getnbasefunctions(cellvalues)
Ke = zeros(n_basefuncs, n_basefuncs)
fe = zeros(n_basefuncs)
f = zeros(ndofs(dh))
assembler = start_assemble(K, f)
for cellid in 1:getncells(dh.grid)
fill!(Ke, 0)
fill!(fe, 0)
coords = getcoordinates(dh.grid, cellid)
reinit!(cellvalues, coords)
for q_point in 1:getnquadpoints(cellvalues)
= getdetJdV(cellvalues, q_point)
for i in 1:n_basefuncs
v = shape_value(cellvalues, q_point, i)
∇v = shape_gradient(cellvalues, q_point, i)
fe[i] += v *
for j in 1:n_basefuncs
∇u = shape_gradient(cellvalues, q_point, j)
Ke[i, j] += (∇v ∇u) *
end
end
end
assemble!(assembler, celldofs(dh,cellid), fe, Ke)
end
return K, f
end

for (dh,u) in zip(dhs,(u1,u2))
push!(dh, :u, 1)
close!(dh)
ch = ConstraintHandler(dh)
add!(ch, dbc)
close!(ch)
update!(ch, 0.0)
K = create_sparsity_pattern(dh);
K, f = doassemble!(cellvalues, K, dh);
apply!(K, f, ch)
sol = K \ f
u .= sol
end

@test Ferrite.ndofs_per_cell(dhs[1]) == Ferrite.ndofs_per_cell(dhs[2])
@test Ferrite.celldofs(dhs[1],3) == Ferrite.celldofs(dhs[2],3)
@test Ferrite.ndofs(dhs[1]) == Ferrite.ndofs(dhs[2])
@test isapprox(u1,u2,atol=1e-8)
end

0 comments on commit b46853e

Please sign in to comment.