Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Rescale compressible Navier-Stokes #1230

Merged
merged 12 commits into from
Oct 7, 2022
Merged

Rescale compressible Navier-Stokes #1230

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
12 commits
Select commit Hold shift + click to select a range
aa9ca7a
remove the nondimensionalization and update the docstring for CNS
andrewwinters5000 Oct 6, 2022
0b93d3e
update TreeMesh elixirs and test values
andrewwinters5000 Oct 6, 2022
b23b9a7
update the elixirs and test values for DGMulti
andrewwinters5000 Oct 6, 2022
df33b79
fix unbalanced parentheses
andrewwinters5000 Oct 6, 2022
356630e
file renaming for consistency
andrewwinters5000 Oct 6, 2022
c3620c6
update file names in the parabolic test script
andrewwinters5000 Oct 6, 2022
605b398
remove the gas constant as a parameter and update comments
andrewwinters5000 Oct 7, 2022
595ef59
rename viscosity to dynamic_viscosity
andrewwinters5000 Oct 7, 2022
f459934
add comment about units in the docstring
andrewwinters5000 Oct 7, 2022
507ebe8
Update src/equations/compressible_navier_stokes_2d.jl
andrewwinters5000 Oct 7, 2022
55aadd8
unify notation to call the viscosity mu
andrewwinters5000 Oct 7, 2022
40ef1c3
fix undefined variable
andrewwinters5000 Oct 7, 2022
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
61 changes: 30 additions & 31 deletions ...ti_2d/elixir_navier_stokes_convergence.jl → ...lti_2d/elixir_navierstokes_convergence.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,14 +4,14 @@ using Trixi
###############################################################################
# semidiscretization of the ideal compressible Navier-Stokes equations

reynolds_number() = 100
prandtl_number() = 0.72
mu() = 0.01

equations = CompressibleEulerEquations2D(1.4)
# Note: If you change the Navier-Stokes parameters here, also change them in the initial condition
# I really do not like this structure but it should work for now
equations_parabolic = CompressibleNavierStokesDiffusion2D(equations, Reynolds=reynolds_number(), Prandtl=prandtl_number(),
Mach_freestream=0.5, gradient_variables=GradientVariablesPrimitive())
equations_parabolic = CompressibleNavierStokesDiffusion2D(equations, mu=mu(), Prandtl=prandtl_number(),
gradient_variables=GradientVariablesPrimitive())

