From 2ce4cad92787bc7f17e99bd27772a2763df8cbe6 Mon Sep 17 00:00:00 2001 From: Hendrik Ranocha Date: Fri, 2 Dec 2022 07:10:40 +0100 Subject: [PATCH] FDSBP constructor --- .../tree_1d_fdsbp/elixir_burgers_basic.jl | 6 ++--- .../elixir_burgers_linear_stability.jl | 6 ++--- .../tree_1d_fdsbp/elixir_euler_convergence.jl | 6 ++--- .../elixir_euler_density_wave.jl | 6 ++--- .../elixir_advection_extended.jl | 6 ++--- .../tree_2d_fdsbp/elixir_euler_convergence.jl | 6 ++--- ...ixir_euler_kelvin_helmholtz_instability.jl | 6 ++--- examples/tree_2d_fdsbp/elixir_euler_vortex.jl | 6 ++--- .../tree_3d_fdsbp/elixir_euler_convergence.jl | 6 ++--- .../elixir_euler_taylor_green_vortex.jl | 6 ++--- src/Trixi.jl | 1 + src/solvers/fdsbp_tree/fdsbp.jl | 22 ++++++++++++++++++- 12 files changed, 52 insertions(+), 31 deletions(-) diff --git a/examples/tree_1d_fdsbp/elixir_burgers_basic.jl b/examples/tree_1d_fdsbp/elixir_burgers_basic.jl index b6d958c66b9..e7008a520a4 100644 --- a/examples/tree_1d_fdsbp/elixir_burgers_basic.jl +++ b/examples/tree_1d_fdsbp/elixir_burgers_basic.jl @@ -17,9 +17,9 @@ D_upw = upwind_operators(SummationByPartsOperators.Mattsson2017, xmin=-1.0, xmax=1.0, N=32) flux_splitting = splitting_lax_friedrichs -solver = DG(D_upw, nothing #= mortar =#, - SurfaceIntegralUpwind(flux_splitting), - VolumeIntegralUpwind(flux_splitting)) +solver = FDSBP(D_upw, + surface_integral=SurfaceIntegralUpwind(flux_splitting), + volume_integral=VolumeIntegralUpwind(flux_splitting)) coordinates_min = 0.0 coordinates_max = 1.0 diff --git a/examples/tree_1d_fdsbp/elixir_burgers_linear_stability.jl b/examples/tree_1d_fdsbp/elixir_burgers_linear_stability.jl index 322421e0954..e5596063eba 100644 --- a/examples/tree_1d_fdsbp/elixir_burgers_linear_stability.jl +++ b/examples/tree_1d_fdsbp/elixir_burgers_linear_stability.jl @@ -20,9 +20,9 @@ D_upw = upwind_operators(SummationByPartsOperators.Mattsson2017, xmin=-1.0, xmax=1.0, N=16) flux_splitting = splitting_lax_friedrichs -solver = DG(D_upw, nothing #= mortar =#, - SurfaceIntegralUpwind(flux_splitting), - VolumeIntegralUpwind(flux_splitting)) +solver = FDSBP(D_upw, + surface_integral=SurfaceIntegralUpwind(flux_splitting), + volume_integral=VolumeIntegralUpwind(flux_splitting)) coordinates_min = -1.0 coordinates_max = 1.0 diff --git a/examples/tree_1d_fdsbp/elixir_euler_convergence.jl b/examples/tree_1d_fdsbp/elixir_euler_convergence.jl index cfd233a0a8c..34717b2c52f 100644 --- a/examples/tree_1d_fdsbp/elixir_euler_convergence.jl +++ b/examples/tree_1d_fdsbp/elixir_euler_convergence.jl @@ -17,9 +17,9 @@ D_upw = upwind_operators(SummationByPartsOperators.Mattsson2017, xmin=-1.0, xmax=1.0, N=32) flux_splitting = splitting_steger_warming -solver = DG(D_upw, nothing #= mortar =#, - SurfaceIntegralUpwind(flux_splitting), - VolumeIntegralUpwind(flux_splitting)) +solver = FDSBP(D_upw, + surface_integral=SurfaceIntegralUpwind(flux_splitting), + volume_integral=VolumeIntegralUpwind(flux_splitting)) coordinates_min = 0.0 coordinates_max = 2.0 diff --git a/examples/tree_1d_fdsbp/elixir_euler_density_wave.jl b/examples/tree_1d_fdsbp/elixir_euler_density_wave.jl index a79fedb2778..a11afaf78ff 100644 --- a/examples/tree_1d_fdsbp/elixir_euler_density_wave.jl +++ b/examples/tree_1d_fdsbp/elixir_euler_density_wave.jl @@ -16,9 +16,9 @@ D_upw = upwind_operators(SummationByPartsOperators.Mattsson2017, xmin=-1.0, xmax=1.0, N=16) flux_splitting = splitting_coirier_vanleer -solver = DG(D_upw, nothing #= mortar =#, - SurfaceIntegralUpwind(flux_splitting), - VolumeIntegralUpwind(flux_splitting)) +solver = FDSBP(D_upw, + surface_integral=SurfaceIntegralUpwind(flux_splitting), + volume_integral=VolumeIntegralUpwind(flux_splitting)) coordinates_min = -1.0 coordinates_max = 1.0 diff --git a/examples/tree_2d_fdsbp/elixir_advection_extended.jl b/examples/tree_2d_fdsbp/elixir_advection_extended.jl index ab9d62b4096..d20894c278c 100644 --- a/examples/tree_2d_fdsbp/elixir_advection_extended.jl +++ b/examples/tree_2d_fdsbp/elixir_advection_extended.jl @@ -15,9 +15,9 @@ initial_condition = initial_condition_convergence_test D_SBP = derivative_operator(SummationByPartsOperators.MattssonNordström2004(), derivative_order=1, accuracy_order=4, xmin=0.0, xmax=1.0, N=100) -solver = DG(D_SBP, nothing #= mortar =#, - SurfaceIntegralStrongForm(flux_lax_friedrichs), - VolumeIntegralStrongForm()) +solver = FDSBP(D_SBP, + surface_integral=SurfaceIntegralStrongForm(flux_lax_friedrichs), + volume_integral=VolumeIntegralStrongForm()) coordinates_min = (-1.0, -1.0) coordinates_max = ( 1.0, 1.0) diff --git a/examples/tree_2d_fdsbp/elixir_euler_convergence.jl b/examples/tree_2d_fdsbp/elixir_euler_convergence.jl index dcdfb35af6d..a5efd90c094 100644 --- a/examples/tree_2d_fdsbp/elixir_euler_convergence.jl +++ b/examples/tree_2d_fdsbp/elixir_euler_convergence.jl @@ -17,9 +17,9 @@ D = upwind_operators(SummationByPartsOperators.Mattsson2017, xmin=-1.0, xmax=1.0, N=16) flux_splitting = splitting_steger_warming -solver = DG(D, nothing #= mortar =#, - SurfaceIntegralUpwind(flux_splitting), - VolumeIntegralUpwind(flux_splitting)) +solver = FDSBP(D_upw, + surface_integral=SurfaceIntegralUpwind(flux_splitting), + volume_integral=VolumeIntegralUpwind(flux_splitting)) coordinates_min = (-1.0, -1.0) coordinates_max = ( 1.0, 1.0) diff --git a/examples/tree_2d_fdsbp/elixir_euler_kelvin_helmholtz_instability.jl b/examples/tree_2d_fdsbp/elixir_euler_kelvin_helmholtz_instability.jl index d11c89875a2..7e86f36fa40 100644 --- a/examples/tree_2d_fdsbp/elixir_euler_kelvin_helmholtz_instability.jl +++ b/examples/tree_2d_fdsbp/elixir_euler_kelvin_helmholtz_instability.jl @@ -30,9 +30,9 @@ D_upw = upwind_operators(SummationByPartsOperators.Mattsson2017, xmin=-1.0, xmax=1.0, N=16) flux_splitting = splitting_vanleer_haenel -solver = DG(D_upw, nothing #= mortar =#, - SurfaceIntegralUpwind(flux_splitting), - VolumeIntegralUpwind(flux_splitting)) +solver = FDSBP(D_upw, + surface_integral=SurfaceIntegralUpwind(flux_splitting), + volume_integral=VolumeIntegralUpwind(flux_splitting)) coordinates_min = (-1.0, -1.0) coordinates_max = ( 1.0, 1.0) diff --git a/examples/tree_2d_fdsbp/elixir_euler_vortex.jl b/examples/tree_2d_fdsbp/elixir_euler_vortex.jl index 9c01d64718c..10f6106f4f7 100644 --- a/examples/tree_2d_fdsbp/elixir_euler_vortex.jl +++ b/examples/tree_2d_fdsbp/elixir_euler_vortex.jl @@ -60,9 +60,9 @@ D_upw = upwind_operators(SummationByPartsOperators.Mattsson2017, xmin=-1.0, xmax=1.0, N=16) flux_splitting = splitting_steger_warming -solver = DG(D_upw, nothing #= mortar =#, - SurfaceIntegralUpwind(flux_splitting), - VolumeIntegralUpwind(flux_splitting)) +solver = FDSBP(D_upw, + surface_integral=SurfaceIntegralUpwind(flux_splitting), + volume_integral=VolumeIntegralUpwind(flux_splitting)) coordinates_min = (-10.0, -10.0) coordinates_max = ( 10.0, 10.0) diff --git a/examples/tree_3d_fdsbp/elixir_euler_convergence.jl b/examples/tree_3d_fdsbp/elixir_euler_convergence.jl index 789929e7b5d..576a07e6aba 100644 --- a/examples/tree_3d_fdsbp/elixir_euler_convergence.jl +++ b/examples/tree_3d_fdsbp/elixir_euler_convergence.jl @@ -17,9 +17,9 @@ D_upw = upwind_operators(SummationByPartsOperators.Mattsson2017, xmin=-1.0, xmax=1.0, N=16) flux_splitting = splitting_steger_warming -solver = DG(D_upw, nothing #= mortar =#, - SurfaceIntegralUpwind(flux_splitting), - VolumeIntegralUpwind(flux_splitting)) +solver = FDSBP(D_upw, + surface_integral=SurfaceIntegralUpwind(flux_splitting), + volume_integral=VolumeIntegralUpwind(flux_splitting)) coordinates_min = (0.0, 0.0, 0.0) coordinates_max = (2.0, 2.0, 2.0) diff --git a/examples/tree_3d_fdsbp/elixir_euler_taylor_green_vortex.jl b/examples/tree_3d_fdsbp/elixir_euler_taylor_green_vortex.jl index 5a8fd8c1520..a138cf94e39 100644 --- a/examples/tree_3d_fdsbp/elixir_euler_taylor_green_vortex.jl +++ b/examples/tree_3d_fdsbp/elixir_euler_taylor_green_vortex.jl @@ -35,9 +35,9 @@ D_upw = upwind_operators(SummationByPartsOperators.Mattsson2017, xmin=-1.0, xmax=1.0, N=16) flux_splitting = splitting_steger_warming -solver = DG(D_upw, nothing #= mortar =#, - SurfaceIntegralUpwind(flux_splitting), - VolumeIntegralUpwind(flux_splitting)) +solver = FDSBP(D_upw, + surface_integral=SurfaceIntegralUpwind(flux_splitting), + volume_integral=VolumeIntegralUpwind(flux_splitting)) coordinates_min = (-1.0, -1.0, -1.0) .* pi coordinates_max = ( 1.0, 1.0, 1.0) .* pi diff --git a/src/Trixi.jl b/src/Trixi.jl index 0de19e38f4b..d3940460063 100644 --- a/src/Trixi.jl +++ b/src/Trixi.jl @@ -189,6 +189,7 @@ export TreeMesh, StructuredMesh, UnstructuredMesh2D, P4estMesh export DG, DGSEM, LobattoLegendreBasis, + FDSBP, VolumeIntegralWeakForm, VolumeIntegralStrongForm, VolumeIntegralFluxDifferencing, VolumeIntegralPureLGLFiniteVolume, diff --git a/src/solvers/fdsbp_tree/fdsbp.jl b/src/solvers/fdsbp_tree/fdsbp.jl index 40cbb6ffe1b..b89d59c1156 100644 --- a/src/solvers/fdsbp_tree/fdsbp.jl +++ b/src/solvers/fdsbp_tree/fdsbp.jl @@ -8,9 +8,29 @@ @muladd begin -# For dispatch +""" + FDSBP(D_SBP; surface_integral, volume_integral) + +Specialization of [`DG`](@ref) methods that uses general summation by parts (SBP) +operators from +[SummationByPartsOperators.jl](https://github.com/ranocha/SummationByPartsOperators.jl). +In particular, this includes classical finite difference (FD) SBP methods. +These methods have the same structure as classical DG methods - local operations +on elements with connectivity through interfaces without imposing any continuity +constraints. + +`D_SBP` is an SBP derivative operator from SummationByPartsOperators.jl. +The other arguments have the same meaning as in [`DG`](@ref) or [`DGSEM`](@ref). + +!!! warning "Experimental implementation (upwind SBP)" + This is an experimental feature and may change in future releases. +""" const FDSBP = DG{Basis} where {Basis<:AbstractDerivativeOperator} +function FDSBP(D_SBP::AbstractDerivativeOperator; surface_integral, volume_integral) + return DG(D_SBP, nothing #= mortar =#, surface_integral, volume_integral) +end + # General interface methods for SummationByPartsOperators.jl and Trixi.jl nnodes(D::AbstractDerivativeOperator) = size(D, 1)