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

Add test for Adami Pressure Extrapolation. #692

Open
wants to merge 3 commits into
base: main
Choose a base branch
from
Open

Add test for Adami Pressure Extrapolation. #692

Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
aac3a4d
Add test for Adami Pressure Extrapolation.
RubberLanding Jan 2, 2025
f450d96
Format file.
RubberLanding Jan 2, 2025
d5865dd
Provide code to debug the malfunctioning tests.
RubberLanding Jan 6, 2025
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
161 changes: 161 additions & 0 deletions test/schemes/boundary/dummy_particles/dummy_particles.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -188,6 +188,167 @@
end
end
end
@testset "Pressure Extrapolation Adami" begin
@testset "Pressure Extrapolation Adami: Constant Pressure" begin
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
@testset "Pressure Extrapolation Adami: Constant Pressure" begin
@testset "Constant Pressure" begin

This is already a subset of "Pressure Extrapolation Adami".

#=
Test whether the pressure is constant and 0.0, if the density of the state equation and in the tank are the same.
=#
particle_spacing = 0.1
n_particles = 4
n_layers = 3
width = particle_spacing * n_particles
height = particle_spacing * n_particles

density = 260
tank1 = RectangularTank(particle_spacing, (width, height), (width, height),
density, n_layers=n_layers,
faces=(true, true, true, false))

smoothing_kernel = SchoenbergCubicSplineKernel{2}()
smoothing_length = 1.2 * particle_spacing
viscosity = ViscosityAdami(nu=1e-6)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
viscosity = ViscosityAdami(nu=1e-6)

I would remove the extra complexity of the viscosity here.

state_equation = StateEquationCole(sound_speed=10, reference_density=257,
RubberLanding marked this conversation as resolved.
Show resolved Hide resolved
exponent=7)

boundary_model = BoundaryModelDummyParticles(tank1.boundary.density,
tank1.boundary.mass,
state_equation=state_equation,
AdamiPressureExtrapolation(),
smoothing_kernel, smoothing_length,
viscosity=viscosity)

boundary_system = BoundarySPHSystem(tank1.boundary, boundary_model)

fluid_system = WeaklyCompressibleSPHSystem(tank1.fluid, SummationDensity(),
state_equation,
smoothing_kernel, smoothing_length)
fluid_system.cache.density .= tank1.fluid.density
v_fluid = zeros(2, TrixiParticles.nparticles(fluid_system))

TrixiParticles.compute_pressure!(fluid_system, v_fluid)

neighborhood_search = TrixiParticles.TrivialNeighborhoodSearch{2}(search_radius=TrixiParticles.compact_support(smoothing_kernel,
smoothing_length),
eachpoint=TrixiParticles.eachparticle(fluid_system))

viscosity = boundary_system.boundary_model.viscosity
TrixiParticles.set_zero!(boundary_model.pressure)
TrixiParticles.reset_cache!(boundary_system.boundary_model.cache,
viscosity)

TrixiParticles.boundary_pressure_extrapolation!(boundary_model, boundary_system,
fluid_system,
tank1.boundary.coordinates,
tank1.fluid.coordinates,
v_fluid,
v_fluid,
neighborhood_search)

@test all(boundary_system.boundary_model.pressure .== boundary_system.boundary_model.pressure[1]) &
all(fluid_system.pressure .== fluid_system.pressure[1])

# TrixiParticles.@autoinfiltrate

#=
Test whether the pressure is constant, if the density of the state equation and in the tank are not the same.
=#
Comment on lines +252 to +254
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
#=
Test whether the pressure is constant, if the density of the state equation and in the tank are not the same.
=#
# Test whether the pressure is constant if the density of the state equation
# and in the tank are not the same.

tank2 = RectangularTank(particle_spacing, (width, height), (width, height),
density, n_layers=n_layers,
faces=(true, true, true, false))
Comment on lines +255 to +257
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is still the same density, or did I miss something?
Shouldn't this be

Suggested change
tank2 = RectangularTank(particle_spacing, (width, height), (width, height),
density, n_layers=n_layers,
faces=(true, true, true, false))
density = 310
tank2 = RectangularTank(particle_spacing, (width, height), (width, height),
density, n_layers=n_layers,
faces=(true, true, true, false))


