Skip to content

Commit

Permalink
Rearrange internal module structure and be more specific with some na…
Browse files Browse the repository at this point in the history
…mes.
  • Loading branch information
termi-official committed Nov 11, 2023
1 parent 16c0cd2 commit 1391b15
Show file tree
Hide file tree
Showing 17 changed files with 324 additions and 293 deletions.
32 changes: 9 additions & 23 deletions src/Thunderbolt.jl
Original file line number Diff line number Diff line change
Expand Up @@ -15,52 +15,38 @@ import Ferrite: vertices, edges, faces, sortedge, sortface

import Krylov: CgSolver

include("collections.jl")
include("utils.jl")

include("mesh/meshes.jl")
include("mesh/coordinate_systems.jl")
include("mesh/tools.jl")
include("mesh/generators.jl")


# TODO where to put these?
abstract type SteadyStateInternalVariable end

include("modeling/coefficients.jl")

include("modeling/boundary_conditions.jl")
# Note that some modules below have an "interface.jl" but this one has only a "common.jl".
# This is simply because there is no modeling interface, but just individual physics modules and couplers.
include("modeling/common.jl")

include("modeling/microstructure.jl")

include("modeling/electrophysiology.jl")
include("modeling/solid_mechanics.jl")
include("modeling/fluid_mechanics.jl")

include("modeling/solid/energies.jl")
include("modeling/solid/contraction.jl")
include("modeling/solid/active.jl")
include("modeling/solid/drivers.jl") # TODO better name. This is basically the quadrature point routine.
include("modeling/multiphysics.jl")

include("modeling/fluid/lumped.jl")

include("modeling/coupler/interface.jl")
include("modeling/coupler/fsi.jl")

include("modeling/problems.jl")
include("modeling/problems.jl") # This is not really "modeling" but a glue layer to translate from model to solver via a discretization

include("solver/interface.jl")
include("solver/operator.jl")
include("solver/newton_raphson.jl")
include("solver/load_stepping.jl")
include("solver/backward_euler.jl")
include("solver/euler.jl")
include("solver/partitioned_solver.jl")
include("solver/operator_splitting.jl")

include("discretization/interface.jl")
include("discretization/fem.jl")


include("io.jl")

# TODO put exports into the individual submodules above!
export
solve,
# Coefficients
Expand Down
140 changes: 5 additions & 135 deletions src/mesh/meshes.jl
Original file line number Diff line number Diff line change
Expand Up @@ -11,143 +11,13 @@
# defining_nodes::NTuple{3,Int}
# end

# TODO we might want to add this to Ferrite (and especially FerriteViz) in one or another way. Maybe traits are better, because they allow more extensibility.
const LinearCellGeometry = Union{Hexahedron, Tetrahedron, Pyramid, Wedge, Triangle, Quadrilateral, Line}

elementtypes(grid::Grid{3,Hexahedron}) = @SVector [Hexahedron]
elementtypes(grid::Grid{3,Tetrahedron}) = @SVector [Tetrahedron]

"""
SimpleMesh3D{C <: AbstractCell, T <: Real} <: AbstractGrid{3}
A grid which also has information abouts its vertices, faces and edges.
"""
struct SimpleMesh3D{C <: AbstractCell, T <: Real} <: AbstractGrid{3}
grid::Grid{3, C, T}
mfaces::OrderedDict{NTuple{3,Int}, Int} # Maps "sortface"-representation to id
medges::OrderedDict{NTuple{2,Int}, Int} # Maps "sortedge"-representation to id
mvertices::OrderedDict{Int, Int} # Maps node to id
number_of_cells_by_type::Dict{DataType, Int}
end

global_edges(mgrid::SimpleMesh3D, cell) = [mgrid.medges[sedge] for sedge first.(sortedge.(edges(cell)))]
global_faces(mgrid::SimpleMesh3D, cell) = [mgrid.mfaces[sface] for sface first.(sortface.(faces(cell)))]
global_vertices(mgrid::SimpleMesh3D, cell) = [mgrid.mvertices[v] for v vertices(cell)]

