Book a Demo!
CoCalc Logo Icon
StoreFeaturesDocsShareSupportNewsAboutPoliciesSign UpSign In
trixi-framework
GitHub Repository: trixi-framework/Trixi.jl
 Path: blob/main/examples/dgmulti_2d/elixir_advection_diffusion.jl
2055 views
1
using OrdinaryDiffEqLowStorageRK

2
using Trixi

3


4
dg = DGMulti(polydeg = 3, element_type = Tri(), approximation_type = Polynomial(),

5
             surface_integral = SurfaceIntegralWeakForm(flux_lax_friedrichs),

6
             volume_integral = VolumeIntegralWeakForm())

7


8
equations = LinearScalarAdvectionEquation2D(1.5, 1.0)

9
equations_parabolic = LaplaceDiffusion2D(5.0e-2, equations)

10


11
initial_condition_zero(x, t, equations::LinearScalarAdvectionEquation2D) = SVector(0.0)

12
initial_condition = initial_condition_zero

13


14
# tag different boundary segments

15
left(x, tol = 50 * eps()) = abs(x[1] + 1) < tol

16
right(x, tol = 50 * eps()) = abs(x[1] - 1) < tol

17
bottom(x, tol = 50 * eps()) = abs(x[2] + 1) < tol

18
top(x, tol = 50 * eps()) = abs(x[2] - 1) < tol

19
is_on_boundary = Dict(:left => left, :right => right, :top => top, :bottom => bottom)

20


21
cells_per_dimension = (16, 16)

22
mesh = DGMultiMesh(dg, cells_per_dimension; is_on_boundary)

23


24
# BC types

25
boundary_condition_left = BoundaryConditionDirichlet((x, t, equations_parabolic) -> SVector(1 +

26
                                                                                            0.1 *

27
                                                                                            x[2]))

28
boundary_condition_zero = BoundaryConditionDirichlet((x, t, equations_parabolic) -> SVector(0.0))

29
boundary_condition_neumann_zero = BoundaryConditionNeumann((x, t, equations_parabolic) -> SVector(0.0))

30


31
# define inviscid boundary conditions

32
boundary_conditions = (; :left => boundary_condition_left,

33
                       :bottom => boundary_condition_zero,

34
                       :top => boundary_condition_do_nothing,

35
                       :right => boundary_condition_do_nothing)

36


37
# define viscous boundary conditions

38
boundary_conditions_parabolic = (; :left => boundary_condition_left,

39
                                 :bottom => boundary_condition_zero,

40
                                 :top => boundary_condition_zero,

41
                                 :right => boundary_condition_neumann_zero)

42


43
semi = SemidiscretizationHyperbolicParabolic(mesh, (equations, equations_parabolic),

44
                                             initial_condition, dg;

45
                                             boundary_conditions = (boundary_conditions,

46
                                                                    boundary_conditions_parabolic))

47


48
tspan = (0.0, 1.5)

49
ode = semidiscretize(semi, tspan)

50


51
summary_callback = SummaryCallback()

52
alive_callback = AliveCallback(alive_interval = 10)

53
analysis_interval = 100

54
analysis_callback = AnalysisCallback(semi, interval = analysis_interval, uEltype = real(dg))

55
save_solution = SaveSolutionCallback(interval = analysis_interval,

56
                                     solution_variables = cons2prim)

57
callbacks = CallbackSet(summary_callback, alive_callback, analysis_callback, save_solution)

58


59
###############################################################################

60
# run the simulation

61


62
time_int_tol = 1e-6

63
sol = solve(ode, RDPK3SpFSAL49(); abstol = time_int_tol, reltol = time_int_tol,

64
            ode_default_options()..., callback = callbacks)

65


66