Skip to content

Remove unused code #9

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
rouckas merged 1 commit into from
Feb 24, 2026
Merged

Remove unused code #9

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
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
102 changes: 0 additions & 102 deletions src/Diagnostic.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -221,105 +221,3 @@ function plt(probe::PSxProbe)
heatmap(y, x, array, yflip = true, legend=false,
title="Nvx", titlefont=font(10), c=:bwr, clim=(-lim, lim))
end



mutable struct Diagnostic
nt::Int64
nv::Int64
nx::Int64
ts::Vector{Float64}
rhos::Array{Float64, 2}
particle_rhos::Array{Float64, 3}
Us::Array{Float64, 2}
Es::Array{Float64, 2}
vs::Array{Float64, 2}
PSs::Array{Float64, 3}
xmax::Float64
vmax::Float64
tmax::Float64
vbins::Vector{Float64}
N::Int64
end

function Diagnostic(nx::Int64, nv::Int64, nt::Int64, nspecies::Int64, xmax::Float64, vmax::Float64, tmax::Float64)
Diagnostic(nt,
nv,
nx,
zeros(nt),
zeros(nt, nx),
zeros(nt, nspecies, nx),
zeros(nt, nx),
zeros(nt, nx),
zeros(nt, nv),
zeros(nt, nx, nv),
xmax,
vmax,
tmax,
LinRange(-vmax-xmax/(nv-1),vmax+xmax/(nv-1),nv+1),
0
)
end

function histogram_v(particles::ParticleEnsemble, diag::Diagnostic, dir::Int64=1)
vs = @view diag.vs[diag.N,:]
for p in particles.coords
v = p.v[dir]
if abs(v) >= diag.vmax
continue
end
vs[Int64(floor(diag.nv * (v / (2 * diag.vmax) + 0.5))) + 1] += particles.q
end
end

function histogram_xv(particles::ParticleEnsemble, pic::PIC, diag::Diagnostic)
PSs = @view diag.PSs[diag.N,:,:]
for p in particles.coords
v = p.v[1]
x = p.r[1]
if abs(v) >= diag.vmax
continue
end
PSs[Int64(floor(x*(pic.nx-1)/pic.xmax+0.5)+1),
Int64(floor(diag.nv * (v / (2 * diag.vmax) + 0.5))) + 1] += particles.q
end
end

function sample(diag::Diagnostic, t::Float64, pic::PIC)
diag.N += 1

for p in pic.particles
# diag.vs[diag.N] += fit(Histogram, p.vx, diag.vbins).weights*p.q
histogram_v(p, diag)
# diag.PSs[diag.N] += fit(Histogram, (p.x, p.vx), (pic.xbins, diag.vbins)).weights*p.q
histogram_xv(p, pic, diag)

end
diag.Us[diag.N,:] .= pic.U
diag.Es[diag.N,:] .= [e[1] for e in pic.E]
diag.rhos[diag.N,:] .= pic.rho
diag.ts[diag.N] = t

end

function plt(diag::Diagnostic, pic::PIC)

function plt_record(x, y, array, title)
lim = maximum(abs.(array))
return heatmap(x, y, array, yflip = true, legend=false,
title=title, titlefont=font(10), c=:bwr, clim=(-lim, lim))
end

p1 = plt_record(pic.x, diag.ts[1:diag.N], diag.rhos[1:diag.N,:], "Charge density")
p2 = plt_record(pic.x, diag.ts[1:diag.N], diag.Us[1:diag.N,:], "Electrical potential")
p3 = plt_record(pic.x, diag.ts[1:diag.N], diag.Es[1:diag.N,:], "Elextrical field")
vx = LinRange(-diag.vmax,diag.vmax,diag.nv)
p4 = plt_record(vx, diag.ts[1:diag.N], diag.vs[1:diag.N,:], "Speed")

return plot(p1, p2, p3, p4, size = (900, 250), layout=(1,4))
end

function plt_PS(diag::Diagnostic, lims=1e-18)
x = dropdims(sum(diag.PSs, dims=1), dims=1)'
return heatmap(x, size=(500, 500), legend=false, c=:bwr, clim=(-lims, lims), title="Phase space")
end
24 changes: 6 additions & 18 deletions src/SimplePIC.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -21,7 +21,7 @@ export init_time, init_monoenergetic, init_thermal
export advance, advance!, energy, scatter

export epsilon_0, e
export PIC, Diagnostic, RhoProbe, EnergyProbe, NxProbe, NvxProbe, EProbe, PSxProbe, UProbe, sample!
export PIC, RhoProbe, EnergyProbe, NxProbe, NvxProbe, EProbe, PSxProbe, UProbe, sample!
export sample, poisson_solve, interpolate, advance, init_leapfrog, solve_init, solve_init_fft, particle_bc

export random_maxwell_v, random_maxwell_vcomponent, random_maxwell_vflux
Expand Down Expand Up @@ -86,12 +86,6 @@ mutable struct PIC{dim, ParticleType, GeometryType, BCType, vdim}
#where ParticleType <: AbstractParticle where FieldType <: AbstractField where BCType <: BoundaryCondition
particles::Vector{ParticleEnsemble{ParticleType}}
interactions::Vector{Interactions{ParticleEnsemble{ParticleType}}}
nx::Int64
xmax::Float64
dx::Float64
x::Vector{Float64}
xbins::Vector{Float64}
dV::Float64
geo::GeometryType
BC::BCType
rho::Array{Float64, dim}
Expand All @@ -105,12 +99,6 @@ function PIC(particles::Vector{ParticleEnsemble{ParticleType}}, interactions::Ve
geo::AbstractGeometry, BC::BCType, epsilon_0::Float64, vdim::Int64) where ParticleType where BCType <: BoundaryCondition
PIC{1, ParticleType, Cartesian1D, BCType, vdim}(particles,
interactions,
geo.nx,
geo.xmax,
geo.dx,
LinRange(0, geo.xmax, geo.nx),
[i-0.5 for i in 0:geo.nx]*geo.xmax/(geo.nx-1),
geo.dV,
geo,
BC,
zeros(geo.nx),
Expand Down Expand Up @@ -239,7 +227,7 @@ end
function sample(pic::PIC)
pic.rho .= pic.rhobg
for p in pic.particles
pic.rho .+= sample_linear(p, pic.nx, pic.dx).*(p.q/pic.dV)
pic.rho .+= sample_linear(p, pic.geo.nx, pic.geo.dx).*(p.q/pic.geo.dV)
end

if typeof(pic.BC) <: BCPeriodic1D
Expand All @@ -265,20 +253,20 @@ end

function interpolate(pic::PIC)
for p in pic.particles
interpolate_linear(p, pic.dx, pic.E)
interpolate_linear(p, pic.geo.dx, pic.E)
end
end

function particle_bc(pic::PIC)
for p in pic.particles
particle_bc(p, pic.xmax, pic.BC)
particle_bc(p, pic.geo.xmax, pic.BC)
end
end

function advance(pic::PIC, dt::Float64, B=SVector(0., 0., 0.))
for (p, inter) in zip(pic.particles, pic.interactions)
advance(p, inter, dt, B)
particle_bc(p, pic.xmax, pic.BC)
particle_bc(p, pic.geo.xmax, pic.BC)
end
end

Expand All @@ -288,7 +276,7 @@ function init_leapfrog(pic::PIC, dt::Float64)
interpolate(pic)
for p in pic.particles
advance_v(p, -dt/2)
particle_bc(p, pic.xmax, pic.BC)
particle_bc(p, pic.geo.xmax, pic.BC)
end
end

Expand Down