num_nodes(mgrid::SimpleMesh3D) = length(mgrid.grid.nodes)
num_faces(mgrid::SimpleMesh3D) = length(mgrid.mfaces)
num_edges(mgrid::SimpleMesh3D) = length(mgrid.medges)
num_vertices(mgrid::SimpleMesh3D) = length(mgrid.mvertices)

elementtypes(::SimpleMesh3D{Tetrahedron}) = @SVector [Tetrahedron]
elementtypes(::SimpleMesh3D{Hexahedron}) = @SVector [Hexahedron]

function to_mesh(grid::Grid{3,C,T}) where {C, T}
mfaces = OrderedDict{NTuple{3,Int}, Int}()
medges = OrderedDict{NTuple{2,Int}, Int}()
mvertices = OrderedDict{Int, Int}()
next_face_idx = 1
next_edge_idx = 1
next_vertex_idx = 1
number_of_cells_by_type = Dict{DataType, Int}()
for cell getcells(grid)
cell_type = typeof(cell)
if haskey(number_of_cells_by_type, cell_type)
number_of_cells_by_type[cell_type] += 1
else
number_of_cells_by_type[cell_type] = 1
end

for v vertices(cell)
if !haskey(mvertices, v)
mvertices[v] = next_vertex_idx
next_vertex_idx += 1
end
end
for e first.(sortedge.(edges(cell)))
if !haskey(medges, e)
medges[e] = next_edge_idx
next_edge_idx += 1
end
end
for f first.(sortface.(faces(cell)))
if !haskey(mfaces, f)
mfaces[f] = next_face_idx
next_face_idx += 1
end
end
end
return SimpleMesh3D(grid, mfaces, medges, mvertices, number_of_cells_by_type)
end


"""
SimpleMesh2D{C <: AbstractCell, T <: Real} <: AbstractGrid{2}
A grid which also has information abouts its vertices, faces and edges.
"""
struct SimpleMesh2D{C <: AbstractCell, T <: Real} <: AbstractGrid{2}
grid::Grid{2, C, T}
mfaces::OrderedDict{NTuple{2,Int}, Int} # Maps "sortface"-representation to id
mvertices::OrderedDict{Int, Int} # Maps node to id
number_of_cells_by_type::Dict{DataType, Int}
end

function to_mesh(grid::Grid{2,C,T}) where {C, T}
mfaces = OrderedDict{NTuple{2,Int}, Int}()
mvertices = OrderedDict{Int, Int}()
next_face_idx = 1
next_vertex_idx = 1
number_of_cells_by_type = Dict{DataType, Int}()
for cell getcells(grid)
cell_type = typeof(cell)
if haskey(number_of_cells_by_type, cell_type)
number_of_cells_by_type[cell_type] += 1
else
number_of_cells_by_type[cell_type] = 1
end

for v vertices(cell)
if !haskey(mvertices, v)
mvertices[v] = next_vertex_idx
next_vertex_idx += 1
end
end
for f first.(sortface.(faces(cell)))
if !haskey(mfaces, f)
mfaces[f] = next_face_idx
next_face_idx += 1
end
end
end
return SimpleMesh2D(grid, mfaces, mvertices, number_of_cells_by_type)
end


global_faces(mgrid::SimpleMesh2D, cell) = [mgrid.mfaces[sface] for sface first.(sortface.(faces(cell)))]
global_vertices(mgrid::SimpleMesh2D, cell) = [mgrid.mvertices[v] for v vertices(cell)]

num_nodes(mgrid::SimpleMesh2D) = length(mgrid.grid.nodes)
num_faces(mgrid::SimpleMesh2D) = length(mgrid.mfaces)
num_vertices(mgrid::SimpleMesh2D) = length(mgrid.mvertices)

elementtypes(::SimpleMesh2D{Triangle}) = @SVector [Triangle]
elementtypes(::SimpleMesh2D{Quadrilateral}) = @SVector [Quadrilateral]