# Create DG solver with polynomial degree = 3 and (local) Lax-Friedrichs/Rusanov flux as surface flux
dg = DGMulti(polydeg = 3, element_type = Tri(), approximation_type = Polynomial(),
Expand Down Expand Up @@ -52,7 +52,7 @@ end
# see also https://github.com/trixi-framework/Trixi.jl/pull/1160
inv_gamma_minus_one = inv(equations.gamma - 1)
Pr = prandtl_number()
Re = reynolds_number()
mu_ = mu()

# Same settings as in `initial_condition`
# Amplitude and shift
Expand Down Expand Up @@ -104,7 +104,6 @@ end

# Some convenience constants
T_const = equations.gamma * inv_gamma_minus_one / Pr
inv_Re = 1.0 / Re
inv_rho_cubed = 1.0 / (rho^3)

# compute the source terms
Expand All @@ -118,50 +117,50 @@ end
+ rho * v1_y * v2
+ rho * v1 * v2_y
# stress tensor from x-direction
- 4.0 / 3.0 * v1_xx * inv_Re
+ 2.0 / 3.0 * v2_xy * inv_Re
- v1_yy * inv_Re
- v2_xy * inv_Re )
- 4.0 / 3.0 * v1_xx * mu_
+ 2.0 / 3.0 * v2_xy * mu_
- v1_yy * mu_
- v2_xy * mu_ )
# y-momentum equation
du3 = ( rho_t * v2 + rho * v2_t + p_y + rho_x * v1 * v2
+ rho * v1_x * v2
+ rho * v1 * v2_x
+ rho_y * v2^2
+ 2.0 * rho * v2 * v2_y
# stress tensor from y-direction
- v1_xy * inv_Re
- v2_xx * inv_Re
- 4.0 / 3.0 * v2_yy * inv_Re
+ 2.0 / 3.0 * v1_xy * inv_Re )
- v1_xy * mu_
- v2_xx * mu_
- 4.0 / 3.0 * v2_yy * mu_
+ 2.0 / 3.0 * v1_xy * mu_ )
# total energy equation
du4 = ( E_t + v1_x * (E + p) + v1 * (E_x + p_x)
+ v2_y * (E + p) + v2 * (E_y + p_y)
# stress tensor and temperature gradient terms from x-direction
- 4.0 / 3.0 * v1_xx * v1 * inv_Re
+ 2.0 / 3.0 * v2_xy * v1 * inv_Re
- 4.0 / 3.0 * v1_x * v1_x * inv_Re
+ 2.0 / 3.0 * v2_y * v1_x * inv_Re
- v1_xy * v2 * inv_Re
- v2_xx * v2 * inv_Re
- v1_y * v2_x * inv_Re
- v2_x * v2_x * inv_Re
- 4.0 / 3.0 * v1_xx * v1 * mu_
+ 2.0 / 3.0 * v2_xy * v1 * mu_
- 4.0 / 3.0 * v1_x * v1_x * mu_
+ 2.0 / 3.0 * v2_y * v1_x * mu_
- v1_xy * v2 * mu_
- v2_xx * v2 * mu_
- v1_y * v2_x * mu_
- v2_x * v2_x * mu_
- T_const * inv_rho_cubed * ( p_xx * rho * rho
- 2.0 * p_x * rho * rho_x
+ 2.0 * p * rho_x * rho_x
- p * rho * rho_xx ) * inv_Re
- p * rho * rho_xx ) * mu_
# stress tensor and temperature gradient terms from y-direction
- v1_yy * v1 * inv_Re
- v2_xy * v1 * inv_Re
- v1_y * v1_y * inv_Re
- v2_x * v1_y * inv_Re
- 4.0 / 3.0 * v2_yy * v2 * inv_Re
+ 2.0 / 3.0 * v1_xy * v2 * inv_Re
- 4.0 / 3.0 * v2_y * v2_y * inv_Re
+ 2.0 / 3.0 * v1_x * v2_y * inv_Re
- v1_yy * v1 * mu_
- v2_xy * v1 * mu_
- v1_y * v1_y * mu_
- v2_x * v1_y * mu_
- 4.0 / 3.0 * v2_yy * v2 * mu_
+ 2.0 / 3.0 * v1_xy * v2 * mu_
- 4.0 / 3.0 * v2_y * v2_y * mu_
+ 2.0 / 3.0 * v1_x * v2_y * mu_
- T_const * inv_rho_cubed * ( p_yy * rho * rho
- 2.0 * p_y * rho * rho_y
+ 2.0 * p * rho_y * rho_y
- p * rho * rho_yy ) * inv_Re )
- p * rho * rho_yy ) * mu_ )

return SVector(du1, du2, du3, du4)
end
Expand Down
11 changes: 5 additions & 6 deletions ...elixir_navier_stokes_lid_driven_cavity.jl → .../elixir_navierstokes_lid_driven_cavity.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5,14 +5,13 @@ using Trixi
# semidiscretization of the ideal compressible Navier-Stokes equations

# TODO: parabolic; unify names of these accessor functions
reynolds_number() = 1000.0
prandtl_number() = 0.72
mach_number() = 0.1
mu() = 0.001

equations = CompressibleEulerEquations2D(1.4)
equations_parabolic = CompressibleNavierStokesDiffusion2D(equations, Reynolds=reynolds_number(),
Prandtl=prandtl_number(),
Mach_freestream=mach_number())
equations_parabolic = CompressibleNavierStokesDiffusion2D(equations, mu=mu(),
Prandtl=prandtl_number())