fluid_system2 = WeaklyCompressibleSPHSystem(tank2.fluid, SummationDensity(),
state_equation,
smoothing_kernel, smoothing_length)
fluid_system2.cache.density .= tank2.fluid.density
v_fluid2 = zeros(2, TrixiParticles.nparticles(fluid_system2))
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You can use v_fluid again, no?

TrixiParticles.compute_pressure!(fluid_system2, v_fluid2)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is only computing the pressure for the fluid system, not for the boundary. I guess it worked because you didn't actually change the density (see my comment above).


@test all(boundary_system.boundary_model.pressure .==
boundary_system.boundary_model.pressure[1]) &
all(fluid_system.pressure .== fluid_system.pressure[1])
end

@testset "Pressure Extrapolation Adami: Linear Pressure" begin
#=
Test whether ...
=#
particle_spacing = 0.1
n_particles = 2
n_layers = 1
width = particle_spacing * n_particles
height = particle_spacing * n_particles
Comment on lines +275 to +279
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why are you not using the same sizes as before?


density = 257

state_equation = StateEquationCole(sound_speed=10, reference_density=257,
exponent=7)

tank1 = RectangularTank(particle_spacing, (width, height), (width, height),
density, acceleration=[0.0, -9.81],
state_equation=state_equation, n_layers=n_layers,
faces=(true, true, true, false))

smoothing_kernel = SchoenbergCubicSplineKernel{2}()
smoothing_length = 1.2 * particle_spacing
viscosity = ViscosityAdami(nu=1e-6)

boundary_model = BoundaryModelDummyParticles(tank1.boundary.density,
tank1.boundary.mass,
state_equation=state_equation,
AdamiPressureExtrapolation(),
smoothing_kernel, smoothing_length,
viscosity=viscosity)

boundary_system = BoundarySPHSystem(tank1.boundary, boundary_model)

fluid_system1 = WeaklyCompressibleSPHSystem(tank1.fluid, SummationDensity(),
state_equation,
smoothing_kernel, smoothing_length)
fluid_system1.cache.density .= tank1.fluid.density
v_fluid1 = zeros(2, TrixiParticles.nparticles(fluid_system1))
TrixiParticles.compute_pressure!(fluid_system1, v_fluid1)

neighborhood_search = TrixiParticles.TrivialNeighborhoodSearch{2}(search_radius=TrixiParticles.compact_support(smoothing_kernel,
smoothing_length),
eachpoint=TrixiParticles.eachparticle(fluid_system1))

viscosity = boundary_system.boundary_model.viscosity

TrixiParticles.set_zero!(boundary_model.pressure)
TrixiParticles.reset_cache!(boundary_system.boundary_model.cache,
viscosity)

TrixiParticles.adami_pressure_extrapolation!(boundary_model, boundary_system,
fluid_system1,
tank1.boundary.coordinates,
tank1.fluid.coordinates,
v_fluid1,
neighborhood_search)

width2 = particle_spacing * (n_particles + 2 * n_layers)
height2 = particle_spacing * (n_particles + n_layers)

tank2 = RectangularTank(particle_spacing, (width2, height2), (width2, height2),
density, acceleration=[0.0, -9.81],
state_equation=state_equation, n_layers=0,
faces=(true, true, true, false))

fluid_system2 = WeaklyCompressibleSPHSystem(tank2.fluid, SummationDensity(),
state_equation,
smoothing_kernel, smoothing_length)

fluid_system2.cache.density .= tank2.fluid.density
v_fluid2 = zeros(2, TrixiParticles.nparticles(fluid_system2))
TrixiParticles.compute_pressure!(fluid_system2, v_fluid2)

# CHECK IF SLICING IS CORRECT
press1 = transpose(reshape(fluid_system1.pressure, (n_particles, n_particles)))
press2 = transpose(reshape(fluid_system2.pressure,
(n_particles + 2 * n_layers, n_particles + n_layers)))[(1 + n_layers):(n_layers + n_particles),
(1 + n_layers):(n_particles + n_layers)]
@test press1 == press2
end
end
end

include("rhs.jl")
Loading