Skip to content

Add multithreading support in the electric field calculation #509

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

Merged
fhagemann merged 2 commits into from
Jul 12, 2025
Merged

Add multithreading support in the electric field calculation #509

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
d69dac4
Add multithreading to electric field calculation
fhagemann Jul 9, 2025
4e2c897
Remove `fieldvector_coordinates` from `get_electric_field_from_potent…
fhagemann Jul 9, 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
23 changes: 11 additions & 12 deletions src/ElectricField/ElectricField.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -49,20 +49,20 @@ Base.convert(T::Type{ElectricField}, x::NamedTuple) = T(x)



function ElectricField(epot::ElectricPotential{T, 3, S}, point_types::PointTypes{T}) where {T, S}
return ElectricField{T, 3, S, typeof(grid.axes)}(get_electric_field_from_potential( epot, point_types ), epot.grid)
function ElectricField(epot::ElectricPotential{T, 3, S}, point_types::PointTypes{T}; use_nthreads::Int = Base.Threads.nthreads()) where {T, S}
return ElectricField{T, 3, S, typeof(grid.axes)}(get_electric_field_from_potential( epot, point_types; use_nthreads ), epot.grid)
end


function get_electric_field_from_potential(epot::ElectricPotential{T, 3, Cylindrical}, point_types::PointTypes{T}, fieldvector_coordinates=:xyz)::ElectricField{T, 3, Cylindrical} where {T <: SSDFloat}
function get_electric_field_from_potential(epot::ElectricPotential{T, 3, Cylindrical}, point_types::PointTypes{T}; use_nthreads::Int = Base.Threads.threadid())::ElectricField{T, 3, Cylindrical} where {T <: SSDFloat}
p = epot.data
axr::Vector{T} = collect(epot.grid.axes[1])
axφ::Vector{T} = collect(epot.grid.axes[2])
axz::Vector{T} = collect(epot.grid.axes[3])

cyclic::T = epot.grid.axes[2].interval.right
ef = Array{SVector{3, T}}(undef, size(p)...)
for iz in 1:size(ef, 3)
@onthreads 1:use_nthreads for iz in workpart(1:size(ef, 3), 1:use_nthreads, Base.Threads.threadid())
for iφ in 1:size(ef, 2)
for ir in 1:size(ef, 1)
### r ###
Expand Down Expand Up @@ -157,9 +157,7 @@ function get_electric_field_from_potential(epot::ElectricPotential{T, 3, Cylindr
end
end
end
if fieldvector_coordinates == :xyz
ef = convert_field_vectors_to_xyz(ef, axφ)
end
Comment on lines -160 to -162
Copy link
Collaborator Author

@fhagemann fhagemann Jul 9, 2025

Choose a reason for hiding this comment

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

Also, I removed fieldvector_coordinates from get_electric_field_from_potential, because I cannot think of a scenario where this should NOT be wanted..

Maybe this function can even be mutating to avoid allocating memory.

ef = convert_field_vectors_to_xyz(ef, axφ)
return ElectricField(ef, point_types.grid)
end

Expand All @@ -177,6 +175,7 @@ function convert_field_vectors_to_xyz(field::Array{SArray{Tuple{3},T,1,3},3}, φ
end



function interpolated_scalarfield(spot::ScalarPotential{T, 3, Cylindrical}) where {T}
@inbounds knots = spot.grid.axes[1].ticks, cat(spot.grid.axes[2].ticks,T(2π),dims=1), spot.grid.axes[3].ticks
ext_data = cat(spot.data, spot.data[:,1:1,:], dims=2)
Expand Down Expand Up @@ -212,17 +211,17 @@ end



function get_electric_field_from_potential(epot::ElectricPotential{T, 3, Cartesian}, point_types::PointTypes{T})::ElectricField{T, 3, Cartesian} where {T <: SSDFloat}
function get_electric_field_from_potential(epot::ElectricPotential{T, 3, Cartesian}, point_types::PointTypes{T}; use_nthreads::Int = Base.Threads.nthreads())::ElectricField{T, 3, Cartesian} where {T <: SSDFloat}
axx::Vector{T} = collect(epot.grid.axes[1])
axy::Vector{T} = collect(epot.grid.axes[2])
axz::Vector{T} = collect(epot.grid.axes[3])
axx_ext::Vector{T} = get_extended_ticks(epot.grid.axes[1])
axy_ext::Vector{T} = get_extended_ticks(epot.grid.axes[2])
axz_ext::Vector{T} = get_extended_ticks(epot.grid.axes[3])
# axx_ext::Vector{T} = get_extended_ticks(epot.grid.axes[1])
# axy_ext::Vector{T} = get_extended_ticks(epot.grid.axes[2])
# axz_ext::Vector{T} = get_extended_ticks(epot.grid.axes[3])

ef::Array{SVector{3, T}} = Array{SVector{3, T}}(undef, size(epot.data))

for ix in eachindex(axx)
@onthreads 1:use_nthreads for ix in workpart(eachindex(axx), 1:use_nthreads, Base.Threads.threadid())
for iy in eachindex(axy)
for iz in eachindex(axz)
if ix - 1 < 1
Expand Down
6 changes: 3 additions & 3 deletions src/Simulation/Simulation.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1167,7 +1167,7 @@ Calculates the [`ElectricField`](@ref) from the [`ElectricPotential`](@ref) stor
!!! note
This method only works if `sim.electric_potential` has already been calculated and is not `missing`.
"""
function calculate_electric_field!(sim::Simulation{T, CS}; n_points_in_φ::Union{Missing, Int} = missing)::Nothing where {T <: SSDFloat, CS}
function calculate_electric_field!(sim::Simulation{T, CS}; n_points_in_φ::Union{Missing, Int} = missing, use_nthreads::Int = Base.Threads.nthreads())::Nothing where {T <: SSDFloat, CS}
@assert !ismissing(sim.electric_potential) "Electric potential has not been calculated yet. Please run `calculate_electric_potential!(sim)` first."
periodicity::T = width(sim.world.intervals[2])
e_pot, point_types = if CS == Cylindrical && periodicity == T(0) # 2D, only one point in φ
Expand All @@ -1181,15 +1181,15 @@ function calculate_electric_field!(sim::Simulation{T, CS}; n_points_in_φ::Union
end
end
get_2π_potential(sim.electric_potential, n_points_in_φ = n_points_in_φ),
get_2π_potential(sim.point_types, n_points_in_φ = n_points_in_φ);
get_2π_potential(sim.point_types, n_points_in_φ = n_points_in_φ)
elseif CS == Cylindrical
get_2π_potential(sim.electric_potential),
get_2π_potential(sim.point_types)
else
sim.electric_potential,
sim.point_types
end
sim.electric_field = get_electric_field_from_potential(e_pot, point_types);
sim.electric_field = get_electric_field_from_potential(e_pot, point_types; use_nthreads)
nothing
end

Expand Down
Loading