# Create DG solver with polynomial degree = 3 and (local) Lax-Friedrichs/Rusanov flux as surface flux
dg = DGMulti(polydeg = 3, element_type = Quad(), approximation_type = GaussSBP(),
Expand All @@ -25,7 +24,7 @@ is_on_boundary = Dict(:top => top, :rest_of_boundary => rest_of_boundary)
mesh = DGMultiMesh(dg, cells_per_dimension=(16, 16); is_on_boundary)

function initial_condition_cavity(x, t, equations::CompressibleEulerEquations2D)
Ma = mach_number()
Ma = 0.1
rho = 1.0
u, v = 0, 0
p = 1.0 / (Ma^2 * equations.gamma)
Expand Down
61 changes: 30 additions & 31 deletions ...dgsem/elixir_navier_stokes_convergence.jl → ..._dgsem/elixir_navierstokes_convergence.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,12 +4,12 @@ using Trixi
###############################################################################
# semidiscretization of the ideal compressible Navier-Stokes equations

reynolds_number() = 100
prandtl_number() = 0.72
mu() = 0.01

equations = CompressibleEulerEquations2D(1.4)
equations_parabolic = CompressibleNavierStokesDiffusion2D(equations, Reynolds=reynolds_number(), Prandtl=prandtl_number(),
Mach_freestream=0.5, gradient_variables=GradientVariablesPrimitive())
equations_parabolic = CompressibleNavierStokesDiffusion2D(equations, mu=mu(), Prandtl=prandtl_number(),
gradient_variables=GradientVariablesPrimitive())

# Create DG solver with polynomial degree = 3 and (local) Lax-Friedrichs/Rusanov flux as surface flux
solver = DGSEM(polydeg=3, surface_flux=flux_lax_friedrichs,
Expand Down Expand Up @@ -54,7 +54,7 @@ end
# see also https://github.com/trixi-framework/Trixi.jl/pull/1160
inv_gamma_minus_one = inv(equations.gamma - 1)
Pr = prandtl_number()
Re = reynolds_number()
mu_ = mu()

# Same settings as in `initial_condition`
# Amplitude and shift
Expand Down Expand Up @@ -106,7 +106,6 @@ end

# Some convenience constants
T_const = equations.gamma * inv_gamma_minus_one / Pr
inv_Re = 1.0 / Re
inv_rho_cubed = 1.0 / (rho^3)

# compute the source terms
Expand All @@ -120,50 +119,50 @@ end
+ rho * v1_y * v2
+ rho * v1 * v2_y
# stress tensor from x-direction
- 4.0 / 3.0 * v1_xx * inv_Re
+ 2.0 / 3.0 * v2_xy * inv_Re
- v1_yy * inv_Re
- v2_xy * inv_Re )
- 4.0 / 3.0 * v1_xx * mu_
+ 2.0 / 3.0 * v2_xy * mu_
- v1_yy * mu_
- v2_xy * mu_ )
# y-momentum equation
du3 = ( rho_t * v2 + rho * v2_t + p_y + rho_x * v1 * v2
+ rho * v1_x * v2
+ rho * v1 * v2_x
+ rho_y * v2^2
+ 2.0 * rho * v2 * v2_y
# stress tensor from y-direction
- v1_xy * inv_Re
- v2_xx * inv_Re
- 4.0 / 3.0 * v2_yy * inv_Re
+ 2.0 / 3.0 * v1_xy * inv_Re )
- v1_xy * mu_
- v2_xx * mu_
- 4.0 / 3.0 * v2_yy * mu_
+ 2.0 / 3.0 * v1_xy * mu_ )
# total energy equation
du4 = ( E_t + v1_x * (E + p) + v1 * (E_x + p_x)
+ v2_y * (E + p) + v2 * (E_y + p_y)
# stress tensor and temperature gradient terms from x-direction
- 4.0 / 3.0 * v1_xx * v1 * inv_Re
+ 2.0 / 3.0 * v2_xy * v1 * inv_Re
- 4.0 / 3.0 * v1_x * v1_x * inv_Re
+ 2.0 / 3.0 * v2_y * v1_x * inv_Re
- v1_xy * v2 * inv_Re
- v2_xx * v2 * inv_Re
- v1_y * v2_x * inv_Re
- v2_x * v2_x * inv_Re
- 4.0 / 3.0 * v1_xx * v1 * mu_
+ 2.0 / 3.0 * v2_xy * v1 * mu_
- 4.0 / 3.0 * v1_x * v1_x * mu_
+ 2.0 / 3.0 * v2_y * v1_x * mu_
- v1_xy * v2 * mu_
- v2_xx * v2 * mu_
- v1_y * v2_x * mu_
- v2_x * v2_x * mu_
- T_const * inv_rho_cubed * ( p_xx * rho * rho
- 2.0 * p_x * rho * rho_x
+ 2.0 * p * rho_x * rho_x
- p * rho * rho_xx ) * inv_Re
- p * rho * rho_xx ) * mu_
# stress tensor and temperature gradient terms from y-direction
- v1_yy * v1 * inv_Re
- v2_xy * v1 * inv_Re
- v1_y * v1_y * inv_Re
- v2_x * v1_y * inv_Re
- 4.0 / 3.0 * v2_yy * v2 * inv_Re
+ 2.0 / 3.0 * v1_xy * v2 * inv_Re
- 4.0 / 3.0 * v2_y * v2_y * inv_Re
+ 2.0 / 3.0 * v1_x * v2_y * inv_Re
- v1_yy * v1 * mu_
- v2_xy * v1 * mu_
- v1_y * v1_y * mu_
- v2_x * v1_y * mu_
- 4.0 / 3.0 * v2_yy * v2 * mu_
+ 2.0 / 3.0 * v1_xy * v2 * mu_
- 4.0 / 3.0 * v2_y * v2_y * mu_
+ 2.0 / 3.0 * v1_x * v2_y * mu_
- T_const * inv_rho_cubed * ( p_yy * rho * rho
- 2.0 * p_y * rho * rho_y
+ 2.0 * p * rho_y * rho_y
- p * rho * rho_yy ) * inv_Re )
- p * rho * rho_yy ) * mu_ )