# Ferrite compat layer for the mesh
@inline Ferrite.getncells(mesh::Union{SimpleMesh2D,SimpleMesh3D}) = Ferrite.getncells(mesh.grid)
@inline Ferrite.getcells(mesh::Union{SimpleMesh2D,SimpleMesh3D}) = Ferrite.getcells(mesh.grid)
@inline Ferrite.getcells(mesh::Union{SimpleMesh2D,SimpleMesh3D}, v::Union{Int, Vector{Int}}) = Ferrite.getcells(mesh.grid, v)
@inline Ferrite.getcells(mesh::Union{SimpleMesh2D,SimpleMesh3D}, setname::String) = Ferrite.getcells(mesh.grid, setname)
@inline Ferrite.getcelltype(mesh::Union{SimpleMesh2D,SimpleMesh3D}) = Ferrite.getcelltype(mesh.grid)
@inline Ferrite.getcelltype(mesh::Union{SimpleMesh2D,SimpleMesh3D}, i::Int) = Ferrite.getcelltype(mesh.grid)

@inline Ferrite.getnnodes(mesh::Union{SimpleMesh2D,SimpleMesh3D}) = Ferrite.getnnodes(mesh.grid)
@inline Ferrite.nnodes_per_cell(mesh::Union{SimpleMesh2D,SimpleMesh3D}, i::Int=1) = Ferrite.nnodes_per_cell(mesh.grid, i)
@inline Ferrite.getnodes(mesh::Union{SimpleMesh2D,SimpleMesh3D}) = Ferrite.getnodes(mesh.grid)
@inline Ferrite.getnodes(mesh::Union{SimpleMesh2D,SimpleMesh3D}, v::Union{Int, Vector{Int}}) = Ferrite.getnodes(mesh.grid, v)
@inline Ferrite.getnodes(mesh::Union{SimpleMesh2D,SimpleMesh3D}, setname::String) = Ferrite.getnodes(mesh.grid, setname)

@inline Ferrite.vtk_grid(filename::AbstractString, mesh::Union{SimpleMesh2D,SimpleMesh3D}; kwargs...) = Ferrite.vtk_grid(filename, mesh.grid, kwargs...)

@inline Ferrite.get_coordinate_type(::SimpleMesh2D{C,T}) where {C,T} = Vec{2,T}
@inline Ferrite.get_coordinate_type(::SimpleMesh3D{C,T}) where {C,T} = Vec{3,T}
include("simple_meshes.jl")
include("coordinate_systems.jl")
include("tools.jl")
include("generators.jl")
136 changes: 136 additions & 0 deletions src/mesh/simple_meshes.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,136 @@

"""
SimpleMesh3D{C <: AbstractCell, T <: Real} <: AbstractGrid{3}
A grid which also has information abouts its vertices, faces and edges.
"""
struct SimpleMesh3D{C <: AbstractCell, T <: Real} <: AbstractGrid{3}
grid::Grid{3, C, T}
mfaces::OrderedDict{NTuple{3,Int}, Int} # Maps "sortface"-representation to id
medges::OrderedDict{NTuple{2,Int}, Int} # Maps "sortedge"-representation to id
mvertices::OrderedDict{Int, Int} # Maps node to id
number_of_cells_by_type::Dict{DataType, Int}
end

global_edges(mgrid::SimpleMesh3D, cell) = [mgrid.medges[sedge] for sedge first.(sortedge.(edges(cell)))]
global_faces(mgrid::SimpleMesh3D, cell) = [mgrid.mfaces[sface] for sface first.(sortface.(faces(cell)))]
global_vertices(mgrid::SimpleMesh3D, cell) = [mgrid.mvertices[v] for v vertices(cell)]

num_nodes(mgrid::SimpleMesh3D) = length(mgrid.grid.nodes)
num_faces(mgrid::SimpleMesh3D) = length(mgrid.mfaces)
num_edges(mgrid::SimpleMesh3D) = length(mgrid.medges)
num_vertices(mgrid::SimpleMesh3D) = length(mgrid.mvertices)

