From 98a5f6c1eab0237ff2be2d66d01d7da720ae9f05 Mon Sep 17 00:00:00 2001 From: Sven Berger Date: Fri, 13 Dec 2024 01:17:24 +0100 Subject: [PATCH 1/9] add boundary system --- test/schemes/fluid/surface_normal_sph.jl | 73 +++++++++++++++++++++--- 1 file changed, 65 insertions(+), 8 deletions(-) diff --git a/test/schemes/fluid/surface_normal_sph.jl b/test/schemes/fluid/surface_normal_sph.jl index 904096000..8cdea7d1c 100644 --- a/test/schemes/fluid/surface_normal_sph.jl +++ b/test/schemes/fluid/surface_normal_sph.jl @@ -1,23 +1,79 @@ +function create_boundary_system(coordinates, particle_spacing, state_equation, kernel, + smoothing_length, NDIMS, walldistance) + # Compute bounding box of fluid particles + xmin = minimum(coordinates[1, :]) + xmax = maximum(coordinates[1, :]) + ymin = minimum(coordinates[2, :]) + ymax = maximum(coordinates[2, :]) + + wall_thickness = 4 * particle_spacing + + if NDIMS == 2 + # In 2D: The wall extends in x-direction and has a thickness in y. + wall_width = xmax - xmin + wall_size = (wall_width, wall_thickness) + wall_coord = (xmin, ymin - walldistance) + elseif NDIMS == 3 + # In 3D: The wall extends in x and z directions, and has thickness in y. + zmin = minimum(coordinates[3, :]) + wall_width_x = xmax - xmin + wall_width_y = ymax - ymin + wall_size = (wall_width_x, wall_width_y, wall_thickness) + wall_coord = (xmin, ymin, zmin - walldistance) + end + + # Create the wall shape + wall = RectangularShape(particle_spacing, + round.(Int, wall_size ./ particle_spacing), + wall_coord, + density=1000.0) + + boundary_model = BoundaryModelDummyParticles(wall.density, + wall.mass, + state_equation=state_equation, + AdamiPressureExtrapolation(), + kernel, + smoothing_length, + correction=nothing) + + boundary_system = BoundarySPHSystem(wall, boundary_model, adhesion_coefficient=0.0) + return boundary_system +end + function create_fluid_system(coordinates, velocity, mass, density, particle_spacing; - NDIMS=2, smoothing_length=1.0) + NDIMS=2, smoothing_length=1.0, wall=false, walldistance=0.0) tspan = (0.0, 0.01) - fluid = InitialCondition(coordinates=coordinates, velocity=velocity, mass=mass, - density=density, particle_spacing=particle_spacing) + fluid = InitialCondition(coordinates=coordinates, + velocity=velocity, + mass=mass, + density=density, + particle_spacing=particle_spacing) state_equation = StateEquationCole(sound_speed=10.0, reference_density=1000.0, exponent=1) - system = WeaklyCompressibleSPHSystem(fluid, SummationDensity(), state_equation, - SchoenbergCubicSplineKernel{NDIMS}(), + kernel = SchoenbergCubicSplineKernel{NDIMS}() + + system = WeaklyCompressibleSPHSystem(fluid, + SummationDensity(), + state_equation, + kernel, smoothing_length; surface_normal_method=ColorfieldSurfaceNormal(), reference_particle_spacing=particle_spacing) - semi = Semidiscretization(system) - ode = semidiscretize(semi, tspan) + if wall + boundary_system = create_boundary_system(coordinates, particle_spacing, + state_equation, kernel, smoothing_length, + NDIMS, walldistance) + semi = Semidiscretization(system, boundary_system) + else + semi = Semidiscretization(system) + end + ode = semidiscretize(semi, tspan) TrixiParticles.update_systems_and_nhs(ode.u0.x..., semi, 0.0) return system, semi, ode @@ -96,7 +152,8 @@ end system, semi, ode = create_fluid_system(coordinates, velocity, mass, density, particle_spacing; NDIMS=NDIMS, - smoothing_length=3.0 * particle_spacing) + smoothing_length=3.0 * particle_spacing, + wall=true, walldistance=2.0) compute_and_test_surface_normals(system, semi, ode; NDIMS=NDIMS) From e05edea8b27287dae8d9d2a663185d9ca6e6a3c0 Mon Sep 17 00:00:00 2001 From: Sven Berger Date: Fri, 13 Dec 2024 01:46:25 +0100 Subject: [PATCH 2/9] update --- test/schemes/fluid/surface_normal_sph.jl | 12 +++++++++--- 1 file changed, 9 insertions(+), 3 deletions(-) diff --git a/test/schemes/fluid/surface_normal_sph.jl b/test/schemes/fluid/surface_normal_sph.jl index 8cdea7d1c..f596da17b 100644 --- a/test/schemes/fluid/surface_normal_sph.jl +++ b/test/schemes/fluid/surface_normal_sph.jl @@ -71,12 +71,13 @@ function create_fluid_system(coordinates, velocity, mass, density, particle_spac semi = Semidiscretization(system, boundary_system) else semi = Semidiscretization(system) + boundary_system = nothing end ode = semidiscretize(semi, tspan) TrixiParticles.update_systems_and_nhs(ode.u0.x..., semi, 0.0) - return system, semi, ode + return system, boundary_system, semi, ode end function compute_and_test_surface_normals(system, semi, ode; NDIMS=2) @@ -126,7 +127,7 @@ end fluid_density = 1000.0 density = fill(fluid_density, nparticles) - system, semi, ode = create_fluid_system(coordinates, velocity, mass, density, + system, bnd_system, semi, ode = create_fluid_system(coordinates, velocity, mass, density, particle_spacing; NDIMS=NDIMS) @@ -149,7 +150,7 @@ end density = sphere_ic.density # To get somewhat accurate normals we increase the smoothing length unrealistically - system, semi, ode = create_fluid_system(coordinates, velocity, mass, density, + system, bnd_system, semi, ode = create_fluid_system(coordinates, velocity, mass, density, particle_spacing; NDIMS=NDIMS, smoothing_length=3.0 * particle_spacing, @@ -186,6 +187,11 @@ end end end + # Boundary system + bnd_color = bnd_system.boundary_model.cache.colorfield_bnd + # this is only true since it assumed that the color is 1 + @test all(bnd_color .>= 0.0) + # Compare computed normals to expected normals for surface particles for i in surface_particles @test isapprox(computed_normals[:, i], expected_normals[:, i], atol=0.04) From 36741d73949be08a30403b02776bc1792837b23c Mon Sep 17 00:00:00 2001 From: Sven Berger Date: Fri, 13 Dec 2024 02:21:30 +0100 Subject: [PATCH 3/9] update --- test/schemes/fluid/surface_normal_sph.jl | 2 -- 1 file changed, 2 deletions(-) diff --git a/test/schemes/fluid/surface_normal_sph.jl b/test/schemes/fluid/surface_normal_sph.jl index f596da17b..7f49690f7 100644 --- a/test/schemes/fluid/surface_normal_sph.jl +++ b/test/schemes/fluid/surface_normal_sph.jl @@ -9,12 +9,10 @@ function create_boundary_system(coordinates, particle_spacing, state_equation, k wall_thickness = 4 * particle_spacing if NDIMS == 2 - # In 2D: The wall extends in x-direction and has a thickness in y. wall_width = xmax - xmin wall_size = (wall_width, wall_thickness) wall_coord = (xmin, ymin - walldistance) elseif NDIMS == 3 - # In 3D: The wall extends in x and z directions, and has thickness in y. zmin = minimum(coordinates[3, :]) wall_width_x = xmax - xmin wall_width_y = ymax - ymin From ba8765f9d7d0c686654f45ae04257e5077ec488c Mon Sep 17 00:00:00 2001 From: Sven Berger Date: Fri, 13 Dec 2024 02:22:57 +0100 Subject: [PATCH 4/9] format --- test/schemes/fluid/surface_normal_sph.jl | 19 +++++++++++-------- 1 file changed, 11 insertions(+), 8 deletions(-) diff --git a/test/schemes/fluid/surface_normal_sph.jl b/test/schemes/fluid/surface_normal_sph.jl index 7f49690f7..e248eaa44 100644 --- a/test/schemes/fluid/surface_normal_sph.jl +++ b/test/schemes/fluid/surface_normal_sph.jl @@ -125,9 +125,10 @@ end fluid_density = 1000.0 density = fill(fluid_density, nparticles) - system, bnd_system, semi, ode = create_fluid_system(coordinates, velocity, mass, density, - particle_spacing; - NDIMS=NDIMS) + system, bnd_system, semi, ode = create_fluid_system(coordinates, velocity, mass, + density, + particle_spacing; + NDIMS=NDIMS) compute_and_test_surface_normals(system, semi, ode; NDIMS=NDIMS) end @@ -148,11 +149,13 @@ end density = sphere_ic.density # To get somewhat accurate normals we increase the smoothing length unrealistically - system, bnd_system, semi, ode = create_fluid_system(coordinates, velocity, mass, density, - particle_spacing; - NDIMS=NDIMS, - smoothing_length=3.0 * particle_spacing, - wall=true, walldistance=2.0) + system, bnd_system, semi, ode = create_fluid_system(coordinates, velocity, mass, + density, + particle_spacing; + NDIMS=NDIMS, + smoothing_length=3.0 * + particle_spacing, + wall=true, walldistance=2.0) compute_and_test_surface_normals(system, semi, ode; NDIMS=NDIMS) From 9936e5d841d26dd06badfa59f3c0f57e5e093db2 Mon Sep 17 00:00:00 2001 From: Sven Berger Date: Fri, 13 Dec 2024 02:52:48 +0100 Subject: [PATCH 5/9] format --- test/schemes/fluid/surface_normal_sph.jl | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/test/schemes/fluid/surface_normal_sph.jl b/test/schemes/fluid/surface_normal_sph.jl index fe6d06662..088e279e3 100644 --- a/test/schemes/fluid/surface_normal_sph.jl +++ b/test/schemes/fluid/surface_normal_sph.jl @@ -142,8 +142,8 @@ end (2, SchoenbergCubicSplineKernel{2}(), 0.25, 3.0, 1.0, (0.0, 0.0)), (2, SchoenbergCubicSplineKernel{2}(), 0.1, 3.5, 1.0, (0.0, 0.0)), (3, SchoenbergCubicSplineKernel{3}(), 0.25, 3.0, 1.0, (0.0, 0.0, 0.0)), - (2, WendlandC2Kernel{2}(), 0.3, 3.0, 1.0, (0.0, 0.0)), - (3, WendlandC2Kernel{3}(), 0.3, 3.0, 1.0, (0.0, 0.0, 0.0)) + (2, WendlandC2Kernel{2}(), 0.3, 3.0, 1.0, (0.0, 0.0)), + (3, WendlandC2Kernel{3}(), 0.3, 3.0, 1.0, (0.0, 0.0, 0.0)) ] for (NDIMS, smoothing_kernel, particle_spacing, smoothing_length_mult, radius, center) in variations @@ -165,7 +165,7 @@ end smoothing_kernel=smoothing_kernel, surface_normal_method=ColorfieldSurfaceNormal(interface_threshold=0.1, ideal_density_threshold=0.9), - wall=true, walldistance=2.0) + wall=true, walldistance=2.0) compute_and_test_surface_normals(system, semi, ode; NDIMS=NDIMS) @@ -197,10 +197,10 @@ end end end - # Boundary system - bnd_color = bnd_system.boundary_model.cache.colorfield_bnd - # this is only true since it assumed that the color is 1 - @test all(bnd_color .>= 0.0) + # Boundary system + bnd_color = bnd_system.boundary_model.cache.colorfield_bnd + # this is only true since it assumed that the color is 1 + @test all(bnd_color .>= 0.0) # Test that computed normals match expected normals for i in surface_particles From 4dc0cf8d0440b04447664abed9cb2c526409e830 Mon Sep 17 00:00:00 2001 From: Sven Berger Date: Fri, 13 Dec 2024 02:55:17 +0100 Subject: [PATCH 6/9] fix merge --- test/schemes/fluid/surface_normal_sph.jl | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/test/schemes/fluid/surface_normal_sph.jl b/test/schemes/fluid/surface_normal_sph.jl index 088e279e3..f470e1b4f 100644 --- a/test/schemes/fluid/surface_normal_sph.jl +++ b/test/schemes/fluid/surface_normal_sph.jl @@ -1,3 +1,4 @@ +include("../../test_util.jl") function create_boundary_system(coordinates, particle_spacing, state_equation, kernel, smoothing_length, NDIMS, walldistance) # Compute bounding box of fluid particles @@ -63,7 +64,7 @@ function create_fluid_system(coordinates, velocity, mass, density, particle_spac if wall boundary_system = create_boundary_system(coordinates, particle_spacing, - state_equation, kernel, smoothing_length, + state_equation, smoothing_kernel, smoothing_length, NDIMS, walldistance) semi = Semidiscretization(system, boundary_system) else From 5e83f0ced9a05789e9753f68e65a0e6461b37e79 Mon Sep 17 00:00:00 2001 From: Sven Berger Date: Fri, 13 Dec 2024 02:55:27 +0100 Subject: [PATCH 7/9] format --- test/schemes/fluid/surface_normal_sph.jl | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/test/schemes/fluid/surface_normal_sph.jl b/test/schemes/fluid/surface_normal_sph.jl index f470e1b4f..2c8a55c7a 100644 --- a/test/schemes/fluid/surface_normal_sph.jl +++ b/test/schemes/fluid/surface_normal_sph.jl @@ -64,7 +64,8 @@ function create_fluid_system(coordinates, velocity, mass, density, particle_spac if wall boundary_system = create_boundary_system(coordinates, particle_spacing, - state_equation, smoothing_kernel, smoothing_length, + state_equation, smoothing_kernel, + smoothing_length, NDIMS, walldistance) semi = Semidiscretization(system, boundary_system) else From 55e6557eb9df6f6f46d98e78cded4f8af3e11a86 Mon Sep 17 00:00:00 2001 From: Sven Berger Date: Fri, 13 Dec 2024 03:00:44 +0100 Subject: [PATCH 8/9] another merge error --- test/schemes/fluid/surface_normal_sph.jl | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/test/schemes/fluid/surface_normal_sph.jl b/test/schemes/fluid/surface_normal_sph.jl index 2c8a55c7a..4b63763cf 100644 --- a/test/schemes/fluid/surface_normal_sph.jl +++ b/test/schemes/fluid/surface_normal_sph.jl @@ -1,4 +1,3 @@ -include("../../test_util.jl") function create_boundary_system(coordinates, particle_spacing, state_equation, kernel, smoothing_length, NDIMS, walldistance) # Compute bounding box of fluid particles @@ -160,7 +159,7 @@ end mass = sphere_ic.mass density = sphere_ic.density - system, semi, ode = create_fluid_system(coordinates, velocity, mass, density, + system, bnd_system, semi, ode = create_fluid_system(coordinates, velocity, mass, density, particle_spacing, nothing; NDIMS=NDIMS, smoothing_length=smoothing_length, From ee5c2da4147c6c7e3610fc29b2628191fff7c430 Mon Sep 17 00:00:00 2001 From: Sven Berger Date: Fri, 13 Dec 2024 03:00:57 +0100 Subject: [PATCH 9/9] format --- test/schemes/fluid/surface_normal_sph.jl | 17 +++++++++-------- 1 file changed, 9 insertions(+), 8 deletions(-) diff --git a/test/schemes/fluid/surface_normal_sph.jl b/test/schemes/fluid/surface_normal_sph.jl index 4b63763cf..fdf861ddb 100644 --- a/test/schemes/fluid/surface_normal_sph.jl +++ b/test/schemes/fluid/surface_normal_sph.jl @@ -159,14 +159,15 @@ end mass = sphere_ic.mass density = sphere_ic.density - system, bnd_system, semi, ode = create_fluid_system(coordinates, velocity, mass, density, - particle_spacing, nothing; - NDIMS=NDIMS, - smoothing_length=smoothing_length, - smoothing_kernel=smoothing_kernel, - surface_normal_method=ColorfieldSurfaceNormal(interface_threshold=0.1, - ideal_density_threshold=0.9), - wall=true, walldistance=2.0) + system, bnd_system, semi, ode = create_fluid_system(coordinates, velocity, mass, + density, + particle_spacing, nothing; + NDIMS=NDIMS, + smoothing_length=smoothing_length, + smoothing_kernel=smoothing_kernel, + surface_normal_method=ColorfieldSurfaceNormal(interface_threshold=0.1, + ideal_density_threshold=0.9), + wall=true, walldistance=2.0) compute_and_test_surface_normals(system, semi, ode; NDIMS=NDIMS)