return SVector(du1, du2, du3, du4)
end
Expand Down
10 changes: 4 additions & 6 deletions ...elixir_navier_stokes_lid_driven_cavity.jl → .../elixir_navierstokes_lid_driven_cavity.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5,14 +5,12 @@ using Trixi
# semidiscretization of the ideal compressible Navier-Stokes equations

# TODO: parabolic; unify names of these accessor functions
reynolds_number() = 1000.0
prandtl_number() = 0.72
mach_number() = 0.1
mu() = 0.001

equations = CompressibleEulerEquations2D(1.4)
equations_parabolic = CompressibleNavierStokesDiffusion2D(equations, Reynolds=reynolds_number(),
Prandtl=prandtl_number(),
Mach_freestream=mach_number())
equations_parabolic = CompressibleNavierStokesDiffusion2D(equations, mu=mu(),
Prandtl=prandtl_number())

# Create DG solver with polynomial degree = 3 and (local) Lax-Friedrichs/Rusanov flux as surface flux
solver = DGSEM(polydeg=3, surface_flux=flux_lax_friedrichs)
Expand All @@ -28,7 +26,7 @@ mesh = TreeMesh(coordinates_min, coordinates_max,


function initial_condition_cavity(x, t, equations::CompressibleEulerEquations2D)
Ma = mach_number()
Ma = 0.1
rho = 1.0
u, v = 0, 0
p = 1.0 / (Ma^2 * equations.gamma)
Expand Down
1 change: 0 additions & 1 deletion src/equations/compressible_euler_2d.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -36,7 +36,6 @@ Here, ``\rho`` is the density, ``v_1``,`v_2` the velocities, ``e`` the specific
p = (\gamma - 1) \left( \rho e - \frac{1}{2} \rho (v_1^2+v_2^2) \right)
```
the pressure.

"""
struct CompressibleEulerEquations2D{RealT<:Real} <: AbstractCompressibleEulerEquations{2, 4}
gamma::RealT # ratio of specific heats
Expand Down
Loading