Using DGMulti FDSBP with parabolic terms in Trixi?

Specific Domains Modelling & Simulations
,
1

Hi, I am comparing the methods in the DGMulti Solver from Trixi with a visquous example but when I try to use the Finite difference methods periodic_derivative_operator(derivative_order=1, accuracy_order=4, xmin=0.0, xmax=1.0, N=50) and give it to the semidiscritisation i get this error at the beginning of the simulation :

DimensionMismatch:
Stacktrace:
[1] _spmatmul!(C::Vector{Float64}, A::SparseArrays.SparseMatrixCSC{Float64, Int64}, B::Vector{Float64}, α::Bool, β::Bool)
@ SparseArrays ~/.julia/juliaup/julia-1.10.4+0.x64.apple.darwin14/share/julia/stdlib/v1.10/SparseArrays/src/linalg.jl:59

Since no example exists for a visquous flow using FD methods in DGMulti, Is there something I am doing wrong ?

2

Here is my code

methodFD = periodic_derivative_operator(derivative_order=1, accuracy_order=4, xmin=0.0, xmax=1.0, N=40)

solver = DGMulti(element_type = Quad(),
             approximation_type = methodFD,
             surface_flux = flux_hll,
             volume_integral = VolumeIntegralWeakForm())

prandtl_number() = 0.72
mu = 0.00005

equations = CompressibleEulerEquations2D(1.4)
equations_parabolic = CompressibleNavierStokesDiffusion2D(equations, mu = mu, Prandtl = prandtl_number(),gradient_variables = GradientVariablesPrimitive())

function initial_condition_BM_vortex(x, t, equations)

  L = 1.0
  x0 = 0.0
  Ma0 = 0.2        # Maximum Mach number
  rho0 = 1.2         # density outside the shear layer
  
  n_vortex=1
  l_ratio = 100.0      # length scale ratio (30.0 default)
  d_ratio = 1.0       # density ratio (4.0 default)
  v_ratio = 0.05       # velocity ratio (0.2 default)
  shear = 200.0        # shear value (15.0 default)

  y0 = L/2
  ya = y0 - L/4
  yb = y0 + L/4
  sigma = L / l_ratio
  u0 = shear*sigma/2
  v0 = v_ratio * u0
  rho1 = d_ratio * rho0

  B = tanh((x[2]-ya)/sigma)/2 - tanh((x[2]-yb)/sigma)/2

  rho = rho0 * (1.0 + B*(d_ratio-1.0))
  u = u0 * (2*B - 1.0)
  v = v0 * cos(2.0 * pi * n_vortex * (x[1]-x0)/L)# + u0

  c0 = u0 / Ma0
  temp0 = c0^2 / (equations.gamma*8.314)
  p0 = rho1 * 8.314 * temp0

    return prim2cons(SVector(rho, u, v, p0), equations)
end

initial_condition = initial_condition_BM_vortex

mesh = DGMultiMesh(solver, coordinates_min = (-0.5, -0.5), coordinates_max = (0.5, 0.5))

semi = SemidiscretizationHyperbolicParabolic(mesh, (equations, equations_parabolic), initial_condition, solver)

tspan = (0.0, 2.0)
ode = semidiscretize(semi, tspan)

summary_callback = SummaryCallback()
alive_callback = AliveCallback(alive_interval = 10)
analysis_interval = 100
analysis_callback = AnalysisCallback(semi, interval = analysis_interval, uEltype = real(solver))
stepsize_callback = StepsizeCallback(cfl = 0.3)
callbacks = CallbackSet(summary_callback, analysis_callback, alive_callback, stepsize_callback)

###############################################################################
# run the simulation

tol = 1.0e-8
sol = solve(ode, RDPK3SpFSAL49(); dt=0.5 * estimate_dt(mesh, solver), abstol = tol, reltol = tol,
save_everystep = true, callback = callbacks);

summary_callback() # print the timer summary
3

I think the combination of FD operators amd viscous terms is iust not supported at the moment - since we don’t have such an example in Trixi.jl right now. Correct me if I’m wrong, @jlchan
That being said, we would be happy to help you to make a PR to fix this if you’re interested

2 Likes
4

Hi ! Thanks for the quick response-
Unfortunatly I am not totally confident in my capabilities on the subject. I am currently an aerospace engineering student working with your package for a project and I am not sure I have the right tools and knowledge to implement something.
I still have trouble with the Julia programming language and not totally familliar with Trixi.jl architecture.
I find it very interesting to work on but I may not be qualified enough…

5

@Kelian_Renoux - yes, @ranocha is correct - the FD operators are not yet supported. However, I’d be happy to help with implementation - I’ve implemented FD viscous terms in an experimental code before using upwind/downwind SBP operators.

1 Like