I am working on a local discontinuous Galerkin discretization of the linearized Navier Stokes equations and comparing it to an interior penalty dG discretization. I have been able to successfully implement an LDG discretization for the momentum equation, but am struggling to implement an LDG for the energy equation. The example I provide to describe the problem is a version where the momentum equation is discretized using the IP method which I know works properly, and the energy equation is discretized using LDG. The energy equation is as follows:
The first order terms and the viscous fluctuation terms are discretized using a Lax-Friedrichs flux. The \dot\omega_T' term is not included as the example flow is non-reactive. I will now present my weak formulation of the second order terms. v is the scalar test function, q is the vector trial function (q = -\bar\alpha\frac{\partial h'}{\partial x_j}) and w is the vector test function. \mathcal{F}^i is the set of all element faces which are not on \partial \Omega. The weak formulation implemented in FEniCS is the sum of the following two equations:
For boundary conditions, parts of the domain have a homogenous Dirichlet condition (\mathcal{F}^D) and others have a homogenous Neumann condition (\mathcal{F}^N). To my understanding, one would weakly impose the homogeneous Dirichlet condition by the addition of the following terms:
I believe the homogeneous Neumann condition should be weakly implemented by the addition of this term:
The included figure shows the results of the simulation. The top half is the solution using the LDG implementation and the bottom half is the solution using the IP method. The areas of the boundary which are marked with a red line indicate a homogenous Dirichlet condition and the portions with a blue line indicate a homogeneous Neumann condition. The condition along the boundary where the two solutions meet is also a homogeneous Neumann condition.
I believe my issue lies in the application of the boundary conditions, as the solution inside the domain seems to be similar between the two formulations, but when it approaches the boundary regions, my LDG formulation doesn’t produce the correct results.
(A zoomed in view showing inconsistencies along both the Dirichlet and Neumann portions of the boundary).
I understand it might be a difficult to provide feedback due to the limited amount of information provided and the lack of code to execute, but any tips, literature references, or other ideas would be greatly appreciated.