From 246f3f587a5cad8432da250d190bb7bbc8fbefb3 Mon Sep 17 00:00:00 2001 From: Jesse Chan Date: Fri, 30 Aug 2024 08:39:23 -0500 Subject: [PATCH 1/2] fix Fmask calculation for SBP Hexes --- src/RefElemData_SBP.jl | 15 +++++++++++++++ test/multidim_sbp_tests.jl | 3 +++ 2 files changed, 18 insertions(+) diff --git a/src/RefElemData_SBP.jl b/src/RefElemData_SBP.jl index 14fba196..a1324bee 100644 --- a/src/RefElemData_SBP.jl +++ b/src/RefElemData_SBP.jl @@ -26,6 +26,21 @@ function RefElemData(elementType::Union{Quad, Hex}, approxType::SBP{TensorProduc Polynomial(TensorProductQuadrature(gauss_lobatto_quad(0, 0, N))), N; kwargs...) + println("Brute force Fmask determination") + # manually determine Fmask so that rd.rf = rd.r[rd.Fmask], ... + new_Fmask = copy(rd.Fmask) + for fid in eachindex(rd.rf) + rstf = SVector(getindex.(rd.rstf, fid)) + for i in eachindex(rd.r) + rst = SVector(getindex.(rd.rst, i)) + if rstf ≈ rst + new_Fmask[fid] = i + break + end + end + end + rd = @set rd.Fmask = new_Fmask; + rd = @set rd.Vf = droptol!(sparse(rd.Vf), tol) rd = @set rd.LIFT = Diagonal(rd.wq) \ (rd.Vf' * Diagonal(rd.wf)) # TODO: make this more efficient with LinearMaps? diff --git a/test/multidim_sbp_tests.jl b/test/multidim_sbp_tests.jl index 2ff50614..d8b2bdee 100644 --- a/test/multidim_sbp_tests.jl +++ b/test/multidim_sbp_tests.jl @@ -47,6 +47,9 @@ @testset "TensorProductLobatto Hex" begin rd = RefElemData(Hex(),SBP(),N) @test propertynames(rd)[1] == :element_type + @test rd.rf == rd.r[rd.Fmask] + @test rd.sf == rd.s[rd.Fmask] + @test rd.tf == rd.t[rd.Fmask] @test rd.t == rd.rst[3] @test rd.tf == rd.rstf[3] @test rd.tq == rd.rstq[3] From 622e4b1f64579e7b00072f9e4ffd67b1855baf79 Mon Sep 17 00:00:00 2001 From: Jesse Chan Date: Fri, 30 Aug 2024 09:04:34 -0500 Subject: [PATCH 2/2] fix comment and remove print statement --- src/RefElemData_SBP.jl | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/src/RefElemData_SBP.jl b/src/RefElemData_SBP.jl index a1324bee..c3e1a848 100644 --- a/src/RefElemData_SBP.jl +++ b/src/RefElemData_SBP.jl @@ -26,8 +26,7 @@ function RefElemData(elementType::Union{Quad, Hex}, approxType::SBP{TensorProduc Polynomial(TensorProductQuadrature(gauss_lobatto_quad(0, 0, N))), N; kwargs...) - println("Brute force Fmask determination") - # manually determine Fmask so that rd.rf = rd.r[rd.Fmask], ... + # brute-force determine Fmask so that rd.rf = rd.r[rd.Fmask], etc. new_Fmask = copy(rd.Fmask) for fid in eachindex(rd.rf) rstf = SVector(getindex.(rd.rstf, fid))