Skip to content

Commit

Permalink
Add numerical support of other real types (Float32) (#1909)
Browse files Browse the repository at this point in the history
* start with src/equations

* revise after the first review

* format src/equations

* another version

* revise after the second review

* apply suggestions from code review - comments

Co-authored-by: Hendrik Ranocha <[email protected]>

* remove TODO

* fix small errors

* complete compressible_euler

* fix error and format

* revise min max

* revise based on new docs

* revise based on new comments

* start sample test

* add more tests

* fix typos and comments

* change code review

* add to CI

* apply inferred macro

Co-authored-by: Hendrik Ranocha <[email protected]>

* complete

* add Trixi

---------

Co-authored-by: Hendrik Ranocha <[email protected]>
  • Loading branch information
huiyuxie and ranocha authored Jun 3, 2024
1 parent 58f9a8a commit 8665300
Show file tree
Hide file tree
Showing 6 changed files with 1,075 additions and 526 deletions.
56 changes: 30 additions & 26 deletions src/equations/acoustic_perturbation_2d.jl
Original file line number Diff line number Diff line change
Expand Up @@ -112,9 +112,9 @@ A constant initial condition where the state variables are zero and the mean flo
Uses the global mean values from `equations`.
"""
function initial_condition_constant(x, t, equations::AcousticPerturbationEquations2D)
v1_prime = 0.0
v2_prime = 0.0
p_prime_scaled = 0.0
v1_prime = 0
v2_prime = 0
p_prime_scaled = 0

return SVector(v1_prime, v2_prime, p_prime_scaled, global_mean_vars(equations)...)
end
Expand All @@ -127,12 +127,13 @@ A smooth initial condition used for convergence tests in combination with
"""
function initial_condition_convergence_test(x, t,
equations::AcousticPerturbationEquations2D)
c = 2.0
A = 0.2
L = 2.0
f = 2.0 / L
a = 1.0
omega = 2 * pi * f
RealT = eltype(x)
a = 1
c = 2
L = 2
f = 2.0f0 / L
A = convert(RealT, 0.2)
omega = 2 * convert(RealT, pi) * f
init = c + A * sin(omega * (x[1] + x[2] - a * t))

v1_prime = init
Expand All @@ -154,12 +155,13 @@ function source_terms_convergence_test(u, x, t,
equations::AcousticPerturbationEquations2D)
v1_mean, v2_mean, c_mean, rho_mean = cons2mean(u, equations)

c = 2.0
A = 0.2
L = 2.0
f = 2.0 / L
a = 1.0
omega = 2 * pi * f
RealT = eltype(u)
a = 1
c = 2
L = 2
f = 2.0f0 / L
A = convert(RealT, 0.2)
omega = 2 * convert(RealT, pi) * f

si, co = sincos(omega * (x[1] + x[2] - a * t))
tmp = v1_mean + v2_mean - a
Expand All @@ -168,7 +170,7 @@ function source_terms_convergence_test(u, x, t,
du3 = A * omega * co * (2 * c_mean^2 * rho_mean + 2 * c * tmp + 2 * A * tmp * si) /
c_mean^2

du4 = du5 = du6 = du7 = 0.0
du4 = du5 = du6 = du7 = 0

return SVector(du1, du2, du3, du4, du5, du6, du7)
end
Expand All @@ -179,8 +181,8 @@ end
A Gaussian pulse in a constant mean flow. Uses the global mean values from `equations`.
"""
function initial_condition_gauss(x, t, equations::AcousticPerturbationEquations2D)
v1_prime = 0.0
v2_prime = 0.0
v1_prime = 0
v2_prime = 0
p_prime = exp(-4 * (x[1]^2 + x[2]^2))

prim = SVector(v1_prime, v2_prime, p_prime, global_mean_vars(equations)...)
Expand Down Expand Up @@ -240,8 +242,8 @@ function boundary_condition_slip_wall(u_inner, normal_direction::AbstractVector,
u_normal = normal[1] * u_inner[1] + normal[2] * u_inner[2]

# create the "external" boundary solution state
u_boundary = SVector(u_inner[1] - 2.0 * u_normal * normal[1],
u_inner[2] - 2.0 * u_normal * normal[2],
u_boundary = SVector(u_inner[1] - 2 * u_normal * normal[1],
u_inner[2] - 2 * u_normal * normal[2],
u_inner[3], cons2mean(u_inner, equations)...)

# calculate the boundary flux
Expand All @@ -257,13 +259,14 @@ end
v1_mean, v2_mean, c_mean, rho_mean = cons2mean(u, equations)

# Calculate flux for conservative state variables
RealT = eltype(u)
if orientation == 1
f1 = v1_mean * v1_prime + v2_mean * v2_prime +
c_mean^2 * p_prime_scaled / rho_mean
f2 = zero(eltype(u))
f2 = zero(RealT)
f3 = rho_mean * v1_prime + v1_mean * p_prime_scaled
else
f1 = zero(eltype(u))
f1 = zero(RealT)
f2 = v1_mean * v1_prime + v2_mean * v2_prime +
c_mean^2 * p_prime_scaled / rho_mean
f3 = rho_mean * v2_prime + v2_mean * p_prime_scaled
Expand All @@ -272,7 +275,7 @@ end
# The rest of the state variables are actually variable coefficients, hence the flux should be
# zero. See https://github.com/trixi-framework/Trixi.jl/issues/358#issuecomment-784828762
# for details.
f4 = f5 = f6 = f7 = zero(eltype(u))
f4 = f5 = f6 = f7 = 0

return SVector(f1, f2, f3, f4, f5, f6, f7)
end
Expand Down Expand Up @@ -313,7 +316,7 @@ end
# The rest of the state variables are actually variable coefficients, hence the flux should be
# zero. See https://github.com/trixi-framework/Trixi.jl/issues/358#issuecomment-784828762
# for details.
f4 = f5 = f6 = f7 = zero(eltype(u))
f4 = f5 = f6 = f7 = 0

return SVector(f1, f2, f3, f4, f5, f6, f7)
end
Expand Down Expand Up @@ -344,8 +347,9 @@ end
equations::AcousticPerturbationEquations2D)
λ = dissipation.max_abs_speed(u_ll, u_rr, orientation_or_normal_direction,
equations)
diss = -0.5 * λ * (u_rr - u_ll)
z = zero(eltype(u_ll))
diss = -0.5f0 * λ * (u_rr - u_ll)
z = 0

return SVector(diss[1], diss[2], diss[3], z, z, z, z)
end

Expand Down
Loading

0 comments on commit 8665300

Please sign in to comment.