From 540e9a6be6a8511d73acd728a69c70c0cd764c57 Mon Sep 17 00:00:00 2001 From: Kai Partmann Date: Tue, 27 Feb 2024 12:56:03 +0100 Subject: [PATCH] Add halo exchange for write fields --- src/core/threads_data_handler.jl | 21 +++- test/integration/material_interface.jl | 129 +++++++++++++++++++++++++ test/integration/test_material.jl | 61 ++++++++++++ 3 files changed, 210 insertions(+), 1 deletion(-) create mode 100644 test/integration/material_interface.jl create mode 100644 test/integration/test_material.jl diff --git a/src/core/threads_data_handler.jl b/src/core/threads_data_handler.jl index 7c41064a..58dff9e1 100644 --- a/src/core/threads_data_handler.jl +++ b/src/core/threads_data_handler.jl @@ -41,7 +41,7 @@ function exchange_read_fields!(dh::ThreadsDataHandler, chunk_id::Int) end function exchange_write_fields!(dh::ThreadsDataHandler, chunk_id::Int) - halo_exchange!(dh, dh.write_halo_exs[chunk_id]) + halo_exchange_add!(dh, dh.write_halo_exs[chunk_id]) return nothing end @@ -54,6 +54,15 @@ function halo_exchange!(dh::ThreadsDataHandler, halo_exs::Vector{HaloExchange}) return nothing end +function halo_exchange_add!(dh::ThreadsDataHandler, halo_exs::Vector{HaloExchange}) + for he in halo_exs + src_field = get_exchange_field(dh.chunks[he.src_chunk_id], he.field) + dest_field = get_exchange_field(dh.chunks[he.dest_chunk_id], he.field) + exchange_add!(dest_field, src_field, he.dest_idxs, he.src_idxs) + end + return nothing +end + @inline function get_exchange_field(b::AbstractBodyChunk, fieldname::Symbol) return getfield(b.store, fieldname)::Matrix{Float64} end @@ -67,3 +76,13 @@ function exchange!(dest::Matrix{T}, src::Matrix{T}, dest_idxs::Vector{Int}, end return nothing end + +function exchange_add!(dest::Matrix{T}, src::Matrix{T}, dest_idxs::Vector{Int}, + src_idxs::Vector{Int}) where {T<:Number} + for i in eachindex(dest_idxs) + for d in axes(dest, 1) + @inbounds dest[d, dest_idxs[i]] += src[d, src_idxs[i]] + end + end + return nothing +end diff --git a/test/integration/material_interface.jl b/test/integration/material_interface.jl new file mode 100644 index 00000000..2feb4743 --- /dev/null +++ b/test/integration/material_interface.jl @@ -0,0 +1,129 @@ + +@testitem "TestMaterial halo exchange" begin + include("test_material.jl") + + position = [0.0 1.0 0.0 0.0 + 0.0 0.0 1.0 0.0 + 0.0 0.0 0.0 1.0] + volume = [1.1, 1.2, 1.3, 1.4] + mat = TestMaterial() + body = Body(mat, position, volume) + material!(body, horizon=2, rho=1, E=1, nu=0.25, Gc=1) + point_set!(body, :a, 1:2) + point_set!(body, :b, 3:4) + velocity_ic!(body, :a, :x, 1.0) + velocity_bc!(t -> t, body, :a, :x) + forcedensity_bc!(t -> t, body, :a, :x) + precrack!(body, :a, :b) + ts = VelocityVerlet(steps=10) + point_decomp = Peridynamics.PointDecomposition(body, 2) + tdh = Peridynamics.ThreadsDataHandler(body, ts, point_decomp) + + b1 = tdh.chunks[1] + @test b1 isa Peridynamics.BodyChunk + @test b1.store.position == position + @test b1.store.displacement == zeros(3, 2) + @test b1.store.velocity == [1.0 1.0; 0.0 0.0; 0.0 0.0] + @test b1.store.velocity_half == zeros(3, 2) + @test b1.store.acceleration == zeros(3, 2) + @test b1.store.b_int == zeros(3, 4) + @test b1.store.b_ext == zeros(3, 2) + @test b1.store.damage ≈ [2/3, 2/3] + @test b1.store.bond_active == [1, 0, 0, 1, 0, 0] + @test b1.store.n_active_bonds == [1, 1] + b2 = tdh.chunks[2] + @test b2 isa Peridynamics.BodyChunk + @test b2.store.position == position[:, [3, 4, 1, 2]] + @test b2.store.displacement == zeros(3, 2) + @test b2.store.velocity == [0.0 0.0; 0.0 0.0; 0.0 0.0] + @test b2.store.velocity_half == zeros(3, 2) + @test b2.store.acceleration == zeros(3, 2) + @test b2.store.b_int == zeros(3, 4) + @test b2.store.b_ext == zeros(3, 2) + @test b2.store.damage ≈ [2/3, 2/3] + @test b2.store.bond_active == [0, 0, 1, 0, 0, 1] + @test b2.store.n_active_bonds == [1, 1] + + randpos = rand(3, 4) + tdh.chunks[2].store.position .= randpos + Peridynamics.exchange_read_fields!(tdh, 1) + + @test b1.store.position[:,1:2] ≈ [0.0 1.0; 0.0 0.0; 0.0 0.0] + @test b1.store.position[:,3:4] ≈ randpos[:,1:2] + @test b1.store.displacement == zeros(3, 2) + @test b1.store.velocity == [1.0 1.0; 0.0 0.0; 0.0 0.0] + @test b1.store.velocity_half == zeros(3, 2) + @test b1.store.acceleration == zeros(3, 2) + @test b1.store.b_int == zeros(3, 4) + @test b1.store.b_ext == zeros(3, 2) + @test b1.store.damage ≈ [2/3, 2/3] + @test b1.store.bond_active == [1, 0, 0, 1, 0, 0] + @test b1.store.n_active_bonds == [1, 1] + + @test b2.store.position ≈ randpos + @test b2.store.displacement == zeros(3, 2) + @test b2.store.velocity == [0.0 0.0; 0.0 0.0; 0.0 0.0] + @test b2.store.velocity_half == zeros(3, 2) + @test b2.store.acceleration == zeros(3, 2) + @test b2.store.b_int == zeros(3, 4) + @test b2.store.b_ext == zeros(3, 2) + @test b2.store.damage ≈ [2/3, 2/3] + @test b2.store.bond_active == [0, 0, 1, 0, 0, 1] + @test b2.store.n_active_bonds == [1, 1] + + randbint = rand(3, 4) + b2.store.b_int .= randbint + b1.store.b_int .+= 1 + Peridynamics.exchange_write_fields!(tdh, 1) + + @test b1.store.position[:,1:2] ≈ [0.0 1.0; 0.0 0.0; 0.0 0.0] + @test b1.store.position[:,3:4] ≈ randpos[:,1:2] + @test b1.store.displacement == zeros(3, 2) + @test b1.store.velocity == [1.0 1.0; 0.0 0.0; 0.0 0.0] + @test b1.store.velocity_half == zeros(3, 2) + @test b1.store.acceleration == zeros(3, 2) + @test b1.store.b_int[:,1:2] ≈ 1 .+ randbint[:,3:4] + @test b1.store.b_int[:,3:4] ≈ ones(3, 2) + @test b1.store.b_ext == zeros(3, 2) + @test b1.store.damage ≈ [2/3, 2/3] + @test b1.store.bond_active == [1, 0, 0, 1, 0, 0] + @test b1.store.n_active_bonds == [1, 1] + + @test b2.store.position ≈ randpos + @test b2.store.displacement == zeros(3, 2) + @test b2.store.velocity == [0.0 0.0; 0.0 0.0; 0.0 0.0] + @test b2.store.velocity_half == zeros(3, 2) + @test b2.store.acceleration == zeros(3, 2) + @test b2.store.b_int ≈ randbint + @test b2.store.b_ext == zeros(3, 2) + @test b2.store.damage ≈ [2/3, 2/3] + @test b2.store.bond_active == [0, 0, 1, 0, 0, 1] + @test b2.store.n_active_bonds == [1, 1] + + Peridynamics.exchange_write_fields!(tdh, 2) + + @test b1.store.position[:,1:2] ≈ [0.0 1.0; 0.0 0.0; 0.0 0.0] + @test b1.store.position[:,3:4] ≈ randpos[:,1:2] + @test b1.store.displacement == zeros(3, 2) + @test b1.store.velocity == [1.0 1.0; 0.0 0.0; 0.0 0.0] + @test b1.store.velocity_half == zeros(3, 2) + @test b1.store.acceleration == zeros(3, 2) + @test b1.store.b_int[:,1:2] ≈ 1 .+ randbint[:,3:4] + @test b1.store.b_int[:,3:4] ≈ ones(3, 2) + @test b1.store.b_ext == zeros(3, 2) + @test b1.store.damage ≈ [2/3, 2/3] + @test b1.store.bond_active == [1, 0, 0, 1, 0, 0] + @test b1.store.n_active_bonds == [1, 1] + + @test b2.store.position ≈ randpos + @test b2.store.displacement == zeros(3, 2) + @test b2.store.velocity == [0.0 0.0; 0.0 0.0; 0.0 0.0] + @test b2.store.velocity_half == zeros(3, 2) + @test b2.store.acceleration == zeros(3, 2) + @test b2.store.b_int[:,1:2] ≈ 1 .+ randbint[:,1:2] + @test b2.store.b_int[:,3:4] ≈ randbint[:,3:4] + @test b2.store.b_ext == zeros(3, 2) + @test b2.store.damage ≈ [2/3, 2/3] + @test b2.store.bond_active == [0, 0, 1, 0, 0, 1] + @test b2.store.n_active_bonds == [1, 1] +end diff --git a/test/integration/test_material.jl b/test/integration/test_material.jl new file mode 100644 index 00000000..622dafe4 --- /dev/null +++ b/test/integration/test_material.jl @@ -0,0 +1,61 @@ +struct TestMaterial <: Peridynamics.AbstractMaterial end +struct TestPointParameters <: Peridynamics.AbstractPointParameters + δ::Float64 + rho::Float64 + E::Float64 + nu::Float64 + G::Float64 + K::Float64 + λ::Float64 + μ::Float64 + Gc::Float64 + εc::Float64 +end +Peridynamics.point_param_type(::TestMaterial) = TestPointParameters +Peridynamics.allowed_material_kwargs(::TestMaterial) = Peridynamics.DEFAULT_POINT_KWARGS +function Peridynamics.get_point_params(::TestMaterial, p::Dict{Symbol,Any}) + δ = Peridynamics.get_horizon(p) + rho = Peridynamics.get_density(p) + E, nu, G, K, λ, μ = Peridynamics.get_elastic_params(p) + Gc, εc = Peridynamics.get_frac_params(p, δ, K) + return TestPointParameters(δ, rho, E, nu, G, K, λ, μ, Gc, εc) +end +Peridynamics.discretization_type(::TestMaterial) = Peridynamics.BondDiscretization +function Peridynamics.init_discretization(body::Body{TestMaterial}, args...) + return Peridynamics.init_bond_discretization(body, args...) +end +struct TestVerletStorage <: Peridynamics.AbstractStorage + position::Matrix{Float64} + displacement::Matrix{Float64} + velocity::Matrix{Float64} + velocity_half::Matrix{Float64} + acceleration::Matrix{Float64} + b_int::Matrix{Float64} + b_ext::Matrix{Float64} + damage::Vector{Float64} + bond_active::Vector{Bool} + n_active_bonds::Vector{Int} +end +function Peridynamics.storage_type(::TestMaterial, ::Peridynamics.VelocityVerlet) + return TestVerletStorage +end +function Peridynamics.init_storage(::Body{TestMaterial}, ::Peridynamics.VelocityVerlet, + bd::Peridynamics.BondDiscretization, + ch::Peridynamics.ChunkHandler) + n_loc_points = length(ch.loc_points) + position = copy(bd.position) + displacement = zeros(3, n_loc_points) + velocity = zeros(3, n_loc_points) + velocity_half = zeros(3, n_loc_points) + acceleration = zeros(3, n_loc_points) + b_int = zeros(3, length(ch.point_ids)) + b_ext = zeros(3, n_loc_points) + damage = zeros(n_loc_points) + bond_active = ones(Bool, length(bd.bonds)) + n_active_bonds = copy(bd.n_neighbors) + return TestVerletStorage(position, displacement, velocity, velocity_half, + acceleration, b_int, b_ext, damage, bond_active, + n_active_bonds) +end +Peridynamics.reads_from_halo(::Type{TestMaterial}) = (:position,) +Peridynamics.writes_to_halo(::Type{TestMaterial}) = (:b_int,)