diff --git a/src/dg1d.jl b/src/dg1d.jl index a6f22d9c8f155ee524207f02c064703372484b80..5b0e5d1378bba7e6bcc0d288e9a7d075c806889f 100644 --- a/src/dg1d.jl +++ b/src/dg1d.jl @@ -43,6 +43,7 @@ module dg1d include("spectralelement.jl") include("subvolumeelement.jl") include("mesh.jl") + include("hybridmesh.jl") # custom evolution include("evolve.jl") diff --git a/src/dg_rhs.jl b/src/dg_rhs.jl index c980efa1f9f5120475c966dffe64a7afff014801..8ab17c9d20a0d9d23fdbac9965d5ba2539b03bf8 100644 --- a/src/dg_rhs.jl +++ b/src/dg_rhs.jl @@ -1,5 +1,6 @@ -export compute_rhs_weak_form!, - compute_rhs_strong_form! +export compute_rhs_weak_form!, + compute_rhs_strong_form!, + compute_rhs_fv_form! """ @@ -55,3 +56,36 @@ function compute_rhs_strong_form!(rhs, f, s, nf_lhs, nf_rhs, mesh::AbstractMesh) rhs .+= s return end + + +""" + compute_rhs_weak_form!(rhs, sv_flag, f, mesh::HybridMesh) + compute_rhs_weak_form!(rhs, sv_flag, f, s, mesh::HybridMesh) +""" +function compute_rhs_fv_form!(rhs, sv_flag, f, hybridmesh::HybridMesh) + @unpack mesh, sv_element = hybridmesh + @unpack invjac = mesh + shape = layout(mesh) + mat_rhs = reshape(view(rhs,:), shape) + mat_f = reshape(view(f,:), shape) + Npts, K = shape + for (k#=cell nr=#, flag) in enumerate(sv_flag) + flag == 0 && continue + for i = 1:Npts-1 + mat_rhs[i,k] = - (mat_f[i+1,k] - mat_f[i,k]) / sv_element.w[i] * invjac[k] + end + end + return +end +function compute_rhs_fv_form!(rhs, sv_flag, f, s, hybridmesh::HybridMesh) + compute_rhs_weak_form!(rhs, tci_flag, f, hybridmesh) + shape = layout(mesh) + mat_rhs = reshape(view(rhs,:), shape) + mat_s = reshape(view(s, :), shape) + for (k#=cell nr=#, flag) in enumerate(sv_flag) + flag == 0 && continue + mat_rhs[:,k] .+= mat_s[:,k] + end + rhs .+= s + return +end diff --git a/src/hybridmesh.jl b/src/hybridmesh.jl new file mode 100644 index 0000000000000000000000000000000000000000..f5545720979360de38f73d73f98a349a2ec5b39b --- /dev/null +++ b/src/hybridmesh.jl @@ -0,0 +1,18 @@ +struct HybridMesh <: AbstractMesh + mesh::Mesh + sv_element::SubVolumeElement + + function HybridMesh(mesh::Mesh) + sv_element = SubVolumeElement(mesh.element) + return new(mesh, sv_element) + end +end + + +function Base.getproperty(hmesh::HybridMesh, p::Symbol) + if p in propertynames(hmesh.mesh) + return getproperty(hmesh.mesh, p) + else + getfield(hmesh.mesh, p) + end +end