elementtypes(::SimpleMesh3D{Tetrahedron}) = @SVector [Tetrahedron]
elementtypes(::SimpleMesh3D{Hexahedron}) = @SVector [Hexahedron]

function to_mesh(grid::Grid{3,C,T}) where {C, T}
mfaces = OrderedDict{NTuple{3,Int}, Int}()
medges = OrderedDict{NTuple{2,Int}, Int}()
mvertices = OrderedDict{Int, Int}()
next_face_idx = 1
next_edge_idx = 1
next_vertex_idx = 1
number_of_cells_by_type = Dict{DataType, Int}()
for cell getcells(grid)
cell_type = typeof(cell)
if haskey(number_of_cells_by_type, cell_type)
number_of_cells_by_type[cell_type] += 1
else
number_of_cells_by_type[cell_type] = 1
end

for v vertices(cell)
if !haskey(mvertices, v)
mvertices[v] = next_vertex_idx
next_vertex_idx += 1
end
end
for e first.(sortedge.(edges(cell)))
if !haskey(medges, e)
medges[e] = next_edge_idx
next_edge_idx += 1
end
end
for f first.(sortface.(faces(cell)))
if !haskey(mfaces, f)
mfaces[f] = next_face_idx
next_face_idx += 1
end
end
end
return SimpleMesh3D(grid, mfaces, medges, mvertices, number_of_cells_by_type)
end


"""
SimpleMesh2D{C <: AbstractCell, T <: Real} <: AbstractGrid{2}
A grid which also has information abouts its vertices, faces and edges.
"""
struct SimpleMesh2D{C <: AbstractCell, T <: Real} <: AbstractGrid{2}
grid::Grid{2, C, T}
mfaces::OrderedDict{NTuple{2,Int}, Int} # Maps "sortface"-representation to id
mvertices::OrderedDict{Int, Int} # Maps node to id
number_of_cells_by_type::Dict{DataType, Int}
end

function to_mesh(grid::Grid{2,C,T}) where {C, T}
mfaces = OrderedDict{NTuple{2,Int}, Int}()
mvertices = OrderedDict{Int, Int}()
next_face_idx = 1
next_vertex_idx = 1
number_of_cells_by_type = Dict{DataType, Int}()
for cell getcells(grid)
cell_type = typeof(cell)
if haskey(number_of_cells_by_type, cell_type)
number_of_cells_by_type[cell_type] += 1
else
number_of_cells_by_type[cell_type] = 1
end

for v vertices(cell)
if !haskey(mvertices, v)
mvertices[v] = next_vertex_idx
next_vertex_idx += 1
end
end
for f first.(sortface.(faces(cell)))
if !haskey(mfaces, f)
mfaces[f] = next_face_idx
next_face_idx += 1
end
end
end
return SimpleMesh2D(grid, mfaces, mvertices, number_of_cells_by_type)
end


global_faces(mgrid::SimpleMesh2D, cell) = [mgrid.mfaces[sface] for sface first.(sortface.(faces(cell)))]
global_vertices(mgrid::SimpleMesh2D, cell) = [mgrid.mvertices[v] for v vertices(cell)]

num_nodes(mgrid::SimpleMesh2D) = length(mgrid.grid.nodes)
num_faces(mgrid::SimpleMesh2D) = length(mgrid.mfaces)
num_vertices(mgrid::SimpleMesh2D) = length(mgrid.mvertices)

elementtypes(::SimpleMesh2D{Triangle}) = @SVector [Triangle]
elementtypes(::SimpleMesh2D{Quadrilateral}) = @SVector [Quadrilateral]

# Ferrite compat layer for the mesh
@inline Ferrite.getncells(mesh::Union{SimpleMesh2D,SimpleMesh3D}) = Ferrite.getncells(mesh.grid)
@inline Ferrite.getcells(mesh::Union{SimpleMesh2D,SimpleMesh3D}) = Ferrite.getcells(mesh.grid)
@inline Ferrite.getcells(mesh::Union{SimpleMesh2D,SimpleMesh3D}, v::Union{Int, Vector{Int}}) = Ferrite.getcells(mesh.grid, v)
@inline Ferrite.getcells(mesh::Union{SimpleMesh2D,SimpleMesh3D}, setname::String) = Ferrite.getcells(mesh.grid, setname)
@inline Ferrite.getcelltype(mesh::Union{SimpleMesh2D,SimpleMesh3D}) = Ferrite.getcelltype(mesh.grid)
@inline Ferrite.getcelltype(mesh::Union{SimpleMesh2D,SimpleMesh3D}, i::Int) = Ferrite.getcelltype(mesh.grid)

@inline Ferrite.getnnodes(mesh::Union{SimpleMesh2D,SimpleMesh3D}) = Ferrite.getnnodes(mesh.grid)
@inline Ferrite.nnodes_per_cell(mesh::Union{SimpleMesh2D,SimpleMesh3D}, i::Int=1) = Ferrite.nnodes_per_cell(mesh.grid, i)
@inline Ferrite.getnodes(mesh::Union{SimpleMesh2D,SimpleMesh3D}) = Ferrite.getnodes(mesh.grid)
@inline Ferrite.getnodes(mesh::Union{SimpleMesh2D,SimpleMesh3D}, v::Union{Int, Vector{Int}}) = Ferrite.getnodes(mesh.grid, v)
@inline Ferrite.getnodes(mesh::Union{SimpleMesh2D,SimpleMesh3D}, setname::String) = Ferrite.getnodes(mesh.grid, setname)

@inline Ferrite.vtk_grid(filename::AbstractString, mesh::Union{SimpleMesh2D,SimpleMesh3D}; kwargs...) = Ferrite.vtk_grid(filename, mesh.grid, kwargs...)

@inline Ferrite.get_coordinate_type(::SimpleMesh2D{C,T}) where {C,T} = Vec{2,T}
@inline Ferrite.get_coordinate_type(::SimpleMesh3D{C,T}) where {C,T} = Vec{3,T}
11 changes: 11 additions & 0 deletions src/modeling/common.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,11 @@
# Common modeling primitives are found here

# TODO actually use this :)
abstract type SteadyStateInternalVariable end

include("core/coefficients.jl")

include("core/boundary_conditions.jl")

include("core/mass.jl")
include("core/diffusion.jl")
File renamed without changes.
File renamed without changes.
38 changes: 38 additions & 0 deletions src/modeling/core/diffusion.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,38 @@
# @doc raw"""
# AssembledDiffusionOperator{MT, DT, CV}

# Assembles the matrix associated to the bilinearform ``a(u,v) = -\int \nabla v(x) \cdot D(x) \nabla u(x) dx`` for a given diffusion tensor ``D(x)`` and ``u,v`` from the same function space.
# """
"""
Represents the integrand of the bilinear form <ϕ,ψ> = -∫ D∇ϕ ⋅ ∇ψ dΩ .
"""
struct BilinearDiffusionIntegrator{CoefficientType}
D::CoefficientType
# coordinate_system
end

struct BilinearDiffusionElementCache{IT <: BilinearDiffusionIntegrator, CV}
integrator::IT
cellvalues::CV
end

function assemble_element!(Kₑ, cell, element_cache::CACHE, time) where {CACHE <: BilinearDiffusionElementCache}
@unpack cellvalues = element_cache
n_basefuncs = getnbasefunctions(cellvalues)

reinit!(cellvalues, cell)

for q_point in 1:getnquadpoints(cellvalues)
ξ = cellvalues.qr.points[q_point]
qp = QuadraturePoint(q_point, ξ)
D_loc = evaluate_coefficient(element_cache.integrator.D, cell, qp, time)
= getdetJdV(cellvalues, q_point)
for i in 1:n_basefuncs
∇Nᵢ = shape_gradient(cellvalues, q_point, i)
for j in 1:n_basefuncs
∇Nⱼ = shape_gradient(cellvalues, q_point, j)
Kₑ[i,j] -= ((D_loc ∇Nᵢ) ∇Nⱼ) *
end
end
end
end
Loading

0 comments on commit 1391b15

Please sign in to comment.